source: git/Singular/LIB/modstd.lib @ 4454fb7

spielwiese
Last change on this file since 4454fb7 was 4454fb7, checked in by Andreas Steenpass <steenpass@…>, 8 years ago
allow modules in Modstd::reduce_parallel
  • Property mode set to 100644
File size: 18.9 KB
Line 
1///////////////////////////////////////////////////////////////////////////////
2version="version modstd.lib 4.0.0.0 May_2014 "; // $Id$
3category="Commutative Algebra";
4info="
5LIBRARY:  modstd.lib      Groebner bases of ideals using modular methods
6
7AUTHORS:  A. Hashemi      Amir.Hashemi@lip6.fr
8          G. Pfister      pfister@mathematik.uni-kl.de
9          H. Schoenemann  hannes@mathematik.uni-kl.de
10          A. Steenpass    steenpass@mathematik.uni-kl.de
11          S. Steidel      steidel@mathematik.uni-kl.de
12
13OVERVIEW:
14  A library for computing Groebner bases of ideals in the polynomial ring over
15  the rational numbers using modular methods.
16
17REFERENCES:
18  E. A. Arnold: Modular algorithms for computing Groebner bases.
19  J. Symb. Comp. 35, 403-419 (2003).
20
21  N. Idrees, G. Pfister, S. Steidel: Parallelization of Modular Algorithms.
22  J. Symb. Comp. 46, 672-684 (2011).
23
24PROCEDURES:
25  modStd(I);    standard basis of I using modular methods
26";
27
28LIB "poly.lib";
29LIB "modular.lib";
30
31proc modStd(ideal I, list #)
32"USAGE:   modStd(I[, exactness]); I ideal, exactness int
33RETURN:   a standard basis of I
34NOTE:     The procedure computes a standard basis of I (over the rational
35          numbers) by using modular methods.
36       @* An optional parameter 'exactness' can be provided.
37          If exactness = 1, the procedure computes a standard basis of I for
38          sure; if exactness = 0, it computes a standard basis of I
39          with high probability.
40SEE ALSO: modular
41EXAMPLE:  example modStd; shows an example"
42{
43    /* read optional parameter */
44    int exactness = 1;
45    if (size(#) > 0) {
46        /* For compatibility, we only test size(#) > 4. This can be changed to
47         * size(#) > 1 in the future. */
48        if (size(#) > 4 || typeof(#[1]) != "int") {
49            ERROR("wrong optional parameter");
50        }
51        exactness = #[1];
52    }
53
54    /* save options */
55    intvec opt = option(get);
56    option(redSB);
57
58    /* choose the right command */
59    string command = "groebner";
60    if (npars(basering) > 0) {
61        command = "Modstd::groebner_norm";
62    }
63
64    /* call modular() */
65    if (exactness) {
66        I = modular(command, list(I), primeTest_std,
67            deleteUnluckyPrimes_std, pTest_std, finalTest_std);
68    }
69    else {
70        I = modular(command, list(I), primeTest_std,
71            deleteUnluckyPrimes_std, pTest_std);
72    }
73
74    /* return the result */
75    attrib(I, "isSB", 1);
76    option(set, opt);
77    return(I);
78}
79example
80{
81    "EXAMPLE:";
82    echo = 2;
83    ring R1 = 0, (x,y,z,t), dp;
84    ideal I = 3x3+x2+1, 11y5+y3+2, 5z4+z2+4;
85    ideal J = modStd(I);
86    J;
87    I = homog(I, t);
88    J = modStd(I);
89    J;
90
91    ring R2 = 0, (x,y,z), ds;
92    ideal I = jacob(x5+y6+z7+xyz);
93    ideal J = modStd(I, 0);
94    J;
95
96    ring R3 = 0, x(1..4), lp;
97    ideal I = cyclic(4);
98    ideal J1 = modStd(I, 1);   // default
99    ideal J2 = modStd(I, 0);
100    size(reduce(J1, J2));
101    size(reduce(J2, J1));
102}
103
104/* compute a normalized GB via groebner() */
105static proc groebner_norm(ideal I)
106{
107    I = simplify(groebner(I), 1);
108    attrib(I, "isSB", 1);
109    return(I);
110}
111
112/* test if the prime p is suitable for the input, i.e. it does not divide
113 * the numerator or denominator of any of the coefficients */
114static proc primeTest_std(int p, alias list args)
115{
116    /* erase zero generators */
117    ideal I = simplify(args[1], 2);
118
119    /* clear denominators and count the terms */
120    ideal J;
121    ideal K;
122    int n = ncols(I);
123    intvec sizes;
124    number cnt;
125    int i;
126    for(i = n; i > 0; i--) {
127        J[i] = cleardenom(I[i]);
128        cnt = leadcoef(J[i])/leadcoef(I[i]);
129        K[i] = numerator(cnt)*var(1)+denominator(cnt);
130    }
131    sizes = size(J[1..n]);
132
133    /* change to characteristic p */
134    def br = basering;
135    list lbr = ringlist(br);
136    if (typeof(lbr[1]) == "int") {
137        lbr[1] = p;
138    }
139    else {
140        lbr[1][1] = p;
141    }
142    def rp = ring(lbr);
143    setring(rp);
144    ideal Jp = fetch(br, J);
145    ideal Kp = fetch(br, K);
146
147    /* test if any coefficient is missing */
148    if (intvec(size(Kp[1..n])) != 2:n) {
149        setring(br);
150        return(0);
151    }
152    if (intvec(size(Jp[1..n])) != sizes) {
153        setring(br);
154        return(0);
155    }
156    setring(br);
157    return(1);
158}
159
160/* find entries in modresults which come from unlucky primes.
161 * For this, sort the entries into categories depending on their leading
162 * ideal and return the indices in all but the biggest category. */
163static proc deleteUnluckyPrimes_std(alias list modresults)
164{
165    int size_modresults = size(modresults);
166
167    /* sort results into categories.
168     * each category is represented by three entries:
169     * - the corresponding leading ideal
170     * - the number of elements
171     * - the indices of the elements
172     */
173    list cat;
174    int size_cat;
175    ideal L;
176    int i;
177    int j;
178    for (i = 1; i <= size_modresults; i++) {
179        L = lead(modresults[i]);
180        attrib(L, "isSB", 1);
181        for (j = 1; j <= size_cat; j++) {
182            if (size(L) == size(cat[j][1])
183                && size(reduce(L, cat[j][1])) == 0
184                && size(reduce(cat[j][1], L)) == 0) {
185                cat[j][2] = cat[j][2]+1;
186                cat[j][3][cat[j][2]] = i;
187                break;
188            }
189        }
190        if (j > size_cat) {
191            size_cat++;
192            cat[size_cat] = list();
193            cat[size_cat][1] = L;
194            cat[size_cat][2] = 1;
195            cat[size_cat][3] = list(i);
196        }
197    }
198
199    /* find the biggest categories */
200    int cat_max = 1;
201    int max = cat[1][2];
202    for (i = 2; i <= size_cat; i++) {
203        if (cat[i][2] > max) {
204            cat_max = i;
205            max = cat[i][2];
206        }
207    }
208
209    /* return all other indices */
210    list unluckyIndices;
211    for (i = 1; i <= size_cat; i++) {
212        if (i != cat_max) {
213            unluckyIndices = unluckyIndices + cat[i][3];
214        }
215    }
216    return(unluckyIndices);
217}
218
219/* test if 'command' applied to 'args' in characteristic p is the same as
220   'result' mapped to characteristic p */
221static proc pTest_std(string command, list args, ideal result, int p)
222{
223    /* change to characteristic p */
224    def br = basering;
225    list lbr = ringlist(br);
226    if (typeof(lbr[1]) == "int") {
227        lbr[1] = p;
228    }
229    else {
230        lbr[1][1] = p;
231    }
232    def rp = ring(lbr);
233    setring(rp);
234    ideal Ip = fetch(br, args)[1];
235    ideal Gp = fetch(br, result);
236    attrib(Gp, "isSB", 1);
237
238    /* test if Ip is in Gp */
239    int i;
240    for (i = ncols(Ip); i > 0; i--) {
241        if (reduce(Ip[i], Gp, 1) != 0) {
242            setring(br);
243            return(0);
244        }
245    }
246
247    /* compute command(args) */
248    execute("Ip = "+command+"(Ip);");
249
250    /* test if Gp is in Ip */
251    for (i = ncols(Gp); i > 0; i--) {
252        if (reduce(Gp[i], Ip, 1) != 0) {
253            setring(br);
254            return(0);
255        }
256    }
257    setring(br);
258    return(1);
259}
260
261/* test if 'result' is a GB of the input ideal */
262static proc finalTest_std(string command, alias list args, ideal result)
263{
264    /* test if args[1] is in result */
265    attrib(result, "isSB", 1);
266    int i;
267    for (i = ncols(args[1]); i > 0; i--) {
268        if (reduce(args[1][i], result, 1) != 0) {
269            return(0);
270        }
271    }
272
273    /* test if result is a GB */
274    ideal G = std(result);
275    if (reduce_parallel(G, result)) {
276        return(0);
277    }
278    return(1);
279}
280
281/* return 1, if I_reduce is _not_ in G_reduce,
282 *        0, otherwise
283 * (same as size(reduce(I_reduce, G_reduce))).
284 * Uses parallelization. */
285static proc reduce_parallel(def I_reduce, def G_reduce)
286{
287    exportto(Modstd, I_reduce);
288    exportto(Modstd, G_reduce);
289    int size_I = ncols(I_reduce);
290    int chunks = Modular::par_range(size_I);
291    intvec range;
292    int i;
293    for (i = chunks; i > 0; i--) {
294        range = Modular::par_range(size_I, i);
295        task t(i) = "Modstd::reduce_task", list(range);
296    }
297    startTasks(t(1..chunks));
298    waitAllTasks(t(1..chunks));
299    int result = 0;
300    for (i = chunks; i > 0; i--) {
301        if (getResult(t(i))) {
302            result = 1;
303            break;
304        }
305    }
306    kill I_reduce;
307    kill G_reduce;
308    return(result);
309}
310
311/* compute a chunk of reductions for reduce_parallel */
312static proc reduce_task(intvec range)
313{
314    int result = 0;
315    int i;
316    for (i = range[1]; i <= range[2]; i++) {
317        if (reduce(I_reduce[i], G_reduce, 1) != 0) {
318            result = 1;
319            break;
320        }
321    }
322    return(result);
323}
324
325////////////////////////////////////////////////////////////////////////////////
326/*
327 * The following procedures are kept for backward compatibility with the old
328 * version of modstd.lib. As of now (May 2014), they are still needed in
329 * modnormal.lib, modwalk.lib, and symodstd.lib. They can be removed here as
330 * soon as they are not longer needed in these libraries.
331 */
332
333LIB "parallel.lib";
334
335static proc mod_init()
336{
337   newstruct("idealPrimeTest", "ideal Ideal");
338}
339
340static proc redFork(ideal I, ideal J, int n)
341{
342   attrib(J,"isSB",1);
343   return(reduce(I,J,1));
344}
345
346proc isIncluded(ideal I, ideal J, list #)
347"USAGE:  isIncluded(I,J); I,J ideals
348RETURN:  1 if J includes I,
349@*       0 if there is an element f in I which does not reduce to 0 w.r.t. J.
350EXAMPLE: example isIncluded; shows an example
351"
352{
353   def R = basering;
354   setring R;
355
356   attrib(J,"isSB",1);
357   int i,j,k;
358
359   if(size(#) > 0)
360   {
361      int n = #[1];
362      if(n >= ncols(I)) { n = ncols(I); }
363      if(n > 1)
364      {
365         for(i = 1; i <= n - 1; i++)
366         {
367            //link l(i) = "MPtcp:fork";
368            link l(i) = "ssi:fork";
369            open(l(i));
370
371            write(l(i), quote(redFork(eval(I[ncols(I)-i]), eval(J), 1)));
372         }
373
374         int t = timer;
375         if(reduce(I[ncols(I)], J, 1) != 0)
376         {
377            for(i = 1; i <= n - 1; i++)
378            {
379               close(l(i));
380            }
381            return(0);
382         }
383         t = timer - t;
384         if(t > 60) { t = 60; }
385         int i_sleep = system("sh", "sleep "+string(t));
386
387         j = ncols(I) - n;
388
389         while(j >= 0)
390         {
391            for(i = 1; i <= n - 1; i++)
392            {
393               if(status(l(i), "read", "ready"))
394               {
395                  if(read(l(i)) != 0)
396                  {
397                     for(i = 1; i <= n - 1; i++)
398                     {
399                        close(l(i));
400                     }
401                     return(0);
402                  }
403                  else
404                  {
405                     if(j >= 1)
406                     {
407                        write(l(i), quote(redFork(eval(I[j]), eval(J), 1)));
408                        j--;
409                     }
410                     else
411                     {
412                        k++;
413                        close(l(i));
414                     }
415                  }
416               }
417            }
418            if(k == n - 1)
419            {
420               j--;
421            }
422            i_sleep = system("sh", "sleep "+string(t));
423         }
424         return(1);
425      }
426   }
427
428   for(i = ncols(I); i >= 1; i--)
429   {
430      if(reduce(I[i],J,1) != 0){ return(0); }
431   }
432   return(1);
433}
434example
435{ "EXAMPLE:"; echo = 2;
436   ring r=0,(x,y,z),dp;
437   ideal I = x+1,x+y+1;
438   ideal J = x+1,y;
439   isIncluded(I,J);
440   isIncluded(J,I);
441   isIncluded(I,J,4);
442
443   ring R = 0, x(1..5), dp;
444   ideal I1 = cyclic(4);
445   ideal I2 = I1,x(5)^2;
446   isIncluded(I1,I2,4);
447}
448
449proc deleteUnluckyPrimes(list T, list L, int ho, list #)
450"USAGE:  deleteUnluckyPrimes(T,L,ho,#); T/L list of polys/primes, ho integer
451RETURN:  lists T,L(,M),lT with T/L(/M) list of polys/primes(/type of #),
452         lT ideal
453NOTE:    - if ho = 1, the polynomials in T are homogeneous, else ho = 0,
454@*       - lT is prevalent, i.e. the most appearing leading ideal in T
455EXAMPLE: example deleteUnluckyPrimes; shows an example
456"
457{
458   ho = ((ho)||(ord_test(basering) == -1));
459   int j,k,c;
460   intvec hl,hc;
461   ideal cT,lT,cK;
462   lT = lead(T[size(T)]);
463   attrib(lT,"isSB",1);
464   if(!ho)
465   {
466      for(j = 1; j < size(T); j++)
467      {
468         cT = lead(T[j]);
469         attrib(cT,"isSB",1);
470         if((size(reduce(cT,lT))!=0)||(size(reduce(lT,cT))!=0))
471         {
472            cK = cT;
473            c++;
474         }
475      }
476      if(c > size(T) div 2){ lT = cK; }
477   }
478   else
479   {
480      hl = hilb(lT,1);
481      for(j = 1; j < size(T); j++)
482      {
483         cT = lead(T[j]);
484         attrib(cT,"isSB",1);
485         hc = hilb(cT,1);
486         if(hl == hc)
487         {
488            for(k = 1; k <= size(lT); k++)
489            {
490               if(lT[k] < cT[k]) { lT = cT; c++; break; }
491               if(lT[k] > cT[k]) { c++; break; }
492            }
493         }
494         else
495         {
496            if(hc < hl){ lT = cT; hl = hilb(lT,1); c++; }
497         }
498      }
499   }
500
501   int addList;
502   if(size(#) > 0) { list M = #; addList = 1; }
503   j = 1;
504   attrib(lT,"isSB",1);
505   while((j <= size(T))&&(c > 0))
506   {
507      cT = lead(T[j]);
508      attrib(cT,"isSB",1);
509      if((size(reduce(cT,lT)) != 0)||(size(reduce(lT,cT)) != 0))
510      {
511         T = delete(T,j);
512         if(j == 1)
513         {
514            L = L[2..size(L)];
515            if(addList == 1) { M = M[2..size(M)]; }
516         }
517         else
518         {
519            if(j == size(L))
520            {
521               L = L[1..size(L)-1];
522               if(addList == 1) { M = M[1..size(M)-1]; }
523            }
524            else
525            {
526               L = L[1..j-1],L[j+1..size(L)];
527               if(addList == 1) { M = M[1..j-1],M[j+1..size(M)]; }
528            }
529         }
530         j--;
531      }
532      j++;
533   }
534
535   for(j = 1; j <= size(L); j++)
536   {
537      L[j] = bigint(L[j]);
538   }
539
540   if(addList == 0) { return(list(T,L,lT)); }
541   if(addList == 1) { return(list(T,L,M,lT)); }
542}
543example
544{ "EXAMPLE:"; echo = 2;
545   list L = 2,3,5,7,11;
546   ring r = 0,(y,x),Dp;
547   ideal I1 = 2y2x,y6;
548   ideal I2 = yx2,y3x,x5,y6;
549   ideal I3 = y2x,x3y,x5,y6;
550   ideal I4 = y2x,11x3y,x5;
551   ideal I5 = y2x,yx3,x5,7y6;
552   list T = I1,I2,I3,I4,I5;
553   deleteUnluckyPrimes(T,L,1);
554   list P = poly(x),poly(x2),poly(x3),poly(x4),poly(x5);
555   deleteUnluckyPrimes(T,L,1,P);
556}
557
558proc primeTest(def II, bigint p)
559{
560   if(typeof(II) == "string")
561   {
562      ideal I = `II`.Ideal;
563   }
564   else
565   {
566      ideal I = II;
567   }
568
569   I = simplify(I, 2);   // erase zero generators
570
571   int i,j;
572   poly f;
573   number cnt;
574   for(i = 1; i <= size(I); i++)
575   {
576      f = cleardenom(I[i]);
577      if(f == 0) { return(0); }
578      cnt = leadcoef(I[i])/leadcoef(f);
579      if((bigint(numerator(cnt)) mod p) == 0) { return(0); }
580      if((bigint(denominator(cnt)) mod p) == 0) { return(0); }
581      for(j = size(f); j > 0; j--)
582      {
583         if((bigint(leadcoef(f[j])) mod p) == 0) { return(0); }
584      }
585   }
586   return(1);
587}
588
589proc primeList(ideal I, int n, list #)
590"USAGE:  primeList(I,n[,ncores]); ( resp. primeList(I,n[,L,ncores]); ) I ideal,
591         n integer
592RETURN:  the intvec of n greatest primes <= 2147483647 (resp. n greatest primes
593         < L[size(L)] union with L) such that none of these primes divides any
594         coefficient occuring in I
595NOTE:    The number of cores to use can be defined by ncores, default is 1.
596EXAMPLE: example primeList; shows an example
597"
598{
599   intvec L;
600   int i,p;
601   int ncores = 1;
602
603//-----------------  Initialize optional parameter ncores  ---------------------
604   if(size(#) > 0)
605   {
606      if(size(#) == 1)
607      {
608         if(typeof(#[1]) == "int")
609         {
610            ncores = #[1];
611            # = list();
612         }
613      }
614      else
615      {
616         ncores = #[2];
617      }
618   }
619
620   if(size(#) == 0)
621   {
622      p = 2147483647;
623      while(!primeTest(I,p))
624      {
625         p = prime(p-1);
626         if(p == 2) { ERROR("no more primes"); }
627      }
628      L[1] = p;
629   }
630   else
631   {
632      L = #[1];
633      p = prime(L[size(L)]-1);
634      while(!primeTest(I,p))
635      {
636         p = prime(p-1);
637         if(p == 2) { ERROR("no more primes"); }
638      }
639      L[size(L)+1] = p;
640   }
641   if(p == 2) { ERROR("no more primes"); }
642   if(ncores == 1)
643   {
644      for(i = 2; i <= n; i++)
645      {
646         p = prime(p-1);
647         while(!primeTest(I,p))
648         {
649            p = prime(p-1);
650            if(p == 2) { ERROR("no more primes"); }
651         }
652         L[size(L)+1] = p;
653      }
654   }
655   else
656   {
657      int neededSize = size(L)+n-1;;
658      list parallelResults;
659      list arguments;
660      int neededPrimes = neededSize-size(L);
661      idealPrimeTest Id;
662      Id.Ideal = I;
663      export(Id);
664      while(neededPrimes > 0)
665      {
666         arguments = list();
667         for(i = ((neededPrimes div ncores)+1-(neededPrimes%ncores == 0))
668            *ncores; i > 0; i--)
669         {
670            p = prime(p-1);
671            if(p == 2) { ERROR("no more primes"); }
672            arguments[i] = list("Id", p);
673         }
674         parallelResults = parallelWaitAll("primeTest", arguments, 0, ncores);
675         for(i = size(arguments); i > 0; i--)
676         {
677            if(parallelResults[i])
678            {
679               L[size(L)+1] = arguments[i][2];
680            }
681         }
682         neededPrimes = neededSize-size(L);
683      }
684      kill Id;
685      if(size(L) > neededSize)
686      {
687         L = L[1..neededSize];
688      }
689   }
690   return(L);
691}
692example
693{ "EXAMPLE:"; echo = 2;
694   ring r = 0,(x,y,z),dp;
695   ideal I = 2147483647x+y, z-181;
696   intvec L = primeList(I,10);
697   size(L);
698   L[1];
699   L[size(L)];
700   L = primeList(I,5,L);
701   size(L);
702   L[size(L)];
703}
704
705////////////////////////////// further examples ////////////////////////////////
706
707/*
708ring r = 0, (x,y,z), lp;
709poly s1 = 5x3y2z+3y3x2z+7xy2z2;
710poly s2 = 3xy2z2+x5+11y2z2;
711poly s3 = 4xyz+7x3+12y3+1;
712poly s4 = 3x3-4y3+yz2;
713ideal i =  s1, s2, s3, s4;
714
715ring r = 0, (x,y,z), lp;
716poly s1 = 2xy4z2+x3y2z-x2y3z+2xyz2+7y3+7;
717poly s2 = 2x2y4z+x2yz2-xy2z2+2x2yz-12x+12y;
718poly s3 = 2y5z+x2y2z-xy3z-xy3+y4+2y2z;
719poly s4 = 3xy4z3+x2y2z-xy3z+4y3z2+3xyz3+4z2-x+y;
720ideal i =  s1, s2, s3, s4;
721
722ring r = 0, (x,y,z), lp;
723poly s1 = 8x2y2 + 5xy3 + 3x3z + x2yz;
724poly s2 = x5 + 2y3z2 + 13y2z3 + 5yz4;
725poly s3 = 8x3 + 12y3 + xz2 + 3;
726poly s4 = 7x2y4 + 18xy3z2 +  y3z3;
727ideal i = s1, s2, s3, s4;
728
729int n = 6;
730ring r = 0,(x(1..n)),lp;
731ideal i = cyclic(n);
732ring s = 0, (x(1..n),t), lp;
733ideal i = imap(r,i);
734i = homog(i,t);
735
736ring r = 0, (x(1..4),s), (dp(4),dp);
737poly s1 = 1 + s^2*x(1)*x(3) + s^8*x(2)*x(3) + s^19*x(1)*x(2)*x(4);
738poly s2 = x(1) + s^8 *x(1)* x(2)* x(3) + s^19* x(2)* x(4);
739poly s3 = x(2) + s^10*x(3)*x(4) + s^11*x(1)*x(4);
740poly s4 = x(3) + s^4*x(1)*x(2) + s^19*x(1)*x(3)*x(4) +s^24*x(2)*x(3)*x(4);
741poly s5 = x(4) + s^31* x(1)* x(2)* x(3)* x(4);
742ideal i =  s1, s2, s3, s4, s5;
743
744ring r = 0, (x,y,z), ds;
745int a = 16;
746int b = 15;
747int c = 4;
748int t = 1;
749poly f = x^a+y^b+z^(3*c)+x^(c+2)*y^(c-1)+x^(c-1)*y^(c-1)*z3
750         +x^(c-2)*y^c*(y2+t*x)^2;
751ideal i = jacob(f);
752
753ring r = 0, (x,y,z), ds;
754int a = 25;
755int b = 25;
756int c = 5;
757int t = 1;
758poly f = x^a+y^b+z^(3*c)+x^(c+2)*y^(c-1)+x^(c-1)*y^(c-1)*z3
759         +x^(c-2)*y^c*(y2+t*x)^2;
760ideal i = jacob(f),f;
761
762ring r = 0, (x,y,z), ds;
763int a = 10;
764poly f = xyz*(x+y+z)^2 +(x+y+z)^3 +x^a+y^a+z^a;
765ideal i = jacob(f);
766
767ring r = 0, (x,y,z), ds;
768int a = 6;
769int b = 8;
770int c = 10;
771int alpha = 5;
772int beta = 5;
773int t = 1;
774poly f = x^a+y^b+z^c+x^alpha*y^(beta-5)+x^(alpha-2)*y^(beta-3)
775         +x^(alpha-3)*y^(beta-4)*z^2+x^(alpha-4)*y^(beta-4)*(y^2+t*x)^2;
776ideal i = jacob(f);
777*/
778
Note: See TracBrowser for help on using the repository browser.