source: git/Singular/LIB/standard.lib @ 7cb9b4

spielwiese
Last change on this file since 7cb9b4 was 7cb9b4, checked in by Yue Ren <ren@…>, 9 years ago
chg: more examples for max and min
  • Property mode set to 100644
File size: 72.3 KB
Line 
1///////////////////////////////////////////////////////////////////////////
2version="version standard.lib 4.0.0.0 Jun_2013 "; // $Id$
3category="Miscellaneous";
4info="
5LIBRARY: standard.lib   Procedures which are always loaded at Start-up
6
7PROCEDURES:
8 stdfglm(ideal[,ord])   standard basis of ideal via fglm [and ordering ord]
9 stdhilb(ideal[,h])     Hilbert driven Groebner basis of ideal
10 groebner(ideal,...)    standard basis using a heuristically chosen method
11 res(ideal/module,[i])  free resolution of ideal or module
12 sprintf(fmt,...)       returns fomatted string
13 fprintf(link,fmt,..)   writes formatted string to link
14 printf(fmt,...)        displays formatted string
15 weightKB(stc,dd,vl)    degree dd part of a kbase w.r.t. some weigths
16 qslimgb(i)             computes a standard basis with slimgb in a qring
17 par2varRing([i])       create a ring making pars to vars, together with i
18 datetime()             return date and time as a string
19 max(i_1,...,i_k)       maximum of i_1, ..., i_k
20 min(i_1,...,i_k)       minimum of i_1, ..., i_k
21
22";
23//AUXILIARY PROCEDURES:
24// hilbRing([i])          ring for computing the (weighted) hilbert series
25// quotientList(L,...)    ringlist for creating a correct quotient ring
26
27//////////////////////////////////////////////////////////////////////////////
28
29proc stdfglm (ideal i, list #)
30"SYNTAX: @code{stdfglm (} ideal_expression @code{)} @*
31         @code{stdfglm (} ideal_expression@code{,} string_expression @code{)}
32TYPE:    ideal
33PURPOSE: computes the standard basis of the ideal in the basering
34         via @code{fglm} from the ordering given as the second argument
35         to the ordering of the basering. If no second argument is given,
36         \"dp\" is used. The standard basis for the given ordering (resp. for
37         \"dp\") is computed via the command groebner except if a further
38         argument \"std\" or \"slimgb\" is given in which case std resp.
39         slimgb is used.
40SEE ALSO: fglm, groebner, std, slimgb, stdhilb
41KEYWORDS: fglm
42EXAMPLE: example stdfglm; shows an example"
43{
44  string os;
45  int s = size(#);
46  def P= basering;
47  string algorithm;
48  int ii;
49  for( ii=1; ii<=s; ii++)
50  {
51    if ( typeof(#[ii])== "string" )
52    {
53      if ( #[ii]=="std" || #[ii]=="slimgb" )
54      {
55        algorithm =  #[ii];
56        # = delete(#,ii);
57        s--;
58        ii--;
59      }
60    }
61  }
62
63  if((s > 0) && (typeof(#[1]) == "string"))
64  {
65    os = #[1];
66    ideal Qideal = ideal(P);
67    int sQ = size(Qideal);
68    int sM = size(minpoly);
69    if ( sM!=0 )
70    {
71      string mpoly = string(minpoly);
72    }
73    if (sQ!=0 )
74    {
75      execute("ring Rfglm=("+charstr(P)+"),("+varstr(P)+"),"+os+";");
76      ideal Qideal = fetch(P,Qideal);
77      qring Pfglm = groebner(Qideal,"std","slimgb");
78    }
79    else
80    {
81      execute("ring Pfglm=("+charstr(P)+"),("+varstr(P)+"),"+os+";");
82    }
83    if ( sM!=0 )
84    {
85      execute("minpoly="+mpoly+";");
86    }
87  }
88  else
89  {
90    list BRlist = ringlist(P);
91    int nvarP = nvars(P);
92    intvec w;                       //for ringweights of basering P
93    int k;
94    for(k=1;  k <= nvarP; k++)
95    {
96      w[k]=deg(var(k));
97    }
98
99    BRlist[3] = list();
100    if( s==0 or (typeof(#[1]) != "string") )
101    {
102      if( w==1 )
103      {
104        BRlist[3][1]=list("dp",w);
105      }
106      else
107      {
108        BRlist[3][1]=list("wp",w);
109      }
110      BRlist[3][2]=list("C",intvec(0));
111      def Pfglm = ring(quotientList(BRlist));
112      setring Pfglm;
113    }
114  }
115  ideal i = fetch(P,i);
116
117  intvec opt = option(get);            //save options
118  option(redSB);
119  if (size(algorithm) > 0)
120  {
121    i = groebner(i,algorithm);
122  }
123  else
124  {
125    i = groebner(i);
126  }
127  option(set,opt);
128  setring P;
129  return (fglm(Pfglm,i));
130}
131example
132{ "EXAMPLE:"; echo = 2;
133   ring r  = 0,(x,y,z),lp;
134   ideal i = y3+x2,x2y+x2,x3-x2,z4-x2-y;
135   stdfglm(i);                   //uses fglm from "dp" (with groebner) to "lp"
136   stdfglm(i,"std");             //uses fglm from "dp" (with std) to "lp"
137
138   ring s  = (0,x),(y,z,u,v),lp;
139   minpoly = x2+1;
140   ideal i = u5-v4,zv-u2,zu3-v3,z2u-v2,z3-uv,yv-zu,yu-z2,yz-v,y2-u,u-xy2;
141   weight(i);
142   stdfglm(i,"(a(2,3,4,5),dp)"); //uses fglm from "(a(2,3,4,5),dp)" to "lp"
143}
144
145/////////////////////////////////////////////////////////////////////////////
146
147proc stdhilb(def i,list #)
148"SYNTAX: @code{stdhilb (} ideal_expression @code{)} @*
149         @code{stdhilb (} module_expression @code{)} @*
150         @code{stdhilb (} ideal_expression, intvec_expression @code{)}@*
151         @code{stdhilb (} module_expression, intvec_expression @code{)}@*
152         @code{stdhilb (} ideal_expression@code{,} list of string_expressions,
153               and intvec_expression @code{)} @*
154TYPE:    type of the first argument
155PURPOSE: Compute a Groebner basis of the ideal/module in the basering by
156         using the Hilbert driven Groebner basis algorithm.
157         If an argument of type string, stating @code{\"std\"} resp. @code{\"slimgb\"},
158         is given, the standard basis computation uses @code{std} or
159         @code{slimgb}, otherwise a heuristically chosen method (default)@*
160         If an optional second argument w of type intvec is given, w is used
161         as variable weights. If w is not given, it is computed as w[i] =
162         deg(var(i)). If the ideal is homogeneous w.r.t. w then the
163         Hilbert series is computed w.r.t. to these weights.
164THEORY:  If the ideal is not homogeneous compute first a Groebner basis
165         of the homogenization [w.r.t. the weights w] of the ideal/module,
166         then the Hilbert function and, finally, a Groebner basis in the
167         original ring by using the computed Hilbert function. If the given
168         w does not coincide with the variable weights of the basering, the
169         result may not be a groebner basis in the original ring.
170NOTE:    'Homogeneous' means weighted homogeneous with respect to the weights
171         w[i] of the variables var(i) of the basering. Parameters are not
172         converted to variables.
173SEE ALSO: stdfglm, std, slimgb, groebner
174KEYWORDS: Hilbert function
175EXAMPLE: example stdhilb;  shows an example"
176{
177
178//--------------------- save data from basering --------------------------
179  def P=basering;
180  int nr;
181  if (typeof(i)=="ideal") { nr=1;}
182  else                    { nr= nrows(i); }    //nr=1 if i is an ideal
183  ideal Qideal = ideal(P);      //defining the quotient ideal if P is a qring
184  int was_qring;                //remembers if basering was a qring
185  int is_homog =homog(i);       //check for homogeneity of i and Qideal
186  if (size(Qideal) > 0)
187  {
188     was_qring = 1;
189  }
190
191  // save ordering of basering P for later use
192  list ord_P = ringlist(P)[3];     //ordering of basering in ringlist
193  string ordstr_P = ordstr(P);     //ordering of basering as string
194  int nvarP = nvars(P);
195
196  //save options:
197  intvec gopt = option(get);
198  int p_opt;
199  string s_opt = option();
200  if (find(s_opt, "prot"))  { p_opt = 1; }
201
202//-------------------- check the given method and weights ---------------------
203//Note: stdhilb is used in elim where it is applied to an elimination ordering
204//a(1..1,0..0),wp(w). In such a ring deg(var(k)=0 for all vars corresponding to
205//0 in a(1..1,0..0), hence we cannot identify w via w[k] = deg(var(k));
206//Therefore hilbstd has the option to give ringweights.
207
208  int k;
209  string method;
210  for (k=1; k<=size(#); k++)
211  {
212    if (typeof(#[k]) == "intvec")
213    {
214        intvec w = #[k];         //given ringweights of basering P
215    }
216    if (typeof(#[k]) == "string")
217    {
218      method = method + "," + #[k];
219    }
220  }
221
222  if ( defined(w)!=voice )
223  {
224    intvec w;
225    for(k=nvarP;  k>=1; k--)
226    {
227       w[k] = deg(var(k));     //compute ring weights
228    }
229  }
230
231
232  if (npars(P) > 0)             //clear denominators of parameters
233  {
234    for( k=ncols(i); k>0; k-- )
235    {
236      i[k]=cleardenom(i[k]);
237    }
238  }
239
240//---------- exclude cases to which stdhilb should no be applied  ----------
241//Note that quotient ideal of qring must be homogeneous too
242
243  int neg=1-attrib (P,"global");
244
245  if( //find(ordstr_P,"s") ||// covered by neg
246     find(ordstr_P,"M") || neg )
247  {
248    // if( defined(hi) && is_homog )
249    // {
250    //  if (p_opt){"std with given Hilbert function in basering";}
251    //  return( std(i,hi,w) );
252    //### here we would need Hibert-Samuel function
253    // }
254
255    if (p_opt)
256    {"//-- stdhilb not implemented, we use std in ring:"; string(P);}
257    return( std(i) );
258  }
259
260//------------------------ change to hilbRing ----------------------------
261//The ground field of P and Philb coincide, Philb has an extra variable
262//@ or @(k). Philb is no qring and the predefined ideal/module Id(1) in
263//Philb is homogeneous (it is the homogenized i w.r.t. @ or @(k))
264//Parameters of P are not converted in Philb
265//Philb has only 1 block dp or wp(w)
266
267  list hiRi = hilbRing(i,w);
268  intvec W = hiRi[2];
269  def Philb = hiRi[1];
270  setring Philb;
271
272//-------- compute Hilbert series of homogenized ideal in Philb ---------
273//There are three cases
274
275  string algorithm;       //possibilities: std, slimgb, stdorslimgb
276  //define algorithm:
277  if( find(method,"std") && !find(method,"slimgb") )
278  {
279    algorithm = "std";
280  }
281  if( find(method,"slimgb") && !find(method,"std") )
282  {
283    algorithm = "slimgb";
284  }
285  if( find(method,"std") && find(method,"slimgb") ||
286    (!find(method,"std") && !find(method,"slimgb")) )
287  {
288    algorithm = "stdorslimgb";
289  }
290
291//### geaendert Dez08: es wird std(Id(1)) statt Id(1) aus Philb nach Phelp
292// weitergegeben fuer hilbertgetriebenen std
293
294  if (( algorithm=="std" || ( algorithm=="stdorslimgb" && char(P)>0 ) )
295  && (defined(hi)!=voice))
296  {
297    if (p_opt) {"compute hilbert series with std in ring " + string(Philb);
298                "weights used for hilbert series:",W;}
299    Id(1) = std(Id(1));
300    intvec hi = hilb( Id(1),1,W );
301  }
302  if (( algorithm=="slimgb" || ( algorithm=="stdorslimgb" && char(P)==0 ) )
303  && (defined(hi)!=voice))
304  {
305    if (p_opt) {"compute hilbert series with slimgb in ring " + string(Philb);
306                "weights used for hilbert series:",W;}
307    Id(1) = qslimgb(Id(1));
308    intvec hi = hilb( Id(1),1,W );
309  }
310
311  //-------------- we need another intermediate ring Phelp ----------------
312  //In Phelp we change only the ordering from Philb (otherwise it coincides
313  //with Philb). Phelp has in addition to P an extra homogenizing variable
314  //with name @ (resp. @(i) if @ and @(1), ..., @(i-1) are defined) with
315  //ordering an extra last block dp(1).
316  //Phelp has the same ordering as P on common variables. In Phelp
317  //a quotient ideal from P is added to the input
318
319  list BRlist = ringlist(Philb);
320  BRlist[3] = list();
321  int so = size(ord_P);
322  if( ord_P[so][1] =="c" || ord_P[so][1] =="C" )
323  {
324    list moduleord = ord_P[so];
325    so = so-1;
326  }
327  for (k=1; k<=so; k++)
328  {
329    BRlist[3][k] = ord_P[k];
330  }
331
332  BRlist[3][so+1] = list("dp",1);
333  w = w,1;
334
335  if( defined(moduleord)==voice )
336  {
337    BRlist[3][so+2] = moduleord;
338  }
339
340//--- change to extended ring Phelp and compute std with hilbert series ----
341  def Phelp = ring(quotientList(BRlist));
342  setring Phelp;
343  def i = imap(Philb, Id(1));
344  kill Philb;
345
346  // compute std with Hilbert series
347  option(redThrough);
348  if (w == 1)
349  {
350    if (p_opt){ "std with hilb in " + string(Phelp);}
351    i = std(i, hi);
352  }
353  else
354  {
355    if(p_opt){"std with weighted hilb in "+string(Phelp);}
356    i = std(i, hi, w);
357  }
358
359//-------------------- go back to original ring ---------------------------
360//The main computation is done. Do not forget to simplfy before maping.
361
362  // subst 1 for homogenizing var
363  if ( p_opt ) { "dehomogenization"; }
364  i = subst(i, var(nvars(basering)), 1);
365
366  if (p_opt) { "simplification"; }
367  i= simplify(i,34);
368
369  setring P;
370  if (p_opt) { "imap to ring "+string(P); }
371  i = imap(Phelp,i);
372  kill Phelp;
373  if( was_qring )
374  {
375    i = NF(i,std(0));
376  }
377  i = simplify(i,34);
378  // compute reduced SB
379  if (find(s_opt, "redSB") > 0)
380  {
381    if (p_opt) { "//interreduction"; }
382    i=interred(i);
383  }
384  attrib(i, "isSB", 1);
385  option(set,gopt);
386  return (i);
387}
388example
389{ "EXAMPLE:"; echo = 2;
390   ring  r = 0,(x,y,z),lp;
391   ideal i = y3+x2,x2y+x2z2,x3-z9,z4-y2-xz;
392   ideal j = stdhilb(i); j;
393
394   ring  r1 = 0,(x,y,z),wp(3,2,1);
395   ideal  i = y3+x2,x2y+x2z2,x3-z9,z4-y2-xz;  //ideal is homogeneous
396   ideal j = stdhilb(i,"std"); j;
397   //this is equivalent to:
398   intvec v = hilb(std(i),1);
399   ideal j1 = std(i,v,intvec(3,2,1)); j1;
400
401   size(NF(j,j1))+size(NF(j1,j));            //j and j1 define the same ideal
402}
403
404///////////////////////////////////////////////////////////////////////////////
405proc quotientList (list RL, list #)
406"SYNTAX: @code{quotientList (} list_expression @code{)} @*
407         @code{quotientList (} list_expression @code{,} string_expression@code{)}
408TYPE:    list
409PURPOSE: define a ringlist, say QL, of the first argument, say RL, which is
410         assumed to be the ringlist of a qring, but where the quotient ideal
411         RL[4] is not a standard basis with respect to the given monomial
412         order in RL[3]. Then QL will be obtained from RL just by replacing
413         RL[4] by a standard of it with respect to this order. RL itself
414         will be returnd if size(RL[4]) <= 1 (in which case it is known to be
415         a standard basis w.r.t. any ordering) or if a second argument
416         \"isSB\" of type string is given.
417NOTE:    the command ring(quotientList(RL)) defines a quotient ring correctly
418         and should be used instead of ring(RL) if the quotient ideal RL[4]
419         is not (or not known to be) a standard basis with respect to the
420         monomial ordering specified in RL[3].
421SEE ALSO: ringlist, ring
422EXAMPLE: example quotientList; shows an example"
423{
424  def P = basering;
425  if( size(#) > 0 )
426  {
427    if ( #[1] == "isSB")
428    {
429      return (RL);
430    }
431  }
432  ideal Qideal  = RL[4];  //##Achtung: falls basering Nullteiler hat, kann
433                           //die SB eines Elements mehrere Elemente enthalten
434  if( size(Qideal) <= 0)
435  {
436    return (RL);
437  }
438
439  RL[4] = ideal(0);
440  def Phelp = ring(RL);
441  setring Phelp;
442  ideal Qideal = groebner(fetch(P,Qideal));
443  setring P;
444  RL[4]=fetch(Phelp,Qideal);
445  return (RL);
446}
447example
448{ "EXAMPLE:"; echo = 2;
449   ring P = 0,(y,z,u,v),lp;
450   ideal i = y+u2+uv3, z+uv3;            //i is an lp-SB but not a dp_SB
451   qring Q = std(i);
452   list LQ = ringlist(Q);
453   LQ[3][1][1]="dp";
454   def Q1 = ring(quotientList(LQ));
455   setring Q1;
456   Q1;
457
458   setring Q;
459   ideal q1 = uv3+z, u2+y-z, yv3-zv3-zu; //q1 is a dp-standard basis
460   LQ[4] = q1;
461   def Q2 = ring(quotientList(LQ,"isSB"));
462   setring Q2;
463   Q2;
464}
465
466///////////////////////////////////////////////////////////////////////////////
467proc par2varRing (list #)
468"USAGE:   par2varRing([l]); l list of ideals/modules [default:l=empty list]
469RETURN:  list, say L, with L[1] a ring where the parameters of the
470         basering have been converted to an additional last block of
471         variables, all of weight 1, and ordering dp.
472         If a list l with l[i] an ideal/module is given, then
473         l[i] + minpoly*freemodule(nrows(l[i])) is mapped to an ideal/module
474         in L[1] with name Id(i).
475         If the basering has no parameters then L[1] is the basering.
476EXAMPLE: example par2varRing; shows an example"
477{
478  def P = basering;
479  int npar = npars(P);  //number of parameters
480  int s = size(#);
481  int ii;
482  if ( npar == 0)
483  {
484    dbprint(printlevel-voice+3,"// ** no parameters, ring was not changed");
485    for( ii = 1; ii <= s; ii++)
486    {
487      def Id(ii) = #[ii];
488      export (Id(ii));
489    }
490    return(list(P));
491  }
492
493  list rlist = ringlist(P);
494  list parlist = rlist[1];
495  rlist[1] = parlist[1];
496
497  string @Minpoly = string(minpoly);     //check for minpoly:
498  int sm = size(minpoly);
499  //now create new ring
500  for( ii = 1; ii <= s; ii++)
501  {
502    def Id(ii) = #[ii];
503  }
504  int nvar = size(rlist[2]);
505  int nblock = size(rlist[3]);
506  int k;
507  for (k=1; k<=npar; k++)
508  {
509    rlist[2][nvar+k] = parlist[2][k];   //change variable list
510  }
511
512  //converted parameters get one block dp. If module ordering was in front
513  //it stays in front, otherwise it will be moved to the end
514  intvec OW = 1:npar;
515  if( rlist[3][nblock][1] =="c" || rlist[3][nblock][1] =="C" )
516  {
517    rlist[3][nblock+1] = rlist[3][nblock];
518    rlist[3][nblock] = list("dp",OW);
519  }
520  else
521  {
522    rlist[3][nblock+1] = list("dp",OW);
523  }
524
525  def Ppar2var = ring(quotientList(rlist));
526  setring Ppar2var;
527  if ( sm == 0 )
528  {
529    for( ii = 1; ii <= s; ii++)
530    {
531      def Id(ii) = imap(P,Id(ii));
532      export (Id(ii));
533    }
534  }
535  else
536  {
537    if( find(option(),"prot") ){"//add minpoly to input";}
538    execute("poly Minpoly = " + @Minpoly + " ;");
539    for( ii = 1; ii <= s; ii++)
540    {
541      def Id(ii) = imap(P,Id(ii));
542      if (typeof(Id(ii))=="module")
543      {
544        Id(ii) = Id(ii),Minpoly*freemodule(nrows(Id(ii)));
545      }
546      else
547      {
548        Id(ii) = Id(ii),Minpoly;
549      }
550      export (Id(ii));
551    }
552  }
553  list Lpar2var = Ppar2var;
554  return(Lpar2var);
555}
556example
557{ "EXAMPLE:"; echo = 2;
558   ring R = (0,x),(y,z,u,v),lp;
559   minpoly = x2+1;
560   ideal i = x3,x2+y+z+u+v,xyzuv-1; i;
561   def P = par2varRing(i)[1]; P;
562   setring(P);
563   Id(1);
564
565   setring R;
566   module m = x3*[1,1,1], (xyzuv-1)*[1,0,1];
567   def Q = par2varRing(m)[1]; Q;
568   setring(Q);
569   print(Id(1));
570}
571
572//////////////////////////////////////////////////////////////////////////////
573proc hilbRing ( list # )
574"USAGE:   hilbRing([w,l]); w = intvec, l = list of ideals/modules
575RETURN:  list, say L: L[1] is a ring and L[2] an intvec
576         L[1] is a ring whith an extra homogenizing variable with name @,
577         resp. @(i) if @ and @(1), ..., @(i-1) are defined.
578         The monomial ordering of L[1] is consists of 1 block: dp if the
579         weights of the variables of the basering, say R, are all 1, resp.
580         wp(w,1) wehre w is either given or the intvec of weights of the
581         variables of R, i.e. w[k]=deg(var(k)).
582         If R is a quotient ring P/Q, then L[1] is not a quotient ring but
583         contains the ideal @Qidealhilb@, the homogenized ideal Q of P.
584         (Parameters of R are not touched).
585         If a list l is given with l[i] an ideal/module, then l[i] is mapped
586         to Id(i), the homogenized l[i]+Q*freemodule(nrows(l[i]) in L[1]
587         (Id(i) = l[i] if l[i] is already homogeneous).
588         L[2] is the intvec (w,1).
589PURPOSE: Prepare a ring for computing the (weighted) hilbert series of
590         an ideal/module with an easy monomial ordering.
591NOTE:    For this purpose we need w[k]=deg(var(k)). However, if the ordering
592         contains an extra weight vector a(v,0..0)) deg(var(k)) returns 0 for
593         k being an index which is 0 in a. Therefore we must compute w
594         beforehand and give it to hilbRing.
595EXAMPLE: example hilbRing; shows an example
596"
597{
598  def P = basering;
599  ideal Qideal = ideal(P);    //defining the quotient ideal if P is a qring
600  if( size(Qideal) != 0 )
601  {
602    int is_qring =1;
603  }
604  list BRlist = ringlist(P);
605  BRlist[4] = ideal(0);      //kill quotient ideal in BRlist
606
607  int nvarP = nvars(P);
608  int s = size(#);
609  int k;
610
611  for(k = 1; k <= s; k++)
612  {
613    if ( typeof(#[k]) == "intvec" )
614    {
615       intvec w = #[k];      //given weights for the variables
616       # = delete (#,k);
617    }
618  }
619
620  s = size(#);
621  for(k = 1; k <= s; k++)
622  {
623     def Id(k) = #[k];
624     int nr(k) = nrows(Id(k));
625  }
626
627  if ( defined(w)!=voice )
628  {
629    intvec w;                   //for ringweights of basering P
630    for(k=1;  k<=nvarP; k++)
631    {
632      w[k]=deg(var(k));        //degree of kth variable
633    }
634  }
635  //--------------------- a homogenizing variable is added ------------------
636  // call it @, resp. @(k) if @(1),...,@(k-1) are defined
637  string homvar;
638  if ( defined(@)==0 )
639  {
640    homvar = "@";
641  }
642  else
643  {
644    k=1;
645    while( defined(@(k)) != 0 )
646    {
647      k++;
648    }
649    homvar = "@("+string(k)+")";
650  }
651  BRlist[2][nvarP+1] = homvar;
652  w[nvarP +1]=1;
653
654  //ordering is set to (dp,C) if weights of all variables are 1
655  //resp. to (wp(w,1),C) where w are the ringweights of basering P
656  //homogenizing var gets weight 1:
657
658  BRlist[3] = list();
659  BRlist[3][2]=list("C",intvec(0));  //put module ordering always last
660  if(w==1)
661  {
662    BRlist[3][1]=list("dp",w);
663  }
664  else
665  {
666    BRlist[3][1]=list("wp",w);
667  }
668
669  //-------------- change ring and get ideal from previous ring ---------------
670  def Philb = ring(quotientList(BRlist));
671  kill BRlist;
672  setring Philb;
673  if( defined(is_qring)==voice )
674  {
675    ideal @Qidealhilb@ =  imap(P,Qideal);
676    if ( ! homog(@Qidealhilb@) )
677    {
678       @Qidealhilb@ =  homog( @Qidealhilb@, `homvar` );
679    }
680    export(@Qidealhilb@);
681
682    if( find(option(),"prot") ){"add quotient ideal to input";}
683
684    for(k = 1; k <= s; k++)
685    { //homogenize if necessary
686      def Id(k) =  imap(P,Id(k));
687      if ( ! homog(Id(k)) )
688      {
689         Id(k) =  homog( imap(P,Id(k)), `homvar` );
690      }
691      if (typeof(Id(k))=="module")
692      {
693        Id(k) =  Id(k),@Qidealhilb@*freemodule(nr(k)) ;
694      }
695      else
696      {
697        Id(k) =  Id(k),@Qidealhilb@ ;
698      }
699      export(Id(k));
700    }
701  }
702  else
703  {
704    for(k = 1; k <= s; k++)
705    { //homogenize if  necessary
706      def Id(k) =  imap(P,Id(k));
707      if ( ! homog(Id(k)) )
708      {
709         Id(k) =  homog( imap(P,Id(k)), `homvar` );
710      }
711      export(Id(k));
712    }
713  }
714  list Lhilb = Philb,w;
715  setring(P); return(Lhilb);
716}
717example
718{ "EXAMPLE:"; echo = 2;
719   ring R = 0,(x,y,z,u,v),lp;
720   ideal i = x+y2+z3,xy+xv+yz+zu+uv,xyzuv-1;
721   intvec w = 6,3,2,1,1;
722   hilbRing(i,w);
723   def P = hilbRing(w,i)[1];
724   setring P;
725   Id(1);
726   hilb(std(Id(1)),1);
727
728   ring S =  0,(x,y,z,u,v),lp;
729   qring T = std(x+y2+z3);
730   ideal i = xy+xv+yz+zu+uv,xyzuv-v5;
731   module m = i*[0,1,1] + (xyzuv-v5)*[1,1,0];
732   def Q = hilbRing(m)[1];  Q;
733   setring Q;
734   print(Id(1));
735}
736
737//////////////////////////////////////////////////////////////////////////////
738proc qslimgb (def i)
739"USAGE:   qslimgb(i); i ideal or module
740RETURN:  same type as input, a standard basis of i computed with slimgb
741NOTE:    Only as long as slimgb does not know qrings qslimgb should be used
742         in case the basering is (possibly) a quotient ring.
743         The quotient ideal is added to the input and slimgb is applied.
744EXAMPLE: example qslimgb; shows an example"
745{
746  def P = basering;
747  ideal Qideal = ideal(P);      //defining the quotient ideal if P is a qring
748  int p_opt;
749  if( find(option(),"prot") )
750  {
751    p_opt=1;
752  }
753  if (size(Qideal) == 0)
754  {
755    if (p_opt) { "slimgb in ring " + string(P); }
756    return(slimgb(i));
757  }
758
759  //case of a qring; since slimgb does not know qrings we
760  //delete the quotient ideal and add it to i
761
762  list BRlist = ringlist(P);
763  BRlist[4] = ideal(0);
764  def Phelp = ring(BRlist);
765  kill BRlist;
766  setring Phelp;
767  // module case:
768  def iq = imap(P,i);
769  iq = iq, imap(P,Qideal)*freemodule(nrows(iq));
770  if (p_opt)
771  {
772    "slimgb in ring " + string(Phelp);
773    "(with quotient ideal added to input)";
774  }
775  iq = slimgb(iq);
776
777  setring P;
778  if (p_opt) { "//imap to original ring"; }
779  i = imap(Phelp,iq);
780  kill Phelp;
781
782  if (find(option(),"redSB") > 0)
783  {
784    if (p_opt) { "//interreduction"; }
785    i=reduce(i,std(0));
786    i=interred(i);
787  }
788  attrib(i, "isSB", 1);
789  return (i);
790}
791example
792{ "EXAMPLE:"; echo = 2;
793   ring R  = (0,v),(x,y,z,u),dp;
794   qring Q = std(x2-y3);
795   ideal i = x+y2,xy+yz+zu+u*v,xyzu*v-1;
796   ideal j = qslimgb(i); j;
797
798   module m = [x+y2,1,0], [1,1,x2+y2+xyz];
799   print(qslimgb(m));
800}
801
802//////////////////////////////////////////////////////////////////////////////
803proc groebner(def i_par, list #)
804"SYNTAX: @code{groebner (} ideal_expression @code{)} @*
805         @code{groebner (} module_expression @code{)} @*
806         @code{groebner (} ideal_expression@code{,} int_expression @code{)} @*
807         @code{groebner (} module_expression@code{,} int_expression @code{)} @*
808         @code{groebner (} ideal_expression@code{,} list of string_expressions
809               @code{)} @*
810         @code{groebner (} ideal_expression@code{,} list of string_expressions
811               and int_expression @code{)} @*
812         @code{groebner (} ideal_expression@code{,} int_expression @code{)}
813TYPE:    type of the first argument
814PURPOSE: computes a standard basis of the first argument @code{I}
815         (ideal or module) by a heuristically chosen method (default)
816         or by a method specified by further arguments of type string.
817         Possible methods are:  @*
818         - the direct methods @code{\"std\"} or @code{\"slimgb\"} without
819           conversion, @*
820         - conversion methods @code{\"hilb\"} or @code{\"fglm\"} where
821           a Groebner basis is first computed with an \"easy\" ordering
822           and then converted to the ordering of the basering by the
823           Hilbert driven Groebner basis computation or by linear algebra.
824           The actual computation of the Groebner basis can be
825           specified by @code{\"std\"} or by @code{\"slimgb\"}
826           (not for all orderings implemented). @*
827         A further string @code{\"par2var\"} converts parameters to an extra
828         block of variables before a Groebner basis computation (and
829         afterwards back).
830         @code{option(prot)} informs about the chosen method.
831NOTE:    If an additional argument, say @code{wait}, of type int is given,
832         then the computation runs for at most @code{wait} seconds.
833         That is, if no result could be computed in @code{wait} seconds,
834         then the computation is interrupted, 0 is returned, a warning
835         message is displayed, and the global variable
836         @code{Standard::groebner_error} is defined.
837         This feature uses MP and hence it is available on UNIX platforms, only.
838HINT:    Since there exists no uniform best method for computing standard
839         bases, and since the difference in performance of a method on
840         different examples can be huge, it is recommended to test, for hard
841         examples, first various methods on a simplified example (e.g. use
842         characteristic 32003 instead of 0 or substitute a subset of
843         parameters/variables by integers, etc.). @*
844SEE ALSO: stdhilb, stdfglm, std, slimgb
845KEYWORDS: time limit on computations; MP, groebner basis computations
846EXAMPLE: example groebner;  shows an example"
847
848{
849//Vorgabe einer Teilmenge aus {hilb,fglm,par2var,std,slimgb}
850//V1: Erste Einstellungen (Jan 2007)
851//V2: Aktuelle Aenderungen (Juni 2008)
852//---------------------------------
853//0. Immer Aufruf von std unabhaengig von der Vorgabe:
854//   gemischte Ordnungen, extra Gewichtsvektor, Matrix Ordnungen
855//   ### Todo: extra Gewichtsvektor sollte nicht immer mit std wirken,
856//   sondern z.B. mit "hilb" arbeiten koennen
857//   ### Todo: es sollte ein Gewichtsvektor mitgegeben werden koennen (oder
858//   berechnet werden), z.B. groebner(I,"hilb",w) oder groebner(I,"withWeights")
859//   wie bei elim in elim.lib
860
861//1. Keine Vorgabe: es wirkt die aktuelle Heuristk:
862//   - Char = p: std
863//V1 - Char = 0: slimgb (im qring wird Quotientenideal zum Input addiert)
864//V2 - Char = 0: std
865//   - 1-Block-Ordnungen/non-commutative: direkt Aufruf von std oder slimgb
866//   - Komplizierte Ordnungen (lp oder > 1 Block): hilb
867//V1 - Parameter werden grundsaetzlich nicht in Variable umgewandelt
868//V2 - Mehr als ein Parmeter wird zu Variable konvertiert
869//   - fglm is keine Heuristik, da sonst vorher dim==0 peprueft werden muss
870
871//2. Vorgabe aus {std,slimgb}: es wird wo immer moeglich das Angegebene
872//   gewaehlt (da slimgb keine Hilbertfunktion kennt, wird std verwendet).
873//   Bei slimgb im qring, wird das Quotientenideal zum Ideal addiert.
874//   Bei Angabe von std zusammen mit slimgb (aequivalent zur Angabe von
875//   keinem von beidem) wirkt obige Heuristik.
876
877//3. Nichtleere Vorgabe aus {hilb,fglm,std,slimgb}:
878//   es wird nur das Angegebene und Moegliche sowie das Notwendige verwendet
879//   und bei Wahlmoeglickeit je nach Heuristik.
880//   Z.B. Vorgabe von {hilb} ist aequivalent zu {hilb,std,slimgb} und es wird
881//   hilb und nach Heuristik std oder slimgb verwendet,
882//   (V1: aber nicht par2var)
883//   bei Vorgabe von {hilb,slimgb} wird hilb und wo moeglich slimgb verwendet.
884
885//4. Bei Vorgabe von {par2var} wird par2var immer mit hilb und nach Heuristik
886//   std oder slimgb verwendet. Zu Variablen konvertierte Parameter haben
887//   extra letzten Block und Gewichte 1.
888
889  def P=basering;
890  if ((typeof(i_par)=="vector")||(typeof(i_par)=="module")||(typeof(i_par)=="matrix")) {module i=i_par;}
891  else {ideal i=i_par; } // int, poly, number, ideal
892  kill i_par;
893// check for integer etc coefficients
894  if (charstr(basering)[1]=="i") // either integer or integer,q
895  {
896    if (find(option(),"prot"))  { "calling std for ideals in ring with ring coefficients"; }
897    return (std(i));
898  }
899
900//----------------------- save the given method ---------------------------
901  string method;                //all given methods as a coma separated string
902  list Method;                  //all given methods as a list
903  int k;
904  for (k=1; k<=size(#); k++)
905  {
906     if (typeof(#[k]) == "int")
907     {
908       int wait = #[k];
909     }
910     if (typeof(#[k]) == "string")
911     {
912       method = method + "," + #[k];
913       Method = Method + list(#[k]);
914     }
915  }
916
917 //======= we have an argument of type int -- try to use MPfork links =======
918  if ( defined(wait) == voice )
919  {
920    if ( system("with", "MP") )
921    {
922        int j = 10;
923        string bs = nameof(basering);
924        link l_fork = "MPtcp:fork";
925        open(l_fork);
926        write(l_fork, quote(system("pid")));
927        int pid = read(l_fork);
928//        write(l_fork, quote(groebner(eval(i))));
929        write(l_fork, quote(groebner(eval(i),eval(Method))));
930//###Fehlermeldung:
931// ***dError: undef. ringorder used
932// occured at:
933
934        // sleep in small intervalls for appr. one second
935        if (wait > 0)
936        {
937          while(j < 1000000)
938          {
939            if (status(l_fork, "read", "ready", j)) {break;}
940            j = j + j;
941          }
942        }
943
944        // sleep in intervalls of one second from now on
945        j = 1;
946        while (j < wait)
947        {
948          if (status(l_fork, "read", "ready", 1000000)) {break;}
949          j = j + 1;
950        }
951
952        if (status(l_fork, "read", "ready"))
953        {
954          def result = read(l_fork);
955          if (bs != nameof(basering))
956          {
957            def PP = basering;
958            setring P;
959            def result = imap(PP, result);
960            kill PP;
961          }
962          if (defined(groebner_error)==1)
963          {
964            kill groebner_error;
965          }
966          kill l_fork;
967        }
968        else
969        {
970          ideal result;
971          if (! defined(groebner_error))
972          {
973            int groebner_error = 1;
974            export groebner_error;
975          }
976          "** groebner did not finish";
977          j = system("sh", "kill " + string(pid));
978        }
979        return (result);
980    }
981    else
982    {
983      "** groebner with a time limit on computation is not supported
984          in this configuration";
985    }
986  }
987
988 //=========== we are still here -- do the actual computation =============
989
990//--------------------- save data from basering ---------------------------
991  string @Minpoly = string(minpoly);      //minimal polynomial
992  int was_minpoly;             //remembers if there was a minpoly in P
993  if (size(minpoly) > 0)
994  {
995     was_minpoly = 1;
996  }
997
998  ideal Qideal = ideal(P);      //defining the quotient ideal if P is a qring
999  int was_qring;                //remembers if basering was a qring
1000  //int is_homog = 1;
1001  if (size(Qideal) > 0)
1002  {
1003     was_qring = 1;
1004     //is_homog = homog(Qideal); //remembers if Qideal was homog (homog(0)=1)
1005  }
1006  list BRlist = ringlist(P);     //ringlist of basering
1007
1008  // save ordering of basering P for later use
1009  list ord_P = BRlist[3];       //should be available in all rings
1010  string ordstr_P = ordstr(P);
1011  int nvars_P = nvars(P);
1012  int npars_P = npars(P);
1013  intvec w;                     //for ringweights of basering P
1014  for(k=1;  k<=nvars_P; k++)
1015  {
1016     w[k]=deg(var(k));
1017  }
1018  int neg=1-attrib (P,"global");
1019
1020  //save options:
1021  intvec opt=option(get);
1022  string s_opt = option();
1023  int p_opt;
1024  if (find(s_opt, "prot"))  { p_opt = 1; }
1025
1026//------------------ cases where std is always used ------------------------
1027//If other methods are not implemented or do not make sense, i.e. for
1028//local or mixed orderings, matrix orderings, extra weight vector
1029//### Todo: extra weight vector should be allowed for e.g. with "hilb"
1030
1031  if(  //( find(ordstr_P,"s") > 0 ) || // covered by neg
1032       ( find(ordstr_P,"M") > 0 )  || ( find(ordstr_P,"a") > 0 )  || neg )
1033  {
1034    if (p_opt) { "std in basering"; }
1035    return(std(i));
1036  }
1037
1038//now we have:
1039//ideal or module, global ordering, no matrix ordering, no extra weight vector
1040//The interesting cases start now.
1041
1042 //------------------ classify the possible settings ---------------------
1043  string algorithm;       //possibilities: std, slimgb, stdorslimgb
1044  string conversion;      //possibilities: hilb, fglm, hilborfglm, no
1045  string partovar;        //possibilities: yes, no
1046  string order;           //possibilities: simple, !simple
1047  string direct;          //possibilities: yes, no
1048
1049  //define algorithm:
1050  if( find(method,"std") && !find(method,"slimgb") )
1051  {
1052    algorithm = "std";
1053  }
1054  if( find(method,"slimgb") && !find(method,"std") )
1055  {
1056    algorithm = "slimgb";
1057  }
1058  if( find(method,"std") && find(method,"slimgb") ||
1059      (!find(method,"std") && !find(method,"slimgb")) )
1060  {
1061    algorithm = "stdorslimgb";
1062  }
1063
1064  //define conversion:
1065  if( find(method,"hilb") && !find(method,"fglm") )
1066  {
1067     conversion = "hilb"; // $Id$
1068  }
1069  if( find(method,"fglm") && !find(method,"hilb") )
1070  {
1071    conversion = "fglm"; // $Id$
1072  }
1073  if( find(method,"fglm") && find(method,"hilb") )
1074  {
1075    conversion = "hilborfglm"; // $Id$
1076  }
1077  if( !find(method,"fglm") && !find(method,"hilb") )
1078  {
1079    conversion = "no"; // $Id$
1080  }
1081
1082  //define partovar:
1083  //if( find(method,"par2var") && npars_P > 0 )   //V1
1084  if( find(method,"par2var") || npars_P > 1 )     //V2
1085  {
1086     partovar = "yes";
1087  }
1088  else
1089  {
1090     partovar = "no";
1091  }
1092
1093  //define order:
1094  if (system("nblocks") <= 2)
1095  {
1096    if ( find(ordstr_P,"M")+find(ordstr_P,"lp")+find(ordstr_P,"rp") <= 0 )
1097    {
1098      order = "simple";
1099    }
1100  }
1101
1102  //define direct:
1103  if ( (order=="simple" && (size(method)==0)) ||
1104       (size(BRlist)>4) ||
1105        (order=="simple" && (method==",par2var" && npars_P==0 )) ||
1106         (conversion=="no" && partovar=="no" &&
1107           (algorithm=="std" || algorithm=="slimgb" ||
1108            (find(method,"std") && find(method,"slimgb")) ) ) )
1109  {
1110    direct = "yes";
1111  }
1112  else
1113  {
1114    direct = "no";
1115  }
1116
1117  //order=="simple" means that the ordering of the variables consists of one
1118  //block which is not a matrix ordering and not a lexicographical ordering.
1119  //(Note:Singular counts always least 2 blocks, one is for module component):
1120  //Call a method "direct" if conversion=="no" && partovar="no" which means
1121  //that we apply std or slimgb dircet in the basering (exception
1122  //as long as slimgb does not know qrings: in a qring of a ring P
1123  //the ideal Qideal is added to the ideal and slimgb is applied in P).
1124  //We apply a direct method if we have a simple monomial ordering, if no
1125  //conversion (fglm or hilb) is specified and if the parameters shall
1126  //not be made to variables
1127  //BRlist (=ringlist of basering) > 4 if the basering is non-commutative
1128//---------------------------- direct methods -----------------------------
1129  if ( direct == "yes" )
1130  {
1131  //if ( algorithm=="std" || (algorithm=="stdorslimgb" && char(P)>0) )   //V1
1132    if ( algorithm=="std" || (algorithm=="stdorslimgb") )                //V2
1133    {
1134      if (p_opt) { "std in " + string(P); }
1135      return(std(i));
1136    }
1137  //if( algorithm=="slimgb" || (algorithm=="stdorslimgb" && char(P)==0)) //V1
1138    if ( algorithm=="slimgb" )                                           //V2
1139    {
1140      return(qslimgb(i));
1141    }
1142  }
1143
1144//--------------------------- indirect methods -----------------------------
1145//indirect methods are methods where a conversion is used with a ring change
1146//We are in the following situation:
1147//direct=="no" (i.e. "hilb" or "fglm" or "par2var" is given)
1148//or no method is given and we have a complicated monomial ordering
1149//V1: "par2var" is not a default strategy, it must be explicitely
1150//given in order to be performed.
1151//V2: "par2var" is a default strategy if there are more than 1 parameters
1152
1153//------------ case where no parameters are made to variables  -------------
1154  if (  partovar == "no" && conversion == "hilb"
1155    || (partovar == "no" && conversion == "fglm" )
1156    || (partovar == "no" && conversion == "hilborfglm" )
1157    || (partovar == "no" && conversion == "no" && direct == "no") )
1158  //last case: heuristic
1159  {
1160    if ( conversion=="fglm" )
1161    {
1162    //if ( algorithm=="std" || (algorithm=="stdorslimgb" && char(P)>0) ) //V1
1163      if ( algorithm=="std" || (algorithm=="stdorslimgb") )              //V2
1164      {
1165        return (stdfglm(i,"std"));
1166      }
1167    //if(algorithm=="slimgb" || (algorithm=="stdorslimgb" && char(P)==0))//V1
1168      if( algorithm=="slimgb" )                                          //V2
1169      {
1170        return (stdfglm(i,"slimgb"));
1171      }
1172    }
1173    else
1174    {
1175    //if ( algorithm=="std" || (algorithm=="stdorslimgb" && char(P)>0) )//V1
1176      if ( algorithm=="std" || (algorithm=="stdorslimgb" ) )            //V2
1177      {
1178        return (stdhilb(i,"std"));
1179      }
1180    //if(algorithm=="slimgb" || (algorithm=="stdorslimgb" && char(P)==0))//V1
1181      if ( algorithm=="slimgb" )                                         //V2
1182      {
1183        return (stdhilb(i,"slimgb"));
1184      }
1185    }
1186  }
1187
1188//------------ case where parameters are made to variables  ----------------
1189//define a ring Phelp via par2varRing in which the parameters are variables
1190
1191  else
1192  {
1193    // reset options
1194    option(none);
1195    // turn on options prot, mem, redSB, intStrategy if previously set
1196    if ( find(s_opt, "prot") )
1197      { option(prot); }
1198    if ( find(s_opt, "mem") )
1199      { option(mem); }
1200    if ( find(s_opt, "redSB") )
1201      { option(redSB); }
1202    if ( find(s_opt, "intStrategy") )
1203      { option(intStrategy); }
1204
1205    //first clear denominators of parameters
1206    if (npars_P > 0)
1207    {
1208      for( k=ncols(i); k>0; k-- )
1209      { i[k]=cleardenom(i[k]); }
1210    }
1211
1212    def Phelp = par2varRing(i)[1];   //minpoly is mapped with i
1213    setring Phelp;
1214    def i = Id(1);
1215    //is_homog = homog(i);
1216
1217    //If parameters are converted to ring variables, they appear in an extra
1218    //block. Therefore we use always hilb for this block ordering:
1219    if ( conversion=="fglm" )
1220    {
1221      i = (stdfglm(i));       //only uesful for 1 parameter with minpoly
1222    }
1223    else
1224    {
1225    //if ( algorithm=="std" || (algorithm=="stdorslimgb" && char(P)>0) )//V1
1226      if ( algorithm=="std" || (algorithm=="stdorslimgb" ))             //V2
1227      {
1228        i = stdhilb(i,"std");
1229      }
1230    //if(algorithm=="slimgb" || (algorithm=="stdorslimgb" && char(P)==0))//V1
1231      if ( algorithm=="slimgb" )                                         //V2
1232      {
1233        i = stdhilb(i,"slimgb");
1234      }
1235    }
1236  }
1237
1238//-------------------- go back to original ring ---------------------------
1239//The main computation is done. However, the SB coming from a ring with
1240//extra variables is in general too big. We simplify it before mapping it
1241//to the basering.
1242
1243  if (p_opt) { "//simplification"; }
1244
1245  if (was_minpoly)
1246  {
1247    execute("ideal Minpoly = " + @Minpoly + ";");
1248    attrib(Minpoly,"isSB",1);
1249    i = simplify(NF(i,Minpoly),2);
1250  }
1251
1252  def Li = lead(i);
1253  setring P;
1254  def Li = imap(Phelp,Li);
1255  Li = simplify(Li,32);
1256  intvec vi;
1257  for (k=1; k<=ncols(Li); k++)
1258  {
1259    vi[k] = Li[k]==0;
1260  }
1261
1262  setring Phelp;
1263  for (k=1;  k<=size(i) ;k++)
1264  {
1265    if(vi[k]==1)
1266    {
1267      i[k]=0;
1268    }
1269  }
1270  i = simplify(i,2);
1271
1272  setring P;
1273  if (p_opt) { "//imap to original ring"; }
1274  i = imap(Phelp,i);
1275  kill Phelp;
1276  i = simplify(i,34);
1277
1278  // clean-up time
1279  option(set, opt);
1280  if (find(s_opt, "redSB") > 0)
1281  {
1282    if (p_opt) { "//interreduction"; }
1283    i=interred(i);
1284  }
1285  attrib(i, "isSB", 1);
1286  return (i);
1287}
1288example
1289{ "EXAMPLE: "; echo=2;
1290  intvec opt = option(get);
1291  option(prot);
1292  ring r  = 0,(a,b,c,d),dp;
1293  ideal i = a+b+c+d,ab+ad+bc+cd,abc+abd+acd+bcd,abcd-1;
1294  groebner(i);
1295
1296  ring s  = 0,(a,b,c,d),lp;
1297  ideal i = imap(r,i);
1298  groebner(i,"hilb");
1299
1300  ring R  = (0,a),(b,c,d),lp;
1301  minpoly = a2+1;
1302  ideal i = a+b+c+d,ab+ad+bc+cd,abc+abd+acd+bcd,d2-c2b2;
1303  groebner(i,"par2var","slimgb");
1304
1305  groebner(i,"fglm");          //computes a reduced standard basis
1306
1307  option(set,opt);
1308}
1309
1310//////////////////////////////////////////////////////////////////////////
1311
1312proc res(list #)
1313"@c we do texinfo here:
1314@cindex resolution, computation of
1315@table @code
1316@item @strong{Syntax:}
1317@code{res (} ideal_expression@code{,} int_expression @code{[,} any_expression @code{])}
1318@*@code{res (} module_expression@code{,} int_expression @code{[,} any_expression @code{])}
1319@item @strong{Type:}
1320resolution
1321@item @strong{Purpose:}
1322computes a (possibly minimal) free resolution of an ideal or module using
1323a heuristically chosen method.
1324@* The second (int) argument (say @code{k}) specifies the length of
1325the resolution. If it is not positive then @code{k} is assumed to be the
1326number of variables of the basering.
1327@* If a third argument is given, the returned resolution is minimized.
1328
1329Depending on the input, the returned resolution is computed using the
1330following methods:
1331@table @asis
1332@item @strong{quotient rings:}
1333@code{nres} (classical method using syzygies) , see @ref{nres}.
1334
1335@item @strong{homogeneous ideals and k=0:}
1336@code{lres} (La'Scala's method), see @ref{lres}.
1337
1338@item @strong{not minimized resolution and (homogeneous input with k not 0, or local rings):}
1339@code{sres} (Schreyer's method), see @ref{sres}.
1340
1341@item @strong{all other inputs:}
1342@code{mres} (classical method), see @ref{mres}.
1343@end table
1344@item @strong{Note:}
1345Accessing single elements of a resolution may require some partial
1346computations to be finished and may therefore take some time.
1347@end table
1348@c ref
1349See also
1350@ref{betti};
1351@ref{ideal};
1352@ref{minres};
1353@ref{module};
1354@ref{mres};
1355@ref{nres};
1356@ref{lres};
1357@ref{hres};
1358@ref{sres};
1359@ref{resolution}.
1360@c ref
1361"
1362{
1363   def P=basering;
1364   if (size(#) < 2)
1365   {
1366     ERROR("res: need at least two arguments: ideal/module, int");
1367   }
1368
1369   def m=#[1]; //the ideal or module
1370   int i=#[2]; //the length of the resolution
1371   if (i< 0) { i=0;}
1372
1373   string varstr_P = varstr(P);
1374
1375   int p_opt;
1376   string s_opt = option();
1377   // set p_opt, if option(prot) is set
1378   if (find(s_opt, "prot"))
1379   {
1380     p_opt = 1;
1381   }
1382
1383   if( (size(ideal(basering)) > 0) || (size(ringlist(P)) > 4) )
1384   {
1385     // the quick hack for qrings - seems to fit most needs
1386     // (lres is not implemented for qrings, sres is not so efficient)
1387     // || non-commutative, since only n/m-res are implemented for NC rings
1388     if (p_opt) { "using nres";}
1389     return(nres(m,i));
1390   }
1391
1392/* if( attrib(basering, "global") == 1 ) // preparations for s_res usage. in testing!
1393   {
1394     if (p_opt) { "using s_res";}
1395     if( !defined(s_res) )
1396     {
1397       def @@v=option(get); option(noloadLib); option(noloadProc); LIB( "schreyer.lib" ); // for s_res
1398       option(set, @@v); kill @@v;
1399     }
1400     resolution re = s_res(m,i);
1401     if(size(#)>2)
1402     {
1403       re=minres(re);
1404     }
1405     return(re);
1406   }*/
1407
1408   if(homog(m)==1)
1409   {
1410      resolution re;
1411      if (((i==0) or (i>=nvars(basering))) && (typeof(m) != "module") && (nvars(basering)>1))
1412      {
1413        //LaScala for the homogeneous case and i == 0
1414        if (p_opt) { "using lres";}
1415        re=lres(m,i);
1416        if(size(#)>2)
1417        {
1418           re=minres(re);
1419        }
1420      }
1421      else
1422      {
1423        if(size(#)>2)
1424        {
1425          if (p_opt) { "using mres";}
1426          re=mres(m,i);
1427        }
1428        else
1429        {
1430          if (p_opt) { "using sres";}
1431          re=sres(std(m),i);
1432        }
1433      }
1434      return(re);
1435   }
1436
1437   //mres for the global non homogeneous case
1438   if(find(ordstr(P),"s")==0)
1439   {
1440      string ri= "ring Phelp ="
1441                  +string(char(P))+",("+varstr_P+"),(dp,C);";
1442      ri = ri + "minpoly = "+string(minpoly) + ";";
1443      execute(ri);
1444      def m=imap(P,m);
1445      if (p_opt) { "using mres in another ring";}
1446      list re=mres(m,i);
1447      setring P;
1448      resolution result=imap(Phelp,re);
1449      if (size(#) > 2) {result = minres(result);}
1450      return(result);
1451   }
1452
1453   //sres for the local case and not minimal resolution
1454   if(size(#)<=2)
1455   {
1456      string ri= "ring Phelp ="
1457                  +string(char(P))+",("+varstr_P+"),(ls,c);";
1458      ri = ri + "minpoly = "+string(minpoly) + ";";
1459      execute(ri);
1460      def m=imap(P,m);
1461      m=std(m);
1462      if (p_opt) { "using sres in another ring";}
1463      list re=sres(m,i);
1464      setring P;
1465      resolution result=imap(Phelp,re);
1466      return(result);
1467   }
1468
1469   //mres for the local case and minimal resolution
1470   string ri= "ring Phelp ="
1471                  +string(char(P))+",("+varstr_P+"),(ls,C);";
1472   ri = ri + "minpoly = "+string(minpoly) + ";";
1473   execute(ri);
1474   def m=imap(P,m);
1475    if (p_opt) { "using mres in another ring";}
1476   list re=mres(m,i);
1477   setring P;
1478   resolution result=imap(Phelp,re);
1479   result = minres(result);
1480   return(result);
1481}
1482example
1483{"EXAMPLE:"; echo = 2;
1484  ring r=0,(x,y,z),dp;
1485  ideal i=xz,yz,x3-y3;
1486  def l=res(i,0); // homogeneous ideal: uses lres
1487  l;
1488  print(betti(l), "betti"); // input to betti may be of type resolution
1489  l[2];         // element access may take some time
1490  i=i,x+1;
1491  l=res(i,0);   // inhomogeneous ideal: uses mres
1492  l;
1493  ring rs=0,(x,y,z),ds;
1494  ideal i=imap(r,i);
1495  def l=res(i,0); // local ring not minimized: uses sres
1496  l;
1497  res(i,0,0);     // local ring and minimized: uses mres
1498}
1499/////////////////////////////////////////////////////////////////////////
1500
1501proc quot (def m1,def m2,list #)
1502"SYNTAX: @code{quot (} module_expression@code{,} module_expression @code{)}
1503         @*@code{quot (} module_expression@code{,} module_expression@code{,}
1504            int_expression @code{)}
1505         @*@code{quot (} ideal_expression@code{,} ideal_expression @code{)}
1506         @*@code{quot (} ideal_expression@code{,} ideal_expression@code{,}
1507            int_expression @code{)}
1508TYPE:    ideal
1509SYNTAX:  @code{quot (} module_expression@code{,} ideal_expression @code{)}
1510TYPE:    module
1511PURPOSE: computes the quotient of the 1st and the 2nd argument.
1512         If a 3rd argument @code{n} is given the @code{n}-th method is used
1513         (@code{n}=1...5).
1514SEE ALSO: quotient
1515EXAMPLE: example quot; shows an example"
1516{
1517  if (((typeof(m1)!="ideal") and (typeof(m1)!="module"))
1518     or ((typeof(m2)!="ideal") and (typeof(m2)!="module")))
1519  {
1520    "USAGE:   quot(m1, m2[, n]); m1, m2 two submodules of k^s,";
1521    "         n (optional) integer (1<= n <=5)";
1522    "RETURN:  the quotient of m1 and m2";
1523    "EXAMPLE: example quot; shows an example";
1524    return();
1525  }
1526  if (typeof(m1)!=typeof(m2))
1527  {
1528    return(quotient(m1,m2));
1529  }
1530  if (size(#)>0)
1531  {
1532    if (typeof(#[1])=="int" )
1533    {
1534      return(quot1(m1,m2,#[1]));
1535    }
1536  }
1537  else
1538  {
1539    return(quot1(m1,m2,2));
1540  }
1541}
1542example
1543{ "EXAMPLE:"; echo = 2;
1544  ring r=181,(x,y,z),(c,ls);
1545  ideal id1=maxideal(4);
1546  ideal id2=x2+xyz,y2-z3y,z3+y5xz;
1547  option(prot);
1548  ideal id3=quotient(id1,id2);
1549  id3;
1550  ideal id4=quot(id1,id2,1);
1551  id4;
1552  ideal id5=quot(id1,id2,2);
1553  id5;
1554}
1555
1556static proc quot1 (module m1, module m2,int n)
1557"USAGE:   quot1(m1, m2, n); m1, m2 two submodules of k^s,
1558         n integer (1<= n <=5)
1559RETURN:  the quotient of m1 and m2
1560EXAMPLE: example quot1; shows an example"
1561{
1562  if (n==1)
1563  {
1564    return(quotient1(m1,m2));
1565  }
1566  else
1567  {
1568    if (n==2)
1569    {
1570      return(quotient2(m1,m2));
1571    }
1572    else
1573    {
1574      if (n==3)
1575      {
1576        return(quotient3(m1,m2));
1577      }
1578      else
1579      {
1580        if (n==4)
1581        {
1582          return(quotient4(m1,m2));
1583        }
1584        else
1585        {
1586          if (n==5)
1587          {
1588            return(quotient5(m1,m2));
1589          }
1590          else
1591          {
1592            return(quotient(m1,m2));
1593          }
1594        }
1595      }
1596    }
1597  }
1598}
1599example
1600{ "EXAMPLE:"; echo = 2;
1601  ring r=181,(x,y,z),(c,ls);
1602  ideal id1=maxideal(4);
1603  ideal id2=x2+xyz,y2-z3y,z3+y5xz;
1604  option(prot);
1605  ideal id6=quotient(id1,id2);
1606  id6;
1607  ideal id7=quot1(id1,id2,1);
1608  id7;
1609  ideal id8=quot1(id1,id2,2);
1610  id8;
1611}
1612
1613static proc quotient0(module a,module b)
1614{
1615  module mm=b+a;
1616  resolution rs=lres(mm,0);
1617  list I=list(rs);
1618  matrix M=I[2];
1619  matrix A[1][nrows(M)]=M[1..nrows(M),1];
1620  ideal i=A;
1621  return (i);
1622}
1623proc quotient1(module a,module b)  //17sec
1624"USAGE:   quotient1(m1, m2); m1, m2 two submodules of k^s,
1625RETURN:  the quotient of m1 and m2"
1626{
1627  int i;
1628  a=std(a);
1629  module dummy;
1630  module B=NF(b,a)+dummy;
1631  ideal re=quotient(a,module(B[1]));
1632  for(i=2;i<=ncols(B);i++)
1633  {
1634     re=intersect1(re,quotient(a,module(B[i])));
1635  }
1636  return(re);
1637}
1638proc quotient2(module a,module b)    //13sec
1639"USAGE:   quotient2(m1, m2); m1, m2 two submodules of k^s,
1640RETURN:  the quotient of m1 and m2"
1641{
1642  a=std(a);
1643  module dummy;
1644  module bb=NF(b,a)+dummy;
1645  int i=ncols(bb);
1646  ideal re=quotient(a,module(bb[i]));
1647  bb[i]=0;
1648  module temp;
1649  module temp1;
1650  module bbb;
1651  int mx;
1652  i=i-1;
1653  while (1)
1654  {
1655    if (i==0) break;
1656    temp = a+bb*re;
1657    temp1 = lead(interred(temp));
1658    mx=ncols(a);
1659    if (ncols(temp1)>ncols(a))
1660    {
1661      mx=ncols(temp1);
1662    }
1663    temp1 = matrix(temp1,1,mx)-matrix(lead(a),1,mx);
1664    temp1 = dummy+temp1;
1665    if (deg(temp1[1])<0) break;
1666    re=intersect1(re,quotient(a,module(bb[i])));
1667    bb[i]=0;
1668    i = i-1;
1669  }
1670  return(re);
1671}
1672proc quotient3(module a,module b)   //89sec
1673"USAGE:   quotient3(m1, m2); m1, m2 two submodules of k^s,
1674         only for global rings
1675RETURN:  the quotient of m1 and m2"
1676{
1677  string s="ring @newr=("+charstr(basering)+
1678           "),("+varstr(basering)+",@t,@w),dp;";
1679  def @newP=basering;
1680  execute(s);
1681  module b=imap(@newP,b);
1682  module a=imap(@newP,a);
1683  int i;
1684  int j=ncols(b);
1685  vector @b;
1686  for(i=1;i<=j;i++)
1687  {
1688     @b=@b+@t^(i-1)*@w^(j-i+1)*b[i];
1689  }
1690  ideal re=quotient(a,module(@b));
1691  setring @newP;
1692  ideal re=imap(@newr,re);
1693  return(re);
1694}
1695proc quotient5(module a,module b)   //89sec
1696"USAGE:   quotient5(m1, m2); m1, m2 two submodules of k^s,
1697         only for global rings
1698RETURN:  the quotient of m1 and m2"
1699{
1700  string s="ring @newr=("+charstr(basering)+
1701           "),("+varstr(basering)+",@t),dp;";
1702  def @newP=basering;
1703  execute(s);
1704  module b=imap(@newP,b);
1705  module a=imap(@newP,a);
1706  int i;
1707  int j=ncols(b);
1708  vector @b;
1709  for(i=1;i<=j;i++)
1710  {
1711     @b=@b+@t^(i-1)*b[i];
1712  }
1713  @b=homog(@b,@w);
1714  ideal re=quotient(a,module(@b));
1715  setring @newP;
1716  ideal re=imap(@newr,re);
1717  return(re);
1718}
1719proc quotient4(module a,module b)   //95sec
1720"USAGE:   quotient4(m1, m2); m1, m2 two submodules of k^s,
1721         only for global rings
1722RETURN:  the quotient of m1 and m2"
1723{
1724  string s="ring @newr=("+charstr(basering)+
1725           "),("+varstr(basering)+",@t),dp;";
1726  def @newP=basering;
1727  execute(s);
1728  module b=imap(@newP,b);
1729  module a=imap(@newP,a);
1730  int i;
1731  vector @b=b[1];
1732  for(i=2;i<=ncols(b);i++)
1733  {
1734     @b=@b+@t^(i-1)*b[i];
1735  }
1736  matrix sy=modulo(@b,a);
1737  ideal re=sy;
1738  setring @newP;
1739  ideal re=imap(@newr,re);
1740  return(re);
1741}
1742static proc intersect1(ideal i,ideal j)
1743{
1744  def R=basering;
1745  execute("ring gnir = ("+charstr(basering)+"),
1746                       ("+varstr(basering)+",@t),(C,dp);");
1747  ideal i=var(nvars(basering))*imap(R,i)+(var(nvars(basering))-1)*imap(R,j);
1748  ideal j=eliminate(i,var(nvars(basering)));
1749  setring R;
1750  map phi=gnir,maxideal(1);
1751  return(phi(j));
1752}
1753
1754//////////////////////////////////////////////////////////////////
1755///
1756/// sprintf, fprintf printf
1757///
1758proc sprintf(string fmt, list #)
1759"SYNTAX:  @code{sprintf (} string_expression @code{[,} any_expressions
1760               @code{] )}
1761RETURN:   string
1762PURPOSE:  @code{sprintf(fmt,...);} performs output formatting. The first
1763          argument is a format control string. Additional arguments may be
1764          required, depending on the content of the control string. A series
1765          of output characters is generated as directed by the control string;
1766          these characters are returned as a string. @*
1767          The control string @code{fmt} is simply text to be copied,
1768          except that the string may contain conversion specifications.@*
1769          Type @code{help print;} for a listing of valid conversion
1770          specifications. As an addition to the conversions of @code{print},
1771          the @code{%n} and @code{%2} conversion specification does not
1772          consume an additional argument, but simply generates a newline
1773          character.
1774NOTE:     If one of the additional arguments is a list, then it should be
1775          wrapped in an additional @code{list()} command, since passing a list
1776          as an argument flattens the list by one level.
1777SEE ALSO: fprintf, printf, print, string
1778EXAMPLE : example sprintf; shows an example
1779"
1780{
1781  int sfmt = size(fmt);
1782  if (sfmt  <= 1)
1783  {
1784    return (fmt);
1785  }
1786  int next, l, nnext;
1787  string ret;
1788  list formats = "%l", "%s", "%2l", "%2s", "%t", "%;", "%p", "%b", "%n", "%2";
1789  while (1)
1790  {
1791    if (size(#) <= 0)
1792    {
1793      return (ret + fmt);
1794    }
1795    nnext = 0;
1796    while (nnext < sfmt)
1797    {
1798      nnext = find(fmt, "%", nnext + 1);
1799      if (nnext == 0)
1800      {
1801        next = 0;
1802        break;
1803      }
1804      l = 1;
1805      while (l <= size(formats))
1806      {
1807        next = find(fmt, formats[l], nnext);
1808        if (next == nnext) break;
1809        l++;
1810      }
1811      if (next == nnext) break;
1812    }
1813    if (next == 0)
1814    {
1815      return (ret + fmt);
1816    }
1817    if (formats[l] != "%2" && formats[l] != "%n")
1818    {
1819      ret = ret + fmt[1, next - 1] + print(#[1], formats[l]);
1820      # = delete(#, 1);
1821    }
1822    else
1823    {
1824      ret = ret + fmt[1, next - 1] + print("", "%2s");
1825    }
1826    if (size(fmt) <= (next + size(formats[l]) - 1))
1827    {
1828      return (ret);
1829    }
1830    fmt = fmt[next + size(formats[l]), size(fmt)-next-size(formats[l]) + 1];
1831  }
1832}
1833example
1834{ "EXAMPLE:"; echo=2;
1835  ring r=0,(x,y,z),dp;
1836  module m=[1,y],[0,x+z];
1837  intmat M=betti(mres(m,0));
1838  list l = r, m, M;
1839  string s = sprintf("s:%s,%n l:%l", 1, 2); s;
1840  s = sprintf("s:%n%s", l); s;
1841  s = sprintf("s:%2%s", list(l)); s;
1842  s = sprintf("2l:%n%2l", list(l)); s;
1843  s = sprintf("%p", list(l)); s;
1844  s = sprintf("%;", list(l)); s;
1845  s = sprintf("%b", M); s;
1846}
1847
1848proc printf(string fmt, list #)
1849"SYNTAX:  @code{printf (} string_expression @code{[,} any_expressions@code{] )}
1850RETURN:   none
1851PURPOSE:  @code{printf(fmt,...);} performs output formatting. The first
1852          argument is a format control string. Additional arguments may be
1853          required, depending on the content of the control string. A series
1854          of output characters is generated as directed by the control string;
1855          these characters are displayed (i.e., printed to standard out). @*
1856          The control string @code{fmt} is simply text to be copied, except
1857          that the string may contain conversion specifications. @*
1858          Type @code{help print;} for a listing of valid conversion
1859          specifications. As an addition to the conversions of @code{print},
1860          the @code{%n} and @code{%2} conversion specification does not
1861          consume an additional argument, but simply generates a newline
1862          character.
1863NOTE:     If one of the additional arguments is a list, then it should be
1864          enclosed once more into a @code{list()} command, since passing a
1865          list as an argument flattens the list by one level.
1866SEE ALSO: sprintf, fprintf, print, string
1867EXAMPLE : example printf; shows an example
1868"
1869{
1870  write("", sprintf(fmt, #));
1871}
1872example
1873{ "EXAMPLE:"; echo=2;
1874  ring r=0,(x,y,z),dp;
1875  module m=[1,y],[0,x+z];
1876  intmat M=betti(mres(m,0));
1877  list l=r,m,matrix(M);
1878  printf("s:%s,l:%l",1,2);
1879  printf("s:%s",l);
1880  printf("s:%s",list(l));
1881  printf("2l:%2l",list(l));
1882  printf("%p",matrix(M));
1883  printf("%;",matrix(M));
1884  printf("%b",M);
1885}
1886
1887
1888proc fprintf(link l, string fmt, list #)
1889"SYNTAX:  @code{fprintf (} link_expression@code{,} string_expression @code{[,}
1890                any_expressions@code{] )}
1891RETURN:   none
1892PURPOSE:  @code{fprintf(l,fmt,...);} performs output formatting.
1893          The second argument is a format control string. Additional
1894          arguments may be required, depending on the content of the
1895          control string. A series of output characters is generated as
1896          directed by the control string; these characters are
1897          written to the link l.
1898          The control string @code{fmt} is simply text to be copied, except
1899          that the string may contain conversion specifications.@*
1900          Type @code{help print;} for a listing of valid conversion
1901          specifications. As an addition to the conversions of @code{print},
1902          the @code{%n} and @code{%2} conversion specification does not
1903          consume an additional argument, but simply generates a newline
1904          character.
1905NOTE:     If one of the additional arguments is a list, then it should be
1906          enclosed once more into a @code{list()} command, since passing
1907          a list as an argument flattens the list by one level.
1908SEE ALSO: sprintf, printf, print, string
1909EXAMPLE : example fprintf; shows an example
1910"
1911{
1912  write(l, sprintf(fmt, #));
1913}
1914example
1915{ "EXAMPLE:"; echo=2;
1916  ring r=0,(x,y,z),dp;
1917  module m=[1,y],[0,x+z];
1918  intmat M=betti(mres(m,0));
1919  list l=r,m,M;
1920  link li="";   // link to stdout
1921  fprintf(li,"s:%s,l:%l",1,2);
1922  fprintf(li,"s:%s",l);
1923  fprintf(li,"s:%s",list(l));
1924  fprintf(li,"2l:%2l",list(l));
1925  fprintf(li,"%p",list(l));
1926  fprintf(li,"%;",list(l));
1927  fprintf(li,"%b",M);
1928}
1929
1930//////////////////////////////////////////////////////////////////////////
1931
1932/*
1933proc minres(list #)
1934{
1935  if (size(#) == 2)
1936  {
1937    if (typeof(#[1]) == "ideal" || typeof(#[1]) == "module")
1938    {
1939      if (typeof(#[2] == "int"))
1940      {
1941        return (res(#[1],#[2],1));
1942      }
1943    }
1944  }
1945
1946  if (typeof(#[1]) == "resolution")
1947  {
1948    return minimizeres(#[1]);
1949  }
1950  else
1951  {
1952    return minimizeres(#);
1953  }
1954
1955}
1956*/
1957///////////////////////////////////////////////////////////////////////////////
1958
1959proc weightKB(def stc, int dd, list wim)
1960"SYNTAX: @code{weightKB (} module_expression@code{,} int_expression @code{,}
1961            list_expression @code{)}@*
1962         @code{weightKB (} ideal_expression@code{,} int_expression@code{,}
1963            list_expression @code{)}
1964RETURN:  the same as the input type of the first argument
1965PURPOSE: If @code{I,d,wim} denotes the three arguments then weightKB
1966         computes the weighted degree- @code{d} part of a vector space basis
1967         (consisting of monomials) of the quotient ring, resp. of the
1968         quotient module, modulo @code{I} w.r.t. weights given by @code{wim}
1969         The information about the weights is given as a list of two intvec:
1970            @code{wim[1]} weights for all variables (positive),
1971            @code{wim[2]} weights for the module generators.
1972NOTE:    This is a generalization of the command @code{kbase} with the same
1973         first two arguments.
1974SEE ALSO: kbase
1975EXAMPLE: example weightKB; shows an example
1976"
1977{
1978  if(checkww(wim)){ERROR("wrong weights";);}
1979  kbclass();
1980  wwtop=wim[1];
1981  stc=interred(lead(stc));
1982  if(typeof(stc)=="ideal")
1983  {
1984    stdtop=stc;
1985    ideal out=widkbase(dd);
1986    delkbclass();
1987    out=simplify(out,2); // delete 0
1988    return(out);
1989  }
1990  list mbase=kbprepare(stc);
1991  module mout;
1992  int im,ii;
1993  if(size(wim)>1){mmtop=wim[2];}
1994  else{mmtop=0;}
1995  for(im=size(mbase);im>0;im--)
1996  {
1997    stdtop=mbase[im];
1998    if(im>size(mmtop)){ii=dd;}
1999    else{ii=dd-mmtop[im];}
2000    mout=mout+widkbase(ii)*gen(im);
2001  }
2002  delkbclass();
2003  mout=simplify(mout,2); // delete 0
2004  return(mout);
2005}
2006example
2007{ "EXAMPLE:"; echo=2;
2008  ring R=0, (x,y), wp(1,2);
2009  weightKB(ideal(0),3,intvec(1,2));
2010}
2011
2012///////////////////////////////////////////////////////////////////////////////
2013
2014proc datetime()
2015"SYNTAX: @code{datetime ()}
2016RETURN:  string
2017PURPOSE: return the curent date and time as a string
2018EXAMPLE: example datetime; shows an example
2019"
2020{
2021  return(read("|: date"));
2022}
2023example
2024{ "EXAMPLE:"; echo=2;
2025  datetime();
2026}
2027
2028///////////////////////////////////////////////////////////////////////////////
2029// construct global values
2030static proc kbclass()
2031{
2032  intvec wwtop,mmtop;
2033  export (wwtop,mmtop);
2034  ideal stdtop,kbtop;
2035  export (stdtop,kbtop);
2036}
2037// delete global values
2038static proc delkbclass()
2039{
2040  kill wwtop,mmtop;
2041  kill stdtop,kbtop;
2042}
2043//  select parts of the modul
2044static proc kbprepare(module mstc)
2045{
2046  list rr;
2047  ideal kk;
2048  int i1,i2;
2049  mstc=transpose(mstc);
2050  for(i1=ncols(mstc);i1>0;i1--)
2051  {
2052    kk=0;
2053    for(i2=nrows(mstc[i1]);i2>0;i2--)
2054    {
2055      kk=kk+mstc[i1][i2];
2056    }
2057    rr[i1]=kk;
2058  }
2059  return(rr);
2060}
2061//  check for weights
2062static proc checkww(list vv)
2063{
2064  if(typeof(vv[1])!="intvec"){return(1);}
2065  intvec ww=vv[1];
2066  int mv=nvars(basering);
2067  if(size(ww)<mv){return(1);}
2068  while(mv>0)
2069  {
2070    if(ww[mv]<=0){return(1);}
2071    mv--;
2072  }
2073  if(size(vv)>1)
2074  {
2075    if(typeof(vv[2])!="intvec"){return(1);}
2076  }
2077  return(0);
2078}
2079///////////////////////////////////////////////////////
2080// The "Caller" for ideals
2081//    dd   - the degree of the result
2082static proc widkbase(int dd)
2083{
2084  if((size(stdtop)==1)&&(deg(stdtop[1])==0)){return(0);}
2085  if(dd<=0)
2086  {
2087    if(dd<0){return(0);}
2088    else{return(1);}
2089  }
2090  int m1,m2;
2091  m1=nvars(basering);
2092  while(wwtop[m1]>dd)
2093  {
2094    m1--;
2095    if(m1==0){return(0);}
2096  }
2097  attrib(stdtop,"isSB",1);
2098  poly mo=1;
2099  if(m1==1)
2100  {
2101    m2=dd div wwtop[1];
2102    if((m2*wwtop[1])==dd)
2103    {
2104      mo=var(1)^m2;
2105      if(reduce(mo,stdtop)==mo){return(mo);}
2106      else{return(0);}
2107    }
2108  }
2109  kbtop=0;
2110  m2=dd;
2111  weightmon(m1-1,m2,mo);
2112  while(m2>=wwtop[m1])
2113  {
2114    m2=m2-wwtop[m1];
2115    mo=var(m1)*mo;
2116    if(m2==0)
2117    {
2118      if((mo!=0) and (reduce(mo,stdtop)==mo))
2119      {
2120        kbtop[ncols(kbtop)+1]=mo;
2121        return(kbtop);
2122      }
2123    }
2124    weightmon(m1-1,m2,mo);
2125  }
2126  return(kbtop);
2127}
2128/////////////////////////////////////////////////////////
2129// the recursive procedure
2130//    va     - number of the variable
2131//    drest  - rest of the degree
2132//    mm     - the candidate
2133static proc weightmon(int va, int drest, poly mm)
2134{
2135  while(wwtop[va]>drest)
2136  {
2137    va--;
2138    if(va==0){return();}
2139  }
2140  int m2;
2141  if(va==1)
2142  {
2143    m2=drest div wwtop[1];
2144    if((m2*wwtop[1])==drest)
2145    {
2146      mm=var(1)^m2*mm;
2147      if ((mm!=0) and (reduce(mm,stdtop)==mm))
2148      {
2149        kbtop[ncols(kbtop)+1]=mm;
2150      }
2151    }
2152    return();
2153  }
2154  m2=drest;
2155  if ((mm!=0) and (reduce(mm,stdtop)==mm))
2156  {
2157    weightmon(va-1,m2,mm);
2158  }
2159  while(m2>=wwtop[va])
2160  {
2161    m2=m2-wwtop[va];
2162    mm=var(va)*mm;
2163    if(m2==0)
2164    {
2165      if ((mm!=0) and (reduce(mm,stdtop)==mm))
2166      {
2167        kbtop[ncols(kbtop)+1]=mm;
2168        return();
2169      }
2170    }
2171     if ((mm!=0) and (reduce(mm,stdtop)==mm))
2172     {
2173       weightmon(va-1,m2,mm);
2174     }
2175  }
2176  return();
2177}
2178example
2179{ "EXAMPLE:"; echo=2;
2180  ring r=0,(x,y,z),dp;
2181  ideal i = x6,y4,xyz;
2182  intvec w = 2,3,6;
2183  weightKB(i, 12, list(w));
2184}
2185
2186///////////////////////////////////////////////////////////////////////////////
2187proc max(def i,list #)
2188"SYNTAX: max (i_1, ..., i_k)
2189TYPE:    same as type of i_1, ..., i_k resp.
2190PURPOSE: returns the maximum for any arguments of a type
2191         for which '>' is defined
2192SEE ALSO: min
2193EXAMPLE: example max; shows an example"
2194{
2195  def maximum = i;
2196  for (int j=1; j<=size(#); j++)
2197  {
2198    if(#[j]>maximum)
2199    {
2200      maximum = #[j];
2201    }
2202  }
2203  return(maximum);
2204}
2205example
2206{ "EXAMPLE:"; echo=2;
2207  // biggest int
2208  max(2,3);
2209  max(1,4,3);
2210  // lexicographically biggest intvec
2211  max(intvec(1,2),intvec(0,1),intvec(1,1));
2212  // polynopmial with biggest leading monomial
2213  ring r = 0,x,dp;
2214  max(x+1,x2+x);
2215}
2216///////////////////////////////////////////////////////////////////////////////
2217proc min(def i,list #)
2218"SYNTAX: min (i_1, ..., i_k)
2219TYPE:    same as type of i_1, ..., i_k resp.
2220PURPOSE: returns the minimum for any arguments of a type
2221         for which '>' is defined
2222SEE ALSO: max
2223EXAMPLE: example min; shows an example"
2224{
2225  def minimum = i;
2226  for (int j=1; j<=size(#); j++)
2227  {
2228    if(#[j]<minimum)
2229    {
2230      minimum = #[j];
2231    }
2232  }
2233  return(minimum);
2234}
2235example
2236{ "EXAMPLE:"; echo=2;
2237  // smallest int
2238  min(2,3);
2239  min(1,4,3);
2240  // lexicographically smallest intvec
2241  min(intvec(1,2),intvec(0,1),intvec(1,1));
2242  // polynopmial with smallest leading monomial
2243  ring r = 0,x,dp;
2244  min(x+1,x2+x);
2245}
2246
2247
2248///////////////////////////////////////////////////////////////////////////////
2249/*
2250                                Versuche:
2251///////////////////////////////////////////////////////////////////////////////
2252proc downsizeSB (def I, list #)
2253"USAGE:   downsizeSB(I [,l]); I ideal, l list of integers [default: l=0]
2254RETURN:  intvec, say v, with v[j] either 1 or 0. We have v[j]=1 if
2255         leadmonom(I[j]) is divisible by some leadmonom(I[k]) or if
2256         leadmonom(i[j]) == leadmonom(i[k]) and l[j] >= l[k], with k!=j.
2257PURPOSE: The procedure is applied in a situation where the standard basis
2258         computation in the basering R is done via a conversion through an
2259         overring Phelp with additional variables and where a direct
2260         imap from Phelp to R is too expensive.
2261         Assume Phelp is created by the procedure @code{par2varRing} or
2262         @code{hilbRing} and IPhelp is a SB in Phelp [ with l[j]=
2263         length(IPhelp(j)) or any other integer reflecting the complexity
2264         of a IPhelp[j] ]. Let I = lead(IPhelp) mapped to R and compute
2265         v = downsizeSB(imap(Phelp,I),l) in R. Then, if Ihelp[j] is deleted
2266         for all j with v[j]=1, we can apply imap to the remaining generators
2267         of Ihelp and still get SB in R  (in general not minimal).
2268EXAMPLE: example downsizeSB; shows an example"
2269{
2270   int k,j;
2271   intvec v,l;
2272   poly M,N,W;
2273   int c=size(I);
2274   if( size(#) != 0 )
2275   {
2276     if ( typeof(#[1]) == "intvec" )
2277     {
2278       l = #[1];
2279     }
2280     else
2281     {
2282        ERROR("// 2nd argument must be an intvec");
2283     }
2284   }
2285
2286   l[c+1]=0;
2287   v[c]=0;
2288
2289   j=0;
2290   while(j<c-1)
2291   {
2292     j++;
2293     M = leadmonom(I[j]);
2294     if( M != 0 )
2295     {
2296        for( k=j+1; k<=c; k++ )
2297        {
2298          N = leadmonom(I[k]);
2299          if( N != 0 )
2300          {
2301             if( (M==N) && (l[j]>l[k]) )
2302             {
2303               I[j]=0;
2304               v[j]=1;
2305               break;
2306             }
2307             if( (M==N) && (l[j]<=l[k]) || N/M != 0 )
2308             {
2309               I[k]=0;
2310               v[k]=1;
2311             }
2312          }
2313        }
2314      }
2315   }
2316   return(v);
2317}
2318example
2319{ "EXAMPLE:"; echo = 2;
2320   ring  r = 0,(x,y,z,t),(dp(3),dp);
2321   ideal i = x+y+z+t,xy+yz+xt+zt,xyz+xyt+xzt+yzt,xyzt-t4;
2322   ideal Id = std(i);
2323   ideal I = lead(Id);  I;
2324   ring S = (0,t),(x,y,z),dp;
2325   downsizeSB(imap(r,I));
2326   //Id[5] can be deleted, we still have a SB of i in the ring S
2327
2328   ring R = (0,x),(y,z,u),lp;
2329   ideal i = x+y+z+u,xy+xu+yz+zu,xyz+xyu+xzu+yzu,xyzu-1;
2330   def Phelp = par2varRing()[1];
2331   setring Phelp;
2332   ideal IPhelp = std(imap(R,i));
2333   ideal I = lead(IPhelp);
2334   setring R;
2335   ideal I = imap(Phelp,I); I;
2336   intvec v = downsizeSB(I); v;
2337}
2338///////////////////////////////////////////////////////////////////////////
2339// PROBLEM: Die Prozedur funktioniert nur fuer Ringe die global bekannt
2340//          sind, also interaktiv, aber nicht aus einer Prozedur.
2341//          Z.B. funktioniert example imapDownsize; nicht
2342
2343proc imapDownsize (string R, string I)
2344"SYNTAX: @code{imapDownsize (} string @code{,} string @code{)} *@
2345         First string must be the string of the name of a ring, second
2346         string must be the string of the name of an object in the ring.
2347TYPE:    same type as the object with name the second string
2348PURPOSE: maps the object given by the second string to the basering.
2349         If R resp. I are the first resp. second string, then
2350         imapDownsize(R,I) is equivalent to simplify(imap(`R`,`I`),34).
2351NOTE:    imapDownsize is usually faster than imap if `I` is large and if
2352         simplify has a great effect, since the procedure maps only those
2353         generators from `I` which are not killed by simplify( - ,34).
2354         This is useful if `I` is a standard bases for a block ordering of
2355         `R` and if some variables from the last block in `R` are mapped
2356         to parameters. Then the returned result is a standard basis in
2357         the basering.
2358SEE ALSO: imap, fetch, map
2359EXAMPLE: example imapDownsize; shows an example"
2360{
2361       def BR = basering;
2362       int k;
2363
2364       setring `R`;
2365       def @leadI@ = lead(`I`);
2366       int s = ncols(@leadI@);
2367       setring BR;
2368       ideal @leadI@ = simplify(imap(`R`,@leadI@),32);
2369       intvec vi;
2370       for (k=1; k<=s; k++)
2371       {
2372         vi[k] = @leadI@[k]==0;
2373       }
2374       kill @leadI@;
2375
2376       setring `R`;
2377       kill @leadI@;
2378       for (k=1;  k<=s; k++)
2379       {
2380           if( vi[k]==1 )
2381           {
2382              `I`[k]=0;
2383           }
2384       }
2385       `I` = simplify(`I`,2);
2386
2387       setring BR;
2388       return(imap(`R`,`I`));
2389}
2390example
2391{ "EXAMPLE:"; echo = 2;
2392   ring  r = 0,(x,y,z,t),(dp(3),dp);
2393   ideal i = x+y+z+t,xy+yz+xt+zt,xyz+xyt+xzt+yzt,xyzt-1;
2394   i = std(i); i;
2395
2396   ring s = (0,t),(x,y,z),dp;
2397   imapDownsize("r","i");     //i[5] is omitted since lead(i[2]) | lead(i[5])
2398}
2399///////////////////////////////////////////////////////////////////////////////
2400//die folgende proc war fuer groebner mit fglm vorgesehen, ist aber zu teuer.
2401//Um die projektive Dimension korrekt zu berechnen, muss man aber teuer
2402//voerher ein SB bzgl. einer Gradordnung berechnen und dann homogenisieren.
2403//Sonst koennen hoeherdimensionale Komponenten in Unendlich entstehen
2404
2405proc projInvariants(ideal i,list #)
2406"SYNTAX: @code{projInvariants (} ideal_expression @code{)} @*
2407         @code{projInvariants (} ideal_expression@code{,} list of string_expres          sions@code{)}
2408TYPE:    list, say L, with L[1] and L[2] of type int and L[3] of type intvec
2409PURPOSE: Computes the (projective) dimension (L[1]), degree (L[2]) and the
2410         first Hilbert series (L[3], as intvec) of the homogenized ideal
2411         in the ring given by the procedure @code{hilbRing} with global
2412         ordering dp (resp. wp if the variables have weights >1)
2413         If an argument of type string @code{\"std\"} resp. @code{\"slimgb\"}
2414         is given, the standard basis computatuion uses @code{std} or
2415         @code{slimgb}, otherwise a heuristically chosen method (default)
2416NOTE:    Homogenized means weighted homogenized with respect to the weights
2417         w[i] of the variables var(i) of the basering. The returned dimension,
2418         degree and Hilbertseries are the respective invariants of the
2419         projective variety defined by the homogenized ideal. The dimension
2420         is equal to the (affine) dimension of the ideal in the basering
2421         (degree and Hilbert series make only sense for homogeneous ideals).
2422SEE ALSO: dim, dmult, hilb
2423KEYWORDS: dimension, degree, Hilbert function
2424EXAMPLE: example projInvariants;  shows an example"
2425{
2426  def P = basering;
2427  int p_opt;
2428  string s_opt = option();
2429  if (find(option(), "prot"))  { p_opt = 1; }
2430
2431//---------------- check method and clear denomintors --------------------
2432  int k;
2433  string method;
2434  for (k=1; k<=size(#); k++)
2435  {
2436     if (typeof(#[k]) == "string")
2437     {
2438       method = method + "," + #[k];
2439     }
2440  }
2441
2442  if (npars(P) > 0)             //clear denominators of parameters
2443  {
2444    for( k=ncols(i); k>0; k-- )
2445    {
2446       i[k]=cleardenom(i[k]);
2447    }
2448  }
2449
2450//------------------------ change to hilbRing ----------------------------
2451     list hiRi = hilbRing(i);
2452     intvec W = hiRi[2];
2453     def Philb = hiRi[1];      //note: Philb is no qring and the predefined
2454     setring Philb;            //ideal Id(1) in Philb is homogeneous
2455     int di, de;               //for dimension, degree
2456     intvec hi;                //for hilbert series
2457
2458//-------- compute Hilbert function of homogenized ideal in Philb ---------
2459//Philb has only 1 block. There are three cases
2460
2461     string algorithm;       //possibilities: std, slimgb, stdorslimgb
2462     //define algorithm:
2463     if( find(method,"std") && !find(method,"slimgb") )
2464     {
2465        algorithm = "std";
2466     }
2467     if( find(method,"slimgb") && !find(method,"std") )
2468     {
2469         algorithm = "slimgb";
2470     }
2471     if( find(method,"std") && find(method,"slimgb") ||
2472         (!find(method,"std") && !find(method,"slimgb")) )
2473     {
2474         algorithm = "stdorslimgb";
2475     }
2476
2477     if ( algorithm=="std" || ( algorithm=="stdorslimgb" && char(P)>0 ) )
2478     {
2479        if (p_opt) {"std in ring " + string(Philb);}
2480        Id(1) = std(Id(1));
2481        di = dim(Id(1))-1;
2482        de = mult(Id(1));
2483        hi = hilb( Id(1),1,W );
2484     }
2485     if ( algorithm=="slimgb" || ( algorithm=="stdorslimgb" && char(P)==0 ) )
2486     {
2487        if (p_opt) {"slimgb in ring " + string(Philb);}
2488        Id(1) = slimgb(Id(1));
2489        di = dim( Id(1) );
2490        if (di > -1)
2491        {
2492           di = di-1;
2493        }
2494        de = mult( Id(1) );
2495        hi = hilb( Id(1),1,W );
2496     }
2497     kill Philb;
2498     list L = di,de,hi;
2499     return(L);
2500}
2501example
2502{ "EXAMPLE:"; echo = 2;
2503   ring r = 32003,(x,y,z),lp;
2504   ideal i = y2-xz,x2-z;
2505   projInvariants(i);
2506
2507   ring R = (0),(x,y,z,u,v),lp;
2508   //minpoly = x2+1;
2509   ideal i = x2+1,x2+y+z+u+v,xyzuv-1;
2510   projInvariants(i);
2511   qring S =std(x2+1);
2512   ideal i = imap(R,i);
2513   projInvariants(i);
2514}
2515
2516*/
2517///////////////////////////////////////////////////////////////////////////////
2518//                           EXAMPLES
2519///////////////////////////////////////////////////////////////////////////////
2520/*
2521example stdfglm;
2522example stdhilb;
2523example groebner;
2524example res;
2525example sprintf;
2526example fprintf;
2527example printf;
2528example weightKB;
2529example qslimgb;
2530example par2varRing;
2531*/
2532static proc mod_init()
2533{
2534  //int pagelength=24;
2535  //exportto(Top,pagelength);
2536}
Note: See TracBrowser for help on using the repository browser.