[BACK]Return to minimal.k CVS log [TXT][DIR] Up to [local] / OpenXM / src / k097 / lib / minimal

Diff for /OpenXM/src/k097/lib/minimal/minimal.k between version 1.5 and 1.25

version 1.5, 2000/05/05 08:13:49 version 1.25, 2000/08/02 05:14:31
Line 1 
Line 1 
 /* $OpenXM: OpenXM/src/k097/lib/minimal/minimal.k,v 1.4 2000/05/04 11:05:20 takayama Exp $ */  /* $OpenXM: OpenXM/src/k097/lib/minimal/minimal.k,v 1.24 2000/08/02 03:23:36 takayama Exp $ */
 #define DEBUG 1  #define DEBUG 1
 /* #define ORDINARY 1 */  Sordinary = false;
 /* If you run this program on openxm version 1.1.2 (FreeBSD),  /* If you run this program on openxm version 1.1.2 (FreeBSD),
    make a symbolic link by the command     make a symbolic link by the command
    ln -s /usr/bin/cpp /lib/cpp     ln -s /usr/bin/cpp /lib/cpp
 */  */
   #define OFFSET 0
   /* #define OFFSET 20*/
 /* Test sequences.  /* Test sequences.
    Use load["minimal.k"];;     Use load["minimal.k"];;
   
Line 32  def load_tower() {
Line 34  def load_tower() {
   if (Boundp("k0-tower.sm1.loaded")) {    if (Boundp("k0-tower.sm1.loaded")) {
   }else{    }else{
     sm1(" [(parse) (k0-tower.sm1) pushfile ] extension ");      sm1(" [(parse) (k0-tower.sm1) pushfile ] extension ");
       sm1(" [(parse) (new.sm1) pushfile ] extension ");
     sm1(" /k0-tower.sm1.loaded 1 def ");      sm1(" /k0-tower.sm1.loaded 1 def ");
   }    }
     sm1(" oxNoX ");
 }  }
 load_tower();  load_tower();
 SonAutoReduce = true;  SonAutoReduce = true;
Line 46  def Reverse(f) {
Line 50  def Reverse(f) {
 def Sgroebner(f) {  def Sgroebner(f) {
    sm1(" [f] groebner /FunctionValue set");     sm1(" [f] groebner /FunctionValue set");
 }  }
   
   def Sinvolutive(f,w) {
     local g,m;
     if (IsArray(f[0])) {
       m = NewArray(Length(f[0]));
     }else{
       m = [0];
     }
     g = Sgroebner(f);
     /* This is a temporary code. */
     sm1(" g 0 get { w m init_w<m>} map /FunctionValue set ");
   }
   
   
   
   def Error(s) {
     sm1(" s error ");
   }
   
   def IsNull(s) {
     if (Stag(s) == 0) return(true);
     else return(false);
   }
   
   def MonomialPart(f) {
     sm1(" [(lmonom) f] gbext /FunctionValue set ");
   }
   
   def Warning(s) {
     Print("Warning: ");
     Println(s);
   }
   def RingOf(f) {
     local r;
     if (IsPolynomial(f)) {
       if (f != Poly("0")) {
         sm1(f," (ring) dc /r set ");
       }else{
         sm1(" [(CurrentRingp)] system_variable /r set ");
       }
     }else{
       Warning("RingOf(f): the argument f must be a polynomial. Return the current ring.");
       sm1(" [(CurrentRingp)] system_variable /r set ");
     }
     return(r);
   }
   
   def Ord_w_m(f,w,m) {
     sm1(" f  w  m ord_w<m> { (universalNumber) dc } map /FunctionValue set ");
   }
   HelpAdd(["Ord_w_m",
   ["Ord_w_m(f,w,m) returns the order of f with respect to w with the shift m.",
    "Note that the order of the ring and the weight w must be the same.",
    "When f is zero, it returns -intInfinity = -999999999.",
    "Example:  Sweyl(\"x,y\",[[\"x\",-1,\"Dx\",1]]); ",
    "          Ord_w_m([x*Dx+1,Dx^2+x^5],[\"x\",-1,\"Dx\",1],[2,0]):"]]);
   
   def Init_w_m(f,w,m) {
     sm1(" f w m init_w<m> /FunctionValue set ");
   }
   HelpAdd(["Init_w_m",
   ["Init_w_m(f,w,m) returns the initial of f with respect to w with the shift m.",
    "Note that the order of the ring and the weight w must be the same.",
    "Example:  Sweyl(\"x,y\",[[\"x\",-1,\"Dx\",1]]); ",
    "          Init_w_m([x*Dx+1,Dx^2+x^5],[\"x\",-1,\"Dx\",1],[2,0]):"]]);
   
   def Max(v) {
     local i,t,n;
     n = Length(v);
     if (n == 0) return(null);
     t = v[0];
     for (i=0; i<n; i++) {
       if (v[i] > t) { t = v[i];}
     }
     return(t);
   }
   HelpAdd(["Max",
   ["Max(v) returns the maximal element in v."]]);
   
   /*  End of standard functions that should be moved to standard libraries. */
 def test0() {  def test0() {
   local f;    local f;
   Sweyl("x,y,z");    Sweyl("x,y,z");
Line 128  sm1(" [(AvoidTheSameRing)] pushEnv 
Line 212  sm1(" [(AvoidTheSameRing)] pushEnv 
       [ [(AvoidTheSameRing) 0] system_variable        [ [(AvoidTheSameRing) 0] system_variable
         [(gbListTower) tower (list) dc] system_variable          [(gbListTower) tower (list) dc] system_variable
       ] pop popEnv ");        ] pop popEnv ");
         /* sm1("(hoge) message show_ring "); */
 }  }
   
 def SresolutionFrameWithTower(g,opt) {  def SresolutionFrameWithTower(g,opt) {
   local gbTower, ans, ff, count, startingGB, opts, skelton,withSkel, autof,    local gbTower, ans, ff, count, startingGB, opts, skelton,withSkel, autof,
         gbasis;          gbasis, nohomog,i,n;
     /* extern Sordinary */
     nohomog = false;
     count = -1;  Sordinary = false; /* default value for options. */
   if (Length(Arglist) >= 2) {    if (Length(Arglist) >= 2) {
     if (IsInteger(opt)) count = opt;      if (IsArray(opt)) {
   }else{        n = Length(opt);
     count = -1;        for (i=0; i<n; i++) {
           if (IsInteger(opt[i])) {
             count = opt[i];
           }
           if (IsString(opt[i])) {
             if (opt[i] == "homogenized") {
               nohomog = true;
             }else if (opt[i] == "Sordinary") {
               Sordinary = true;
             }else{
               Println("Warning: unknown option");
               Println(opt);
             }
           }
         }
       } else if (IsNull(opt)){
       } else {
         Println("Warning: option should be given by an array.");
         Println(opt);
         Println("--------------------------------------------");
       }
   }    }
   
   sm1(" setupEnvForResolution ");    sm1(" setupEnvForResolution ");
Line 148  def SresolutionFrameWithTower(g,opt) {
Line 256  def SresolutionFrameWithTower(g,opt) {
   */    */
   
   sm1(" (mmLarger) (matrix) switch_function ");    sm1(" (mmLarger) (matrix) switch_function ");
   g = Map(g,"Shomogenize");    if (! nohomog) {
       Println("Automatic homogenization.");
       g = Map(g,"Shomogenize");
     }else{
       Println("No automatic homogenization.");
     }
   if (SonAutoReduce) {    if (SonAutoReduce) {
     sm1("[ (AutoReduce) ] system_variable /autof set ");      sm1("[ (AutoReduce) ] system_variable /autof set ");
     sm1("[ (AutoReduce) 1 ] system_variable ");      sm1("[ (AutoReduce) 1 ] system_variable ");
Line 188  def SresolutionFrameWithTower(g,opt) {
Line 301  def SresolutionFrameWithTower(g,opt) {
 }  }
 HelpAdd(["SresolutionFrameWithTower",  HelpAdd(["SresolutionFrameWithTower",
 ["It returs [resolution of the initial, gbTower, skelton, gbasis]",  ["It returs [resolution of the initial, gbTower, skelton, gbasis]",
    "option: \"homogenized\" (no automatic homogenization) ",
  "Example: Sweyl(\"x,y\");",   "Example: Sweyl(\"x,y\");",
  "         a=SresolutionFrameWithTower([x^3,x*y,y^3-1]);"]]);   "         a=SresolutionFrameWithTower([x^3,x*y,y^3-1]);"]]);
   
 def SresolutionFrame(f,opt) {  def SresolutionFrame(f,opt) {
   local ans;    local ans;
   ans = SresolutionFrameWithTower(f);    ans = SresolutionFrameWithTower(f,opt);
   return(ans[0]);    return(ans[0]);
 }  }
 /* ---------------------------- */  /* ---------------------------- */
Line 287  def Sres0FrameWithSkelton(g) {
Line 401  def Sres0FrameWithSkelton(g) {
   
   
 def StotalDegree(f) {  def StotalDegree(f) {
   sm1(" [(grade) f] gbext (universalNumber) dc /FunctionValue set ");    local d0;
     sm1(" [(grade) f] gbext (universalNumber) dc /d0 set ");
     /* Print("degree of "); Print(f); Print(" is "); Println(d0); */
     return(d0);
 }  }
   
   HelpAdd(["Sord_w",
   ["Sord_w(f,w) returns the w-order of f",
    "Example: Sord_w(x^2*Dx*Dy,[x,-1,Dx,1]):"]]);
 /* Sord_w(x^2*Dx*Dy,[x,-1,Dx,1]); */  /* Sord_w(x^2*Dx*Dy,[x,-1,Dx,1]); */
 def Sord_w(f,w) {  def Sord_w(f,w) {
   local neww,i,n;    local neww,i,n;
Line 336  def test_SinitOfArray() {
Line 456  def test_SinitOfArray() {
   
 /* f is assumed to be a monomial with toes. */  /* f is assumed to be a monomial with toes. */
 def Sdegree(f,tower,level) {  def Sdegree(f,tower,level) {
   local i;    local i,ww, wd;
     /* extern WeightOfSweyl; */
     ww = WeightOfSweyl;
   f = Init(f);    f = Init(f);
   if (level <= 1) return(StotalDegree(f));    if (level <= 1) return(StotalDegree(f));
   i = Degree(f,es);    i = Degree(f,es);
   return(StotalDegree(f)+Sdegree(tower[level-2,i],tower,level-1));    return(StotalDegree(f)+Sdegree(tower[level-2,i],tower,level-1));
   
 }  }
   
 def SgenerateTable(tower) {  def SgenerateTable(tower) {
   local height, n,i,j, ans, ans_at_each_floor;    local height, n,i,j, ans, ans_at_each_floor;
   
     /*
     Print("SgenerateTable: tower=");Println(tower);
     sm1(" print_switch_status "); */
   height = Length(tower);    height = Length(tower);
   ans = NewArray(height);    ans = NewArray(height);
   for (i=0; i<height; i++) {    for (i=0; i<height; i++) {
     n = Length(tower[i]);      n = Length(tower[i]);
     ans_at_each_floor=NewArray(n);      ans_at_each_floor=NewArray(n);
     for (j=0; j<n; j++) {      for (j=0; j<n; j++) {
       ans_at_each_floor[j] = Sdegree(tower[i,j],tower,i+1)-(i+1);        ans_at_each_floor[j] = Sdegree(tower[i,j],tower,i+1)-(i+1)
                               + OFFSET;
       /* Println([i,j,ans_at_each_floor[j]]); */        /* Println([i,j,ans_at_each_floor[j]]); */
     }      }
     ans[i] = ans_at_each_floor;      ans[i] = ans_at_each_floor;
Line 419  def SmaxOfStrategy(a) {
Line 547  def SmaxOfStrategy(a) {
 }  }
   
   
 def SlaScala(g) {  def SlaScala(g,opt) {
   local rf, tower, reductionTable, skel, redundantTable, bases,    local rf, tower, reductionTable, skel, redundantTable, bases,
         strategy, maxOfStrategy, height, level, n, i,          strategy, maxOfStrategy, height, level, n, i,
         freeRes,place, f, reducer,pos, redundant_seq,bettiTable,freeResV,ww,          freeRes,place, f, reducer,pos, redundant_seq,bettiTable,freeResV,ww,
Line 427  def SlaScala(g) {
Line 555  def SlaScala(g) {
         reductionTable_tmp;          reductionTable_tmp;
   /* extern WeightOfSweyl; */    /* extern WeightOfSweyl; */
   ww = WeightOfSweyl;    ww = WeightOfSweyl;
   Print("WeghtOfSweyl="); Println(WeightOfSweyl);    Print("WeightOfSweyl="); Println(WeightOfSweyl);
   rf = SresolutionFrameWithTower(g);    rf = SresolutionFrameWithTower(g,opt);
     Print("rf="); sm1_pmat(rf);
   redundant_seq = 1;   redundant_seq_ordinary = 1;    redundant_seq = 1;   redundant_seq_ordinary = 1;
   tower = rf[1];    tower = rf[1];
   
     Println("Generating reduction table which gives an order of reduction.");
     Print("WeghtOfSweyl="); Println(WeightOfSweyl);
     Print("tower"); Println(tower);
   reductionTable = SgenerateTable(tower);    reductionTable = SgenerateTable(tower);
     Print("reductionTable="); sm1_pmat(reductionTable);
   
   skel = rf[2];    skel = rf[2];
   redundantTable = SnewArrayOfFormat(rf[1]);    redundantTable = SnewArrayOfFormat(rf[1]);
   redundantTable_ordinary = SnewArrayOfFormat(rf[1]);    redundantTable_ordinary = SnewArrayOfFormat(rf[1]);
Line 452  def SlaScala(g) {
Line 587  def SlaScala(g) {
         Println([level,i]);          Println([level,i]);
         reductionTable_tmp[i] = -200000;          reductionTable_tmp[i] = -200000;
         if (reductionTable[level,i] == strategy) {          if (reductionTable[level,i] == strategy) {
            Print("Processing "); Print([level,i]);             Print("Processing [level,i]= "); Print([level,i]);
            Print("   Strategy = "); Println(strategy);             Print("   Strategy = "); Println(strategy);
            if (level == 0) {             if (level == 0) {
              if (IsNull(redundantTable[level,i])) {               if (IsNull(redundantTable[level,i])) {
Line 473  def SlaScala(g) {
Line 608  def SlaScala(g) {
                   place = f[3];                    place = f[3];
                   /* (level-1, place) is the place for f[0],                    /* (level-1, place) is the place for f[0],
                      which is a newly obtained  GB. */                       which is a newly obtained  GB. */
 #ifdef ORDINARY  if (Sordinary) {
                   redundantTable[level-1,place] = redundant_seq;                    redundantTable[level-1,place] = redundant_seq;
                   redundant_seq++;                    redundant_seq++;
 #else  }else{
                   if (f[4] > f[5]) {                    if (f[4] > f[5]) {
                     /* Zero in the gr-module */                      /* Zero in the gr-module */
                     Print("v-degree of [org,remainder] = ");                      Print("v-degree of [org,remainder] = ");
Line 487  def SlaScala(g) {
Line 622  def SlaScala(g) {
                     redundantTable[level-1,place] = redundant_seq;                      redundantTable[level-1,place] = redundant_seq;
                     redundant_seq++;                      redundant_seq++;
                   }                    }
 #endif  }
                   redundantTable_ordinary[level-1,place]                    redundantTable_ordinary[level-1,place]
                      =redundant_seq_ordinary;                       =redundant_seq_ordinary;
                   redundant_seq_ordinary++;                    redundant_seq_ordinary++;
Line 520  def SlaScala(g) {
Line 655  def SlaScala(g) {
     bases = Sbases_to_vec(bases,bettiTable[i]);      bases = Sbases_to_vec(bases,bettiTable[i]);
     freeResV[i] = bases;      freeResV[i] = bases;
   }    }
   return([freeResV, redundantTable,reducer,bettiTable,redundantTable_ordinary]);    return([freeResV, redundantTable,reducer,bettiTable,redundantTable_ordinary,rf]);
 }  }
   
 def SthereIs(reductionTable_tmp,strategy) {  def SthereIs(reductionTable_tmp,strategy) {
Line 613  def SunitOfFormat(pos,forms) {
Line 748  def SunitOfFormat(pos,forms) {
   return(ans);    return(ans);
 }  }
   
 def Error(s) {  
   sm1(" s error ");  
 }  
   
 def IsNull(s) {  
   if (Stag(s) == 0) return(true);  
   else return(false);  
 }  
   
 def StowerOf(tower,level) {  def StowerOf(tower,level) {
   local ans,i;    local ans,i;
   ans = [ ];    ans = [ ];
Line 642  def Sspolynomial(f,g) {
Line 769  def Sspolynomial(f,g) {
   sm1("f g spol /FunctionValue set");    sm1("f g spol /FunctionValue set");
 }  }
   
 def MonomialPart(f) {  
   sm1(" [(lmonom) f] gbext /FunctionValue set ");  
 }  
   
   /* WARNING:
     When you use SwhereInTower, you have to change gbList
     as below. Ofcourse, you should restrore the gbList
     SsetTower(StowerOf(tower,level));
     pos = SwhereInTower(syzHead,tower[level]);
   */
 def SwhereInTower(f,tower) {  def SwhereInTower(f,tower) {
   local i,n,p,q;    local i,n,p,q;
   if (f == Poly("0")) return(-1);    if (f == Poly("0")) return(-1);
Line 682  def SpairAndReduction(skel,level,ii,freeRes,tower,ww) 
Line 812  def SpairAndReduction(skel,level,ii,freeRes,tower,ww) 
   
   tower2 = StowerOf(tower,level-1);    tower2 = StowerOf(tower,level-1);
   SsetTower(tower2);    SsetTower(tower2);
     Println(["level=",level]);
     Println(["tower2=",tower2]);
   /** sm1(" show_ring ");   */    /** sm1(" show_ring ");   */
   
   gi = Stoes_vec(bases[i]);    gi = Stoes_vec(bases[i]);
Line 715  def SpairAndReduction(skel,level,ii,freeRes,tower,ww) 
Line 847  def SpairAndReduction(skel,level,ii,freeRes,tower,ww) 
   sj = sj*tmp[1]+t_syz[j];    sj = sj*tmp[1]+t_syz[j];
   t_syz[i] = si;    t_syz[i] = si;
   t_syz[j] = sj;    t_syz[j] = sj;
   
     SsetTower(StowerOf(tower,level));
   pos = SwhereInTower(syzHead,tower[level]);    pos = SwhereInTower(syzHead,tower[level]);
   
     SsetTower(StowerOf(tower,level-1));
   pos2 = SwhereInTower(tmp[0],tower[level-1]);    pos2 = SwhereInTower(tmp[0],tower[level-1]);
   ans = [tmp[0],t_syz,pos,pos2,vdeg,vdeg_reduced];    ans = [tmp[0],t_syz,pos,pos2,vdeg,vdeg_reduced];
   /* pos is the place to put syzygy at level. */    /* pos is the place to put syzygy at level. */
Line 753  def Sreduction(f,myset) {
Line 889  def Sreduction(f,myset) {
   return([tmp[0],tmp[1],t_syz]);    return([tmp[0],tmp[1],t_syz]);
 }  }
   
 def Warning(s) {  
   Print("Warning: ");  
   Println(s);  
 }  
 def RingOf(f) {  
   local r;  
   if (IsPolynomial(f)) {  
     if (f != Poly("0")) {  
       sm1(f," (ring) dc /r set ");  
     }else{  
       sm1(" [(CurrentRingp)] system_variable /r set ");  
     }  
   }else{  
     Warning("RingOf(f): the argument f must be a polynomial. Return the current ring.");  
     sm1(" [(CurrentRingp)] system_variable /r set ");  
   }  
   return(r);  
 }  
   
 def Sfrom_es(f,size) {  def Sfrom_es(f,size) {
   local c,ans, i, d, myes, myee, j,n,r,ans2;    local c,ans, i, d, myes, myee, j,n,r,ans2;
Line 828  def Sbases_to_vec(bases,size) {
Line 946  def Sbases_to_vec(bases,size) {
   return(newbases);    return(newbases);
 }  }
   
 def Sminimal(g) {  HelpAdd(["Sminimal",
   ["It constructs the V-minimal free resolution by LaScala's algorithm",
    "option: \"homogenized\" (no automatic homogenization ",
    "      : \"Sordinary\"   (no (u,v)-minimal resolution)",
    "Options should be given as an array.",
    "Example:  Sweyl(\"x,y\",[[\"x\",-1,\"y\",-1,\"Dx\",1,\"Dy\",1]]);",
    "          v=[[2*x*Dx + 3*y*Dy+6, 0],",
    "             [3*x^2*Dy + 2*y*Dx, 0],",
    "             [0,  x^2+y^2],",
    "             [0,  x*y]];",
    "         a=Sminimal(v);",
    "         Sweyl(\"x,y\",[[\"x\",-1,\"y\",-1,\"Dx\",1,\"Dy\",1]]);",
    "         b = ReParse(a[0]); sm1_pmat(b); ",
    "         IsExact_h(b,[x,y]):",
    "Note:  a[0] is the V-minimal resolution. a[3] is the Schreyer resolution."]]);
   
   def Sminimal(g,opt) {
   local r, freeRes, redundantTable, reducer, maxLevel,    local r, freeRes, redundantTable, reducer, maxLevel,
         minRes, seq, maxSeq, level, betti, q, bases, dr,          minRes, seq, maxSeq, level, betti, q, bases, dr,
         betti_levelplus, newbases, i, j,qq;          betti_levelplus, newbases, i, j,qq, tminRes,bettiTable, ansSminimal;
   r = SlaScala(g);    if (Length(Arglist) < 2) {
        opt = null;
     }
     /* Sordinary is set in SlaScala(g,opt) --> SresolutionFrameWithTower */
   
     ScheckIfSchreyer("Sminimal:0");
     r = SlaScala(g,opt);
   /* Should I turn off the tower?? */    /* Should I turn off the tower?? */
     ScheckIfSchreyer("Sminimal:1");
   freeRes = r[0];    freeRes = r[0];
   redundantTable = r[1];    redundantTable = r[1];
   reducer = r[2];    reducer = r[2];
     bettiTable = SbettiTable(redundantTable);
     Println("Betti numbers are ------");
     sm1_pmat(bettiTable);
   minRes = SnewArrayOfFormat(freeRes);    minRes = SnewArrayOfFormat(freeRes);
   seq = 0;    seq = 0;
   maxSeq = SgetMaxSeq(redundantTable);    maxSeq = SgetMaxSeq(redundantTable);
Line 889  def Sminimal(g) {
Line 1033  def Sminimal(g) {
       }        }
     }      }
    }     }
    return([Stetris(minRes,redundantTable),     tminRes = Stetris(minRes,redundantTable);
           [ minRes, redundantTable, reducer,r[3],r[4]],r[0]]);     ansSminimal = [SpruneZeroRow(tminRes), tminRes,
                     [ minRes, redundantTable, reducer,r[3],r[4]],r[0],r[5]];
      Println("------------ Note -----------------------------");
      Println("To get shift vectors, use Reparse and SgetShifts(resmat,w)");
      Println("To get initial of the complex, use Reparse and Sinit_w(resmat,w)");
      Println("0: minimal resolution, 3: Schreyer resolution ");
      Println("------------ Resolution Summary  --------------");
      Print("Betti numbers : ");
      Println(Map(ansSminimal[0],"Length"));
      Print("Betti numbers of the Schreyer frame: ");
      Println(Map(ansSminimal[3],"Length"));
      Println("-----------------------------------------------");
   
      sm1(" restoreEnvAfterResolution ");
   
      return(ansSminimal);
   /* r[4] is the redundantTable_ordinary */    /* r[4] is the redundantTable_ordinary */
   /* r[0] is the freeResolution */    /* r[0] is the freeResolution */
     /* r[5] is the skelton */
 }  }
   
   
Line 1024  def Sannfs(f,v) {
Line 1184  def Sannfs(f,v) {
 def Sannfs2(f) {  def Sannfs2(f) {
   local p,pp;    local p,pp;
   p = Sannfs(f,"x,y");    p = Sannfs(f,"x,y");
 /*    sm1(" p 0 get { [(x) (y) (Dx) (Dy)] laplace0 } map /p set ");
   Sweyl("x,y",[["x",1,"y",1,"Dx",1,"Dy",1,"h",1],    Sweyl("x,y",[["x",-1,"y",-1,"Dx",1,"Dy",1]]);
                ["x",-1,"y",-1,"Dx",1,"Dy",1]]); */    pp = Map(p,"Spoly");
   Sweyl("x,y",[["x",-1,"y",-1,"Dx",1,"Dy",1]]);    return(Sminimal(pp));
   pp = Map(p[0],"Spoly");  
   return(Sminimal(pp));  
 }  }
   
   HelpAdd(["Sannfs2",
   ["Sannfs2(f) constructs the V-minimal free resolution for the weight (-1,1)",
    "of the Laplace transform of the annihilating ideal of the polynomial f in x,y.",
    "See also Sminimal, Sannfs3.",
    "Example: a=Sannfs2(\"x^3-y^2\");",
    "         b=a[0]; sm1_pmat(b);",
    "         b[1]*b[0]:",
    "Example: a=Sannfs2(\"x*y*(x-y)*(x+y)\");",
    "         b=a[0]; sm1_pmat(b);",
    "         b[1]*b[0]:"
   ]]);
   /* Some samples.
     The betti numbers of most examples are 2,1. (0-th and 1-th).
     a=Sannfs2("x*y*(x+y-1)"); ==> The betti numbers are 3, 2.
     a=Sannfs2("x^3-y^2-x");
     a=Sannfs2("x*y*(x-y)");
   */
   
   
 def Sannfs3(f) {  def Sannfs3(f) {
   local p,pp;    local p,pp;
   p = Sannfs(f,"x,y,z");    p = Sannfs(f,"x,y,z");
     sm1(" p 0 get { [(x) (y) (z) (Dx) (Dy) (Dz)] laplace0 } map /p set ");
   Sweyl("x,y,z",[["x",-1,"y",-1,"z",-1,"Dx",1,"Dy",1,"Dz",1]]);    Sweyl("x,y,z",[["x",-1,"y",-1,"z",-1,"Dx",1,"Dy",1,"Dz",1]]);
   pp = Map(p[0],"Spoly");    pp = Map(p,"Spoly");
   return(Sminimal(pp));    return(Sminimal(pp));
 }  }
   
 /*  HelpAdd(["Sannfs3",
   The betti numbers of most examples are 2,1. (0-th and 1-th).  ["Sannfs3(f) constructs the V-minimal free resolution for the weight (-1,1)",
   a=Sannfs2("x*y*(x+y-1)"); ==> The betti numbers are 3, 2.   "of the Laplace transform of the annihilating ideal of the polynomial f in x,y,z.",
   a=Sannfs2("x^3-y^2-x");    : it causes an error. It should be fixed.   "See also Sminimal, Sannfs2.",
   a=Sannfs2("x*y*(x-y)");    : it causes an error. It should be fixed.   "Example: a=Sannfs3(\"x^3-y^2*z^2\");",
    "         b=a[0]; sm1_pmat(b);",
    "         b[1]*b[0]: b[2]*b[1]:"]]);
   
 */  
   
   
   /* Sannfs2("x*y*(x-y)*(x+y)"); is a test problem */
   /* x y (x+y-1)(x-2),  x^3-y^2, x^3 - y^2 z^2,
      x y z (x+y+z-1) seems to be interesting, because the first syzygy
     contains 1.
   */
   
 /*  The below is under construction. */  def CopyArray(m) {
 def Sschreyer(g) {    local ans,i,n;
   local rf, tower, reductionTable, skel, redundantTable, bases,    if (IsArray(m)) {
         strategy, maxOfStrategy, height, level, n, i,       n = Length(m);
         freeRes,place, f, reducer,pos, redundant_seq,bettiTable,freeResV,ww,       ans = NewArray(n);
         redundantTable_ordinary, redundant_seq_ordinary,       for (i=0; i<n; i++) {
         reductionTable_tmp,c2,ii,nn;         ans[i] = CopyArray(m[i]);
   /* extern WeightOfSweyl; */       }
   ww = WeightOfSweyl;       return(ans);
   Print("WeghtOfSweyl="); Println(WeightOfSweyl);    }else{
   rf = SresolutionFrameWithTower(g);       return(m);
   redundant_seq = 1;   redundant_seq_ordinary = 1;    }
   tower = rf[1];  }
   reductionTable = SgenerateTable(tower);  HelpAdd(["CopyArray",
   skel = rf[2];  ["It duplicates the argument array recursively.",
   redundantTable = SnewArrayOfFormat(rf[1]);   "Example: m=[1,[2,3]];",
   redundantTable_ordinary = SnewArrayOfFormat(rf[1]);   "         a=CopyArray(m); a[1] = \"Hello\";",
   reducer = SnewArrayOfFormat(rf[1]);   "         Println(m); Println(a);"]]);
   freeRes = SnewArrayOfFormat(rf[1]);  
   bettiTable = SsetBettiTable(rf[1],g);  
   
   height = Length(reductionTable);  def IsZeroVector(m) {
   for (level = 0; level < height; level++) {    local n,i;
       n = Length(reductionTable[level]);    n = Length(m);
       for (i=0; i<n; i++) {    for (i=0; i<n; i++) {
            Println([level,i]);      if (!IsZero(m[i])) {
            Print("Processing "); Print([level,i]);        return(false);
            if (level == 0) {      }
              if (IsNull(redundantTable[level,i])) {    }
                bases = freeRes[level];    return(true);
                /* Println(["At floor : GB=",i,bases,tower[0,i]]); */  }
                pos = SwhereInGB(tower[0,i],rf[3,0]);  
                bases[i] = rf[3,0,pos];  
                /* redundantTable[level,i] = 0;  
                redundantTable_ordinary[level,i] = 0; */  
                freeRes[level] = bases;  
                /* Println(["GB=",i,bases,tower[0,i]]); */  
              }  
            }else{ /* level >= 1 */  
              if (IsNull(redundantTable[level,i])) {  
                bases = freeRes[level];  
                f = SpairAndReduction2(skel,level,i,freeRes,tower,  
                                       ww,redundantTable);  
                if (f[0] != Poly("0")) {  
                   place = f[3];  
                   /* (level-1, place) is the place for f[0],  
                      which is a newly obtained  GB. */  
 #ifdef ORDINARY  
                   redundantTable[level-1,place] = redundant_seq;  
                   redundant_seq++;  
 #else  
                   if (f[4] > f[5]) {  
                     /* Zero in the gr-module */  
                     Print("v-degree of [org,remainder] = ");  
                     Println([f[4],f[5]]);  
                     Print("[level,i] = "); Println([level,i]);  
                     redundantTable[level-1,place] = 0;  
                   }else{  
                     redundantTable[level-1,place] = redundant_seq;  
                     redundant_seq++;  
                   }  
 #endif  
                   redundantTable_ordinary[level-1,place]  
                      =redundant_seq_ordinary;  
                   redundant_seq_ordinary++;  
                   bases[i] = SunitOfFormat(place,f[1])-f[1];  /* syzygy */  
                   /* redundantTable[level,i] = 0;  
                   redundantTable_ordinary[level,i] = 0; */  
                   /* i must be equal to f[2], I think. Double check. */  
   
                   /* Correction Of Constant */  def SpruneZeroRow(res) {
                   c2 = f[6];    local minRes, n,i,j,m, base,base2,newbase,newbase2, newMinRes;
                   nn = Length(bases);  
                   for (ii=0; ii<nn;ii++) {  
                      if (ii != place) {  
                        bases[ii] = bases[ii]*c2;  
                      }  
                   }  
   
                   freeRes[level] = bases;    minRes = CopyArray(res);
                   /* bases = freeRes[level-1];    n = Length(minRes);
                      bases[place] = f[0];  
                      freeRes[level-1] = bases;  It is already set. */  
                   reducer[level-1,place] = f[1];  
                }else{  
                   /* redundantTable[level,i] = 0; */  
                   bases = freeRes[level];  
                   bases[i] = f[1];  /* Put the syzygy. */  
                   freeRes[level] = bases;  
                }  
              }  /* end of level >= 1 */  
           }  
     } /* i loop */  
   } /* level loop */  
   n = Length(freeRes);  
   freeResV = SnewArrayOfFormat(freeRes);  
   for (i=0; i<n; i++) {    for (i=0; i<n; i++) {
     bases = freeRes[i];      base = minRes[i];
     bases = Sbases_to_vec(bases,bettiTable[i]);      m = Length(base);
     freeResV[i] = bases;      if (i != n-1) {
         base2 = minRes[i+1];
         base2 = Transpose(base2);
       }
       newbase = [ ];
       newbase2 = [ ];
       for (j=0; j<m; j++) {
         if (!IsZeroVector(base[j])) {
           newbase = Append(newbase,base[j]);
           if (i != n-1) {
             newbase2 = Append(newbase2,base2[j]);
           }
         }
       }
       minRes[i] = newbase;
       if (i != n-1) {
         if (newbase2 == [ ]) {
           minRes[i+1] = [ ];
         }else{
           minRes[i+1] = Transpose(newbase2);
         }
       }
   }    }
   return([freeResV, redundantTable,reducer,bettiTable,redundantTable_ordinary]);  
     newMinRes = [ ];
     n = Length(minRes);
     i = 0;
     while (i < n ) {
       base = minRes[i];
       if (base == [ ]) {
         i = n; /* break; */
       }else{
         newMinRes = Append(newMinRes,base);
       }
       i++;
     }
     return(newMinRes);
 }  }
   
 def SpairAndReduction2(skel,level,ii,freeRes,tower,ww,redundantTable) {  def testAnnfs2(f) {
   local i, j, myindex, p, bases, tower2, gi, gj,    local a,i,n;
        si, sj, tmp, t_syz, pos, ans, ssp, syzHead,pos2,    a = Sannfs2(f);
        vdeg,vdeg_reduced,n,c2;    b=a[0];
   Println("SpairAndReduction2:");    n = Length(b);
     Println("------ V-minimal free resolution -----");
     sm1_pmat(b);
     Println("----- Is it complex?  ---------------");
     for (i=0; i<n-1; i++) {
       Println(b[i+1]*b[i]);
     }
     return(a);
   }
   def testAnnfs3(f) {
     local a,i,n;
     a = Sannfs3(f);
     b=a[0];
     n = Length(b);
     Println("------ V-minimal free resolution -----");
     sm1_pmat(b);
     Println("----- Is it complex?  ---------------");
     for (i=0; i<n-1; i++) {
       Println(b[i+1]*b[i]);
     }
     return(a);
   }
   
   if (level < 1) Error("level should be >= 1 in SpairAndReduction.");  def ToString_array(p) {
   p = skel[level,ii];    local ans;
   myindex = p[0];    if (IsArray(p)) {
   i = myindex[0]; j = myindex[1];      ans = Map(p,"ToString_array");
   bases = freeRes[level-1];    }else{
   Println(["p and bases ",p,bases]);      ans = ToString(p);
   if (IsNull(bases[i]) || IsNull(bases[j])) {  
     Println([level,i,j,bases[i],bases[j]]);  
     Error("level, i, j : bases[i], bases[j]  must not be NULL.");  
   }    }
     return(ans);
   }
   
   tower2 = StowerOf(tower,level-1);  /* sm1_res_div([[x],[y]],[[x^2],[x*y],[y^2]],[x,y]): */
   SsetTower(tower2);  
   /** sm1(" show_ring ");   */  
   
   gi = Stoes_vec(bases[i]);  def sm1_res_div(I,J,V) {
   gj = Stoes_vec(bases[j]);    I = ToString_array(I);
     J = ToString_array(J);
     V = ToString_array(V);
     sm1(" [[ I J]  V ] res*div /FunctionValue set ");
   }
   
   ssp = Sspolynomial(gi,gj);  /* It has not yet been working */
   si = ssp[0,0];  def sm1_res_kernel_image(m,n,v) {
   sj = ssp[0,1];    m = ToString_array(m);
   syzHead = si*es^i;    n = ToString_array(n);
   /* This will be the head term, I think. But, double check. */    v = ToString_array(v);
   Println([si*es^i,sj*es^j]);    sm1(" [m n v] res-kernel-image /FunctionValue set ");
   }
   def Skernel(m,v) {
     m = ToString_array(m);
     v = ToString_array(v);
     sm1(" [ m v ] syz /FunctionValue set ");
   }
   
   Print("[gi, gj] = "); Println([gi,gj]);  
   sm1(" [(Homogenize)] system_variable message ");  
   Print("Reduce the element "); Println(si*gi+sj*gj);  
   Print("by  "); Println(bases);  
   
   tmp = Sreduction(si*gi+sj*gj, bases);  def sm1_gb(f,v) {
     f =ToString_array(f);
     v = ToString_array(v);
     sm1(" [f v] gb /FunctionValue set ");
   }
   
   Print("result is "); Println(tmp);  
   t_syz = tmp[2];  
   si = si*tmp[1]+t_syz[i];  
   sj = sj*tmp[1]+t_syz[j];  
   t_syz[i] = si;  
   t_syz[j] = sj;  
   
   c2 = null;  def SisComplex(a) {
   /* tmp[0] must be zero */    local n,i,j,k,b,p,q;
   n = Length(t_syz);    n = Length(a);
   for (i=0; i<n; i++) {    for (i=0; i<n-1; i++) {
      if (IsConstant(t_syz[i])) {      if (Length(a[i+1]) != 0) {
        if (IsNull(redundantTable[level-1,i])) {        b = a[i+1]*a[i];
          /* i must equal to pos2 below. */        p = Length(b); q = Length(b[0]);
          c2 = -t_syz[i];        for (j=0; j<p; j++) {
          tmp[0] = freeRes[level-1,i];          for (k=0; k<q; k++) {
          t_syz[i] = 0;            if (!IsZero(b[j,k])) {
          /* break; does not work. Use */               Print("Is is not complex at ");
          i = n;               Println([i,j,k]);
        }               return(false);
      }            }
           }
         }
       }
   }    }
     return(true);
   }
   
   /* This is essential part for V-minimal resolution. */  def IsExact_h(c,v) {
   /* vdeg = SvDegree(si*gi+sj*gj,tower,level-1,ww); */    local a;
   vdeg = SvDegree(si*gi,tower,level-1,ww);    v = ToString_array(v);
   vdeg_reduced = SvDegree(tmp[0],tower,level-1,ww);    a = [c,v];
   Print("vdegree of the original = "); Println(vdeg);    sm1(a," isExact_h /FunctionValue set ");
   Print("vdegree of the remainder = "); Println(vdeg_reduced);  }
   HelpAdd(["IsExact_h",
   ["IsExact_h(complex,var): bool",
    "It checks the given complex is exact or not in D<h> (homogenized Weyl algebra)",
    "cf. ReParse"
   ]]);
   
   pos = SwhereInTower(syzHead,tower[level]);  def IsSameIdeal_h(ii,jj,v) {
   pos2 = SwhereInTower(tmp[0],tower[level-1]);    local a;
   ans = [tmp[0],t_syz,pos,pos2,vdeg,vdeg_reduced,c2];    v = ToString_array(v);
   /* pos is the place to put syzygy at level. */    a = [ii,jj,v];
   /* pos2 is the place to put a new GB at level-1. */    sm1(a," isSameIdeal_h /FunctionValue set ");
   Println(ans);  }
   HelpAdd(["IsSameIdeal_h",
   ["IsSameIdeal_h(ii,jj,var): bool",
    "It checks the given ideals are the same or not in D<h> (homogenized Weyl algebra)",
    "cf. ReParse"
   ]]);
   
   def ReParse(a) {
     local c;
     if (IsArray(a)) {
       c = Map(a,"ReParse");
     }else{
       sm1(a," toString . /c set");
     }
     return(c);
   }
   HelpAdd(["ReParse",
   ["Reparse(obj): obj",
    "It parses the given object in the current ring.",
    "Outputs from SlaScala, Sschreyer may cause a trouble in other functions,",
    "because it uses the Schreyer order.",
    "In this case, ReParse the outputs from these functions.",
    "cf. IsExaxt_h"
   ]]);
   
   def ScheckIfSchreyer(s) {
     local ss;
     sm1(" (report) (grade) switch_function /ss set ");
     if (ss != "module1v") {
        Print("ScheckIfSchreyer: from "); Println(s);
        Error("grade is not module1v");
     }
     /*
     sm1(" (report) (mmLarger) switch_function /ss set ");
     if (ss != "tower") {
        Print("ScheckIfSchreyer: from "); Println(s);
        Error("mmLarger is not tower");
     }
     */
     sm1(" [(Schreyer)] system_variable (universalNumber) dc /ss set ");
     if (ss != 1) {
        Print("ScheckIfSchreyer: from "); Println(s);
        Error("Schreyer order is not set.");
     }
     /* More check will be necessary. */
     return(true);
   }
   
   def SgetShift(mat,w,m) {
     local omat;
     sm1(" mat { w m ord_w<m> {(universalNumber) dc}map } map /omat set");
     return(Map(omat,"Max"));
   }
   HelpAdd(["SgetShift",
   ["SgetShift(mat,w,m) returns the shift vector of mat with respect to w with the shift m.",
    "Note that the order of the ring and the weight w must be the same.",
    "Example:  Sweyl(\"x,y\",[[\"x\",-1,\"Dx\",1]]); ",
    "          SgetShift([[x*Dx+1,Dx^2+x^5],[Poly(\"0\"),x],[x,x]],[\"x\",-1,\"Dx\",1],[2,0]):"]]);
   
   def SgetShifts(resmat,w) {
     local i,n,ans,m0;
     n = Length(resmat);
     ans = NewArray(n);
     m0 = NewArray(Length(resmat[0,0]));
     ans[0] = m0;
     for (i=0; i<n-1; i++) {
       ans[i+1] = SgetShift(resmat[i],w,m0);
       m0 = ans[i+1];
     }
     return(ans);
   }
   HelpAdd(["SgetShifts",
   ["SgetShifts(resmat,w) returns the shift vectors of the resolution resmat",
    " with respect to w with the shift m.",
    "Note that the order of the ring and the weight w must be the same.",
    "Zero row is not allowed.",
    "Example:   a=Sannfs2(\"x^3-y^2\");",
    "           b=a[0]; w = [\"x\",-1,\"y\",-1,\"Dx\",1,\"Dy\",1];",
    "           Sweyl(\"x,y\",[w]); b = Reparse(b);",
    "           SgetShifts(b,w):"]]);
   
   def Sinit_w(resmat,w) {
     local shifts,ans,n,i,m,mat,j;
     shifts = SgetShifts(resmat,w);
     n = Length(resmat);
     ans = NewArray(n);
     for (i=0; i<n; i++) {
       m = shifts[i];
       mat = ScopyArray(resmat[i]);
       for (j=0; j<Length(mat); j++) {
         mat[j] = Init_w_m(mat[j],w,m);
       }
       ans[i] = mat;
     }
     return(ans);
   }
   HelpAdd(["Sinit_w",
   ["Sinit_w(resmat,w) returns the initial of the complex resmat with respect to the weight w.",
    "Example:   a=Sannfs2(\"x^3-y^2\");",
    "           b=a[0]; w = [\"x\",-1,\"y\",-1,\"Dx\",1,\"Dy\",1];",
    "           Sweyl(\"x,y\",[w]); b = Reparse(b);",
    "           c=Sinit_w(b,w); c:"
   ]]);
   
   /* This method does not work, because we have zero rows.
      Think about it later. */
   def SbettiTable(rtable) {
     local ans,i,j,pp;
     ans = SnewArrayOfFormat(rtable);
     for (i=0; i<Length(rtable); i++) {
       pp = 0;
       for (j=0; j<Length(rtable[i]); j++) {
          if (rtable[i,j] != 0) {pp = pp+1;}
       }
       ans[i] = pp;
     }
   return(ans);    return(ans);
 }  }
   
Legend:
Removed from v.1.5  
changed lines
  Added in v.1.25
FreeBSD-CVSweb <freebsd-cvsweb@FreeBSD.org>