%% $OpenXM: OpenXM/src/kan96xx/Doc/sst.sm1,v 1.1 2000/02/09 02:30:39 takayama Exp $
%% Sample implementations of small algorithms in the book SST.
%% gbhg3/Int/sst.sm1 --> kan96xx/lib
[(parse) (oxasir.sm1) pushfile] extension

%% 1999, 9/21  saturation, intersection are moved to SSkan/lib/complex.sm1

/stdpair0 {
  /arg1 set
  [/in-stdpair0 /ff /mm /vlist /sigma /sigma* /srule
   /msigma /satMsigma ] pushVariables
  [(CurrentRingp) (KanGBmessage)] pushEnv
  [
     /ff arg1 def
     /mm ff 0 get def
     /sigma ff 1 get {toString} map def
     /vlist [ff 2 get to_records pop ] def
     
     sigma vlist complement /sigma* set
     [vlist from_records ring_of_polynomials 0] define_ring 
     sigma { [ 2 1 roll . (1).] } map /srule set
     mm {toString . srule replace} map /msigma set
     [(KanGBmessage) 0] system_variable    
     [msigma] groebner_sugar 0 get /msigma set

     [msigma sigma* sigma* from_records] saturation /satMsigma set
     [satMsigma msigma sigma*] /arg1 set
  ] pop
  popEnv
  popVariables
  arg1
} def
(stdpair0 ) message

[(stdpair0)
[ $[M sigma vlist] stdpair0 [M_sigma M sigma*]$
  $Example: [[(x2^2) (x2 x4) (x2 x3) (x1 x4^2)]  $
  $          [(x1) (x3)] (x1,x2,x3,x4)] stdpair0 :: $
]] putUsages

/stdpair {
  /arg1 set
  [/in-stdpair /ffs /vlist /mm /ass /i /pp /vlist0 /ass*] pushVariables
  [(CurrentRingp) (KanGBmessage)] pushEnv
  [
     /ffs arg1 def
     ffs 0 get /mm set
     ffs 1 get /vlist0 set

     [vlist0 to_records pop] /vlist set
     [mm vlist] primadec /ass set
     ass { 1 get } map /ass set
     %%(Associated primes are : ) messagen ass pmat 
     ass { {toString} map vlist complement } map /ass* set

     /pp ass length def
    [
     0 1 pp 1 sub {
       /i set
       %%(stdpair0 of : ) messagen
       %%[mm ass* i get vlist0] pmat
       [mm ass* i get vlist0] stdpair0 
       [ass* i get] join
     } for
    ] /arg1 set
  ] pop
  popEnv
  popVariables
  arg1
} def
[(stdpair)
[ $ [[ A  B  sigma^c  sigma] [  ] ....] $
  $   Monomials in A/B give u part of (x^u,sigma). $
  $Example 1: [[(x2^2) (x2 x4) (x2 x3) (x1 x4^2)]  $
  $            (x1,x2,x3,x4)] stdpair pmat $
  $Example 2: [[(x2 x3) (x1^3 x4^2)]  $
  $            (x1,x2,x3,x4,x5)] stdpair pmat $
  $Example 3: [[(x4^2) (x3^7) (x2 x3^6) (x2^2)]  $
  $            (x1,x2,x3,x4,x5)] stdpair pmat $
]] putUsages
(stpair ) message

/stien {
  (See sec 5.6, Consider submodule generated by x5. Or, annihilating ideal of x5 which can be obtained by the syzygy computation. See the paper by Stienstra.) message
  %%  (x1+x2+x3+x4+x5+1)  does not work!  
  /ff1
   [(x1+x2+x3+x4+x5) (x1+x2-x4) (x2+x3-x4) (x1 x3) (x2 x4)]
  def
  [(x1,x2,x3,x4,x5) ring_of_polynomials
    [[(x1) 1] [(x2) 1] [(x3) 1] [(x4) 1]] weight_vector 0] define_ring
  ff1 {. } map /ff1 set
  [ff1] groebner_sugar 0 get /ff2 set  ff2 pmat
  [[(x5).] ff1 join [(needSyz)]] groebner_sugar /ff3 set
  ff3 2 get /ff4.syz set
  ff4.syz [[(x5).] ff1 join] transpose mul pmat  %% check!
  ff4.syz { 0 get dup (0). eq { pop } { } ifelse } map /ann-x5 set ann-x5 pmat
  [ff1 ann-x5 join] groebner_sugar 0 get pmat
} def
(stien ) messagen

/stien2 {
  (See sec 5.6, Consider submodule generated by x5+1. Or, annihilating ideal of x5+1 which can be obtained by the syzygy computation. See the paper by Stienstra.) message
  /ff1
   [(x1+x2+x3+x4+x5+1) (x1+x2-x4) (x2+x3-x4) (x1 x3) (x2 x4)]
  def
  [(x1,x2,x3,x4,x5) ring_of_polynomials
    [[(x1) 1] [(x2) 1] [(x3) 1] [(x4) 1]] weight_vector 0] define_ring
  ff1 {. } map /ff1 set
  [ff1] groebner_sugar 0 get /ff2 set  ff2 pmat
  [[(x5+1).] ff1 join [(needSyz)]] groebner_sugar /ff3 set
  ff3 2 get /ff4.syz set
  ff4.syz [[(x5+1).] ff1 join] transpose mul pmat  %% check!
  ff4.syz { 0 get dup (0). eq { pop } { } ifelse } map /ann-x5 set ann-x5 pmat
  [ff1 ann-x5 join] groebner_sugar 0 get /ff-stien set ff-stien pmat
  [[(3 x5 + 3 + 8 x4). (x5 + 1 + 4 x3). (x5 + 1 + 8 x2). (x5 + 1 + 4 x1).
   (x5^2 + 2 x5 + 1).]] groebner_sugar 0 get /ff-int set ff-int pmat
  [(Homogenize) 0] system_variable
  ff-stien { ff-int reduction 0 get } map pmat
  ff-int { ff-stien reduction 0 get } map pmat
  [(Homogenize) 1] system_variable

} def
(stien2 ) messagen

/stien3 {
  (See sec 5.6,) message
  [(x1,x2,y1,y2,z,zz) ring_of_polynomials
   [[(zz) 1]] weight_vector
  0] define_ring
  [(x1 z +x1 x2 + x2 z + y1 y2 + z^2). 
   (x1 y1 - z^2).
   (x2 y2 - z^2).
   (z zz -1).
  ]
  /ff1 set
  [ff1] groebner_sugar 0 get /ff2 set  

  ff2 [(zz)] eliminatev /ff2 set
  ff2 pmat
  (--------------------------- see stien3.mac --) message
  (  ) message
  [(x1 z  + x2 z + y1 y2 + z^2). 
   (x1 y1 - z^2).
   (x2 y2 - z^2).
   (z zz -1).
  ]
  /ff1a set
  [ff1a] groebner_sugar 0 get /ff2a set  

  ff2a [(zz)] eliminatev /ff2a set
  ff2a pmat

} def
(stien3 ) messagen

/stien3a {
  (See sec 5.6, msri/book/stien3.mac ) message
  [(a,b,A,B,z,zz) ring_of_polynomials
   [[(zz) 1]] weight_vector
  0] define_ring
  [(B A z + B z^2+ B a z + B a^2 + A z^2+ z^3 + a z^2+ b A z+ b z^2+ b^2 A).
   (a A - z^2).
   (b B - z^2).
   (z zz -1).
  ]
  /ff1 set
  [ff1] groebner_sugar 0 get /ff2 set  

  ff2 [(zz)] eliminatev /ff2 set
  ff2 pmat
  (--------------------------- see stien3a.mac , ff2--) message
  1 1 ff2 length { /ii set
    ff2 ii 1 sub get /ff22 set
    (poly f) messagen ii messagen ( ) messagen ff22 message  
  } for

  (---------------------------------------------) message
  [(B A z + 2 B z^2+ 3 B a z + 4 B a^2 + 5 A z^2+ z^3 + 6 a z^2+ 7 b A z+ 8 b z^2+ 9 b^2 A).
   (a A - z^2).
   (b B - z^2).
   (z zz -1).
  ]
  /ff1 set
  [ff1] groebner_sugar 0 get /ff2 set  

  ff2 [(zz)] eliminatev /ff2 set
  ff2 pmat
  (--------------------------- see stien3a.mac , ff2--) message
  1 1 ff2 length { /ii set
    ff2 ii 1 sub get /ff22 set
    (poly f) messagen ii messagen ( ) messagen ff22 message  
  } for
} def
(stien3a ) message


(  ) message

/foo1 {
            [(x1 y1 + x2 y2 + x3 y3 + x4 y4) 
             (x2 y2 + x4 y4) (x3 y3 + x4 y4) (y1 y4 - y2 y3)
             (y1 u1 -1) (y2 u2 -1) (y3 u3 -1) (y4 u4-1)] /ff set
            [(x1,x2,x3,x4,y1,y2,y3,y4,u1,u2,u3,u4) ring_of_polynomials
             [[(y1) 1 (y2) 1 (y3) 1 (y4) 1 (u1) 1 (u2) 1 (u3) 1 (u4) 1]]
             weight_vector 0] define_ring
            [ff {.} map] groebner_sugar 0 get
            [(y1) (y2) (y3) (y4) (u1) (u2) (u3) (u4)] eliminatev ::
} def