source: git/Singular/LIB/standard.lib @ 334c21f

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