Changeset 6cd6983 in git for Singular/LIB/bfct.lib


Ignore:
Timestamp:
Dec 1, 2008, 9:58:20 PM (15 years ago)
Author:
Viktor Levandovskyy <levandov@…>
Branches:
(u'spielwiese', '8e0ad00ce244dfd0756200662572aef8402f13d5')
Children:
a6a5d9ad89ba42bbc129bacba382bb6f463aae65
Parents:
8c66fb0e0dc173c2eb5aff1e988308a8527bc045
Message:
*levandov: improvements, bugfixes, docu


git-svn-id: file:///usr/local/Singular/svn/trunk@11212 2c84dea3-7e68-4137-9b89-c4e89433aadc
File:
1 edited

Legend:

Unmodified
Added
Removed
  • Singular/LIB/bfct.lib

    r8c66fb r6cd6983  
    11//////////////////////////////////////////////////////////////////////////////
    2 version="$Id: bfct.lib,v 1.5 2008-10-06 17:04:26 Singular Exp $";
     2version="$Id: bfct.lib,v 1.6 2008-12-01 20:58:20 levandov Exp $";
    33category="Noncommutative";
    44info="
     
    1818MAIN PROCEDURES:
    1919
    20 bfct(f[,s,t,v]);           compute the global Bernstein-Sato polynomial of a given poly
    21 bfctsyz(f[,r,s,t,u,v]);    compute the global Bernstein-Sato polynomial of a given poly
    22 bfctonestep(f[,s,t]);      compute the global Bernstein-Sato polynomial of a given poly
    23 bfctideal(I,w[,s,t]);      compute the global b-function of a given ideal w.r.t. a given weight
    24 minpol(f,I);               compute the minimal polynomial of the endormorphism in basering modulo ideal given by a poly
    25 minpolsyz(f,I[,p,s,t]);    compute the minimal polynomial of the endormorphism in basering modulo ideal given by a poly
    26 linreduce(f,I[,s]);        reduce a poly by linear reductions of its leading term
    27 ncsolve(I[,s]);            find and compute a linear dependency of the elements of an ideal
     20bfct(f[,s,t,v]);            compute the global Bernstein-Sato polynomial of a given poly
     21bfctsyz(f[,r,s,t,u,v]);     compute the global Bernstein-Sato polynomial of a given poly
     22bfctann(f[,s]);             compute the global Bernstein-Sato polynomial of a given poly
     23bfctonestep(f[,s,t]);       compute the global Bernstein-Sato polynomial of a given poly
     24bfctideal(I,w[,s,t]);       compute the global b-function of a given ideal w.r.t. a given weight
     25pintersect(f,I);            compute the intersection of the ideal generated by a given poly with a given ideal
     26pintersectsyz(f,I[,p,s,t]); compute the intersection of the ideal generated by a given poly with a given ideal
     27linreduce(f,I[,s]);         reduce a poly by linear reductions only
     28ncsolve(I[,s]);             find and compute a linear dependency of the elements of an ideal
    2829
    2930AUXILIARY PROCEDURES:
     
    3233isin(l,i);       check whether an element is a member of a list
    3334scalarprod(v,w); compute the standard scalar product of two intvecs
     35vec2poly(v[,i]); convert a coefficient vector to a poly
    3436
    3537SEE ALSO: dmod_lib, dmodapp_lib, gmssing_lib
     
    3840
    3941LIB "qhmoduli.lib"; // for Max
    40 LIB "dmodapp.lib"; // for initialideal etc
     42LIB "dmodapp.lib";  // for initialideal etc
    4143
    4244
     
    4850  example isin;
    4951  example scalarprod;
     52  example vec2poly;
    5053  "MAIN PROCEDURES:";
    5154  example bfct;
    5255  example bfctsyz;
     56  example bfctann;
    5357  example bfctonestep;
    5458  example bfctideal;
    55   example minpol;
    56   example minpolsyz;
     59  example pintersect;
     60  example pintersectsyz;
    5761  example linreduce;
    5862  example ncsolve;
     
    156160proc isin (list l, i)
    157161"USAGE:  isin(l,i);  l a list, i an argument of any type
    158 RETURN:  1 if i is a member of l, or 0 otherwise
     162RETURN:  an int, the position of the first appearance of i in l, or 0 if i is not a member of l
    159163PURPOSE: check whether the second argument is a member of a list
    160164EXAMPLE: example isin; shows an example
     
    166170    if (l[j]==i)
    167171    {
    168       return(1);
     172      return(j);
    169173      break;
    170174    }
     
    216220
    217221proc linreduce(poly f, ideal I, list #)
    218 "USAGE:  linreduce(f, I [,s]);  f a poly, I an ideal, s an optional int
    219 RETURN:  a poly obtained by linear reductions of the leading term of the given poly with an ideal
    220 PURPOSE: reduce a poly only by linear reductions of its leading term
    221 NOTE:    If s<>0, a list consisting of the reduced poly and the vector of the used
     222"USAGE:  linreduce(f, I [,s,t]);  f a poly, I an ideal, s,t optional ints
     223RETURN:  a poly obtained by linear reductions of the given poly with the given ideal
     224PURPOSE: reduce a poly only by linear reductions
     225NOTE:    If s<>0, a list consisting of the reduced poly and the coefficient vector of the used
    222226@*       reductions is returned.
     227@*       If t<>0, only leading monomials are reduced, otherwise, and by default, all monomials
     228@*       are reduced, if possible.
    223229EXAMPLE: example linreduce; shows examples
    224230"
    225231{
     232  int ppl = printlevel - voice + 2;
    226233  int remembercoeffs = 0; // default
     234  int redlm          = 0; // default
    227235  if (size(#)>0)
    228236  {
     
    231239      remembercoeffs = #[1];
    232240    }
    233   }
    234   int i;
     241    if (size(#)>1)
     242    {
     243      if (typeof(#[2])=="int" || typeof(#[2])=="number")
     244      {
     245        redlm = #[2];
     246      }
     247    }
     248  }
     249  int i,j,k;
    235250  int sI = ncols(I);
    236251  ideal lmI,lcI;
     
    241256  }
    242257  vector v;
    243   poly lm,c;
    244   int reduction;
    245   lm = leadmonom(f);
    246   reduction = 1;
    247   while (reduction == 1) // while there was a reduction
    248   {
    249     reduction = 0;
    250     for (i=sI;i>=1;i--)
    251     {
    252       if (lm <> 0 && lm == lmI[i])
    253       {
    254         c = leadcoef(f)/lcI[i];
    255         f = f - c*I[i];
    256         lm = leadmonom(f);
    257         reduction = 1;
    258         if (remembercoeffs <> 0)
     258  poly c;
     259  int reduction = 1;
     260  if (redlm == 0)
     261  {
     262    ideal monomf;
     263    for (k=1; k<=size(f); k++)
     264    {
     265      monomf[k] = normalize(f[k]);
     266    }
     267    while (reduction == 1) // while there was a reduction
     268    {
     269      reduction = 0;
     270      for (i=sI; i>=1; i--)
     271      {
     272        dbprint(ppl,"testing ideal entry:",i);
     273        for (j=1; j<=size(f); j++)
    259274        {
    260           v = v - c * gen(i);
     275          if (monomf[j] == lmI[i])
     276          {
     277            c = leadcoef(f[j])/lcI[i];
     278            f = f - c*I[i];
     279            dbprint(ppl,"reducing poly to ",f);
     280            monomf = 0;
     281            for (k=1; k<=size(f); k++)
     282            {
     283              monomf[k] = normalize(f[k]);
     284            }
     285            reduction = 1;
     286            if (remembercoeffs <> 0)
     287            {
     288              v = v - c * gen(i);
     289            }
     290            break;
     291          }
     292        }
     293        if (reduction == 1)
     294        {
     295          break;
     296        }
     297      }
     298    }
     299  }
     300  else // reduce only leading monomials
     301  {
     302    poly lm = leadmonom(f);
     303    while (reduction == 1) // while there was a reduction
     304    {
     305      reduction = 0;
     306      for (i=sI;i>=1;i--)
     307      {
     308        if (lm <> 0 && lm == lmI[i])
     309        {
     310          c = leadcoef(f)/lcI[i];
     311          f = f - c*I[i];
     312          lm = leadmonom(f);
     313          reduction = 1;
     314          if (remembercoeffs <> 0)
     315          {
     316            v = v - c * gen(i);
     317          }
    261318        }
    262319      }
     
    279336  ideal I = 1,y,xy;
    280337  poly f = 5xy+7y+3;
    281   poly g = 5y+7x+3;
    282   linreduce(f,I);
     338  poly g = 7x+5y+3;
    283339  linreduce(g,I);
     340  linreduce(g,I,0,1);
    284341  linreduce(f,I,1);
    285342}
     
    386443}
    387444
    388 proc minpol (poly s, ideal I)
    389 "USAGE:  minpol(f, I);  f a poly, I an ideal
    390 RETURN:  coefficient vector of the minimal polynomial of the endomorphism of basering modulo I defined by f
    391 PURPOSE: compute the minimal polynomial
    392 NOTE:    If f does not define an endomorphism, this proc will not terminate.
     445proc pintersect (poly s, ideal I)
     446"USAGE:  pintersect(f, I);  f a poly, I an ideal
     447RETURN:  coefficient vector of the monic generator of the intersection of the ideal generated by f with I
     448PURPOSE: compute the intersection of an ideal with a principal ideal defined by f
     449NOTE:    If the intersection is zero, this proc might not terminate.
    393450@*       I should be given as standard basis.
    394451@*       If printlevel=1, progress debug messages will be printed,
    395452@*       if printlevel>=2, all the debug messages will be printed.
    396 EXAMPLE: example minpol; shows examples
     453EXAMPLE: example pintersect; shows examples
    397454"
    398455{
     
    400457  attrib(I,"isSB",1); // set attribute for suppressing NF messages
    401458  int ppl = printlevel-voice+2;
     459  // ---case 1: I = basering---
     460  if (size(I) == 1)
     461  {
     462    if (simplify(I,1) == ideal(1))
     463    {
     464      return(gen(2)); // = s
     465    }
     466  }
    402467  def save = basering;
     468  int n = nvars(save);
    403469  int i,j,k;
     470  // ---case 2: intersection is zero---
     471  intvec degs = leadexp(s);
     472  intvec possdegbounds;
     473  list degI;
     474  i = 1;
     475  while (i <= ncols(I))
     476  {
     477    if (i == ncols(I)+1)
     478    {
     479      break;
     480    }
     481    degI[i] = leadexp(I[i]);
     482    for (j=1; j<=n; j++)
     483    {
     484      if (degs[j] == 0)
     485      {
     486        if (degI[i][j] <> 0)
     487        {
     488          break;
     489        }
     490      }
     491      if (j == n)
     492      {
     493        k++;
     494        possdegbounds[k] = Max(degI[i]);
     495      }
     496    }
     497    i++;
     498  }
     499  int degbound = Min(possdegbounds);
     500  dbprint(ppl,"a lower bound for the degree of the insection is:");
     501  dbprint(ppl,degbound);
     502  if (degbound == 0) // lm(s) does not appear in lm(I)
     503  {
     504    return(vector(0));
     505  }
     506  // ---case 3: intersection is non-trivial---
     507  ideal redNI = 1;
    404508  vector v;
    405509  list l,ll;
    406510  l[1] = vector(0);
    407   poly toNF, tobracket, newNF, rednewNF;
    408   ideal NI = 1;
     511  poly toNF,tobracket,newNF,rednewNF,oldNF,secNF;
    409512  i = 1;
    410   ideal redNI = 1;
    411   newNF = NF(s,I);
    412   dbprint(ppl+1,"minpol starts...");
    413   dbprint(ppl+1,"with ideal I=", I);
     513  dbprint(ppl+1,"pintersect starts...");
    414514  while (1)
    415515  {
     
    417517    if (i>1)
    418518    {
    419       tobracket = s^(i-1)-NI[i];
    420       if (tobracket==0) // bracket doesn't like zeros
     519      oldNF = newNF;
     520      tobracket = s^(i-1) - oldNF;
     521      if (tobracket==0) // todo bug in bracket?
    421522      {
    422523        toNF = 0;
     
    424525      else
    425526      {
    426         toNF = bracket(tobracket,NI[2]);
    427       }
    428       newNF = NF(toNF+NI[i]*NI[2],I);  // = NF(s^i,I)
     527        toNF = bracket(tobracket,secNF);
     528      }
     529      newNF = NF(toNF+oldNF*secNF,I);  // = NF(s^i,I)
     530    }
     531    else
     532    {
     533      newNF = NF(s,I);
     534      secNF = newNF;
    429535    }
    430536    ll = linreduce(newNF,redNI,1);
    431537    rednewNF = ll[1];
    432538    l[i+1] = ll[2];
    433     dbprint(ppl+1,"newNF is:", newNF);
    434     dbprint(ppl+1,"rednewNF is:", rednewNF);
     539    dbprint(ppl,"newNF is:", newNF);
     540    dbprint(ppl,"rednewNF is:", rednewNF);
    435541    if (rednewNF != 0) // no linear dependency
    436542    {
    437       NI[i+1] = newNF;
    438543      redNI[i+1] = rednewNF;
    439       dbprint(ppl+1,"NI is:", NI);
    440       dbprint(ppl+1,"redNI is:", redNI);
    441544      i++;
    442545    }
    443546    else // there is a linear dependency, hence we are done
    444547    {
    445       dbprint(ppl+1,"the degree of the minimal polynomial is:", i);
     548      dbprint(ppl+1,"the degree of the generator of the intersection is:", i);
    446549      break;
    447550    }
    448551  }
    449552  dbprint(ppl,"used linear reductions:", l);
    450   // we obtain the coefficients of the minimal polynomial by the used reductions:
     553  // we obtain the coefficients of the generator of the intersection by the used reductions:
    451554  ring @R = 0,(a(1..i+1)),dp;
    452555  setring @R;
     
    461564    }
    462565  }
    463   for(j=i;j>=1;j--)
     566  for (j=i;j>=1;j--)
    464567  {
    465568    C[i+1]= subst(C[i+1],a(j),a(j)+C[j]);
     
    474577  v = imap(@R,v);
    475578  kill @R;
    476   dbprint(ppl+1,"minpol finished");
     579  dbprint(ppl+1,"pintersect finished");
    477580  return(v);
    478581}
     
    480583{
    481584  "EXAMPLE:"; echo = 2;
    482   printlevel = 0;
    483585  ring r = 0,(x,y),dp;
    484586  poly f = x^2+y^3+x*y^2;
     
    486588  setring D;
    487589  inF;
    488   poly s = t*Dt;
    489   minpol(s,inF);
    490 }
    491 
    492 proc minpolsyz (poly s, ideal II, list #)
    493 "USAGE:  minpolsyz(f, I [,p,s,t]);  f a poly, I an ideal, p, t optial ints, p a prime number
    494 RETURN:  coefficient vector of the minimal polynomial of the endomorphism of basering modulo I defined by f
    495 PURPOSE: compute the minimal polynomial
    496 NOTE:    If f does not define an endomorphism, this proc will not terminate.
     590  pintersect(t*Dt,inF);
     591}
     592
     593proc pintersectsyz (poly s, ideal II, list #)
     594"USAGE:  pintersectsyz(f, I [,p,s,t]);  f a poly, I an ideal, p, t optial ints, p a prime number
     595RETURN:  coefficient vector of the monic generator of the intersection of the ideal generated by f with I
     596PURPOSE: compute the intersection of an ideal with a principal ideal defined by f
     597NOTE:    If the intersection is zero, this proc might not terminate.
    497598@*       I should be given as standard basis.
    498 @*       If p>0 is given, the proc computes the minimal polynomial in char p first and
    499 @*       then only searches for a minimal polynomial of the obtained degree in the basering.
     599@*       If p>0 is given, this proc computes the generator of the intersection in char p first and
     600@*       then only searches for a generator of the obtained degree in the basering.
    500601@*       Otherwise, it searched for all degrees.
    501602@*       This is done by computing syzygies.
     
    506607@*       If printlevel=1, progress debug messages will be printed,
    507608@*       if printlevel>=2, all the debug messages will be printed.
    508 EXAMPLE: example minpolsyz; shows examples
     609EXAMPLE: example pintersectsyz; shows examples
    509610"
    510611{
     
    566667  }
    567668  i = 1;
    568   dbprint(ppl+1,"minpolynomial starts...");
     669  dbprint(ppl+1,"pintersectsyz starts...");
    569670  dbprint(ppl+1,"with ideal I=", I);
    570671  while (1)
     
    639740    i++;
    640741  }
    641   dbprint(ppl+1,"minpol finished");
     742  dbprint(ppl+1,"pintersectsyz finished");
    642743  return(v);
    643744}
     
    645746{
    646747  "EXAMPLE:"; echo = 2;
    647   printlevel = 0;
    648748  ring r = 0,(x,y),dp;
    649749  poly f = x^2+y^3+x*y^2;
     
    652752  inF;
    653753  poly s = t*Dt;
    654   minpolsyz(s,inF);
     754  pintersectsyz(s,inF);
    655755  int p = prime(20000);
    656   minpolsyz(s,inF,p,0,0);
     756  pintersectsyz(s,inF,p,0,0);
     757}
     758
     759proc vec2poly (list #)
     760"USAGE:  vec2poly(v [,i]);  v a vector or an intvec, i an optional int
     761RETURN:  a poly with coefficient vector v
     762PURPOSE: convert a coefficient vector to a poly
     763NOTE:    If i>0 is given, the returned poly is an element of K[var(i)],
     764@*       otherwise, and by default, @code{i=1} is used.
     765@*       The first entry of v is the coefficient of 1.
     766EXAMPLE: example vec2poly; shows examples
     767"
     768{
     769  def save = basering;
     770  int i,ringvar;
     771  ringvar = 1; // default
     772  if (size(#) > 0)
     773  {
     774    if (typeof(#[1])=="vector" || typeof(#[1])=="intvec")
     775    {
     776      def v = #[1];
     777    }
     778    else
     779    {
     780      ERROR("wrong input: expected vector/intvec expression");
     781    }
     782    if (size(#) > 1)
     783    {
     784      if (typeof(#[2])=="int" || typeof(#[2])=="number")
     785      {
     786        ringvar = int(#[2]);
     787      }
     788    }
     789  }
     790  if (ringvar > nvars(save))
     791  {
     792    ERROR("var out of range");
     793  }
     794  poly p;
     795  for (i=1; i<=nrows(v); i++)
     796  {
     797    p = p + v[i]*(var(ringvar))^(i-1);
     798  }
     799  return(p);
     800}
     801example
     802{
     803  "EXAMPLE:"; echo = 2;
     804  ring r = 0,(x,y),dp;
     805  vector v = gen(1) + 3*gen(3) + 22/9*gen(4);
     806  intvec iv = 3,2,1;
     807  vec2poly(v,2);
     808  vec2poly(iv);
    657809}
    658810
     
    705857}
    706858
    707 static proc bfctengine (poly f, int whichengine, int methodord, int methodminpol, int minpolchar, int modengine, intvec u0)
     859static proc bfctengine (poly f, int inorann, int whichengine, int methodord, int methodpintersect, int pintersectchar, int modengine, intvec u0)
    708860{
    709861  int ppl = printlevel - voice +2;
     
    711863  def save = basering;
    712864  int n = nvars(save);
    713   def DD = initialmalgrange(f,whichengine,methodord,1,u0);
    714   setring DD;
    715   ideal inI = inF;
    716   kill inF;
    717   poly s = t*Dt;
     865  if (inorann == 0) // bfct using initial ideal
     866  {
     867    def D = initialmalgrange(f,whichengine,methodord,1,u0);
     868    setring D;
     869    ideal J = inF;
     870    kill inF;
     871    poly s = t*Dt;
     872  }
     873  else // bfct using Ann(f^s)
     874  {
     875    def D = SannfsBFCT(f,whichengine);
     876    setring D;
     877    ideal J = LD;
     878    kill LD;
     879  }
    718880  vector b;
    719881  // try it modular
    720   if (methodminpol <> 0) // minpolsyz
    721   {
    722     if (minpolchar == 0) // minpolsyz::modular
     882  if (methodpintersect <> 0) // pintersectsyz
     883  {
     884    if (pintersectchar == 0) // pintersectsyz::modular
    723885    {
    724886      int lb = 30000;
     
    734896          usedprimes = usedprimes,q;
    735897          dbprint(ppl,"used prime is: "+string(q));
    736           b = minpolsyz(s,inI,q,whichengine,modengine);
     898          b = pintersectsyz(s,J,q,whichengine,modengine);
    737899        }
    738900        i++;
    739901      }
    740902    }
    741     else // minpolsyz::non-modular
    742     {
    743       b = minpolsyz(s,inI,0,whichengine);
    744     }
    745   }
    746   else // minpol: linreduce
    747   {
    748     b = minpol(s,inI);
     903    else // pintersectsyz::non-modular
     904    {
     905      b = pintersectsyz(s,J,0,whichengine);
     906    }
     907  }
     908  else // pintersect: linreduce
     909  {
     910    b = pintersect(s,J);
    749911  }
    750912  setring save;
    751   vector b = imap(DD,b);
    752   list l = listofroots(b,1);
     913  vector b = imap(D,b);
     914  if (inorann == 0)
     915  {
     916    list l = listofroots(b,1);
     917  }
     918  else
     919  {
     920    list l = listofroots(b,0);
     921  }
    753922  return(l);
    754923}
     
    758927RETURN:  list of roots of the Bernstein-Sato polynomial bs(f) and their multiplicies
    759928PURPOSE: compute the global Bernstein-Sato polynomial for a hypersurface, defined by f, according to the algorithm by Masayuki Noro
    760 NOTE:    In this proc, a system of linear equations is solved by linear reductions.
     929NOTE:    In this proc, the initial Malgrange ideal is computed.
     930@*       Further, a system of linear equations is solved by linear reductions.
    761931@*       If s<>0, @code{std} is used for Groebner basis computations,
    762932@*       otherwise, and by default, @code{slimgb} is used.
     
    801971    }
    802972  }
    803   list b = bfctengine(f,whichengine,methodord,0,0,0,u0);
     973  list b = bfctengine(f,0,whichengine,methodord,0,0,0,u0);
    804974  return(b);
    805975}
     
    818988RETURN:  list of roots of the Bernstein-Sato polynomial bs(f) and its multiplicies
    819989PURPOSE: compute the global Bernstein-Sato polynomial for a hypersurface, defined by f, according to the algorithm by Masayuki Noro
    820 NOTE:    In this proc, a system of linear equations is solved by computing syzygies.
     990NOTE:    In this proc, the initial Malgrange ideal is computed.
     991@*       Further, a system of linear equations is solved by computing syzygies.
    821992@*       If r<>0, @code{std} is used for Groebner basis computations in characteristic 0,
    822993@*       otherwise, and by default, @code{slimgb} is used.
    823994@*       If s<>0, a matrix ordering is used for Groebner basis computations,
    824995@*       otherwise, and by default, a block ordering is used.
    825 @*       If t<>0, the minimal polynomial computation is solely performed over charasteristic 0,
     996@*       If t<>0, the computation of the intersection is solely performed over charasteristic 0,
    826997@*       otherwise and by default, a modular method is used.
    827998@*       If u<>0 and by default, @code{std} is used for Groebner basis computations in characteristic >0,
     
    8391010  // one for the engine used for Groebner basis computations in char 0,
    8401011  // one for  M() ordering or its realization
    841   // one for a modular method when computing the minimal polynomial
     1012  // one for a modular method when computing the intersection
    8421013  // and one for the engine used for Groebner basis computations in char >0
    8431014  // in # can also be the optional weight vector
     
    8461017  int whichengine = 0; // default
    8471018  int methodord   = 0; // default
    848   int minpolchar  = 0; // default
     1019  int pintersectchar  = 0; // default
    8491020  int modengine   = 1; // default
    8501021  intvec u0       = 0; // default
     
    8651036        if (typeof(#[3])=="int" || typeof(#[3])=="number")
    8661037        {
    867           minpolchar = int(#[3]);
     1038          pintersectchar = int(#[3]);
    8681039        }
    8691040        if (size(#)>3)
     
    8841055    }
    8851056  }
    886   list b = bfctengine(f,whichengine,methodord,1,minpolchar,modengine,u0);
     1057  list b = bfctengine(f,0,whichengine,methodord,1,pintersectchar,modengine,u0);
    8871058  return(b);
    8881059}
     
    9381109    s = s + w[i]*var(i)*var(n+i);
    9391110  }
    940   vector b = minpol(s,J);
     1111  vector b = pintersect(s,J);
    9411112  list l = listofroots(b,0);
    9421113  return(l);
     
    9941165  setring DDh;
    9951166  dbprint(ppl, "the initial ideal:", string(matrix(inF)));
    996   inF = nselect(inF,3..2*n+4);
    997   inF = nselect(inF,1);
     1167  intvec tonselect = 1;
     1168  for (i=3; i<=2*n+4; i++)
     1169  {
     1170    tonselect = tonselect,i;
     1171  }
     1172  inF = nselect(inF,tonselect);
    9981173  dbprint(ppl, "generators containing only s:", string(matrix(inF)));
    9991174  inF = engine(inF, whichengine); // is now a principal ideal;
     
    10131188  bfctonestep(f);
    10141189  bfctonestep(f,1,1);
     1190}
     1191
     1192proc bfctann (poly f, list #)
     1193"USAGE:  bfctann(f [,r]);  f a poly, r an optional int
     1194RETURN:  list of roots of the Bernstein-Sato polynomial bs(f) and their multiplicies
     1195PURPOSE: compute the global Bernstein-Sato polynomial for a hypersurface, defined by f
     1196NOTE:    In this proc, ann(f^s) is computed.
     1197@*       If r<>0, @code{std} is used for Groebner basis computations,
     1198@*       otherwise, and by default, @code{slimgb} is used.
     1199@*       If printlevel=1, progress debug messages will be printed,
     1200@*       if printlevel>=2, all the debug messages will be printed.
     1201EXAMPLE: example bfctann; shows examples
     1202"
     1203{
     1204  def save = basering;
     1205  int ppl = printlevel - voice + 2;
     1206  int whichengine = 0; // default
     1207  if (size(#)>0)
     1208  {
     1209    if (typeof(#[1])=="int" || typeof(#[1])=="number")
     1210    {
     1211      whichengine = int(#[1]);
     1212    }
     1213  }
     1214  list b = bfctengine(f,1,whichengine,0,1,0,0,0);
     1215  return(b);
     1216}
     1217example
     1218{
     1219  "EXAMPLE:"; echo = 2;
     1220  ring r = 0,(x,y),dp;
     1221  poly f = x^2+y^3+x*y^2;
     1222  bfctann(f);
    10151223}
    10161224
     
    10371245  bfct(xyzreiffen45);
    10381246}
     1247
Note: See TracChangeset for help on using the changeset viewer.