Changeset 7f3ad4 in git for Singular/LIB/bfct.lib


Ignore:
Timestamp:
Jan 14, 2009, 5:07:05 PM (15 years ago)
Author:
Hans Schönemann <hannes@…>
Branches:
(u'spielwiese', '5b153614cbc72bfa198d75b1e9e33dab2645d9fe')
Children:
0721816437af5ddabc83aa203a12d9b58b42a33c
Parents:
95edd5641377e851d4a1d4e986853687991d0895
Message:
*hannes: format


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

Legend:

Unmodified
Added
Removed
  • Singular/LIB/bfct.lib

    r95edd5 r7f3ad4  
    11//////////////////////////////////////////////////////////////////////////////
    2 version="$Id: bfct.lib,v 1.10 2009-01-12 19:31:23 Singular Exp $";
     2version="$Id: bfct.lib,v 1.11 2009-01-14 16:07:03 Singular Exp $";
    33category="Noncommutative";
    44info="
     
    3535AUXILIARY PROCEDURES:
    3636
    37 allPositive(v);  checks whether all entries of an intvec are positive 
     37allPositive(v);  checks whether all entries of an intvec are positive
    3838scalarProd(v,w); computes the standard scalar product of two intvecs
    3939vec2poly(v[,i]); constructs an univariate poly with given coefficients
     
    7979EXAMPLE: example gradedWeyl; shows examples
    8080NOTE:    u[i] is the weight of x(i), v[i] the weight of D(i).
    81 @*       u+v has to be a non-negative intvec. 
     81@*       u+v has to be a non-negative intvec.
    8282"
    8383{
     
    8888  {
    8989    ERROR("weight vectors have wrong dimension");
    90   } 
     90  }
    9191  intvec uv,gr;
    9292  uv = u+v;
     
    9797      if (uv[i]==0)
    9898      {
    99         gr[i] = 0;
     99        gr[i] = 0;
    100100      }
    101101      else
    102102      {
    103         gr[i] = 1;
     103        gr[i] = 1;
    104104      }
    105105    }
     
    114114  for (i=1; i<=n; i++)
    115115  {
    116     if (gr[i] == 1) 
     116    if (gr[i] == 1)
    117117    {
    118118       l2[n+i] = "xi("+string(i)+")";
     
    231231RETURN:  ideal/list, linear reductum (over field) of f by elements from I
    232232PURPOSE: reduce a poly only by linear reductions (no monomial multiplications)
    233 NOTE:    If s<>0, a list consisting of the reduced ideal and the coefficient 
     233NOTE:    If s<>0, a list consisting of the reduced ideal and the coefficient
    234234@*       vectors of the used reductions given as module is returned.
    235235@*       Otherwise (and by default), only the reduced ideal is returned.
     
    255255      if (typeof(#[2])=="int" || typeof(#[2])=="number")
    256256      {
    257         redtail = #[2];
     257        redtail = #[2];
    258258      }
    259259    }
     
    273273      for (i=1; i<=sI; i++)
    274274      {
    275         if (I[i] == 0)
    276         {
    277           j++;
    278           J[j] = 0;
    279           ordJ[j] = -1;
    280           M[j] = gen(i);
    281         }
    282         else
    283         {
    284           M[i+sZeros-j] = gen(lJ[2][i-j]+j);
    285         }
     275        if (I[i] == 0)
     276        {
     277          j++;
     278          J[j] = 0;
     279          ordJ[j] = -1;
     280          M[j] = gen(i);
     281        }
     282        else
     283        {
     284          M[i+sZeros-j] = gen(lJ[2][i-j]+j);
     285        }
    286286      }
    287287    }
     
    290290      for (i=1; i<=sZeros; i++)
    291291      {
    292         J[i] = 0;
    293         ordJ[i] = -1;
     292        J[i] = 0;
     293        ordJ[i] = -1;
    294294      }
    295295    }
     
    327327      for (j=sZeros+1; j<i; j++)
    328328      {
    329         if (lm == 0) { break; }
    330         if (ordlm > maxordJ) { break; }
    331         if (ordlm == ordJ[j])
    332         {           
    333           if (lm > maxlmJ) { break; }
    334           if (lm == lmJ[j])
    335           {
    336             dbprint(ppl,"reducing " + string(redpoly));
    337             dbprint(ppl," with " + string(J[j]));
    338             c = leadcoef(redpoly)/leadcoef(J[j]);
    339             redpoly = redpoly - c*J[j];
    340             dbprint(ppl," to " + string(redpoly));
    341             lm = leadmonom(redpoly);
    342             ordlm = ord(lm);
    343             if (remembercoeffs <> 0) { M[i] = M[i] - c * M[j]; }
    344             reduction = 1;
    345           }
    346         }
     329        if (lm == 0) { break; }
     330        if (ordlm > maxordJ) { break; }
     331        if (ordlm == ordJ[j])
     332        {
     333          if (lm > maxlmJ) { break; }
     334          if (lm == lmJ[j])
     335          {
     336            dbprint(ppl,"reducing " + string(redpoly));
     337            dbprint(ppl," with " + string(J[j]));
     338            c = leadcoef(redpoly)/leadcoef(J[j]);
     339            redpoly = redpoly - c*J[j];
     340            dbprint(ppl," to " + string(redpoly));
     341            lm = leadmonom(redpoly);
     342            ordlm = ord(lm);
     343            if (remembercoeffs <> 0) { M[i] = M[i] - c * M[j]; }
     344            reduction = 1;
     345          }
     346        }
    347347      }
    348348    }
     
    354354    lmJ = insertGenerator(lmJ,lm,j);
    355355    ordJ = insertGenerator(ordJ,poly(ordlm),j);
    356     if (remembercoeffs <> 0) 
     356    if (remembercoeffs <> 0)
    357357    {
    358358      v = M[i];
     
    373373      for (j=i+1; j<=sI; j++)
    374374      {
    375         for (k=2; k<=size(J[j]); k++) // run over all terms in J[j]
    376         {
    377           if (ordJ[i] == ord(J[j][k]))
    378           {   
    379             if (lm == normalize(J[j][k]))
    380             {
    381               c = leadcoef(J[j][k])/leadcoef(J[i]);
    382               dbprint(ppl,"reducing " + string(J[j]));
    383               dbprint(ppl," with " + string(J[i]));
    384               J[j] = J[j] - c*J[i];
    385               dbprint(ppl," to " + string(J[j]));
    386               if (remembercoeffs <> 0) { M[j] = M[j] - c * M[i]; }
    387             }
    388           }
    389         }
     375        for (k=2; k<=size(J[j]); k++) // run over all terms in J[j]
     376        {
     377          if (ordJ[i] == ord(J[j][k]))
     378          {
     379            if (lm == normalize(J[j][k]))
     380            {
     381              c = leadcoef(J[j][k])/leadcoef(J[i]);
     382              dbprint(ppl,"reducing " + string(J[j]));
     383              dbprint(ppl," with " + string(J[i]));
     384              J[j] = J[j] - c*J[i];
     385              dbprint(ppl," to " + string(J[j]));
     386              if (remembercoeffs <> 0) { M[j] = M[j] - c * M[i]; }
     387            }
     388          }
     389        }
    390390      }
    391391    }
     
    410410RETURN:  poly/list, linear reductum (over field) of f by elements from I
    411411PURPOSE: reduce a poly only by linear reductions (no monomial multiplications)
    412 NOTE:    If s<>0, a list consisting of the reduced poly and the coefficient 
     412NOTE:    If s<>0, a list consisting of the reduced poly and the coefficient
    413413@*       vector of the used reductions is returned, otherwise (and by default)
    414414@*       only reduced poly is returned.
     
    433433      if (typeof(#[2])=="int" || typeof(#[2])=="number")
    434434      {
    435         redtail = #[2];
     435        redtail = #[2];
    436436      }
    437437      if (size(#)>2)
    438438      {
    439         if (typeof(#[3])=="int" || typeof(#[3])=="number")
     439        if (typeof(#[3])=="int" || typeof(#[3])=="number")
    440440        {
    441           prepareideal = #[3];
    442         }
     441          prepareideal = #[3];
     442        }
    443443      }
    444444    }
     
    469469      for (i=1; i<=sI; i++)
    470470      {
    471         M[i] = gen(i);
     471        M[i] = gen(i);
    472472      }
    473473    }
     
    495495      if (ordf == ordI[i])
    496496      {
    497         if (lm == lmI[i])
    498         {
    499           c = leadcoef(f)/lcI[i];
    500           f = f - c*I[i];
    501           lm = leadmonom(f);
    502           ordf = ord(lm);
    503           if (remembercoeffs <> 0)
    504           {
    505             v = v - c * M[i];
    506           }
    507         }
     497        if (lm == lmI[i])
     498        {
     499          c = leadcoef(f)/lcI[i];
     500          f = f - c*I[i];
     501          lm = leadmonom(f);
     502          ordf = ord(lm);
     503          if (remembercoeffs <> 0)
     504          {
     505            v = v - c * M[i];
     506          }
     507        }
    508508      }
    509509    }
     
    518518      for (j=1; j<=size(f); j++)
    519519      {
    520         if (ord(f[j]) == ordI[i])
    521         {
    522           if (normalize(f[j]) == lmI[i])
    523           {
    524             c = leadcoef(f[j])/lcI[i];
    525             f = f - c*I[i];
    526             dbprint(ppl,"reducing poly to ",f);
    527             if (remembercoeffs <> 0)
    528             {
    529               v = v - c * M[i];
    530             }
    531           }
    532         }
     520        if (ord(f[j]) == ordI[i])
     521        {
     522          if (normalize(f[j]) == lmI[i])
     523          {
     524            c = leadcoef(f[j])/lcI[i];
     525            f = f - c*I[i];
     526            dbprint(ppl,"reducing poly to ",f);
     527            if (remembercoeffs <> 0)
     528            {
     529              v = v - c * M[i];
     530            }
     531          }
     532        }
    533533      }
    534534    }
     
    556556  f = x3 + y2 + x2 + y + x;
    557557  I = x3 - y3, y3 - x2, x3 - y2, x2 - y, y2-x;
    558   list l = linReduce(f, I, 1); 
     558  list l = linReduce(f, I, 1);
    559559  l;
    560560  module M = I;
     
    601601      if (p <> 0)
    602602      {
    603         whichengine = 1;
     603        whichengine = 1;
    604604      }
    605605    }
     
    668668proc pIntersect (poly s, ideal I)
    669669"USAGE:  pIntersect(f, I);  f a poly, I an ideal
    670 RETURN:  vector, coefficient vector of the monic polynomial 
     670RETURN:  vector, coefficient vector of the monic polynomial
    671671PURPOSE: compute the intersection of ideal I with the subalgebra K[f]
    672672ASSUME:  I is given as Groebner basis.
     
    712712      if (degs[j] == 0)
    713713      {
    714         if (degI[i][j] <> 0)
    715         {
    716           break;
    717         }
     714        if (degI[i][j] <> 0)
     715        {
     716          break;
     717        }
    718718      }
    719719      if (j == n)
    720720      {
    721         k++;
    722         possdegbounds[k] = Max(degI[i]);
     721        k++;
     722        possdegbounds[k] = Max(degI[i]);
    723723      }
    724724    }
     
    749749      if (tobracket==0) // todo bug in bracket?
    750750      {
    751         toNF = 0;
     751        toNF = 0;
    752752      }
    753753      else
    754754      {
    755         toNF = bracket(tobracket,secNF);
     755        toNF = bracket(tobracket,secNF);
    756756      }
    757757      newNF = NF(toNF+oldNF*secNF,I);  // = NF(s^i,I)
     
    802802    v = v + m[j,1]*gen(j);
    803803  }
    804   setring save; 
     804  setring save;
    805805  v = imap(@R,v);
    806806  kill @R;
     
    864864      if (size(#)>2)
    865865      {
    866         if (typeof(#[3])=="int" || typeof(#[3])=="number")
    867         {
    868           modengine = int(#[3]);
    869         }
     866        if (typeof(#[3])=="int" || typeof(#[3])=="number")
     867        {
     868          modengine = int(#[3]);
     869        }
    870870      }
    871871    }
     
    909909      if (tobracket!=0)
    910910      {
    911         toNF = bracket(tobracket,NI[2]) + NI[i]*NI[2];
     911        toNF = bracket(tobracket,NI[2]) + NI[i]*NI[2];
    912912      }
    913913      else
    914914      {
    915         toNF = NI[i]*NI[2];
     915        toNF = NI[i]*NI[2];
    916916      }
    917917      newNF =  NF(toNF,I);
     
    927927      if (v!=0) // there is a modular solution
    928928      {
    929         dbprint(ppl,"got solution in char ",solveincharp," of degree " ,i);
    930         setring save;
    931         v = linSyzSolve(NI,whichengine);
    932         if (v==0)
    933         {
    934           break;
    935         }
     929        dbprint(ppl,"got solution in char ",solveincharp," of degree " ,i);
     930        setring save;
     931        v = linSyzSolve(NI,whichengine);
     932        if (v==0)
     933        {
     934          break;
     935        }
    936936      }
    937937      else // no modular solution
    938938      {
    939         setring save;
    940         v = 0;
     939        setring save;
     940        v = 0;
    941941      }
    942942    }
     
    965965      {
    966966        dbprint(ppl,"linSyzSolve: got solution!");
    967         // "got solution!";
     967        // "got solution!";
    968968        break;
    969969      }
     
    10151015      if (typeof(#[2])=="int" || typeof(#[2])=="number")
    10161016      {
    1017         ringvar = int(#[2]);
     1017        ringvar = int(#[2]);
    10181018      }
    10191019    }
     
    11211121      while (b == 0)
    11221122      {
    1123         dbprint(ppl,"number of run in the loop: "+string(i));
    1124         int q = prime(random(lb,ub));
    1125         if (findFirst(usedprimes,q)==0) // if q was not already used
     1123        dbprint(ppl,"number of run in the loop: "+string(i));
     1124        int q = prime(random(lb,ub));
     1125        if (findFirst(usedprimes,q)==0) // if q was not already used
    11261126        {
    1127           usedprimes = usedprimes,q;
    1128           dbprint(ppl,"used prime is: "+string(q));
    1129           b = pIntersectSyz(s,J,q,whichengine,modengine);
    1130         }
    1131         i++;
    1132       }
    1133     }
    1134     else // pIntersectSyz::non-modular 
     1127          usedprimes = usedprimes,q;
     1128          dbprint(ppl,"used prime is: "+string(q));
     1129          b = pIntersectSyz(s,J,q,whichengine,modengine);
     1130        }
     1131        i++;
     1132      }
     1133    }
     1134    else // pIntersectSyz::non-modular
    11351135    {
    11361136      b = pIntersectSyz(s,J,0,whichengine);
     
    11581158RETURN:  list of ideal and intvec
    11591159PURPOSE: computes the roots of the Bernstein-Sato polynomial b(s)
    1160 @*       for the hypersurface defined by f. 
     1160@*       for the hypersurface defined by f.
    11611161ASSUME:  The basering is a commutative polynomial ring in char 0.
    1162 BACKGROUND: In this proc, the initial Malgrange ideal is computed according to the algorithm 
     1162BACKGROUND: In this proc, the initial Malgrange ideal is computed according to the algorithm
    11631163@*       by Masayuki Noro and then a system of linear equations is solved by linear reductions.
    1164 NOTE:    In the output list, the ideal contains all the roots 
     1164NOTE:    In the output list, the ideal contains all the roots
    11651165@*       and the intvec their multiplicities.
    11661166@*       If s<>0, @code{std} is used for GB computations,
    1167 @*       otherwise, and by default, @code{slimgb} is used. 
     1167@*       otherwise, and by default, @code{slimgb} is used.
    11681168@*       If t<>0, a matrix ordering is used for Groebner basis computations,
    11691169@*       otherwise, and by default, a block ordering is used.
    1170 @*       If v is a positive weight vector, v is used for homogenization computations, 
     1170@*       If v is a positive weight vector, v is used for homogenization computations,
    11711171@*       otherwise and by default, no weights are used.
    11721172DISPLAY: If printlevel=1, progress debug messages will be printed,
     
    11991199      if (size(#)>2)
    12001200      {
    1201         if (typeof(#[3])=="intvec" && size(#[3])==n && allPositive(#[3])==1)
     1201        if (typeof(#[3])=="intvec" && size(#[3])==n && allPositive(#[3])==1)
    12021202        {
    1203           u0 = #[3];
    1204         }
     1203          u0 = #[3];
     1204        }
    12051205      }
    12061206    }
     
    12221222"USAGE:  bfctSyz(f [,r,s,t,u,v]);  f a poly, r,s,t,u optional ints, v an optional intvec
    12231223RETURN:  list of ideal and intvec
    1224 PURPOSE: computes the roots of the Bernstein-Sato polynomial b(s) 
     1224PURPOSE: computes the roots of the Bernstein-Sato polynomial b(s)
    12251225@*       for the hypersurface defined by f
    12261226ASSUME:  The basering is a commutative polynomial ring in char 0.
    1227 BACKGROUND: In this proc, the initial Malgrange ideal is computed according to the algorithm 
     1227BACKGROUND: In this proc, the initial Malgrange ideal is computed according to the algorithm
    12281228@*       by Masayuki Noro and then a system of linear equations is solved by computing syzygies.
    1229 NOTE:    In the output list, the ideal contains all the roots 
     1229NOTE:    In the output list, the ideal contains all the roots
    12301230@*       and the intvec their multiplicities.
    1231 @*       If r<>0, @code{std} is used for GB computations in characteristic 0, 
    1232 @*       otherwise, and by default, @code{slimgb} is used. 
    1233 @*       If s<>0, a matrix ordering is used for GB computations, otherwise, 
     1231@*       If r<>0, @code{std} is used for GB computations in characteristic 0,
     1232@*       otherwise, and by default, @code{slimgb} is used.
     1233@*       If s<>0, a matrix ordering is used for GB computations, otherwise,
    12341234@*       and by default, a block ordering is used.
    1235 @*       If t<>0, the computation of the intersection is solely performed over 
     1235@*       If t<>0, the computation of the intersection is solely performed over
    12361236@*       charasteristic 0, otherwise and by default, a modular method is used.
    1237 @*       If u<>0 and by default, @code{std} is used for GB computations in 
    1238 @*       characteristic >0, otherwise, @code{slimgb} is used. 
    1239 @*       If v is a positive weight vector, v is used for homogenization 
     1237@*       If u<>0 and by default, @code{std} is used for GB computations in
     1238@*       characteristic >0, otherwise, @code{slimgb} is used.
     1239@*       If v is a positive weight vector, v is used for homogenization
    12401240@*       computations, otherwise and by default, no weights are used.
    12411241DISPLAY: If printlevel=1, progress debug messages will be printed,
     
    12751275        if (typeof(#[3])=="int" || typeof(#[3])=="number")
    12761276        {
    1277           pIntersectchar = int(#[3]);
    1278         }
    1279         if (size(#)>3)
     1277          pIntersectchar = int(#[3]);
     1278        }
     1279        if (size(#)>3)
    12801280        {
    1281           if (typeof(#[4])=="int" || typeof(#[4])=="number")
     1281          if (typeof(#[4])=="int" || typeof(#[4])=="number")
    12821282          {
    1283             modengine = int(#[4]);
    1284           }
    1285           if (size(#)>4)
    1286           {
    1287             if (typeof(#[5])=="intvec" && size(#[5])==n && allPositive(#[5])==1)
    1288             {
    1289               u0 = #[5];
    1290             }
    1291           }
    1292         }
     1283            modengine = int(#[4]);
     1284          }
     1285          if (size(#)>4)
     1286          {
     1287            if (typeof(#[5])=="intvec" && size(#[5])==n && allPositive(#[5])==1)
     1288            {
     1289              u0 = #[5];
     1290            }
     1291          }
     1292        }
    12931293      }
    12941294    }
     
    13151315BACKGROUND:  In this proc, the initial ideal of I is computed according to the algorithm by
    13161316@*       Masayuki Noro and then a system of linear equations is solved by linear reductions.
    1317 NOTE:    In the output list, the ideal contains all the roots 
     1317NOTE:    In the output list, the ideal contains all the roots
    13181318@*       and the intvec their multiplicities.
    13191319@*       If s<>0, @code{std} is used for GB computations in characteristic 0,
    1320 @*       otherwise, and by default, @code{slimgb} is used. 
     1320@*       otherwise, and by default, @code{slimgb} is used.
    13211321@*       If t<>0, a matrix ordering is used for GB computations, otherwise,
    13221322@*       and by default, a block ordering is used.
     
    13741374"USAGE:  bfctOneGB(f [,s,t]);  f a poly, s,t optional ints
    13751375RETURN:  list of ideal and intvec
    1376 PURPOSE: computes the roots of the Bernstein-Sato polynomial b(s) for the 
     1376PURPOSE: computes the roots of the Bernstein-Sato polynomial b(s) for the
    13771377@*       hypersurface defined by f, using only one GB computation
    13781378ASSUME:  The basering is a commutative polynomial ring in char 0.
    1379 BACKGROUND: In this proc, the initial Malgrange ideal is computed based on the 
     1379BACKGROUND: In this proc, the initial Malgrange ideal is computed based on the
    13801380@*       algorithm by Masayuki Noro and combined with an elimination ordering.
    1381 NOTE:    In the output list, the ideal contains all the roots 
     1381NOTE:    In the output list, the ideal contains all the roots
    13821382@*       and the intvec their multiplicities.
    13831383@*       If s<>0, @code{std} is used for the GB computation, otherwise,
    1384 @*       and by default, @code{slimgb} is used. 
     1384@*       and by default, @code{slimgb} is used.
    13851385@*       If t<>0, a matrix ordering is used for GB computations,
    13861386@*       otherwise, and by default, a block ordering is used.
     
    14921492"USAGE:  bfctAnn(f [,r]);  f a poly, r an optional int
    14931493RETURN:  list of ideal and intvec
    1494 PURPOSE: computes the roots of the Bernstein-Sato polynomial b(s) 
     1494PURPOSE: computes the roots of the Bernstein-Sato polynomial b(s)
    14951495@*       for the hypersurface defined by f
    14961496ASSUME:  The basering is a commutative polynomial ring in char 0.
    14971497BACKGROUND: In this proc, ann(f^s) is computed and then a system of linear
    14981498@*       equations is solved by linear reductions.
    1499 NOTE:    In the output list, the ideal contains all the roots 
     1499NOTE:    In the output list, the ideal contains all the roots
    15001500@*       and the intvec their multiplicities.
    15011501@*       If r<>0, @code{std} is used for GB computations,
    1502 @*       otherwise, and by default, @code{slimgb} is used. 
     1502@*       otherwise, and by default, @code{slimgb} is used.
    15031503DISPLAY: If printlevel=1, progress debug messages will be printed,
    15041504@*       if printlevel>=2, all the debug messages will be printed.
Note: See TracChangeset for help on using the changeset viewer.