% $OpenXM: OpenXM/src/kan96xx/Doc/dhecart.sm1,v 1.4 2005/06/19 08:29:02 takayama Exp $
% Stdbasis via the double homogenization:  dx x = x dx + h H
% Homogenize=3
(ecart_loaded) boundp { }
{ [(parse) (ecart.sm1) pushfile] extension } ifelse

/dh.begin {
  [(Homogenize) (AutoReduce) (KanGBmessage)] pushEnv /dh.saved.env set
  [(Homogenize) 3] system_variable
  dh.autoReduce { [(AutoReduce) 1] system_variable } { } ifelse
} def

/dh.end {
  dh.saved.env popEnv
  [(Homogenize) 1] system_variable
} def

/dh.dehomogenize {
  dehomogenize
} def

% Global environmental variables
/dh.gb.verbose 1 def
/dh.autoHomogenize 1 def
/dh.autoReduce 1 def
/dh.needSyz 0 def

/dh.message {
   (dh.ecart: ) messagen  message
} def
/dh.messagen {
   (dh.ecart: ) messagen  messagen
} def

%%test
% [(x,y) ring_of_differential_operators [[(Dx) 1]] ecart.weight_vector 0] define_ring ;   dh.begin ;
% [[(x Dx + 1). homogenize]] groebner ::

%%test
% [ [(x Dx + y Dy + 1) (x Dx y Dy -1)] (x,y) [[(x) -1 (y) -1]]] dh.gb pmat
%  --> It is not an admissible order.  
% [ [(x Dx + y Dy + 1) (x Dx y Dy -1)] (x,y) [[(Dx) 1 (Dy) 1 (x) -1 (y) -1] [(Dx) 1 (Dy) 1] [(x) -1 (y) -1]]] dh.gb pmat
   
/dh.gb {
  /arg1 set
  [/in-dh.gb /aa /typev /setarg /f /v 
   /gg /wv /vec /ans /rr /mm
   /env2 /ans.gb /groebnerOptions
  ] pushVariables 
  [(CurrentRingp) (KanGBmessage)] pushEnv
  [
    /aa arg1 def
    aa isArray { } { ( << array >> dh.gb) error } ifelse
    aa getAttributeList configureGroebnerOption /groebnerOptions set
    /setarg 0 def
    /wv 0 def

    aa { tag } map /typev set
    typev [ ArrayP ] eq
    {  /f aa 0 get def
       /v gb.v def
       /setarg 1 def
    } { } ifelse
    typev [ArrayP StringP] eq
    {  /f aa 0 get def
       /v aa 1 get def
       /setarg 1 def
    } { } ifelse
    typev [ArrayP RingP] eq
    {  /f aa 0 get def
       /v aa 1 get def
       /setarg 1 def
    } { } ifelse
    typev [ArrayP ArrayP] eq
    {  /f aa 0 get def
       /v aa 1 get from_records def
       /setarg 1 def
    } { } ifelse
    typev [ArrayP StringP ArrayP] eq
    {  /f aa 0 get def
       /v aa 1 get def
       /wv aa 2 get def
       /setarg 1 def
    } { } ifelse
    typev [ArrayP ArrayP ArrayP] eq
    {  /f aa 0 get def
       /v aa 1 get from_records def
       /wv aa 2 get def
       /setarg 1 def
    } { } ifelse

    /env1 getOptions def

    setarg { } { (dh.gb : Argument mismatch) error } ifelse
    
    [(KanGBmessage) dh.gb.verbose ] system_variable

    %%% Start of the preprocess
    v tag RingP eq {
       /rr v def 
    }{
      f getRing /rr set
    } ifelse
    %% To the normal form : matrix expression.
    f gb.toMatrixOfString /f set
    /mm gb.itWasMatrix def

    rr tag 0 eq {
      %% Define the ring.
      v isInteger {
        (Error in dh.gb: Specify variables) error
      } {  } ifelse
      %% wv is set when parsing the arguments.
      wv isInteger {        
        (Give a weight vector) error
      }{
         [v ring_of_differential_operators 
          wv ecart.weight_vector
          gb.characteristic
          ] define_ring
      } ifelse
    } {
      %% Use the ring structre given by the input.
      v isInteger not {
        gb.warning {
         (Warning : the given ring definition is not used.) message
        } { } ifelse
      } {  } ifelse
      rr ring_def
      /wv rr gb.getWeight def
    } ifelse
    %%% Enf of the preprocess

    dh.begin

    [v] ecart.checkOrder 

    groebnerOptions gb.options mergeGroebnerOptions /groebnerOptions set
    gb.verbose { (groebnerOptions = ) messagen groebnerOptions message } { } ifelse


   dh.autoHomogenize not {
% No automatic hH-homogenization.
       f { {. } map} map /f set   
   } {
% Automatic hH-homogenization
      (dh.gb : Input polynomial or vectors are automatically homogenized) dh.message
       f { {. } map} map /f set   
       f { { [[@@@.Hsymbol . (1).] [@@@.hsymbol . (1).] ] replace } map } map /f set
       f { { homogenize } map } map /f set
       f dh.message
   } ifelse

   dh.needSyz {
     [f [(needSyz)] groebnerOptions join ] groebner /gg set
   } {  
     [f groebnerOptions] groebner 0 get /gg set
   } ifelse


    dh.needSyz {
      mm {
       gg 0 get { /tmp set [mm tmp] toVectors } map /ans.gb set
      } { /ans.gb gg 0 get def } ifelse
      /ans [gg 2 get , ans.gb , gg 1 get , f ] def
%      ans pmat ;
    } { 
      wv isInteger {
        /ans [gg gg {init} map] def
      }{
%% Get the initial ideal
        /ans [gg gg {wv 0 get weightv init} map] def
      }ifelse

      %% Postprocess : recover the matrix expression.
      mm {
        ans { /tmp set [mm tmp] toVectors } map
        /ans set
      }{ }
      ifelse
    } ifelse

    dh.end
    ans gg getAttributeList setAttributeList /ans set

    ans getRing (oxRingStructure) dc /dh.gb.oxRingStructure set
    %% 
    env1 restoreOptions  %% degreeShift changes "grade"

    /arg1 ans def
  ] pop
  popEnv
  popVariables
  arg1
} def

[(dh.gb)
 [(a dh.gb b)
  (array a; array b;)
  $b : [g ii];  array g; array in; g is a standard (Grobner) basis of f$
  (             in the ring of differential operators.)
  (The computation is done in the doubly homogenized Weyl algebra.)
  (Dx x = x Dx + h H)
   $            ii is the initial ideal in case of w is given or <<a>> belongs$
   $            to a ring. In the other cases, it returns the initial monominal.$
  (a : [f ];    array f;  f is a set of generators of an ideal in a ring.)
  (a : [f v];   array f; string v;  v is the variables. )
  (a : [f v w]; array f; string v; array of array w; w is the weight matirx.)
  (  )
  (Globals:   dh.autoHomogenize dh.gb.verbose dh.needSyz dh.gb.oxRingStructure)
  (cf. dh.begin dh.end dh.message dh.messagen)
  ( )
  $Example 1: [ [( (x Dx)^2 + (y Dy)^2 -1) ( x y Dx Dy -1)] (x,y) $
  $             [ [ (Dx) 1 ] [(x) -1 (y) -1 (Dx) 1 (Dy) 1]] ] dh.gb pmat ; $
  (Example 2: )
  $ [ [(2 x Dx + 3 y Dy+6) (2 y Dx + 3 x^2 Dy)] (x,y) $
  $   [[(x) -1 (Dx) 1 (y) -1 (Dy) 1]]] dh.gb  /ff set ff pmat ;$
  (To set the current ring to the ring in which ff belongs )
  (      ff getRing ring_def  )
  ( )
  (To set the current ring to the ring in which ff belongs )
  (      ff getRing ring_def  )
  ( )
  (Data:  dh.p1, dh.p2, dh.p3 )
  (In order to get a standard basis of the test data, type in dh.test.p1, ...)
  (  )
  (Example 3: )
  $ /gb.verbose 1 def $
  $ [ [(2 x Dx + 3 y Dy+6) (2 y Dx + 3 x^2 Dy)] (x,y) $
  $   [[(x) -1 (Dx) 1 (y) -1 (Dy) 1]]] $
  $ [(gbCheck) 1] setAttributeList dh.gb getAttributeList :: $
  ( )
]] putUsages

%Test input.
%misc-2003/09/oaku/b.sm1, Granger-Oaku-Takayama, Tangent cone algorithm ...
/dh.p1 {
 [
  [(t-(x^3 - y^2 z^2 - w^2))
   (Dx + (3 x^2 ) Dt)
   (Dy - (2 y z^2) Dt)
   (Dz - (2 y^2 z) Dt)
   (Dw - (2 w ) Dt)
  ]
  [ [(t) -1 (Dt) 1]
    [(Dt) 1 (Dx) 1 (Dy) 1 (Dz) 1 (Dw) 1]
    [(t) -1 (x) -1 (y) -1 (z) -1 (w) -1]]
 ]
} def
/dh.test.p1 {
  [(KanGBmessage) 1] system_variable
  { [dh.p1 0 get (x,y,z,t,w)  dh.p1 1 get] dh.gb } timer
} def

%misc-2003/09/oaku/ob.sm1,
% fw2 [(x) (y) (z) (w)] fw_delta
%  > 30min, degree 25.
/dh.p2 {
 [
  [   (-w^8-z^4-y^3*w-x^3+t)  (3*x^2*Dt+Dx)  (3*y^2*w*Dt+Dy)  (4*z^3*Dt+Dz)  
       (8*w^7*Dt+y^3*Dt+Dw) ]
  [ [(t) -1 (Dt) 1]
    [(Dt) 1 (Dx) 1 (Dy) 1 (Dz) 1 (Dw) 1]
    [(t) -1 (x) -1 (y) -1 (z) -1 (w) -1]]
 ]
} def
/dh.test.p2 {
  [(KanGBmessage) 1] system_variable
  { [dh.p2 0 get (x,y,z,t,w)  dh.p2 1 get] dh.gb } timer
} def

%misc-2003/09/oaku/
% x^3 + (x+1)*y*z,  x^3+x*y*z is easy, but it is difficult in ecart.
/dh.p3 {
 [
  [ $-x^3-x*y*z-y*z+t$ , $3*x^2*Dt+y*z*Dt+Dx$ , $x*z*Dt+z*Dt+Dy$ , 
    $x*y*Dt+y*Dt+Dz$ ]
  [ [(t) -1 (Dt) 1]
    [(Dt) 1 (Dx) 1 (Dy) 1 (Dz) 1]
    [(t) -1 (x) -1 (y) -1 (z) -1]]
 ]
} def
/dh.test.p3 {
  [(KanGBmessage) 1] system_variable
  { [dh.p3 0 get (x,y,z,t)  dh.p3 1 get] dh.gb } timer
} def