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

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