source: git/Singular/LIB/standard.lib @ 1ebcb0

spielwiese
Last change on this file since 1ebcb0 was 1ebcb0, checked in by Hans Schoenemann <hannes@…>, 9 years ago
library cleanup: moved min, max to standard.lib, removed varNum (->rvar)
  • Property mode set to 100644
File size: 71.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,j)               maximum of i and j
20 min(i,j)               minimum of i and j
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(homog(m)==1)
1393   {
1394      resolution re;
1395      if (((i==0) or (i>=nvars(basering))) && (typeof(m) != "module") && (nvars(basering)>1))
1396      {
1397        //LaScala for the homogeneous case and i == 0
1398        if (p_opt) { "using lres";}
1399        re=lres(m,i);
1400        if(size(#)>2)
1401        {
1402           re=minres(re);
1403        }
1404      }
1405      else
1406      {
1407        if(size(#)>2)
1408        {
1409          if (p_opt) { "using mres";}
1410          re=mres(m,i);
1411        }
1412        else
1413        {
1414          if (p_opt) { "using sres";}
1415          re=sres(std(m),i);
1416        }
1417      }
1418      return(re);
1419   }
1420
1421   //mres for the global non homogeneous case
1422   if(find(ordstr(P),"s")==0)
1423   {
1424      string ri= "ring Phelp ="
1425                  +string(char(P))+",("+varstr_P+"),(dp,C);";
1426      ri = ri + "minpoly = "+string(minpoly) + ";";
1427      execute(ri);
1428      def m=imap(P,m);
1429      if (p_opt) { "using mres in another ring";}
1430      list re=mres(m,i);
1431      setring P;
1432      resolution result=imap(Phelp,re);
1433      if (size(#) > 2) {result = minres(result);}
1434      return(result);
1435   }
1436
1437   //sres for the local case and not minimal resolution
1438   if(size(#)<=2)
1439   {
1440      string ri= "ring Phelp ="
1441                  +string(char(P))+",("+varstr_P+"),(ls,c);";
1442      ri = ri + "minpoly = "+string(minpoly) + ";";
1443      execute(ri);
1444      def m=imap(P,m);
1445      m=std(m);
1446      if (p_opt) { "using sres in another ring";}
1447      list re=sres(m,i);
1448      setring P;
1449      resolution result=imap(Phelp,re);
1450      return(result);
1451   }
1452
1453   //mres for the local case and minimal resolution
1454   string ri= "ring Phelp ="
1455                  +string(char(P))+",("+varstr_P+"),(ls,C);";
1456   ri = ri + "minpoly = "+string(minpoly) + ";";
1457   execute(ri);
1458   def m=imap(P,m);
1459    if (p_opt) { "using mres in another ring";}
1460   list re=mres(m,i);
1461   setring P;
1462   resolution result=imap(Phelp,re);
1463   result = minres(result);
1464   return(result);
1465}
1466example
1467{"EXAMPLE:"; echo = 2;
1468  ring r=0,(x,y,z),dp;
1469  ideal i=xz,yz,x3-y3;
1470  def l=res(i,0); // homogeneous ideal: uses lres
1471  l;
1472  print(betti(l), "betti"); // input to betti may be of type resolution
1473  l[2];         // element access may take some time
1474  i=i,x+1;
1475  l=res(i,0);   // inhomogeneous ideal: uses mres
1476  l;
1477  ring rs=0,(x,y,z),ds;
1478  ideal i=imap(r,i);
1479  def l=res(i,0); // local ring not minimized: uses sres
1480  l;
1481  res(i,0,0);     // local ring and minimized: uses mres
1482}
1483/////////////////////////////////////////////////////////////////////////
1484
1485proc quot (def m1,def m2,list #)
1486"SYNTAX: @code{quot (} module_expression@code{,} module_expression @code{)}
1487         @*@code{quot (} module_expression@code{,} module_expression@code{,}
1488            int_expression @code{)}
1489         @*@code{quot (} ideal_expression@code{,} ideal_expression @code{)}
1490         @*@code{quot (} ideal_expression@code{,} ideal_expression@code{,}
1491            int_expression @code{)}
1492TYPE:    ideal
1493SYNTAX:  @code{quot (} module_expression@code{,} ideal_expression @code{)}
1494TYPE:    module
1495PURPOSE: computes the quotient of the 1st and the 2nd argument.
1496         If a 3rd argument @code{n} is given the @code{n}-th method is used
1497         (@code{n}=1...5).
1498SEE ALSO: quotient
1499EXAMPLE: example quot; shows an example"
1500{
1501  if (((typeof(m1)!="ideal") and (typeof(m1)!="module"))
1502     or ((typeof(m2)!="ideal") and (typeof(m2)!="module")))
1503  {
1504    "USAGE:   quot(m1, m2[, n]); m1, m2 two submodules of k^s,";
1505    "         n (optional) integer (1<= n <=5)";
1506    "RETURN:  the quotient of m1 and m2";
1507    "EXAMPLE: example quot; shows an example";
1508    return();
1509  }
1510  if (typeof(m1)!=typeof(m2))
1511  {
1512    return(quotient(m1,m2));
1513  }
1514  if (size(#)>0)
1515  {
1516    if (typeof(#[1])=="int" )
1517    {
1518      return(quot1(m1,m2,#[1]));
1519    }
1520  }
1521  else
1522  {
1523    return(quot1(m1,m2,2));
1524  }
1525}
1526example
1527{ "EXAMPLE:"; echo = 2;
1528  ring r=181,(x,y,z),(c,ls);
1529  ideal id1=maxideal(4);
1530  ideal id2=x2+xyz,y2-z3y,z3+y5xz;
1531  option(prot);
1532  ideal id3=quotient(id1,id2);
1533  id3;
1534  ideal id4=quot(id1,id2,1);
1535  id4;
1536  ideal id5=quot(id1,id2,2);
1537  id5;
1538}
1539
1540static proc quot1 (module m1, module m2,int n)
1541"USAGE:   quot1(m1, m2, n); m1, m2 two submodules of k^s,
1542         n integer (1<= n <=5)
1543RETURN:  the quotient of m1 and m2
1544EXAMPLE: example quot1; shows an example"
1545{
1546  if (n==1)
1547  {
1548    return(quotient1(m1,m2));
1549  }
1550  else
1551  {
1552    if (n==2)
1553    {
1554      return(quotient2(m1,m2));
1555    }
1556    else
1557    {
1558      if (n==3)
1559      {
1560        return(quotient3(m1,m2));
1561      }
1562      else
1563      {
1564        if (n==4)
1565        {
1566          return(quotient4(m1,m2));
1567        }
1568        else
1569        {
1570          if (n==5)
1571          {
1572            return(quotient5(m1,m2));
1573          }
1574          else
1575          {
1576            return(quotient(m1,m2));
1577          }
1578        }
1579      }
1580    }
1581  }
1582}
1583example
1584{ "EXAMPLE:"; echo = 2;
1585  ring r=181,(x,y,z),(c,ls);
1586  ideal id1=maxideal(4);
1587  ideal id2=x2+xyz,y2-z3y,z3+y5xz;
1588  option(prot);
1589  ideal id6=quotient(id1,id2);
1590  id6;
1591  ideal id7=quot1(id1,id2,1);
1592  id7;
1593  ideal id8=quot1(id1,id2,2);
1594  id8;
1595}
1596
1597static proc quotient0(module a,module b)
1598{
1599  module mm=b+a;
1600  resolution rs=lres(mm,0);
1601  list I=list(rs);
1602  matrix M=I[2];
1603  matrix A[1][nrows(M)]=M[1..nrows(M),1];
1604  ideal i=A;
1605  return (i);
1606}
1607proc quotient1(module a,module b)  //17sec
1608"USAGE:   quotient1(m1, m2); m1, m2 two submodules of k^s,
1609RETURN:  the quotient of m1 and m2"
1610{
1611  int i;
1612  a=std(a);
1613  module dummy;
1614  module B=NF(b,a)+dummy;
1615  ideal re=quotient(a,module(B[1]));
1616  for(i=2;i<=ncols(B);i++)
1617  {
1618     re=intersect1(re,quotient(a,module(B[i])));
1619  }
1620  return(re);
1621}
1622proc quotient2(module a,module b)    //13sec
1623"USAGE:   quotient2(m1, m2); m1, m2 two submodules of k^s,
1624RETURN:  the quotient of m1 and m2"
1625{
1626  a=std(a);
1627  module dummy;
1628  module bb=NF(b,a)+dummy;
1629  int i=ncols(bb);
1630  ideal re=quotient(a,module(bb[i]));
1631  bb[i]=0;
1632  module temp;
1633  module temp1;
1634  module bbb;
1635  int mx;
1636  i=i-1;
1637  while (1)
1638  {
1639    if (i==0) break;
1640    temp = a+bb*re;
1641    temp1 = lead(interred(temp));
1642    mx=ncols(a);
1643    if (ncols(temp1)>ncols(a))
1644    {
1645      mx=ncols(temp1);
1646    }
1647    temp1 = matrix(temp1,1,mx)-matrix(lead(a),1,mx);
1648    temp1 = dummy+temp1;
1649    if (deg(temp1[1])<0) break;
1650    re=intersect1(re,quotient(a,module(bb[i])));
1651    bb[i]=0;
1652    i = i-1;
1653  }
1654  return(re);
1655}
1656proc quotient3(module a,module b)   //89sec
1657"USAGE:   quotient3(m1, m2); m1, m2 two submodules of k^s,
1658         only for global rings
1659RETURN:  the quotient of m1 and m2"
1660{
1661  string s="ring @newr=("+charstr(basering)+
1662           "),("+varstr(basering)+",@t,@w),dp;";
1663  def @newP=basering;
1664  execute(s);
1665  module b=imap(@newP,b);
1666  module a=imap(@newP,a);
1667  int i;
1668  int j=ncols(b);
1669  vector @b;
1670  for(i=1;i<=j;i++)
1671  {
1672     @b=@b+@t^(i-1)*@w^(j-i+1)*b[i];
1673  }
1674  ideal re=quotient(a,module(@b));
1675  setring @newP;
1676  ideal re=imap(@newr,re);
1677  return(re);
1678}
1679proc quotient5(module a,module b)   //89sec
1680"USAGE:   quotient5(m1, m2); m1, m2 two submodules of k^s,
1681         only for global rings
1682RETURN:  the quotient of m1 and m2"
1683{
1684  string s="ring @newr=("+charstr(basering)+
1685           "),("+varstr(basering)+",@t),dp;";
1686  def @newP=basering;
1687  execute(s);
1688  module b=imap(@newP,b);
1689  module a=imap(@newP,a);
1690  int i;
1691  int j=ncols(b);
1692  vector @b;
1693  for(i=1;i<=j;i++)
1694  {
1695     @b=@b+@t^(i-1)*b[i];
1696  }
1697  @b=homog(@b,@w);
1698  ideal re=quotient(a,module(@b));
1699  setring @newP;
1700  ideal re=imap(@newr,re);
1701  return(re);
1702}
1703proc quotient4(module a,module b)   //95sec
1704"USAGE:   quotient4(m1, m2); m1, m2 two submodules of k^s,
1705         only for global rings
1706RETURN:  the quotient of m1 and m2"
1707{
1708  string s="ring @newr=("+charstr(basering)+
1709           "),("+varstr(basering)+",@t),dp;";
1710  def @newP=basering;
1711  execute(s);
1712  module b=imap(@newP,b);
1713  module a=imap(@newP,a);
1714  int i;
1715  vector @b=b[1];
1716  for(i=2;i<=ncols(b);i++)
1717  {
1718     @b=@b+@t^(i-1)*b[i];
1719  }
1720  matrix sy=modulo(@b,a);
1721  ideal re=sy;
1722  setring @newP;
1723  ideal re=imap(@newr,re);
1724  return(re);
1725}
1726static proc intersect1(ideal i,ideal j)
1727{
1728  def R=basering;
1729  execute("ring gnir = ("+charstr(basering)+"),
1730                       ("+varstr(basering)+",@t),(C,dp);");
1731  ideal i=var(nvars(basering))*imap(R,i)+(var(nvars(basering))-1)*imap(R,j);
1732  ideal j=eliminate(i,var(nvars(basering)));
1733  setring R;
1734  map phi=gnir,maxideal(1);
1735  return(phi(j));
1736}
1737
1738//////////////////////////////////////////////////////////////////
1739///
1740/// sprintf, fprintf printf
1741///
1742proc sprintf(string fmt, list #)
1743"SYNTAX:  @code{sprintf (} string_expression @code{[,} any_expressions
1744               @code{] )}
1745RETURN:   string
1746PURPOSE:  @code{sprintf(fmt,...);} performs output formatting. The first
1747          argument is a format control string. Additional arguments may be
1748          required, depending on the content of the control string. A series
1749          of output characters is generated as directed by the control string;
1750          these characters are returned as a string. @*
1751          The control string @code{fmt} is simply text to be copied,
1752          except that the string may contain conversion specifications.@*
1753          Type @code{help print;} for a listing of valid conversion
1754          specifications. As an addition to the conversions of @code{print},
1755          the @code{%n} and @code{%2} conversion specification does not
1756          consume an additional argument, but simply generates a newline
1757          character.
1758NOTE:     If one of the additional arguments is a list, then it should be
1759          wrapped in an additional @code{list()} command, since passing a list
1760          as an argument flattens the list by one level.
1761SEE ALSO: fprintf, printf, print, string
1762EXAMPLE : example sprintf; shows an example
1763"
1764{
1765  int sfmt = size(fmt);
1766  if (sfmt  <= 1)
1767  {
1768    return (fmt);
1769  }
1770  int next, l, nnext;
1771  string ret;
1772  list formats = "%l", "%s", "%2l", "%2s", "%t", "%;", "%p", "%b", "%n", "%2";
1773  while (1)
1774  {
1775    if (size(#) <= 0)
1776    {
1777      return (ret + fmt);
1778    }
1779    nnext = 0;
1780    while (nnext < sfmt)
1781    {
1782      nnext = find(fmt, "%", nnext + 1);
1783      if (nnext == 0)
1784      {
1785        next = 0;
1786        break;
1787      }
1788      l = 1;
1789      while (l <= size(formats))
1790      {
1791        next = find(fmt, formats[l], nnext);
1792        if (next == nnext) break;
1793        l++;
1794      }
1795      if (next == nnext) break;
1796    }
1797    if (next == 0)
1798    {
1799      return (ret + fmt);
1800    }
1801    if (formats[l] != "%2" && formats[l] != "%n")
1802    {
1803      ret = ret + fmt[1, next - 1] + print(#[1], formats[l]);
1804      # = delete(#, 1);
1805    }
1806    else
1807    {
1808      ret = ret + fmt[1, next - 1] + print("", "%2s");
1809    }
1810    if (size(fmt) <= (next + size(formats[l]) - 1))
1811    {
1812      return (ret);
1813    }
1814    fmt = fmt[next + size(formats[l]), size(fmt)-next-size(formats[l]) + 1];
1815  }
1816}
1817example
1818{ "EXAMPLE:"; echo=2;
1819  ring r=0,(x,y,z),dp;
1820  module m=[1,y],[0,x+z];
1821  intmat M=betti(mres(m,0));
1822  list l = r, m, M;
1823  string s = sprintf("s:%s,%n l:%l", 1, 2); s;
1824  s = sprintf("s:%n%s", l); s;
1825  s = sprintf("s:%2%s", list(l)); s;
1826  s = sprintf("2l:%n%2l", list(l)); s;
1827  s = sprintf("%p", list(l)); s;
1828  s = sprintf("%;", list(l)); s;
1829  s = sprintf("%b", M); s;
1830}
1831
1832proc printf(string fmt, list #)
1833"SYNTAX:  @code{printf (} string_expression @code{[,} any_expressions@code{] )}
1834RETURN:   none
1835PURPOSE:  @code{printf(fmt,...);} performs output formatting. The first
1836          argument is a format control string. Additional arguments may be
1837          required, depending on the content of the control string. A series
1838          of output characters is generated as directed by the control string;
1839          these characters are displayed (i.e., printed to standard out). @*
1840          The control string @code{fmt} is simply text to be copied, except
1841          that the string may contain conversion specifications. @*
1842          Type @code{help print;} for a listing of valid conversion
1843          specifications. As an addition to the conversions of @code{print},
1844          the @code{%n} and @code{%2} conversion specification does not
1845          consume an additional argument, but simply generates a newline
1846          character.
1847NOTE:     If one of the additional arguments is a list, then it should be
1848          enclosed once more into a @code{list()} command, since passing a
1849          list as an argument flattens the list by one level.
1850SEE ALSO: sprintf, fprintf, print, string
1851EXAMPLE : example printf; shows an example
1852"
1853{
1854  write("", sprintf(fmt, #));
1855}
1856example
1857{ "EXAMPLE:"; echo=2;
1858  ring r=0,(x,y,z),dp;
1859  module m=[1,y],[0,x+z];
1860  intmat M=betti(mres(m,0));
1861  list l=r,m,matrix(M);
1862  printf("s:%s,l:%l",1,2);
1863  printf("s:%s",l);
1864  printf("s:%s",list(l));
1865  printf("2l:%2l",list(l));
1866  printf("%p",matrix(M));
1867  printf("%;",matrix(M));
1868  printf("%b",M);
1869}
1870
1871
1872proc fprintf(link l, string fmt, list #)
1873"SYNTAX:  @code{fprintf (} link_expression@code{,} string_expression @code{[,}
1874                any_expressions@code{] )}
1875RETURN:   none
1876PURPOSE:  @code{fprintf(l,fmt,...);} performs output formatting.
1877          The second argument is a format control string. Additional
1878          arguments may be required, depending on the content of the
1879          control string. A series of output characters is generated as
1880          directed by the control string; these characters are
1881          written to the link l.
1882          The control string @code{fmt} is simply text to be copied, except
1883          that the string may contain conversion specifications.@*
1884          Type @code{help print;} for a listing of valid conversion
1885          specifications. As an addition to the conversions of @code{print},
1886          the @code{%n} and @code{%2} conversion specification does not
1887          consume an additional argument, but simply generates a newline
1888          character.
1889NOTE:     If one of the additional arguments is a list, then it should be
1890          enclosed once more into a @code{list()} command, since passing
1891          a list as an argument flattens the list by one level.
1892SEE ALSO: sprintf, printf, print, string
1893EXAMPLE : example fprintf; shows an example
1894"
1895{
1896  write(l, sprintf(fmt, #));
1897}
1898example
1899{ "EXAMPLE:"; echo=2;
1900  ring r=0,(x,y,z),dp;
1901  module m=[1,y],[0,x+z];
1902  intmat M=betti(mres(m,0));
1903  list l=r,m,M;
1904  link li="";   // link to stdout
1905  fprintf(li,"s:%s,l:%l",1,2);
1906  fprintf(li,"s:%s",l);
1907  fprintf(li,"s:%s",list(l));
1908  fprintf(li,"2l:%2l",list(l));
1909  fprintf(li,"%p",list(l));
1910  fprintf(li,"%;",list(l));
1911  fprintf(li,"%b",M);
1912}
1913
1914//////////////////////////////////////////////////////////////////////////
1915
1916/*
1917proc minres(list #)
1918{
1919  if (size(#) == 2)
1920  {
1921    if (typeof(#[1]) == "ideal" || typeof(#[1]) == "module")
1922    {
1923      if (typeof(#[2] == "int"))
1924      {
1925        return (res(#[1],#[2],1));
1926      }
1927    }
1928  }
1929
1930  if (typeof(#[1]) == "resolution")
1931  {
1932    return minimizeres(#[1]);
1933  }
1934  else
1935  {
1936    return minimizeres(#);
1937  }
1938
1939}
1940*/
1941///////////////////////////////////////////////////////////////////////////////
1942
1943proc weightKB(def stc, int dd, list wim)
1944"SYNTAX: @code{weightKB (} module_expression@code{,} int_expression @code{,}
1945            list_expression @code{)}@*
1946         @code{weightKB (} ideal_expression@code{,} int_expression@code{,}
1947            list_expression @code{)}
1948RETURN:  the same as the input type of the first argument
1949PURPOSE: If @code{I,d,wim} denotes the three arguments then weightKB
1950         computes the weighted degree- @code{d} part of a vector space basis
1951         (consisting of monomials) of the quotient ring, resp. of the
1952         quotient module, modulo @code{I} w.r.t. weights given by @code{wim}
1953         The information about the weights is given as a list of two intvec:
1954            @code{wim[1]} weights for all variables (positive),
1955            @code{wim[2]} weights for the module generators.
1956NOTE:    This is a generalization of the command @code{kbase} with the same
1957         first two arguments.
1958SEE ALSO: kbase
1959EXAMPLE: example weightKB; shows an example
1960"
1961{
1962  if(checkww(wim)){ERROR("wrong weights";);}
1963  kbclass();
1964  wwtop=wim[1];
1965  stc=interred(lead(stc));
1966  if(typeof(stc)=="ideal")
1967  {
1968    stdtop=stc;
1969    ideal out=widkbase(dd);
1970    delkbclass();
1971    out=simplify(out,2); // delete 0
1972    return(out);
1973  }
1974  list mbase=kbprepare(stc);
1975  module mout;
1976  int im,ii;
1977  if(size(wim)>1){mmtop=wim[2];}
1978  else{mmtop=0;}
1979  for(im=size(mbase);im>0;im--)
1980  {
1981    stdtop=mbase[im];
1982    if(im>size(mmtop)){ii=dd;}
1983    else{ii=dd-mmtop[im];}
1984    mout=mout+widkbase(ii)*gen(im);
1985  }
1986  delkbclass();
1987  mout=simplify(mout,2); // delete 0
1988  return(mout);
1989}
1990example
1991{ "EXAMPLE:"; echo=2;
1992  ring R=0, (x,y), wp(1,2);
1993  weightKB(ideal(0),3,intvec(1,2));
1994}
1995
1996///////////////////////////////////////////////////////////////////////////////
1997
1998proc datetime()
1999"SYNTAX: @code{datetime ()}
2000RETURN:  string
2001PURPOSE: return the curent date and time as a string
2002EXAMPLE: example datetime; shows an example
2003"
2004{
2005  return(read("|: date"));
2006}
2007example
2008{ "EXAMPLE:"; echo=2;
2009  datetime();
2010}
2011
2012///////////////////////////////////////////////////////////////////////////////
2013// construct global values
2014static proc kbclass()
2015{
2016  intvec wwtop,mmtop;
2017  export (wwtop,mmtop);
2018  ideal stdtop,kbtop;
2019  export (stdtop,kbtop);
2020}
2021// delete global values
2022static proc delkbclass()
2023{
2024  kill wwtop,mmtop;
2025  kill stdtop,kbtop;
2026}
2027//  select parts of the modul
2028static proc kbprepare(module mstc)
2029{
2030  list rr;
2031  ideal kk;
2032  int i1,i2;
2033  mstc=transpose(mstc);
2034  for(i1=ncols(mstc);i1>0;i1--)
2035  {
2036    kk=0;
2037    for(i2=nrows(mstc[i1]);i2>0;i2--)
2038    {
2039      kk=kk+mstc[i1][i2];
2040    }
2041    rr[i1]=kk;
2042  }
2043  return(rr);
2044}
2045//  check for weights
2046static proc checkww(list vv)
2047{
2048  if(typeof(vv[1])!="intvec"){return(1);}
2049  intvec ww=vv[1];
2050  int mv=nvars(basering);
2051  if(size(ww)<mv){return(1);}
2052  while(mv>0)
2053  {
2054    if(ww[mv]<=0){return(1);}
2055    mv--;
2056  }
2057  if(size(vv)>1)
2058  {
2059    if(typeof(vv[2])!="intvec"){return(1);}
2060  }
2061  return(0);
2062}
2063///////////////////////////////////////////////////////
2064// The "Caller" for ideals
2065//    dd   - the degree of the result
2066static proc widkbase(int dd)
2067{
2068  if((size(stdtop)==1)&&(deg(stdtop[1])==0)){return(0);}
2069  if(dd<=0)
2070  {
2071    if(dd<0){return(0);}
2072    else{return(1);}
2073  }
2074  int m1,m2;
2075  m1=nvars(basering);
2076  while(wwtop[m1]>dd)
2077  {
2078    m1--;
2079    if(m1==0){return(0);}
2080  }
2081  attrib(stdtop,"isSB",1);
2082  poly mo=1;
2083  if(m1==1)
2084  {
2085    m2=dd div wwtop[1];
2086    if((m2*wwtop[1])==dd)
2087    {
2088      mo=var(1)^m2;
2089      if(reduce(mo,stdtop)==mo){return(mo);}
2090      else{return(0);}
2091    }
2092  }
2093  kbtop=0;
2094  m2=dd;
2095  weightmon(m1-1,m2,mo);
2096  while(m2>=wwtop[m1])
2097  {
2098    m2=m2-wwtop[m1];
2099    mo=var(m1)*mo;
2100    if(m2==0)
2101    {
2102      if((mo!=0) and (reduce(mo,stdtop)==mo))
2103      {
2104        kbtop[ncols(kbtop)+1]=mo;
2105        return(kbtop);
2106      }
2107    }
2108    weightmon(m1-1,m2,mo);
2109  }
2110  return(kbtop);
2111}
2112/////////////////////////////////////////////////////////
2113// the recursive procedure
2114//    va     - number of the variable
2115//    drest  - rest of the degree
2116//    mm     - the candidate
2117static proc weightmon(int va, int drest, poly mm)
2118{
2119  while(wwtop[va]>drest)
2120  {
2121    va--;
2122    if(va==0){return();}
2123  }
2124  int m2;
2125  if(va==1)
2126  {
2127    m2=drest div wwtop[1];
2128    if((m2*wwtop[1])==drest)
2129    {
2130      mm=var(1)^m2*mm;
2131      if ((mm!=0) and (reduce(mm,stdtop)==mm))
2132      {
2133        kbtop[ncols(kbtop)+1]=mm;
2134      }
2135    }
2136    return();
2137  }
2138  m2=drest;
2139  if ((mm!=0) and (reduce(mm,stdtop)==mm))
2140  {
2141    weightmon(va-1,m2,mm);
2142  }
2143  while(m2>=wwtop[va])
2144  {
2145    m2=m2-wwtop[va];
2146    mm=var(va)*mm;
2147    if(m2==0)
2148    {
2149      if ((mm!=0) and (reduce(mm,stdtop)==mm))
2150      {
2151        kbtop[ncols(kbtop)+1]=mm;
2152        return();
2153      }
2154    }
2155     if ((mm!=0) and (reduce(mm,stdtop)==mm))
2156     {
2157       weightmon(va-1,m2,mm);
2158     }
2159  }
2160  return();
2161}
2162example
2163{ "EXAMPLE:"; echo=2;
2164  ring r=0,(x,y,z),dp;
2165  ideal i = x6,y4,xyz;
2166  intvec w = 2,3,6;
2167  weightKB(i, 12, list(w));
2168}
2169
2170///////////////////////////////////////////////////////////////////////////////
2171proc max(def i,def j)
2172"SYNTAX: max (i, j)
2173TYPE:    same as type of i resp. j
2174PURPOSE: returns the maximum for any 2 arguments of a type
2175         for which '>' is defined
2176SEE ALSO: min
2177EXAMPLE: example max; shows an example"
2178{
2179  if(i>j){return(i);}
2180  return(j);
2181}
2182example
2183{ "EXAMPLE:"; echo=2;
2184  max(2,3);
2185  max(4,3);
2186}
2187///////////////////////////////////////////////////////////////////////////////
2188proc min(def i,def j)
2189"SYNTAX: min (i, j)
2190TYPE:    same as type of i resp. j
2191PURPOSE: returns the minimum for any 2 arguments of a type
2192         for which '>' is defined
2193SEE ALSO: max
2194EXAMPLE: example min; shows an example"
2195{
2196  if(i>j){return(j);}
2197  return(i);
2198}
2199example
2200{ "EXAMPLE:"; echo=2;
2201  min(2,3);
2202  min(4,3);
2203}
2204
2205
2206///////////////////////////////////////////////////////////////////////////////
2207/*
2208                                Versuche:
2209///////////////////////////////////////////////////////////////////////////////
2210proc downsizeSB (def I, list #)
2211"USAGE:   downsizeSB(I [,l]); I ideal, l list of integers [default: l=0]
2212RETURN:  intvec, say v, with v[j] either 1 or 0. We have v[j]=1 if
2213         leadmonom(I[j]) is divisible by some leadmonom(I[k]) or if
2214         leadmonom(i[j]) == leadmonom(i[k]) and l[j] >= l[k], with k!=j.
2215PURPOSE: The procedure is applied in a situation where the standard basis
2216         computation in the basering R is done via a conversion through an
2217         overring Phelp with additional variables and where a direct
2218         imap from Phelp to R is too expensive.
2219         Assume Phelp is created by the procedure @code{par2varRing} or
2220         @code{hilbRing} and IPhelp is a SB in Phelp [ with l[j]=
2221         length(IPhelp(j)) or any other integer reflecting the complexity
2222         of a IPhelp[j] ]. Let I = lead(IPhelp) mapped to R and compute
2223         v = downsizeSB(imap(Phelp,I),l) in R. Then, if Ihelp[j] is deleted
2224         for all j with v[j]=1, we can apply imap to the remaining generators
2225         of Ihelp and still get SB in R  (in general not minimal).
2226EXAMPLE: example downsizeSB; shows an example"
2227{
2228   int k,j;
2229   intvec v,l;
2230   poly M,N,W;
2231   int c=size(I);
2232   if( size(#) != 0 )
2233   {
2234     if ( typeof(#[1]) == "intvec" )
2235     {
2236       l = #[1];
2237     }
2238     else
2239     {
2240        ERROR("// 2nd argument must be an intvec");
2241     }
2242   }
2243
2244   l[c+1]=0;
2245   v[c]=0;
2246
2247   j=0;
2248   while(j<c-1)
2249   {
2250     j++;
2251     M = leadmonom(I[j]);
2252     if( M != 0 )
2253     {
2254        for( k=j+1; k<=c; k++ )
2255        {
2256          N = leadmonom(I[k]);
2257          if( N != 0 )
2258          {
2259             if( (M==N) && (l[j]>l[k]) )
2260             {
2261               I[j]=0;
2262               v[j]=1;
2263               break;
2264             }
2265             if( (M==N) && (l[j]<=l[k]) || N/M != 0 )
2266             {
2267               I[k]=0;
2268               v[k]=1;
2269             }
2270          }
2271        }
2272      }
2273   }
2274   return(v);
2275}
2276example
2277{ "EXAMPLE:"; echo = 2;
2278   ring  r = 0,(x,y,z,t),(dp(3),dp);
2279   ideal i = x+y+z+t,xy+yz+xt+zt,xyz+xyt+xzt+yzt,xyzt-t4;
2280   ideal Id = std(i);
2281   ideal I = lead(Id);  I;
2282   ring S = (0,t),(x,y,z),dp;
2283   downsizeSB(imap(r,I));
2284   //Id[5] can be deleted, we still have a SB of i in the ring S
2285
2286   ring R = (0,x),(y,z,u),lp;
2287   ideal i = x+y+z+u,xy+xu+yz+zu,xyz+xyu+xzu+yzu,xyzu-1;
2288   def Phelp = par2varRing()[1];
2289   setring Phelp;
2290   ideal IPhelp = std(imap(R,i));
2291   ideal I = lead(IPhelp);
2292   setring R;
2293   ideal I = imap(Phelp,I); I;
2294   intvec v = downsizeSB(I); v;
2295}
2296///////////////////////////////////////////////////////////////////////////
2297// PROBLEM: Die Prozedur funktioniert nur fuer Ringe die global bekannt
2298//          sind, also interaktiv, aber nicht aus einer Prozedur.
2299//          Z.B. funktioniert example imapDownsize; nicht
2300
2301proc imapDownsize (string R, string I)
2302"SYNTAX: @code{imapDownsize (} string @code{,} string @code{)} *@
2303         First string must be the string of the name of a ring, second
2304         string must be the string of the name of an object in the ring.
2305TYPE:    same type as the object with name the second string
2306PURPOSE: maps the object given by the second string to the basering.
2307         If R resp. I are the first resp. second string, then
2308         imapDownsize(R,I) is equivalent to simplify(imap(`R`,`I`),34).
2309NOTE:    imapDownsize is usually faster than imap if `I` is large and if
2310         simplify has a great effect, since the procedure maps only those
2311         generators from `I` which are not killed by simplify( - ,34).
2312         This is useful if `I` is a standard bases for a block ordering of
2313         `R` and if some variables from the last block in `R` are mapped
2314         to parameters. Then the returned result is a standard basis in
2315         the basering.
2316SEE ALSO: imap, fetch, map
2317EXAMPLE: example imapDownsize; shows an example"
2318{
2319       def BR = basering;
2320       int k;
2321
2322       setring `R`;
2323       def @leadI@ = lead(`I`);
2324       int s = ncols(@leadI@);
2325       setring BR;
2326       ideal @leadI@ = simplify(imap(`R`,@leadI@),32);
2327       intvec vi;
2328       for (k=1; k<=s; k++)
2329       {
2330         vi[k] = @leadI@[k]==0;
2331       }
2332       kill @leadI@;
2333
2334       setring `R`;
2335       kill @leadI@;
2336       for (k=1;  k<=s; k++)
2337       {
2338           if( vi[k]==1 )
2339           {
2340              `I`[k]=0;
2341           }
2342       }
2343       `I` = simplify(`I`,2);
2344
2345       setring BR;
2346       return(imap(`R`,`I`));
2347}
2348example
2349{ "EXAMPLE:"; echo = 2;
2350   ring  r = 0,(x,y,z,t),(dp(3),dp);
2351   ideal i = x+y+z+t,xy+yz+xt+zt,xyz+xyt+xzt+yzt,xyzt-1;
2352   i = std(i); i;
2353
2354   ring s = (0,t),(x,y,z),dp;
2355   imapDownsize("r","i");     //i[5] is omitted since lead(i[2]) | lead(i[5])
2356}
2357///////////////////////////////////////////////////////////////////////////////
2358//die folgende proc war fuer groebner mit fglm vorgesehen, ist aber zu teuer.
2359//Um die projektive Dimension korrekt zu berechnen, muss man aber teuer
2360//voerher ein SB bzgl. einer Gradordnung berechnen und dann homogenisieren.
2361//Sonst koennen hoeherdimensionale Komponenten in Unendlich entstehen
2362
2363proc projInvariants(ideal i,list #)
2364"SYNTAX: @code{projInvariants (} ideal_expression @code{)} @*
2365         @code{projInvariants (} ideal_expression@code{,} list of string_expres          sions@code{)}
2366TYPE:    list, say L, with L[1] and L[2] of type int and L[3] of type intvec
2367PURPOSE: Computes the (projective) dimension (L[1]), degree (L[2]) and the
2368         first Hilbert series (L[3], as intvec) of the homogenized ideal
2369         in the ring given by the procedure @code{hilbRing} with global
2370         ordering dp (resp. wp if the variables have weights >1)
2371         If an argument of type string @code{\"std\"} resp. @code{\"slimgb\"}
2372         is given, the standard basis computatuion uses @code{std} or
2373         @code{slimgb}, otherwise a heuristically chosen method (default)
2374NOTE:    Homogenized means weighted homogenized with respect to the weights
2375         w[i] of the variables var(i) of the basering. The returned dimension,
2376         degree and Hilbertseries are the respective invariants of the
2377         projective variety defined by the homogenized ideal. The dimension
2378         is equal to the (affine) dimension of the ideal in the basering
2379         (degree and Hilbert series make only sense for homogeneous ideals).
2380SEE ALSO: dim, dmult, hilb
2381KEYWORDS: dimension, degree, Hilbert function
2382EXAMPLE: example projInvariants;  shows an example"
2383{
2384  def P = basering;
2385  int p_opt;
2386  string s_opt = option();
2387  if (find(option(), "prot"))  { p_opt = 1; }
2388
2389//---------------- check method and clear denomintors --------------------
2390  int k;
2391  string method;
2392  for (k=1; k<=size(#); k++)
2393  {
2394     if (typeof(#[k]) == "string")
2395     {
2396       method = method + "," + #[k];
2397     }
2398  }
2399
2400  if (npars(P) > 0)             //clear denominators of parameters
2401  {
2402    for( k=ncols(i); k>0; k-- )
2403    {
2404       i[k]=cleardenom(i[k]);
2405    }
2406  }
2407
2408//------------------------ change to hilbRing ----------------------------
2409     list hiRi = hilbRing(i);
2410     intvec W = hiRi[2];
2411     def Philb = hiRi[1];      //note: Philb is no qring and the predefined
2412     setring Philb;            //ideal Id(1) in Philb is homogeneous
2413     int di, de;               //for dimension, degree
2414     intvec hi;                //for hilbert series
2415
2416//-------- compute Hilbert function of homogenized ideal in Philb ---------
2417//Philb has only 1 block. There are three cases
2418
2419     string algorithm;       //possibilities: std, slimgb, stdorslimgb
2420     //define algorithm:
2421     if( find(method,"std") && !find(method,"slimgb") )
2422     {
2423        algorithm = "std";
2424     }
2425     if( find(method,"slimgb") && !find(method,"std") )
2426     {
2427         algorithm = "slimgb";
2428     }
2429     if( find(method,"std") && find(method,"slimgb") ||
2430         (!find(method,"std") && !find(method,"slimgb")) )
2431     {
2432         algorithm = "stdorslimgb";
2433     }
2434
2435     if ( algorithm=="std" || ( algorithm=="stdorslimgb" && char(P)>0 ) )
2436     {
2437        if (p_opt) {"std in ring " + string(Philb);}
2438        Id(1) = std(Id(1));
2439        di = dim(Id(1))-1;
2440        de = mult(Id(1));
2441        hi = hilb( Id(1),1,W );
2442     }
2443     if ( algorithm=="slimgb" || ( algorithm=="stdorslimgb" && char(P)==0 ) )
2444     {
2445        if (p_opt) {"slimgb in ring " + string(Philb);}
2446        Id(1) = slimgb(Id(1));
2447        di = dim( Id(1) );
2448        if (di > -1)
2449        {
2450           di = di-1;
2451        }
2452        de = mult( Id(1) );
2453        hi = hilb( Id(1),1,W );
2454     }
2455     kill Philb;
2456     list L = di,de,hi;
2457     return(L);
2458}
2459example
2460{ "EXAMPLE:"; echo = 2;
2461   ring r = 32003,(x,y,z),lp;
2462   ideal i = y2-xz,x2-z;
2463   projInvariants(i);
2464
2465   ring R = (0),(x,y,z,u,v),lp;
2466   //minpoly = x2+1;
2467   ideal i = x2+1,x2+y+z+u+v,xyzuv-1;
2468   projInvariants(i);
2469   qring S =std(x2+1);
2470   ideal i = imap(R,i);
2471   projInvariants(i);
2472}
2473
2474*/
2475///////////////////////////////////////////////////////////////////////////////
2476//                           EXAMPLES
2477///////////////////////////////////////////////////////////////////////////////
2478/*
2479example stdfglm;
2480example stdhilb;
2481example groebner;
2482example res;
2483example sprintf;
2484example fprintf;
2485example printf;
2486example weightKB;
2487example qslimgb;
2488example par2varRing;
2489*/
2490static proc mod_init()
2491{
2492  //int pagelength=24;
2493  //exportto(Top,pagelength);
2494}
Note: See TracBrowser for help on using the repository browser.