source: git/Singular/LIB/standard.lib @ c99fd4

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