source: git/Singular/LIB/teachstd.lib @ 2c957af

spielwiese
Last change on this file since 2c957af was 2c957af, checked in by Hans Schönemann <hannes@…>, 18 years ago
*hannes: Stanislav Bulygins changes git-svn-id: file:///usr/local/Singular/svn/trunk@9345 2c84dea3-7e68-4137-9b89-c4e89433aadc
  • Property mode set to 100644
File size: 22.8 KB
Line 
1// $Id:
2//GMG, last modified 28.9.01
3///////////////////////////////////////////////////////////////////////////////
4version="$Id: teachstd.lib,v 1.6 2006-07-25 12:02:44 Singular Exp $";
5category="Teaching";
6info="
7LIBRARY:  teachstd.lib   Procedures for teaching standard bases
8AUTHOR:                  G.-M. Greuel, greuel@mathematik.uni-kl.de
9
10NOTE: The library is intended to be used for teaching purposes, but not
11      for serious computations. Sufficiently high printlevel allows to
12      control each step, thus illustrating the algorithms at work.
13      The procedures are implemented exactly as described in the book
14      'A SINGULAR Introduction to Commutative Algebra' by G.-M. Greuel and
15      G. Pfister (Springer 2002).
16
17PROCEDURES:
18 ecart(f);              ecart of f
19 tail(f);               tail of f
20 sameComponent(f,g);    test for same module component of lead(f) and lead(g)
21 leadmonomial(f);       leading monomial as poly (also for vectors)
22 monomialLcm(m,n);      lcm of monomials m and n as poly (also for vectors)
23 spoly(f[,1]);          s-polynomial of f [symmetric form]
24 minEcart(T,h);         element g from T of minimal ecart s.t. LM(g)|LM(h)
25 NFMora(i);             normal form of i w.r.t Mora algorithm
26 prodcrit(f,g);         test for product criterion
27 chaincrit(f,g,h);      test for chain criterion
28 pairset(G);            pairs form G neither satifying prodcrit nor chaincrit
29 updatePairs(P,S,h);    pairset P enlarded by not useless pairs (h,f), f in S
30 standard(id);          standard basis of ideal/module
31 localstd(id);          local standard basis of id using Lazard's method
32
33              [parameters in square brackets are optional]
34";
35
36LIB "poly.lib";
37
38///////////////////////////////////////////////////////////////////////////////
39proc ecart(f)
40"USAGE:   ecart(f); f poly or vector
41RETURN:  the ecart e of f of type int
42EXAMPLE: example ecart; shows an example
43"
44{
45  int e = maxdeg1(f)-maxdeg1(lead(f));
46  return(e);
47}
48example
49{ "EXAMPLE:"; echo = 2;
50   ring r=0,(x,y,z),ls;
51   ecart((y+z+x+xyz)**2);
52   ring s=0,(x,y,z),dp;
53   ecart((y+z+x+xyz)**2);
54}
55///////////////////////////////////////////////////////////////////////////////
56proc leadmonomial(f)
57"USAGE:   leadmonomial(f); f poly or vector
58RETURN:  the leading monomial of f of type poly
59NOTE:    if f is of type poly, leadmonomial(f)=leadmonom(f), if f is of type
60         vector and if leadmonom(f)=m*gen(i) then leadmonomial(f)=m
61EXAMPLE: example leadmonomial; shows an example
62"
63{
64  int e;
65  poly m;
66  if(typeof(f) == "vector")
67  {
68    e=leadexp(f)[nvars(basering)+1];
69    m=leadmonom(f)[e,1];
70  }
71  if(typeof(f) == "poly")
72  {
73    m=leadmonom(f);
74  }
75  return(m);
76}
77example
78{ "EXAMPLE:"; echo = 2;
79   ring s=0,(x,y,z),(c,dp);
80   leadmonomial((y+z+x+xyz)**2);
81   leadmonomial([(y+z+x+xyz)**2,xyz5]);
82}
83///////////////////////////////////////////////////////////////////////////////
84proc tail(f)
85"USAGE:   tail(f); f poly or vector
86RETURN:  f-lead(f), the tail of f of type poly
87EXAMPLE: example tail; shows an example
88"
89{
90  def t = f-lead(f);
91  return(t);
92}
93example
94{ "EXAMPLE:"; echo = 2;
95   ring r=0,(x,y,z),ls;
96   tail((y+z+x+xyz)**2);
97   ring s=0,(x,y,z),dp;
98   tail((y+z+x+xyz)**2);
99}
100///////////////////////////////////////////////////////////////////////////////
101proc sameComponent(f,g)
102"USAGE:   sameComponent(f,g); f,g poly or vector
103RETURN:  1 if f and g are of type poly or if f and g are of type vector and
104         their leading monomials involve the same module component,
105         0 if not
106EXAMPLE: example sameComponent; shows an example
107"
108{
109  if(typeof(f) != typeof(g))
110  {
111    ERROR("** arguments must be of same type");
112  }
113  if(typeof(f) == "vector")
114  {
115     if( leadexp(f)[nvars(basering)+1] !=  leadexp(g)[nvars(basering)+1] )
116     {
117       return(0);
118     }
119  }
120  return(1);
121}
122example
123{ "EXAMPLE:"; echo = 2;
124   ring r=0,(x,y,z),dp;
125   sameComponent([y+z+x,xyz],[z2,xyz]);
126   sameComponent([y+z+x,xyz],[z4,xyz]);
127   sameComponent(y+z+x+xyz, xy+z5);
128}
129///////////////////////////////////////////////////////////////////////////////
130proc monomialLcm(m,n)
131"USAGE:   monomialLcm(m,n); m,n of type poly or vector
132RETURN:  least common multiple of leading monomials of m and n, of type poly
133NOTE:    if m = (x1...xr)^(a1,...,ar)*gen(i) (gen(i)=1 if m is of type poly)
134         and n = (x1...xr)^(b1,...,br)*gen(j), then the proc returns
135         (x1,...,xr)^(max(a1,b1),...,max(ar,br)) if i=j and 0 if i!=j.
136EXAMPLE: example monomialLcm; shows an example
137"
138{
139  if(typeof(n) != typeof(m))
140  {
141    ERROR("** arguments must be of same type");
142  }
143
144  poly u ;
145  if(sameComponent(m,n) == 0)   //leading term of vectors involve
146  {                             //different module components
147     return(u);
148  }
149
150  intvec v = leadexp(m);        //now start to compute lcm
151  intvec w = leadexp(n);
152  u=1;
153  int i;
154  for (i=1; i<=nvars(basering); i++)
155  {
156      if(v[i]>=w[i])
157      {
158         u = u*var(i)**v[i];
159      }
160      else
161      {
162         u = u*var(i)**w[i];
163      }
164   }
165   return(u);
166}
167example
168{ "EXAMPLE:"; echo = 2;
169   ring r=0,(x,y,z),ds;
170   monomialLcm(xy2,yz3);
171   monomialLcm([xy2,xz],[yz3]);
172   monomialLcm([xy2,xz3],[yz3]);
173}
174///////////////////////////////////////////////////////////////////////////////
175proc spoly(f,g,list #)
176"USAGE:   spoly(f,g[,s]); f,g poly or vector, s int
177RETURN:  the s-polynomial of f and g, of type poly or vector
178         if s!=0 the symmetric s-polynomial (without division) is returned
179EXAMPLE: example spoly; shows an example
180"
181{
182  if(typeof(f) != typeof(g))
183  {
184    ERROR("** arguments must be of same type");
185  }
186  if(size(#) == 0)
187  {
188     #[1]=0;
189  }
190
191  int e;
192  poly o = monomialLcm(f,g);
193  if( o == 0)                    //can only happen, if vectors f and g involve
194  {                              //different module components
195    vector sp;
196    return(sp);
197  }
198
199  poly m=leadmonomial(f);         //compute the leading monomial as poly
200  poly n=leadmonomial(g);
201
202  if (#[1]==0)                    //the asymmetric s-poly
203  {
204    def sp = (o/m)*f - (leadcoef(f)/leadcoef(g))*(o/n)*g;
205  }
206  else                            //the symmetric s-poly, avoiding division
207  {
208    def sp = leadcoef(g)*(o/m)*f - leadcoef(f)*(o/n)*g;
209  }
210  return(sp);
211}
212example
213{ "EXAMPLE:"; echo = 2;
214   ring r=0,(x,y,z),ls;
215   spoly(2x2+x2y,3y3+xyz);
216   ring s=0,(x,y,z),(c,dp);
217   spoly(2x2+x2y,3y3+xyz);
218   spoly(2x2+x2y,3y3+xyz,1);             //symmetric s-poly without division
219   spoly([5x2+x2y,z5],[x2,y3,y4]);       //s-poly for vectors
220}
221///////////////////////////////////////////////////////////////////////////////
222proc minEcart(T,h)
223"USAGE:   minEcart(T,h); T ideal or module, h poly or vector
224RETURN:  element g from T such that leadmonom(g) divides leadmonom(h)
225         ecart(g) is minimal with this property (if T != 0);
226         return 0 if T is 0 or h = 0
227EXAMPLE: example minEcart; shows an example
228"
229{
230  int i,k,e;
231  if (size(T)==0 or h==0 )       //trivial cases
232  {
233     h = 0;
234     return(h);
235  }
236//---- check whether some element in T has the same module component as h ----
237  int v;
238  intvec w;
239  T = simplify(T,2);
240
241  if (typeof(h) == "vector")
242  {
243     e=1;
244     v = leadexp(h)[nvars(basering)+1];
245     for( i=1; i<=size(T); i++)
246     {
247       w[i]=leadexp(T[i])[nvars(basering)+1];
248       if(v == w[i])
249       {
250         e=0;         //some element in T involves the same component as h
251       }
252     }
253  }
254  if ( e == 1 )       //no element in T involves the same component as h
255  {
256     h = 0;
257     return(h);
258  }
259
260  if (typeof(h) == "poly")   //for polys v=w[i] for all i
261  {
262    v = 1;
263    w[size(T)]=0;
264    w=w+1;
265  }
266
267//------ check whether for some g in T leadmonom(g) divides leadmonom(h) -----
268  def g = T[1];
269  poly f = leadmonomial(h);
270
271  for( i=1; i<=size(T); i++)
272  {
273    if( f/leadmonomial(T[i]) != 0 and v==w[i] )
274    {
275       g=T[i];
276       k=i;
277       break;
278    }
279  }
280  if (k == 0)        //no leadmonom(g) divides leadmonom(h)
281  {
282     g = 0;
283     return(g);
284  }
285//--look for T[i] with minimal ecart s.t.leadmonom(T[i]) divides leadmonom(h)--
286
287  for( i=k+1; i<=size(T); i++)
288  {
289    if( f/leadmonomial(T[i]) != 0 and v==w[i] )
290    {
291       if (ecart(T[i]) < ecart(g))
292       {
293         g=T[i];
294        }
295     }
296  }
297  return(g);
298}
299example
300{ "EXAMPLE:"; echo = 2;
301   ring R=0,(x,y,z),dp;
302   ideal T = x2y+x2,y3+xyz,xyz2+z4;
303   poly h = x2y2z2+x5+yx3+z6;
304   minEcart(T,h);"";
305   ring S=0,(x,y,z),(c,ds);
306   module T = [x2+x2y,y2],[y3+xyz,x3-z3],[x3y+z4,0,x2];
307   vector h = [x3y+x5+x2y2z2+z6,x3];
308   minEcart(T,h);
309}
310///////////////////////////////////////////////////////////////////////////////
311proc NFMora(f,G,list #)
312"USAGE:   NFMora(f,G[,s]); f poly or vector,G ideal or module, s int
313RETURN:  the Mora normal form of f w.r.t. G, same type as f
314         if s!=0 the symmetric s-polynomial (without division) is used
315NOTE:    Show comments if printlevel > 0, pauses computation if printlevel > 1
316EXAMPLE: example NFMora; shows an example
317"
318{
319  if(size(#) == 0)
320  {
321     #[1]=0;
322  }
323
324  int y = printlevel - voice + 2;
325  if( f==0 or size(G) ==0 )
326  {
327     if (y>0)
328     {
329        "// 1st or 2nd argument 0";
330     }
331     return(f);
332  }
333
334  int i,e;
335  def h = f;
336  def T = G;
337// -------------------- start with f to be reduced by G --------------------
338  if (y>0)
339  {"";
340   "// Input for NFMora is (f,T):";
341   "// f:"; f;
342   "// T:"; T;
343   "// Reduce f with T, eventually enlarging T for local ordering";
344  }
345
346// ----------------------- enter the reduction loop ------------------------
347  def g = minEcart(T,h);
348  while (h!=0 and g!=0)
349  {
350     if (y>0)
351     {  "";
352        "//                Reduction-step in NFMora:",i;
353        "// h = (f after",i,"reductions) reduction with g from T:";
354        "// g = element of minimal ecart in T s.t. LM(g)|LM(h):";
355        "// h:";h;
356        "// g:";g;
357     }
358     if (y>1)
359     {
360        pause("press <return> to continue");
361        "// T, set used for reduction:"; T;
362        pause("press <return> to continue");
363     }
364     e=0;
365     if( ecart(g) > ecart(h) )
366     {
367        T=T,h;  e=1;
368     }
369     if (y>0 )
370     {
371        "// T-set enlarged for next reduction? (yes/no = 1/0):  ", e;
372        if( e==1 )
373        {
374          "// T-set for next reduction got enlarged by h:";
375          "// h:";h;
376          if (y>1)
377          {
378             pause("press <return> to continue");
379          }
380        }
381     }
382     h = spoly(h,g,#[1]);
383     g = minEcart(T,h);
384     i=i+1;
385  }
386  if(y>0)
387  { "";
388    "// normal form is:";
389  }
390  return(h);
391}
392example
393{ "EXAMPLE:"; echo = 2;
394   ring r=0,(x,y,z),dp;
395   poly f = x2y2z2+x5+yx3+z6-3y3;
396   ideal G = x2y+x2,y3+xyz,xyz2+z6;
397   NFMora(f,G);"";
398   ring s=0,(x,y,z),ds;
399   poly f = x3y+x5+x2y2z2+z6;
400   ideal G = x2+x2y,y3+xyz,x3y2+z4;
401   NFMora(f,G);"";
402   vector v = [f,x2+x2y];
403   module M = [x2+x2y,f],[y3+xyz,y3],[x3y2+z4,z2];
404   NFMora(v,M);
405}
406///////////////////////////////////////////////////////////////////////////////
407proc prodcrit(f,g)
408"USAGE:   prodcrit(f,g); f,g poly or vector
409RETURN:  1 if product criterion applies in the same module component,
410         2 if lead(f) and lead(g) involve different components, 0 else
411NOTE:    if product criterion applies we can delete (f,g) from pairset
412EXAMPLE: example prodcrit; shows an example
413"
414{
415  if(typeof(f)=="poly")    //product criterion for polynomials
416  {
417    if( monomialLcm(f,g)==leadmonom(f)*leadmonom(g))
418    {
419      return(1);
420    }
421    return(0);
422  }
423// ------------------- product criterion for modules ---------------------
424  if(sameComponent(f,g)==1)
425  {
426     if( monomialLcm(f,g)==leadmonomial(f)*leadmonomial(g) )
427     {
428         int c = leadexp(f)[nvars(basering)+1];  //component involving lead(f)
429         if((f-f[c]*gen(c))-(g-g[c]*gen(c))==0)  //other components are 0
430         {
431           return(1);
432         }
433     }
434     return(0);
435  }
436  return(2);
437}
438example
439{ "EXAMPLE:"; echo = 2;
440   ring r=0,(x,y,z),dp;
441   poly f = y3z3+x5+yx3+z6;
442   poly g = x5+yx3;
443   prodcrit(f,g);
444   vector v = x3z2*gen(1)+x3y*gen(1)+x2y*gen(2);
445   vector w = y4*gen(1)+y3*gen(2)+xyz*gen(1);
446   prodcrit(v,w);
447}
448///////////////////////////////////////////////////////////////////////////////
449proc chaincrit(f,g,h)
450"USAGE:   chaincrit(f,g,h); f,g,h poly or module
451RETURN:  1 if chain criterion applies, 0 else
452NOTE:    if chain criterion applies to f,g,h we can delete (g,h) from pairset
453EXAMPLE: example chaincrit; shows an example
454"
455{
456  if(sameComponent(f,g) and sameComponent(f,h))
457  {
458    if( monomialLcm(g,h)/leadmonomial(f) !=0 )
459    {
460      return(1);
461    }
462  }
463  return(0);
464}
465example
466{ "EXAMPLE:"; echo = 2;
467   ring r=0,(x,y,z),dp;
468   poly f = x2y2z2+x5+yx3+z6;
469   poly g = x5+yx3;
470   poly h = y2z5+x5+yx3;
471   chaincrit(f,g,h);
472   vector u = [x2y3-z2,x2y];
473   vector v = [x2y2+z2,x2-y2];
474   vector w = [x2y4+z3,x2+y2];
475   chaincrit(u,v,w);
476}
477///////////////////////////////////////////////////////////////////////////////
478proc pairset(G)
479"USAGE:   pairset(G); G ideal or module
480RETURN:  list L,
481         L[1] = the pairset of G as list (not containing pairs for
482         which the product or the chain criterion applies)
483         L[2] = intvec v, v[1]= # product criterion, v[2]= # chain criterion
484EXAMPLE: example pairset; shows an example
485"
486{
487   int i,j,k,s,c,ccrit,pcrit,pr;
488   int y = printlevel - voice + 2;
489   G = simplify(G,10);
490   def g = G;
491   ideal pair;
492   list P,I;      //P=pairlist of G, I=list of corresponding indices of pairs
493   for (i=1; i<=size(G); i++)
494   {
495      for(j = i+1; j<=size(G); j++)
496      {
497         pr = prodcrit(G[i],G[j]);       //check first product criterion
498         if( pr != 0 )
499         {
500            pcrit=pcrit+(pr==1);
501         }
502         else
503         {
504            s = size(P);                //now check chain criterion
505            for(k=1; k<=s; k++)
506            {
507              if( i==I[k][2] )
508              {
509                if ( chaincrit(P[k][1],P[k][2],G[j]) )
510                {                     //need not include (G[i],G[j]) in P
511                  c=1; ccrit=ccrit+1;
512                  break;
513                }
514              }
515              if( j==I[k][1] and c==0 )
516              {
517                "########### enter pairset2 #############";
518                if ( chaincrit(G[i],P[k][1],P[k][2]) )
519                {                     //can delete P[k]=(P[k][1],P[k][2])
520                  ccrit=ccrit+1;
521                  P = delete(P,k);
522                  s = s-1;
523                }
524              }
525            }
526            if ( c==0 )
527            {
528               g = G[i],G[j];
529               P[s+1]=g;
530               I[s+1]=intvec(i,j);
531            }
532            c=0;
533         }
534      }
535   }
536   if (y>0)
537   {  "";
538      "// product criterion:",pcrit," chain criterion:",ccrit;
539   }
540   intvec v = pcrit,ccrit;
541   P=P,v;
542   return(P);
543}
544example
545{ "EXAMPLE:"; echo = 2;
546   ring r=0,(x,y,z),dp;
547   ideal G = x2y+x2,y3+xyz,xyz2+z4;
548   pairset(G);"";
549   module T = [x2y3-z2,x2y],[x2y2+z2,x2-y2],[x2y4+z3,x2+y2];
550   pairset(T);
551}
552///////////////////////////////////////////////////////////////////////////////
553proc updatePairs(P,S,h)
554"USAGE:   updatePairs(P,S,h); P list, S ideal or module, h poly or vector
555         P a list of pairs of polys or vectors (obtained from pairset)
556RETURN:  list Q,
557         Q[1] = the pairset P enlarged  by all pairs (f,h), f from S,
558         without pairs for which the product or the chain criterion applies
559         Q[2] = intvec v, v[1]= # product criterion, v[2]= # chain criterion
560EXAMPLE: example updatePairs; shows an example
561"
562{
563   int i,j,k,s,r,c,ccrit,pcrit,pr;
564   int y = printlevel - voice + 2;
565   ideal pair;
566   list Q = P;           //Q will become enlarged pairset
567   s = size(P);
568   r = size(Q);          //r will grow with Q
569   list R;
570   def g = S;            //give g the correct type ideal/module
571   for (i=1; i<=size(S); i++)
572   {
573      pr =  prodcrit(h,S[i]);
574      if( pr != 0 )     //check product criterion
575      {
576        pcrit=pcrit+(pr==1);       //count product criterion in same component
577      }
578      else
579      {                 //prodcrit did not apply, check for chaincrit
580         r=size(Q);
581         for(k=1; k<=r; k++)
582         {
583           if( Q[k][2]==S[i] )    //S[i]=Q[k][2]
584           {
585             if( chaincrit(Q[k][1],S[i],h) )
586             {                   //can forget (S[i],h)
587                c=1; ccrit=ccrit+1;
588                break;
589             }
590           }
591         }
592         if ( c==0 )
593         {
594            g = S[i],h;        //add pair (S[i],h)
595            Q[r+1] = g;
596         }
597         c=0;
598      }
599   }
600   if (y>0)
601   {  "";;
602      "// product criterion:",pcrit," chain criterion:",ccrit;
603   }
604   intvec v = pcrit,ccrit;
605   Q = Q,v;
606   return(Q);
607}
608example
609{ "EXAMPLE:"; echo = 2;
610   ring R1=0,(x,y,z),(c,dp);
611   ideal S = x2y+x2,y3+xyz;
612   poly h = x2y+xyz;
613   list P = pairset(S)[1];
614   P;"";
615   updatePairs(P,S,h);"";
616   module T = [x2y3-z2,x2y],[x2y4+z3,x2+y2];
617   P = pairset(T)[1];
618   P;"";
619   updatePairs(P,T,[x2+x2y,y3+xyz]);
620}
621///////////////////////////////////////////////////////////////////////////////
622proc standard(id, list #)
623"USAGE:   standard(i[,s]); id ideal or module, s int
624RETURN:  a standard basis of id, using generalized Mora's algorithm
625         which is Buchberger's algorithm for global monomial orderings.
626         If s!=0 the symmetric s-polynomial (without division) is used
627NOTE:    Show comments if printlevel > 0, pauses computation if printlevel > 1
628EXAMPLE: example standard; shows an example
629"
630{
631  if(size(#) == 0)
632  {
633     #[1]=0;
634  }
635
636   def S = id;                 //S will become the standard basis of id
637   def h = S[1];
638   int i,z;
639   int y = printlevel - voice + 2;
640   if(y>0)
641   {  "";
642      "// the set S, to become a standard basis:"; S;
643      if(y>1)
644      {
645         "// create pairset, i.e. pairs from S,";
646         "// after application of product and chain criterion";
647      }
648   }
649   list P = pairset(S);        //create pairset of S=id
650   intvec v = P[2];
651   P = P[1];
652//-------------------------- Main loop in SB lgorithm ----------------------
653   while (size(P) !=0)
654   {  z=z+1;
655      if(y>0)
656      {  "";
657         "//              Enter NFMora for next pair, count",z;
658         "// size of partial standard basis S:    (",size(S),")";
659         "// number of pairs of S after updating: (",size(P),")";
660         if(y>1)
661         {
662            "// 1st pair of new pairset:"; P[1];
663            "// set T=S used for reduction:";S;
664            "// apply NFMora to (spoly,S), spoly = spoly(1st pair)";
665         }
666      }
667      //-------------------- apply NFMora = Mora's normal form -------------
668      h = spoly(P[1][1],P[1][2],#[1]);
669      if(y>1)
670      {
671         "// spoly:";h;
672      }
673      h = NFMora(h,S,#[1]);
674      if(h==0)           //normal form is 0
675      {
676         if(y==1)
677         {
678            "// pair has reduced to 0";
679         }
680          if(y>1)
681         { h;"";
682           pause("press <return> to continue");
683         }
684     }
685      P = delete(P,1);   //spoly of pair reduced to 0, pair can be deleted
686      //--- spoly of pair did not reduce to 0, update S and paiset of S ----
687      if( h != 0)
688      {
689         if(y==1)
690         {
691            "// ** new spoly in degree **:", deg(h);
692         }
693         if(y>1)
694         { h;"";
695           pause("press <return> to continue");
696           "// update pairset";
697         }
698        P=updatePairs(P,S,h);    //update P (=paisert of S)
699        v=v+P[2];                //with useful pairs (g,h), g from S
700         P=P[1];
701         S=S,h;                  //update S, will become the standard basis
702      }
703   }
704//------------------------------ finished ---------------------------------
705   if( find(option(),"prot") or y>0 )
706   {  "";                //show how often prodcrit and chaincrit applied
707      "// product criterion:",v[1]," chain criterion:",v[2];
708      "";
709      "// Final standard basis:";
710   }
711   return(S);
712}
713example
714{ "EXAMPLE:"; echo = 2;
715   ring r=0,(x,y,z),dp;
716   ideal G = x2y+x2,y3+xyz,xyz2+z4;
717   standard(G);"";
718   ring s=0,(x,y,z),(c,ds);
719   ideal G = 2x2+x2y,y3+xyz,3x3y+z4;
720   standard(G);"";
721   standard(G,1);"";           //use symmetric s-poly without division
722   module M = [2x2,x3y+z4],[3y3+xyz,y3],[5z4,z2];
723   standard(M);
724}
725///////////////////////////////////////////////////////////////////////////////
726proc localstd (id)
727"USAGE:   localstd (id); id = ideal
728RETURN:  A standard basis for a local degree ordering, using Lazard's method.
729NOTE:    The procedure homogenizes id w.r.t. a new 1st variable local@t,
730         computes a SB w.r.t. (dp(1),dp) and substitutes local@t by 1.
731         Hence the result is a SB with respect to an ordering which sorts
732         first w.r.t. the subdegree of the original variables and then refines
733         it with dp. This is the local degree ordering ds.
734         localstd may be used in order to avoid cancellation of units and thus
735         to be able to use option(contentSB) also for local orderings.
736EXAMPLE: example localstd; shows an example
737"
738{
739  int ii;
740  def bas = basering;
741  execute("ring  @r_locstd
742     =("+charstr(bas)+"),(local@t,"+varstr(bas)+"),(dp(1),dp);");
743  ideal id = imap(bas,id);
744  ideal hid = homog(id,local@t);
745  hid = std(hid);
746  hid = subst(hid,local@t,1);
747  setring bas;
748  def hid = imap(@r_locstd,hid);
749  attrib(hid,"isSB",1);
750  kill @r_locstd;
751  return(hid);
752}
753example
754{ "EXAMPLE:"; echo = 2;
755   ring R = 0,(x,y,z),ds;
756   ideal i  = xyz+z5,2x2+y3+z7,3z5+y5;
757   localstd(i);
758}
759///////////////////////////////////////////////////////////////////////////////
760
761/*
762// some examples:
763LIB"teachstd.lib";
764option(prot); printlevel=3;
765ring r0 = 0,(t,x,y),lp;
766ideal i = x-t2,y-t3;
767standard(i);
768
769printlevel=1;
770standard(i);
771
772option(prot); printlevel =1;
773ring r1 = (0,a,b),(x,y,z),(c,ds);
774module M = [ax2,bx3y+z4],[a3y3+xyz,by3],[5az4,(a+b)*z2];
775module N1= std(M);
776module N2 = standard(M,1);
777NF(lead(N2),lead(N1));
778NF(lead(N1),lead(N2));rom T
779ring r2 = 0,(x,y,z),dp;
780ideal I = x2y+x2,y3+xyz,xyz2+z4;
781option(prot);
782int tt = timer;
783ideal J = standard(I);
784timer -tt;                 //4sec, product criterion: 9  chain criterion: 6
785ideal J1 = std(I);
786NF(lead(J),lead(J1));
787NF(lead(J1),lead(J));
788
789ring r3 = 0,(x,y,z),ds;
790poly f = x3*y4+z5+xyz;
791ideal I = f,jacob(f);
792option(prot);
793int tt = timer;
794ideal J = standard(I);
795timer -tt;                 //3sec, product criterion: 1  chain criterion: 3
796ideal J1 = std(I);
797NF(lead(J),lead(J1));
798NF(lead(J1),lead(J));
799
800//Becker:
801ring r4 = 32003,(x,y,z),lp;
802ideal I = x3-1, y3-1,
803-27x3-243x2y+27x2z-729xy2+162xyz-9xz2-729y3+243y2z-27yz2+z3-27;
804option(prot);
805int tt = timer;
806ideal J = standard(I);
807timer -tt;                 //201sec, product criterion: 42  chain criterion: 33
808ideal J1 = std(I);
809NF(lead(J),lead(J1));
810NF(lead(J1),lead(J));
811
812//Alex
813ring r5 = 32003,(x,y,z,t),dp;
814ideal I =
8152t3xy2z+x2ty+2x2y,
8162tz+y3x2t+z2t3y2x;
817option(prot); printlevel =1;
818ideal J1 = std(I);
819int tt = timer;
820ideal J = standard(I);
821timer -tt;                 //15sec product criterion: 0  chain criterion: 12
822NF(lead(J),lead(J1));
823NF(lead(J1),lead(J));
824
825ring r6 = 32003,(x,y,z,t),dp; //is already SB for ds, for dp too long
826ideal I=
8279x8+y7t3z4+5x4y2t2+2xy2z3t2,
8289y8+7xy6t+2x5y4t2+2x2yz3t2,
8299z8+3x2y3z2t4;
830option(prot);
831int tt = timer;
832ideal J = standard(I);
833timer -tt;                 //0sec, product criterion: 3  chain criterion: 0
834ideal J1 = std(I);
835NF(lead(J),lead(J1));
836NF(lead(J1),lead(J));
837
838
839*/
840
Note: See TracBrowser for help on using the repository browser.