Changeset 6e3023a in git


Ignore:
Timestamp:
Dec 18, 2013, 2:57:02 PM (9 years ago)
Author:
Hans Schoenemann <hannes@…>
Branches:
(u'jengelh-datetime', 'ceac47cbc86fe4a15902392bdbb9bd2ae0ea02c6')(u'spielwiese', 'a800fe4b3e9d37a38c5a10cc0ae9dfa0c15a4ee6')
Children:
3a9e92a717b194d4a06a429b452258c0269579076235be46a7d664256bb99b1bd6340cd8ad800e2c
Parents:
74c446222122656dc7859919b51f3806768e7fd331b00db0e867fa73d4415e16969423de4337f6f7
Message:
Merge pull request #437 from ederc/sig-comparisons

Improvements for SBA
Files:
30 edited
1 moved

Legend:

Unmodified
Added
Removed
  • .gitignore

    r31b00db r6e3023a  
    88*.lo
    99*.la
     10*.gcno
     11*.gcda
    1012*~
    1113.libs
  • Singular/LIB/COPYING

    r31b00db r6e3023a  
    134134paramet.lib     Thomas Keilen                   keilen@mathematik.uni-kl.de
    135135perron.lib      Oleksandr Motsak                motsak at mathematik dot uni minus kl dot de
    136 phindex.lib 
     136phindex.lib
    137137pointid.lib     Stefan Steidel                  steidel@mathematik.uni-kl.de
    138138poly.lib        Olaf Bachmann                   obachman@mathematik.uni-kl.de
    139139                Gert-Martin Greuel              greuel@mathematik.uni-kl.de
    140140                Anne Fruehbis-Krueger           anne@mathematik.uni-kl.de
    141 oldpolymake.lib Thomas Markwig                  keilen@mathematik.uni-kl.de
     141polymake.lib    Thomas Markwig                  keilen@mathematik.uni-kl.de
    142142presolve.lib    Gert-Martin Greuel              greuel@mathematik.uni-kl.de
    143143primdec.lib     Gerhard Pfister                 pfister@mathematik.uni-kl.de
  • Singular/LIB/divisors.lib

    r31b00db r6e3023a  
    335335ASSUME:  A is a divisor on X.
    336336RETURN:  a list with a basis of the space of global sections of D.
    337 THEORY:  We assume that the qring of X satisfies the S2-condition. We compute sat((f*J) : I) /f
     337THEORY:  We assume that the qring of X satisfies the S2-condition and that X is smooth. We compute sat((f*J) : I) /f
    338338         where D = (I)-(J).
    339339KEYWORDS: divisors
     
    370370"
    371371{
    372   return( deg(std(A.num))-deg(std(A.den)));
     372  return( mult(std(A.num))-mult(std(A.den)));
    373373}
    374374example
  • Singular/LIB/gitfan.lib

    r31b00db r6e3023a  
    11///////////////////////////////////////////////////////////////////
    2 version="version gitfan.lib 4.0.0.0 Jun_2013 "; // $Id$
     2version="version gitfan.lib 4.0.0.0 Jun_2013 ";  // $Id$
    33category="Algebraic Geometry";
    44info="
     
    2727
    2828LIB "parallel.lib"; // for parallelWaitAll
    29 LIB "general.lib";
    3029
    3130////////////////////////////////////////////////////
    32 
    33 proc int2face(int n, int r)
     31proc mod_init()
     32{
     33  LIB"gfanlib.so";
     34}
     35
     36static proc int2face(int n, int r)
    3437{
    3538  int k = r-1;
     
    4043    while(2^k > n0){
    4144      k--;
     45      //v[size(v)+1] = 0;
    4246    }
    4347
     
    5256/////////////////////////////////
    5357
    54 proc isAface(ideal a, intvec gam0, int n)
     58proc isAface(ideal a, intvec gam0)
    5559"USAGE:  isAface(a,gam0); a: ideal, gam0:intvec
    5660PURPOSE: Checks whether the face of the positive orthant,
     
    6266"
    6367{
     68  poly pz;
     69
    6470  // special case: gam0 is the zero-cone:
    6571  if (size(gam0) == 1 and gam0[1] == 0){
    66     poly pz; ideal G;
     72    ideal G;
    6773
    6874    // is an a-face if and only if RL0(0,...,0) = const
     
    8793  }
    8894
    89   // ring is too big: Switch to KK[T_i | e_i\in gam0] instead:
    90   intvec w=ringlist(basering)[3][1][2];
    91   intvec w0;
     95
     96  // ring is too big: Switch to KK[T_i; e_i\in gam0] instead:
    9297  string initNewRing = "ring Rgam0 = 0,(";
    9398
    9499  for (int i=1; i<size(gam0); i++){
    95100    initNewRing = initNewRing + string(var(gam0[i])) + ",";
    96     w0[i]=w[gam0[i]];
    97   }
    98 
    99   w0 = w0,w[gam0[i]],1;
    100   initNewRing = initNewRing + string(var(gam0[size(gam0)])) + ",U),Wp("+string(w0)+");";
     101  }
     102
     103  initNewRing = initNewRing + string(var(gam0[size(gam0)])) + "),dp;";
    101104  def R = basering;
    102105  execute(initNewRing);
    103106
    104107  ideal agam0 = imap(R,a);
    105 
    106   for (i = 1; i<=size(agam0); i=i+1)
    107   {
    108     if (size(agam0[i]) == 1)
    109     { return(0,0); }
    110   }
    111108
    112109  poly p = var(1); // first entry of g; p = prod T_i with i element of g
     
    117114
    118115  agam0 = agam0, p - 1; // rad-membership
    119   agam0 = timeStd(agam0,5);
    120   // "timestd finished after: "+string(timer-t);
    121   // int timeout = 0;
    122   if (attrib(agam0,"isSB") < 1)
    123   {
    124     kill agam0; kill Rgam0; kill initNewRing; kill w; kill w0; setring R; kill R;
    125     return(0,1);
    126     // // "timestd failed in "+string(gam0)+", falling back to saturation!";
    127     // setring R; kill Rgam0;
    128 
    129     // initNewRing = "ring Rgam0 = 0,(";
    130 
    131     // w0 = 0;
    132     // for (i=1; i<size(gam0); i++){
    133     //   initNewRing = initNewRing + string(var(gam0[i])) + ",";
    134     //   w0[i]=w[gam0[i]];
    135     // }
    136 
    137     // w0 = w0,w[gam0[i]];
    138     // initNewRing = initNewRing + string(var(gam0[size(gam0)])) + "),Wp("+string(w0)+");";
    139     // execute(initNewRing);
    140 
    141     // ideal G = imap(R,a);
    142 
    143     // timeout = 1;
    144     // int t=rtimer;
    145     // for(int k=nvars(basering); k>0; k--) { G=sat(G,var(k))[1]; }
    146     // t = (rtimer - t) div 1000;
    147     // string(n)+": saturation successful after "+string(t)+" with: "+string(G<>1);
    148   }
     116  ideal G = std(agam0);
    149117
    150118  // does G contain 1?, i.e. is G = 1?
    151   if(agam0 <> 1) {
    152     kill agam0; kill Rgam0; kill initNewRing; kill w; kill w0; setring R; kill R;
    153     return(1,0); // true
    154   }
    155 
    156   kill agam0; kill Rgam0; kill initNewRing; kill w; kill w0; setring R; kill R;
    157   return(0,0); // false
     119  if(G <> 1) {
     120    return(1); // true
     121  }
     122
     123  return(0); // false
    158124}
    159125example
     
    172138}
    173139
    174 
    175 proc isAfaceNonZero(ideal a, intvec gam0)
    176 {
    177   string initNewRing = "ring Rgam0 = 0,(";
    178   for (int i=1; i<size(gam0); i++)
    179     { initNewRing = initNewRing + string(var(gam0[i])) + ","; }
    180   initNewRing = initNewRing + string(var(gam0[size(gam0)])) + "),dp;";
    181   def R = basering;
    182   execute(initNewRing);
    183 
    184   ideal agam0 = imap(R,a);
    185 
    186   for ( i = 1; i<=size(agam0); i=i+1)
    187     { if (size(agam0[i]) == 1) { return(0); } }
    188 
    189   poly p = var(1);
    190   for ( i = 2; i <= nvars(basering); i++ )
    191     { p = p * var(i); }
    192 
    193   agam0 = agam0, p - 1;
    194   ideal G = std(agam0);
    195 
    196   if(G <> 1)
    197     { return(1); }
    198 
    199   return(0);
    200 }
    201140////////////////////////////////////////////////////
    202141
     
    206145  list AF;
    207146
    208   for(i = start; i <= end; i=i+1){
     147  for(i = start; i <= end; i++){
    209148    if(i < 2^r){
    210       string(i)+"    "+string(size(AF));
    211149      gam0 = int2face(i,r);
    212150
     
    227165
    228166proc afaces(ideal a, list #)
    229 "USAGE:  afaces(a, b, c); a: ideal, b: int, c: int
     167"USAGE:  afaces(a, b, c); a: ideal, d: int, c: int
    230168PURPOSE: Returns a list of all a-faces (represented by intvecs).
    231169         Moreover, it is possible to specify a dimensional bound b,
     
    272210
    273211  // do remaining ones:
    274   for(i = ncores * step +1; i < 2^r; i=i+1){
     212  for(i = ncores * step +1; i < 2^r; i++){
    275213    "another one needed";
    276214    gam0 = int2face(i,r);
     
    285223  }
    286224
    287   int l;
    288225  // read out afaces of out into AF:
    289   for(l = 1; l <= size(out); l++){
    290     AF = AF + out[l];
     226  for(i = 1; i <= size(out); i++){
     227    AF = AF + out[i];
    291228  }
    292229
  • Singular/LIB/grobcov.lib

    r31b00db r6e3023a  
    113113                   Can be used without calling presiouly setglobalrings(),
    114114
    115 locus(G);     Special routine for determining the locus of points
     115locus(G):     Special routine for determining the locus of points
    116116                   of  objects. Given a parametric ideal J with
    117117                   parameters (a_1,..a_m) and variables (x_1,..,xn),
     
    131131                   Can be used without calling presiouly setglobalrings(),
    132132
    133 locusdg(G);  Is a special routine for computing the locus in dinamical geometry.
     133locusdg(G):  Is a special routine for computing the locus in dinamical geometry.
    134134                    It detects automatically a possible point that is to be avoided by the mover,
    135135                    whose coordinates must be the last two coordinates in the definition of the ring.
     
    139139                    Can be used without calling presiouly setglobalrings(),
    140140
    141 locusto(L);  Transforms the output of locus to a string that
     141locusto(L):  Transforms the output of locus to a string that
    142142                  can be reed from different computational systems.
    143143                  Can be used without calling presiouly setglobalrings(),
    144144
    145 addcons(L); Given a disjoint set of locally closed subsets in P-representation,
     145addcons(L): Given a disjoint set of locally closed subsets in P-representation,
    146146                  it returns the canonical P-representation of the constructible set
    147147                  formed by the union of them.
     
    12221222proc PtoCrep(list L)
    12231223{
     1224  int te=0;
    12241225  def RR=basering;
     1226  if(defined(@P)){te=1;}
     1227  else{te=0; setglobalrings();}
    12251228  setring(@P);
    12261229  def Lp=imap(RR,L);
     
    12401243  def La=list(ida,idb);
    12411244  setring(RR);
    1242   return(imap(@P,La));
     1245  def LL=imap(@P,La);
     1246  if(te==0){kill @P; kill @R; kill @RP;}
     1247  return(LL);
    12431248}
    12441249
     
    38953900  int i; int j; int k; int l;
    38963901  intvec v; intvec v1; intvec v0;
    3897   ideal J; list K;
     3902  ideal J; list K; list L1;
     3903  for(i=1;i<=size(L);i++)
     3904  {
     3905    if(equalideals(L[i][1][1],ideal(1))==0)
     3906    {L1[size(L1)+1]=L[i];}
     3907  }
     3908  L=L1;
    38983909  for(i=1;i<=size(L);i++)
    38993910  {
  • Singular/LIB/modstd.lib

    r31b00db r6e3023a  
    3333static proc mod_init()
    3434{
    35    newstruct("ideal_primeTest", "ideal Ideal");
     35   newstruct("idealPrimeTest", "ideal Ideal");
    3636}
    3737
     
    487487      list arguments;
    488488      int neededPrimes = neededSize-size(L);
    489       ideal_primeTest Id;
     489      idealPrimeTest Id;
    490490      Id.Ideal = I;
    491491      export(Id);
  • Singular/LIB/normal.lib

    r31b00db r6e3023a  
    44454445        printlevel = printlevel + 1;
    44464446        ideal JDefault = 0; // Now it uses the default J;
    4447         list nor1 = normalM(Id1, decomp, withDelta, denomOption, JDefault, JDefault)[1];
    4448         list nor2 = normalM(Id2, decomp, withDelta, denomOption, JDefault, JDefault)[1];
     4447        list nor1 = normalM(Id1, decomp, withDelta, denomOption, JDefault, JDefault);
     4448        list nor2 = normalM(Id2, decomp, withDelta, denomOption, JDefault, JDefault);
    44494449        printlevel = printlevel - 1;
    44504450        option(set,save_opt);
    4451         return(list(nor1, nor2));
     4451        list res = nor1 + nor2;
     4452        return(res);
    44524453      }
    44534454    }
     
    45684569
    45694570      ideal JDefault = 0;  // Now it uses the default J;
    4570       list nor1 = normalM(Id1, decomp, withDelta, denomOption, JDefault, JDefault)[1];
    4571       list nor2 = normalM(Id2, decomp, withDelta, denomOption, JDefault, JDefault)[1];
     4571      list nor1 = normalM(Id1, decomp, withDelta, denomOption, JDefault, JDefault);
     4572      list nor2 = normalM(Id2, decomp, withDelta, denomOption, JDefault, JDefault);
    45724573      printlevel = printlevel - 1;
    45734574      option(set,save_opt);
    4574       return(list(nor1, nor2));
     4575      list res = nor1 + nor2;
     4576      return(res);
    45754577    }
    45764578  }
  • Singular/LIB/polymake.lib

    r31b00db r6e3023a  
    1 ////
    2 version="version oldpolymake.lib 4.0.0.0 Jun_2013 "; // $Id$
     1version="version polymake.lib 4.0.0.0 Jun_2013 ";
    32category="Tropical Geometry";
    43info="
    5 LIBRARY:   oldpolymake.lib  Computations with polytopes and fans,
    6                             interface to polymake and TOPCOM
     4LIBRARY:   polymake.lib    Computations with polytopes and fans,
     5                           interface to polymake and TOPCOM
    76AUTHOR:    Thomas Markwig,  email: keilen@mathematik.uni-kl.de
    87
     
    109   Most procedures will not work unless polymake or topcom is installed and
    1110   if so, they will only work with the operating system LINUX!
    12    For more detailed information see the following note or consult the
     11   For more detailed information see IMPORTANT NOTE respectively consult the
    1312   help string of the procedures.
    1413
    15 NOTE:
     14   The conventions used in this library for polytopes and fans, e.g. the
     15   length and labeling of their vertices resp. rays, differs from the conventions
     16   used in polymake and thus from the conventions used in the polymake
     17   extension polymake.so of Singular. We recommend to use the newer polymake.so
     18   whenever possible.
     19
     20IMPORTANT NOTE:
    1621   Even though this is a Singular library for computing polytopes and fans
    1722   such as the Newton polytope or the Groebner fan of a polynomial, most of
     
    1924@* - polymake by Ewgenij Gawrilow, TU Berlin and Michael Joswig, TU Darmstadt
    2025@*   (see http://www.math.tu-berlin.de/polymake/),
    21 @* respectively (only in the procedure triangularions) by the program
    22 @* - topcom by Joerg Rambau, Universitaet Bayreuth (see http://www.uni-bayreuth.de/
    23      departments/wirtschaftsmathematik/rambau/TOPCOM);
    24 @*   This library should rather be seen as an interface which allows to use a
     26@* respectively (only in the procedure triangulations) by the program
     27@* - topcom by Joerg Rambau, Universitaet Bayreuth (see
     28@*   http://www.uni-bayreuth.de/departments/wirtschaftsmathematik/rambau/TOPCOM);
     29@*   this library should rather be seen as an interface which allows to use a
    2530     (very limited) number of options which polymake respectively topcom offers
    2631     to compute with polytopes and fans and to make the results available in
     
    3237     independent of both, polymake and topcom.
    3338
    34 PROCEDURES:
    35   polymakePolytope()   computes the vertices of a polytope using polymake
    36   newtonPolytopeP()     computes the Newton polytope of a polynomial
    37   newtonPolytopeLP()   computes the lattice points of the Newton polytope
    38   normalFan()          computes the normal fan of a polytope
    39   groebnerFan()        computes the Groebner fan of a polynomial
    40   intmatToPolymake()   transforms an integer matrix into polymake format
    41   polymakeToIntmat()   transforms polymake output into an integer matrix
    42 
    43   triangulations()     computes all triangulations of a marked polytope
    44   secondaryPolytope()  computes the secondary polytope of a marked polytope
    45 
    46   secondaryFan()       computes the secondary fan of a marked polytope
    47 
    48   cycleLength()      computes the cycleLength of cycle
    49   splitPolygon()     splits a marked polygon into vertices, facets, interior points
    50   eta()              computes the eta-vector of a triangulation
    51   findOrientedBoundary()    computes the boundary of a convex hull
    52   cyclePoints()      computes lattice points connected to some lattice point
    53   latticeArea()      computes the lattice area of a polygon
    54   picksFormula()     computes the ingrediants of Pick's formula for a polygon
    55   ellipticNF()       computes the normal form of an elliptic polygon
    56   ellipticNFDB()     displays the 16 normal forms of elliptic polygons
    57 
    58   polymakeKeepTmpFiles()  determines if the files created in /tmp should be kept
     39PROCEDURES USING POLYMAKE:
     40  polymakePolytope()  computes the vertices of a polytope using polymake
     41  newtonPolytopeP()    computes the Newton polytope of a polynomial
     42  newtonPolytopeLP()  computes the lattice points of the Newton polytope
     43  normalFanL()         computes the normal fan of a polytope
     44  groebnerFan()       computes the Groebner fan of a polynomial
     45
     46PROCEDURES USING TOPCOM:
     47  triangulations()    computes all triangulations of a marked polytope
     48  secondaryPolytope() computes the secondary polytope of a marked polytope
     49
     50PROCEDURES USING POLYMAKE AND TOPCOM:
     51  secondaryFan()      computes the secondary fan of a marked polytope
     52
     53PROCEDURES CONERNED WITH PLANAR POLYGONS:
     54  cycleLength()    computes the cycleLength of cycle
     55  splitPolygon()   splits a marked polygon into vertices, facets, interior points
     56  eta()            computes the eta-vector of a triangulation
     57  findOrientedBoundary()  computes the boundary of a convex hull
     58  cyclePoints()    computes lattice points connected to some lattice point
     59  latticeArea()    computes the lattice area of a polygon
     60  picksFormula()   computes the ingrediants of Pick's formula for a polygon
     61  ellipticNF()     computes the normal form of an elliptic polygon
     62  ellipticNFDB()   displays the 16 normal forms of elliptic polygons
    5963
    6064KEYWORDS:    polytope; fan; secondary fan; secondary polytope; polymake;
     
    8589////////////////////////////////////////////////////////////////////////////////
    8690
     91static proc mod_init ()
     92{
     93  LIB "gfanlib.so";
     94  LIB "polymake.so";
     95}
    8796
    8897/////////////////////////////////////////////////////////////////////////////
     
    9099/////////////////////////////////////////////////////////////////////////////
    91100
    92 proc polymakePolytope (intmat polytop,list #)
    93 "USAGE:  polymakePolytope(polytope[,#]);   polytope list, # string
    94 ASSUME:  each row of polytope gives the coordinates of a lattice point of a
     101proc polymakePolytope (intmat points)
     102"USAGE:  polymakePolytope(points);   polytope intmat
     103ASSUME:  each row of points gives the coordinates of a lattice point of a
    95104         polytope with their affine coordinates as given by the output of
    96105         secondaryPolytope
    97106PURPOSE: the procedure calls polymake to compute the vertices of the polytope
    98107         as well as its dimension and information on its facets
    99 RETURN:  list L with four entries
     108RETURN:  list, L with four entries
    100109@*            L[1] : an integer matrix whose rows are the coordinates of vertices
    101110                     of the polytope
    102111@*            L[2] : the dimension of the polytope
    103 @*            L[3] : a list whose i-th entry explains to which vertices the
     112@*            L[3] : a list whose ith entry explains to which vertices the
    104113                     ith vertex of the Newton polytope is connected
    105114                     -- i.e. L[3][i] is an integer vector and an entry k in
    106115                        there means that the vertex L[1][i] is connected to the
    107116                        vertex L[1][k]
    108 @*            L[4] : an integer matrix whose rows mulitplied by
     117@*            L[4] : an matrix of type bigintmat whose rows mulitplied by
    109118                     (1,var(1),...,var(nvar)) give a linear system of equations
    110119                     describing the affine hull of the polytope,
     
    118127         convention which starts indexing its vertices by zero while we start
    119128         with one !
    120 @*    -  the procedure creates the file  /tmp/polytope.polymake which contains
    121          the polytope in polymake format; if you wish to use this for further
    122          computations with polymake, you have to use the procedure
    123          polymakeKeepTmpFiles in before
    124 @*    -  moreover, the procedure creates the file /tmp/polytope.output which
    125          it deletes again before ending
    126 @*    -  it is possible to provide an optional second argument a string
    127          which then will be used instead of 'polytope' in the name of the
    128          polymake output file
    129129EXAMPLE: example polymakePolytope;   shows an example"
    130130{
    131   // the header for the file secendarypolytope.polymake
    132   string sp="_application polytope
    133 _version 2.2
    134 _type RationalPolytope
    135 
    136 POINTS
    137 ";
     131  // add a first column to polytope as homogenising coordinate
     132  points=intmatAddFirstColumn(points,"points");
     133  polytope polytop=polytopeViaPoints(points);
     134  list graph=vertexAdjacencyGraph(polytop)[2];
    138135  int i,j;
    139   // set the name for the polymake output file
    140   if (size(#)>0)
    141   {
    142     if (typeof(#[1])=="string")
    143     {
    144       string dateiname=#[1];
    145     }
    146     else
    147     {
    148       string dateiname="polytope";
    149     }
    150   }
    151   else
    152   {
    153     string dateiname="polytope";
    154   }
    155   // create the lattice point list for polymake
    156   sp=sp+intmatToPolymake(polytop,"points");
    157   // initialise dateiname.polymake and compute the vertices
    158   write(":w /tmp/"+dateiname+".polymake",sp);
    159   system("sh","cd /tmp; polymake "+dateiname+".polymake VERTICES > "+dateiname+".output");
    160   string vertices=read("/tmp/"+dateiname+".output");
    161   system("sh","/bin/rm /tmp/"+dateiname+".output");
    162   intmat np=polymakeToIntmat(vertices,"affine");
    163   // compute the dimension
    164   system("sh","cd /tmp; polymake "+dateiname+".polymake DIM > "+dateiname+".output");
    165   string pdim=read("/tmp/"+dateiname+".output");
    166   system("sh","/bin/rm /tmp/"+dateiname+".output");
    167   pdim=pdim[5,size(pdim)-6];
    168   execute("int nd="+pdim+";");
    169   // compute the vertex-edge graph
    170   system("sh","cd /tmp; polymake "+dateiname+".polymake GRAPH > "+dateiname+".output");
    171   string vertexedgegraph=read("/tmp/"+dateiname+".output");
    172   system("sh","/bin/rm /tmp/"+dateiname+".output");
    173   vertexedgegraph=vertexedgegraph[7,size(vertexedgegraph)-8];
    174   string newveg;
    175   for (i=1;i<=size(vertexedgegraph);i++)
    176   {
    177     if (vertexedgegraph[i]=="{")
    178     {
    179       newveg=newveg+"intvec(";
    180     }
    181     else
    182     {
    183       if (vertexedgegraph[i]=="}")
    184       {
    185         newveg=newveg+"),";
    186       }
    187       else
    188       {
    189         if (vertexedgegraph[i]==" ")
    190         {
    191           newveg=newveg+",";
    192         }
    193         else
    194         {
    195           newveg=newveg+vertexedgegraph[i];
    196         }
    197       }
    198     }
    199   }
    200   newveg=newveg[1,size(newveg)-1];
    201   execute("list nveg="+newveg+";");
    202   // raise each entry in nveg by one
    203   for (i=1;i<=size(nveg);i++)
    204   {
    205     for (j=1;j<=size(nveg[i]);j++)
    206     {
    207       nveg[i][j]=nveg[i][j]+1;
    208     }
    209   }
    210   // compute the affine hull
    211   system("sh","cd /tmp; polymake "+dateiname+".polymake AFFINE_HULL > "+dateiname+".output");
    212   string equations=read("/tmp/"+dateiname+".output");
    213   system("sh","/bin/rm /tmp/"+dateiname+".output");
    214   if (size(equations)>14)
    215   {
    216     intmat neq=polymakeToIntmat(equations,"cleardenom");
    217   }
    218   else
    219   {
    220     intmat neq[1][ncols(polytop)+1];
    221   }
    222   // delete the tmp-files, if polymakekeeptmpfiles is not set
    223   if (defined(polymakekeeptmpfiles)==0)
    224   {
    225     system("sh","/bin/rm /tmp/"+dateiname+".polymake");
    226   }
    227   // return the files
    228   return(list(np,nd,nveg,neq));
     136  for (i=1;i<=size(graph);i++)
     137  {
     138    for (j=1;j<=size(graph[i]);j++)
     139    {
     140      graph[i][j]=graph[i][j]+1;
     141    }
     142  }
     143  return(list(intmatcoldelete(vertices(polytop),1),dimension(polytop),graph,equations(polytop)));
    229144}
    230145example
     
    247162   ring r=0,x(1..4),dp;
    248163   matrix M[5][1]=1,x(1),x(2),x(3),x(4);
    249    np[4]*M;
     164   intmat(np[4])*M;
    250165}
    251166
    252167/////////////////////////////////////////////////////////////////////////////
    253168
    254 proc newtonPolytopeP (poly f,list #)
    255 "USAGE: newtonPolytopeP(f[,#]);  f poly, # string
    256 RETURN: list L with four entries
     169proc newtonPolytopeP (poly f)
     170"USAGE: newtonPolytopeP(f);  f poly
     171RETURN: list, L with four entries
    257172@*            L[1] : an integer matrix whose rows are the coordinates of vertices
    258173                     of the Newton polytope of f
     
    263178                        there means that the vertex L[1][i] is
    264179                        connected to the vertex L[1][k]
    265 @*            L[4] : an integer matrix whose rows mulitplied by
     180@*            L[4] : an matrix of type bigintmat whose rows mulitplied by
    266181                     (1,var(1),...,var(nvar)) give a linear system of equations
    267182                     describing the affine hull of the Newton polytope, i.e. the
     
    269184NOTE: -  if we replace the first column of L[4] by zeros, i.e. if we move
    270185         the affine hull to the origin, then we get the equations for the
    271          orthogonal comploment of the linearity space of the normal fan dual
     186         orthogonal complement of the linearity space of the normal fan dual
    272187         to the Newton polytope, i.e. we get the EQUATIONS that
    273188         we need as input for polymake when computing the normal fan
     
    275190         TU Berlin and Michael Joswig, so it only works if polymake is installed;
    276191         see http://www.math.tu-berlin.de/polymake/
    277 @*    -  the procedure creates the file  /tmp/newtonPolytope.polymake which
    278          contains the polytope in polymake format and which can be used for
    279          further computations with polymake
    280 @*    -  moreover, the procedure creates the file /tmp/newtonPolytope.output
    281          and deletes it again before ending
    282 @*    -  it is possible to give as an optional second argument a string
    283          which then will be used instead of 'newtonPolytope' in the name of
    284          the polymake output file
    285192EXAMPLE: example newtonPolytopeP;   shows an example"
    286193{
     
    296203    f=f-lead(f);
    297204  }
    298   if (size(#)==0)
    299   {
    300     #[1]="newtonPolytope";
    301   }
    302205  // call polymakePolytope with exponents
    303   return(polymakePolytope(exponents,#));
     206  return(polymakePolytope(exponents));
    304207}
    305208example
     
    330233   // the Newton polytope is contained in the affine space given
    331234   //     by the equations
    332    np[4]*M;
     235   intmat(np[4])*M;
    333236}
    334237
     
    339242RETURN: list, the exponent vectors of the monomials occuring in f,
    340243              i.e. the lattice points of the Newton polytope of f
    341 EXAMPLE: example normalFan;   shows an example"
     244EXAMPLE: example newtonPolytopeLP;   shows an example"
    342245{
    343246  list np;
     
    363266/////////////////////////////////////////////////////////////////////////////
    364267
    365 proc normalFan (intmat vertices,intmat affinehull,list graph,int er,list #)
    366 "USAGE:  normalFan (vert,aff,graph,rays,[,#]);   vert,aff intmat,  graph list, rays int, # string
     268proc normalFanL (def vertices, def affinehull,list graph,int er,list #)
     269"USAGE:  normalFanL (vert,aff,graph,rays,[,#]);   vert,aff intmat,  graph list, rays int, # string
    367270ASSUME:  - vert is an integer matrix whose rows are the coordinate of
    368271           the vertices of a convex lattice polytope;
     
    396299@*       - in the optional argument # it is possible to hand over other names
    397300           for the variables to be used -- be careful, the format must be correct
    398            which is not tested, e.g. if you want the variable names to be
    399            u00,u10,u01,u11 then you must hand over the string 'u11,u10,u01,u11'
    400 EXAMPLE: example normalFan;   shows an example"
    401 {
     301           and that is not tested, e.g. if you want the variable names to be
     302           u00,u10,u01,u11 then you must hand over the string u11,u10,u01,u11
     303EXAMPLE: example normalFanL;   shows an example"
     304{
     305  if (typeof(affinehull) != "intmat" && typeof (affinehull) != "bigintmat")
     306  {
     307    ERROR("normalFanL: input affinehull has to be either intmat or bigintmat");
     308    list L;
     309    return (L);
     310  }
    402311  list ineq; // stores the inequalities of the cones
    403312  int i,j,k;
     
    431340    list pp;  // contain strings representing the inequalities
    432341              // describing the normal cone
    433     intmat ie[size(graph[i])][ncols(vertices)]; // contains the inequalities
    434                                                 // as rows
     342    if (typeof (vertices) == "intmat")
     343    {
     344      intmat ie[size(graph[i])][ncols(vertices)]; // contains the inequalities
     345    }                                             // as rows
     346    if (typeof (vertices) == "bigintmat")
     347    {
     348      bigintmat ie[size(graph[i])][ncols(vertices)]; // contains the inequalities
     349    }                                                // as rows
    435350    // consider all the vertices to which the ith vertex in the
    436351    // polytope is connected by an edge
     
    469384  }
    470385  // remove the first column of affine hull to compute the linearity space
    471   intmat linearity=intmatcoldelete(affinehull,1);
     386  bigintmat linearity[1][ncols(vertices)];
     387  if (nrows(affinehull)>0)
     388  {
     389    linearity=intmatcoldelete(affinehull,1);
     390  }
    472391  //////////////////////////////////////////////////////////////////
    473392  // Compute next the extreme rays of the cones
     
    476395  {
    477396    list extremerays; // keeps the result
    478     string polymake; // keeps polymake output
    479     // the header for ineq.polymake
    480     string head="_application polytope
    481 _version 2.2
    482 _type RationalPolytope
    483 
    484 INEQUALITIES
    485 ";
    486     // the tail for both polymake files
    487     string tail="
    488 EQUATIONS
    489 ";
    490     tail=tail+intmatToPolymake(linearity,"rays");
    491     string ungleichungen; // keeps the inequalities for the polymake code
     397    cone kegel;
     398    bigintmat linearspan=intmatAddFirstColumn(linearity,"rays");
    492399    intmat M; // the matrix keeping the inequalities
    493     // create the file ineq.output
    494     write(":w /tmp/ineq.output","");
    495     int dimension; // keeps the dimension of the intersection the
    496                    // bad cones with the u11tobeseencone
    497400    for (i=1;i<=size(ineq);i++)
    498401    {
    499       i,". Cone of ",nrows(vertices); // indicate how many
    500                                       // vertices have been dealt with
    501       ungleichungen=intmatToPolymake(ineq[i][1],"rays");
    502       // write the inequalities to ineq.polymake and call polymake
    503       write(":w /tmp/ineq.polymake",head+ungleichungen+tail);
    504       ungleichungen=""; // clear ungleichungen
    505       system("sh","cd /tmp; /bin/rm ineq.output; polymake ineq.polymake VERTICES > ineq.output");
    506       // read the result of polymake
    507       polymake=read("/tmp/ineq.output");
    508       intmat VERT=polymakeToIntmat(polymake,"affine");
    509       extremerays[i]=VERT;
    510       kill VERT;
     402      kegel=coneViaInequalities(intmatAddFirstColumn(ineq[i][1],"rays"),linearspan);
     403      extremerays[i]=intmatcoldelete(rays(kegel),1);
    511404    }
    512405    for (i=1;i<=size(ineq);i++)
     
    514407      ineq[i]=ineq[i]+list(extremerays[i]);
    515408    }
    516   }
    517   // delete the tmp-files, if polymakekeeptmpfiles is not set
    518   if (defined(polymakekeeptmpfiles)==0)
    519   {
    520     system("sh","/bin/rm /tmp/ineq.polymake");
    521     system("sh","/bin/rm /tmp/ineq.output");
    522409  }
    523410  // get the linearity space
     
    534421   list np=newtonPolytopeP(f);
    535422   // the Groebner fan of f, i.e. the normal fan of the Newton polytope
    536    list gf=normalFan(np[1],np[4],np[3],1,"x,y,z");
     423   list gf=normalFanL(np[1],np[4],np[3],1,"x,y,z");
    537424   // the number of cones in the Groebner fan of f is:
    538425   size(gf[1]);
     
    549436/////////////////////////////////////////////////////////////////////////////
    550437
    551 proc groebnerFan (poly f,list #)
    552 "USAGE:  groebnerFan(f[,#]);  f poly, # string
     438proc groebnerFan (poly f)
     439"USAGE:  groebnerFan(f);  f poly
    553440RETURN:  list, the ith entry of L[1] contains information about the ith cone
    554441               in the Groebner fan dual to the ith vertex in the Newton
     
    564451                      in each cone
    565452@*             L[3] = the Newton polytope of f in the format of the procedure
    566                       newtonPolytope
    567 @*             L[4] = integer matrix where each row represents the exponet
     453                      newtonPolytopeP
     454@*             L[4] = integer matrix where each row represents the exponent
    568455                      vector of one monomial occuring in the input polynomial
    569456NOTE: - if you have already computed the Newton polytope of f then you might want
    570         to use the procedure normalFan instead in order to avoid doing costly
     457        to use the procedure normalFanL instead in order to avoid doing costly
    571458        computation twice
    572459@*    - the procedure calls for its computation polymake by Ewgenij Gawrilow,
    573460        TU Berlin and Michael Joswig, so it only works if polymake is installed;
    574461        see http://www.math.tu-berlin.de/polymake/
    575 @*    - the procedure creates the file  /tmp/newtonPolytope.polymake which
    576         contains the Newton polytope of f in polymake format and which can
    577         be used for further computations with polymake
    578 @*    - it is possible to give as an optional second argument as string which
    579         then will be used instead of 'newtonPolytope' in the name of the
    580         polymake output file
    581462EXAMPLE: example groebnerFan;   shows an example"
    582463{
     
    591472    f=f-lead(f);
    592473  }
    593   if (size(#)==0)
    594   {
    595     #[1]="newtonPolytope";
    596   }
    597474  // call polymakePolytope with exponents
    598   list newtonp=polymakePolytope(exponents,"newtonPolytope");
     475  list newtonp=polymakePolytope(exponents);
    599476  // get the variables as string
    600477  string variablen;
     
    604481  }
    605482  variablen=variablen[1,size(variablen)-1];
    606   // call normalFan in order to compute the Groebner fan
    607   list gf=normalFan(newtonp[1],newtonp[4],newtonp[3],1,variablen);
     483  // call normalFanL in order to compute the Groebner fan
     484  list gf=normalFanL(newtonp[1],newtonp[4],newtonp[3],1,variablen);
    608485  // append newtonp to gf
    609486  gf[3]=newtonp;
     
    642519
    643520
    644 /////////////////////////////////////////////////////////////////////////////
    645 
    646 proc intmatToPolymake (intmat M,string art)
    647 "USAGE:  intmatToPolymake(M,art);  M intmat, art string
    648 ASSUME:  - M is an integer matrix which should be transformed into polymake
    649            format;
    650 @*       - art is one of the following strings:
    651 @*           + 'rays'   : indicating that a first column of 0's should be added
    652 @*           + 'points' : indicating that a first column of 1's should be added
    653 RETURN:  string, the matrix is transformed in a string and a first column has
    654                  been added
    655 EXAMPLE: example intmatToPolymake;   shows an example"
    656 {
    657   if (art=="rays")
    658   {
    659     string anf="0 ";
    660   }
    661   else
    662   {
    663     string anf="1 ";
    664   }
    665   string sp;
    666   int i,j;
    667   // create the lattice point list for polymake
    668   for (i=1;i<=nrows(M);i++)
    669   {
    670     sp=sp+anf;
    671     for (j=1;j<=ncols(M);j++)
    672     {
    673       sp=sp+string(M[i,j])+" ";
    674       if (j==ncols(M))
    675       {
    676         sp=sp+"
    677 ";
    678       }
    679     }
    680   }
    681   return(sp);
    682 }
    683 example
    684 {
    685    "EXAMPLE:";
    686    echo=2;
    687    intmat M[3][4]=1,2,3,4,5,6,7,8,9,10,11,12;
    688    intmatToPolymake(M,"rays");
    689    intmatToPolymake(M,"points");
    690 }
    691 
    692 /////////////////////////////////////////////////////////////////////////////
    693 
    694 proc polymakeToIntmat (string pm,string art)
    695 "USAGE:  polymakeToIntmat(pm,art);  pm, art string
    696 ASSUME:  pm is the result of calling polymake with one 'argument' like
    697          VERTICES, AFFINE_HULL, etc., so that the first row of the string is
    698          the name of the corresponding 'argument' and the further rows contain
    699          the result which consists of vectors either over the integers
    700          or over the rationals
    701 RETURN:  intmat, the rows of the matrix are basically the vectors in pm, starting
    702                  from the second row, where each row has been multiplied with the
    703                  lowest common multiple of the denominators of its entries as if
    704                  it is an integer matrix; moreover, if art=='affine', then
    705                  the first column is omitted since we only want affine
    706                  coordinates
    707 EXAMPLE: example polymakeToIntmat;   shows an example"
    708 {
    709   // we need a line break
    710   string zeilenumbruch="
    711 ";
    712   // remove the 'argment' name, i.e. the first row of pm
    713   while (pm[1]!=zeilenumbruch)
    714   {
    715     pm=stringdelete(pm,1);
    716   }
    717   pm=stringdelete(pm,1);
    718   // find out how many entries each vector has, namely one more
    719   // than 'spaces' in a row
    720   int i=1;
    721   int s=1;
    722   int z=1;
    723   while (pm[i]!=zeilenumbruch)
    724   {
    725     if (pm[i]==" ")
    726     {
    727       s++;
    728     }
    729     i++;
    730   }
    731   // if we want to have affine coordinates
    732   if (art=="affine")
    733   {
    734     s--; // then there is one column less
    735     // and the entry of the first column (in the first row) has to be removed
    736     while (pm[1]!=" ")
    737     {
    738       pm=stringdelete(pm,1);
    739     }
    740     pm=stringdelete(pm,1);
    741   }
    742   // we add two line breaks at the end in order to have this as
    743   // a stopping criterion
    744   pm=pm+zeilenumbruch+zeilenumbruch;
    745   // we now have to work through each row
    746   for (i=1;i<=size(pm);i++)
    747   {
    748     // if there are two consecutive line breaks we are done
    749     if ((pm[i]==zeilenumbruch) and (pm[i+1]==zeilenumbruch))
    750     {
    751       i=size(pm)+1;
    752     }
    753     else
    754     {
    755       // a line break has to be replaced by a comma
    756       if (pm[i]==zeilenumbruch)
    757       {
    758         z++;
    759         pm[i]=",";
    760         // if we want to have affine coordinates,
    761         // then we have to delete the first entry in each row
    762         if (art=="affine")
    763         {
    764          while (pm[i+1]!=" ")
    765          {
    766            pm=stringdelete(pm,i+1);
    767          }
    768          pm=stringdelete(pm,i+1);
    769         }
    770       }
    771       // a space has to be replaced by a comma
    772       if (pm[i]==" ")
    773       {
    774         pm[i]=",";
    775       }
    776     }
    777   }
    778   // if we have introduced superflous commata at the end, they should be removed
    779   while (pm[size(pm)]==",")
    780   {
    781     pm=stringdelete(pm,size(pm));
    782   }
    783   // since the matrix could be over the rationals,
    784   // we need a ring with rational coefficients
    785   ring zwischering=0,x,lp;
    786   // create the matrix with the elements of pm as entries
    787   execute("matrix mm["+string(z)+"]["+string(s)+"]="+pm+";");
    788   // transform this into an integer matrix
    789   matrix M[1][ncols(mm)]; // takes a row of mm
    790   int cm;  // takes a lowest common multiple
    791   // multiply each row by an integer such that its entries are integers
    792   for (int j=1;j<=nrows(mm);j++)
    793   {
    794     M=mm[j,1..ncols(mm)];
    795     cm=commondenominator(M);
    796     for (i=1;i<=ncols(mm);i++)
    797     {
    798       mm[j,i]=cm*mm[j,i];
    799     }
    800   }
    801   // transform the matrix mm into an integer matrix
    802   execute("intmat im["+string(z)+"]["+string(s)+"]="+string(mm)+";");
    803   return(im);
    804 }
    805 example
    806 {
    807    "EXAMPLE:";
    808    echo=2;
    809    // this is the usual output of some polymake computation
    810    string pm="VERTICES
    811 0 1 3 5/3 1/3 -1 -23/3 -1/3 5/3 1/3 1
    812 0 1 3 -23/3 5/3 1 5/3 1/3 1/3 -1/3 -1
    813 0 1 1 1/3 -1/3 -1 5/3 1/3 -23/3 5/3 3
    814 0 1 1 5/3 -23/3 3 1/3 5/3 -1/3 1/3 -1
    815 0 1 -1 1/3 5/3 3 -1/3 -23/3 1/3 5/3 1
    816 0 1 -1 -1/3 1/3 1 1/3 5/3 5/3 -23/3 3
    817 0 1 -1 1 3 -5 -1 3 -1 1 -1
    818 0 1 -1 -1 -1 -1 1 1 3 3 -5
    819 0 1 -5 3 1 -1 3 -1 1 -1 -1
    820 
    821 ";
    822    intmat PM=polymakeToIntmat(pm,"affine");
    823    // note that the first column has been removed, since we asked for
    824    // affine coordinates, and the denominators have been cleared
    825    print(PM);
    826 }
    827 
    828521///////////////////////////////////////////////////////////////////////////////
    829522/// PROCEDURES USING TOPCOM
    830523///////////////////////////////////////////////////////////////////////////////
    831524
    832 proc triangulations (list polygon)
    833 "USAGE:  triangulations(polygon); list polygon
     525proc triangulations (list polygon,list #)
     526"USAGE:  triangulations(polygon[,#]); list polygon, list #
    834527ASSUME:  polygon is a list of integer vectors of the same size representing
    835528         the affine coordinates of the lattice points
     
    846539       this  procedure; see
    847540@*     http://www.uni-bayreuth.de/departments/wirtschaftsmathematik/rambau/TOPCOM
     541@*   - if you only want to have the regular triangulations the procedure should
     542       be called with the string 'regular' as optional argument
    848543@*   - the procedure creates the files /tmp/triangulationsinput and
    849544       /tmp/triangulationsoutput;
     
    851546       output containing the triangulations of corresponding to points in the
    852547       format of points2triangs; if you wish to use this for further
    853        computations with topcom, you have to use the procedure
    854        polymakeKeepTmpFiles in before
     548       computations with topcom, you have to call the procedure with the
     549       string 'keepfiles' as optional argument
    855550@*   - note that an integer i in an integer vector representing a triangle
    856551       refers to the ith lattice point, i.e. polygon[i]; this convention is
     
    860555{
    861556  int i,j;
     557  // check for optional arguments
     558  int regular,keepfiles;
     559  if (size(#)>0)
     560  {
     561    for (i=1;i<=size(#);i++)
     562    {
     563      if (typeof(#[i])=="string")
     564      {
     565        if (#[i]=="keepfiles")
     566        {
     567          keepfiles=1;
     568        }
     569        if (#[i]=="regular")
     570        {
     571          regular=1;
     572        }
     573      }
     574    }
     575  }
    862576  // prepare the input for points2triangs by writing the input polygon in the
    863577  // necessary format
     
    875589  write(":w /tmp/triangulationsinput",spi);
    876590  // call points2triangs
    877   system("sh","cd /tmp; points2triangs < triangulationsinput > triangulationsoutput");
     591  if (regular==1) // compute only regular triangulations
     592  {
     593    system("sh","cd /tmp; points2triangs --regular < triangulationsinput > triangulationsoutput");
     594  }
     595  else // compute all triangulations
     596  {
     597    system("sh","cd /tmp; points2triangs < triangulationsinput > triangulationsoutput");
     598  }
    878599  string p2t=read("/tmp/triangulationsoutput"); // takes result of points2triangs
    879   // delete the tmp-files, if polymakekeeptmpfiles is not set
    880   if (defined(polymakekeeptmpfiles)==0)
     600  // delete the tmp-files, if no second argument is given
     601  if (keepfiles==0)
    881602  {
    882603    system("sh","cd /tmp; rm -f triangulationsinput; rm -f triangulationsoutput");
     
    953674            else
    954675            {
    955               np2t=np2t+p2t[i];
     676              if (p2t[i]=="[")
     677              {
     678                // in Topcom version 17.4 (and maybe also in earlier versions) the list
     679                // of triangulations is indexed starting with index 0, in Singular
     680                // we have to start with index 1
     681                np2t=np2t+p2t[i]+"1+";
     682              }
     683              else
     684              {
     685                np2t=np2t+p2t[i];
     686              }
    956687            }
    957688          }
     
    962693  list T;
    963694  execute(np2t);
     695  // depending on the version of Topcom, the list T has or has not an entry T[1]
     696  // if it has none, the entry should be removed
     697  while (typeof(T[1])=="none")
     698  {
     699    T=delete(T,1);
     700  }
    964701  // raise each index by one
    965702  for (i=1;i<=size(T);i++)
     
    1000737RETURN:  list, say L, such that:
    1001738@*             L[1] = intmat, each row gives the affine coordinates of a lattice
    1002                       point in the secondary polytope given by the marked polytope
    1003                       corresponding to polygon
     739                      point in the secondary polytope given by the marked
     740                      polytope corresponding to polygon
    1004741@*             L[2] = the list of corresponding triangulations
    1005742NOTE: if the triangluations are not handed over as optional argument the
     
    1114851       see http://www.math.tu-berlin.de/polymake/
    1115852@*   - in the optional argument # it is possible to hand over other names for
    1116        the variables to be used -- be careful, the format must be correct
    1117        which is not tested, e.g. if you want the variable names to be
     853       the variables to be used -- be careful, the format must be correct and
     854       that is not tested, e.g. if you want the variable names to be
    1118855       u00,u10,u01,u11 then you must hand over the string 'u11,u10,u01,u11'
    1119856@*   - if the triangluations are not handed over as optional argument the
     
    1135872  list sp=secondaryPolytope(polygon,triang);
    1136873  list spp=polymakePolytope(sp[1]);
    1137   list sf=normalFan(spp[1],spp[4],spp[3],1);
     874  list sf=normalFanL(spp[1],spp[4],spp[3],1);
    1138875  return(list(sf[1],sf[2],spp,triang));
    1139876}
     
    13781115@*       triang is a list of integer vectors all of size three describing a
    13791116         triangulation of the polygon described by polygon; if an entry of
    1380          triang is the vector (i,j,k) then the triangle is built from the vertices
     1117         triang is the vector (i,j,k) then the triangle is built by the vertices
    13811118         with indices i, j and k
    13821119RETURN:  intvec, the integer vector eta describing that vertex of the Newton
     
    19481685  // points and
    19491686  // the number b of boundary points are connected by the formula: A=b+2g-2
    1950   return(list(area,bdpts,(area-bdpts+2)/2));
     1687  return(list(area,bdpts,(area-bdpts+2) div 2));
    19511688}
    19521689example
     
    19731710"USAGE:  ellipticNF(polygon);   polygon list
    19741711ASSUME:  polygon is a list of integer vectors in the plane such that their
    1975          convex hull C has precisely one interior lattice point, i.e. C is the
     1712         convex hull C has precisely one interior lattice point; i.e. C is the
    19761713         Newton polygon of an elliptic curve
    19771714PURPOSE: compute the normal form of the polygon with respect to the unimodular
     
    19911728  int i;            // index
    19921729  intvec edge;      // stores the vector of an edge
    1993   intvec boundary;  // stores lattice lengths of the edges of the Newton cylce
     1730  intvec boundary;  // stores lattice lengths of the edges of the Newton cycle
    19941731  // find the vertices of the Newton cycle and order it clockwise
    19951732  list pg=findOrientedBoundary(polygon)[2];
     
    23022039
    23032040/////////////////////////////////////////////////////////////////////////////////
    2304 /// AUXILARY PROCEDURES
    2305 /////////////////////////////////////////////////////////////////////////////////
    2306 
    2307 proc polymakeKeepTmpFiles (int i)
    2308 "USAGE:  polymakeKeepTmpFiles(int i);   i int
    2309 PURPOSE: some procedures create files in the directory /tmp which are used for
    2310          computations with polymake respectively topcom; these will be removed
    2311          when the corresponding procedure is left; however, it might be
    2312          desireable to keep them for further computations with either polymake or
    2313          topcom; this can be achieved by this procedure; call the procedure as:
    2314 @*       - polymakeKeepTmpFiles(1); - then the files will be kept
    2315 @*       - polymakeKeepTmpFiles(0); - then files will be removed in the future
    2316 RETURN:  none"
    2317 {
    2318   if ((i==1) and (defined(polymakekeeptmpfiles)==0))
    2319   {
    2320     int polymakekeeptmpfiles;
    2321     export(polymakekeeptmpfiles);
    2322   }
    2323   if (i!=1)
    2324   {
    2325     if (defined(polymakekeeptmpfiles))
    2326     {
    2327       kill polymakekeeptmpfiles;
    2328     }
    2329   }
    2330 }
    2331 
    2332 
    2333 /////////////////////////////////////////////////////////////////////////////////
    23342041/////////////////////////////////////////////////////////////////////////////////
    23352042/// AUXILARY PROCEDURES, WHICH ARE DECLARED STATIC
     
    23642071}
    23652072
    2366 static proc intmatcoldelete (intmat w,int i)
     2073static proc intmatcoldelete (def w,int i)
    23672074"USAGE:      intmatcoldelete(w,i); w intmat, i int
    23682075RETURN:      intmat, the integer matrix w with the ith comlumn deleted
    2369 NOTE:        the procedure is called by intmatsort and normalFan"
    2370 {
    2371   if ((i<1) or (i>ncols(w)) or (ncols(w)==1))
    2372   {
    2373     return(w);
    2374   }
    2375   if (i==1)
    2376   {
    2377     intmat M[nrows(w)][ncols(w)-1]=w[1..nrows(w),2..ncols(w)];
    2378     return(M);
    2379   }
    2380   if (i==ncols(w))
    2381   {
    2382     intmat M[nrows(w)][ncols(w)-1]=w[1..nrows(w),1..ncols(w)-1];
    2383     return(M);
    2384   }
    2385   else
    2386   {
    2387     intmat M[nrows(w)][i-1]=w[1..nrows(w),1..i-1];
    2388     intmat N[nrows(w)][ncols(w)-i]=w[1..nrows(w),i+1..ncols(w)];
    2389     return(intmatconcat(M,N));
     2076NOTE:        the procedure is called by intmatsort and normalFanL"
     2077{
     2078  if (typeof(w)=="intmat")
     2079  {
     2080    if ((i<1) or (i>ncols(w)) or (ncols(w)==1))
     2081    {
     2082      return(w);
     2083    }
     2084    if (i==1)
     2085    {
     2086      intmat M[nrows(w)][ncols(w)-1]=w[1..nrows(w),2..ncols(w)];
     2087      return(M);
     2088    }
     2089    if (i==ncols(w))
     2090    {
     2091      intmat M[nrows(w)][ncols(w)-1]=w[1..nrows(w),1..ncols(w)-1];
     2092      return(M);
     2093    }
     2094    else
     2095    {
     2096      intmat M[nrows(w)][i-1]=w[1..nrows(w),1..i-1];
     2097      intmat N[nrows(w)][ncols(w)-i]=w[1..nrows(w),i+1..ncols(w)];
     2098      return(intmatconcat(M,N));
     2099    }
     2100  }
     2101  if (typeof(w)=="bigintmat")
     2102  {
     2103    if ((i<1) or (i>ncols(w)) or (ncols(w)==1))
     2104    {
     2105      return(w);
     2106    }
     2107    if (i==1)
     2108    {
     2109      bigintmat M[nrows(w)][ncols(w)-1]=w[1..nrows(w),2..ncols(w)];
     2110      return(M);
     2111    }
     2112    if (i==ncols(w))
     2113    {
     2114      bigintmat M[nrows(w)][ncols(w)-1]=w[1..nrows(w),1..ncols(w)-1];
     2115      return(M);
     2116    }
     2117    else
     2118    {
     2119      bigintmat MN[nrows(w)][ncols(w)-1];
     2120      MN[1..nrows(w),1..i-1]=w[1..nrows(w),1..i-1];
     2121      MN[1..nrows(w),i..ncols(w)-1]=w[1..nrows(w),i+1..ncols(w)];
     2122      return(MN);
     2123    }
     2124  } else
     2125  {
     2126    ERROR("intmatcoldelete: input matrix has to be of type intmat or bigintmat");
     2127    intmat M; return(M);
    23902128  }
    23912129}
     
    26832421  return(list(coord,coords,latex));
    26842422}
     2423
     2424static proc intmatAddFirstColumn (def M,string art)
     2425"USAGE:  intmatAddFirstColumn(M,art);  M intmat, art string
     2426ASSUME:  - M is an integer matrix where a first column of 0's or 1's should be added
     2427@*       - art is one of the following strings:
     2428@*         + 'rays'   : indicating that a first column of 0's should be added
     2429@*         + 'points' : indicating that a first column of 1's should be added
     2430RETURN:  intmat, a first column has been added to the matrix"
     2431{
     2432  if (typeof (M) == "intmat")
     2433  {
     2434    intmat N[nrows(M)][ncols(M)+1];
     2435    int i,j;
     2436    for (i=1;i<=nrows(M);i++)
     2437    {
     2438      if (art=="rays")
     2439      {
     2440        N[i,1]=0;
     2441      }
     2442      else
     2443      {
     2444        N[i,1]=1;
     2445      }
     2446      for (j=1;j<=ncols(M);j++)
     2447      {
     2448        N[i,j+1]=M[i,j];
     2449      }
     2450    }
     2451    return(N);
     2452  }
     2453  if (typeof (M) == "bigintmat")
     2454  {
     2455    bigintmat N[nrows(M)][ncols(M)+1];
     2456    int i,j;
     2457    for (i=1;i<=nrows(M);i++)
     2458    {
     2459      if (art=="rays")
     2460      {
     2461        N[i,1]=0;
     2462      }
     2463      else
     2464      {
     2465        N[i,1]=1;
     2466      }
     2467      for (j=1;j<=ncols(M);j++)
     2468      {
     2469        N[i,j+1]=M[i,j];
     2470      }
     2471    }
     2472    return(N);
     2473  }
     2474  else
     2475  {
     2476    ERROR ("intmatAddFirstColumn: input matrix has to be either intmat or bigintmat");
     2477    intmat N;
     2478    return (N);
     2479  }
     2480}
     2481
     2482
     2483
  • Singular/LIB/primdec.lib

    r31b00db r6e3023a  
    709709    return(l);
    710710  }
     711  def op = option(get);
    711712  def P=basering;
    712713  int i,j,k,m,q,d,o;
    713   int n=nvars(basering);
     714  int n = nvars(basering);
    714715  ideal s,t,u,sact;
    715716  poly ni;
     
    769770        {
    770771          setring RL;
     772
    771773          option(redSB);
    772774          t=std(t);
     
    785787          rp=delete(rp,1);
    786788          keep1=keep1+rp;
     789
    787790          option(noredSB);
    788791        }
     
    803806  }
    804807  l=l+keep1;
     808  option(set,op);
    805809  return(l);
    806810}
     
    18981902   list @res,empty;
    18991903   ideal ser;
    1900    option(redSB);
     1904   def op = option( get );
     1905   option( redSB );
    19011906   list @pr=facstd(i);
    19021907   //if(size(@pr)==1)
     
    19231928//      }
    19241929//   }
    1925     option(noredSB);
     1930   // option( noredSB );
     1931   option( set, op );
    19261932   int j,k,odim,ndim,count;
    19271933   attrib(@pr[1],"isSB",1);
     
    20962102  }
    20972103
    2098   if(dim(i) == -1){setring P0;return(ideal(1));}
    2099   if((dim(i) == 0) && (npars(P) == 0))
     2104  if( dim(i) == -1 )
     2105  { 
     2106    option( set,op );
     2107    setring P0;
     2108    return( ideal(1) );
     2109  }
     2110  if( (dim(i) == 0 ) && ( npars(P) == 0) )
    21002111  {
    21012112    int di = vdim(i);
     
    25462557   int n=nvars(R);
    25472558
     2559   def op = option(get);
     2560
    25482561//---Anfang Provisorium
    25492562   if((size(i)==2) && (w==2))
    25502563   {
    2551       option(redSB);
    2552       ideal J=std(i);
    2553       option(noredSB);
    2554       if((size(J)==2)&&(deg(J[1])==1))
     2564      option( redSB );
     2565      ideal J = std(i);
     2566      option( noredSB );
     2567      if ((size(J)==2)&&(deg(J[1])==1))
    25552568      {
    25562569         ideal keep;
     
    25772590            resu[j]=L;
    25782591         }
     2592         option( set, op );
    25792593         return(resu);
    25802594      }
     
    26522666      re=convList(pr);
    26532667   }
    2654    return(re);
     2668   option( set, op );
     2669   return( re );
    26552670}
    26562671///////////////////////////////////////////////////////////////////////////////
     
    70097024      primary=resu;
    70107025    }
     7026    option(set,op);
    70117027    if (intersectOption == "intersect")
    70127028    {
  • Singular/LIB/tropical.lib

    r31b00db r6e3023a  
    1 ///////////////////////////////////////////////////////////
    2 version="version tropical.lib 4.0.0.0 Jun_2013 "; // $Id$
     1//
     2version="version tropical.lib 4.0.0.0 Dec_2013 ";
    33category="Tropical Geometry";
    44info="
     
    204204LIB "elim.lib";
    205205LIB "linalg.lib";
    206 LIB "oldpolymake.lib";
     206LIB "polymake.lib";
    207207LIB "primdec.lib";
    208208LIB "absfact.lib";
     
    650650          // pass first to a ring where a and @a
    651651          // are variables in order to use maps
    652           string @mp= string(minpoly);
     652          poly mp=minpoly;
    653653          ring INTERRING=char(LIFTRING),(t,@a,a),dp;
    654           execute("poly mp=" + @mp + ";");
     654          poly mp=imap(LIFTRING,mp);
    655655          ideal LIFT=imap(LIFTRING,LIFT);
    656656          kill LIFTRING;
     
    10121012  {
    10131013    def GLOBALRING=basering;
    1014     string @mp= string(minpoly);
     1014    number mp=minpoly;
    10151015    execute("ring LOCALRING=("+charstr(basering)+"),("+varstr(basering)+"),ds;");
    1016     execute("minpoly= " + @mp + ";");
    10171016    poly f=imap(GLOBALRING,f);
     1017    minpoly=imap(GLOBALRING,mp);
    10181018  }
    10191019  // check if a substitution is necessary
     
    10491049    if (minpoly!=0)
    10501050    {
    1051       string @mp= string(minpoly);
     1051      poly mp=minpoly;
    10521052      def OLDRING=basering;
    10531053      execute("ring NEWRING=0,("+varstr(basering)+","+parstr(basering)+"),ds;");
    1054       execute("ideal I = " + @mp + ";");
    1055       I = I,imap(OLDRING,f);
     1054      ideal I=imap(OLDRING,mp),imap(OLDRING,f);
    10561055    }
    10571056    else
     
    10651064      w=NewtP[jj]-NewtP[jj+1];
    10661065      ggteiler=gcd(w[1],w[2]);
    1067       zw=w[1]/ggteiler;
    1068       w[1]=w[2]/ggteiler;
     1066      zw=w[1] div ggteiler;
     1067      w[1]=w[2] div ggteiler;
    10691068      w[2]=zw;
    10701069      // if we have introduced a third variable for the parameter, then w needs a third component
     
    11341133    else
    11351134    {
    1136       string @mp= string(minpoly);
     1135      poly mp=minpoly;
    11371136      ring degreering=0,a,dp;
    1138       execute("poly mp=" + @mp + ";");
     1137      poly mp=imap(countring,mp);
    11391138      numberofbranchesfound=numberofbranchesfound+deg(mp);
    11401139      setring countring;
     
    11581157   ring r=0,(x,y),ds;
    11591158   poly f=x2-y4+x5y7;
    1160    puiseuxExpansion(puiseuxExpansion(f,3));
     1159   puiseuxExpansion(f,3,"subst");
     1160   displayPuiseuxExpansion(puiseuxExpansion(f,3));
    11611161}
    11621162
     
    15401540      ggt=gcd(normalvector[1],normalvector[2]);   // the gcd of the entries
    15411541      zw=normalvector[2];    // create the outward pointing normal vector
    1542       normalvector[2]=-normalvector[1]/ggt;
    1543       normalvector[1]=zw/ggt;
     1542      normalvector[2]=-normalvector[1] div ggt;
     1543      normalvector[1]=zw div ggt;
    15441544      if (size(#)==0) // we are computing w.r.t. minimum
    15451545      {
     
    27182718        if(k<>j)  // except the unit vector of the non-zero component
    27192719        {
    2720           intvec u; u[size(w)]=0; u[k]=1;
    2721           O[r,1..size(w)]=u; r=r+1;
     2720          intvec u;
     2721          u[size(w)]=0;
     2722          u[k]=1;
     2723          O[r,1..size(w)]=u;
     2724          r=r+1;
     2725          kill u;
    27222726        }
    27232727      }
     
    47084712        wneu=NewtP[1]-NewtP[2];
    47094713        int ggteiler=gcd(wneu[1],wneu[2]);
    4710         wneu[1]=-wneu[1]/ggteiler;
    4711         wneu[2]=wneu[2]/ggteiler;
     4714        wneu[1]=-wneu[1] div ggteiler;
     4715        wneu[2]=wneu[2] div ggteiler;
    47124716        if (wneu[1]>0)
    47134717        {
     
    48074811  for (jj=1;jj<=anzahlvariablen+numberdeletedvariables-1;jj++)
    48084812  {
    4809     PARA[jj]=(PARA[jj]+a[jj+1])*t^(tw[jj+1]*tweight/ww[1]);
     4813    PARA[jj]=(PARA[jj]+a[jj+1])*t^(tw[jj+1]*tweight div ww[1]);
    48104814  }
    48114815  // if we have reached the stop-level, i.e. either
     
    53535357    if (koeffizienten[j-ord(p)+1]!=0)
    53545358    {
    5355       if ((j-(N*wj)/w1)==0)
     5359      if ((j-(N*wj) div w1)==0)
    53565360      {
    53575361        dp=dp+displaycoef(koeffizienten[j-ord(p)+1]);
    53585362      }
    5359       if (((j-(N*wj)/w1)==1) and (N!=1))
     5363      if (((j-(N*wj) div w1)==1) and (N!=1))
    53605364      {
    53615365        dp=dp+displaycoef(koeffizienten[j-ord(p)+1])+"*t^(1/"+string(N)+")";
    53625366      }
    5363       if (((j-(N*wj)/w1)==1) and (N==1))
     5367      if (((j-(N*wj) div w1)==1) and (N==1))
    53645368      {
    53655369        dp=dp+displaycoef(koeffizienten[j-ord(p)+1])+"*t";
    53665370      }
    5367       if (((j-(N*wj)/w1)==-1) and (N==1))
     5371      if (((j-(N*wj) div w1)==-1) and (N==1))
    53685372      {
    53695373        dp=dp+displaycoef(koeffizienten[j-ord(p)+1])+"*1/t";
    53705374      }
    5371       if (((j-(N*wj)/w1)==-1) and (N!=1))
     5375      if (((j-(N*wj) div w1)==-1) and (N!=1))
    53725376      {
    53735377        dp=dp+displaycoef(koeffizienten[j-ord(p)+1])+"*1/t^(1/"+string(N)+")";
    53745378      }
    5375       if (((j-(N*wj)/w1)>1) and (N==1))
    5376       {
    5377         dp=dp+displaycoef(koeffizienten[j-ord(p)+1])+"*t^"+string(j-(N*wj)/w1);
    5378       }
    5379       if (((j-(N*wj)/w1)>1) and (N!=1))
    5380       {
    5381         dp=dp+displaycoef(koeffizienten[j-ord(p)+1])+"*t^("+string(j-(N*wj)/w1)+"/"+string(N)+")";
    5382       }
    5383       if (((j-(N*wj)/w1)<-1) and (N==1))
     5379      if (((j-(N*wj) div w1)>1) and (N==1))
     5380      {
     5381        dp=dp+displaycoef(koeffizienten[j-ord(p)+1])+"*t^"+string(j-(N*wj) div w1);
     5382      }
     5383      if (((j-(N*wj) div w1)>1) and (N!=1))
     5384      {
     5385        dp=dp+displaycoef(koeffizienten[j-ord(p)+1])+"*t^("+string(j-(N*wj) div w1)+"/"+string(N)+")";
     5386      }
     5387      if (((j-(N*wj) div w1)<-1) and (N==1))
    53845388      {
    53855389        dp=dp+displaycoef(koeffizienten[j-ord(p)+1])+"*1/t^"+string(wj-j);
    53865390      }
    5387       if (((j-(N*wj)/w1)<-1) and (N!=1))
    5388       {
    5389         dp=dp+displaycoef(koeffizienten[j-ord(p)+1])+"*1/t^("+string((N*wj)/w1-j)+"/"+string(N)+")";
     5391      if (((j-(N*wj) div w1)<-1) and (N!=1))
     5392      {
     5393        dp=dp+displaycoef(koeffizienten[j-ord(p)+1])+"*1/t^("+string((N*wj) div w1-j)+"/"+string(N)+")";
    53905394      }
    53915395      if (j<deg(p))
     
    54955499    {
    54965500      setring LIFTRing;
    5497       string @mp= string(minpoly);
     5501      poly mp=minpoly;
    54985502      execute("ring TESTRing=("+charstr(LIFTRing)+"),("+varstr(BASERING)+"),dp;");
    5499       execute("minpoly= " + @mp + ";");
     5503      minpoly=number(imap(LIFTRing,mp));
    55005504      ideal i=imap(BASERING,i);
    55015505    }
     
    57695773              wneu=NewtP[jjj]-NewtP[jjj+1];
    57705774              int ggteiler=gcd(wneu[1],wneu[2]);
    5771               wneu[1]=-wneu[1]/ggteiler;
    5772               wneu[2]=wneu[2]/ggteiler;
     5775              wneu[1]=-wneu[1] div ggteiler;
     5776              wneu[2]=wneu[2] div ggteiler;
    57735777              if (wneu[1]>0)
    57745778              {
     
    58695873        {
    58705874          execute("ring PARARing=("+string(char(basering))+",@a),t,ls;");
    5871           minpoly=number(imap(PREGFANRING,m)); // ?
     5875          minpoly=number(imap(PREGFANRING,m));
    58725876        }
    58735877        else
     
    59085912        for (jjj=1;jjj<=anzahlvariablen+numberdeletedvariables-1;jjj++)
    59095913        {
    5910           PARA[jjj]=(PARA[jjj]+zeros[size(zeros)][jjj+1])*t^(ww[jjj+1]*tweight/ww[1]);
     5914          PARA[jjj]=(PARA[jjj]+zeros[size(zeros)][jjj+1])*t^(ww[jjj+1]*tweight div ww[1]);
    59115915        }
    59125916        // delete the last entry in zero, since that one has
     
    60336037    {
    60346038      I=subst(i,var(lastvar),0);
    6035       while (subst(I[1],t,0)==0)
     6039      while ((I[1]!=0) and (subst(I[1],t,0)==0))
    60366040      {
    60376041        I[1]=I[1]/t;
     
    78147818    poly denom=delta*hn^5;
    78157819    poly ggt=gcd(num,denom);
    7816     num=num/ggt;
    7817     denom=denom/ggt;
     7820    num=num div ggt;
     7821    denom=denom div ggt;
    78187822    setring BASERING;
    78197823    poly num=imap(TRING,num);
  • Singular/ipid.cc

    r31b00db r6e3023a  
    433433  else if ((IDTYP(h)==RING_CMD)||(IDTYP(h)==QRING_CMD))
    434434    rKill(h);
    435   else
     435  else if (IDDATA(h)!=NULL)
    436436    s_internalDelete(IDTYP(h),IDDATA(h),r);
    437437  //  general  -------------------------------------------------------------
  • Singular/singular-libs

    r31b00db r6e3023a  
    2323        paraplanecurves.lib phindex.lib \
    2424        pointid.lib poly.lib \
    25         oldpolymake.lib presolve.lib primdec.lib primdecint.lib \
     25        presolve.lib primdec.lib primdecint.lib \
    2626        primitiv.lib qhmoduli.lib random.lib realclassify.lib \
    2727        realrad.lib reesclos.lib resbinomial.lib \
     
    3232        surf.lib surfacesignature.lib surfex.lib \
    3333        teachstd.lib toric.lib triang.lib \
    34         tropical.lib weierstr.lib zeroset.lib
     34        weierstr.lib zeroset.lib
    3535
    3636# libs in beta testing:
     
    4141# for new libs
    4242
    43 SLIB1 = classifyci.lib classify_aeq.lib ffsolve.lib decomp.lib template.lib\
     43SLIB1 = classifyci.lib classify_aeq.lib \
     44        divisors.lib \
     45        ffsolve.lib decomp.lib template.lib \
    4446        findifs.lib finitediff.lib \
     47        gitfan.lib \
    4548        locnormal.lib modnormal.lib \
    46         JMBTest.lib JMSConst.lib multigrading.lib parallel.lib \
     49        JMBTest.lib JMSConst.lib multigrading.lib\
     50        orbitparam.lib parallel.lib polymake.lib\
    4751        realizationMatroids.lib resources.lib ringgb.lib \
    4852        schreyer.lib symodstd.lib derham.lib polybori.lib ellipticcovers.lib \
    49         schubert.lib tasks.lib
     53        schubert.lib tasks.lib tropical.lib
    5054
    5155PLIBS = bfun.lib central.lib dmod.lib dmodapp.lib dmodvar.lib fpadim.lib \
     
    5357        ncfactor.lib nctools.lib perron.lib qmatrix.lib \
    5458        ncall.lib
    55        
    5659
     60
  • Singular/subexpr.cc

    r31b00db r6e3023a  
    455455void s_internalDelete(const int t,  void *d, const ring r)
    456456{
     457  assume(d!=NULL);
    457458  switch (t)
    458459  {
  • Tst/Short/bug_newstruct.res.gz.uu

    r31b00db r6e3023a  
    11begin 640 bug_newstruct.res.gz
    2 M'XL(")/1P5$``V)U9U]N97=S=')U8W0N<F5S`)6236O"0!"&[_D5+\%#`A+=
    3 MJJUMZI86+P$1P=*K^=K(0K*1_6CIO^_&F(W7GG9FWOEXAMGCYS;9`R`4N^0#
    4 MOE8ZJGGNQ][QICQ0V.")"ZZ#,/:Z%Y0B-^>38#]*2U/HR%J1TIEV50L*9R\C
    5 MS&9PR5#:5)535]$H!?ZEY-]<M=*?^C57-M<T329*-471"@:=\=H/1[C'NS%/
    6 M$:XENU%>T[ZLB#WKO2'MO!1<P8B255RPLA>8E*U$6Q1&2E:""U@W9U4K&?IF
    7 MM4W&^J5O81NF;L@SQ4"-@YM-YA2':,#?C%"$=$*WR&;`*C(!T6ID2O&SP"39
    8 M?[WODNTDF(?0+1K6Y,S25="_%W9=,EBN%N&_T`FQ[&ZNHR?VO!;Z&G[M.E//
    9 =,=_\(?/NI,2>M/L4W<F-"D@83[P_D)E)>TH"````
     2M'XL("+JRJ5(``V)U9U]N97=S=')U8W0N<F5S`)62RVK#,!!%]_Z*B\G"AN!$
     3M;=JF=:/2DHTAA$!*MXD?<A#84I#DEOY]I3B6L^U*,W/G<8;1_G.=;0$0BDWV
     4M@=!HDS2\"--@?U7N*&SPP`4W49P&[@6E*+K30;`?;517FL1:B3:Y\57W%-Y>
     5M))C-X).A35?77GU(1BD*SQ7_YEJJ<!HV7-O<KFUS4>DI2BD83,Z;,![A'F_&
     6M/"6XE&Q&>4G[LC(-K/>&H_..X!J=J%C-!:MZ@2DE%619=DJQ"ES`N@6KI6+H
     7MFS4V&<N7OH5M>/1#GBD&:NS\;#*GV"4#_FJ$(L0);I'5@%7F`D(:Y%KSD\`D
     8MVWZ];[+U))K',!(M:PMFZ6J8WS.[+!DM'A;QO]`)L>Q^KJ<G]KP6^A)^=9UI
     9>X)FO_I!Y<U)B3^H^A3MYIR,2IY/@#R.]3(I*`@``
    1010`
    1111end
  • Tst/Short/polymake.tst

    r31b00db r6e3023a  
    33LIB "tst.lib";
    44tst_init();
    5 LIB "oldpolymake.lib"; // TODO: update test result!
     5LIB "polymake.lib";
    66///////////////////////////////////////////////////////////////////////////
    77// A) Test for Procedures using Polymake
    88///////////////////////////////////////////////////////////////////////////
    99example polymakePolytope;
    10 example newtonPolytope;
     10example newtonPolytopeP;
    1111example newtonPolytopeLP;
    12 example normalFan;
     12example normalFanL;
    1313example groebnerFan;
    14 example intmatToPolymake;
    15 example polymakeToIntmat;
    1614///////////////////////////////////////////////////////////////////////////
    1715// B) Test for Procedures using Topcom
  • Tst/Short/sba_s.res.gz.uu

    r74c4462 r6e3023a  
    1 begin 640 sba_s.res.gz
    2 M'XL("-A/A5(``W-B85]S+G)E<P#M6UM/&T<4?N=7C%"DVO$E<Y_9.""EA#:H
    3 M4:@P>4BB-"+$%&-\D6U40I7_WMF+9\:[9S:[6F@5-3P`GLLWYWQSSG?.KF!X
    4 M^N+H-4*([*,G3]"CH\^/=H;9$.TC^SOKQ].K3V<?5_WU:HUZ/70]G_V)UJ/5
    5 M>H4NYLMXSJ[FR>J;U6B%EN=?SJ_'YVAT>S9=7)N!B^5\BHZ&P^<'/T4:+<X6
    6 MHZ7=)[P391^].OH9[9K3^M?C3[L#.Z/VD1G\.)Z-UZWV8"?^B?;W,^-FH[_Z
    7 MJ_79VJ[6^PXS2C!W%_/K+]N8!+M5A/318CD_1Y.S]>IF>39KC6=K=->V\X:7
    8 MO^T'0XSYBE?,T!ZZ<X@BF5B.#4DG9H91C%FW==O"_7YKUFZWNY\7;K%*43Z/
    9 MSJ[1V*S.CC8+W9HH69-\K=8)[J?1[/S2K-[][?GI\,W)\]YN)YWQ]U&2[!O=
    10 M+N;+=>O$FV#^Q-B;R$P?K6^6LU9RB#=I3/UJ/WCDTBA'VV6.-X8]WA@)\,98
    11 M.6_=RRWFF,@Q=SF?SO]L.?ZZGNU,U:*PL]M[Z4*$10$>.0GPR%D)CUQX/'+I
    12 M>.0JXS%-G'ST<>VQR*,`BX*4L=@CN?@3+,=B>K8?1D*$N3MX>_#JZ`",/J$"
    13 MK(DHP)HD):Q)YK$FN6--BFW6\L$GI4>;5`':9/0MVG+AIP@8?I8^/_H4J\/@
    14 M=O`I$:!1J0"-*BJA41./1DT=C9IE-([.YVGDS=J.-\V=/FO9SRB\,PX8;ARX
    15 M[H<XO-L.NPC'*^TGND&<=J^Z$[>*]QW%;C0Y/U9RY`$F1R_,J=8<@I-#3(5J
    16 M3??P8/IL[ZY'!]-.Q_E%,$UO)MZX>'S;FK;-M[NV@^#I_,3,8S<JT]$8^FIO
    17 MVB&#JQC<_-@&UYMK3^!1!]VV)O$!5[E3"-XLG'0Z;C2S[:L=X)ZQJ(=:YN3V
    18 MP-)(2&;6^+V9^!`O<EC:;74$T>1<BT_3:WAOKO0#,@N-C2;P$U-[OKF4;[AV
    19 ME%!9X)KDN*:ZGR-CZC#9MBF,]KV5SF+&,Q.WW6/2AB=AEO5BLAT>'+M,NXO9
    20 MV^SBR?F%9"*<>N,GWKA+",(3W\%\([%H6Z]X9/.-".SEVV4^X4@LX_9W%D@X
    21 MDDIS(.&V!8ND>FP_16#*D52#<RE')"NF')&BF'*IOI:F7*:T^91[?&E!,FG-
    22 M)=U&1+^1=$I43#JEH*3+Q-/=F29`TCV^_(-YB:=9,/&T@!)/JZUHUU%IXIG3
    23 MJ-T:D6+R1>R;R1>)DN2+MLV)(C_Y[+U03*#DHYC94*6XI%<H)-]6H:-8@1E(
    24 M<01F("4N0VC:BX,92(G7;E'BVBT:M]YQ`DS<4M/4Q@\TK<F>N=[)LSTYB(.C
    25 M;><CUTK05#QS3FXZMXD[/U75U<*L6U_DC:->A;L86-,R-36`YM&F->Z2+O;V
    26 MZ/RDTS&:ZJB;Q/[.5%7]27]G8LID?&T2OWO1/>DFIKII&7;WTO<WU=^`OZG,
    27 M%OSEM,1?SDO\Y;+$7ZY+_!6XU%]!`7_M<Z'GK^`E_@H)^BL*5^A9+7&)O[)P
    28 MA?Y.7N*OE*7^2EWB[]8%*USBL**@PZIPAY[9JA#MGMFJ<(?>3EV(=F^GIJ4.
    29 M:PXXG#2_GJ]IKQOP-6UW"[Y&A?OS+(X*D>Y9'!7NS]]9B'1_IR[SE6%(J]+&
    30 MPSG+<(E2,0PJ%<,E2L5PB5(Q4J)4C)0H%2.E2L7B-OCK3O9<Q7;>'9X<HY/#
    31 M%V\.3H^.7P^?FJUXY_3X=S3L;8WRG>'1NT-T_`OZ]:FM8<P.#M^^VPSS_S7F
    32 MAMC>R^_&Y.\%<_,FJGK0*@B=-[/X(3#E?XQIF:T1M3^XK889]_75(Y9!R!I"
    33 MEO\6)H$P(PA30)CB03`35FM$ZP]>JV%FY0M,!1"<T<HLJ,J@%*16-0.M8:E^
    34 M"$NU)1>.6Q@?E,2&]!(P(AK26\/2ZO36L%3;*E9=QRDHY-5M?A!0D`@)LAL@
    35 MXD%`+;LU5/<[XA<$;<IO'="XHE6/7`)">Z,>-JE>TIK#,K#+%R`L!F$)!,LA
    36 M6`$^/01@:<IP'>W%(#P8::2I^M:!A3D&*U%SCD&-"'*<E3BP.0$C3DBP++,(
    37 M/""J_F"I09KKX"H(5RO07D&;XO)ZN!NB:ZBQD/IAJ([`&&E.-07M;4XUK')!
    38 MW$WEJQ[4&M:[INK<'!94#J;!T"/@:Q@"UB@85\`"&L*U/->(:0V7@,9,-X:%
    39 M&8$SL#G3DM3"C2MA]6B6\$LZD`R.J\."!8N!#Z8,;KK@_I"!9%`"`5/X21?N
    40 M!^!^F8*UD*;O*`*Q#);:"$QO#K\-!>T&837\5@VD@X.-08!G\5`\1Z!N!'G.
    41 M2B&8AAI\S:@"RB3!F%;5D:6`7P@0`8:UA*4#K%I$J``O`FP=!=PZPMA,<;`F
    42 M4A["WI`.QS=,NU:P3-T#[8$V]1YHESB@)_=!N\8U:=_4R!IO*2EH/@.['3-:
    43 M^3%<@7)5!S<"<8F&.0&!*5C&`L`4['."P);K.N^9,5PL&[,-UY[F;%,,"M8]
    44 ML,WA:AD"CJMEC:B&'^\TW/'@ZLU.@&@%`FL0U[@()3N&VW?X,0FN[Q0L\(00
    45 M^,$.?--DA#ZA.O`B!*S(7,#O+#!8ZPG<$</0C,,7"=H>@8$=X)LHN,%LSC<E
    46 G<#*&^-[\Y86R?UW'XK\XCO]7)_Z/G)M5B[0'CW;^`87ADYAR-```
     1begin 644 sba_s.res.gz
     2M'XL("%K7JE(``W-B85]S+G)E<P"5F&UOVS80Q]_G4Q!&@,FQK8J/HN#%0)85
     3M:+!B`Y+NQ5!TA>,XM>-'6`J6>LAW'R7*Y)&2F-4O9(K'(T\_W?U)^^[3KS>_
     4M(X3P!+U[A\YO'L[/[NHN$B/3IG%ISN^G7_.XR`LT&J'U;OL-%?.\R-'C[E#:
     5MS&A6C7[.YSDZS+[/ULL9FK],-_NUZG@\[#;HYN[NZOJG3*+]=#\_&#\.5A0Q
     6M^GCS"^JIU>+U\KXW-I9T@E3GU^5V643]\5GYC2:3.KCM_)\X+Z:%&2TG=LZL
     7MFK.WWZV_NW/BQ([".$;[PVZ&5M,B?SY,M]%R6Z!CW]@5EW_-C0*C/N6(+;I$
     8M1SLCKPR'I8)TJRR4)`D=1B]1$L?1MM_O#Q_V=G"J9WF83]=HJ4;72ZN!=DQ6
     9MC:D^>5'->S_?SA9J=.^WJT]W?]Y>C7H#;8%^!%=^\Y?][E!$M\!`H6$)#'7H
     10M\^+YL(VJ18!1A?IJ;@!<DGG8%AXWF@!N%'=PHS3,;;APR%'ND5OL-KMOD>4W
     11M!+'3](<0#GJC#S9%:-;!D>$.CHP&.#(..#)A.;*TYJ@+Q\\^)@%%EG50Y#A$
     12M<82]_./4HZC7AFG$>3>[Z[^N/]Y<MV8?3SNH\:R#FL`!:H(":H)9:H*[U/SD
     13M$P)@$VD'-I&]A<U+OQ2WII_!![,OI3]"T$V^E'=@3-,.C&D6P"@QP"B)Q2AI
     14MC7$^V^G,V_8M-\FL/DL1UPB/Z@$4&SNYC+L8'MVTRY)RI+DCIQDWPZ?ARHYB
     15ML45L>ZOU2R5'8,)JZ;U:U82#DVH1M4-%F\MDO/GY\C@BX\U@8)\+)T2_F=)Q
     16M?_$2;?KJ<NS;*9BVKY0]L;U"]Y93/UUN!GC\5$ZNOMS)Y>FU5].C`7J)5N4"
     17M3]XJ.#D-7`T&MK>.[=5T,!`L&J%(K=P?&XP8UV$M/RO#EW*0G4M:5PN(5.N:
     18M^8E^#9_5*_V"U$`5HTK\*M01#)>P$VN+A(@&:^RQ)C+V8&SLG-0-A9(8C+01
     19M4U:'Z#X>%28],374F\7V_OH/6VG'DM[)BU7K-XH),P+Z;T&_+0C,JF=OK3=<
     20MBK9Y*I:9>L,\`?6V\`L.ES)NVK2CX+"6YHZ"<P4+:STV=UEKR6&MP5[)84&;
     21M)8<%;Y:<UM=@R=5*ZY?<Q<),4DNK5W0G$7VCZ%+^/XLN3=N*KA9/^\XD;BFZ
     22MB\7?%!2>I)V%)WE;X<G4R7:9!0M/K4:,:X:;Q9?1-XLOXX'BR]QPL@P6GWDO
     23M),%MQ4<2:E*5)(&S0J/XG(V.)&EK!9(D:ZU`@FV%$'T6;ZU`@L%QBV![W"+E
     24MT;LL@)4=J@ZUY0^::'4IQF@U(>HR&O6-.;,G":*UTWO&T\%M99?7HIKOU;CB
     25MT8^-@`WN<6PBJ\543:A^V43+(1TFP$?Z1BMC1,NH-1+HJ445&J$G<XW8\12^
     26M$7IZ`270DR6^$7AJ<5TMUTIKAH_#VV&%QYI9-^(%9*P%N(,QDZV,>1)@S$F`
     27M,6<!QEP$&',98"R2`&-!`HP%"S`6(LA8R!;&YN<O8)PF`<8I:66<L@#CM)'D
     28M(.I4!AC+1I(#3TD"C&4CR:&G"#"6C20'GOHTV\E8'V\[&#N)K`^\'9#UN;<!
     29M.6OH@0V;)HTLQ\#8T`/HV<ARZ-G0`^C9R'+@B1MZ`#QQ(\NA)PM!IOKPZT&N
     30M?LNLP!RRFR_5BN[SI:2A!2!BTLAP$+$OXPY?7\8=OKZ,.WQ]&7?X^C+N\/5E
     31MW.&K9;R3+VO;[_39%0!F@=V.LM;=CK+`;D=98+>COH8[@'T-=P#[&NX`]C7<
     32M`>QKN`/8UW`'L`CN=K14\=>S^N\`81JC#V>G/ZB$;:G>\C@E]%7=U<.Y:0`_
     33M;ENU']=7Z\=,`_@QVZK]F+Y:/VH:P(_:5NU']57=G1Z6FW,8+7>H\E_=\K_;
     340YSQ2Q,[/_@-@&(*'G!8`````
    4735`
    4836end
  • Tst/Short/sba_s.stat

    r74c4462 r6e3023a  
    1 1 >> tst_memory_0 :: 1384468440:3.1.3.sw, 32 bit:spielwiese:x86_64-Darwin:schmetterhemd.mathematik.uni-kl.de:19521252
    2 1 >> tst_memory_1 :: 1384468440:3.1.3.sw, 32 bit:spielwiese:x86_64-Darwin:schmetterhemd.mathematik.uni-kl.de:20971520
    3 1 >> tst_memory_2 :: 1384468440:3.1.3.sw, 32 bit:spielwiese:x86_64-Darwin:schmetterhemd.mathematik.uni-kl.de:21269452
    4 1 >> tst_timer_1 :: 1384468440:3.1.3.sw, 32 bit:spielwiese:x86_64-Darwin:schmetterhemd.mathematik.uni-kl.de:38
     11 >> tst_memory_0 :: 1385654086:3.1.3.sw, 64  1385654090:3.1.3.sw, 64  1385655331:3.1.3.sw, 64  1386104581:3.1.3.sw, 64  1386927962:3.1.3.sw, 64 bit:spielwiese:x86_64-Linux:schosch:245535672
     21 >> tst_memory_1 :: 1385654086:3.1.3.sw, 64  1385654090:3.1.3.sw, 64  1385655331:3.1.3.sw, 64  1386104581:3.1.3.sw, 64  1386927962:3.1.3.sw, 64 bit:spielwiese:x86_64-Linux:schosch:249974784
     31 >> tst_memory_2 :: 1385654086:3.1.3.sw, 64  1385654090:3.1.3.sw, 64  1385655331:3.1.3.sw, 64  1386104581:3.1.3.sw, 64  1386927962:3.1.3.sw, 64 bit:spielwiese:x86_64-Linux:schosch:250400768
     41 >> tst_timer_1 :: 1385654086:3.1.3.sw, 64  1385654090:3.1.3.sw, 64  1385655331:3.1.3.sw, 64  1386104581:3.1.3.sw, 64  1386927962:3.1.3.sw, 64 bit:spielwiese:x86_64-Linux:schosch:524
  • Tst/Short/sba_s.tst

    r74c4462 r6e3023a  
    124124
    125125int k;
    126 for (k=3; k<=6; k++)
     126for (k=6; k>2; k--)
    127127{
    128128  string bench = cyclicn(k);
    129129  sprintf(bench);
    130130  ideal f;
     131  f = sba(i,3,0);
     132  f = sba(i,3,1);
     133  f = sba(i,2,0);
     134  f = sba(i,2,1);
    131135  f = sba(i,1,0);
    132136  f = sba(i,1,1);
     
    137141  sprintf(bench);
    138142  ideal f;
     143  f = sba(i,3,0);
     144  f = sba(i,3,1);
     145  f = sba(i,2,0);
     146  f = sba(i,2,1);
    139147  f = sba(i,1,0);
    140148  f = sba(i,1,1);
     
    145153  sprintf(bench);
    146154  ideal f;
     155  f = sba(i,3,0);
     156  f = sba(i,3,1);
     157  f = sba(i,2,0);
     158  f = sba(i,2,1);
    147159  f = sba(i,1,0);
    148160  f = sba(i,1,1);
     
    153165  sprintf(bench);
    154166  ideal f;
     167  f = sba(i,3,0);
     168  f = sba(i,3,1);
     169  f = sba(i,2,0);
     170  f = sba(i,2,1);
    155171  f = sba(i,1,0);
    156172  f = sba(i,1,1);
     
    161177  sprintf(bench);
    162178  ideal f;
     179  f = sba(i,3,0);
     180  f = sba(i,3,1);
     181  f = sba(i,2,0);
     182  f = sba(i,2,1);
    163183  f = sba(i,1,0);
    164184  f = sba(i,1,1);
     
    169189  sprintf(bench);
    170190  ideal f;
     191  f = sba(i,3,0);
     192  f = sba(i,3,1);
     193  f = sba(i,2,0);
     194  f = sba(i,2,1);
    171195  f = sba(i,1,0);
    172196  f = sba(i,1,1);
  • debian/changelog

    r31b00db r6e3023a  
    1 singular (3.1.3.sw-1) unstable; urgency=low
     1singular (4.0.0-1) unstable; urgency=low
    22
    33  * Initial release.
  • doc/COPYING.texi

    r31b00db r6e3023a  
    131131@copyright{} Winfried Bruns and Bogdan Ichim
    132132@* @uref{http://www.mathematik.uni-osnabrueck.de/normaliz/}
    133 @item polymake (used by oldpolymake.lib, @pxref{oldpolymake_lib})
     133@item polymake (used by polymake.lib, @pxref{polymake_lib})
    134134@copyright{} Ewgenij Gawrilow and Michael Joswig
    135135@* @uref{http://www.polymake.de/}
     
    142142@copyright{} Oliver Labs (2001-2008), Stephan Holzer (2004-2005)
    143143@* @uref{http://surfex.AlgebraicSurface.net}
    144 @item TOPCOM (used by oldpolymake.lib, @pxref{oldpolymake_lib})
     144@item TOPCOM (used by polymake.lib, @pxref{polymake_lib})
    145145@copyright{} J@"org Rambau
    146146@* @uref{http://www.rambau.wm.uni-bayreuth.de/TOPCOM/}
  • doc/Makefile.bsp

    r31b00db r6e3023a  
    6363CHKSUM_DB       = ${DOC_SUBDIR}/chksum
    6464# which tags to avoid:
    65 DOC2TEX_EXAMPLE_EXCLUSIONS = -exclude oldpolymake
     65DOC2TEX_EXAMPLE_EXCLUSIONS = -exclude polymake
    6666# which tags to avoid:
    67 TAG             = oldpolymake
     67TAG             = polymake
    6868DOC2TEX         = ${PERL} ./doc2tex.pl -docdir ${DOC_SUBDIR} \
    6969                  -Singular ${SINGULAR} -verbose ${VERBOSE} -make ${MAKE} \
  • doc/NEWS.texi

    r31b00db r6e3023a  
    3434@item
    3535name spaces (@nref{package})
     36@item
     37algebraic/transcendental extension rewritten
    3638@end itemize
    3739Besides theses internal changes, @sc{Singular} version 4 offers many new
     
    4244@itemize
    4345@item @nref{sba}: an F5 like Groebner base computation
     46@item @nref{ASSUME}: support for debugging libraries
    4447@end itemize
    4548
     
    4750@itemize
    4851@item @nref{classifyci_lib}: Isolated complete intersection singularities in characteristic  0
    49 @item @xref{schubert_lib}:  Proceduces for Intersection Theory
     52@item @nref{derham_lib}: Computation of deRham cohomology
    5053@end itemize
    5154
     
    5457@item
    5558improved grobcov.lib (@nref{grobcov_lib})
     59@item
     60fixed normal.lib (@nref{normal_lib})
    5661@end itemize
    5762
  • kernel/kInline.h

    r74c4462 r6e3023a  
    11561156// dummy function for function pointer strat->rewCrit being usable in all
    11571157// possible choices for criteria
    1158 KINLINE BOOLEAN arriRewDummy(poly /*sig*/, unsigned long /*not_sevSig*/, kStrategy /*strat*/, int /*start=0*/)
     1158KINLINE BOOLEAN arriRewDummy(poly /*sig*/, unsigned long /*not_sevSig*/, poly /*lm*/, kStrategy /*strat*/, int /*start=0*/)
    11591159{
    11601160  return FALSE;
  • kernel/kspoly.cc

    r74c4462 r6e3023a  
    201201   * TODO:
    202202   * --------------------------------------------
    203    * if strat->incremental
     203   * if strat->sbaOrder == 1
    204204   * Since we are subdividing lower index and
    205205   * current index reductions it is enough to
     
    207207   * for a check. This should speed-up checking
    208208   * a lot!
    209    * if !strat->incremental
     209   * if !strat->sbaOrder == 0
    210210   * We are not subdividing lower and current index
    211211   * due to the fact that we are using the induced
     
    223223  if (!PW->is_sigsafe)
    224224  {
    225     poly f1 = p_Copy(PR->GetLmCurrRing(),currRing);
     225    poly ftmp;
     226    ring rtmp;
     227    PR->SetLmCurrRing();
     228    poly f1 = pCopy(PR->GetLmCurrRing());
     229    //PR->GetLm(ftmp,rtmp);
     230    //poly f1 = k_LmInit_tailRing_2_currRing(ftmp,rtmp);
    226231    poly f2 = PW->GetLmCurrRing();
    227232    poly sigMult = pCopy(PW->sig);   // copy signature of reducer
     
    235240    printf("--------------\n");
    236241#endif
    237     sigMult = pp_Mult_qq(f1,sigMult,currRing);
     242    sigMult = p_Mult_q(f1,sigMult,currRing);
    238243//#if 1
    239244#ifdef DEBUGF5
     
    253258
    254259#endif
    255     pDelete(&f1);
     260    //pDelete(&f1);
    256261    pDelete(&sigMult);
    257262    // go on with the computations only if the signature of p2 is greater than the
     
    262267      return 3;
    263268    }
    264     PW->is_sigsafe  = TRUE;
     269    //PW->is_sigsafe  = TRUE;
    265270  }
    266271  PR->is_redundant = FALSE;
  • kernel/kstd1.cc

    r74c4462 r6e3023a  
    14221422  // standard basis computations
    14231423  strat->red          = redSig;
    1424   //strat->incremental  = TRUE;
     1424  //strat->sbaOrder  = 1;
    14251425  strat->currIdx      = 1;
    14261426}
     
    22532253}
    22542254
    2255 ideal kSba(ideal F, ideal Q, tHomog h,intvec ** w, int incremental, int arri, intvec *hilb,int syzComp,
     2255ideal kSba(ideal F, ideal Q, tHomog h,intvec ** w, int sbaOrder, int arri, intvec *hilb,int syzComp,
    22562256          int newIdeal, intvec *vw)
    22572257{
     
    22632263  BOOLEAN delete_w=(w==NULL);
    22642264  kStrategy strat=new skStrategy;
    2265   if (incremental!=0)
    2266   {
    2267     strat->incremental = TRUE;
    2268   }
    2269   else
    2270   {
    2271     strat->incremental = FALSE;
    2272   }
     2265  strat->sbaOrder = sbaOrder;
    22732266  if (arri!=0)
    22742267  {
    22752268    strat->rewCrit1 = arriRewDummy;
    22762269    strat->rewCrit2 = arriRewCriterion;
     2270    strat->rewCrit3 = arriRewCriterionPre;
    22772271  }
    22782272  else
     
    22802274    strat->rewCrit1 = faugereRewCriterion;
    22812275    strat->rewCrit2 = faugereRewCriterion;
     2276    strat->rewCrit3 = faugereRewCriterion;
    22822277  }
    22832278
  • kernel/kstd2.cc

    r74c4462 r6e3023a  
    4848#endif
    4949
    50 #define F5C       0
     50#define F5C       1
    5151#if F5C
    5252  #define F5CTAILRED 1
    5353#endif
    5454
    55 #define SBA_PRODUCT_CRITERION         0
    56 #define SBA_PRINT_ZERO_REDUCTIONS     1
    57 #define SBA_PRINT_REDUCTION_STEPS     1
    58 #define SBA_PRINT_SIZE_G              1
    59 #define SBA_PRINT_SIZE_SYZ            1
    60 #define SBA_PRINT_PRODUCT_CRITERION   0
    61 
     55#define SBA_INTERRED_START                  0
     56#define SBA_TAIL_RED                        1
     57#define SBA_PRODUCT_CRITERION               0
     58#define SBA_PRINT_ZERO_REDUCTIONS           0
     59#define SBA_PRINT_REDUCTION_STEPS           0
     60#define SBA_PRINT_OPERATIONS                0
     61#define SBA_PRINT_SIZE_G                    0
     62#define SBA_PRINT_SIZE_SYZ                  0
     63#define SBA_PRINT_PRODUCT_CRITERION         0
     64
     65// counts sba's reduction steps
     66#if SBA_PRINT_REDUCTION_STEPS
    6267long sba_reduction_steps;
     68long sba_interreduction_steps;
     69#endif
     70#if SBA_PRINT_OPERATIONS
     71long sba_operations;
     72long sba_interreduction_operations;
     73#endif
     74
    6375/***********************************************
    6476 * SBA stuff -- done
     
    466478
    467479    ksReducePoly(h, &(strat->T[ii]), NULL, NULL, strat);
     480#if SBA_PRINT_REDUCTION_STEPS
     481    sba_interreduction_steps++;
     482#endif
     483#if SBA_PRINT_OPERATIONS
     484    sba_interreduction_operations  +=  pLength(strat->T[ii].p);
     485#endif
    468486
    469487#ifdef KDEBUG
     
    513531    }
    514532  }
     533}
     534
     535KINLINE int ksReducePolyTailSig(LObject* PR, TObject* PW, LObject* Red)
     536{
     537  BOOLEAN ret;
     538  number coef;
     539
     540  assume(PR->GetLmCurrRing() != PW->GetLmCurrRing());
     541  Red->HeadNormalize();
     542  /*
     543  printf("------------------------\n");
     544  pWrite(Red->GetLmCurrRing());
     545  */
     546  ret = ksReducePolySig(Red, PW, 1, NULL, &coef);
     547
     548
     549  if (!ret)
     550  {
     551    if (! n_IsOne(coef, currRing->cf))
     552    {
     553      PR->Mult_nn(coef);
     554      // HANNES: mark for Normalize
     555    }
     556    n_Delete(&coef, currRing->cf);
     557  }
     558  return ret;
    515559}
    516560
     
    623667    sigSafe = ksReducePolySig(h, &(strat->T[ii]), strat->S_2_R[ii], NULL, NULL, strat);
    624668#if SBA_PRINT_REDUCTION_STEPS
    625     sba_reduction_steps++;
     669    if (sigSafe != 3)
     670      sba_reduction_steps++;
     671#endif
     672#if SBA_PRINT_OPERATIONS
     673    if (sigSafe != 3)
     674      sba_operations  +=  pLength(strat->T[ii].p);
    626675#endif
    627676    // if reduction has taken place, i.e. the reduction was sig-safe
     
    685734    }
    686735  }
     736}
     737
     738// tail reduction for SBA
     739poly redtailSba (LObject* L, int pos, kStrategy strat, BOOLEAN withT, BOOLEAN normalize)
     740{
     741#define REDTAIL_CANONICALIZE 100
     742  strat->redTailChange=FALSE;
     743  if (strat->noTailReduction) return L->GetLmCurrRing();
     744  poly h, p;
     745  p = h = L->GetLmTailRing();
     746  if ((h==NULL) || (pNext(h)==NULL))
     747    return L->GetLmCurrRing();
     748
     749  TObject* With;
     750  // placeholder in case strat->tl < 0
     751  TObject  With_s(strat->tailRing);
     752
     753  LObject Ln(pNext(h), strat->tailRing);
     754  Ln.sig      = L->sig;
     755  Ln.sevSig   = L->sevSig;
     756  Ln.pLength  = L->GetpLength() - 1;
     757
     758  pNext(h) = NULL;
     759  if (L->p != NULL) pNext(L->p) = NULL;
     760  L->pLength = 1;
     761
     762  Ln.PrepareRed(strat->use_buckets);
     763
     764  int cnt=REDTAIL_CANONICALIZE;
     765  while(!Ln.IsNull())
     766  {
     767    loop
     768    {
     769      Ln.SetShortExpVector();
     770      if (withT)
     771      {
     772        int j;
     773        j = kFindDivisibleByInT(strat->T, strat->sevT, strat->tl, &Ln);
     774        if (j < 0) break;
     775        With = &(strat->T[j]);
     776      }
     777      else
     778      {
     779        With = kFindDivisibleByInS(strat, pos, &Ln, &With_s);
     780        if (With == NULL) break;
     781      }
     782      cnt--;
     783      if (cnt==0)
     784      {
     785        cnt=REDTAIL_CANONICALIZE;
     786        /*poly tmp=*/Ln.CanonicalizeP();
     787        if (normalize)
     788        {
     789          Ln.Normalize();
     790          //pNormalize(tmp);
     791          //if (TEST_OPT_PROT) { PrintS("n"); mflush(); }
     792        }
     793      }
     794      if (normalize && (!TEST_OPT_INTSTRATEGY) && (!nIsOne(pGetCoeff(With->p))))
     795      {
     796        With->pNorm();
     797      }
     798      strat->redTailChange=TRUE;
     799      int ret = ksReducePolyTailSig(L, With, &Ln);
     800#if SBA_PRINT_REDUCTION_STEPS
     801      if (ret != 3)
     802        sba_reduction_steps++;
     803#endif
     804#if SBA_PRINT_OPERATIONS
     805      if (ret != 3)
     806        sba_operations  +=  pLength(With->p);
     807#endif
     808      if (ret)
     809      {
     810        // reducing the tail would violate the exp bound
     811        //  set a flag and hope for a retry (in bba)
     812        strat->completeReduce_retry=TRUE;
     813        if ((Ln.p != NULL) && (Ln.t_p != NULL)) Ln.p=NULL;
     814        do
     815        {
     816          pNext(h) = Ln.LmExtractAndIter();
     817          pIter(h);
     818          L->pLength++;
     819        } while (!Ln.IsNull());
     820        goto all_done;
     821      }
     822      if (Ln.IsNull()) goto all_done;
     823      if (! withT) With_s.Init(currRing);
     824    }
     825    pNext(h) = Ln.LmExtractAndIter();
     826    pIter(h);
     827    pNormalize(h);
     828    L->pLength++;
     829  }
     830
     831  all_done:
     832  Ln.Delete();
     833  if (L->p != NULL) pNext(L->p) = pNext(p);
     834
     835  if (strat->redTailChange)
     836  {
     837    L->length = 0;
     838  }
     839
     840  //if (TEST_OPT_PROT) { PrintS("N"); mflush(); }
     841  //L->Normalize(); // HANNES: should have a test
     842  assume(kTest_L(L));
     843  return L->GetLmCurrRing();
    687844}
    688845
     
    772929
    773930    ksReducePoly(h, &(strat->T[ii]), NULL, NULL, strat);
     931#if SBA_PRINT_REDUCTION_STEPS
     932    sba_interreduction_steps++;
     933#endif
     934#if SBA_PRINT_OPERATIONS
     935    sba_interreduction_operations  +=  pLength(strat->T[ii].p);
     936#endif
    774937
    775938#ifdef KDEBUG
     
    9421105    number coef;
    9431106    ksReducePoly(h,&(strat->T[ii]),strat->kNoetherTail(),&coef,strat);
     1107#if SBA_PRINT_REDUCTION_STEPS
     1108    sba_interreduction_steps++;
     1109#endif
     1110#if SBA_PRINT_OPERATIONS
     1111    sba_interreduction_operations  +=  pLength(strat->T[ii].p);
     1112#endif
    9441113#ifdef KDEBUG
    9451114    if (TEST_OPT_DEBUG)
     
    15091678  //    induced Schreyer order.
    15101679  // The corresponding orders are computed in sbaRing(), depending
    1511   // on the flag strat->incremental
    1512   long zeroreductions     = 0;
    1513   long product_criterion  = 0;
    1514   long size_g             = 0;
    1515   long size_syz           = 0;
     1680  // on the flag strat->sbaOrder
     1681#if SBA_PRINT_ZERO_REDUCTIONS
     1682  long zeroreductions           = 0;
     1683#endif
     1684#if SBA_PRINT_PRODUCT_CRITERION
     1685  long product_criterion        = 0;
     1686#endif
     1687#if SBA_PRINT_SIZE_G
     1688  long size_g                   = 0;
     1689#endif
     1690#if SBA_PRINT_SIZE_SYZ
     1691  long size_syz                 = 0;
     1692#endif
    15161693  // global variable
    1517   sba_reduction_steps     = 0;
    1518 
    1519   ideal F = F0;
     1694#if SBA_PRINT_REDUCTION_STEPS
     1695  sba_reduction_steps           = 0;
     1696  sba_interreduction_steps      = 0;
     1697#endif
     1698#if SBA_PRINT_OPERATIONS
     1699  sba_operations                = 0;
     1700  sba_interreduction_operations = 0;
     1701#endif
     1702
     1703  ideal F1 = F0;
    15201704  ring sRing, currRingOld;
    15211705  currRingOld  = currRing;
    1522   if (strat->incremental)
     1706  if (strat->sbaOrder == 1 || strat->sbaOrder == 3)
    15231707  {
    15241708    sRing = sbaRing(strat);
     
    15261710    {
    15271711      rChangeCurrRing (sRing);
    1528       F = idrMoveR (F0, currRingOld, currRing);
    1529     }
    1530   }
    1531 #if 0
     1712      F1 = idrMoveR (F0, currRingOld, currRing);
     1713    }
     1714  }
     1715  // sort ideal F
     1716  ideal F       = idInit(IDELEMS(F1),F1->rank);
     1717  intvec *sort  = idSort(F1);
     1718  for (int i=0; i<sort->length();++i)
     1719    F->m[i] = F1->m[(*sort)[i]-1];
     1720#if SBA_INTERRED_START
     1721  F = kInterRed(F,NULL);
     1722#endif
     1723#if F5DEBUG
    15321724  printf("SBA COMPUTATIONS DONE IN THE FOLLOWING RING:\n");
    15331725  rWrite (currRing);
     1726  printf("ordSgn = %d\n",currRing->OrdSgn);
    15341727  printf("\n");
    15351728#endif
     
    15421735  int hilbeledeg=1,hilbcount=0,minimcnt=0;
    15431736  LObject L;
    1544   // BOOLEAN withT     = FALSE;
     1737  BOOLEAN withT     = TRUE;
    15451738  strat->max_lower_index = 0;
    15461739
     
    16051798    #endif
    16061799    if (strat->Ll== 0) strat->interpt=TRUE;
     1800    /*
    16071801    if (TEST_OPT_DEGBOUND
    16081802        && ((strat->honey && (strat->L[strat->Ll].ecart+currRing->pFDeg(strat->L[strat->Ll].p,currRing)>Kstd1_deg))
    16091803            || ((!strat->honey) && (currRing->pFDeg(strat->L[strat->Ll].p,currRing)>Kstd1_deg))))
    16101804    {
    1611       /*
    1612        *stops computation if
    1613        * 24 IN test and the degree +ecart of L[strat->Ll] is bigger then
    1614        *a predefined number Kstd1_deg
    1615        */
     1805     
     1806       //stops computation if
     1807       // 24 IN test and the degree +ecart of L[strat->Ll] is bigger then
     1808       //a predefined number Kstd1_deg
    16161809      while ((strat->Ll >= 0)
    16171810        && (strat->L[strat->Ll].p1!=NULL) && (strat->L[strat->Ll].p2!=NULL)
     
    16231816      else strat->noClearS=TRUE;
    16241817    }
    1625     if (strat->incremental && pGetComp(strat->L[strat->Ll].sig) != strat->currIdx)
     1818    */
     1819    if (strat->sbaOrder == 1 && pGetComp(strat->L[strat->Ll].sig) != strat->currIdx)
    16261820    {
    16271821      strat->currIdx  = pGetComp(strat->L[strat->Ll].sig);
     
    16301824      // 2. generation of new principal syzygy rules for syzCriterion
    16311825      f5c ( strat, olddeg, minimcnt, hilbeledeg, hilbcount, srmax,
    1632             lrmax, reduc, Q, w, hilb );
     1826          lrmax, reduc, Q, w, hilb );
    16331827#endif
    16341828      // initialize new syzygy rules for the next iteration step
    16351829      initSyzRules(strat);
     1830
    16361831    }
    16371832    /*********************************************************************
    1638      * interrreduction step is done, we can go on with the next iteration
    1639      * step of the signature-based algorithm
    1640      ********************************************************************/
     1833      * interrreduction step is done, we can go on with the next iteration
     1834      * step of the signature-based algorithm
     1835      ********************************************************************/
    16411836    /* picks the last element from the lazyset L */
    16421837    strat->P = strat->L[strat->Ll];
    16431838    strat->Ll--;
    1644 //#if 1
     1839    /* reduction of the element choosen from L */
     1840
     1841    if (!strat->rewCrit2(strat->P.sig, ~strat->P.sevSig, strat->P.GetLmCurrRing(), strat, strat->P.checked+1)) {
     1842      //#if 1
    16451843#ifdef DEBUGF5
    1646     Print("SIG OF NEXT PAIR TO HANDLE IN SIG-BASED ALGORITHM\n");
    1647     Print("-------------------------------------------------\n");
    1648     pWrite(strat->P.sig);
    1649     pWrite(pHead(strat->P.p));
    1650     pWrite(pHead(strat->P.p1));
    1651     pWrite(pHead(strat->P.p2));
    1652     Print("-------------------------------------------------\n");
    1653 #endif
    1654     if (pNext(strat->P.p) == strat->tail)
    1655     {
    1656       // deletes the short spoly
     1844      Print("SIG OF NEXT PAIR TO HANDLE IN SIG-BASED ALGORITHM\n");
     1845      Print("-------------------------------------------------\n");
     1846      pWrite(strat->P.sig);
     1847      pWrite(pHead(strat->P.p));
     1848      pWrite(pHead(strat->P.p1));
     1849      pWrite(pHead(strat->P.p2));
     1850      Print("-------------------------------------------------\n");
     1851#endif
     1852      if (pNext(strat->P.p) == strat->tail)
     1853      {
     1854        // deletes the short spoly
     1855        /*
    16571856#ifdef HAVE_RINGS
    1658       if (rField_is_Ring(currRing))
    1659         pLmDelete(strat->P.p);
     1857        if (rField_is_Ring(currRing))
     1858          pLmDelete(strat->P.p);
     1859        else
     1860#endif
     1861          pLmFree(strat->P.p);
     1862*/
     1863          // TODO: needs some masking
     1864          // TODO: masking needs to vanish once the signature
     1865          //       sutff is completely implemented
     1866          strat->P.p = NULL;
     1867        poly m1 = NULL, m2 = NULL;
     1868
     1869        // check that spoly creation is ok
     1870        while (strat->tailRing != currRing &&
     1871            !kCheckSpolyCreation(&(strat->P), strat, m1, m2))
     1872        {
     1873          assume(m1 == NULL && m2 == NULL);
     1874          // if not, change to a ring where exponents are at least
     1875          // large enough
     1876          if (!kStratChangeTailRing(strat))
     1877          {
     1878            WerrorS("OVERFLOW...");
     1879            break;
     1880          }
     1881        }
     1882        // create the real one
     1883        ksCreateSpoly(&(strat->P), NULL, strat->use_buckets,
     1884            strat->tailRing, m1, m2, strat->R);
     1885
     1886      }
     1887      else if (strat->P.p1 == NULL)
     1888      {
     1889        if (strat->minim > 0)
     1890          strat->P.p2=p_Copy(strat->P.p, currRing, strat->tailRing);
     1891        // for input polys, prepare reduction
     1892        strat->P.PrepareRed(strat->use_buckets);
     1893      }
     1894      if (strat->P.p == NULL && strat->P.t_p == NULL)
     1895      {
     1896        red_result = 0;
     1897      }
    16601898      else
    1661 #endif
    1662         pLmFree(strat->P.p);
    1663 
    1664       // TODO: needs some masking
    1665       // TODO: masking needs to vanish once the signature
    1666       //       sutff is completely implemented
    1667       strat->P.p = NULL;
    1668       poly m1 = NULL, m2 = NULL;
    1669 
    1670       // check that spoly creation is ok
    1671       while (strat->tailRing != currRing &&
    1672              !kCheckSpolyCreation(&(strat->P), strat, m1, m2))
    1673       {
    1674         assume(m1 == NULL && m2 == NULL);
    1675         // if not, change to a ring where exponents are at least
    1676         // large enough
    1677         if (!kStratChangeTailRing(strat))
    1678         {
    1679           WerrorS("OVERFLOW...");
    1680           break;
    1681         }
    1682       }
    1683       // create the real one
    1684       ksCreateSpoly(&(strat->P), NULL, strat->use_buckets,
    1685                     strat->tailRing, m1, m2, strat->R);
    1686 
    1687     }
    1688     else if (strat->P.p1 == NULL)
    1689     {
    1690       if (strat->minim > 0)
    1691         strat->P.p2=p_Copy(strat->P.p, currRing, strat->tailRing);
    1692       // for input polys, prepare reduction
    1693       strat->P.PrepareRed(strat->use_buckets);
    1694     }
    1695 
    1696     if (strat->P.p == NULL && strat->P.t_p == NULL)
    1697     {
    1698       red_result = 0;
    1699     }
    1700     else
    1701     {
    1702       if (TEST_OPT_PROT)
    1703         message((strat->honey ? strat->P.ecart : 0) + strat->P.pFDeg(),
    1704                 &olddeg,&reduc,strat, red_result);
    1705 
    1706 //#if 1
     1899      {
     1900        //#if 1
    17071901#ifdef DEBUGF5
    1708       Print("Poly before red: ");
    1709       pWrite(strat->P.p);
    1710 #endif
    1711       /* reduction of the element choosen from L */
    1712       if (!strat->rewCrit2(strat->P.sig, ~strat->P.sevSig, strat, strat->P.checked+1)) {
     1902        Print("Poly before red: ");
     1903        pWrite(pHead(strat->P.p));
     1904        pWrite(strat->P.sig);
     1905#endif
    17131906#if SBA_PRODUCT_CRITERION
    1714         if (strat->P.checked == 3) {
     1907        if (strat->P.prod_crit) {
     1908#if SBA_PRINT_PRODUCT_CRITERION
    17151909          product_criterion++;
    1716           enterSyz(strat->P, strat);
     1910#endif
     1911          int pos = posInSyz(strat, strat->P.sig);
     1912          enterSyz(strat->P, strat, pos);
    17171913          if (strat->P.lcm!=NULL)
    17181914            pLmFree(strat->P.lcm);
     
    17221918        }
    17231919#else
    1724       red_result = strat->red(&strat->P,strat);
    1725 #endif
    1726       } else {
    1727         if (strat->P.lcm!=NULL)
    1728           pLmFree(strat->P.lcm);
    1729         red_result = 2;
    1730       }
    1731       if (errorreported)  break;
    1732     }
     1920        red_result = strat->red(&strat->P,strat);
     1921#endif
     1922      }
     1923    } else {
     1924      /*
     1925      if (strat->P.lcm != NULL)
     1926        pLmFree(strat->P.lcm);
     1927        */
     1928      red_result = 2;
     1929    }
     1930    if (errorreported)  break;
     1931
     1932//#if 1
     1933#ifdef DEBUGF5
     1934    if (red_result != 0) {
     1935        Print("Poly after red: ");
     1936        pWrite(pHead(strat->P.p));
     1937        pWrite(strat->P.GetLmCurrRing());
     1938        pWrite(strat->P.sig);
     1939        printf("%d\n",red_result);
     1940    }
     1941#endif
    17331942
    17341943    if (strat->overflow)
    17351944    {
    17361945        if (!kStratChangeTailRing(strat)) { Werror("OVERFLOW.."); break;}
    1737     }
    1738     if (strat->incremental)
    1739     {
    1740       for (int jj = 0; jj<strat->tl+1; jj++)
    1741       {
    1742         if (pGetComp(strat->T[jj].sig) == strat->currIdx)
    1743         {
    1744           strat->T[jj].is_sigsafe = FALSE;
    1745         }
    1746       }
    1747     }
    1748     else
    1749     {
    1750       for (int jj = 0; jj<strat->tl+1; jj++)
    1751       {
    1752         strat->T[jj].is_sigsafe = FALSE;
    1753       }
    17541946    }
    17551947
     
    17851977      // in the ring case we cannot expect LC(f) = 1,
    17861978      // therefore we call pContent instead of pNorm
    1787       /*
    1788       if ((TEST_OPT_INTSTRATEGY) || (rField_is_Ring(currRing)))
    1789       {
    1790         strat->P.pCleardenom();
    1791         if ((TEST_OPT_REDSB)||(TEST_OPT_REDTAIL))
    1792         {
    1793           strat->P.p = redtailBba(&(strat->P),pos-1,strat, withT);
     1979#if SBA_TAIL_RED
     1980      if (strat->sbaOrder != 2) {
     1981        if ((TEST_OPT_INTSTRATEGY) || (rField_is_Ring(currRing)))
     1982        {
    17941983          strat->P.pCleardenom();
    1795         }
    1796       }
    1797       else
    1798       {
    1799         strat->P.pNorm();
    1800         if ((TEST_OPT_REDSB)||(TEST_OPT_REDTAIL))
    1801           strat->P.p = redtailBba(&(strat->P),pos-1,strat, withT);
    1802       }
    1803       */
     1984          if ((TEST_OPT_REDSB)||(TEST_OPT_REDTAIL))
     1985          {
     1986            strat->P.p = redtailSba(&(strat->P),pos-1,strat, withT);
     1987            strat->P.pCleardenom();
     1988          }
     1989        }
     1990        else
     1991        {
     1992          strat->P.pNorm();
     1993          if ((TEST_OPT_REDSB)||(TEST_OPT_REDTAIL))
     1994            strat->P.p = redtailSba(&(strat->P),pos-1,strat, withT);
     1995        }
     1996      }
     1997#endif
     1998
     1999    // remove sigsafe label since it is no longer valid for the next element to
     2000    // be reduced
     2001    if (strat->sbaOrder == 1)
     2002    {
     2003      for (int jj = 0; jj<strat->tl+1; jj++)
     2004      {
     2005        if (pGetComp(strat->T[jj].sig) == strat->currIdx)
     2006        {
     2007          strat->T[jj].is_sigsafe = FALSE;
     2008        }
     2009      }
     2010    }
     2011    else
     2012    {
     2013      for (int jj = 0; jj<strat->tl+1; jj++)
     2014      {
     2015        strat->T[jj].is_sigsafe = FALSE;
     2016      }
     2017    }
    18042018#ifdef KDEBUG
    18052019      if (TEST_OPT_DEBUG){PrintS("new s:");strat->P.wrp();PrintLn();}
     
    18332047      // enter into S, L, and T
    18342048      //if ((!TEST_OPT_IDLIFT) || (pGetComp(strat->P.p) <= strat->syzComp))
    1835       if(!strat->incremental)
    1836       {
    1837         BOOLEAN overwrite = TRUE;
     2049      enterT(strat->P, strat);
     2050      strat->T[strat->tl].is_sigsafe = FALSE;
     2051      /*
     2052      printf("hier\n");
     2053      pWrite(strat->P.GetLmCurrRing());
     2054      pWrite(strat->P.sig);
     2055      */
     2056#ifdef HAVE_RINGS
     2057      if (rField_is_Ring(currRing))
     2058        superenterpairs(strat->P.p,strat->sl,strat->P.ecart,pos,strat, strat->tl);
     2059      else
     2060#endif
     2061        enterpairsSig(strat->P.p,strat->P.sig,strat->sl+1,strat->sl,strat->P.ecart,pos,strat, strat->tl);
     2062      // posInS only depends on the leading term
     2063      strat->enterS(strat->P, pos, strat, strat->tl);
     2064      if(strat->sbaOrder != 1)
     2065      {
     2066        BOOLEAN overwrite = FALSE;
    18382067        for (int tk=0; tk<strat->sl+1; tk++)
    18392068        {
     
    18552084
    18562085          strat->P.sevSig = pGetShortExpVector (strat->P.sig);
     2086          int i;
     2087          LObject Q;
    18572088          for(int ps=0;ps<strat->sl+1;ps++)
    18582089          {
    1859             int i = strat->syzl;
    18602090
    18612091            strat->newt = TRUE;
     
    18692099              strat->syzmax += setmaxTinc;
    18702100            }
    1871             strat->syz[i] = pCopy(strat->P.sig);
     2101            Q.sig = pCopy(strat->P.sig);
    18722102            // add LM(F->m[i]) to the signature to get a Schreyer order
    18732103            // without changing the underlying polynomial ring at all
    1874             p_ExpVectorAdd (strat->syz[i],strat->S[ps],currRing);
     2104            if (strat->sbaOrder == 0)
     2105              p_ExpVectorAdd (Q.sig,strat->S[ps],currRing);
    18752106            // since p_Add_q() destroys all input
    18762107            // data we need to recreate help
     
    18812112            // the corresponding principal syzygy
    18822113            // => we do not need to compute the "real" syzygy completely
    1883             poly help = pCopy(strat->sig[ps]);
     2114            poly help = p_Copy(strat->sig[ps],currRing);
    18842115            p_ExpVectorAdd (help,strat->P.p,currRing);
    1885             strat->syz[i] = p_Add_q(strat->syz[i],help,currRing);
     2116            Q.sig = p_Add_q(Q.sig,help,currRing);
    18862117            //printf("%d. SYZ  ",i+1);
    18872118            //pWrite(strat->syz[i]);
    1888             strat->sevSyz[i] = p_GetShortExpVector(strat->syz[i],currRing);
    1889             strat->syzl++;
     2119            Q.sevSig = p_GetShortExpVector(Q.sig,currRing);
     2120            i = posInSyz(strat, Q.sig);
     2121            enterSyz(Q, strat, i);
    18902122          }
    18912123        }
    18922124      }
    1893       enterT(strat->P, strat);
    1894       strat->T[strat->tl].is_sigsafe = FALSE;
    1895 #ifdef HAVE_RINGS
    1896       if (rField_is_Ring(currRing))
    1897         superenterpairs(strat->P.p,strat->sl,strat->P.ecart,pos,strat, strat->tl);
    1898       else
    1899 #endif
    1900         enterpairsSig(strat->P.p,strat->P.sig,strat->sl+1,strat->sl,strat->P.ecart,pos,strat, strat->tl);
    1901       // posInS only depends on the leading term
    1902       strat->enterS(strat->P, pos, strat, strat->tl);
     2125      // deg - idx - lp/rp
     2126      // => we need to add syzygies with indices > pGetComp(strat->P.sig)
     2127      if(strat->sbaOrder == 0 || strat->sbaOrder == 3)
     2128      {
     2129        int cmp     = pGetComp(strat->P.sig);
     2130        int max_cmp = IDELEMS(F);
     2131        int* vv = (int*)omAlloc((currRing->N+1)*sizeof(int));
     2132        pGetExpV (strat->P.p,vv);
     2133        LObject Q;
     2134        int pos;
     2135        int idx = p_GetComp(strat->P.sig,currRing);
     2136        //printf("++ -- adding syzygies -- ++\n");
     2137        // if new element is the first one in this index
     2138        if (strat->currIdx < idx) {
     2139          for (int i=0; i<strat->sl; ++i) {
     2140            Q.sig = p_Copy(strat->P.sig,currRing);
     2141            p_ExpVectorAdd(Q.sig,strat->S[i],currRing);
     2142            poly help = p_Copy(strat->sig[i],currRing);
     2143            p_ExpVectorAdd(help,strat->P.p,currRing);
     2144            Q.sig = p_Add_q(Q.sig,help,currRing);
     2145            //pWrite(Q.sig);
     2146            pos = posInSyz(strat, Q.sig);
     2147            enterSyz(Q, strat, pos);
     2148          }
     2149          strat->currIdx = idx;
     2150        } else {
     2151          // if the element is not the first one in the given index we build all
     2152          // possible syzygies with elements of higher index
     2153          for (int i=cmp+1; i<=max_cmp; ++i) {
     2154            pos = -1;
     2155            for (int j=0; j<strat->sl; ++j) {
     2156              if (p_GetComp(strat->sig[j],currRing) == i) {
     2157                pos = j;
     2158                break;
     2159              }
     2160            }
     2161            if (pos != -1) {
     2162              Q.sig = p_One(currRing);
     2163              p_SetExpV(Q.sig, vv, currRing);
     2164              // F->m[i-1] corresponds to index i
     2165              p_ExpVectorAdd(Q.sig,F->m[i-1],currRing);
     2166              p_SetComp(Q.sig, i, currRing);
     2167              poly help = p_Copy(strat->P.sig,currRing);
     2168              p_ExpVectorAdd(help,strat->S[pos],currRing);
     2169              Q.sig = p_Add_q(Q.sig,help,currRing);
     2170              if (strat->sbaOrder == 0) {
     2171                if (p_LmCmp(Q.sig,strat->syz[strat->syzl-1],currRing) == -currRing->OrdSgn) {
     2172                  pos = posInSyz(strat, Q.sig);
     2173                  enterSyz(Q, strat, pos);
     2174                }
     2175              } else {
     2176                pos = posInSyz(strat, Q.sig);
     2177                enterSyz(Q, strat, pos);
     2178              }
     2179            }
     2180          }
     2181          //printf("++ -- done adding syzygies -- ++\n");
     2182        }
     2183      }
    19032184//#if 1
    19042185#if DEBUGF50
     
    19452226      // pair was not detected by the rewritten criterion in strat->red = redSig
    19462227      if (red_result!=2) {
     2228#if SBA_PRINT_ZERO_REDUCTIONS
    19472229        zeroreductions++;
    1948         enterSyz(strat->P,strat);
     2230#endif
     2231        int pos = posInSyz(strat, strat->P.sig);
     2232        enterSyz(strat->P, strat, pos);
    19492233//#if 1
    19502234#ifdef DEBUGF5
    19512235        Print("ADDING STUFF TO SYZ :  ");
    1952         pWrite(strat->P.p);
     2236        //pWrite(strat->P.p);
    19532237        pWrite(strat->P.sig);
    19542238#endif
     
    20122296#endif
    20132297#if SBA_PRINT_SIZE_SYZ
    2014   size_syz = strat->syzl+1;
    2015 #endif
    2016 
     2298  // that is correct, syzl is counting one too far
     2299  size_syz = strat->syzl;
     2300#endif
    20172301  exitSba(strat);
    20182302//  if (TEST_OPT_WEIGHTM)
     
    20432327  }
    20442328#endif
    2045   if (strat->incremental && sRing!=currRingOld)
     2329  if ((strat->sbaOrder == 1 || strat->sbaOrder == 3) && sRing!=currRingOld)
    20462330  {
    20472331    rChangeCurrRing (currRingOld);
    2048     F0          = idrMoveR (F, sRing, currRing);
     2332    F0          = idrMoveR (F1, sRing, currRing);
    20492333    strat->Shdl = idrMoveR_NoSort (strat->Shdl, sRing, currRing);
    20502334    rDelete (sRing);
     
    20632347#endif
    20642348#if SBA_PRINT_ZERO_REDUCTIONS
    2065   printf("ZERO REDUCTIONS:   %ld\n",zeroreductions);
     2349  printf("----------------------------------------------------------\n");
     2350  printf("ZERO REDUCTIONS:            %ld\n",zeroreductions);
     2351  zeroreductions  = 0;
    20662352#endif
    20672353#if SBA_PRINT_REDUCTION_STEPS
    2068   printf("TOP S-REDUCTIONS:  %ld\n",sba_reduction_steps);
     2354  printf("----------------------------------------------------------\n");
     2355  printf("S-REDUCTIONS:               %ld\n",sba_reduction_steps);
     2356#endif
     2357#if SBA_PRINT_OPERATIONS
     2358  printf("OPERATIONS:                 %ld\n",sba_operations);
     2359#endif
     2360#if SBA_PRINT_REDUCTION_STEPS
     2361  printf("- - - - - - - - - - - - - - - - - - - - - - - - - - - - - \n");
     2362  printf("INTERREDUCTIONS:            %ld\n",sba_interreduction_steps);
     2363#endif
     2364#if SBA_PRINT_OPERATIONS
     2365  printf("INTERREDUCTION OPERATIONS:  %ld\n",sba_interreduction_operations);
     2366#endif
     2367#if SBA_PRINT_REDUCTION_STEPS
     2368  printf("- - - - - - - - - - - - - - - - - - - - - - - - - - - - - \n");
     2369  printf("ALL REDUCTIONS:             %ld\n",sba_reduction_steps+sba_interreduction_steps);
     2370  sba_interreduction_steps  = 0;
     2371  sba_reduction_steps       = 0;
     2372#endif
     2373#if SBA_PRINT_OPERATIONS
     2374  printf("ALL OPERATIONS:             %ld\n",sba_operations+sba_interreduction_operations);
     2375  sba_interreduction_operations = 0;
     2376  sba_operations                = 0;
    20692377#endif
    20702378#if SBA_PRINT_SIZE_G
    2071   printf("SIZE OF G:         %ld\n",size_g);
     2379  printf("----------------------------------------------------------\n");
     2380  printf("SIZE OF G:                  %ld\n",size_g);
     2381  size_g  = 0;
    20722382#endif
    20732383#if SBA_PRINT_SIZE_SYZ
    2074   printf("SIZE OF SYZ:       %ld\n",size_syz);
     2384  printf("SIZE OF SYZ:                %ld\n",size_syz);
     2385  printf("----------------------------------------------------------\n");
     2386  size_syz  = 0;
    20752387#endif
    20762388#if SBA_PRINT_PRODUCT_CRITERION
    2077   printf("PRODUCT CRITERIA:  %ld\n",product_criterion);
    2078 #endif
    2079   zeroreductions      = 0;
    2080   size_g              = 0;
    2081   size_syz            = 0;
    2082   product_criterion   = 0;
    2083   sba_reduction_steps = 0;
     2389  printf("PRODUCT CRITERIA:           %ld\n",product_criterion);
     2390  product_criterion = 0;
     2391#endif
    20842392  return (strat->Shdl);
    20852393}
     
    24042712      // therefore we call pContent instead of pNorm
    24052713#if F5CTAILRED
     2714      BOOLEAN withT = TRUE;
    24062715      if ((TEST_OPT_INTSTRATEGY) || (rField_is_Ring(currRing)))
    24072716      {
  • kernel/kutil.cc

    r74c4462 r6e3023a  
    987987  memmove(&(strat->sevSig[i]),&(strat->sevSig[i+1]),(strat->sl - i)*sizeof(unsigned long));
    988988  memmove(&(strat->S_2_R[i]),&(strat->S_2_R[i+1]),(strat->sl - i)*sizeof(int));
    989   memmove(&(strat->fromS[i]),&(strat->fromS[i+1]),(strat->sl - i)*sizeof(int));
    990989#else
    991990  int j;
     
    998997    strat->sevSig[j] = strat->sevSig[j+1];
    999998    strat->S_2_R[j] = strat->S_2_R[j+1];
    1000     strat->fromS[j] = strat->fromS[j+1];
    1001999  }
    10021000#endif
     
    10451043#endif
    10461044      pLmFree(set[j].lcm);
     1045  }
     1046  if (set[j].sig!=NULL)
     1047  {
     1048#ifdef HAVE_RINGS
     1049    if (pGetCoeff(set[j].sig) != NULL)
     1050      pLmDelete(set[j].sig);
     1051    else
     1052#endif
     1053      pLmFree(set[j].sig);
    10471054  }
    10481055  if (set[j].p!=NULL)
     
    17791786              // the corresponding signatures for criteria checks
    17801787  LObject  Lp;
    1781   // poly last;
    17821788  poly pSigMult = p_Copy(pSig,currRing);
    17831789  poly sSigMult = p_Copy(strat->sig[i],currRing);
     
    18061812  Print("P1  ");
    18071813  pWrite(pHead(p));
    1808   Print("FROM: %d\n", from);
    18091814  Print("P2  ");
    18101815  pWrite(pHead(strat->S[i]));
    1811   Print("FROM: %d\n", strat->fromS[i]);
    18121816  Print("M1  ");
    18131817  pWrite(m1);
     
    18271831  pWrite(sSigMult);
    18281832  Print("----------------\n");
    1829 #endif
    1830   // testing by syzCrit = F5 Criterion
    1831   // testing by rewCrit1 = Rewritten Criterion
    1832   if  ( strat->syzCrit(pSigMult,pSigMultNegSev,strat) ||
    1833         strat->syzCrit(sSigMult,sSigMultNegSev,strat)
    1834         || strat->rewCrit1(sSigMult,sSigMultNegSev,strat,i+1)
    1835       )
    1836   {
    1837     pDelete(&pSigMult);
    1838     pDelete(&sSigMult);
    1839     strat->cp++;
    1840     pLmFree(Lp.lcm);
    1841     Lp.lcm=NULL;
    1842     pDelete (&m1);
    1843     pDelete (&m2);
    1844     return;
    1845   }
    1846   // in any case Lp is checked up to the next strat->P which is added
    1847   // to S right after this critical pair creation.
    1848   // NOTE: this even holds if the 2nd generator gives the bigger signature
    1849   //       moreover, this improves rewCriterion,
    1850   //       i.e. strat->checked > strat->from if and only if the 2nd generator
    1851   //       gives the bigger signature.
    1852   Lp.checked = strat->sl+1;
     1833  Lp.checked  = 0;
     1834#endif
    18531835  int sigCmp = p_LmCmp(pSigMult,sSigMult,currRing);
    18541836//#if 1
     
    18621844    // printf("!!!!   EQUAL SIGS   !!!!\n");
    18631845    // pSig = sSig, delete element due to Rewritten Criterion
    1864     strat->cp++;
    18651846    pDelete(&pSigMult);
    18661847    pDelete(&sSigMult);
     
    18711852    return;
    18721853  }
    1873   // at this point it is clear that the pair will be added to L, since it has
    1874   // passed all tests up to now
    1875 
    1876   // store from which element this pair comes from for further tests
    1877   Lp.from = strat->sl+1;
    1878   if(sigCmp==currRing->OrdSgn)
    1879   {
    1880     // pSig > sSig
    1881     pDelete (&sSigMult);
    1882     Lp.sig    = pSigMult;
    1883     Lp.sevSig = ~pSigMultNegSev;
    1884   }
    1885   else
    1886   {
    1887     // pSig < sSig
    1888     pDelete (&pSigMult);
    1889     Lp.sig    = sSigMult;
    1890     Lp.sevSig = ~sSigMultNegSev;
    1891   }
    1892 // adds buchberger's first criterion
    1893   if (pLmCmp(m2,pHead(p)) == 0) {
    1894     Lp.checked  = 3; // 3 == Product Criterion
    1895 #if 0
    1896     enterSyz(Lp, strat);
     1854  // testing by syzCrit = F5 Criterion
     1855  // testing by rewCrit1 = Rewritten Criterion
     1856  // NOTE: Arri's Rewritten Criterion is tested below, we need Lp.p for it!
     1857  if  ( strat->syzCrit(pSigMult,pSigMultNegSev,strat) ||
     1858        strat->syzCrit(sSigMult,sSigMultNegSev,strat)
     1859        || strat->rewCrit1(sSigMult,sSigMultNegSev,Lp.lcm,strat,i+1)
     1860      )
     1861  {
     1862    pDelete(&pSigMult);
     1863    pDelete(&sSigMult);
     1864    pLmFree(Lp.lcm);
    18971865    Lp.lcm=NULL;
    18981866    pDelete (&m1);
    18991867    pDelete (&m2);
    19001868    return;
    1901 #endif
    1902   }
    1903   pDelete (&m1);
    1904   pDelete (&m2);
    1905 #if DEBUGF5
    1906   printf("SIGNATURE OF PAIR:  ");
    1907   pWrite(Lp.sig);
    1908 #endif
     1869  }
    19091870  /*
    19101871  *the pair (S[i],p) enters B if the spoly != 0
     
    19851946      }
    19861947  }
     1948  // store from which element this pair comes from for further tests
     1949  //Lp.from = strat->sl+1;
     1950  if(sigCmp==currRing->OrdSgn)
     1951  {
     1952    // pSig > sSig
     1953    pDelete (&sSigMult);
     1954    Lp.sig    = pSigMult;
     1955    Lp.sevSig = ~pSigMultNegSev;
     1956  }
     1957  else
     1958  {
     1959    // pSig < sSig
     1960    pDelete (&pSigMult);
     1961    Lp.sig    = sSigMult;
     1962    Lp.sevSig = ~sSigMultNegSev;
     1963  }
    19871964  if (Lp.p == NULL)
    19881965  {
    1989     /*- the case that the s-poly is 0 -*/
    1990     if (strat->pairtest==NULL) initPairtest(strat);
    1991     strat->pairtest[i] = TRUE;/*- hint for spoly(S^[i],p)=0 -*/
    1992     strat->pairtest[strat->sl+1] = TRUE;
    1993     /*hint for spoly(S[i],p) == 0 for some i,0 <= i <= sl*/
    1994     /*
    1995     *suppose we have (s,r),(r,p),(s,p) and spoly(s,p) == 0 and (r,p) is
    1996     *still in B (i.e. lcm(r,p) == lcm(s,p) or the leading term of s does not
    1997     *devide lcm(r,p)). In the last case (s,r) can be canceled if the leading
    1998     *term of p devides the lcm(s,r)
    1999     *(this canceling should be done here because
    2000     *the case lcm(s,p) == lcm(s,r) is not covered in chainCrit)
    2001     *the first case is handeled in chainCrit
    2002     */
    20031966    if (Lp.lcm!=NULL) pLmFree(Lp.lcm);
     1967    int pos = posInSyz(strat, Lp.sig);
     1968    enterSyz(Lp, strat, pos);
    20041969  }
    20051970  else
    20061971  {
     1972    // testing by rewCrit3 = Arris Rewritten Criterion (for F5 nothing happens!)
     1973    if (strat->rewCrit3(Lp.sig,~Lp.sevSig,Lp.p,strat,strat->sl+1)) {
     1974      pLmFree(Lp.lcm);
     1975      pDelete(&Lp.sig);
     1976      Lp.lcm=NULL;
     1977      pDelete (&m1);
     1978      pDelete (&m2);
     1979      return;
     1980    }
     1981    // in any case Lp is checked up to the next strat->P which is added
     1982    // to S right after this critical pair creation.
     1983    // NOTE: this even holds if the 2nd generator gives the bigger signature
     1984    //       moreover, this improves rewCriterion,
     1985    //       i.e. strat->checked > strat->from if and only if the 2nd generator
     1986    //       gives the bigger signature.
     1987    Lp.checked = strat->sl+1;
     1988    // at this point it is clear that the pair will be added to L, since it has
     1989    // passed all tests up to now
     1990
     1991  // adds buchberger's first criterion
     1992    if (pLmCmp(m2,pHead(p)) == 0) {
     1993      Lp.prod_crit = TRUE; // Product Criterion
     1994#if 0
     1995      int pos = posInSyz(strat, Lp.sig);
     1996      enterSyz(Lp, strat, pos);
     1997      Lp.lcm=NULL;
     1998      pDelete (&m1);
     1999      pDelete (&m2);
     2000      return;
     2001#endif
     2002    }
     2003    pDelete (&m1);
     2004    pDelete (&m2);
     2005#if DEBUGF5
     2006    printf("SIGNATURE OF PAIR:  ");
     2007    pWrite(Lp.sig);
     2008#endif
    20072009    /*- the pair (S[i],p) enters B -*/
    20082010    Lp.p1 = strat->S[i];
     
    46314633}
    46324634
     4635// for sba, sorting syzygies
     4636int posInSyz (const kStrategy strat, poly sig)
     4637{
     4638if (strat->syzl==0) return 0;
     4639if (pLmCmp(strat->syz[strat->syzl-1],sig) != currRing->OrdSgn)
     4640  return strat->syzl;
     4641int i;
     4642int an = 0;
     4643int en= strat->syzl-1;
     4644loop
     4645{
     4646  if (an >= en-1)
     4647  {
     4648    if (pLmCmp(strat->syz[an],sig) != currRing->OrdSgn) return en;
     4649    return an;
     4650  }
     4651  i=(an+en) / 2;
     4652  if (pLmCmp(strat->syz[i],sig) != currRing->OrdSgn) an=i;
     4653  else                                      en=i;
     4654  /*aend. fuer lazy == in !=- machen */
     4655}
     4656}
     4657
    46334658/*2
    46344659*
     
    50655090  for (int k=0; k<strat->syzl; k++)
    50665091  {
     5092    //printf("-%d",k);
    50675093//#if 1
    50685094#ifdef DEBUGF5
    5069     Print("checking with: %d --  ",k);
     5095    Print("checking with: %d / %d --  \n",k,strat->syzl);
    50705096    pWrite(pHead(strat->syz[k]));
    50715097#endif
     
    50765102      printf("DELETE!\n");
    50775103#endif
     5104      //printf("- T -\n\n");
    50785105      return TRUE;
    50795106    }
    50805107  }
     5108  //printf("- F -\n\n");
    50815109  return FALSE;
    50825110}
     
    50895117//#if 1
    50905118#ifdef DEBUGF5
    5091   Print("syzygy criterion checks:  ");
     5119  Print("--- syzygy criterion checks:  ");
    50925120  pWrite(sig);
    50935121#endif
     
    51125140    for (int k=min; k<max; k++)
    51135141    {
    5114 #ifdef DEBUGF5
     5142#ifdef F5DEBUG
    51155143      printf("COMP %d/%d - MIN %d - MAX %d - SYZL %ld\n",comp,strat->currIdx,min,max,strat->syzl);
    51165144      Print("checking with: %d --  ",k);
     
    51275155 * REWRITTEN CRITERION for signature-based standard basis algorithms
    51285156 */
    5129 BOOLEAN faugereRewCriterion(poly sig, unsigned long not_sevSig, kStrategy strat, int start=0)
     5157BOOLEAN faugereRewCriterion(poly sig, unsigned long not_sevSig, poly /*lm*/, kStrategy strat, int start=0)
    51305158{
    51315159  //printf("Faugere Rewritten Criterion\n");
     
    51355163  pWrite(sig);
    51365164#endif
    5137   //for(int k = start; k<strat->sl+1; k++)
    51385165  for(int k = strat->sl; k>start; k--)
    51395166  {
     
    51455172#endif
    51465173    if (p_LmShortDivisibleBy(strat->sig[k], strat->sevSig[k], sig, not_sevSig, currRing))
    5147     //if (p_LmEqual(strat->sig[k], sig, currRing))
    51485174    {
    51495175//#if 1
     
    51535179      return TRUE;
    51545180    }
     5181    //k--;
    51555182  }
    51565183#ifdef DEBUGF5
     
    51845211//        critical pair. In this situation we can discard the critical pair
    51855212//        completely.
    5186 BOOLEAN arriRewCriterion(poly /*sig*/, unsigned long /*not_sevSig*/, kStrategy strat, int /*start=0*/)
    5187 {
    5188   //printf("Arri Rewritten Criterion\n");
    5189   while (strat->Ll > 0 && pLmEqual(strat->L[strat->Ll].sig,strat->P.sig))
    5190   {
    5191     // deletes the short spoly
    5192 #ifdef HAVE_RINGS
    5193     if (rField_is_Ring(currRing))
    5194       pLmDelete(strat->L[strat->Ll].p);
    5195     else
    5196 #endif
    5197       pLmFree(strat->L[strat->Ll].p);
    5198 
    5199     // TODO: needs some masking
    5200     // TODO: masking needs to vanish once the signature
    5201     //       sutff is completely implemented
    5202     strat->L[strat->Ll].p = NULL;
    5203     poly m1 = NULL, m2 = NULL;
    5204 
    5205     // check that spoly creation is ok
    5206     while (strat->tailRing != currRing &&
    5207           !kCheckSpolyCreation(&(strat->L[strat->Ll]), strat, m1, m2))
    5208     {
    5209       assume(m1 == NULL && m2 == NULL);
    5210       // if not, change to a ring where exponents are at least
    5211       // large enough
    5212       if (!kStratChangeTailRing(strat))
    5213       {
    5214         WerrorS("OVERFLOW...");
    5215         break;
    5216       }
    5217     }
    5218     // create the real one
    5219     ksCreateSpoly(&(strat->L[strat->Ll]), NULL, strat->use_buckets,
    5220                   strat->tailRing, m1, m2, strat->R);
    5221     if (strat->P.GetLmCurrRing() == NULL)
    5222     {
    5223       deleteInL(strat->L,&strat->Ll,strat->Ll,strat);
    5224     }
    5225     if (strat->L[strat->Ll].GetLmCurrRing() == NULL)
    5226     {
    5227       strat->P.Delete();
    5228       strat->P = strat->L[strat->Ll];
    5229       strat->Ll--;
    5230     }
    5231 
    5232     if ((strat->P.GetLmCurrRing() != NULL)
    5233     && (strat->L[strat->Ll].GetLmCurrRing() != NULL))
    5234     {
    5235       if (pLmCmp(strat->P.GetLmCurrRing(),strat->L[strat->Ll].GetLmCurrRing()) == -1)
    5236       {
    5237         deleteInL(strat->L,&strat->Ll,strat->Ll,strat);
    5238       }
    5239       else
    5240       {
    5241         strat->P.Delete();
    5242         strat->P = strat->L[strat->Ll];
    5243         strat->Ll--;
    5244       }
    5245     }
    5246   }
     5213BOOLEAN arriRewCriterion(poly /*sig*/, unsigned long /*not_sevSig*/, poly /*lm*/, kStrategy strat, int start=0)
     5214{
     5215  poly p1 = pOne();
     5216  poly p2 = pOne();
     5217  for (int ii=strat->sl; ii>start; ii--)
     5218  {
     5219    if (p_LmShortDivisibleBy(strat->sig[ii], strat->sevSig[ii], strat->P.sig, ~strat->P.sevSig, currRing))
     5220    {
     5221      p_ExpVectorSum(p1,strat->P.sig,strat->S[ii],currRing);
     5222      p_ExpVectorSum(p2,strat->sig[ii],strat->P.p,currRing);
     5223      if (!(pLmCmp(p1,p2) == 1))
     5224      {
     5225        pDelete(&p1);
     5226        pDelete(&p2);
     5227        return TRUE;
     5228      }
     5229    }
     5230  }
     5231  pDelete(&p1);
     5232  pDelete(&p2);
     5233  return FALSE;
     5234}
     5235
     5236BOOLEAN arriRewCriterionPre(poly sig, unsigned long not_sevSig, poly lm, kStrategy strat, int /*start=0*/)
     5237{
     5238  int found = -1;
     5239  for (int i=strat->Bl; i>-1; i--) {
     5240    if (pLmEqual(strat->B[i].sig,sig)) {
     5241      found = i;
     5242      break;
     5243    }
     5244  }
     5245  if (found != -1) {
     5246    if (pLmCmp(lm,strat->B[found].GetLmCurrRing()) == -1) {
     5247      deleteInL(strat->B,&strat->Bl,found,strat);
     5248    } else {
     5249      return TRUE;
     5250    }
     5251  }
     5252  poly p1 = pOne();
     5253  poly p2 = pOne();
    52475254  for (int ii=strat->sl; ii>-1; ii--)
    52485255  {
    5249     if (p_LmShortDivisibleBy(strat->sig[ii], strat->sevSig[ii], strat->P.sig, ~strat->P.sevSig, currRing))
    5250     {
    5251       if (!(pLmCmp(ppMult_mm(strat->P.sig,pHead(strat->S[ii])),ppMult_mm(strat->sig[ii],strat->P.GetLmCurrRing())) == 1))
    5252       {
    5253         strat->P.Delete();
     5256    if (p_LmShortDivisibleBy(strat->sig[ii], strat->sevSig[ii], sig, not_sevSig, currRing))
     5257    {
     5258      p_ExpVectorSum(p1,sig,strat->S[ii],currRing);
     5259      p_ExpVectorSum(p2,strat->sig[ii],lm,currRing);
     5260      if (!(pLmCmp(p1,p2) == 1))
     5261      {
     5262        pDelete(&p1);
     5263        pDelete(&p2);
    52545264        return TRUE;
    52555265      }
    52565266    }
    52575267  }
     5268  pDelete(&p1);
     5269  pDelete(&p2);
    52585270  return FALSE;
    52595271}
     
    59295941  else i=setmaxT;
    59305942  strat->ecartS =   initec(i);
    5931   strat->fromS  =   initec(i);
    59325943  strat->sevS   =   initsevS(i);
    59335944  strat->sevSig =   initsevS(i);
     
    59375948  strat->S      =   strat->Shdl->m;
    59385949  strat->sig    =   (poly *)omAlloc0(i*sizeof(poly));
    5939   if (!strat->incremental)
     5950  if (strat->sbaOrder != 1)
    59405951  {
    59415952    strat->syz    = (poly *)omAlloc0(i*sizeof(poly));
     
    59996010      // => we can keep the underlying monomial order and get a Schreyer
    60006011      //    order without any bigger overhead
    6001       if (!strat->incremental)
     6012      if (strat->sbaOrder == 0 || strat->sbaOrder == 3)
    60026013      {
    60036014        p_ExpVectorAdd (h.sig,F->m[i],currRing);
     
    60366047      }
    60376048      /*
    6038       if (!strat->incremental)
     6049      if (strat->sbaOrder != 1)
    60396050      {
    60406051        for(j=0;j<i;j++)
     
    61036114    strat->sevSyz     = initsevS(ps);
    61046115    strat->syz        = (poly *)omAlloc(ps*sizeof(poly));
    6105     strat->syzl       = strat->syzmax = ps;
     6116    strat->syzmax     = ps;
     6117    strat->syzl       = 0;
    61066118    strat->syzidxmax  = comp;
    61076119#if defined(DEBUGF5) || defined(DEBUGF51)
     
    61416153        strat->syzIdx[j]  = ctr;
    61426154        j++;
     6155        LObject Q;
     6156        int pos;
    61436157        for (k = 0; k<i; k++)
    61446158        {
    6145           poly p          = pOne();
    6146           pLcm(strat->S[k],strat->S[i],p);
    6147           strat->syz[ctr] = p;
    6148           p_SetCompP (strat->syz[ctr], comp, currRing);
    6149           poly q          = p_Copy(p, currRing);
     6159          Q.sig          = pOne();
     6160          p_ExpVectorCopy(Q.sig,strat->S[k],currRing);
     6161          p_SetCompP (Q.sig, comp, currRing);
     6162          poly q          = p_One(currRing);
     6163          p_ExpVectorCopy(q,strat->S[i],currRing);
    61506164          q               = p_Neg (q, currRing);
    61516165          p_SetCompP (q, p_GetComp(strat->sig[k], currRing), currRing);
    6152           strat->syz[ctr] = p_Add_q (strat->syz[ctr], q, currRing);
    6153 #if defined(DEBUGF5) || defined(DEBUGF51)
    6154           pWrite(strat->syz[ctr]);
    6155 #endif
    6156           strat->sevSyz[ctr] = p_GetShortExpVector(strat->syz[ctr],currRing);
     6166          Q.sig = p_Add_q (Q.sig, q, currRing);
     6167          Q.sevSig  = p_GetShortExpVector(Q.sig,currRing);
     6168          pos = posInSyz(strat, Q.sig);
     6169          enterSyz(Q, strat, pos);
    61576170          ctr++;
    61586171        }
     
    61806193    }
    61816194    strat->syzIdx[j]  = ctr;
     6195    LObject Q;
     6196    int pos;
    61826197    for (k = 0; k<strat->sl+1; k++)
    61836198    {
    6184       strat->syz[ctr] = p_Copy (pHead(strat->S[k]), currRing);
    6185       p_SetCompP (strat->syz[ctr], comp, currRing);
    6186       poly q          = p_Copy (pHead(strat->L[strat->Ll].p), currRing);
     6199      Q.sig          = pOne();
     6200      p_ExpVectorCopy(Q.sig,strat->S[k],currRing);
     6201      p_SetCompP (Q.sig, comp, currRing);
     6202      poly q          = p_One(currRing);
     6203      p_ExpVectorCopy(q,strat->L[strat->Ll].p,currRing);
    61876204      q               = p_Neg (q, currRing);
    61886205      p_SetCompP (q, p_GetComp(strat->sig[k], currRing), currRing);
    6189       strat->syz[ctr] = p_Add_q (strat->syz[ctr], q, currRing);
    6190 //#if 1
    6191 #if DEBUGF5 || DEBUGF51
    6192       printf("..");
    6193       pWrite(strat->syz[ctr]);
    6194 #endif
    6195       strat->sevSyz[ctr] = p_GetShortExpVector(strat->syz[ctr],currRing);
     6206      Q.sig = p_Add_q (Q.sig, q, currRing);
     6207      Q.sevSig = p_GetShortExpVector(Q.sig,currRing);
     6208      pos = posInSyz(strat, Q.sig);
     6209      enterSyz(Q, strat, pos);
    61966210      ctr++;
    61976211    }
     
    61996213#ifdef DEBUGF5
    62006214    Print("Principal syzygies:\n");
     6215    printf("syzl   %d\n",strat->syzl);
     6216    printf("syzmax %d\n",strat->syzmax);
     6217    printf("ps     %d\n",ps);
    62016218    Print("--------------------------------\n");
    6202     for(i=0;i<=ps-1;i++)
    6203     {
     6219    for(i=0;i<=strat->syzl-1;i++)
     6220    {
     6221      printf("%d - ",i);
    62046222      pWrite(strat->syz[i]);
     6223    }
     6224    for(i=0;i<strat->currIdx;i++)
     6225    {
     6226      printf("%d - %d\n",i,strat->syzIdx[i]);
    62056227    }
    62066228    Print("--------------------------------\n");
     
    63676389  else i=setmaxT;
    63686390  i=((i+IDELEMS(F)+IDELEMS(P)+15)/16)*16;
    6369   strat->fromS=initec(i);
    63706391  strat->sevS=initsevS(i);
    63716392  strat->sevSig=initsevS(i);
     
    70247045                                          (IDELEMS(strat->Shdl)+setmaxTinc)
    70257046                                                  *sizeof(int));
    7026     strat->fromS = (intset)omReallocSize(strat->fromS,
    7027                                           IDELEMS(strat->Shdl)*sizeof(int),
    7028                                           (IDELEMS(strat->Shdl)+setmaxTinc)
    7029                                                   *sizeof(int));
    70307047    strat->S_2_R = (int*) omRealloc0Size(strat->S_2_R,
    70317048                                         IDELEMS(strat->Shdl)*sizeof(int),
     
    70647081    memmove(&(strat->ecartS[atS+1]), &(strat->ecartS[atS]),
    70657082            (strat->sl - atS + 1)*sizeof(int));
    7066     memmove(&(strat->fromS[atS+1]), &(strat->fromS[atS]),
    7067             (strat->sl - atS + 1)*sizeof(int));
    70687083    memmove(&(strat->sevS[atS+1]), &(strat->sevS[atS]),
    70697084            (strat->sl - atS + 1)*sizeof(unsigned long));
     
    70817096      strat->S[i] = strat->S[i-1];
    70827097      strat->ecartS[i] = strat->ecartS[i-1];
    7083       strat->fromS[i] = strat->fromS[i-1];
    70847098      strat->sevS[i] = strat->sevS[i-1];
    70857099      strat->S_2_R[i] = strat->S_2_R[i-1];
     
    71287142  }
    71297143  strat->ecartS[atS] = p.ecart;
    7130   strat->fromS[atS] = p.from;
    71317144  strat->S_2_R[atS] = atR;
    71327145  strat->sl++;
     
    72147227}
    72157228
     7229
    72167230/*2
    72177231* puts signature p.sig to the set syz
    72187232*/
    7219 void enterSyz(LObject p, kStrategy strat)
    7220 {
    7221   int i = strat->syzl;
    7222 
     7233void enterSyz(LObject p, kStrategy strat, int atT)
     7234{
     7235  int i;
    72237236  strat->newt = TRUE;
    7224   if (strat->syzl == strat->syzmax)
     7237  if (strat->syzl == strat->syzmax-1)
    72257238  {
    72267239    pEnlargeSet(&strat->syz,strat->syzmax,setmaxTinc);
     
    72317244    strat->syzmax += setmaxTinc;
    72327245  }
    7233   strat->syz[i] = p.sig;
    7234   strat->sevSyz[i] = p.sevSig;
     7246  if (atT < strat->syzl)
     7247  {
     7248#ifdef ENTER_USE_MEMMOVE
     7249    memmove(&(strat->syz[atT+1]), &(strat->syz[atT]),
     7250            (strat->syzl-atT+1)*sizeof(poly));
     7251    memmove(&(strat->sevSyz[atT+1]), &(strat->sevSyz[atT]),
     7252            (strat->syzl-atT+1)*sizeof(unsigned long));
     7253#endif
     7254    for (i=strat->syzl; i>=atT+1; i--)
     7255    {
     7256#ifndef ENTER_USE_MEMMOVE
     7257      strat->syz[i] = strat->syz[i-1];
     7258      strat->sevSyz[i] = strat->sevSyz[i-1];
     7259#endif
     7260    }
     7261  }
     7262  //i = strat->syzl;
     7263  i = atT;
     7264  strat->syz[atT] = p.sig;
     7265  strat->sevSyz[atT] = p.sevSig;
    72357266  strat->syzl++;
    7236 #ifdef DEBUGF5
    7237   Print("last element in strat->syz: %d--%d  ",i+1,strat->syzmax);
    7238   pWrite(strat->syz[i]);
     7267#if F5DEBUG
     7268  Print("element in strat->syz: %d--%d  ",atT+1,strat->syzmax);
     7269  pWrite(strat->syz[atT]);
    72397270#endif
    72407271  // recheck pairs in strat->L with new rule and delete correspondingly
     
    72427273  while (cc>-1)
    72437274  {
    7244     if (p_LmShortDivisibleBy( strat->syz[strat->syzl-1], strat->sevSyz[strat->syzl-1],
     7275    if (p_LmShortDivisibleBy( strat->syz[atT], strat->sevSyz[atT],
    72457276                              strat->L[cc].sig, ~strat->L[cc].sevSig, currRing))
    72467277    {
     
    72497280    cc--;
    72507281  }
    7251 
     7282//#if 1
     7283#ifdef DEBUGF5
     7284    Print("--- Syzygies ---\n");
     7285    printf("syzl   %d\n",strat->syzl);
     7286    printf("syzmax %d\n",strat->syzmax);
     7287    Print("--------------------------------\n");
     7288    for(i=0;i<=strat->syzl-1;i++)
     7289    {
     7290      printf("%d - ",i);
     7291      pWrite(strat->syz[i]);
     7292    }
     7293    Print("--------------------------------\n");
     7294#endif
    72527295}
    72537296
     
    73417384   *****************************************/
    73427385  //strat->rewCrit1     = faugereRewCriterion;
    7343   if (strat->incremental)
     7386  if (strat->sbaOrder == 1)
    73447387  {
    73457388    strat->syzCrit  = syzCriterionInc;
     
    77607803  omFreeSize(strat->sevT, (strat->tmax)*sizeof(unsigned long));
    77617804  omFreeSize(strat->ecartS,IDELEMS(strat->Shdl)*sizeof(int));
    7762   omFreeSize(strat->fromS,IDELEMS(strat->Shdl)*sizeof(int));
    77637805  omFreeSize((ADDRESS)strat->sevS,IDELEMS(strat->Shdl)*sizeof(unsigned long));
    77647806  omFreeSize((ADDRESS)strat->sevSig,IDELEMS(strat->Shdl)*sizeof(unsigned long));
    77657807  omFreeSize((ADDRESS)strat->syz,(strat->syzmax)*sizeof(poly));
    77667808  omFreeSize((ADDRESS)strat->sevSyz,(strat->syzmax)*sizeof(unsigned long));
    7767   if (strat->incremental)
     7809  if (strat->sbaOrder == 1)
    77687810  {
    77697811    omFreeSize(strat->syzIdx,(strat->syzidxmax)*sizeof(int));
     
    82548296{
    82558297  int n = rBlocks(r); // Including trailing zero!
    8256   // if incremental => use (C,monomial order from r)
    8257   if (strat->incremental)
     8298  // if sbaOrder == 1 => use (C,monomial order from r)
     8299  if (strat->sbaOrder == 1)
    82588300  {
    82598301    if (r->order[0] == ringorder_C || r->order[0] == ringorder_c)
     
    82618303      return r;
    82628304    }
    8263     ring res = rCopy0(r, FALSE, TRUE);
    8264     for (int i=1; i<n-1; i++)
    8265     {
    8266       res->order[i] = res->order[i-1];
    8267       res->block0[i] = res->block0[i-1];
    8268       res->block1[i] = res->block1[i-1];
    8269       res->wvhdl[i] = res->wvhdl[i-1];
     8305    ring res = rCopy0(r, TRUE, FALSE);
     8306    res->order  = (int *)omAlloc0((n+1)*sizeof(int));
     8307    res->block0 = (int *)omAlloc0((n+1)*sizeof(int));
     8308    res->block1 = (int *)omAlloc0((n+1)*sizeof(int));
     8309    int **wvhdl = (int **)omAlloc0((n+1)*sizeof(int*));
     8310    res->wvhdl  = wvhdl;
     8311    for (int i=1; i<n; i++)
     8312    {
     8313      res->order[i]   = r->order[i-1];
     8314      res->block0[i]  = r->block0[i-1];
     8315      res->block1[i]  = r->block1[i-1];
     8316      res->wvhdl[i]   = r->wvhdl[i-1];
    82708317    }
    82718318
    82728319    // new 1st block
    82738320    res->order[0]   = ringorder_C; // Prefix
    8274     res->block0[0]  = 1;
    8275     res->block1[0]  = res->N;
    8276     //res->wvhdl[j]   = NULL;
    8277     // res->order [j] = 0; // The End!
     8321    // removes useless secondary component order if defined in old ring
     8322    for (int i=rBlocks(res); i>0; --i) {
     8323      if (res->order[i] == ringorder_C || res->order[i] == ringorder_c) {
     8324        res->order[i] = 0;
     8325        break;
     8326      }
     8327    }
    82788328    rComplete(res, 1);
    82798329#ifdef HAVE_PLURAL
     
    82978347    return (res);
    82988348  }
    8299 
    8300   // not incremental => use Schreyer order
     8349  // if sbaOrder == 3 => degree - position - ring order
     8350  if (strat->sbaOrder == 3)
     8351  {
     8352    ring res = rCopy0(r, TRUE, FALSE);
     8353    res->order  = (int *)omAlloc0((n+2)*sizeof(int));
     8354    res->block0 = (int *)omAlloc0((n+2)*sizeof(int));
     8355    res->block1 = (int *)omAlloc0((n+2)*sizeof(int));
     8356    int **wvhdl = (int **)omAlloc0((n+2)*sizeof(int*));
     8357    res->wvhdl  = wvhdl;
     8358    for (int i=2; i<n+2; i++)
     8359    {
     8360      res->order[i]   = r->order[i-2];
     8361      res->block0[i]  = r->block0[i-2];
     8362      res->block1[i]  = r->block1[i-2];
     8363      res->wvhdl[i]   = r->wvhdl[i-2];
     8364    }
     8365
     8366    // new 1st block
     8367    res->order[0]   = ringorder_a; // Prefix
     8368    res->block0[0]  = 1;
     8369    res->wvhdl[0]   = (int *)omAlloc(res->N*sizeof(int));
     8370    for (int i=0; i<res->N; ++i)
     8371      res->wvhdl[0][i]  = 1;
     8372    res->block1[0]  = si_min(res->N, rVar(res));
     8373    // new 2nd block
     8374    res->order[1]   = ringorder_C; // Prefix
     8375    res->wvhdl[1]   = NULL;
     8376    // removes useless secondary component order if defined in old ring
     8377    for (int i=rBlocks(res); i>0; --i) {
     8378      if (res->order[i] == ringorder_C || res->order[i] == ringorder_c) {
     8379        res->order[i] = 0;
     8380        break;
     8381      }
     8382    }
     8383    rComplete(res, 1);
     8384#ifdef HAVE_PLURAL
     8385    if (rIsPluralRing(r))
     8386    {
     8387      if ( nc_rComplete(r, res, false) ) // no qideal!
     8388      {
     8389#ifndef NDEBUG
     8390        WarnS("error in nc_rComplete");
     8391#endif
     8392        // cleanup?
     8393
     8394        //      rDelete(res);
     8395        //      return r;
     8396
     8397        // just go on..
     8398      }
     8399    }
     8400#endif
     8401    strat->tailRing = res;
     8402    return (res);
     8403  }
     8404
     8405  // not sbaOrder == 1 => use Schreyer order
    83018406  // this is done by a trick when initializing the signatures
    83028407  // in initSLSba():
  • kernel/kutil.h

    r74c4462 r6e3023a  
    180180public:
    181181  unsigned long sev;
    182   unsigned long from; // from which polynomial it comes from
    183             // this is important for signature-based
    184             // algorithms
     182  unsigned long from; // index in sig up to which the correspongin LObject was already checked
    185183  unsigned long checked; // this is the index of S up to which
    186184                      // the corresponding LObject was already checked in
     
    188186                      // reduction process it is enough to start a second
    189187                      // rewritten criterion check from checked+1 onwards
    190                       // NOTE: If checked = 3 then the corresponding pair is
     188  BOOLEAN prod_crit;
     189                      // NOTE: If prod_crit = TRUE then the corresponding pair is
    191190                      // detected by Buchberger's Product Criterion and can be
    192191                      // deleted
     
    289288  void (*chainCrit) (poly p,int ecart,kStrategy strat);
    290289  BOOLEAN (*syzCrit) (poly sig, unsigned long not_sevSig, kStrategy strat);
    291   BOOLEAN (*rewCrit1) (poly sig, unsigned long not_sevSig, kStrategy strat, int start /*= 0*/);
    292   BOOLEAN (*rewCrit2) (poly sig, unsigned long not_sevSig, kStrategy strat, int start /*= 0*/);
     290  BOOLEAN (*rewCrit1) (poly sig, unsigned long not_sevSig, poly lm, kStrategy strat, int start /*= 0*/);
     291  BOOLEAN (*rewCrit2) (poly sig, unsigned long not_sevSig, poly lm, kStrategy strat, int start /*= 0*/);
     292  BOOLEAN (*rewCrit3) (poly sig, unsigned long not_sevSig, poly lm, kStrategy strat, int start /*= 0*/);
    293293  pFDegProc pOrigFDeg;
    294294  pLDegProc pOrigLDeg;
     
    310310                // syzygy of component i comes up
    311311                // important for signature-based algorithms
    312   BOOLEAN incremental;
     312  unsigned sbaOrder;
    313313  unsigned long currIdx;
    314314  int max_lower_index;
     
    446446int posInLSig (const LSet set, const int length,
    447447               LObject* L,const kStrategy strat);
     448int posInSyz (const kStrategy strat, const poly sig);
    448449int posInL0 (const LSet set, const int length,
    449450             LObject* L,const kStrategy strat);
     
    467468poly redtailBba (LObject *L, int pos,kStrategy strat,
    468469                 BOOLEAN withT = FALSE,BOOLEAN normalize=FALSE);
     470poly redtailSba (LObject *L, int pos,kStrategy strat,
     471                 BOOLEAN withT = FALSE,BOOLEAN normalize=FALSE);
    469472poly redtailBba (TObject *T, int pos,kStrategy strat);
    470473poly redtail (poly p,int pos,kStrategy strat);
     
    518521void initSyzRules (kStrategy strat);
    519522void updateS(BOOLEAN toT,kStrategy strat);
    520 void enterSyz (LObject p,kStrategy strat);
     523void enterSyz (LObject p,kStrategy strat, int atT);
    521524void enterT (LObject p,kStrategy strat, int atT = -1);
    522525void cancelunit (LObject* p,BOOLEAN inNF=FALSE);
     
    542545BOOLEAN syzCriterion(poly sig, unsigned long not_sevSig, kStrategy strat);
    543546BOOLEAN syzCriterionInc(poly sig, unsigned long not_sevSig, kStrategy strat);
    544 KINLINE BOOLEAN arriRewDummy(poly sig, unsigned long not_sevSig, kStrategy strat, int start);
    545 BOOLEAN arriRewCriterion(poly sig, unsigned long not_sevSig, kStrategy strat, int start);
    546 BOOLEAN faugereRewCriterion(poly sig, unsigned long not_sevSig, kStrategy strat, int start);
     547KINLINE BOOLEAN arriRewDummy(poly sig, unsigned long not_sevSig, poly lm, kStrategy strat, int start);
     548BOOLEAN arriRewCriterion(poly sig, unsigned long not_sevSig, poly lm, kStrategy strat, int start);
     549BOOLEAN arriRewCriterionPre(poly sig, unsigned long not_sevSig, poly lm, kStrategy strat, int start);
     550BOOLEAN faugereRewCriterion(poly sig, unsigned long not_sevSig, poly lm, kStrategy strat, int start);
    547551BOOLEAN findMinLMPair(poly sig, unsigned long not_sevSig, kStrategy strat, int start);
    548552// returns index of p in TSet, or -1 if not found
Note: See TracChangeset for help on using the changeset viewer.