Changeset a5ab27 in git


Ignore:
Timestamp:
Aug 8, 2014, 11:31:25 AM (10 years ago)
Author:
Hans Schoenemann <hannes@…>
Branches:
(u'spielwiese', '60097d763e0b541617a3b864b5310c523edaf81d')
Children:
30994898b960fef81ba9392060ceb124019cd069
Parents:
fde66a9b59619d7999ad28b0c18422dcddba436a
git-author:
Hans Schoenemann <hannes@mathematik.uni-kl.de>2014-08-08 11:31:25+02:00
git-committer:
Hans Schoenemann <hannes@mathematik.uni-kl.de>2014-08-08 12:12:17+02:00
Message:
Revert "Groebner Walk, Sagbi Walk"

This reverts commit ffd4a9600207a9670cfe0cb81450accbaaec26d7.

Conflicts:
	Singular/LIB/modwalk.lib
	Singular/LIB/rwalk.lib
	Singular/LIB/swalk.lib
	Singular/walk.cc
	Singular/walk.h

To many bugs: - syntax error (mainly: include paths)
              - changed interface without changed libs (grwalk.lib)
	      - several failing tests (in Manual/*, Buch/*, etc.)
Location:
Singular
Files:
6 edited

Legend:

Unmodified
Added
Removed
  • Singular/LIB/modwalk.lib

    • Property mode changed from 100755 to 100644
    rfde66a9 ra5ab27  
    3131proc modpWalk(def II, int p, int variant, list #)
    3232"USAGE:  modpWalk(I,p,#); I ideal, p integer, variant integer
    33 ASSUME:  If size(#) > 0, then
     33ASSUME:  If size(#) > 0, then 
    3434           #[1] is an intvec describing the current weight vector
    3535           #[2] is an intvec describing the target weight vector
    3636RETURN:  ideal - a standard basis of I mod p, integer - p
    3737NOTE:    The procedure computes a standard basis of the ideal I modulo p and
    38          fetches the result to the basering.
     38         fetches the result to the basering. 
    3939EXAMPLE: example modpWalk; shows an example
    4040"
     
    8282
    8383//-------------------------  make i homogeneous  -----------------------------
    84 /*  if(!mixedTest() && !h)
     84
     85  if(!hasMixedOrdering() && !h)
    8586  {
    8687    if(!((find(ordstr_R0, "M") > 0) || (find(ordstr_R0, "a") > 0) || neg))
     
    8990      {
    9091        list rl@r = ringlist(@r);
    91         nvar@r = nvars(@r);
     92        nvar@r = nvars(@r);     
    9293        intvec w;
    9394        for(k = 1; k <= nvar@r; k++)
     
    105106    }
    106107  }
    107 */
     108
    108109//-------------------------  compute a standard basis mod p  -----------------------------
    109110
     
    146147    if(size(#) == 2)
    147148    {
    148       i=fwalk(i,#);
     149      trwalk(i,radius,pert_deg,#);
    149150    }
    150151    else
    151152    {
    152       i=fwalk(i);
    153     }
    154   }
    155   if(variant == 5)
    156   {
    157     if(size(#) == 2)
    158     {
    159      i=prwalk(i,radius,pert_deg,pert_deg,#);
    160     }
    161     else
    162     {
    163       i=prwalk(i,radius,pert_deg,pert_deg);
    164     }
    165   }
    166   if(variant == 6)
    167   {
    168     if(size(#) == 2)
    169     {
    170       i=pwalk(i,pert_deg,pert_deg,#);
    171     }
    172     else
    173     {
    174       i=pwalk(i,pert_deg,pert_deg);
    175     }
    176   }
    177   /*
    178   if(!mixedTest() && !h)
     153      trwalk(i,radius,pert_deg);
     154    }
     155  }
     156  if(!hasMixedOrdering() && !h)
    179157  {
    180158    if(!((find(ordstr_R0, "M") > 0) || (find(ordstr_R0, "a") > 0) || neg))
     
    190168      }
    191169    }
    192   }*/
     170  }
    193171  setring R0;
    194172  return(list(fetch(@r,i),p));
     
    202180  intvec a = 2,1,3,4;
    203181  intvec b = 1,9,1,1;
    204   ring ra = 0,x(1..4),(a(a),lp);
     182  ring ra = 0,(w,x,y,z),(a(a),lp);
    205183  ideal I = std(cyclic(4));
    206   ring rb = 0,x(1..4),(a(b),lp);
     184  ring rb = 0,(w,x,y,z),(a(b),lp);
    207185  ideal I = imap(ra,I);
    208186  modpWalk(I,p,1,a,b);
     
    236214proc modWalk(def II, int variant, list #)
    237215"USAGE:  modWalk(II); II ideal or list(ideal,int)
    238 ASSUME:  If variant =
    239 @*       - 1 the Random Walk algorithm with radius II[2] is applied
    240            to II[1] if II = list(ideal, int). It is applied to II with radius 2
    241            if II is an ideal.
    242 @*       - 2, the Groebner Walk algorithm is applied to II[1] or to II, respectively.
    243 @*       - 3, the Fractal Walk algorithm with random element is applied to II[1] or II,
    244            respectively.
    245 @*       - 4, the Fractal Walk algorithm is applied to II[1] or II, respectively.
    246 @*       - 5, the Perturbation Walk algorithm with random element is applied to II[1]
    247            or II, respectively, with radius II[3] and perturbation degree II[2].
    248 @*       - 6, the Perturbation Walk algorithm is applied to II[1] or II, respectively,
    249            with perturbation degree II[3].
     216ASSUME:  If variant = 1 the random walk algorithm with radius II[2] is applied
     217         to II[1] if II = list(ideal, int). It is applied to II with radius 2
     218         if II is an ideal. If variant = 2, the Groebner walk algorithm is
     219         applied to II[1] or to II, respectively.
    250220         If size(#) > 0, then # contains either 1, 2 or 4 integers such that
    251221@*       - #[1] is the number of available processors for the computation,
     
    300270        II[3] = 2;
    301271      }
     272   
    302273    }
    303274    else
     
    377348          ERROR("Unexpected type of input.");
    378349        }
    379       }
     350      } 
    380351    }
    381352    if(size(#) == 3)
     
    497468//-------------------------  Save current options  -----------------------------
    498469  intvec opt = option(get);
    499   //option(redSB);
     470  option(redSB);
    500471
    501472//--------------------  Initialize the list of primes  -------------------------
     
    556527    }
    557528    H = chinrem(T1,T2);
    558     J = parallelFarey(H,N,n1);
    559     //J=farey(H,N);
     529    //J = parallelFarey(H,N,n1);
     530    J=farey(H,N);
    560531    if(printlevel >= 10)
    561532    {
     
    567538
    568539    tt = timer; rt = rtimer;
    569     pTest = pTestSB(I,J,L,variant);
    570     //pTest = primeTestSB(I,J,L,variant);
     540    //pTest = pTestSB(I,J,L,variant);
     541    pTest = primeTestSB(I,J,L,variant);
    571542    if(printlevel >= 10)
    572543    {
     
    622593    tt = timer; rt = rtimer;
    623594    L = primeList(I,n3,L,n1);
     595L;
    624596    if(printlevel >= 10)
    625597    {
     
    672644  ideal I=-x+y2z-z,xz+1,x2+y2-1;
    673645  // I is a standard basis in dp
    674   ideal J = modWalk(I,1);
    675   J;
    676 /*  modWalk(I,2,2,0);
     646  modWalk(I,1);
     647  modWalk(I,2,2,0);
    677648  modWalk(I,3,system("cpu"),0);
    678649  std(I);
     
    683654  ideal i=fetch(r0,i0);
    684655  modWalk(i,1,system("cpu"),0);
    685   modWalk(i,3);*/
     656  modWalk(i,3);
    686657}
    687658
     
    738709  }
    739710  results_chinrem = parallelWaitAll("chinrem",arguments_chinrem);
    740     for(j=1; j <= size(results_chinrem); j++)
     711    for(j=1; j <= size(results_chinrem); j++) 
    741712    {
    742713      J = results_chinrem[j];
  • Singular/LIB/rwalk.lib

    • Property mode changed from 100755 to 100644
    rfde66a9 ra5ab27  
    1 //
    2 version="version rwalk.lib 4.0.0.0 Aug_2014 ";
     1/////////////////////////////////////////////////
     2version="version rwalk.lib 4.0.0.0 Jun_2013 "; // $Id$
    33category="Commutative Algebra";
    44
     
    1414 * Argument string for Random Walk *
    1515 ***********************************/
    16 proc OrderStringalp_NP(string Wpal,list #)
     16static proc OrderStringalp_NP(string Wpal,list #)
    1717{
    1818  int n= nvars(basering);
     
    3535        if(Wpal == "al"){
    3636          order_str = "(a("+string(#[1])+"),lp("+string(n) + "),C)";
    37         }
    38         if(Wpal == "M"){
    39           order_str = "(M("+string(#[1])+"),C)";
    4037        }
    4138        else {
     
    6461           order_str = "(a("+string(#[1])+"),lp("+string(n) + "),C)";
    6562         }
    66          if(Wpal == "M"){
    67           order_str = "(M("+string(#[1])+"),C)";
    68          }
    6963         else {
    7064           order_str = "(Wp("+string(#[1])+"),C)";
     
    7973             order_str = "(a("+string(#[1])+"),lp("+string(n) + "),C)";
    8074           }
    81            if(Wpal == "M"){
    82              order_str = "(M("+string(#[1])+"),C)";
    83            }
    8475           else {
    8576             order_str = "(Wp("+string(#[1])+"),C)";
    86            //  order_str = "(a("+string(#[1])+"),C)";
    8777           }
    8878         }
     
    10595             order_str = "(a("+string(#[1])+"),lp("+string(n) + "),C)";
    10696           }
    107            if(Wpal == "M"){
    108              order_str = "(M("+string(#[1])+"),C)";
    109            }
    11097           else {
    11198             order_str = "(Wp("+string(#[1])+"),C)";
     
    135122}
    136123
    137 /****************
    138  * Random Walk *
    139  ****************/
     124/************************************
     125 * Random Walk with perturbation    *
     126 ************************************/
    140127proc rwalk(ideal Go, int radius, int pert_deg, list #)
    141128"SYNTAX: rwalk(ideal i, int radius);
    142          if size(#)>0 then rwalk(ideal i, int radius, intvec v, intvec w);
     129         rwalk(ideal i, int radius, intvec v, intvec w);
    143130TYPE:    ideal
    144131PURPOSE: compute the standard basis of the ideal, calculated via
     
    151138{
    152139//--------------------  Initialize parameters  ------------------------
    153 int n= nvars(basering);
    154 list OSCTW = OrderStringalp_NP("al",#);
    155 if(size(#)>1)
    156   {
    157   if(size(#[2]) == n*n)
    158     {
    159     OSCTW= OrderStringalp_NP("M", #);
    160     }
    161   }
    162 else
    163   {
    164   OSCTW= OrderStringalp_NP("al", #);
    165   }
     140list OSCTW = OrderStringalp_NP("al", #);
    166141string ord_str = OSCTW[2];
    167142intvec curr_weight = OSCTW[3]; // original weight vector
    168143intvec target_weight = OSCTW[4]; // target weight vector
    169144kill OSCTW;
     145option(redSB);
    170146
    171147//--------------------  Initialize parameters  ------------------------
    172148def xR = basering;
    173 execute("ring ostR = "+charstr(xR)+",("+varstr(xR)+"),"+ord_str+";");
     149execute("ring ostR = ("+charstr(xR)+"),("+varstr(xR)+"),"+ord_str+";");
    174150def old_ring = basering;
    175151
    176152ideal G = fetch(xR, Go);
    177 G = system("Mrwalk", G, curr_weight, target_weight, radius, pert_deg, basering);
     153G = system("Mrwalk", G, curr_weight, target_weight, radius, pert_deg);
    178154
    179155setring xR;
     
    194170  int perturb_deg = 2;
    195171  rwalk(I,radius,perturb_deg);
    196 }
    197 
    198 /*****************************************
    199  * Perturbation Walk with random element *
    200  *****************************************/
    201 proc prwalk(ideal Go, int radius, int o_pert_deg, int t_pert_deg, list #)
    202 "SYNTAX: rwalk(ideal i, int radius);
    203          if size(#)>0 then rwalk(ideal i, int radius, intvec v, intvec w);
    204 TYPE:    ideal
    205 PURPOSE: compute the standard basis of the ideal, calculated via
    206          the Random walk algorithm  from the ordering
    207          \"(a(v),lp)\", \"dp\" or \"Dp\"
    208          to the ordering  \"(a(w),lp)\" or \"(a(1,0,...,0),lp)\".
    209 SEE ALSO: std, stdfglm, groebner, gwalk, pwalk, fwalk, twalk, awalk1, awalk2
    210 KEYWORDS: Groebner walk
    211 EXAMPLE: example rwalk; shows an example"
    212 {
    213 //--------------------  Initialize parameters  ------------------------
    214 list OSCTW = OrderStringalp_NP("al", #);
    215 string ord_str = OSCTW[2];
    216 intvec curr_weight = OSCTW[3]; // original weight vector
    217 intvec target_weight = OSCTW[4]; // target weight vector
    218 kill OSCTW;
    219 
    220 //--------------------  Initialize parameters  ------------------------
    221 def xR = basering;
    222 execute("ring ostR = ("+charstr(xR)+"),("+varstr(xR)+"),"+ord_str+";");
    223 def old_ring = basering;
    224 
    225 ideal G = fetch(xR, Go);
    226 G = system("Mprwalk", G, curr_weight, target_weight, radius, o_pert_deg, t_pert_deg, basering);
    227 
    228 setring xR;
    229 kill Go;
    230 
    231 keepring basering;
    232 ideal result = fetch(old_ring, G);
    233 attrib(result,"isSB",1);
    234 return (result);
    235 }
    236 example
    237 {
    238   "EXAMPLE:"; echo = 2;
    239   // compute a Groebner basis of I w.r.t. lp.
    240   ring r = 32003,(z,y,x), lp;
    241   ideal I = y3+xyz+y2z+xz3, 3+xy+x2y+y2z;
    242   int radius = 1;
    243   int o_perturb_deg = 2;
    244   int t_perturb_deg = 2;
    245   prwalk(I,radius,o_perturb_deg,t_perturb_deg);
    246172}
    247173
     
    264190
    265191   string ord_str =   OSCTW[2];
    266    intvec curr_weight   =   OSCTW[3]; /* current weight vector */
     192   intvec curr_weight   =   OSCTW[3]; /* original weight vector */
    267193   intvec target_weight =   OSCTW[4]; /* target weight vector */
    268194   kill OSCTW;
     195   option(redSB);
    269196   def xR = basering;
    270197
     
    273200   //print("//** help ring = " + string(basering));
    274201   ideal G = fetch(xR, Go);
    275    int pert_deg = 2;
    276    G = system("Mfrwalk", G, curr_weight, target_weight, pert_deg, radius, basering);
     202
     203   G = system("Mfrwalk", G, curr_weight, target_weight, radius);
    277204
    278205   setring xR;
     
    315242  kill L;
    316243
     244  option(redSB);
    317245  def xR = basering;
    318246
  • Singular/LIB/swalk.lib

    • Property mode changed from 100755 to 100644
    rfde66a9 ra5ab27  
    2323
    2424LIB "sagbi.lib";
    25 LIB "crypto.lib";
    26 LIB "atkins.lib";
    27 LIB "random.lib";
    2825//////////////////////////////////////////////////////////////////////////////
    2926proc swalk(ideal Gox, list #)
     
    4239  intvec itarget_weight =   OSCTW[4]; /* terget weight vector */
    4340  kill OSCTW;
    44       option(redSB);
     41  option(redSB);
    4542  def xR = basering;
    4643  list rl=ringlist(xR);
    47   if(size(#) == 0)
    48   {
    49     rl[3][1][1]="dp";
    50   }
     44  rl[3][1][1]="dp";
    5145  def ostR=ring(rl);
    5246  setring ostR;
     
    5751  vector curr_weight=changeTypeInt(icurr_weight);
    5852  vector target_weight=changeTypeInt(itarget_weight);
    59   curr_weight;
    6053  ideal Gold;
    6154  list l;
     
    8679              ideal Gnew=Convert(Gold);
    8780              Gnew=interreduceSd(Gnew);
    88          }
    89     }
    90    setring xR;
    91    ideal result = fetch(old_ring, Gnew);
    92    kill old_ring;
    93    attrib(result,"isSB",1);
    94    return (result);
    95 }
    96 example
    97 {
    98   "EXAMPLE:";echo = 2;
    99   ring r = 0,(x,y), lp;
    100   ideal I =x2,y2,xy+y,2xy2+y3;
    101   swalk(I);
    102 }
    103 //////////////////////////////////////////////////////////////////////////////
    104 proc rswalk(ideal Gox, int weight_rad, int pdeg, list #)
    105 "USAGE:  rswalk(i[,v,w]); i ideal, v,w int vectors
    106 RETURN: The sagbi basis of the subalgebra defined by the generators of i,
    107         calculated via the Sagbi walk algorithm from the ordering dp to lp
    108         if v,w are not given (resp. from the ordering (a(v),lp) to the
    109         ordering (a(w),lp) if v and w are given).
    110 EXAMPLE: example swalk; shows an example
    111 "
    112 {
    113   /* we use ring with ordering (a(...),lp,C) */
    114   list OSCTW    = OrderStringalp_NP("al", #);//"dp"
    115   string ord_str =   OSCTW[2];
    116   intvec icurr_weight   =   OSCTW[3]; /* original weight vector */
    117   intvec itarget_weight =   OSCTW[4]; /* terget weight vector */
    118   kill OSCTW;
    119   option(redSB);
    120   def xR = basering;
    121   list rl=ringlist(xR);
    122   rl[3][1][1]="dp";
    123   def ostR=ring(rl);
    124   setring ostR;
    125   def new_ring = basering;
    126   ideal Gnew = fetch(xR, Gox);
    127   Gnew=sagbi(Gnew,1);
    128   Gnew=interreduceSd(Gnew);
    129   vector curr_weight=changeTypeInt(icurr_weight);
    130   vector target_weight=changeTypeInt(itarget_weight);
    131   ideal Gold;
    132   list l;
    133   intvec v;
    134   int n=0;
    135   while(n==0)
    136     {
    137        Gold=Gnew;
    138        def old_ring=new_ring;
    139        setring old_ring;
    140 
    141        kill new_ring;
    142        if(curr_weight==target_weight){n=1;}
    143        else {
    144               l=collectDiffExpo(Gold);
    145               vector new_weight=RandomNextWeight(Gold, l, curr_weight, target_weight, weight_rad, pdeg);
    146               vector w=cleardenom(new_weight);
    147               v=changeType(w);
    148               list p= ringlist(old_ring);
    149               p[3][1][2]= v;
    150               def new_ring=ring(p);
    151               setring new_ring;
    152               ideal Gold=fetch(old_ring,Gold);
    153               vector curr_weight=fetch(old_ring,new_weight);
    154               vector target_weight=fetch(old_ring,target_weight);
    155               kill old_ring;
    156               ideal Gnew=Convert(Gold);
    157               Gnew=interreduceSd(Gnew);
    15881           }
    15982    }
     
    16891  ring r = 0,(x,y), lp;
    16992  ideal I =x2,y2,xy+y,2xy2+y3;
    170   swalk(I,7);
    171 }
     93  swalk(I);
     94}
     95
    17296//////////////////////////////////////////////////////////////////////////////
    17397static proc inprod(vector v,vector w)
     
    287211
    288212//////////////////////////////////////////////////////////////////////////////
    289 static proc test_in_cone(vector w, list l)
    290 {
    291  int i,j,k;
    292  vector v;
    293  poly n;
    294  number a;
    295  for(i=1;i<=size(l);i++)
    296  {
    297     for(j=1;j<=size(l[i]);j++)
    298     {
    299         v=0;
    300         for(k=1;k<=size(l[i][j]);k++)
    301         {
    302             v=v+l[i][j][k]*gen(k);
    303         }
    304         n = inprod(w,v);
    305         a = leadcoef(n);
    306         if(a<0)
    307         {
    308            return(0);
    309         }
    310     }
    311  }
    312  return(1);
    313 }
    314 
    315 
    316 //////////////////////////////////////////////////////////////////////////////
    317213static proc last( vector c, vector t,list l)
    318214"USAGE: last(c,t,l); c,t vectors, l list
     
    325221 number ul=1;
    326222 int i,j,k;
    327  number a,b,z,u;
    328  poly n,q;
     223 number u;
    329224 vector v;
    330225 for(i=1;i<=size(l);i++)
     
    337232            v=v+l[i][j][k]*gen(k);
    338233        }
    339         n= inprod(c,v);
    340         q= inprod(t,v);
    341         a=leadcoef(n);
    342         b=leadcoef(q);
    343         z=a-b;
     234        poly n= inprod(c,v);
     235        poly q= inprod(t,v);
     236        number a=leadcoef(n);
     237        number b=leadcoef(q);
     238        number z=a-b;
    344239        if(b<0)
    345             {
     240        {
    346241            u=a/z;
    347             if(u<ul)
    348             {
    349                ul=u;
    350              }
    351         }
     242            if(u<ul) {ul=u;}
     243        }
     244        kill a,b,z,n,q ;
    352245    }
    353  }~
    354  kill a,b,z,n,q,u;
     246 }
    355247 return(ul);
    356248}
     
    448340   Lift(In,InG,Gold);
    449341}
    450 //////////////////////////////////////////////////////////////////////////////
    451 static proc PertVectors(ideal Gold, vector target_weight, int pdeg)
    452 {
    453 int nV = nvars(basering);
    454 int nG = size(Gold);
    455 int i;
    456 number ntemp, maxAi, maxA;
    457 if(pdeg > nV || pdeg <= 0)
    458   {
    459     intvec v_null=0;
    460     return v_null;
    461   }
    462 if(pdeg == 1)
    463   {
    464     return target_weight;
    465   }
    466 maxAi=0;
    467 for(i=1; i<=nV; i++)
    468   {
    469     ntemp = leadcoef(inprod(target_weight,gen(i)));
    470     if(ntemp < 0)
    471       {
    472         ntemp = -ntemp;
    473       }
    474     if(maxAi < ntemp)
    475       {
    476         maxAi = ntemp;
    477       }
    478   }
    479 maxA = maxAi+pdeg-1;
    480 number epsilon = maxA*deg(Gold)+1;
    481 vector pert_weight = epsilon^(pdeg-1)*target_weight;
    482 for(i=2; i<=pdeg; i++)
    483   {
    484     pert_weight = pert_weight + epsilon^(pdeg-i)*gen(i);
    485   }
    486 return(pert_weight);
    487 }
    488 
    489 
    490 //////////////////////////////////////////////////////////////////////////////
    491 static proc RandomNextWeight(ideal Gold, list L, vector curr_weight, vector target_weight,int weight_rad, int pdeg)
    492 "USAGE: RandomNextWeight(Gold, L, curr_weight, target_weight);
    493 RETURN: Intermediate next weight vector
    494 EXAMPLE: example RandomNextWeight; shows an example
    495 "
    496 {
    497   int i,n1,n2,n3;
    498   number norm, weight_norm;
    499   def Rold = basering;
    500   int nV = nvars(basering);
    501   number ulast=last(curr_weight, target_weight, L);
    502   vector new_weight=(1-ulast)*curr_weight+ulast*target_weight;
    503   vector w1=cleardenom(new_weight);
    504   intvec v1=changeType(w1);
    505   list p= ringlist(Rold);
    506   p[3][1][2]= v1;
    507   def new_ring=ring(p);
    508   setring new_ring;
    509   ideal Gold = fetch(Rold, Gold);
    510   n1=size(Initial(Gold));
    511   setring Rold;
    512   intvec next_weight;
    513   kill new_ring;
    514   while(1)
    515   {
    516     weight_norm = 0;
    517     while(weight_norm == 0)
    518     {
    519       for(i=1; i<=nV; i++)
    520       {
    521         next_weight[i] = random(1,10000)-5000;
    522         weight_norm = weight_norm + next_weight[i]^2;
    523       }
    524       norm = 0;
    525       while(norm^2 < weight_norm)
    526       {
    527          norm=norm+1;
    528       }
    529       weight_norm = 1+norm;
    530     }
    531     new_weight = 0;
    532     for(i=1; i<=nV;i++)
    533     {
    534       if(next_weight[i] < 0)
    535       {
    536         new_weight = new_weight + (1 + round(weight_rad*leadcoef(next_weight[i])/weight_norm))*gen(i);
    537       }
    538       else
    539       {
    540         new_weight = new_weight + ( round(weight_rad*leadcoef(next_weight[i])/weight_norm))*gen(i);
    541       }
    542     }
    543     new_weight = new_weight + curr_weight;
    544     if(test_in_cone(new_weight, L)==1)
    545     {
    546       break;
    547     }
    548   }
    549   kill next_weight;
    550   kill norm;
    551   vector w2=cleardenom(new_weight);
    552   intvec v2=changeType(w2);
    553   p[3][1][2]= v2;
    554   def new_ring=ring(p);
    555   setring new_ring;
    556   ideal Gold = fetch(Rold, Gold);
    557   n2=size(Initial(Gold));
    558   setring Rold;
    559   kill new_ring;
    560 
    561   vector w3=cleardenom(PertVectors(Gold,target_weight,pdeg));
    562   intvec v3=changeType(w3);
    563   p[3][1][2]= v1;
    564   def new_ring=ring(p);
    565   setring new_ring;
    566   ideal Gold = fetch(Rold, Gold);
    567   n3=size(Initial(Gold));
    568   setring Rold;
    569   kill new_ring;
    570   kill p;
    571 
    572   if(n2<n1)
    573   {
    574     if(n3<n2)
    575     {
    576       // n3<n2<n1
    577       return(w3);
    578     }
    579     else
    580     {
    581       // n2<n1 und n2<=n3
    582       return(w2);
    583     }
    584   }
    585   else
    586   {
    587     if(n3<n1)
    588     {
    589       //n3<n1<=n2
    590       return(w3);
    591     }
    592     else
    593     {
    594       // n1<=n3 und n1<=n2
    595       return(w1);
    596     }
    597   }
    598 }
    599342
    600343//////////////////////////////////////////////////////////////////////////////
     
    606349{
    607350 ideal In=Initial(Gold);
    608  ideal InG=sagbi(In); //+In;
     351 ideal InG=sagbi(In,1)+In;
    609352 ideal Gnew=Lift(In,InG,Gold);
    610353 return(Gnew);
     
    637380     for(i=1;i<=size(M);i++)
    638381     { B[i]=M[i];}
    639      f=sagbiNF(f,B,1)[1];
     382     f=sagbiNF(f,B,1);
    640383     J=J+f;
    641384  }
     
    855598
    856599///////////////////////////////////////////////////////////////////////////////*
    857 /*Further examples
     600Further examples
    858601ring r=0,(x,y,z),lp;
    859602
  • Singular/extra.cc

    • Property mode changed from 100755 to 100644
    rfde66a9 ra5ab27  
    88#define HAVE_WALK 1
    99
    10 #ifdef HAVE_CONFIG_H
    11 #include "singularconfig.h"
    12 #endif /* HAVE_CONFIG_H */
     10
     11
     12
    1313#include <kernel/mod2.h>
    1414#include <misc/auxiliary.h>
    15 
    16 #define SI_DONT_HAVE_GLOBAL_VARS
     15#include <misc/sirandom.h>
     16
    1717#include <factory/factory.h>
    1818
     
    5555
    5656
     57#include <resources/feResource.h>
    5758#include <polys/monomials/ring.h>
    5859#include <kernel/polys.h>
     
    6869#include <kernel/fast_mult.h>
    6970#include <kernel/digitech.h>
    70 #include <kernel/stairc.h>
    71 #include <kernel/febase.h>
     71#include <kernel/GBEngine/stairc.h>
    7272#include <kernel/ideals.h>
    73 #include <kernel/kstd1.h>
    74 #include <kernel/syz.h>
    75 #include <kernel/kutil.h>
    76 
    77 #include <kernel/shiftgb.h>
    78 #include <kernel/linearAlgebra.h>
     73#include <kernel/GBEngine/kstd1.h>
     74#include <kernel/GBEngine/syz.h>
     75#include <kernel/GBEngine/kutil.h>
     76
     77#include <kernel/GBEngine/shiftgb.h>
     78#include <kernel/linear_algebra/linearAlgebra.h>
     79
     80#include <kernel/combinatorics/hutil.h>
    7981
    8082// for tests of t-rep-GB
    81 #include <kernel/tgb.h>
    82 
     83#include <kernel/GBEngine/tgb.h>
     84
     85#include <kernel/linear_algebra/minpoly.h>
     86
     87#include <numeric/mpr_base.h>
    8388
    8489#include "tok.h"
     
    9297#include "distrib.h"
    9398
    94 #include "minpoly.h"
    9599#include "misc_ip.h"
    96100
     
    105109
    106110#ifdef HAVE_RINGS
    107 #include <kernel/ringgb.h>
     111#include <kernel/GBEngine/ringgb.h>
    108112#endif
    109113
    110114#ifdef HAVE_F5
    111 #include <kernel/f5gb.h>
     115#include <kernel/GBEngine/f5gb.h>
    112116#endif
    113117
     
    118122
    119123#ifdef HAVE_SPECTRUM
    120 #include <kernel/spectrum.h>
     124#include <kernel/spectrum/spectrum.h>
    121125#endif
    122126
     
    125129#include <polys/nc/ncSAMult.h> // for CMultiplier etc classes
    126130#include <polys/nc/sca.h>
    127 #include <kernel/nc.h>
     131#include <kernel/GBEngine/nc.h>
    128132#include "ipconv.h"
    129133#ifdef HAVE_RATGRING
    130 #include <kernel/ratgring.h>
     134#include <kernel/GBEngine/ratgring.h>
    131135#endif
    132136#endif
     
    147151
    148152#include <polys/clapconv.h>
    149 #include <kernel/kstdfac.h>
     153#include <kernel/GBEngine/kstdfac.h>
    150154
    151155#include <polys/clapsing.h>
     
    185189static BOOLEAN jjEXTENDED_SYSTEM(leftv res, leftv h);
    186190#endif
    187 
    188 extern BOOLEAN jjJanetBasis(leftv res, leftv v);
    189191
    190192#ifdef ix86_Win  /* PySingular initialized? */
     
    19061908        if (h == NULL || h->Typ() != IDEAL_CMD ||
    19071909            h->next == NULL || h->next->Typ() != INTVEC_CMD ||
    1908             h->next->next == NULL || h->next->next->Typ() != INTVEC_CMD ||
    1909             h->next->next->next == NULL || h->next->next->next->Typ() != RING_CMD)
    1910         {
    1911           Werror("system(\"Mwalk\", ideal, intvec, intvec, ring) expected");
    1912           return TRUE;
    1913         }
    1914 
    1915         if ((((intvec*) h->next->Data())->length() != currRing->N &&
    1916             ((intvec*) h->next->next->Data())->length() != currRing->N ) &&
    1917             (((intvec*) h->next->Data())->length() != (currRing->N)*(currRing->N) &&
    1918             ((intvec*) h->next->next->Data())->length() != (currRing->N)*(currRing->N) ))
    1919         {
    1920           Werror("system(\"Mwalk\" ...) intvecs not of length %d or %d\n",
    1921                  currRing->N,(currRing->N)*(currRing->N));
    1922           return TRUE;
    1923         }
     1910            h->next->next == NULL || h->next->next->Typ() != INTVEC_CMD)
     1911        {
     1912          Werror("system(\"Mwalk\", ideal, intvec, intvec) expected");
     1913          return TRUE;
     1914        }
     1915
     1916        if (((intvec*) h->next->Data())->length() != currRing->N &&
     1917            ((intvec*) h->next->next->Data())->length() != currRing->N )
     1918        {
     1919          Werror("system(\"Mwalk\" ...) intvecs not of length %d\n",
     1920                 currRing->N);
     1921          return TRUE;
     1922        }
    19241923        ideal arg1 = (ideal) h->Data();
    19251924        intvec* arg2 = (intvec*) h->next->Data();
    19261925        intvec* arg3   =  (intvec*) h->next->next->Data();
    1927         ring arg4 = (ring) h->next->next->next->Data();
    1928 
    1929         ideal result = (ideal) Mwalk(arg1, arg2, arg3, arg4);
     1926
     1927
     1928        ideal result = (ideal) Mwalk(arg1, arg2, arg3);
    19301929
    19311930        res->rtyp = IDEAL_CMD;
     
    19761975      {
    19771976        if (h == NULL || h->Typ() != IDEAL_CMD ||
    1978             h->next == NULL || h->next->Typ() != INTVEC_CMD ||
    1979             h->next->next == NULL || h->next->next->Typ() != INTVEC_CMD ||
    1980             h->next->next->next == NULL || h->next->next->next->Typ() != INT_CMD ||
    1981             h->next->next->next->next == NULL || h->next->next->next->next->Typ() != INT_CMD||
    1982             h->next->next->next->next->next == NULL || h->next->next->next->next->next->Typ() != RING_CMD)
    1983         {
    1984           Werror("system(\"Mpwalk\", ideal, intvec, intvec, int, int, ring) expected");
    1985           return TRUE;
    1986         }
    1987 
    1988         if (((intvec*) h->next->Data())->length() != currRing->N &&
    1989             ((intvec*) h->next->next->Data())->length()!=currRing->N)
    1990         {
    1991           Werror("system(\"Mpwalk\" ...) intvecs not of length %d\n",currRing->N);
     1977            h->next == NULL || h->next->Typ() != INT_CMD ||
     1978            h->next->next == NULL || h->next->next->Typ() != INT_CMD ||
     1979            h->next->next->next == NULL ||
     1980              h->next->next->next->Typ() != INTVEC_CMD ||
     1981            h->next->next->next->next == NULL ||
     1982              h->next->next->next->next->Typ() != INTVEC_CMD||
     1983            h->next->next->next->next->next == NULL ||
     1984              h->next->next->next->next->next->Typ() != INT_CMD)
     1985        {
     1986          Werror("system(\"Mpwalk\", ideal, int, int, intvec, intvec, int) expected");
     1987          return TRUE;
     1988        }
     1989
     1990        if (((intvec*) h->next->next->next->Data())->length() != currRing->N &&
     1991            ((intvec*) h->next->next->next->next->Data())->length()!=currRing->N)
     1992        {
     1993          Werror("system(\"Mpwalk\" ...) intvecs not of length %d\n",
     1994                 currRing->N);
    19921995          return TRUE;
    19931996        }
    19941997        ideal arg1 = (ideal) h->Data();
    1995         intvec* arg2 = (intvec*) h->next->Data();
    1996         intvec* arg3 = (intvec*) h->next->next->Data();
    1997         int arg4 = (int) (long) h->next->next->next->Data();
    1998         int arg5 = (int) (long) h->next->next->next->next->Data();
    1999         ring arg6 = (ring) h->next->next->next->next->next->Data();
     1998        int arg2 = (int) ((long)(h->next->Data()));
     1999        int arg3 = (int) ((long)(h->next->next->Data()));
     2000        intvec* arg4 = (intvec*) h->next->next->next->Data();
     2001        intvec* arg5   =  (intvec*) h->next->next->next->next->Data();
     2002        int arg6   =  (int) ((long)(h->next->next->next->next->next->Data()));
     2003
    20002004
    20012005        ideal result = (ideal) Mpwalk(arg1, arg2, arg3, arg4, arg5, arg6);
     
    20082012      else
    20092013      if (strcmp(sys_cmd, "Mrwalk") == 0)
    2010       {
     2014      { // Random Walk
    20112015        if (h == NULL || h->Typ() != IDEAL_CMD ||
    20122016            h->next == NULL || h->next->Typ() != INTVEC_CMD ||
    20132017            h->next->next == NULL || h->next->next->Typ() != INTVEC_CMD ||
    20142018            h->next->next->next == NULL || h->next->next->next->Typ() != INT_CMD ||
    2015             h->next->next->next->next == NULL || h->next->next->next->next->Typ() != INT_CMD ||
    2016             h->next->next->next->next->next == NULL || h->next->next->next->next->next->Typ() != RING_CMD)
    2017         {
    2018           Werror("system(\"Mrwalk\", ideal, intvec, intvec, int, int, ring) expected");
    2019           return TRUE;
    2020         }
    2021 
    2022         if ((((intvec*) h->next->Data())->length() != currRing->N &&
    2023             ((intvec*) h->next->next->Data())->length() != currRing->N ) &&
    2024             (((intvec*) h->next->Data())->length() != (currRing->N)*(currRing->N) &&
    2025             ((intvec*) h->next->next->Data())->length() != (currRing->N)*(currRing->N) ))
    2026         {
    2027           Werror("system(\"Mrwalk\" ...) intvecs not of length %d or %d\n",
    2028                  currRing->N,(currRing->N)*(currRing->N));
     2019            h->next->next->next->next == NULL || h->next->next->next->next->Typ() != INT_CMD)
     2020        {
     2021          Werror("system(\"Mrwalk\", ideal, intvec, intvec, int, int) expected");
     2022          return TRUE;
     2023        }
     2024
     2025        if (((intvec*) h->next->Data())->length() != currRing->N &&
     2026            ((intvec*) h->next->next->Data())->length() != currRing->N )
     2027        {
     2028          Werror("system(\"Mrwalk\" ...) intvecs not of length %d\n",
     2029                 currRing->N);
    20292030          return TRUE;
    20302031        }
     
    20342035        int arg4 = (int)(long) h->next->next->next->Data();
    20352036        int arg5 = (int)(long) h->next->next->next->next->Data();
    2036         ring arg6 = (ring) h->next->next->next->next->next->Data();
    2037 
    2038 
    2039         ideal result = (ideal) Mrwalk(arg1, arg2, arg3, arg4, arg5, arg6);
     2037
     2038
     2039        ideal result = (ideal) Mrwalk(arg1, arg2, arg3, arg4, arg5);
    20402040
    20412041        res->rtyp = IDEAL_CMD;
     
    21182118        if (h == NULL || h->Typ() != IDEAL_CMD ||
    21192119            h->next == NULL || h->next->Typ() != INTVEC_CMD ||
    2120             h->next->next == NULL || h->next->next->Typ() != INTVEC_CMD ||
    2121             h->next->next->next == NULL || h->next->next->next->Typ() != INT_CMD ||
    2122             h->next->next->next->next == NULL || h->next->next->next->next->Typ() != RING_CMD)
    2123         {
    2124           Werror("system(\"Mfwalk\", ideal, intvec, intvec, int, ring) expected");
    2125           return TRUE;
    2126         }
     2120            h->next->next == NULL || h->next->next->Typ() != INTVEC_CMD)
     2121        {
     2122          Werror("system(\"Mfwalk\", ideal, intvec, intvec) expected");
     2123          return TRUE;
     2124        }
     2125
    21272126        if (((intvec*) h->next->Data())->length() != currRing->N &&
    2128             ((intvec*) h->next->next->Data())->length() != currRing->N)
    2129         {
    2130           Werror("system(\"Mfwalk\" ...) intvecs not of length %d\n",currRing->N);
    2131           return TRUE;
    2132         }
    2133         ideal arg1      =       (ideal) h->Data();
    2134         intvec* arg2    =       (intvec*) h->next->Data();
    2135         intvec* arg3    =       (intvec*) h->next->next->Data();
    2136         int arg4        =       (int)(long) h->next->next->next->Data();
    2137         ring arg5       =       (ring) h->next->next->next->next->Data();
    2138         ideal result    =       (ideal) Mfwalk(arg1, arg2, arg3, arg4, arg5);
     2127            ((intvec*) h->next->next->Data())->length() != currRing->N )
     2128        {
     2129          Werror("system(\"Mfwalk\" ...) intvecs not of length %d\n",
     2130                 currRing->N);
     2131          return TRUE;
     2132        }
     2133        ideal arg1 = (ideal) h->Data();
     2134        intvec* arg2 = (intvec*) h->next->Data();
     2135        intvec* arg3   =  (intvec*) h->next->next->Data();
     2136
     2137        ideal result = (ideal) Mfwalk(arg1, arg2, arg3);
    21392138
    21402139        res->rtyp = IDEAL_CMD;
     
    21452144      else
    21462145      if (strcmp(sys_cmd, "Mfrwalk") == 0)
     2146      {
     2147        if (h == NULL || h->Typ() != IDEAL_CMD ||
     2148            h->next == NULL || h->next->Typ() != INTVEC_CMD ||
     2149            h->next->next == NULL || h->next->next->Typ() != INTVEC_CMD ||
     2150            h->next->next->next == NULL || h->next->next->next->Typ() != INT_CMD)
     2151        {
     2152          Werror("system(\"Mfrwalk\", ideal, intvec, intvec, int) expected");
     2153          return TRUE;
     2154        }
     2155
     2156        if (((intvec*) h->next->Data())->length() != currRing->N &&
     2157            ((intvec*) h->next->next->Data())->length() != currRing->N )
     2158        {
     2159          Werror("system(\"Mfrwalk\" ...) intvecs not of length %d\n",
     2160                 currRing->N);
     2161          return TRUE;
     2162        }
     2163        ideal arg1 = (ideal) h->Data();
     2164        intvec* arg2 = (intvec*) h->next->Data();
     2165        intvec* arg3 = (intvec*) h->next->next->Data();
     2166        int arg4 = (int)(long) h->next->next->next->Data();
     2167
     2168        ideal result = (ideal) Mfrwalk(arg1, arg2, arg3, arg4);
     2169
     2170        res->rtyp = IDEAL_CMD;
     2171        res->data =  result;
     2172
     2173        return FALSE;
     2174      }
     2175      else
     2176
     2177  #ifdef TRAN_Orig
     2178      if (strcmp(sys_cmd, "TranMImprovwalk") == 0)
     2179      {
     2180        if (h == NULL || h->Typ() != IDEAL_CMD ||
     2181            h->next == NULL || h->next->Typ() != INTVEC_CMD ||
     2182            h->next->next == NULL || h->next->next->Typ() != INTVEC_CMD)
     2183        {
     2184          Werror("system(\"TranMImprovwalk\", ideal, intvec, intvec) expected");
     2185          return TRUE;
     2186        }
     2187
     2188        if (((intvec*) h->next->Data())->length() != currRing->N &&
     2189            ((intvec*) h->next->next->Data())->length() != currRing->N )
     2190        {
     2191          Werror("system(\"TranMImprovwalk\" ...) intvecs not of length %d\n",
     2192                 currRing->N);
     2193          return TRUE;
     2194        }
     2195        ideal arg1 = (ideal) h->Data();
     2196        intvec* arg2 = (intvec*) h->next->Data();
     2197        intvec* arg3   =  (intvec*) h->next->next->Data();
     2198
     2199
     2200        ideal result = (ideal) TranMImprovwalk(arg1, arg2, arg3);
     2201
     2202        res->rtyp = IDEAL_CMD;
     2203        res->data =  result;
     2204
     2205        return FALSE;
     2206      }
     2207      else
     2208  #endif
     2209      if (strcmp(sys_cmd, "MAltwalk2") == 0)
     2210        {
     2211        if (h == NULL || h->Typ() != IDEAL_CMD ||
     2212            h->next == NULL || h->next->Typ() != INTVEC_CMD ||
     2213            h->next->next == NULL || h->next->next->Typ() != INTVEC_CMD)
     2214        {
     2215          Werror("system(\"MAltwalk2\", ideal, intvec, intvec) expected");
     2216          return TRUE;
     2217        }
     2218
     2219        if (((intvec*) h->next->Data())->length() != currRing->N &&
     2220            ((intvec*) h->next->next->Data())->length() != currRing->N )
     2221        {
     2222          Werror("system(\"MAltwalk2\" ...) intvecs not of length %d\n",
     2223                 currRing->N);
     2224          return TRUE;
     2225        }
     2226        ideal arg1 = (ideal) h->Data();
     2227        intvec* arg2 = (intvec*) h->next->Data();
     2228        intvec* arg3   =  (intvec*) h->next->next->Data();
     2229
     2230
     2231        ideal result = (ideal) MAltwalk2(arg1, arg2, arg3);
     2232
     2233        res->rtyp = IDEAL_CMD;
     2234        res->data =  result;
     2235
     2236        return FALSE;
     2237      }
     2238      else
     2239      if (strcmp(sys_cmd, "TranMImprovwalk") == 0)
     2240      {
     2241        if (h == NULL || h->Typ() != IDEAL_CMD ||
     2242            h->next == NULL || h->next->Typ() != INTVEC_CMD ||
     2243            h->next->next == NULL || h->next->next->Typ() != INTVEC_CMD||
     2244            h->next->next->next == NULL || h->next->next->next->Typ() != INT_CMD)
     2245        {
     2246          Werror("system(\"TranMImprovwalk\", ideal, intvec, intvec, int) expected");
     2247          return TRUE;
     2248        }
     2249
     2250        if (((intvec*) h->next->Data())->length() != currRing->N &&
     2251            ((intvec*) h->next->next->Data())->length() != currRing->N )
     2252        {
     2253          Werror("system(\"TranMImprovwalk\" ...) intvecs not of length %d\n",
     2254                 currRing->N);
     2255          return TRUE;
     2256        }
     2257        ideal arg1 = (ideal) h->Data();
     2258        intvec* arg2 = (intvec*) h->next->Data();
     2259        intvec* arg3   =  (intvec*) h->next->next->Data();
     2260        int arg4   =  (int) ((long)(h->next->next->next->Data()));
     2261
     2262        ideal result = (ideal) TranMImprovwalk(arg1, arg2, arg3, arg4);
     2263
     2264        res->rtyp = IDEAL_CMD;
     2265        res->data =  result;
     2266
     2267        return FALSE;
     2268      }
     2269      else
     2270      if (strcmp(sys_cmd, "TranMrImprovwalk") == 0)
    21472271      {
    21482272        if (h == NULL || h->Typ() != IDEAL_CMD ||
     
    21502274            h->next->next == NULL || h->next->next->Typ() != INTVEC_CMD ||
    21512275            h->next->next->next == NULL || h->next->next->next->Typ() != INT_CMD ||
    2152             h->next->next->next->next == NULL || h->next->next->next->next->Typ() != INT_CMD ||
    2153             h->next->next->next->next->next == NULL || h->next->next->next->next->next->Typ() != RING_CMD)
    2154         {
    2155           Werror("system(\"Mfrwalk\", ideal, intvec, intvec, int, int, ring) expected");
    2156           return TRUE;
    2157         }
     2276            h->next->next->next == NULL || h->next->next->next->next->Typ() != INT_CMD ||
     2277            h->next->next->next == NULL || h->next->next->next->next->next->Typ() != INT_CMD)
     2278        {
     2279          Werror("system(\"TranMrImprovwalk\", ideal, intvec, intvec) expected");
     2280          return TRUE;
     2281        }
     2282
    21582283        if (((intvec*) h->next->Data())->length() != currRing->N &&
    2159             ((intvec*) h->next->next->Data())->length() != currRing->N)
    2160         {
    2161           Werror("system(\"Mfrwalk\" ...) intvecs not of length %d\n",currRing->N);
     2284            ((intvec*) h->next->next->Data())->length() != currRing->N )
     2285        {
     2286          Werror("system(\"TranMrImprovwalk\" ...) intvecs not of length %d\n", currRing->N);
    21622287          return TRUE;
    21632288        }
     
    21672292        int arg4 = (int)(long) h->next->next->next->Data();
    21682293        int arg5 = (int)(long) h->next->next->next->next->Data();
    2169         ring arg6 = (ring) h->next->next->next->next->next->Data();
    2170 
    2171         ideal result = (ideal) Mfrwalk(arg1, arg2, arg3, arg4, arg5, arg6);
     2294        int arg6 = (int)(long) h->next->next->next->next->next->Data();
     2295
     2296        ideal result = (ideal) TranMrImprovwalk(arg1, arg2, arg3, arg4, arg5, arg6);
    21722297
    21732298        res->rtyp = IDEAL_CMD;
     
    21772302      }
    21782303      else
    2179       if (strcmp(sys_cmd, "Mprwalk") == 0)
    2180       {
    2181         if (h == NULL || h->Typ() != IDEAL_CMD ||
    2182             h->next == NULL || h->next->Typ() != INTVEC_CMD ||
    2183             h->next->next == NULL || h->next->next->Typ() != INTVEC_CMD ||
    2184             h->next->next->next == NULL || h->next->next->next->Typ() != INT_CMD ||
    2185             h->next->next->next->next == NULL || h->next->next->next->next->Typ() != INT_CMD ||
    2186             h->next->next->next->next->next == NULL || h->next->next->next->next->next->Typ() != INT_CMD ||
    2187             h->next->next->next->next->next->next == NULL || h->next->next->next->next->next->next->Typ() != RING_CMD)
    2188         {
    2189           Werror("system(\"Mprwalk\", ideal, intvec, intvec, int, int, int, ring) expected");
    2190           return TRUE;
    2191         }
    2192         if (((intvec*) h->next->Data())->length() != currRing->N &&
    2193             ((intvec*) h->next->next->Data())->length() != currRing->N )
    2194         {
    2195           Werror("system(\"Mrwalk\" ...) intvecs not of length %d\n",
    2196                  currRing->N);
    2197           return TRUE;
    2198         }
    2199         ideal arg1 = (ideal) h->Data();
    2200         intvec* arg2 = (intvec*) h->next->Data();
    2201         intvec* arg3 =  (intvec*) h->next->next->Data();
    2202         int arg4 = (int)(long) h->next->next->next->Data();
    2203         int arg5 = (int)(long) h->next->next->next->next->Data();
    2204         int arg6 = (int)(long) h->next->next->next->next->next->Data();
    2205         ring arg7 = (ring) h->next->next->next->next->next->next->Data();
    2206 
    2207 
    2208         ideal result = (ideal) Mprwalk(arg1, arg2, arg3, arg4, arg5, arg6, arg7);
    2209 
    2210         res->rtyp = IDEAL_CMD;
    2211         res->data =  result;
    2212 
    2213         return FALSE;
    2214       }
    2215       else
    2216 
    2217   #ifdef TRAN_Orig
    2218       if (strcmp(sys_cmd, "TranMImprovwalk") == 0)
    2219       {
    2220         if (h == NULL || h->Typ() != IDEAL_CMD ||
    2221             h->next == NULL || h->next->Typ() != INTVEC_CMD ||
    2222             h->next->next == NULL || h->next->next->Typ() != INTVEC_CMD)
    2223         {
    2224           Werror("system(\"TranMImprovwalk\", ideal, intvec, intvec) expected");
    2225           return TRUE;
    2226         }
    2227 
    2228         if (((intvec*) h->next->Data())->length() != currRing->N &&
    2229             ((intvec*) h->next->next->Data())->length() != currRing->N )
    2230         {
    2231           Werror("system(\"TranMImprovwalk\" ...) intvecs not of length %d\n",
    2232                  currRing->N);
    2233           return TRUE;
    2234         }
    2235         ideal arg1 = (ideal) h->Data();
    2236         intvec* arg2 = (intvec*) h->next->Data();
    2237         intvec* arg3   =  (intvec*) h->next->next->Data();
    2238 
    2239 
    2240         ideal result = (ideal) TranMImprovwalk(arg1, arg2, arg3);
    2241 
    2242         res->rtyp = IDEAL_CMD;
    2243         res->data =  result;
    2244 
    2245         return FALSE;
    2246       }
    2247       else
    2248   #endif
    2249       if (strcmp(sys_cmd, "MAltwalk2") == 0)
    2250         {
    2251         if (h == NULL || h->Typ() != IDEAL_CMD ||
    2252             h->next == NULL || h->next->Typ() != INTVEC_CMD ||
    2253             h->next->next == NULL || h->next->next->Typ() != INTVEC_CMD)
    2254         {
    2255           Werror("system(\"MAltwalk2\", ideal, intvec, intvec) expected");
    2256           return TRUE;
    2257         }
    2258 
    2259         if (((intvec*) h->next->Data())->length() != currRing->N &&
    2260             ((intvec*) h->next->next->Data())->length() != currRing->N )
    2261         {
    2262           Werror("system(\"MAltwalk2\" ...) intvecs not of length %d\n",
    2263                  currRing->N);
    2264           return TRUE;
    2265         }
    2266         ideal arg1 = (ideal) h->Data();
    2267         intvec* arg2 = (intvec*) h->next->Data();
    2268         intvec* arg3   =  (intvec*) h->next->next->Data();
    2269 
    2270 
    2271         ideal result = (ideal) MAltwalk2(arg1, arg2, arg3);
    2272 
    2273         res->rtyp = IDEAL_CMD;
    2274         res->data =  result;
    2275 
    2276         return FALSE;
    2277       }
    2278       else
    2279       if (strcmp(sys_cmd, "TranMImprovwalk") == 0)
    2280       {
    2281         if (h == NULL || h->Typ() != IDEAL_CMD ||
    2282             h->next == NULL || h->next->Typ() != INTVEC_CMD ||
    2283             h->next->next == NULL || h->next->next->Typ() != INTVEC_CMD||
    2284             h->next->next->next == NULL || h->next->next->next->Typ() != INT_CMD)
    2285         {
    2286           Werror("system(\"TranMImprovwalk\", ideal, intvec, intvec, int) expected");
    2287           return TRUE;
    2288         }
    2289 
    2290         if (((intvec*) h->next->Data())->length() != currRing->N &&
    2291             ((intvec*) h->next->next->Data())->length() != currRing->N )
    2292         {
    2293           Werror("system(\"TranMImprovwalk\" ...) intvecs not of length %d\n",
    2294                  currRing->N);
    2295           return TRUE;
    2296         }
    2297         ideal arg1 = (ideal) h->Data();
    2298         intvec* arg2 = (intvec*) h->next->Data();
    2299         intvec* arg3   =  (intvec*) h->next->next->Data();
    2300         int arg4   =  (int) ((long)(h->next->next->next->Data()));
    2301 
    2302         ideal result = (ideal) TranMImprovwalk(arg1, arg2, arg3, arg4);
    2303 
    2304         res->rtyp = IDEAL_CMD;
    2305         res->data =  result;
    2306 
    2307         return FALSE;
    2308       }
    2309       else
    2310 #endif
     2304
     2305  #endif
    23112306  /*================= Extended system call ========================*/
    23122307     {
     
    23242319#ifdef HAVE_EXTENDED_SYSTEM
    23252320  // You can put your own system calls here
    2326 #  include <kernel/fglmcomb.cc>
    2327 #  include <kernel/fglm.h>
     2321#  include <kernel/fglm/fglmcomb.cc>
     2322#  include <kernel/fglm/fglm.h>
    23282323#  ifdef HAVE_NEWTON
    23292324#    include <hc_newton.h>
     
    23312326#  include <polys/mod_raw.h>
    23322327#  include <polys/monomials/ring.h>
    2333 #  include <kernel/shiftgb.h>
     2328#  include <kernel/GBEngine/shiftgb.h>
    23342329
    23352330static BOOLEAN jjEXTENDED_SYSTEM(leftv res, leftv h)
     
    26872682//       }
    26882683//       else
    2689   /*==================== isSqrFree =============================*/
    2690       if(strcmp(sys_cmd,"isSqrFree")==0)
    2691       {
    2692         if ((h!=NULL) &&(h->Typ()==POLY_CMD))
    2693         {
    2694           res->rtyp=INT_CMD;
    2695           res->data=(void *)(long) singclap_isSqrFree((poly)h->Data(), currRing);
    2696           return FALSE;
    2697         }
    2698         else
    2699           WerrorS("poly expected");
    2700       }
    2701       else
    27022684  /*==================== pDivStat =============================*/
    27032685  #if defined(PDEBUG) || defined(PDIV_DEBUG)
     
    31603142      else
    31613143  #endif
     3144  /*==================== Roune Hilb  =================*/
     3145       if (strcmp(sys_cmd, "hilbroune") == 0)
     3146       {
     3147         ideal I;
     3148         if ((h!=NULL) && (h->Typ()==IDEAL_CMD))
     3149         {
     3150           I=(ideal)h->CopyD();
     3151           slicehilb(I);
     3152         }
     3153         else return TRUE;
     3154         return FALSE;
     3155       }
    31623156  /*==================== minor =================*/
    31633157      if (strcmp(sys_cmd, "minor")==0)
     
    35973591        {
    35983592#ifdef HAVE_PLURAL
    3599           Print("NTL_0:%d (use NTL for gcd of polynomials in char 0)\n",isOn(SW_USE_NTL_GCD_0));
    3600           Print("NTL_p:%d (use NTL for gcd of polynomials in char p)\n",isOn(SW_USE_NTL_GCD_P));
    36013593          Print("EZGCD:%d (use EZGCD for gcd of polynomials in char 0)\n",isOn(SW_USE_EZGCD));
    36023594          Print("EZGCD_P:%d (use EZGCD_P for gcd of polynomials in char p)\n",isOn(SW_USE_EZGCD_P));
     
    36143606          char *s=(char *)h->Data();
    36153607#ifdef HAVE_PLURAL
    3616           if (strcmp(s,"NTL_0")==0) { if (d) On(SW_USE_NTL_GCD_0); else Off(SW_USE_NTL_GCD_0); } else
    3617           if (strcmp(s,"NTL_p")==0) { if (d) On(SW_USE_NTL_GCD_P); else Off(SW_USE_NTL_GCD_P); } else
    36183608          if (strcmp(s,"EZGCD")==0) { if (d) On(SW_USE_EZGCD); else Off(SW_USE_EZGCD); } else
    36193609          if (strcmp(s,"EZGCD_P")==0) { if (d) On(SW_USE_EZGCD_P); else Off(SW_USE_EZGCD_P); } else
     
    38783868  if (strcmp(sys_cmd,"reservedLink")==0)
    38793869  {
    3880     si_link ssiCommandLink();
     3870    extern si_link ssiCommandLink();
    38813871    res->rtyp=LINK_CMD;
    38823872    si_link p=ssiCommandLink();
  • Singular/walk.cc

    • Property mode changed from 100755 to 100644
    rfde66a9 ra5ab27  
    1818//#define CHECK_IDEAL
    1919//#define CHECK_IDEAL_MWALK
    20 //#define CHECK_IDEAL_MFRWALK
    2120
    2221//#define NEXT_VECTORS_CC
     
    3635/* includes */
    3736
    38 #include <stdio.h>
    39 #include <stdlib.h>
    40 // === Zeit & System (Holger Croeni ===
    41 #include <time.h>
    42 #include <sys/time.h>
    43 #include <math.h>
    44 #include <sys/stat.h>
    45 #include <unistd.h>
    46 #include <float.h>
    47 #include <misc/mylimits.h>
    48 #include <sys/types.h>
    49 
    50 
    51 //#include "config.h"
    5237#include <kernel/mod2.h>
    5338#include <misc/intvec.h>
     
    8873#include <polys/monomials/ring.h>
    8974//#include <polys/ext_fields/longalg.h>
    90 #ifdef HAVE_FACTORY
    9175#include <polys/clapsing.h>
    92 #endif
    9376
    9477#include <coeffs/mpr_complex.h>
     78
     79#include <stdio.h>
     80// === Zeit & System (Holger Croeni ===
     81#include <time.h>
     82#include <sys/time.h>
     83#include <math.h>
     84#include <sys/stat.h>
     85#include <unistd.h>
     86#include <float.h>
     87#include <misc/mylimits.h>
     88#include <sys/types.h>
    9589
    9690int nstep;
     
    124118}
    125119
     120/************************************
     121 * construct the set s from F u {P} *
     122 ************************************/
     123// unused
     124#if 0
     125static void initSSpecialCC (ideal F, ideal Q, ideal P,kStrategy strat)
     126{
     127  int   i,pos;
     128
     129  if (Q!=NULL) i=((IDELEMS(Q)+(setmaxTinc-1))/setmaxTinc)*setmaxTinc;
     130  else i=setmaxT;
     131
     132  strat->ecartS=initec(i);
     133  strat->sevS=initsevS(i);
     134  strat->S_2_R=initS_2_R(i);
     135  strat->fromQ=NULL;
     136  strat->Shdl=idInit(i,F->rank);
     137  strat->S=strat->Shdl->m;
     138
     139  // - put polys into S -
     140  if (Q!=NULL)
     141  {
     142    strat->fromQ=initec(i);
     143    memset(strat->fromQ,0,i*sizeof(int));
     144    for (i=0; i<IDELEMS(Q); i++)
     145    {
     146      if (Q->m[i]!=NULL)
     147      {
     148        LObject h;
     149        h.p = pCopy(Q->m[i]);
     150        //if (TEST_OPT_INTSTRATEGY)
     151        //{
     152        //  //pContent(h.p);
     153        //  h.pCleardenom(); // also does a pContent
     154        //}
     155        //else
     156        //{
     157        //  h.pNorm();
     158        //}
     159        strat->initEcart(&h);
     160        if (rHasLocalOrMixedOrdering_currRing())
     161        {
     162          deleteHC(&h,strat);
     163        }
     164        if (h.p!=NULL)
     165        {
     166          if (strat->sl==-1)
     167            pos =0;
     168          else
     169          {
     170            pos = posInS(strat,strat->sl,h.p,h.ecart);
     171          }
     172          h.sev = pGetShortExpVector(h.p);
     173          h.SetpFDeg();
     174          strat->enterS(h,pos,strat, strat->tl+1);
     175          enterT(h, strat);
     176          strat->fromQ[pos]=1;
     177        }
     178      }
     179    }
     180  }
     181  //- put polys into S -
     182  for (i=0; i<IDELEMS(F); i++)
     183  {
     184    if (F->m[i]!=NULL)
     185    {
     186      LObject h;
     187      h.p = pCopy(F->m[i]);
     188      if (rHasGlobalOrdering(currRing))
     189      {
     190        //h.p=redtailBba(h.p,strat->sl,strat);
     191        h.p=redtailBba(h.p,strat->sl,strat);
     192      }
     193      else
     194      {
     195        deleteHC(&h,strat);
     196      }
     197      strat->initEcart(&h);
     198      if (h.p!=NULL)
     199      {
     200        if (strat->sl==-1)
     201          pos =0;
     202        else
     203          pos = posInS(strat,strat->sl,h.p,h.ecart);
     204        h.sev = pGetShortExpVector(h.p);
     205        strat->enterS(h,pos,strat, strat->tl+1);
     206        h.length = pLength(h.p);
     207        h.SetpFDeg();
     208        enterT(h,strat);
     209      }
     210    }
     211  }
     212#ifdef INITSSPECIAL
     213  for (i=0; i<IDELEMS(P); i++)
     214  {
     215    if (P->m[i]!=NULL)
     216    {
     217      LObject h;
     218      h.p=pCopy(P->m[i]);
     219      strat->initEcart(&h);
     220      h.length = pLength(h.p);
     221      if (TEST_OPT_INTSTRATEGY)
     222      {
     223        h.pCleardenom();
     224      }
     225      else
     226      {
     227        h.pNorm();
     228      }
     229      if(strat->sl>=0)
     230      {
     231        if (rHasGlobalOrdering(currRing))
     232        {
     233          h.p=redBba(h.p,strat->sl,strat);
     234          if (h.p!=NULL)
     235            h.p=redtailBba(h.p,strat->sl,strat);
     236        }
     237        else
     238        {
     239          h.p=redMora(h.p,strat->sl,strat);
     240          strat->initEcart(&h);
     241        }
     242        if(h.p!=NULL)
     243        {
     244          if (TEST_OPT_INTSTRATEGY)
     245          {
     246            h.pCleardenom();
     247          }
     248          else
     249          {
     250            h.is_normalized = 0;
     251            h.pNorm();
     252          }
     253          h.sev = pGetShortExpVector(h.p);
     254          h.SetpFDeg();
     255          pos = posInS(strat->S,strat->sl,h.p,h.ecart);
     256          enterpairsSpecial(h.p,strat->sl,h.ecart,pos,strat,strat->tl+1);
     257          strat->enterS(h,pos,strat, strat->tl+1);
     258          enterT(h,strat);
     259        }
     260      }
     261      else
     262      {
     263        h.sev = pGetShortExpVector(h.p);
     264        h.SetpFDeg();
     265        strat->enterS(h,0,strat, strat->tl+1);
     266        enterT(h,strat);
     267      }
     268    }
     269  }
     270#endif
     271}
     272#endif
     273
    126274/*****************
    127275 *interreduce F  *
     
    132280  kStrategy strat = new skStrategy;
    133281
     282//  if (TEST_OPT_PROT)
     283//  {
     284//    writeTime("start InterRed:");
     285//    mflush();
     286//  }
    134287  //strat->syzComp     = 0;
    135288  strat->kHEdgeFound = (currRing->ppNoether) != NULL;
     
    159312  initS(F,Q,strat);
    160313
     314  /*
     315  timetmp=clock();//22.01.02
     316  initSSpecialCC(F,Q,NULL,strat);
     317  tininitS=tininitS+clock()-timetmp;//22.01.02
     318  */
    161319  if(TEST_OPT_REDSB)
    162320  {
     
    190348    strat->fromQ = NULL;
    191349  }
     350//  if (TEST_OPT_PROT)
     351//  {
     352//    writeTime("end Interred:");
     353//    mflush();
     354//  }
    192355  ideal shdl=strat->Shdl;
    193356  idSkipZeroes(shdl);
     
    197360}
    198361
    199 #ifdef TIME_TEST
     362//unused
     363#if 0
    200364static void TimeString(clock_t tinput, clock_t tostd, clock_t tif,clock_t tstd,
    201365                       clock_t tlf,clock_t tred, clock_t tnw, int step)
     
    237401#endif
    238402
    239 #ifdef TIME_TEST
     403//unused
     404#if 0
    240405static void TimeStringFractal(clock_t tinput, clock_t tostd, clock_t tif,clock_t tstd,
    241406                       clock_t textra, clock_t tlf,clock_t tred, clock_t tnw)
     
    266431#endif
    267432
     433#ifdef CHECK_IDEAL_MWALK
    268434static void idString(ideal L, const char* st)
    269435{
     
    277443  Print(" %s;", pString(L->m[nL-1]));
    278444}
    279 
     445#endif
     446
     447#if defined(CHECK_IDEAL_MWALK) || defined(ENDWALKS)
    280448static void headidString(ideal L, char* st)
    281449{
     
    289457  Print(" %s;", pString(pHead(L->m[nL-1])));
    290458}
    291 
     459#endif
     460
     461#if defined(CHECK_IDEAL_MWALK) || defined(ENDWALKS)
    292462static void idElements(ideal L, char* st)
    293463{
     
    301471  }
    302472  int j, nsame;
     473  // int  nk=0;
    303474  for(i=0; i<nL; i++)
    304475  {
     
    326497  omFree(K);
    327498}
     499#endif
     500
    328501
    329502static void ivString(intvec* iv, const char* ch)
     
    339512}
    340513
     514//unused
     515//#if 0
    341516static void MivString(intvec* iva, intvec* ivb, intvec* ivc)
    342517{
     
    361536  Print("%d)", (*ivc)[nV]);
    362537}
     538//#endif
    363539
    364540/********************************************************************
     
    401577  mpz_clear(g);
    402578}
     579
     580//unused
     581#if 0
     582static int isVectorNeg(intvec* omega)
     583{
     584  int i;
     585
     586  for(i=omega->length(); i>=0; i--)
     587  {
     588    if((*omega)[i]<0)
     589    {
     590      return 1;
     591    }
     592  }
     593  return 0;
     594}
     595#endif
    403596
    404597/********************************************************************
     
    434627    {
    435628      PrintLn();
    436       PrintS("\n//** MLmWeightedDegree: OVERFLOW: ");
     629      PrintS("\n// ** OVERFLOW in \"MwalkInitialForm\": ");
    437630      mpz_out_str( stdout, 10, zsum);
    438       PrintS(" is greater than 2147483647 (max. integer representation).\n");
     631      PrintS(" is greater than 2147483647 (max. integer representation)");
    439632      Overflow_Error = TRUE;
    440633    }
     
    571764  if(G->m[0] == NULL)
    572765  {
    573     idString(G,"G");
    574     WerrorS("//** test_w_in_ConeCC: the result may be wrong, i.e. equal zero.\n");
     766    PrintS("//** the result may be WRONG, i.e. 0!!\n");
     767    return 0;
    575768  }
    576769
     
    762955}
    763956
    764 /*****************************************************************************
    765 * create a weight matrix order as intvec of an extra weight vector (a(iv),lp)*
    766 ******************************************************************************/
    767 intvec* MivMatrixOrderRefine(intvec* iv, intvec* iw)
    768 {
    769   assume(iv->length() == iw->length());
    770   int i, nR = iv->length();
    771 
    772   intvec* ivm = new intvec(nR*nR);
    773 
    774   for(i=0; i<nR; i++)
    775   {
    776     (*ivm)[i] = (*iv)[i];
    777     (*ivm)[i+nR] = (*iw)[i];
    778   }
    779   for(i=2; i<nR; i++)
    780   {
    781     (*ivm)[i*nR+i-2] = 1;
    782   }
    783   return ivm;
    784 }
    785 
    786957/*******************************
    787958 * return intvec = (1, ..., 1) *
     
    810981}
    811982
    812 /**********************************************
    813  * Look for the smallest absolut value in vec *
    814  **********************************************/
    815 static int MivAbsMax(intvec* vec)
    816 {
    817   int i,k;
    818   if((*vec)[0] < 0)
    819   {
    820     k = -(*vec)[0];
    821   }
    822   else
    823   {
    824     k = (*vec)[0];
    825   }
    826   for(i=1; i < (vec->length()); i++)
    827   {
    828     if((*vec)[i] < 0)
    829     {
    830       if(-(*vec)[i] > k)
    831       {
    832         k = -(*vec)[i];
    833       }
    834     }
    835     else
    836     {
    837       if((*vec)[i] > k)
    838       {
    839         k = (*vec)[i];
    840       }
    841     }
    842   }
    843   return k;
    844 }
     983//unused
     984/*****************************************************************************
     985 * print the max total degree and the max coefficient of G                   *
     986 *****************************************************************************/
     987#if 0
     988static void checkComplexity(ideal G, char* cG)
     989{
     990  int nV = currRing->N;
     991  int nG = IDELEMS(G);
     992  intvec* ivUnit = Mivdp(nV);
     993  int i, tmpdeg, maxdeg=0;
     994  number tmpcoeff , maxcoeff=currRing->cf->nNULL;
     995  poly p;
     996  for(i=nG-1; i>=0; i--)
     997  {
     998    tmpdeg = MwalkWeightDegree(G->m[i], ivUnit);
     999    if(tmpdeg > maxdeg )
     1000    {
     1001      maxdeg = tmpdeg;
     1002    }
     1003  }
     1004
     1005  for(i=nG-1; i>=0; i--)
     1006  {
     1007    p = pCopy(G->m[i]);
     1008    while(p != NULL)
     1009    {
     1010      //tmpcoeff = pGetCoeff(pHead(p));
     1011      tmpcoeff = pGetCoeff(p);
     1012      if(nGreater(tmpcoeff,maxcoeff))
     1013      {
     1014         maxcoeff = nCopy(tmpcoeff);
     1015      }
     1016      pIter(p);
     1017    }
     1018    pDelete(&p);
     1019  }
     1020  p = pNSet(maxcoeff);
     1021  char* pStr = pString(p);
     1022  delete ivUnit;
     1023  Print("// max total degree of %s = %d\n",cG, maxdeg);
     1024  Print("// max coefficient of %s  = %s", cG, pStr);//ing(p));
     1025  Print(" which consists of %d digits", (int)strlen(pStr));
     1026  PrintLn();
     1027}
     1028#endif
     1029
    8451030/*****************************************************************************
    8461031* If target_ord = intmat(A1, ..., An) then calculate the perturbation        *
     
    8561041intvec* MPertVectors(ideal G, intvec* ivtarget, int pdeg)
    8571042{
    858   int nV = currRing->N, i=0, j=0, nG = IDELEMS(G);
     1043  // ivtarget is a matrix order of a degree reverse lex. order
     1044  int nV = currRing->N;
     1045  //assume(pdeg <= nV && pdeg >= 0);
     1046
     1047  int i, j, nG = IDELEMS(G);
    8591048  intvec* v_null =  new intvec(nV);
     1049
    8601050
    8611051  // Check that the perturbed degree is valid
    8621052  if(pdeg > nV || pdeg <= 0)
    8631053  {
    864     WerrorS("\n//** MPertVectors: The perturbed degree is wrong.\n");
     1054    WerrorS("//** The perturbed degree is wrong!!");
    8651055    return v_null;
    8661056  }
     
    8721062  }
    8731063  mpz_t *pert_vector = (mpz_t*)omAlloc(nV*sizeof(mpz_t));
     1064  //mpz_t *pert_vector1 = (mpz_t*)omAlloc(nV*sizeof(mpz_t));
    8741065
    8751066  for(i=0; i<nV; i++)
    8761067  {
    8771068    mpz_init_set_si(pert_vector[i], (*ivtarget)[i]);
    878   }
    879   // Compute max1 = Max(A_2)+Max(A_3)+...+Max(A_pdeg), where the A_i is the i-th row of the matrix target_ord.
     1069   // mpz_init_set_si(pert_vector1[i], (*ivtarget)[i]);
     1070  }
     1071  // Calculate max1 = Max(A2)+Max(A3)+...+Max(Apdeg),
     1072  // where the Ai are the i-te rows of the matrix target_ord.
    8801073  int ntemp, maxAi, maxA=0;
    8811074  for(i=1; i<pdeg; i++)
     
    9011094  }
    9021095
    903   // Compute inveps = 1/eps, where 1/eps > totaldeg(p)*max1 for all p in G.
     1096  // Calculate inveps = 1/eps, where 1/eps > totaldeg(p)*max1 for all p in G.
     1097
     1098  intvec* ivUnit = Mivdp(nV);
     1099
    9041100  mpz_t tot_deg; mpz_init(tot_deg);
    9051101  mpz_t maxdeg; mpz_init(maxdeg);
    9061102  mpz_t inveps; mpz_init(inveps);
    9071103
     1104
    9081105  for(i=nG-1; i>=0; i--)
    9091106  {
    910      mpz_set_ui(maxdeg, MwalkWeightDegree(G->m[i], Mivdp(nV)));
     1107     mpz_set_ui(maxdeg, MwalkWeightDegree(G->m[i], ivUnit));
    9111108     if (mpz_cmp(maxdeg,  tot_deg) > 0 )
    9121109     {
     
    9141111     }
    9151112  }
     1113
     1114  delete ivUnit;
    9161115  mpz_mul_ui(inveps, tot_deg, maxA);
    9171116  mpz_add_ui(inveps, inveps, 1);
     
    9221121  if(mpz_cmp_ui(inveps, pdeg)>0 && pdeg > 3)
    9231122  {
     1123    //  Print("\n// choose the\"small\" inverse epsilon := %d / %d = ", mpz_get_si(inveps), pdeg);
    9241124    mpz_fdiv_q_ui(inveps, inveps, pdeg);
     1125    // mpz_out_str(stdout, 10, inveps);
    9251126  }
    9261127#else
    927   PrintS("\n//** MPertVectors: the \"big\" inverse epsilon: ");
     1128  // PrintS("\n// the \"big\" inverse epsilon: ");
    9281129  mpz_out_str(stdout, 10, inveps);
    9291130#endif
    9301131
    931   // pert(A1) = inveps^(pdeg-1)*A1 + inveps^(pdeg-2)*A2+...+A_pdeg, pert_vector := A1
    932   for(i=1; i<pdeg; i++)
     1132  // pert(A1) = inveps^(pdeg-1)*A1 + inveps^(pdeg-2)*A2+...+A_pdeg,
     1133  // pert_vector := A1
     1134  for( i=1; i < pdeg; i++ )
    9331135  {
    9341136    for(j=0; j<nV; j++)
     
    9651167
    9661168  intvec *pert_vector1= new intvec(nV);
    967 
     1169  j = 0;
    9681170  for(i=0; i<nV; i++)
    9691171  {
    9701172    (* pert_vector1)[i] = mpz_get_si(pert_vector[i]);
    971   }
    972   while(MivAbsMax(pert_vector1)>100000)
    973   {
    974     j = 0;
     1173    (* pert_vector1)[i] = 0.1*(* pert_vector1)[i];
     1174    (* pert_vector1)[i] = floor((* pert_vector1)[i] + 0.5);
     1175    if((* pert_vector1)[i] == 0)
     1176    {
     1177      j++;
     1178    }
     1179  }
     1180  if(j > nV - 1)
     1181  {
     1182    // Print("\n//  MPertVectors: geaenderter vector gleich Null! \n");
     1183    delete pert_vector1;
     1184    goto CHECK_OVERFLOW;
     1185  }
     1186
     1187// check that the perturbed weight vector lies in the Groebner cone
     1188  if(test_w_in_ConeCC(G,pert_vector1) != 0)
     1189  {
     1190    // Print("\n//  MPertVectors: geaenderter vector liegt in Groebnerkegel! \n");
    9751191    for(i=0; i<nV; i++)
    9761192    {
    977       (* pert_vector1)[i] = 0.1*(* pert_vector1)[i];
    978       (* pert_vector1)[i] = floor((* pert_vector1)[i] + 0.5);
    979       if((* pert_vector1)[i] == 0)
    980       {
    981         j++;
    982       }
    983     }
    984     if(j > nV - 1)
    985     {
    986       break;
    987     }
    988 
    989 // check that the perturbed weight vector lies in the Groebner cone
    990     if(test_w_in_ConeCC(G,pert_vector1) != 0)
    991     {
    992       for(i=0; i<nV; i++)
    993       {
    994         mpz_set_si(pert_vector[i], (*pert_vector1)[i]);
    995       }
    996     }
    997     else
    998     {
    999       break;
    1000     }
     1193      mpz_set_si(pert_vector[i], (*pert_vector1)[i]);
     1194    }
     1195  }
     1196  else
     1197  {
     1198    //Print("\n// MpertVectors: geaenderter vector liegt nicht in Groebnerkegel! \n");
    10011199  }
    10021200  delete pert_vector1;
    10031201
    1004   //CHECK_OVERFLOW:
     1202  CHECK_OVERFLOW:
    10051203  intvec* result = new intvec(nV);
    10061204
    1007 // 2147483647 is max. integer representation in SINGULAR
     1205  /* 2147483647 is max. integer representation in SINGULAR */
    10081206  mpz_t sing_int;
    10091207  mpz_init_set_ui(sing_int,  2147483647);
    10101208
     1209  int ntrue=0;
    10111210  for(i=0; i<nV; i++)
    10121211  {
     
    10141213    if(mpz_cmp(pert_vector[i], sing_int)>=0)
    10151214    {
     1215      ntrue++;
    10161216      if(Overflow_Error == FALSE)
    10171217      {
    10181218        Overflow_Error = TRUE;
    1019         PrintS("\n//** MPertvectors: OVERFLOW: ");
     1219        PrintS("\n// ** OVERFLOW in \"MPertvectors\": ");
    10201220        mpz_out_str( stdout, 10, pert_vector[i]);
    1021         PrintS(" is greater than 2147483647 (max. integer representation).\n");
    1022         Print("\n//** MPertvectors: So vector[%d] := %d may be wrong.\n", i+1, (*result)[i]);
    1023       }
    1024     }
     1221        PrintS(" is greater than 2147483647 (max. integer representation)");
     1222        Print("\n//  So vector[%d] := %d is wrong!!", i+1, (*result)[i]);
     1223      }
     1224    }
     1225  }
     1226
     1227  if(Overflow_Error == TRUE)
     1228  {
     1229    ivString(result, "pert_vector");
     1230    Print("\n// %d element(s) of it is overflow!!", ntrue);
    10251231  }
    10261232
     
    10281234  mpz_clear(sing_int);
    10291235  omFree(pert_vector);
    1030 
     1236  //omFree(pert_vector1);
    10311237  mpz_clear(tot_deg);
    10321238  mpz_clear(maxdeg);
     
    10631269  if(pdeg > nV || pdeg <= 0)
    10641270  {
    1065     WerrorS("\n//** MPertVectorslp: The perturbed degree is wrong.\n");
     1271    WerrorS("//** The perturbed degree is wrong!!");
    10661272    return pert_vector;
    10671273  }
     
    11161322  // Print(" %d", inveps);
    11171323#else
    1118   PrintS("\n//** MPertVectorslp: the \"big\" inverse epsilon %d.\n", inveps);
     1324  PrintS("\n// the \"big\" inverse epsilon %d", inveps);
    11191325#endif
    11201326
     
    12281434  return(ivM);
    12291435}
     1436
     1437//unused
     1438#if 0
     1439static intvec* MatrixOrderdp(int nV)
     1440{
     1441  int i;
     1442  intvec* ivM = new intvec(nV*nV);
     1443
     1444  for(i=0; i<nV; i++)
     1445  {
     1446    (*ivM)[i] = 1;
     1447  }
     1448  for(i=1; i<nV; i++)
     1449  {
     1450    (*ivM)[(i+1)*nV - i] = -1;
     1451  }
     1452  return(ivM);
     1453}
     1454#endif
    12301455
    12311456intvec* MivUnit(int nV)
     
    13511576  BOOLEAN nflow = FALSE;
    13521577
    1353   // compute gcd
     1578  // computes gcd
    13541579  mpz_set(ztmp, pert_vector[0]);
    13551580  for(i=0; i<niv; i++)
     
    13811606  if(j > nV - 1)
    13821607  {
     1608    // Print("\n//  MfPertwalk: geaenderter vector gleich Null! \n");
    13831609    delete result1;
    13841610    goto CHECK_OVERFLOW;
    13851611  }
    13861612
    1387   // check whether or not the perturbed weight vector lies in the Groebner cone
     1613// check that the perturbed weight vector lies in the Groebner cone
    13881614  if(test_w_in_ConeCC(G,result1) != 0)
    13891615  {
     1616    // Print("\n//  MfPertwalk: geaenderter vector liegt in Groebnerkegel! \n");
    13901617    delete result;
    13911618    result = result1;
     
    13981625  {
    13991626    delete result1;
     1627    // Print("\n// Mfpertwalk: geaenderter vector liegt nicht in Groebnerkegel! \n");
    14001628  }
    14011629
     
    14121640        Overflow_Error = TRUE;
    14131641        Print("\n// Xlev = %d and the %d-th element is", Xnlev,  i+1);
    1414         PrintS("\n//** Mfpertvector: OVERFLOW: ");
     1642        PrintS("\n// ** OVERFLOW in \"Mfpertvector\": ");
    14151643        mpz_out_str( stdout, 10, pert_vector[i]);
    1416         PrintS(" is greater than 2147483647 (max. integer representation).\n");
    1417         Print("\n//** Mfpertvector: So vector[%d] := %d may be wrong.\n", i+1, (*result)[i]);
     1644        PrintS(" is greater than 2147483647 (max. integer representation)");
     1645        Print("\n//  So vector[%d] := %d is wrong!!", i+1, (*result)[i]);
    14181646      }
    14191647    }
     
    14781706}
    14791707
    1480 
    1481 /*********************************************************
    1482 * Let Gomega be the inital form ideal of G.              *
    1483 * Let be M = {m_1,...,m_s} be the reduced GB of Gomega.  *
    1484 * Hence, mi=sum(h_ij.in_g_j) for suitable hij, i=1,...,s *
    1485 * Compute F = {f_1,...,f_s} with f_i=sum h_ij.g_j        *
    1486 **********************************************************/
     1708/*********************************************************************
     1709 * G is a red. Groebner basis w.r.t. <_1                             *
     1710 * Gomega is an initial form ideal of <G> w.r.t. a weight vector w   *
     1711 * M is a subideal of <Gomega> and M selft is a red. Groebner basis  *
     1712 *    of the ideal <Gomega> w.r.t. <_w                               *
     1713 * Let m_i = h1.gw1 + ... + hs.gws for each m_i in M; gwi in Gomega  *
     1714 * return F with n(F) = n(M) and f_i = h1.g1 + ... + hs.gs for each i*
     1715 ********************************************************************/
    14871716static ideal MLifttwoIdeal(ideal Gw, ideal M, ideal G)
    14881717{
    1489   ideal M_temp = idLift(Gw, M, NULL, FALSE, TRUE, TRUE, NULL);
    1490   idSkipZeroes(M_temp);
    1491   int i, j, nM = IDELEMS(M_temp);
     1718  ideal Mtmp = idLift(Gw, M, NULL, FALSE, TRUE, TRUE, NULL);
     1719
     1720  // If Gw is a GB, then isSB = TRUE, otherwise FALSE
     1721  // So, it is better, if one tests whether Gw is a GB
     1722  // in ideals.cc:
     1723  // idLift (ideal mod, ideal submod,ideal * rest, BOOLEAN goodShape,
     1724  //           BOOLEAN isSB,BOOLEAN divide,matrix * unit)
     1725
     1726  // Let be Mtmp = {m1,...,ms}, where mi=sum hij.in_gj, for all i=1,...,s
     1727  // We compute F = {f1,...,fs}, where fi=sum hij.gj
     1728  int i, j, nM = IDELEMS(Mtmp);
    14921729  ideal idpol, idLG;
    14931730  ideal F = idInit(nM, 1);
     1731
    14941732  for(i=0; i<nM; i++)
    14951733  {
    1496      //idpol = idVec2Ideal(M_temp->m[i]);
    1497      //idLG = MidMult(idpol, G);
    1498      //idpol = NULL;
    1499      idLG = MidMult(idVec2Ideal(M_temp->m[i]), G);
    1500     // PrintS("uups");
     1734     idpol = idVec2Ideal(Mtmp->m[i]);
     1735     idLG = MidMult(idpol, G);
     1736     idpol = NULL;
    15011737     F->m[i] = NULL;
    15021738     for(j=IDELEMS(idLG)-1; j>=0; j--)
     
    15041740       F->m[i] = pAdd(F->m[i], idLG->m[j]);
    15051741       idLG->m[j]=NULL;
    1506        idSkipZeroes(idLG);
    15071742     }
    15081743     idDelete(&idLG);
    15091744  }
    1510   //PrintS("aaps");
    1511   //idDelete(&idpol);
    1512   idDelete(&M_temp);
    1513   //F=kStd(F,NULL,testHomog,NULL,NULL,0,NULL,NULL);
    1514   return(F);
     1745  idDelete(&Mtmp);
     1746  return F;
    15151747}
    15161748
     
    15731805}
    15741806
    1575 /*****************************************************************************
    1576  * compute the gcd of the entries of the vectors curr_weight and diff_weight *
    1577  *****************************************************************************/
    1578 static int simplify_gcd(intvec* curr_weight, intvec* diff_weight)
    1579 {
    1580   int j;
    1581   int nRing = currRing->N;
    1582   int gcd_tmp = (*curr_weight)[0];
    1583   for (j=1; j<nRing; j++)
    1584   {
    1585     gcd_tmp = gcd(gcd_tmp, (*curr_weight)[j]);
    1586     if(gcd_tmp == 1)
    1587     {
    1588       break;
    1589     }
    1590   }
    1591   if(gcd_tmp != 1)
    1592   {
    1593     for (j=0; j<nRing; j++)
    1594     {
    1595     gcd_tmp = gcd(gcd_tmp, (*diff_weight)[j]);
    1596     if(gcd_tmp == 1)
    1597       {
    1598         break;
    1599       }
    1600     }
    1601   }
    1602   return gcd_tmp;
    1603 }
    1604 
     1807/**********************************************
     1808 * Look for the smallest absolut value in vec *
     1809 **********************************************/
     1810static int MivAbsMax(intvec* vec)
     1811{
     1812  int i,k;
     1813  if((*vec)[0] < 0)
     1814  {
     1815    k = -(*vec)[0];
     1816  }
     1817  else
     1818  {
     1819    k = (*vec)[0];
     1820  }
     1821  for(i=1; i < (vec->length()); i++)
     1822  {
     1823    if((*vec)[i] < 0)
     1824    {
     1825      if(-(*vec)[i] > k)
     1826      {
     1827        k = -(*vec)[i];
     1828      }
     1829    }
     1830    else
     1831    {
     1832      if((*vec)[i] > k)
     1833      {
     1834        k = (*vec)[i];
     1835      }
     1836    }
     1837  }
     1838  return k;
     1839}
    16051840
    16061841/**********************************************************************
     
    16141849  Overflow_Error = FALSE;
    16151850
    1616   assume(currRing != NULL && curr_weight != NULL && target_weight != NULL && G != NULL);
     1851  assume(currRing != NULL && curr_weight != NULL &&
     1852         target_weight != NULL && G != NULL);
    16171853
    16181854  int nRing = currRing->N;
    1619   int nG = IDELEMS(G);
    1620   int checkRed, j;
     1855  int checkRed, j, kkk, nG = IDELEMS(G);
    16211856  intvec* ivtemp;
    16221857
     
    16351870  mpz_set_si(sing_int,  2147483647);
    16361871
     1872  mpz_t sing_int_half;
     1873  mpz_init(sing_int_half);
     1874  mpz_set_si(sing_int_half,  3*(1073741824/2));
     1875
    16371876  mpz_t deg_w0_p1, deg_d0_p1;
    16381877  mpz_init(deg_w0_p1);
     
    16451884  mpz_t t_null;
    16461885  mpz_init(t_null);
    1647   mpz_set_si(t_null,0);
    16481886
    16491887  mpz_t ggt;
     
    16531891  mpz_init(dcw);
    16541892
     1893  //int tn0, tn1, tz1, ncmp, gcd_tmp, ntmp;
    16551894  int gcd_tmp;
    16561895  intvec* diff_weight = MivSub(target_weight, curr_weight);
     1896
    16571897  intvec* diff_weight1 = MivSub(target_weight, curr_weight);
    16581898  poly g;
    1659 #ifdef  NEXT_VECTORS_CC
    1660   ivString(diff_weight, "//** MwalkNextWeightCC: diff_weight");
    1661 #endif
     1899  //poly g, gw;
    16621900  for (j=0; j<nG; j++)
    16631901  {
     
    16681906      mpz_set_si(deg_w0_p1, MivDotProduct(ivtemp, curr_weight));
    16691907      mpz_set_si(deg_d0_p1, MivDotProduct(ivtemp, diff_weight));
    1670 #ifdef  NEXT_VECTORS_CC
    1671       ivString(ivtemp, "//** MwalkNextWeightCC: ivtemp");
    1672       PrintS("\n//** MwalkNextWeightCC: weighted degree w. r. t. curr_weight: ");
    1673       mpz_out_str( stdout, 10, deg_w0_p1);
    1674       PrintS("\n//** MwalkNextWeightCC: weighted degree w. r. t. diff_weight: ");
    1675       mpz_out_str( stdout, 10, deg_d0_p1);
    1676 #endif
    16771908      delete ivtemp;
    16781909
     
    16831914        mpz_set_si(MwWd, MivDotProduct(ivtemp, curr_weight));
    16841915        mpz_sub(s_zaehler, deg_w0_p1, MwWd);
    1685 #ifdef NEXT_VECTORS_CC
    1686         PrintS("\n//** MwalkNextWeightCC: Grad bzgl. curr_weight: ");
    1687         mpz_out_str( stdout, 10, MwWd);
    1688         PrintS("\n//** MwalkNextWeightCC: s_zaehler: ");
    1689         mpz_out_str( stdout, 10, s_zaehler);
    1690 #endif
    1691         if(mpz_cmp(s_zaehler,t_null) != 0)
     1916
     1917        if(mpz_cmp(s_zaehler, t_null) != 0)
    16921918        {
    16931919          mpz_set_si(MwWd, MivDotProduct(ivtemp, diff_weight));
    16941920          mpz_sub(s_nenner, MwWd, deg_d0_p1);
    1695 #ifdef NEXT_VECTORS_CC
    1696           PrintS("\n//** MwalkNextWeightCC: Grad bzgl. diff_weight: ");
    1697           mpz_out_str( stdout, 10, MwWd);
    1698           PrintS("\n//** MwalkNextWeightCC: s_nenner: ");
    1699           mpz_out_str( stdout, 10, s_nenner);
    1700 #endif
     1921
    17011922          // check for 0 < s <= 1
    1702           if((mpz_cmp_si(s_zaehler,0) > 0 && mpz_cmp(s_nenner,s_zaehler)>=0) ||
    1703              (mpz_cmp_si(s_zaehler,0) < 0 && mpz_cmp(s_nenner,s_zaehler)<=0))
     1923          if( (mpz_cmp(s_zaehler,t_null) > 0 &&
     1924               mpz_cmp(s_nenner, s_zaehler)>=0) ||
     1925              (mpz_cmp(s_zaehler, t_null) < 0 &&
     1926               mpz_cmp(s_nenner, s_zaehler)<=0))
    17041927          {
    17051928            // make both positive
    1706             if(mpz_cmp_si(s_zaehler,0) < 0)
     1929            if (mpz_cmp(s_zaehler, t_null) < 0)
    17071930            {
    17081931              mpz_neg(s_zaehler, s_zaehler);
     
    17131936            cancel(s_zaehler, s_nenner);
    17141937
    1715             if(mpz_cmp_si(t_nenner,0) != 0)
     1938            if(mpz_cmp(t_nenner, t_null) != 0)
    17161939            {
    17171940              mpz_mul(sztn, s_zaehler, t_nenner);
    17181941              mpz_mul(sntz, s_nenner, t_zaehler);
    1719 #ifdef NEXT_VECTORS_CC
    1720               PrintS("\n//** MwalkNextWeightCC: sntz: ");
    1721               mpz_out_str( stdout, 10, sntz);
    1722               PrintS("\n//** MwalkNextWeightCC: sztn: ");
    1723               mpz_out_str( stdout, 10, sztn);
    1724 #endif
     1942
    17251943              if(mpz_cmp(sztn,sntz) < 0)
    17261944              {
    1727                 mpz_add(t_nenner,s_nenner,t_null);
    1728                 mpz_add(t_zaehler,s_zaehler,t_null);
    1729 #ifdef NEXT_VECTORS_CC
    1730                 PrintS("\n//** MwalkNextWeightCC: t_nenner: ");
    1731                 mpz_out_str( stdout, 10, t_nenner);
    1732                 PrintS("\n//** MwalkNextWeightCC: t_zaehler: ");
    1733                 mpz_out_str( stdout, 10, t_zaehler);
    1734 #endif
     1945                mpz_add(t_nenner, t_null, s_nenner);
     1946                mpz_add(t_zaehler,t_null, s_zaehler);
    17351947              }
    17361948            }
    17371949            else
    17381950            {
    1739               mpz_add(t_nenner,s_nenner,t_null);
    1740               mpz_add(t_zaehler,s_zaehler,t_null);
    1741 #ifdef NEXT_VECTORS_CC
    1742               PrintS("\n//** MwalkNextWeightCC: t_nenner: ");
    1743               mpz_out_str( stdout, 10, t_nenner);
    1744               PrintS("\n//** MwalkNextWeightCC: t_zaehler: ");
    1745               mpz_out_str( stdout, 10, t_zaehler);
    1746 #endif
     1951              mpz_add(t_nenner, t_null, s_nenner);
     1952              mpz_add(t_zaehler,t_null, s_zaehler);
    17471953            }
    17481954          }
     
    17531959    }
    17541960  }
     1961//Print("\n// Alloc Size = %d \n", nRing*sizeof(mpz_t));
    17551962  mpz_t *vec=(mpz_t*)omAlloc(nRing*sizeof(mpz_t));
    17561963
    17571964
    1758   // there is no 0<t<=1, so define the next weight vector to be equal to the current weight vector
     1965  // there is no 0<t<1 and define the next weight vector that is equal to the current weight vector
    17591966  if(mpz_cmp(t_nenner, t_null) == 0)
    17601967  {
    1761 #ifdef NEXT_VECTORS_CC
    1762     Print("\n//** MwalkNextWeightCC: t_nenner equals zero.\n");
    1763 #endif
     1968    #ifndef SING_NDEBUG
     1969    Print("\n//MwalkNextWeightCC: t_nenner ist Null!");
     1970    #endif
    17641971    delete diff_weight;
    1765     diff_weight = ivCopy(curr_weight);
     1972    diff_weight = ivCopy(curr_weight);//take memory
    17661973    goto FINISH;
    17671974  }
    17681975
    1769   // t = 1, so define the next weight vector to be equal to the target vector
     1976  // define the target vector as the next weight vector, if t = 1
    17701977  if(mpz_cmp_si(t_nenner, 1)==0 && mpz_cmp_si(t_zaehler,1)==0)
    17711978  {
    17721979    delete diff_weight;
    1773     diff_weight = ivCopy(target_weight);
     1980    diff_weight = ivCopy(target_weight); //this takes memory
    17741981    goto FINISH;
    17751982  }
    17761983
     1984   checkRed = 0;
     1985
     1986  SIMPLIFY_GCD:
     1987
    17771988  // simplify the vectors curr_weight and diff_weight (C-int)
    1778   gcd_tmp = simplify_gcd(curr_weight,diff_weight);
     1989  gcd_tmp = (*curr_weight)[0];
     1990
     1991  for (j=1; j<nRing; j++)
     1992  {
     1993    gcd_tmp = gcd(gcd_tmp, (*curr_weight)[j]);
     1994    if(gcd_tmp == 1)
     1995    {
     1996      break;
     1997    }
     1998  }
    17791999  if(gcd_tmp != 1)
    17802000  {
    17812001    for (j=0; j<nRing; j++)
    17822002    {
    1783     (*curr_weight)[j] =  (*curr_weight)[j]/gcd_tmp;
    1784     (*diff_weight)[j] =  (*diff_weight)[j]/gcd_tmp;
    1785     }
     2003      gcd_tmp = gcd(gcd_tmp, (*diff_weight)[j]);
     2004      if(gcd_tmp == 1)
     2005      {
     2006        break;
     2007      }
     2008    }
     2009  }
     2010  if(gcd_tmp != 1)
     2011  {
     2012    for (j=0; j<nRing; j++)
     2013    {
     2014      (*curr_weight)[j] =  (*curr_weight)[j]/gcd_tmp;
     2015      (*diff_weight)[j] =  (*diff_weight)[j]/gcd_tmp;
     2016    }
     2017  }
     2018  if(checkRed > 0)
     2019  {
     2020    for (j=0; j<nRing; j++)
     2021    {
     2022      mpz_set_si(vec[j], (*diff_weight)[j]);
     2023    }
     2024    goto TEST_OVERFLOW;
    17862025  }
    17872026
    17882027#ifdef  NEXT_VECTORS_CC
    1789   Print("\n//** MwalkNextWeightCC: gcd of the weight vectors (current and target) = %d", gcd_tmp);
    1790   ivString(curr_weight, "new curr_weight");
    1791   ivString(diff_weight, "new diff_weight");
    1792   PrintS("\n//** MwalkNextWeightCC: t_zaehler: ");
    1793   mpz_out_str( stdout, 10, t_zaehler);
    1794   PrintS(", t_nenner: ");
    1795   mpz_out_str( stdout, 10, t_nenner);
    1796 #endif
     2028  Print("\n// gcd of the weight vectors (current and target) = %d", gcd_tmp);
     2029  ivString(curr_weight, "new cw");
     2030  ivString(diff_weight, "new dw");
     2031
     2032  PrintS("\n// t_zaehler: ");  mpz_out_str( stdout, 10, t_zaehler);
     2033  PrintS(", t_nenner: ");  mpz_out_str( stdout, 10, t_nenner);
     2034#endif
     2035
     2036  //  BOOLEAN isdwpos;
    17972037
    17982038  // construct a new weight vector
     
    18012041    mpz_set_si(dcw, (*curr_weight)[j]);
    18022042    mpz_mul(s_nenner, t_nenner, dcw);
    1803     mpz_mul_si(s_zaehler, t_zaehler, (*diff_weight)[j]);
     2043
     2044    if( (*diff_weight)[j]>0)
     2045    {
     2046      mpz_mul_ui(s_zaehler, t_zaehler, (*diff_weight)[j]);
     2047    }
     2048    else
     2049    {
     2050      mpz_mul_ui(s_zaehler, t_zaehler, -(*diff_weight)[j]);
     2051      mpz_neg(s_zaehler, s_zaehler);
     2052    }
    18042053    mpz_add(sntz, s_nenner, s_zaehler);
    18052054    mpz_init_set(vec[j], sntz);
    18062055
    18072056#ifdef NEXT_VECTORS_CC
    1808     Print("\n//** MwalkNextWeightCC:  j = %d ==> ", j);
     2057    Print("\n//   j = %d ==> ", j);
    18092058    PrintS("(");
    18102059    mpz_out_str( stdout, 10, t_nenner);
     
    18332082
    18342083#ifdef  NEXT_VECTORS_CC
    1835   PrintS("\n//** MwalkNextWeightCC: gcd of elements of the vector: ");
     2084  PrintS("\n// gcd of elements of the vector: ");
    18362085  mpz_out_str( stdout, 10, ggt);
    18372086#endif
    18382087
    1839 
    1840   // check whether vec[j]>2147483647, reduce it
    1841 
     2088/**********************************************************************
     2089* construct a new weight vector and check whether vec[j] is overflow, *
     2090* i.e. vec[j] > 2^31.                                                 *
     2091* If vec[j] doesn't overflow, define a weight vector. Otherwise,      *
     2092* report that overflow appears. In the second case, test whether the  *
     2093* the correctness of the new vector plays an important role           *
     2094**********************************************************************/
     2095  kkk=0;
     2096  for(j=0; j<nRing; j++)
     2097    {
     2098    if(mpz_cmp(vec[j], sing_int_half) >= 0)
     2099      {
     2100      goto REDUCTION;
     2101      }
     2102    }
     2103  checkRed = 1;
    18422104  for (j=0; j<nRing; j++)
     2105    {
     2106      (*diff_weight)[j] = mpz_get_si(vec[j]);
     2107    }
     2108  goto SIMPLIFY_GCD;
     2109
     2110  REDUCTION:
     2111  for (j=0; j<nRing; j++)
    18432112  {
    18442113    (*diff_weight)[j] = mpz_get_si(vec[j]);
    18452114  }
    1846   while(MivAbsMax(diff_weight) > 100000000)
     2115  while(MivAbsMax(diff_weight) >= 5)
    18472116  {
    18482117    for (j=0; j<nRing; j++)
     
    18502119      if(mpz_cmp_si(ggt,1)==0)
    18512120      {
    1852         (*diff_weight1)[j] = floor(0.001*(*diff_weight)[j] + 0.5);
    1853 #ifdef  NEXT_VECTORS_CC
    1854         Print("\n//** MwalkNextWeightCC: vector[%d] = %d \n",j+1, (*diff_weight1)[j]);
    1855 #endif
     2121        (*diff_weight1)[j] = floor(0.1*(*diff_weight)[j] + 0.5);
     2122        // Print("\n//  vector[%d] = %d \n",j+1, (*diff_weight1)[j]);
    18562123      }
    18572124      else
    18582125      {
    18592126        mpz_divexact(vec[j], vec[j], ggt);
    1860         (*diff_weight1)[j] = floor(0.001*(*diff_weight)[j] + 0.5);
    1861 #ifdef  NEXT_VECTORS_CC
    1862         Print("\n//** MwalkNextWeightCC: vector[%d] = %d \n",j+1, (*diff_weight1)[j]);
    1863 #endif
    1864       }
    1865     }
     2127        (*diff_weight1)[j] = floor(0.1*(*diff_weight)[j] + 0.5);
     2128        // Print("\n//  vector[%d] = %d \n",j+1, (*diff_weight1)[j]);
     2129      }
     2130/*
     2131      if((*diff_weight1)[j] == 0)
     2132      {
     2133        kkk = kkk + 1;
     2134      }
     2135*/
     2136    }
     2137
     2138
     2139/*
     2140    if(kkk > nRing - 1)
     2141    {
     2142      // diff_weight was reduced to zero
     2143      // Print("\n // MwalkNextWeightCC: geaenderter Vector gleich Null! \n");
     2144      goto TEST_OVERFLOW;
     2145    }
     2146*/
    18662147
    18672148    if(test_w_in_ConeCC(G,diff_weight1) != 0)
    18682149    {
    1869 #ifdef  NEXT_VECTORS_CC
    1870       Print("\n//** MwalkNextWeightCC: changed weight vector in correct cone.\n");
    1871 #endif
     2150      Print("\n// MwalkNextWeightCC: geaenderter vector liegt in Groebnerkegel! \n");
    18722151      for (j=0; j<nRing; j++)
    18732152      {
    18742153        (*diff_weight)[j] = (*diff_weight1)[j];
    18752154      }
    1876 
    1877       // simplify the vectors curr_weight and diff_weight (C-int)
    1878       gcd_tmp = simplify_gcd(curr_weight,diff_weight);
    1879       if(gcd_tmp != 1)
    1880       {
    1881         for (j=0; j<nRing; j++)
    1882         {
    1883           (*curr_weight)[j] =  (*curr_weight)[j]/gcd_tmp;
    1884           (*diff_weight)[j] =  (*diff_weight)[j]/gcd_tmp;
    1885           mpz_set_si(vec[j],(*diff_weight)[j]);
    1886         }
     2155      if(MivAbsMax(diff_weight) < 5)
     2156      {
     2157        checkRed = 1;
     2158        goto SIMPLIFY_GCD;
    18872159      }
    18882160    }
    18892161    else
    18902162    {
    1891 #ifdef  NEXT_VECTORS_CC
    1892       Print("\n//** MwalkNextWeightCC: changed weight vector not in correct cone.\n");
    1893 #endif
    1894 
    1895       for (j=0; j<nRing; j++)
    1896       {
    1897         (*diff_weight)[j] = mpz_get_si(vec[j]);
    1898       }
     2163      // Print("\n// MwalkNextWeightCC: geaenderter vector liegt nicht in Groebnerkegel! \n");
    18992164      break;
    19002165    }
    19012166  }
    19022167
     2168 TEST_OVERFLOW:
     2169
    19032170  for (j=0; j<nRing; j++)
    19042171  {
    1905     if((mpz_cmp(vec[j], sing_int)>=0) && ((*diff_weight)[j] > 2147483647))
     2172    if(mpz_cmp(vec[j], sing_int)>=0)
    19062173    {
    19072174      if(Overflow_Error == FALSE)
    19082175      {
    19092176        Overflow_Error = TRUE;
    1910         PrintS("\n//** MwalkNextWeightCC: OVERFLOW: ");
     2177        PrintS("\n// ** OVERFLOW in \"MwalkNextWeightCC\": ");
    19112178        mpz_out_str( stdout, 10, vec[j]);
    1912         PrintS(" is greater than 2147483647 (max. integer representation)");
    1913         Print("\n//  So vector[%d] may be wrong.\n",j+1/*,gmp: vec[j]*/);
     2179        PrintS(" is greater than 2147483647 (max. integer representation)\n");
     2180        //Print("//  So vector[%d] := %d is wrong!!\n",j+1, vec[j]);// vec[j] is mpz_t
    19142181      }
    19152182    }
     
    19172184
    19182185 FINISH:
    1919   delete diff_weight1;
    1920   mpz_clear(t_zaehler);
    1921   mpz_clear(t_nenner);
    1922   mpz_clear(s_zaehler);
    1923   mpz_clear(s_nenner);
    1924   mpz_clear(sntz);
    1925   mpz_clear(sztn);
    1926   mpz_clear(temp);
    1927   mpz_clear(MwWd);
    1928   mpz_clear(deg_w0_p1);
    1929   mpz_clear(deg_d0_p1);
    1930   mpz_clear(ggt);
    1931   omFree(vec);
    1932   mpz_clear(sing_int);
    1933   mpz_clear(dcw);
    1934   mpz_clear(t_null);
     2186   delete diff_weight1;
     2187   mpz_clear(t_zaehler);
     2188   mpz_clear(t_nenner);
     2189   mpz_clear(s_zaehler);
     2190   mpz_clear(s_nenner);
     2191   mpz_clear(sntz);
     2192   mpz_clear(sztn);
     2193   mpz_clear(temp);
     2194   mpz_clear(MwWd);
     2195   mpz_clear(deg_w0_p1);
     2196   mpz_clear(deg_d0_p1);
     2197   mpz_clear(ggt);
     2198   omFree(vec);
     2199   mpz_clear(sing_int_half);
     2200   mpz_clear(sing_int);
     2201   mpz_clear(dcw);
     2202   mpz_clear(t_null);
     2203
     2204
    19352205
    19362206  if(Overflow_Error == FALSE)
     
    19382208    Overflow_Error = nError;
    19392209  }
    1940   rComplete(currRing);
    1941   for(j=0; j<IDELEMS(G); j++)
    1942   {
    1943     poly p=G->m[j];
     2210 rComplete(currRing);
     2211  for(kkk=0; kkk<IDELEMS(G);kkk++)
     2212  {
     2213    poly p=G->m[kkk];
    19442214    while(p!=NULL)
    19452215    {
     
    19812251}
    19822252
    1983 /*****************************************************
    1984  * define and execute a new ring with ordering (M,C) *
    1985  *****************************************************/
    1986 static ring VMatrDefault(intvec* va)
     2253/**************************************************************
     2254 * define and execute a new ring which order is (a(vb),a(va),lp,C) *
     2255 * ************************************************************/
     2256static void VMrHomogeneous(intvec* va, intvec* vb)
    19872257{
    19882258
     
    20022272  r->cf  = currRing->cf;
    20032273  r->N   = currRing->N;
    2004 
    20052274  int nb = 4;
     2275
    20062276
    20072277  //names
     
    20142284  }
    20152285
    2016   /*weights: entries for 3 blocks: NULL Made:???*/
    2017   r->wvhdl = (int **)omAlloc0(nb * sizeof(int_ptr));
    2018   r->wvhdl[0] = (int*) omAlloc(nv*nv*sizeof(int));
    2019   r->wvhdl[1] =NULL; // (int*) omAlloc(nv*sizeof(int));
    2020   r->wvhdl[2]=NULL;
    2021   r->wvhdl[3]=NULL;
    2022   for(i=0; i<nv*nv; i++)
    2023     r->wvhdl[0][i] = (*va)[i];
    2024 
    2025   /* order: a,lp,C,0 */
    2026   r->order = (int *) omAlloc(nb * sizeof(int *));
    2027   r->block0 = (int *)omAlloc0(nb * sizeof(int *));
    2028   r->block1 = (int *)omAlloc0(nb * sizeof(int *));
    2029 
    2030   // ringorder a for the first block: var 1..nv
    2031   r->order[0]  = ringorder_M;
    2032   r->block0[0] = 1;
    2033   r->block1[0] = nv;
    2034 
    2035   // ringorder C for the second block
    2036   r->order[1]  = ringorder_C;
    2037   r->block0[1] = 1;
    2038   r->block1[1] = nv;
    2039 
    2040 
    2041 // ringorder C for the third block: var 1..nv
    2042   r->order[2]  = ringorder_C;
    2043   r->block0[2] = 1;
    2044   r->block1[2] = nv;
    2045 
    2046 
    2047   // the last block: everything is 0
    2048   r->order[3]  = 0;
    2049 
    2050   // polynomial ring
    2051   r->OrdSgn    = 1;
    2052 
    2053   // complete ring intializations
    2054 
    2055   rComplete(r);
    2056 
    2057   //rChangeCurrRing(r);
    2058   return r;
    2059 }
    2060 
    2061 /***********************************************************
    2062  * define and execute a new ring with ordering (a(vb),M,C) *
    2063  ***********************************************************/
    2064 static ring VMatrRefine(intvec* va, intvec* vb)
    2065 {
    2066 
    2067   if ((currRing->ppNoether)!=NULL)
    2068   {
    2069     pDelete(&(currRing->ppNoether));
    2070   }
    2071   if (((sLastPrinted.rtyp>BEGIN_RING) && (sLastPrinted.rtyp<END_RING)) ||
    2072       ((sLastPrinted.rtyp==LIST_CMD)&&(lRingDependend((lists)sLastPrinted.data))))
    2073   {
    2074     sLastPrinted.CleanUp();
    2075   }
    2076 
    2077   ring r = (ring) omAlloc0Bin(sip_sring_bin);
    2078   int i, nv = currRing->N;
    2079   int nvs = nv*nv;
    2080   r->cf  = currRing->cf;
    2081   r->N   = currRing->N;
    2082 
    2083   int nb = 4;
    2084 
    2085   //names
    2086   char* Q; // In order to avoid the corrupted memory, do not change.
    2087   r->names = (char **) omAlloc0(nv * sizeof(char_ptr));
    2088   for(i=0; i<nv; i++)
    2089   {
    2090     Q = currRing->names[i];
    2091     r->names[i]  = omStrDup(Q);
    2092   }
    2093 
    2094   /*weights: entries for 3 blocks: NULL Made:???*/
     2286  //weights: entries for 3 blocks: NULL Made:???
    20952287  r->wvhdl = (int **)omAlloc0(nb * sizeof(int_ptr));
    20962288  r->wvhdl[0] = (int*) omAlloc(nv*sizeof(int));
    2097   r->wvhdl[1] = (int*) omAlloc(nvs*sizeof(int));
    2098   r->wvhdl[2]=NULL;
    2099   r->wvhdl[3]=NULL;
    2100   for(i=0; i<nvs; i++)
    2101   {
    2102     r->wvhdl[1][i] = (*va)[i];
    2103   }
    2104   for(i=0; i<nv; i++)
    2105   {
    2106     r->wvhdl[0][i] = (*vb)[i];
    2107   }
    2108   /* order: a,lp,C,0 */
     2289  r->wvhdl[1] = (int*) omAlloc((nv-1)*sizeof(int));
     2290
     2291  for(i=0; i<nv-1; i++)
     2292  {
     2293    r->wvhdl[1][i] = (*vb)[i];
     2294    r->wvhdl[0][i] = (*va)[i];
     2295  }
     2296  r->wvhdl[0][nv] = (*va)[nv];
     2297
     2298  // order: (1..1),a,lp,C
    21092299  r->order = (int *) omAlloc(nb * sizeof(int *));
    21102300  r->block0 = (int *)omAlloc0(nb * sizeof(int *));
     
    21162306  r->block1[0] = nv;
    21172307
    2118   // ringorder M for the second block: var 1..nv
    2119   r->order[1]  = ringorder_M;
    2120   r->block0[1] = 1;
     2308 // ringorder a for the second block: var 2..nv
     2309  r->order[1]  = ringorder_a;
     2310  r->block0[1] = 2;
    21212311  r->block1[1] = nv;
    21222312
    2123   // ringorder C for the third block: var 1..nv
    2124   r->order[2]  = ringorder_C;
    2125   r->block0[2] = 1;
    2126   r->block1[2] = nv;
    2127 
    2128 
    2129   // the last block: everything is 0
    2130   r->order[3]  = 0;
    2131 
    2132   // polynomial ring
    2133   r->OrdSgn    = 1;
    2134 
    2135   // complete ring intializations
    2136 
    2137   rComplete(r);
    2138 
    2139   //rChangeCurrRing(r);
    2140   return r;
    2141 }
    2142 
    2143 /****************************************************************
    2144  * define and execute a new ring with ordering (a(va),Wp(vb),C) *
    2145  * **************************************************************/
    2146 
    2147 static ring VMrRefine(intvec* va, intvec* vb)
    2148 {
    2149 
    2150   if ((currRing->ppNoether)!=NULL)
    2151   {
    2152     pDelete(&(currRing->ppNoether));
    2153   }
    2154   if (((sLastPrinted.rtyp>BEGIN_RING) && (sLastPrinted.rtyp<END_RING)) ||
    2155       ((sLastPrinted.rtyp==LIST_CMD)&&(lRingDependend((lists)sLastPrinted.data))))
    2156   {
    2157     sLastPrinted.CleanUp();
    2158   }
    2159 
    2160   ring r = (ring) omAlloc0Bin(sip_sring_bin);
    2161   int i, nv = currRing->N;
    2162 
    2163   r->cf  = currRing->cf;
    2164   r->N   = currRing->N;
    2165   //int nb = nBlocks(currRing)  + 1;
    2166   int nb = 4;
    2167 
    2168   //names
    2169   char* Q; // In order to avoid the corrupted memory, do not change.
    2170   r->names = (char **) omAlloc0(nv * sizeof(char_ptr));
    2171   for(i=0; i<nv; i++)
    2172   {
    2173     Q = currRing->names[i];
    2174     r->names[i]  = omStrDup(Q);
    2175   }
    2176 
    2177   //weights: entries for 3 blocks: NULL Made:???
    2178   r->wvhdl = (int **)omAlloc0(nb * sizeof(int_ptr));
    2179   r->wvhdl[0] = (int*) omAlloc(nv*sizeof(int));
    2180   r->wvhdl[1] = (int*) omAlloc(nv*sizeof(int));
    2181 
    2182   for(i=0; i<nv; i++)
    2183   {
    2184     r->wvhdl[0][i] = (*va)[i];
    2185     r->wvhdl[1][i] = (*vb)[i];
    2186   }
    2187   r->wvhdl[2]=NULL;
    2188   r->wvhdl[3]=NULL;
    2189 
    2190   // order: (1..1),a,lp,C
    2191   r->order = (int *) omAlloc(nb * sizeof(int *));
    2192   r->block0 = (int *)omAlloc0(nb * sizeof(int *));
    2193   r->block1 = (int *)omAlloc0(nb * sizeof(int *));
    2194 
    2195   // ringorder a for the first block: var 1..nv
    2196   r->order[0]  = ringorder_a;
    2197   r->block0[0] = 1;
    2198   r->block1[0] = nv;
    2199 
    2200  // ringorder Wp for the second block: var 1..nv
    2201   r->order[1]  = ringorder_Wp;
    2202   r->block0[1] = 1;
    2203   r->block1[1] = nv;
    2204 
    2205   // ringorder lp for the third block: var 1..nv
    2206   r->order[2]  = ringorder_C;
    2207   r->block0[2] = 1;
     2313  // ringorder lp for the third block: var 2..nv
     2314  r->order[2]  = ringorder_lp;
     2315  r->block0[2] = 2;
    22082316  r->block1[2] = nv;
    22092317
     
    22122320  // especially, by ring syz_ring=rCurrRingAssure_SyzComp();
    22132321  // therefore, nb  must be (nBlocks(currRing)  + 1)
    2214   r->order[3]  = 0;
     2322  r->order[3]  = ringorder_C;
    22152323
    22162324  // polynomial ring
     
    22212329  rComplete(r);
    22222330
    2223   //rChangeCurrRing(r);
    2224   return r;
    2225 }
     2331  rChangeCurrRing(r);
     2332}
     2333
    22262334
    22272335/**************************************************************
    22282336 * define and execute a new ring which order is (a(va),lp,C)  *
    22292337 * ************************************************************/
    2230 //static void VMrDefault(intvec* va)
    2231 static ring VMrDefault(intvec* va)
     2338static void VMrDefault(intvec* va)
    22322339{
    22332340
     
    22962403  rComplete(r);
    22972404
    2298   //rChangeCurrRing(r);
    2299   return r;
     2405  rChangeCurrRing(r);
    23002406}
    23012407
     
    24862592  r->OrdSgn    = 1;
    24872593
     2594
     2595//   if (rParameter(currRing)!=NULL)
     2596//   {
     2597//     r->cf->extRing->qideal->m[0]=p_Copy(currRing->cf->extRing->qideal->m[0], currRing->cf->extRing);
     2598//     int l=rPar(currRing);
     2599//     r->cf->extRing->names=(char **)omAlloc(l*sizeof(char_ptr));
     2600//
     2601//     for(i=l-1;i>=0;i--)
     2602//     {
     2603//       rParameter(r)[i]=omStrDup(rParameter(currRing)[i]);
     2604//     }
     2605//   }
     2606
    24882607  // complete ring intializations
    24892608
     
    25002619}
    25012620
    2502 
     2621//unused
     2622/**************************************************************
     2623 * check wheather one or more components of a vector are zero *
     2624 **************************************************************/
     2625//#if 0
     2626static int isNolVector(intvec* hilb)
     2627{
     2628  int i;
     2629  for(i=hilb->length()-1; i>=0; i--)
     2630  {
     2631    if((* hilb)[i]==0)
     2632    {
     2633      return 1;
     2634    }
     2635  }
     2636  return 0;
     2637}
     2638//#endif
    25032639
    25042640/******************************  Februar 2002  ****************************
     
    27892925  for(i=IDELEMS(G)-1; i>=0; i--)
    27902926  {
     2927#if 0
     2928    if(pLength(G->m[i])>2)
     2929    {
     2930      return 1;
     2931    }
     2932#else
    27912933    if((G->m[i]!=NULL) /* len >=0 */
    27922934       && (G->m[i]->next!=NULL) /* len >=1 */
     
    27962938       )
    27972939    {
    2798       return 1;
    2799     }
     2940    return 1;
     2941    }
     2942#endif
    28002943  }
    28012944  return 0;
     
    28703013  for(i=nG-1; i>=0; i--)
    28713014  {
     3015#if 0
     3016    poly t;
     3017    if((t=pSub(pCopy(H0->m[i]), pCopy(H1->m[i]))) != NULL)
     3018    {
     3019      pDelete(&t);
     3020      return 0;
     3021    }
     3022    pDelete(&t);
     3023#else
    28723024    if(!pEqualPolys(H0->m[i],H1->m[i]))
    28733025    {
    28743026      return 0;
    28753027    }
     3028#endif
    28763029  }
    28773030  return 1;
     
    29033056  return result;
    29043057}
    2905 
    2906 
    2907 
     3058#endif
     3059
     3060//unused
    29083061/************************************************************************
    29093062 * perturb the weight vector iva w.r.t. the ideal G.                    *
     
    29133066 *  d := (2*maxdeg*maxdeg + (nV+1)*maxdeg)*m;                           *
    29143067 ************************************************************************/
     3068#if 0
    29153069static intvec* TranPertVector(ideal G, intvec* iva)
    29163070{
     
    30533207  return repr_vector;
    30543208}
    3055 
    3056 
    3057 
    3058 
     3209#endif
     3210
     3211//unused
     3212#if 0
    30593213static intvec* TranPertVector_lp(ideal G)
    30603214{
     
    31493303  return repr_vector;
    31503304}
    3151 
    3152 
     3305#endif
     3306
     3307//unused
     3308#if 0
    31533309static intvec* RepresentationMatrix_Dp(ideal G, intvec* M)
    31543310{
     
    32483404 *****************************************************************************/
    32493405
     3406//unused
     3407#if 0
     3408static int testnegintvec(intvec* v)
     3409{
     3410  int n = v->length();
     3411  int i;
     3412  for(i=0; i<n; i++)
     3413  {
     3414    if((*v)[i]<0)
     3415    {
     3416      return(1);
     3417    }
     3418  }
     3419  return(0);
     3420}
     3421#endif
     3422
     3423// npwinc = 0, if curr_weight doesn't stay in the correct Groebner cone
    32503424static ideal Rec_LastGB(ideal G, intvec* curr_weight,
    32513425                        intvec* orig_target_weight, int tp_deg, int npwinc)
     
    36133787    nstep ++;
    36143788    to = clock();
    3615     // compute an initial form ideal of <G> w.r.t. "curr_vector"
     3789    /* compute an initial form ideal of <G> w.r.t. "curr_vector" */
    36163790    Gomega = MwalkInitialForm(G, curr_weight);
    36173791    xtif=xtif+clock()-to;
    3618 
     3792#if 0
     3793    if(Overflow_Error == TRUE)
     3794    {
     3795      for(i=nV-1; i>=0; i--)
     3796        (*curr_weight)[i] = (*extra_curr_weight)[i];
     3797      delete extra_curr_weight;
     3798      goto LAST_GB_ALT2;
     3799    }
     3800#endif
    36193801    oldRing = currRing;
    36203802
    3621     // define a new ring that its ordering is "(a(curr_weight),lp)
     3803    /* define a new ring that its ordering is "(a(curr_weight),lp) */
    36223804    if (rParameter(currRing) != NULL)
    36233805    {
     
    36723854    if(Overflow_Error == TRUE)
    36733855    {
     3856      /*
     3857        ivString(next_weight, "omega");
     3858        PrintS("\n// ** The weight vector does NOT stay in Cone!!\n");
     3859      */
    36743860#ifdef TEST_OVERFLOW
    36753861      goto  TEST_OVERFLOW_OI;
    36763862#endif
     3863
    36773864      newRing = currRing;
    36783865      if (rParameter(currRing) != NULL)
     
    37353922  TimeStringFractal(xftinput, tostd, xtif, xtstd, xtextra,xtlift, xtred,xtnw);
    37363923
    3737   //Print("\n// pSetm_Error = (%d)", ErrorCheck());
     3924  Print("\n// pSetm_Error = (%d)", ErrorCheck());
    37383925  //Print("\n// Overflow_Error? (%d)", nOverflow_Error);
    37393926  Print("\n// Awalk2 took %d steps!!", nstep);
     
    37693956{
    37703957  int i, weight_norm;
    3771   //int randCount=0;
    37723958  int nV = currRing->N;
    37733959  intvec* next_weight2;
    37743960  intvec* next_weight22 = new intvec(nV);
    3775   intvec* result = new intvec(nV);
    3776 
    3777   //compute a perturbed next weight vector "next_weight1"
    3778   intvec* next_weight1 = MkInterRedNextWeight(MPertVectors(G,MivMatrixOrderRefine(curr_weight,target_weight),pert_deg),target_weight,G);
    3779   //compute a random next weight vector "next_weight2"
    3780   while(1)
    3781   {
    3782     weight_norm = 0;
    3783     while(weight_norm == 0)
    3784     {
    3785       for(i=0; i<nV; i++)
    3786       {
    3787         (*next_weight22)[i] = rand() % 60000 - 30000;
    3788         weight_norm = weight_norm + (*next_weight22)[i]*(*next_weight22)[i];
    3789       }
    3790       weight_norm = 1 + floor(sqrt(weight_norm));
    3791     }
    3792     for(i=0; i<nV; i++)
    3793     {
    3794       if((*next_weight22)[i] < 0)
    3795       {
    3796         (*next_weight22)[i] = 1 + (*curr_weight)[i] + floor(weight_rad*(*next_weight22)[i]/weight_norm);
    3797       }
    3798       else
    3799       {
    3800         (*next_weight22)[i] = (*curr_weight)[i] + floor(weight_rad*(*next_weight22)[i]/weight_norm);
    3801       }
    3802     }
    3803     if(test_w_in_ConeCC(G, next_weight22) == 1)
    3804     {
    3805       next_weight2 = MkInterRedNextWeight(next_weight22,target_weight,G);
    3806       delete next_weight22;
    3807       break;
    3808     }
    3809   }
    3810   // compute "usual" next weight vector
    38113961  intvec* next_weight = MwalkNextWeightCC(curr_weight,target_weight, G);
    3812   ideal G_test = MwalkInitialForm(G, next_weight);
    3813   ideal G_test2 = MwalkInitialForm(G, next_weight2);
    3814 
    3815   // compare next weights
    3816   if(Overflow_Error == FALSE)
    3817   {
     3962  if(MivComp(next_weight, target_weight) == 1)
     3963  {
     3964    return(next_weight);
     3965  }
     3966  else
     3967  {
     3968    //compute a perturbed next weight vector "next_weight1"
     3969    intvec* next_weight1 = MkInterRedNextWeight(MPertVectors(G, MivMatrixOrder(curr_weight), pert_deg), target_weight, G);
     3970    //Print("\n // size of next_weight1 = %d", sizeof((*next_weight1)));
     3971
     3972    //compute a random next weight vector "next_weight2"
     3973    while(1)
     3974    {
     3975      weight_norm = 0;
     3976      while(weight_norm == 0)
     3977      {
     3978        for(i=0; i<nV; i++)
     3979        {
     3980          //Print("\n//  next_weight[%d]  = %d", i, (*next_weight)[i]);
     3981          (*next_weight22)[i] = rand() % 60000 - 30000;
     3982          weight_norm = weight_norm + (*next_weight22)[i]*(*next_weight22)[i];
     3983        }
     3984        weight_norm = 1 + floor(sqrt(weight_norm));
     3985      }
     3986
     3987      for(i=nV-1; i>=0; i--)
     3988      {
     3989        if((*next_weight22)[i] < 0)
     3990        {
     3991          (*next_weight22)[i] = 1 + (*curr_weight)[i] + floor(weight_rad*(*next_weight22)[i]/weight_norm);
     3992        }
     3993        else
     3994        {
     3995          (*next_weight22)[i] = (*curr_weight)[i] + floor(weight_rad*(*next_weight22)[i]/weight_norm);
     3996        }
     3997      //Print("\n//  next_weight22[%d]  = %d", i, (*next_weight22)[i]);
     3998      }
     3999
     4000      if(test_w_in_ConeCC(G, next_weight22) == 1)
     4001      {
     4002        //Print("\n//MWalkRandomNextWeight: next_weight2 im Kegel\n");
     4003        next_weight2 = MkInterRedNextWeight(next_weight22, target_weight, G);
     4004        delete next_weight22;
     4005        break;
     4006      }
     4007    }
     4008    intvec* result = new intvec(nV);
     4009    ideal G_test = MwalkInitialForm(G, next_weight);
    38184010    ideal G_test1 = MwalkInitialForm(G, next_weight1);
     4011    ideal G_test2 = MwalkInitialForm(G, next_weight2);
     4012
     4013    // compare next_weights
    38194014    if(IDELEMS(G_test1) < IDELEMS(G_test))
    38204015    {
    3821       if(IDELEMS(G_test2) < IDELEMS(G_test1))
    3822       {
    3823         // |G_test2| < |G_test1| < |G_test|
     4016      if(IDELEMS(G_test2) <= IDELEMS(G_test1)) // |G_test2| <= |G_test1| < |G_test|
     4017      {
    38244018        for(i=0; i<nV; i++)
    38254019        {
     
    38274021        }
    38284022      }
    3829       else
    3830       {
    3831         // |G_test1| < |G_test|, |G_test1| <= |G_test2|
     4023      else // |G_test1| < |G_test|, |G_test1| < |G_test2|
     4024      {
    38324025        for(i=0; i<nV; i++)
    38334026        {
     
    38384031    else
    38394032    {
    3840       if(IDELEMS(G_test2) < IDELEMS(G_test)) // |G_test2| < |G_test| <= |G_test1|
     4033      if(IDELEMS(G_test2) <= IDELEMS(G_test)) // |G_test2| <= |G_test| <= |G_test1|
    38414034      {
    38424035        for(i=0; i<nV; i++)
     
    38454038        }
    38464039      }
    3847       else
    3848       {
    3849         // |G_test| < |G_test1|, |G_test| <= |G_test2|
     4040      else // |G_test| <= |G_test1|, |G_test| < |G_test2|
     4041      {
    38504042        for(i=0; i<nV; i++)
    38514043        {
     
    38544046      }
    38554047    }
    3856     idDelete(&G_test1);
    3857   }
    3858   else
    3859   {
    3860     Overflow_Error = FALSE;
    3861     if(IDELEMS(G_test2) < IDELEMS(G_test))
    3862     {
    3863       for(i=1; i<nV; i++)
    3864       {
    3865         (*result)[i] = (*next_weight2)[i];
    3866       }
    3867     }
    3868     else
    3869     {
    3870       for(i=0; i<nV; i++)
    3871       {
    3872         (*result)[i] = (*next_weight)[i];
    3873       }
    3874     }
    3875   }
    3876   idDelete(&G_test);
    3877   idDelete(&G_test2);
    3878   if(test_w_in_ConeCC(G, result) == 1)
    3879   {
    3880     delete next_weight2;
    38814048    delete next_weight;
    38824049    delete next_weight1;
    3883     return result;
    3884   }
    3885   else
    3886   {
    3887     delete result;
    3888     delete next_weight2;
    3889     delete next_weight1;
    3890     return next_weight;
     4050    idDelete(&G_test);
     4051    idDelete(&G_test1);
     4052    idDelete(&G_test2);
     4053    if(test_w_in_ConeCC(G, result) == 1)
     4054    {
     4055      delete next_weight2;
     4056      return result;
     4057    }
     4058    else
     4059    {
     4060      delete result;
     4061      return next_weight2;
     4062    }
    38914063  }
    38924064}
     
    38944066
    38954067/***************************************************************************
    3896  * The procedure REC_GB_Mwalk computes a GB for <G> w.r.t. the weight order*
     4068 * The procedur REC_GB_Mwalk computes a GB for <G> w.r.t. the weight order *
    38974069 * otw, where G is a reduced GB w.r.t. the weight order cw.                *
    38984070 * The new procedur Mwalk calls REC_GB.                                    *
     
    39454117    if(test_G_GB_walk(H0_tmp,H1)==1)
    39464118    {
    3947       Print("\n//REC_GB_Mwalk: input in %d-th recursive is a GB!\n",tp_deg);
     4119      //Print("\n//REC_GB_Mwalk: input in %d-th recursive is a GB!\n",tp_deg);
    39484120      idDelete(&H0_tmp);
    39494121      idDelete(&H1);
     
    39794151      goto NEXT_STEP;
    39804152    }
    3981     Print("\n//REC_GB_Mwalk: Entering the %d-th step in the %d-th recursive:\n",nwalk,tp_deg);
     4153    //Print("\n//REC_GB_Mwalk: Entering the %d-th step in the %d-th recursive:\n",nwalk,tp_deg);
    39824154    to = clock();
    39834155    // compute an initial form ideal of <G> w.r.t. "curr_vector"
     
    40644236    if(Overflow_Error == TRUE)
    40654237    {
    4066       PrintS("\n//REC_GB_Mwalk: The computed vector does NOT stay in the correct cone!!\n");
     4238      //PrintS("\n//REC_GB_Mwalk: The computed vector does NOT stay in the correct cone!!\n");
    40674239      nnwinC = 0;
    40684240      if(tp_deg == nV)
     
    41794351}
    41804352
    4181 /*******************************
    4182  * THE GROEBNER WALK ALGORITHM *
    4183  *******************************/
    4184 ideal Mwalk(ideal Go, intvec* orig_M, intvec* target_M, ring baseRing)
    4185 {
    4186   BITSET save1 = si_opt_1; // save current options
    4187   si_opt_1 &= (~Sy_bit(OPT_REDSB)); // no reduced Groebner basis
    4188   si_opt_1 &= (~Sy_bit(OPT_REDTAIL)); // not tail reductions
    4189   //si_opt_1|=(Sy_bit(OPT_REDTAIL)|Sy_bit(OPT_REDSB));
     4353
     4354// THE NEW GROEBNER WALK ALGORITHM
     4355// Groebnerwalk with a recursive "second" alternative GW, called REC_GB_Mwalk that only computes the last reduced GB
     4356ideal Mwalk(ideal Go, intvec* curr_weight, intvec* target_weight)
     4357{
    41904358  Set_Error(FALSE);
    41914359  Overflow_Error = FALSE;
    4192 #ifdef TIME_TEST
     4360  //Print("// pSetm_Error = (%d)", ErrorCheck());
     4361
    41934362  clock_t tinput, tostd, tif=0, tstd=0, tlift=0, tred=0, tnw=0;
    41944363  xtif=0; xtstd=0; xtlift=0; xtred=0; xtnw=0;
    41954364  tinput = clock();
    41964365  clock_t tim;
    4197 #endif
     4366  nstep=0;
     4367  int i;
     4368  int nV = currRing->N;
     4369  int nwalk=0;
     4370  int endwalks=0;
     4371
     4372  ideal Gomega, M, F, Gomega1, Gomega2, M1, F1, G;
     4373  //ideal G1;
     4374  //ring endRing;
     4375  ring newRing, oldRing;
     4376  intvec* ivNull = new intvec(nV);
     4377  intvec* exivlp = Mivlp(nV);
     4378#ifndef BUCHBERGER_ALG
     4379  intvec* hilb_func;
     4380#endif
     4381  intvec* tmp_weight = new intvec(nV);
     4382  for(i=nV-1; i>=0; i--)
     4383    (*tmp_weight)[i] = (*curr_weight)[i];
     4384
     4385   // to avoid (1,0,...,0) as the target vector
     4386  intvec* last_omega = new intvec(nV);
     4387  for(i=nV-1; i>0; i--)
     4388    (*last_omega)[i] = 1;
     4389  (*last_omega)[0] = 10000;
     4390
     4391  ring XXRing = currRing;
     4392
     4393  to = clock();
     4394  // the monomial ordering of this current ring would be "dp"
     4395  G = MstdCC(Go);
     4396  tostd = clock()-to;
     4397
     4398  if(currRing->order[0] == ringorder_a)
     4399    goto NEXT_VECTOR;
     4400
     4401  while(1)
     4402  {
     4403    nwalk ++;
     4404    nstep ++;
     4405    to = clock();
     4406    // compute an initial form ideal of <G> w.r.t. "curr_vector"
     4407    Gomega = MwalkInitialForm(G, curr_weight);
     4408    tif = tif + clock()-to;
     4409    oldRing = currRing;
     4410
     4411    if(endwalks == 1)
     4412    {
     4413      /* compute a reduced Groebner basis of Gomega w.r.t. >>_cw by
     4414         the recursive changed perturbation walk alg. */
     4415      tim = clock();
     4416      /*
     4417        Print("\n// **** Grï¿œbnerwalk took %d steps and ", nwalk);
     4418        PrintS("\n// **** call the rec. Pert. Walk to compute a red GB of:");
     4419        idElements(Gomega, "G_omega");
     4420      */
     4421
     4422      if(MivSame(exivlp, target_weight)==1)
     4423        M = REC_GB_Mwalk(idCopy(Gomega), tmp_weight, curr_weight, 2,1);
     4424      else
     4425        goto NORMAL_GW;
     4426      /*
     4427        Print("\n//  time for the last std(Gw)  = %.2f sec",
     4428        ((double) (clock()-tim)/1000000));
     4429        PrintS("\n// ***************************************************\n");
     4430      */
     4431#ifdef CHECK_IDEAL_MWALK
     4432      idElements(Gomega, "G_omega");
     4433      headidString(Gomega, "Gw");
     4434      idElements(M, "M");
     4435      //headidString(M, "M");
     4436#endif
     4437      to = clock();
     4438      F = MLifttwoIdeal(Gomega, M, G);
     4439      xtlift = xtlift + clock() - to;
     4440
     4441      idDelete(&Gomega);
     4442      idDelete(&M);
     4443      idDelete(&G);
     4444
     4445      oldRing = currRing;
     4446
     4447      /* create a new ring newRing */
     4448       if (rParameter(currRing) != NULL)
     4449       {
     4450         DefRingPar(curr_weight);
     4451       }
     4452       else
     4453       {
     4454         VMrDefault(curr_weight);
     4455       }
     4456      newRing = currRing;
     4457      F1 = idrMoveR(F, oldRing,currRing);
     4458    }
     4459    else
     4460    {
     4461    NORMAL_GW:
     4462#ifndef  BUCHBERGER_ALG
     4463      if(isNolVector(curr_weight) == 0)
     4464      {
     4465        hilb_func = hFirstSeries(Gomega,NULL,NULL,curr_weight,currRing);
     4466      }
     4467      else
     4468      {
     4469        hilb_func = hFirstSeries(Gomega,NULL,NULL,last_omega,currRing);
     4470      }
     4471#endif // BUCHBERGER_ALG
     4472
     4473      // define a new ring that its ordering is "(a(curr_weight),lp)
     4474      if (rParameter(currRing) != NULL)
     4475      {
     4476        DefRingPar(curr_weight);
     4477      }
     4478      else
     4479      {
     4480        VMrDefault(curr_weight);
     4481      }
     4482      newRing = currRing;
     4483      Gomega1 = idrMoveR(Gomega, oldRing,currRing);
     4484
     4485      to = clock();
     4486      // compute a reduced Groebner basis of <Gomega> w.r.t. "newRing"
     4487#ifdef  BUCHBERGER_ALG
     4488      M = MstdhomCC(Gomega1);
     4489#else
     4490      M=kStd(Gomega1,NULL,isHomog,NULL,hilb_func,0,NULL,curr_weight);
     4491      delete hilb_func;
     4492#endif // BUCHBERGER_ALG
     4493      tstd = tstd + clock() - to;
     4494
     4495      // change the ring to oldRing
     4496      rChangeCurrRing(oldRing);
     4497      M1 =  idrMoveR(M, newRing,currRing);
     4498      Gomega2 =  idrMoveR(Gomega1, newRing,currRing);
     4499
     4500      to = clock();
     4501      // compute a representation of the generators of submod (M) with respect to those of mod (Gomega). Gomega is a reduced Groebner basis w.r.t. the current ring.
     4502      F = MLifttwoIdeal(Gomega2, M1, G);
     4503      tlift = tlift + clock() - to;
     4504
     4505      idDelete(&M1);
     4506      idDelete(&Gomega2);
     4507      idDelete(&G);
     4508
     4509      // change the ring to newRing
     4510      rChangeCurrRing(newRing);
     4511      F1 = idrMoveR(F, oldRing,currRing);
     4512    }
     4513
     4514    to = clock();
     4515    // reduce the Groebner basis <G> w.r.t. new ring
     4516    G = kInterRedCC(F1, NULL);
     4517    if(endwalks != 1)
     4518    {
     4519      tred = tred + clock() - to;
     4520    }
     4521    else
     4522    {
     4523      xtred = xtred + clock() - to;
     4524    }
     4525    idDelete(&F1);
     4526    if(endwalks == 1)
     4527    {
     4528      break;
     4529    }
     4530  NEXT_VECTOR:
     4531    to = clock();
     4532    // compute a next weight vector
     4533    intvec* next_weight = MkInterRedNextWeight(curr_weight,target_weight,G);
     4534    tnw = tnw + clock() - to;
     4535#ifdef PRINT_VECTORS
     4536    MivString(curr_weight, target_weight, next_weight);
     4537#endif
     4538
     4539    //if(test_w_in_ConeCC(G, next_weight) != 1)
     4540    if(Overflow_Error == TRUE)
     4541    {
     4542      newRing = currRing;
     4543      PrintS("\n// ** The computed vector does NOT stay in Cone!!\n");
     4544
     4545      if (rParameter(currRing) != NULL)
     4546      {
     4547        DefRingPar(target_weight);
     4548      }
     4549      else
     4550      {
     4551        VMrDefault(target_weight);
     4552      }
     4553      F1 = idrMoveR(G, newRing,currRing);
     4554      G = MstdCC(F1);
     4555      idDelete(&F1);
     4556
     4557      newRing = currRing;
     4558      break;
     4559    }
     4560
     4561    if(MivComp(next_weight, ivNull) == 1)
     4562    {
     4563      newRing = currRing;
     4564      delete next_weight;
     4565      break;
     4566    }
     4567    if(MivComp(next_weight, target_weight) == 1)
     4568    {
     4569      endwalks = 1;
     4570    }
     4571    for(i=nV-1; i>=0; i--)
     4572    {
     4573      (*tmp_weight)[i] = (*curr_weight)[i];
     4574      (*curr_weight)[i] = (*next_weight)[i];
     4575    }
     4576    delete next_weight;
     4577  }
     4578  rChangeCurrRing(XXRing);
     4579  G = idrMoveR(G, newRing,currRing);
     4580
     4581  delete tmp_weight;
     4582  delete ivNull;
     4583  delete exivlp;
     4584
     4585#ifdef TIME_TEST
     4586  TimeString(tinput, tostd, tif, tstd, tlift, tred, tnw, nstep);
     4587
     4588  Print("\n// pSetm_Error = (%d)", ErrorCheck());
     4589  Print("\n// Overflow_Error? (%d)\n", Overflow_Error);
     4590#endif
     4591  return(G);
     4592}
     4593
     4594// 07.11.2012
     4595// THE RANDOM WALK ALGORITHM
     4596ideal Mrwalk(ideal Go, intvec* curr_weight, intvec* target_weight, int weight_rad, int pert_deg)
     4597{
     4598  Set_Error(FALSE);
     4599  Overflow_Error = FALSE;
     4600  //Print("// pSetm_Error = (%d)", ErrorCheck());
     4601
     4602  clock_t tinput, tostd, tif=0, tstd=0, tlift=0, tred=0, tnw=0;
     4603  xtif=0; xtstd=0; xtlift=0; xtred=0; xtnw=0;
     4604  tinput = clock();
     4605  clock_t tim;
    41984606  nstep=0;
    41994607  int i,nwalk,endwalks = 0;
    4200   int nV = baseRing->N;
    4201 
    4202   ideal Gomega, M, F, Gomega1, Gomega2, M1; //, F1;
    4203   ring newRing;
    4204   ring XXRing = baseRing;
     4608  int nV = currRing->N;
     4609
     4610  ideal Gomega, M, F, Gomega1, Gomega2, M1, F1, G;
     4611  //ideal G1;
     4612  //ring endRing;
     4613  ring newRing, oldRing;
    42054614  intvec* ivNull = new intvec(nV);
    4206   intvec* curr_weight = new intvec(nV);
    4207   intvec* target_weight = new intvec(nV);
    42084615  intvec* exivlp = Mivlp(nV);
    42094616  intvec* tmp_weight = new intvec(nV);
    4210   for(i=0; i<nV; i++)
    4211   {
    4212     (*tmp_weight)[i] = (*target_M)[i];
    4213   }
    4214   for(i=0; i<nV; i++)
    4215   {
    4216     (*curr_weight)[i] = (*orig_M)[i];
    4217     (*target_weight)[i] = (*target_M)[i];
     4617  for(i=nV-1; i>=0; i--)
     4618  {
     4619    (*tmp_weight)[i] = (*curr_weight)[i];
    42184620  }
    42194621#ifndef BUCHBERGER_ALG
     
    42274629  (*last_omega)[0] = 10000;
    42284630#endif
    4229   rComplete(currRing);
    4230 #ifdef CHECK_IDEAL_MWALK
    4231     idString(Go,"Go");
    4232 #endif
    4233 #ifdef TIME_TEST
     4631  ring XXRing = currRing;
     4632
    42344633  to = clock();
    4235 #endif
    4236      if(orig_M->length() == nV)
    4237       {
    4238         newRing = VMrDefault(curr_weight); // define a new ring with ordering "(a(curr_weight),lp)
    4239       }
    4240       else
    4241       {
    4242         newRing = VMatrDefault(orig_M);
    4243       }
    4244   rChangeCurrRing(newRing);
    4245   ideal G = MstdCC(idrMoveR(Go,baseRing,currRing));
    4246   baseRing = currRing;
    4247 #ifdef TIME_TEST
     4634  G = MstdCC(Go);
    42484635  tostd = clock()-to;
    4249 #endif
    4250 
    4251   nwalk = 0;
     4636
     4637  //if(currRing->order[0] == ringorder_a)
     4638  //  goto NEXT_VECTOR;
     4639nwalk = 0;
    42524640  while(1)
    42534641  {
    42544642    nwalk ++;
    42554643    nstep ++;
    4256 #ifdef TIME_TEST
    42574644    to = clock();
    4258 #endif
    42594645#ifdef CHECK_IDEAL_MWALK
    42604646    idString(G,"G");
    42614647#endif
    4262     Gomega = MwalkInitialForm(G, curr_weight); // compute an initial form ideal of <G> w.r.t. "curr_vector"
    4263 #ifdef TIME_TEST
    4264     tif = tif + clock()-to; //time for computing initial form ideal
    4265 #endif
     4648    Gomega = MwalkInitialForm(G, curr_weight);// compute an initial form ideal of <G> w.r.t. "curr_vector"
     4649    tif = tif + clock()-to;
     4650    oldRing = currRing;
    42664651#ifdef CHECK_IDEAL_MWALK
    4267     idString(Gomega,"Gomega");
    4268 #endif
     4652    idString(Gomega,"G_omega");
     4653#endif
     4654    if(endwalks == 1)
     4655    {
     4656      // compute a reduced Groebner basis of Gomega w.r.t. >>_cw by the recursive changed perturbation walk alg.
     4657      tim = clock();
     4658      Print("\n//Mrwalk: Groebnerwalk took %d steps.\n", nwalk);
     4659      //PrintS("\n//Mrwalk: call the rec. Pert. Walk to compute a red GB of:");
     4660      //idString(Gomega, "G_omega");
     4661      M = REC_GB_Mwalk(idCopy(Gomega), tmp_weight, curr_weight,pert_deg,1);
     4662      Print("\n//Mrwalk: time for the last std(Gw)  = %.2f sec\n",((double) (clock()-tim)/1000000));
     4663
     4664#ifdef CHECK_IDEAL_MWALK
     4665      idString(Gomega, "G_omega");
     4666      idString(M, "M");
     4667#endif
     4668      to = clock();
     4669      F = MLifttwoIdeal(Gomega, M, G);
     4670      xtlift = xtlift + clock() - to;
     4671
     4672      idDelete(&Gomega);
     4673      idDelete(&M);
     4674      idDelete(&G);
     4675
     4676      oldRing = currRing;
     4677      VMrDefault(curr_weight); //define a new ring with ordering "(a(curr_weight),lp)
     4678
     4679       /*if (rParameter(currRing) != NULL)
     4680       {
     4681         DefRingPar(curr_weight);
     4682       }
     4683       else
     4684       {
     4685         VMrDefault(curr_weight);
     4686       }*/
     4687      newRing = currRing;
     4688      F1 = idrMoveR(F, oldRing,currRing);
     4689    }
     4690    else
     4691    {
    42694692#ifndef  BUCHBERGER_ALG
    4270     if(isNolVector(curr_weight) == 0)
    4271     {
    4272       hilb_func = hFirstSeries(Gomega,NULL,NULL,curr_weight,currRing);
    4273     }
    4274     else
    4275     {
    4276       hilb_func = hFirstSeries(Gomega,NULL,NULL,last_omega,currRing);
    4277     }
    4278 #endif
    4279     if(nwalk == 1)
    4280     {
    4281       if(orig_M->length() == nV)
    4282       {
    4283         newRing = VMrDefault(curr_weight); // define a new ring with ordering "(a(curr_weight),lp)
     4693      if(isNolVector(curr_weight) == 0)
     4694      {
     4695        hilb_func = hFirstSeries(Gomega,NULL,NULL,curr_weight,currRing);
    42844696      }
    42854697      else
    42864698      {
    4287         newRing = VMatrDefault(orig_M);
    4288       }
    4289     }
    4290     else
    4291     {
    4292      if(target_M->length() == nV)
    4293      {
    4294        newRing = VMrRefine(curr_weight,target_weight); //define a new ring with ordering "(a(curr_weight),Wp(target_weight))"
    4295      }
    4296      else
    4297      {
    4298        newRing = VMatrRefine(target_M,curr_weight);
    4299      }
    4300     }
    4301     rChangeCurrRing(newRing);
    4302     Gomega1 = idrMoveR(Gomega, baseRing,currRing);
    4303     idDelete(&Gomega);
    4304     // compute a reduced Groebner basis of <Gomega> w.r.t. "newRing"
    4305 #ifdef TIME_TEST
     4699        hilb_func = hFirstSeries(Gomega,NULL,NULL,last_omega,currRing);
     4700      }
     4701#endif
     4702/*
     4703      if (rParameter(currRing) != NULL)
     4704      {
     4705        DefRingPar(curr_weight);
     4706      }
     4707      else
     4708      {
     4709        VMrDefault(curr_weight);
     4710      }*/
     4711
     4712      VMrDefault(curr_weight); //define a new ring with ordering "(a(curr_weight),lp)
     4713      newRing = currRing;
     4714      Gomega1 = idrMoveR(Gomega, oldRing,currRing);
     4715#ifdef CHECK_IDEAL_MWALK
     4716      idString(Gomega1, "G_omega1");
     4717#endif
     4718      // compute a reduced Groebner basis of <Gomega> w.r.t. "newRing"
     4719      to = clock();
     4720#ifndef  BUCHBERGER_ALG
     4721      M=kStd(Gomega1,NULL,isHomog,NULL,hilb_func,0,NULL,curr_weight);
     4722      delete hilb_func;
     4723#else
     4724      //M = MstdhomCC(Gomega1);
     4725      //M = MstdCC(Gomega1);
     4726      //M = kStd(Gomega1, NULL, testHomog,NULL, NULL,0,0,curr_weight);
     4727      M = kStd(Gomega1,NULL,testHomog,NULL,NULL,0,0,NULL);
     4728#endif
     4729      tstd = tstd + clock() - to;
     4730#ifdef CHECK_IDEAL_MWALK
     4731      idString(M, "M");
     4732#endif
     4733      //change the ring to oldRing
     4734      rChangeCurrRing(oldRing);
     4735      M1 =  idrMoveR(M, newRing,currRing);
     4736#ifdef CHECK_IDEAL_MWALK
     4737      idString(M1, "M1");
     4738#endif
     4739      Gomega2 =  idrMoveR(Gomega1, newRing,currRing);
     4740
     4741      to = clock();
     4742      // compute a representation of the generators of submod (M) with respect to those of mod (Gomega), where Gomega is a reduced Groebner basis w.r.t. the current ring
     4743      F = MLifttwoIdeal(Gomega2, M1, G);
     4744#ifdef CHECK_IDEAL_MWALK
     4745      idString(F, "F");
     4746#endif
     4747      tlift = tlift + clock() - to;
     4748
     4749      idDelete(&M1);
     4750      idDelete(&Gomega2);
     4751      idDelete(&G);
     4752
     4753      rChangeCurrRing(newRing); // change the ring to newRing
     4754      F1 = idrMoveR(F,oldRing,currRing);
     4755    }
     4756
    43064757    to = clock();
    4307 #endif
    4308 #ifndef  BUCHBERGER_ALG
    4309     M=kStd(Gomega1,NULL,isHomog,NULL,hilb_func,0,NULL,curr_weight);
    4310     delete hilb_func;
    4311 #else
    4312     M = kStd(Gomega1,NULL,testHomog,NULL,NULL,0,0,NULL);
    4313 #endif
    4314 #ifdef TIME_TEST
    4315     tstd = tstd + clock() - to;
    4316 #endif
    4317     idSkipZeroes(M);
    4318 #ifdef CHECK_IDEAL_MWALK
    4319     PrintS("\n//** Mwalk: computed M.\n");
    4320     idString(M, "M");
    4321 #endif
    4322     //change the ring to baseRing
    4323     rChangeCurrRing(baseRing);
    4324     M1 =  idrMoveR(M, newRing,currRing);
    4325     idDelete(&M);
    4326     Gomega2 = idrMoveR(Gomega1, newRing,currRing);
    4327     idDelete(&Gomega1);
    4328 #ifdef TIME_TEST
    4329     to = clock();
    4330 #endif
    4331     // compute a representation of the generators of submod (M) with respect to those of mod (Gomega), where Gomega is a reduced Groebner basis w.r.t. the current ring
    4332     F = MLifttwoIdeal(Gomega2, M1, G);
    4333 #ifdef TIME_TEST
    4334     tlift = tlift + clock() - to;
    4335 #endif
    4336 #ifdef CHECK_IDEAL_MWALK
    4337     idString(F, "F");
    4338 #endif
    4339     idDelete(&Gomega2);
    4340     idDelete(&M1);
    4341     rChangeCurrRing(newRing); // change the ring to newRing
    4342     G = idrMoveR(F,baseRing,currRing);
    4343     idDelete(&F);
    4344     baseRing = currRing;
    4345 #ifdef TIME_TEST
    4346     to = clock();
    4347 #endif
    4348     //G = kStd(F1,NULL,testHomog,NULL,NULL,0,0,NULL);
    4349 #ifdef TIME_TEST
    4350     tstd = tstd + clock() - to;
    4351 #endif
    4352     idSkipZeroes(G);
     4758    G = kInterRedCC(F1, NULL); //reduce the Groebner basis <G> w.r.t. new ring
    43534759#ifdef CHECK_IDEAL_MWALK
    43544760    idString(G, "G");
    43554761#endif
    4356 #ifdef TIME_TEST
     4762    idDelete(&F1);
     4763    if(endwalks == 1)
     4764    {
     4765      xtred = xtred + clock() - to;
     4766      break;
     4767    }
     4768    else
     4769    {
     4770      tred = tred + clock() - to;
     4771    }
     4772
     4773  //NEXT_VECTOR:
    43574774    to = clock();
    4358 #endif
    4359     intvec* next_weight = MwalkNextWeightCC(curr_weight,target_weight,G);
    4360 #ifdef TIME_TEST
    4361     tnw = tnw + clock() - to;
    4362 #endif
    4363 #ifdef PRINT_VECTORS
    4364     MivString(curr_weight, target_weight, next_weight);
    4365 #endif
    4366     if(MivComp(next_weight, ivNull) == 1 || test_w_in_ConeCC(G, target_weight) == 1 || MivComp(target_weight,curr_weight) == 1 || MivComp(next_weight,curr_weight) == 1)
    4367     {
    4368 #ifdef CHECK_IDEAL_MWALK
    4369       PrintS("\n//** Mwalk: entering last cone.\n");
    4370 #endif
    4371       Gomega = MwalkInitialForm(G, curr_weight); // compute an initial form ideal of <G> w.r.t. "curr_vector"
    4372       if(target_M->length() == nV)
    4373       {
    4374         newRing = VMrDefault(target_weight); // define a new ring with ordering "(a(curr_weight),lp)
     4775    //intvec* next_weight = MkInterRedNextWeight(curr_weight,target_weight,G);
     4776   intvec* next_weight = MWalkRandomNextWeight(G, curr_weight, target_weight, weight_rad, pert_deg);
     4777
     4778  tnw = tnw + clock() - to;
     4779//#ifdef PRINT_VECTORS
     4780  MivString(curr_weight, target_weight, next_weight);
     4781//#endif
     4782    if(Overflow_Error == TRUE)
     4783    {
     4784      newRing = currRing;
     4785      PrintS("\n//Mrwalk: The computed vector does NOT stay in Cone!!\n");
     4786
     4787      if (rParameter(currRing) != NULL)
     4788      {
     4789        DefRingPar(target_weight);
    43754790      }
    43764791      else
    43774792      {
    4378         newRing = VMatrDefault(target_M);
    4379       }
    4380       rChangeCurrRing(newRing);
    4381       Gomega1 = idrMoveR(Gomega, baseRing,currRing);
    4382       idDelete(&Gomega);
    4383 #ifdef CHECK_IDEAL_MWALK
    4384       idString(Gomega1, "Gomega");
    4385 #endif
    4386       M = kStd(Gomega1,NULL,testHomog,NULL,NULL,0,0,NULL);
    4387 #ifdef CHECK_IDEAL_MWALK
    4388       idString(M,"M");
    4389 #endif
    4390       rChangeCurrRing(baseRing);
    4391       M1 =  idrMoveR(M, newRing,currRing);
    4392       idDelete(&M);
    4393       Gomega2 = idrMoveR(Gomega1, newRing,currRing);
    4394       idDelete(&Gomega1);
    4395       F = MLifttwoIdeal(Gomega2, M1, G);
    4396 #ifdef CHECK_IDEAL_MWALK
    4397       idString(F,"F");
    4398 #endif
    4399       idDelete(&Gomega2);
    4400       idDelete(&M1);
    4401       rChangeCurrRing(newRing); // change the ring to newRing
    4402       G = idrMoveR(F,baseRing,currRing);
    4403       idDelete(&F);
    4404       baseRing = currRing;
    4405       si_opt_1 = save1; //set original options, e. g. option(RedSB)
    4406       idSkipZeroes(G);
    4407 #ifdef TIME_TEST
    4408       to = clock();
    4409 #endif
    4410       if(si_opt_1 == (Sy_bit(OPT_REDSB)))
    4411       {
    4412         G = kInterRedCC(G,NULL); //reduce the Groebner basis <G> w.r.t. currRing, if option(redSB) is set
    4413       }
    4414 #ifdef TIME_TEST
    4415       tred = tred + clock() - to;
    4416 #endif
    4417       idSkipZeroes(G);
     4793        VMrDefault(target_weight);  //define a new ring with ordering "(a(curr_weight),lp)
     4794      }
     4795      F1 = idrMoveR(G, newRing,currRing);
     4796      G = MstdCC(F1);
     4797      idDelete(&F1);
     4798
     4799      newRing = currRing;
     4800      break;
     4801    }
     4802
     4803    if(MivComp(next_weight, ivNull) == 1)
     4804    {
     4805      newRing = currRing;
    44184806      delete next_weight;
    44194807      break;
    4420 #ifdef CHECK_IDEAL_MWALK
    4421       PrintS("\n//** Mwalk: last cone.\n");
    4422 #endif
    4423     }
    4424 #ifdef CHECK_IDEAL_MWALK
    4425     PrintS("\n//** Mwalk: update weight vectors.\n");
    4426 #endif
     4808    }
     4809    if(MivComp(next_weight, target_weight) == 1)
     4810    {
     4811      endwalks = 1;
     4812    }
    44274813    for(i=nV-1; i>=0; i--)
    44284814    {
     
    44334819  }
    44344820  rChangeCurrRing(XXRing);
    4435   ideal result = idrMoveR(G,baseRing,currRing);
    4436   idDelete(&G);
    4437 /*#ifdef CHECK_IDEAL_MWALK
    4438   pDelete(&p);
    4439 #endif*/
     4821  G = idrMoveR(G, newRing,currRing);
     4822
    44404823  delete tmp_weight;
    44414824  delete ivNull;
     
    44454828#endif
    44464829#ifdef TIME_TEST
    4447   Print("\n//** Mwalk: Groebner Walk took %d steps.\n", nstep);
    44484830  TimeString(tinput, tostd, tif, tstd, tlift, tred, tnw, nstep);
    4449   Print("\n//** Mwalk: Ergebnis.\n");
    4450   //Print("\n// pSetm_Error = (%d)", ErrorCheck());
    4451   //Print("\n// Overflow_Error? (%d)\n", Overflow_Error);
    4452 #endif
    4453   return(result);
    4454 }
    4455 
    4456 /*****************************
    4457  * THE RANDOM WALK ALGORITHM *
    4458  *****************************/
    4459 ideal Mrwalk(ideal Go, intvec* orig_M, intvec* target_M, int weight_rad, int pert_deg, ring baseRing)
    4460 {
    4461   BITSET save1 = si_opt_1; // save current options
    4462   si_opt_1 &= (~Sy_bit(OPT_REDSB)); // no reduced Groebner basis
    4463   si_opt_1 &= (~Sy_bit(OPT_REDTAIL)); // not tail reductions
    4464   //si_opt_1|=(Sy_bit(OPT_REDTAIL)|Sy_bit(OPT_REDSB));
    4465   Set_Error(FALSE);
     4831
     4832  Print("\n// pSetm_Error = (%d)", ErrorCheck());
     4833  Print("\n// Overflow_Error? (%d)\n", Overflow_Error);
     4834#endif
     4835  return(G);
     4836}
     4837
     4838//unused
     4839#if 0
     4840ideal Mwalk_tst(ideal Go, intvec* curr_weight, intvec* target_weight)
     4841{
     4842  //clock_t tinput=clock();
     4843  //idString(Go,"Ginp");
     4844  int i, nV = currRing->N;
     4845  int nwalk=0, endwalks=0;
     4846
     4847  ideal Gomega, M, F, Gomega1, Gomega2, M1, F1, G;
     4848  // ideal G1; ring endRing;
     4849  ring newRing, oldRing;
     4850  intvec* ivNull = new intvec(nV);
     4851  ring XXRing = currRing;
     4852
     4853  intvec* tmp_weight = new intvec(nV);
     4854  for(i=nV-1; i>=0; i--)
     4855  {
     4856    (*tmp_weight)[i] = (*curr_weight)[i];
     4857  }
     4858  /* the monomial ordering of this current ring would be "dp" */
     4859  G = MstdCC(Go);
     4860#ifndef BUCHBERGER_ALG
     4861  intvec* hilb_func;
     4862#endif
     4863  /* to avoid (1,0,...,0) as the target vector */
     4864  intvec* last_omega = new intvec(nV);
     4865  for(i=nV-1; i>0; i--)
     4866    (*last_omega)[i] = 1;
     4867  (*last_omega)[0] = 10000;
     4868
     4869  while(1)
     4870  {
     4871    nwalk ++;
     4872    //Print("\n// Entering the %d-th step:", nwalk);
     4873    //Print("\n// ring r[%d] = %s;", nwalk, rString(currRing));
     4874    idString(G,"G");
     4875    /* compute an initial form ideal of <G> w.r.t. "curr_vector" */
     4876    Gomega = MwalkInitialForm(G, curr_weight);
     4877    //ivString(curr_weight, "omega");
     4878    idString(Gomega,"Gw");
     4879
     4880#ifndef  BUCHBERGER_ALG
     4881    if(isNolVector(curr_weight) == 0)
     4882      hilb_func = hFirstSeries(Gomega,NULL,NULL,curr_weight,currRing);
     4883    else
     4884      hilb_func = hFirstSeries(Gomega,NULL,NULL,last_omega,currRing);
     4885#endif // BUCHBERGER_ALG
     4886
     4887
     4888    oldRing = currRing;
     4889
     4890    /* define a new ring that its ordering is "(a(curr_weight),lp) */
     4891    VMrDefault(curr_weight);
     4892    newRing = currRing;
     4893
     4894    Gomega1 = idrMoveR(Gomega, oldRing,currRing);
     4895
     4896    /* compute a reduced Groebner basis of <Gomega> w.r.t. "newRing" */
     4897#ifdef  BUCHBERGER_ALG
     4898    M = MstdhomCC(Gomega1);
     4899#else
     4900    M=kStd(Gomega1,NULL,isHomog,NULL,hilb_func,0,NULL,curr_weight);
     4901    delete hilb_func;
     4902#endif // BUCHBERGER_ALG
     4903
     4904    idString(M,"M");
     4905
     4906      /* change the ring to oldRing */
     4907    rChangeCurrRing(oldRing);
     4908    M1 =  idrMoveR(M, newRing,currRing);
     4909    Gomega2 =  idrMoveR(Gomega1, newRing,currRing);
     4910
     4911      /* compute a representation of the generators of submod (M)
     4912         with respect to those of mod (Gomega).
     4913         Gomega is a reduced Groebner basis w.r.t. the current ring */
     4914    F = MLifttwoIdeal(Gomega2, M1, G);
     4915    idDelete(&M1);
     4916    idDelete(&Gomega2);
     4917    idDelete(&G);
     4918    idString(F,"F");
     4919
     4920    /* change the ring to newRing */
     4921    rChangeCurrRing(newRing);
     4922    F1 = idrMoveR(F, oldRing,currRing);
     4923
     4924    /* reduce the Groebner basis <G> w.r.t. new ring */
     4925    G = kInterRedCC(F1, NULL);
     4926    //idSkipZeroes(G);//done by kInterRed
     4927    idDelete(&F1);
     4928    idString(G,"G");
     4929    if(endwalks == 1)
     4930      break;
     4931
     4932    /* compute a next weight vector */
     4933    intvec* next_weight = MkInterRedNextWeight(curr_weight,target_weight,G);
     4934#ifdef PRINT_VECTORS
     4935    MivString(curr_weight, target_weight, next_weight);
     4936#endif
     4937
     4938    if(MivComp(next_weight, ivNull) == 1)
     4939    {
     4940      delete next_weight;
     4941      break;
     4942    }
     4943    if(MivComp(next_weight, target_weight) == 1)
     4944      endwalks = 1;
     4945
     4946    for(i=nV-1; i>=0; i--)
     4947      (*tmp_weight)[i] = (*curr_weight)[i];
     4948
     4949    /* 06.11.01 to free the memory: NOT Changed!!*/
     4950    for(i=nV-1; i>=0; i--)
     4951      (*curr_weight)[i] = (*next_weight)[i];
     4952    delete next_weight;
     4953  }
     4954  rChangeCurrRing(XXRing);
     4955  G = idrMoveR(G, newRing,currRing);
     4956
     4957  delete tmp_weight;
     4958  delete ivNull;
     4959  PrintLn();
     4960  return(G);
     4961}
     4962#endif
     4963
     4964/**************************************************************/
     4965/*     Implementation of the perturbation walk algorithm      */
     4966/**************************************************************/
     4967/* If the perturbed target weight vector or an intermediate weight vector
     4968   doesn't stay in the correct Groebner cone, we have only
     4969   a reduced Groebner basis for the given ideal with respect to
     4970   a monomial order which differs to the given order.
     4971   Then we have to compute the wanted  reduced Groebner basis for it.
     4972   For this, we can use
     4973   1) the improved Buchberger algorithm or
     4974   2) the changed perturbation walk algorithm with a decreased degree.
     4975*/
     4976// use kStd, if nP = 0, else call LastGB
     4977ideal Mpwalk(ideal Go, int op_deg, int tp_deg,intvec* curr_weight,
     4978             intvec* target_weight, int nP)
     4979{
     4980  Set_Error(FALSE  );
    44664981  Overflow_Error = FALSE;
    4467 #ifdef TIME_TEST
    4468   clock_t tinput, tostd, tif=0, tstd=0, tlift=0, tred=0, tnw=0;
     4982  //Print("// pSetm_Error = (%d)", ErrorCheck());
     4983
     4984  clock_t  tinput, tostd, tif=0, tstd=0, tlift=0, tred=0, tnw=0;
     4985  xtextra=0;
    44694986  xtif=0; xtstd=0; xtlift=0; xtred=0; xtnw=0;
    44704987  tinput = clock();
     4988
    44714989  clock_t tim;
    4472 #endif
    4473   nstep=0;
    4474   int i,nwalk,endwalks = 0;
    4475   int nV = baseRing->N;
    4476 
    4477   ideal Gomega, M, F, Gomega1, Gomega2, M1; //, F1;
    4478   ring newRing;
    4479   ring XXRing = baseRing;
    4480   intvec* ivNull = new intvec(nV);
    4481   intvec* curr_weight = new intvec(nV);
    4482   intvec* target_weight = new intvec(nV);
    4483   intvec* exivlp = Mivlp(nV);
    4484   intvec* tmp_weight = new intvec(nV);
    4485 #ifdef CHECK_IDEAL_MWALK
    4486   poly p;
    4487 #endif
    4488   for(i=0; i<nV; i++)
    4489   {
    4490     (*tmp_weight)[i] = (*target_M)[i];
    4491   }
    4492   for(i=0; i<nV; i++)
    4493   {
    4494     (*curr_weight)[i] = (*orig_M)[i];
    4495     (*target_weight)[i] = (*target_M)[i];
    4496   }
    4497 #ifndef BUCHBERGER_ALG
    4498   intvec* hilb_func;
    4499    // to avoid (1,0,...,0) as the target vector
    4500   intvec* last_omega = new intvec(nV);
    4501   for(i=nV-1; i>0; i--)
    4502   {
    4503     (*last_omega)[i] = 1;
    4504   }
    4505   (*last_omega)[0] = 10000;
    4506 #endif
    4507   rComplete(currRing);
    4508 #ifdef CHECK_IDEAL_MWALK
    4509   for(i=0; i<IDELEMS(Go); i++)
    4510   {
    4511     poly p=Go->m[i];
    4512     while(p!=NULL)
    4513     {
    4514       p_Setm(p,currRing);
    4515       pIter(p);
    4516     }
    4517     pDelete(&p);
    4518   }
    4519     idString(Go,"Go");
    4520 #endif
    4521 #ifdef TIME_TEST
    4522   to = clock();
    4523 #endif
    4524      if(orig_M->length() == nV)
    4525       {
    4526         newRing = VMrDefault(curr_weight); // define a new ring with ordering "(a(curr_weight),lp)
    4527       }
    4528       else
    4529       {
    4530         newRing = VMatrDefault(orig_M);
    4531       }
    4532   rChangeCurrRing(newRing);
    4533   ideal G = MstdCC(idrMoveR(Go,baseRing,currRing));
    4534   baseRing = currRing;
    4535 #ifdef TIME_TEST
    4536   tostd = clock()-to;
    4537 #endif
    4538 
    4539   nwalk = 0;
    4540   while(1)
    4541   {
    4542     nwalk ++;
    4543     nstep ++;
    4544 #ifdef TIME_TEST
    4545     to = clock();
    4546 #endif
    4547 #ifdef CHECK_IDEAL_MWALK
    4548     idString(G,"G");
    4549 #endif
    4550     Gomega = MwalkInitialForm(G, curr_weight); // compute an initial form ideal of <G> w.r.t. "curr_vector"
    4551 #ifdef TIME_TEST
    4552     tif = tif + clock()-to; //time for computing initial form ideal
    4553 #endif
    4554 #ifdef CHECK_IDEAL_MWALK
    4555     idString(Gomega,"Gomega");
    4556 #endif
    4557 #ifndef  BUCHBERGER_ALG
    4558     if(isNolVector(curr_weight) == 0)
    4559     {
    4560       hilb_func = hFirstSeries(Gomega,NULL,NULL,curr_weight,currRing);
    4561     }
    4562     else
    4563     {
    4564       hilb_func = hFirstSeries(Gomega,NULL,NULL,last_omega,currRing);
    4565     }
    4566 #endif
    4567     if(nwalk == 1)
    4568     {
    4569       if(orig_M->length() == nV)
    4570       {
    4571         newRing = VMrDefault(curr_weight); // define a new ring with ordering "(a(curr_weight),lp)
    4572       }
    4573       else
    4574       {
    4575         newRing = VMatrDefault(orig_M);
    4576       }
    4577     }
    4578     else
    4579     {
    4580      if(target_M->length() == nV)
    4581      {
    4582        newRing = VMrRefine(curr_weight,target_weight); //define a new ring with ordering "(a(curr_weight),Wp(target_weight))"
    4583      }
    4584      else
    4585      {
    4586        newRing = VMatrRefine(target_M,curr_weight);
    4587      }
    4588     }
    4589     rChangeCurrRing(newRing);
    4590     Gomega1 = idrMoveR(Gomega, baseRing,currRing);
    4591     idDelete(&Gomega);
    4592     // compute a reduced Groebner basis of <Gomega> w.r.t. "newRing"
    4593 #ifdef TIME_TEST
    4594     to = clock();
    4595 #endif
    4596 #ifndef  BUCHBERGER_ALG
    4597     M=kStd(Gomega1,NULL,isHomog,NULL,hilb_func,0,NULL,curr_weight);
    4598     delete hilb_func;
    4599 #else
    4600     M = kStd(Gomega1,NULL,testHomog,NULL,NULL,0,0,NULL);
    4601 #endif
    4602 #ifdef TIME_TEST
    4603     tstd = tstd + clock() - to;
    4604 #endif
    4605     idSkipZeroes(M);
    4606 #ifdef CHECK_IDEAL_MWALK
    4607     idString(M, "M");
    4608 #endif
    4609     //change the ring to baseRing
    4610     rChangeCurrRing(baseRing);
    4611     M1 =  idrMoveR(M, newRing,currRing);
    4612     idDelete(&M);
    4613     Gomega2 = idrMoveR(Gomega1, newRing,currRing);
    4614     idDelete(&Gomega1);
    4615 #ifdef TIME_TEST
    4616     to = clock();
    4617 #endif
    4618     // compute a representation of the generators of submod (M) with respect to those of mod (Gomega), where Gomega is a reduced Groebner basis w.r.t. the current ring
    4619     F = MLifttwoIdeal(Gomega2, M1, G);
    4620 #ifdef CHECK_IDEAL_MWALK
    4621     idString(F,"F");
    4622 #endif
    4623 #ifdef TIME_TEST
    4624     tlift = tlift + clock() - to;
    4625 #endif
    4626     idDelete(&Gomega2);
    4627     idDelete(&M1);
    4628     rChangeCurrRing(newRing); // change the ring to newRing
    4629     G = idrMoveR(F,baseRing,currRing);
    4630     idDelete(&F);
    4631     baseRing = currRing;
    4632     idSkipZeroes(G);
    4633 #ifdef CHECK_IDEAL_MWALK
    4634     idString(G, "G");
    4635 #endif
    4636 #ifdef TIME_TEST
    4637     to = clock();
    4638 #endif
    4639     intvec* next_weight = MWalkRandomNextWeight(G, curr_weight, target_weight, weight_rad, pert_deg);
    4640 #ifdef TIME_TEST
    4641     tnw = tnw + clock() - to;
    4642 #endif
    4643 #ifdef PRINT_VECTORS
    4644     MivString(curr_weight, target_weight, next_weight);
    4645 #endif
    4646 
    4647     if(MivComp(next_weight, ivNull) == 1 || test_w_in_ConeCC(G, target_weight) == 1 || MivComp(target_weight,curr_weight) == 1 || MivComp(next_weight,curr_weight) == 1)
    4648     {
    4649 #ifdef CHECK_IDEAL_MWALK
    4650       PrintS("\n//** Mrwalk: entering last cone.\n");
    4651 #endif
    4652       Gomega = MwalkInitialForm(G, curr_weight); // compute an initial form ideal of <G> w.r.t. "curr_vector"
    4653       if(target_M->length() == nV)
    4654       {
    4655         newRing = VMrDefault(target_weight); // define a new ring with ordering "(a(curr_weight),lp)
    4656       }
    4657       else
    4658       {
    4659         newRing = VMatrDefault(target_M);
    4660       }
    4661       rChangeCurrRing(newRing);
    4662       Gomega1 = idrMoveR(Gomega, baseRing,currRing);
    4663       idDelete(&Gomega);
    4664       M = kStd(Gomega1,NULL,testHomog,NULL,NULL,0,0,NULL);
    4665       rChangeCurrRing(baseRing);
    4666       M1 =  idrMoveR(M, newRing,currRing);
    4667       idDelete(&M);
    4668       Gomega2 = idrMoveR(Gomega1, newRing,currRing);
    4669       idDelete(&Gomega1);
    4670       F = MLifttwoIdeal(Gomega2, M1, G);
    4671       idDelete(&Gomega2);
    4672       idDelete(&M1);
    4673       rChangeCurrRing(newRing); // change the ring to newRing
    4674       G = idrMoveR(F,baseRing,currRing);
    4675       idDelete(&F);
    4676       baseRing = currRing;
    4677       si_opt_1 = save1; //set original options, e. g. option(RedSB)
    4678       idSkipZeroes(G);
    4679 #ifdef TIME_TEST
    4680       to = clock();
    4681 #endif
    4682       if(si_opt_1 == (Sy_bit(OPT_REDSB)))
    4683       {
    4684         G = kInterRedCC(G,NULL); //reduce the Groebner basis <G> w.r.t. currRing, if option(redSB) is set
    4685       }
    4686 #ifdef TIME_TEST
    4687       tred = tred + clock() - to;
    4688 #endif
    4689       idSkipZeroes(G);
    4690       delete next_weight;
    4691       break;
    4692     }
    4693     for(i=nV-1; i>=0; i--)
    4694     {
    4695       (*tmp_weight)[i] = (*curr_weight)[i];
    4696       (*curr_weight)[i] = (*next_weight)[i];
    4697     }
    4698     delete next_weight;
    4699   }
    4700   rChangeCurrRing(XXRing);
    4701   ideal result = idrMoveR(G,baseRing,currRing);
    4702   idDelete(&G);
    4703 #ifdef CHECK_IDEAL_MWALK
    4704   pDelete(&p);
    4705 #endif
    4706   delete tmp_weight;
    4707   delete ivNull;
    4708   delete exivlp;
    4709 #ifndef BUCHBERGER_ALG
    4710   delete last_omega;
    4711 #endif
    4712 #ifdef TIME_TEST
    4713   Print("\n//** Mrwalk: Groebner Walk took %d steps.\n", nstep);
    4714   TimeString(tinput, tostd, tif, tstd, tlift, tred, tnw, nstep);
    4715 #endif
    4716   return(result);
    4717 }
    4718 
    4719 /**************************************************************************
    4720  * Implementation of the perturbation walk algorithm                      *
    4721  *                                                                        *
    4722  * If the perturbed target weight vector or an intermediate weight vector *
    4723  * doesn't stay in the correct Groebner cone, we have only                *
    4724  * a reduced Groebner basis for the given ideal with respect to           *
    4725  * a monomial order which differs to the given order.                     *
    4726  * Then we have to compute the wanted  reduced Groebner basis for it.     *
    4727  * For this, we can use                                                   *
    4728  * 1) the improved Buchberger algorithm or                                *
    4729  * 2) the changed perturbation walk algorithm with a decreased degree.    *
    4730  **************************************************************************/
    4731 
    4732 /***********************************
    4733  * THE PERTURBATION WALK ALGORITHM *
    4734  ***********************************/
    4735 ideal Mpwalk(ideal Go, intvec* curr_weight, intvec* target_weight, int op_deg, int tp_deg, ring baseRing)
    4736 {
    4737   BITSET save1 = si_opt_1; // save current options
    4738   si_opt_1 &= (~Sy_bit(OPT_REDSB)); // no reduced Groebner basis
    4739   Set_Error(FALSE);
    4740   Overflow_Error = FALSE;
    4741 #ifdef TIME_TEST
    4742   clock_t tinput=0, tostd=0, tif=0, tstd=0, tlift=0, tred=0, tnw=0;
    4743   xtif=0; xtstd=0; xtlift=0; xtred=0; xtnw=0;
    4744   tinput = clock();
    4745   clock_t tim;
    4746 #endif
    4747   int i,nwalk,nV = baseRing->N;
    4748   ideal G, Gomega, M, F, Gomega1, Gomega2, M1;
    4749   ring newRing;
    4750   ring XXRing = baseRing;
     4990
     4991  nstep = 0;
     4992  int i, ntwC=1, ntestw=1,  nV = currRing->N;
     4993  int endwalks=0;
     4994
     4995  ideal Gomega, M, F, G, Gomega1, Gomega2, M1,F1,Eresult,ssG;
     4996  ring newRing, oldRing, TargetRing;
     4997  intvec* iv_M_dp;
     4998  intvec* iv_M_lp;
    47514999  intvec* exivlp = Mivlp(nV);
    47525000  intvec* orig_target = target_weight;
    47535001  intvec* pert_target_vector = target_weight;
    47545002  intvec* ivNull = new intvec(nV);
    4755   intvec* tmp_weight = new intvec(nV);
    4756 #ifdef CHECK_IDEAL_MWALK
    4757   poly p;
    4758 #endif
    4759   for(i=0; i<nV; i++)
    4760   {
    4761     (*tmp_weight)[i] = (*curr_weight)[i];
    4762   }
     5003  intvec* iv_dp = MivUnit(nV);// define (1,1,...,1)
    47635004#ifndef BUCHBERGER_ALG
    47645005  intvec* hilb_func;
    4765    // to avoid (1,0,...,0) as the target vector
     5006#endif
     5007  intvec* next_weight;
     5008
     5009  // to avoid (1,0,...,0) as the target vector
    47665010  intvec* last_omega = new intvec(nV);
    4767   for(i=0 i<nV; i++)
    4768   {
     5011  for(i=nV-1; i>0; i--)
    47695012    (*last_omega)[i] = 1;
    4770   }
    47715013  (*last_omega)[0] = 10000;
    4772 #endif
    4773   idString(Go,"Go");
    4774   baseRing = currRing;
    4775   newRing = VMrDefault(curr_weight);
    4776   rChangeCurrRing(newRing);
    4777   G = idrMoveR(Go,baseRing,currRing);
    4778 #ifdef TIME_TEST
     5014
     5015  ring XXRing = currRing;
     5016
     5017
    47795018  to = clock();
    4780 #endif
    4781   G = kStd(G,NULL,testHomog,NULL,NULL,0,0,NULL);
    4782   idSkipZeroes(G);
    4783 #ifdef TIME_TEST
    4784   tostd = tostd + to - clock();
    4785 #endif
    4786 #ifdef CHECK_IDEAL_MWALK
    4787   idString(G,"G");
    4788 #endif
    4789   if(op_deg >1)
    4790   {
    4791     if(MivComp(curr_weight,MivUnit(nV)) == 1) //ring order is "dp"
    4792     {
    4793       curr_weight = MPertVectors(G, MivMatrixOrderdp(nV), op_deg);
    4794     }
    4795     else //ring order is not "dp"
    4796     {
    4797       curr_weight = MPertVectors(G, MivMatrixOrder(curr_weight), op_deg);
    4798     }
    4799   }
    4800   baseRing = currRing;
     5019  /* perturbs the original vector */
     5020  if(MivComp(curr_weight, iv_dp) == 1) //rOrdStr(currRing) := "dp"
     5021  {
     5022    G = MstdCC(Go);
     5023    tostd = clock()-to;
     5024    if(op_deg != 1){
     5025      iv_M_dp = MivMatrixOrderdp(nV);
     5026      //ivString(iv_M_dp, "iv_M_dp");
     5027      curr_weight = MPertVectors(G, iv_M_dp, op_deg);
     5028    }
     5029  }
     5030  else
     5031  {
     5032    //define ring order := (a(curr_weight),lp);
     5033    if (rParameter(currRing) != NULL)
     5034      DefRingPar(curr_weight);
     5035    else
     5036      VMrDefault(curr_weight);
     5037
     5038    G = idrMoveR(Go, XXRing,currRing);
     5039    G = MstdCC(G);
     5040    tostd = clock()-to;
     5041    if(op_deg != 1){
     5042      iv_M_dp = MivMatrixOrder(curr_weight);
     5043      curr_weight = MPertVectors(G, iv_M_dp, op_deg);
     5044    }
     5045  }
     5046  delete iv_dp;
     5047  if(op_deg != 1) delete iv_M_dp;
     5048
     5049  ring HelpRing = currRing;
     5050
     5051  /* perturbs the target weight vector */
    48015052  if(tp_deg > 1 && tp_deg <= nV)
    48025053  {
     5054    if (rParameter(currRing) != NULL)
     5055      DefRingPar(target_weight);
     5056    else
     5057      VMrDefault(target_weight);
     5058
     5059    TargetRing = currRing;
     5060    ssG = idrMoveR(G,HelpRing,currRing);
     5061    if(MivSame(target_weight, exivlp) == 1)
     5062    {
     5063      iv_M_lp = MivMatrixOrderlp(nV);
     5064      //ivString(iv_M_lp, "iv_M_lp");
     5065      //target_weight = MPertVectorslp(ssG, iv_M_lp, tp_deg);
     5066      target_weight = MPertVectors(ssG, iv_M_lp, tp_deg);
     5067    }
     5068    else
     5069    {
     5070      iv_M_lp = MivMatrixOrder(target_weight);
     5071      //target_weight = MPertVectorslp(ssG, iv_M_lp, tp_deg);
     5072      target_weight = MPertVectors(ssG, iv_M_lp, tp_deg);
     5073    }
     5074    delete iv_M_lp;
    48035075    pert_target_vector = target_weight;
    4804   }
    4805 #ifdef TIME_TEST
    4806   Print("\n//** Mpwalk: perturbation walk with degrees (%d,%d):",op_deg,tp_deg);
    4807   ivString(curr_weight, "new curr_weight");
    4808   ivString(target_weight, "new target_weight");
    4809 #endif
    4810   nwalk = 0;
     5076    rChangeCurrRing(HelpRing);
     5077    G = idrMoveR(ssG, TargetRing,currRing);
     5078  }
     5079  /*
     5080    Print("\n// Perturbationwalkalg. vom Gradpaar (%d,%d):",op_deg,tp_deg);
     5081    ivString(curr_weight, "new sigma");
     5082    ivString(target_weight, "new tau");
     5083  */
    48115084  while(1)
    48125085  {
    4813     nwalk ++;
    4814 #ifdef CHECK_IDEAL_MWALK
    4815     idString(G,"G");
    4816 #endif
    4817 #ifdef TIME_TEST
     5086    nstep ++;
    48185087    to = clock();
    4819 #endif
    4820     Gomega = MwalkInitialForm(G, curr_weight); // compute an initial form ideal of <G> w.r.t. "curr_vector"
    4821 #ifdef TIME_TEST
    4822     tif = tif + clock()-to; //time for computing initial form ideal
    4823 #endif
    4824 #ifdef CHECK_IDEAL_MWALK
    4825     idString(Gomega,"G_omega");
    4826 #endif
     5088    /* compute an initial form ideal of <G> w.r.t. the weight vector
     5089       "curr_weight" */
     5090    Gomega = MwalkInitialForm(G, curr_weight);
     5091
     5092
     5093#ifdef ENDWALKS
     5094    if(endwalks == 1){
     5095      Print("\n// ring r%d = %s;\n", nstep, rString(currRing));
     5096      idElements(G, "G");
     5097      // idElements(Gomega, "Gw");
     5098      headidString(G, "G");
     5099      //headidString(Gomega, "Gw");
     5100    }
     5101#endif
     5102
     5103    tif = tif + clock()-to;
     5104
    48275105#ifndef  BUCHBERGER_ALG
    48285106    if(isNolVector(curr_weight) == 0)
    4829     {
    48305107      hilb_func = hFirstSeries(Gomega,NULL,NULL,curr_weight,currRing);
    4831     }
    48325108    else
    4833     {
    48345109      hilb_func = hFirstSeries(Gomega,NULL,NULL,last_omega,currRing);
    4835     }
    4836 #endif
    4837     if(nwalk == 1)
    4838     {
    4839       newRing = VMrDefault(curr_weight); // define a new ring with ordering "(a(curr_weight),lp)
    4840     }
     5110#endif // BUCHBERGER_ALG
     5111
     5112    oldRing = currRing;
     5113
     5114    // define a new ring with ordering "(a(curr_weight),lp)
     5115    if (rParameter(currRing) != NULL)
     5116      DefRingPar(curr_weight);
    48415117    else
    4842     {
    4843       newRing = VMrRefine(curr_weight,target_weight); //define a new ring with ordering "(a(curr_weight),Wp(target_weight))"
    4844     }
    4845     rChangeCurrRing(newRing);
    4846     Gomega1 = idrMoveR(Gomega, baseRing,currRing);
    4847     idDelete(&Gomega);
    4848 #ifdef CHECK_IDEAL_MWALK
    4849     idString(Gomega1, "Gomega1");
    4850 #endif
    4851     // compute a Groebner basis of <Gomega> w.r.t. "newRing"
    4852 #ifdef TIME_TEST
     5118      VMrDefault(curr_weight);
     5119
     5120    newRing = currRing;
     5121    Gomega1 = idrMoveR(Gomega, oldRing,currRing);
     5122
     5123#ifdef ENDWALKS
     5124    if(endwalks==1)
     5125    {
     5126      Print("\n// ring r%d = %s;\n", nstep, rString(currRing));
     5127      idElements(Gomega1, "Gw");
     5128      headidString(Gomega1, "headGw");
     5129      PrintS("\n// compute a rGB of Gw:\n");
     5130
     5131#ifndef  BUCHBERGER_ALG
     5132      ivString(hilb_func, "w");
     5133#endif
     5134    }
     5135#endif
     5136
     5137    tim = clock();
    48535138    to = clock();
    4854 #endif
    4855 #ifndef  BUCHBERGER_ALG
     5139    /* compute a reduced Groebner basis of <Gomega> w.r.t. "newRing" */
     5140#ifdef  BUCHBERGER_ALG
     5141    M = MstdhomCC(Gomega1);
     5142#else
    48565143    M=kStd(Gomega1,NULL,isHomog,NULL,hilb_func,0,NULL,curr_weight);
    48575144    delete hilb_func;
    4858 #else
    4859     M = kStd(Gomega1,NULL,testHomog,NULL,NULL,0,0,NULL);
    4860 #endif
    4861     idSkipZeroes(M);
    4862 #ifdef TIME_TEST
    4863     tstd = tstd + clock() - to;
    4864 #endif
    4865 #ifdef CHECK_IDEAL_MWALK
    4866     idString(M, "M");
    4867 #endif
    4868     //change the ring to baseRing
    4869     rChangeCurrRing(baseRing);
     5145#endif // BUCHBERGER_ALG
     5146
     5147    if(endwalks == 1){
     5148      xtstd = xtstd+clock()-to;
     5149#ifdef ENDWALKS
     5150      Print("\n// time for the last std(Gw)  = %.2f sec\n",
     5151            ((double) clock())/1000000 -((double)tim) /1000000);
     5152#endif
     5153    }
     5154    else
     5155      tstd=tstd+clock()-to;
     5156
     5157    /* change the ring to oldRing */
     5158    rChangeCurrRing(oldRing);
    48705159    M1 =  idrMoveR(M, newRing,currRing);
    4871     idDelete(&M);
    4872 #ifdef CHECK_IDEAL_MWALK
    4873     idString(M1, "M1");
    4874 #endif
    4875     Gomega2 = idrMoveR(Gomega1, newRing,currRing);
    4876     idDelete(&Gomega1);
    4877 #ifdef CHECK_IDEAL_MWALK
    4878     idString(Gomega2, "G_omega2");
    4879 #endif
    4880     to = clock();
    4881     // compute a representation of the generators of submod (M) with respect to those of mod (Gomega), where Gomega is a reduced Groebner basis w.r.t. the current ring
     5160    Gomega2 =  idrMoveR(Gomega1, newRing,currRing);
     5161
     5162    //if(endwalks==1)  PrintS("\n// Lifting is working:..");
     5163
     5164    to=clock();
     5165    /* compute a representation of the generators of submod (M)
     5166       with respect to those of mod (Gomega).
     5167       Gomega is a reduced Groebner basis w.r.t. the current ring */
    48825168    F = MLifttwoIdeal(Gomega2, M1, G);
    4883     idSkipZeroes(F);
    4884 #ifdef TIME_TEST
    4885     tlift = tlift + clock() - to;
    4886 #endif
    4887 #ifdef CHECK_IDEAL_MWALK
    4888     idString(F,"F");
    4889 #endif
    4890     rChangeCurrRing(newRing); // change the ring to newRing
    4891     G = idrMoveR(F,baseRing,currRing);
    4892     idDelete(&F);
    4893     baseRing = currRing; // set baseRing equal to newRing
    4894 #ifdef CHECK_IDEAL_MWALK
    4895     idString(G,"G");
    4896 #endif
    4897 #ifdef TIME_TEST
    4898     to = clock();
    4899 #endif
    4900     intvec* next_weight = MwalkNextWeightCC(curr_weight,target_weight,G);
    4901 #ifdef TIME_TEST
    4902     tnw = tnw + clock() - to;
    4903 #endif
     5169    if(endwalks != 1)
     5170      tlift = tlift+clock()-to;
     5171    else
     5172      xtlift=clock()-to;
     5173
     5174    idDelete(&M1);
     5175    idDelete(&Gomega2);
     5176    idDelete(&G);
     5177
     5178    /* change the ring to newRing */
     5179    rChangeCurrRing(newRing);
     5180    F1 = idrMoveR(F, oldRing,currRing);
     5181
     5182    //if(endwalks==1)PrintS("\n// InterRed is working now:");
     5183
     5184    to=clock();
     5185    /* reduce the Groebner basis <G> w.r.t. new ring */
     5186    G = kInterRedCC(F1, NULL);
     5187    if(endwalks != 1)
     5188      tred = tred+clock()-to;
     5189    else
     5190      xtred=clock()-to;
     5191
     5192    idDelete(&F1);
     5193
     5194    if(endwalks == 1)
     5195      break;
     5196
     5197    to=clock();
     5198    /* compute a next weight vector */
     5199    next_weight = MkInterRedNextWeight(curr_weight,target_weight, G);
     5200    tnw=tnw+clock()-to;
    49045201#ifdef PRINT_VECTORS
    49055202    MivString(curr_weight, target_weight, next_weight);
    49065203#endif
     5204
    49075205    if(Overflow_Error == TRUE)
    49085206    {
    4909       PrintS("\n//**Mpwalk: OVERFLOW! The computed vector does not stay in cone, the result may be wrong.\n");
     5207      ntwC = 0;
     5208      //ntestomega = 1;
     5209      //Print("\n// ring r%d = %s;\n", nstep, rString(currRing));
     5210      //idElements(G, "G");
    49105211      delete next_weight;
     5212      goto FINISH_160302;
     5213    }
     5214    if(MivComp(next_weight, ivNull) == 1){
     5215      newRing = currRing;
     5216      delete next_weight;
     5217      //Print("\n// ring r%d = %s;\n", nstep, rString(currRing));
    49115218      break;
    49125219    }
    4913     if(test_w_in_ConeCC(G,target_weight) == 1 || MivComp(next_weight, ivNull) == 1)
    4914     {
    4915       delete next_weight;
    4916       break;
    4917     }
    4918     //update tmp_weight and curr_weight
     5220    if(MivComp(next_weight, target_weight) == 1)
     5221      endwalks = 1;
     5222
    49195223    for(i=nV-1; i>=0; i--)
    4920     {
    4921       (*tmp_weight)[i] = (*curr_weight)[i];
    49225224      (*curr_weight)[i] = (*next_weight)[i];
    4923     }
     5225
    49245226    delete next_weight;
    4925   } //end of while-loop
    4926   Print("\n// Mpwalk took %d steps. Ring= %s;\n", nwalk, rString(currRing));
    4927   idSkipZeroes(G);
    4928   si_opt_1 = save1; //set original options, e. g. option(RedSB)
    4929   baseRing = currRing;
    4930   rChangeCurrRing(XXRing);
    4931   ideal Res = idrMoveR(G,baseRing,currRing);
    4932   delete tmp_weight;
     5227  }//while
     5228
     5229  if(tp_deg != 1)
     5230  {
     5231  FINISH_160302:
     5232    if(MivSame(orig_target, exivlp) == 1)
     5233      if (rParameter(currRing) != NULL)
     5234        DefRingParlp();
     5235      else
     5236        VMrDefaultlp();
     5237    else
     5238      if (rParameter(currRing) != NULL)
     5239        DefRingPar(orig_target);
     5240      else
     5241        VMrDefault(orig_target);
     5242
     5243    TargetRing=currRing;
     5244    F1 = idrMoveR(G, newRing,currRing);
     5245#ifdef CHECK_IDEAL
     5246      headidString(G, "G");
     5247#endif
     5248
     5249
     5250    // check whether the pertubed target vector stays in the correct cone
     5251    if(ntwC != 0){
     5252      ntestw = test_w_in_ConeCC(F1, pert_target_vector);
     5253    }
     5254
     5255    if( ntestw != 1 || ntwC == 0)
     5256    {
     5257      /*
     5258      if(ntestw != 1){
     5259        ivString(pert_target_vector, "tau");
     5260        PrintS("\n// ** perturbed target vector doesn't stay in cone!!");
     5261        Print("\n// ring r%d = %s;\n", nstep, rString(currRing));
     5262        idElements(F1, "G");
     5263      }
     5264      */
     5265      // LastGB is "better" than the kStd subroutine
     5266      to=clock();
     5267      ideal eF1;
     5268      if(nP == 0 || tp_deg == 1 || MivSame(orig_target, exivlp) != 1){
     5269        // PrintS("\n// ** calls \"std\" to compute a GB");
     5270        eF1 = MstdCC(F1);
     5271        idDelete(&F1);
     5272      }
     5273      else {
     5274        // PrintS("\n// ** calls \"LastGB\" to compute a GB");
     5275        rChangeCurrRing(newRing);
     5276        ideal F2 = idrMoveR(F1, TargetRing,currRing);
     5277        eF1 = LastGB(F2, curr_weight, tp_deg-1);
     5278        F2=NULL;
     5279      }
     5280      xtextra=clock()-to;
     5281      ring exTargetRing = currRing;
     5282
     5283      rChangeCurrRing(XXRing);
     5284      Eresult = idrMoveR(eF1, exTargetRing,currRing);
     5285    }
     5286    else{
     5287      rChangeCurrRing(XXRing);
     5288      Eresult = idrMoveR(F1, TargetRing,currRing);
     5289    }
     5290  }
     5291  else {
     5292    rChangeCurrRing(XXRing);
     5293    Eresult = idrMoveR(G, newRing,currRing);
     5294  }
    49335295  delete ivNull;
     5296  if(tp_deg != 1)
     5297    delete target_weight;
     5298
     5299  if(op_deg != 1 )
     5300    delete curr_weight;
     5301
    49345302  delete exivlp;
    4935 #ifndef BUCHBERGER_ALG
    49365303  delete last_omega;
    4937 #endif
     5304
    49385305#ifdef TIME_TEST
    4939   TimeString(tinput, tostd, tif, tstd, tlift, tred, tnw, nstep);
    4940 #endif
    4941 
    4942   return(Res);
    4943 }
    4944 
    4945 /*******************************************************
    4946  * THE PERTURBATION WALK ALGORITHM WITH RANDOM ELEMENT *
    4947  *******************************************************/
    4948 ideal Mprwalk(ideal Go, intvec* curr_weight, intvec* target_weight, int weight_rad, int op_deg, int tp_deg, ring baseRing)
    4949 {
    4950   BITSET save1 = si_opt_1; // save current options
    4951   si_opt_1 &= (~Sy_bit(OPT_REDSB)); // no reduced Groebner basis
    4952   Set_Error(FALSE);
    4953   Overflow_Error = FALSE;
    4954 #ifdef TIME_TEST
    4955   clock_t tinput=0, tostd=0, tif=0, tstd=0, tlift=0, tred=0, tnw=0;
    4956   xtif=0; xtstd=0; xtlift=0; xtred=0; xtnw=0;
    4957   tinput = clock();
    4958   clock_t tim;
    4959 #endif
    4960   int i,nwalk,nV = baseRing->N;
    4961 
    4962   ideal G, Gomega, M, F, Gomega1, Gomega2, M1;
    4963   ring newRing;
    4964   ring XXRing = baseRing;
    4965   intvec* exivlp = Mivlp(nV);
    4966   intvec* orig_target = target_weight;
    4967   intvec* pert_target_vector = target_weight;
    4968   intvec* ivNull = new intvec(nV);
    4969   intvec* tmp_weight = new intvec(nV);
    4970 #ifdef CHECK_IDEAL_MWALK
    4971   poly p;
    4972 #endif
    4973   for(i=0; i<nV; i++)
    4974   {
    4975     (*tmp_weight)[i] = (*curr_weight)[i];
    4976   }
    4977 #ifndef BUCHBERGER_ALG
    4978   intvec* hilb_func;
    4979    // to avoid (1,0,...,0) as the target vector
    4980   intvec* last_omega = new intvec(nV);
    4981   for(i=0 i<nV; i++)
    4982   {
    4983     (*last_omega)[i] = 1;
    4984   }
    4985   (*last_omega)[0] = 10000;
    4986 #endif
    4987   baseRing = currRing;
    4988   newRing = VMrDefault(curr_weight);
    4989   rChangeCurrRing(newRing);
    4990   G = idrMoveR(Go,baseRing,currRing);
    4991 #ifdef TIME_TEST
    4992   to = clock();
    4993 #endif
    4994   G = kStd(G,NULL,testHomog,NULL,NULL,0,0,NULL);
    4995   idSkipZeroes(G);
    4996 #ifdef TIME_TEST
    4997   tostd = tostd + to - clock();
    4998 #endif
    4999 #ifdef CHECK_IDEAL_MWALK
    5000   idString(G,"G");
    5001 #endif
    5002   if(op_deg >1)
    5003   {
    5004     if(MivComp(curr_weight,MivUnit(nV)) == 1) //ring order is "dp"
    5005     {
    5006       curr_weight = MPertVectors(G, MivMatrixOrderdp(nV), op_deg);
    5007     }
    5008     else //ring order is not "dp"
    5009     {
    5010       curr_weight = MPertVectors(G, MivMatrixOrder(curr_weight), op_deg);
    5011     }
    5012   }
    5013   baseRing = currRing;
    5014   if(tp_deg > 1 && tp_deg <= nV)
    5015   {
    5016     pert_target_vector = target_weight;
    5017   }
    5018 #ifdef CHECK_IDEAL_MWALK
    5019   ivString(curr_weight, "new curr_weight");
    5020   ivString(target_weight, "new target_weight");
    5021 #endif
    5022   nwalk = 0;
    5023   while(1)
    5024   {
    5025     nwalk ++;
    5026 #ifdef TIME_TEST
    5027     to = clock();
    5028 #endif
    5029     Gomega = MwalkInitialForm(G, curr_weight); // compute an initial form ideal of <G> w.r.t. "curr_vector"
    5030 #ifdef TIME_TEST
    5031     tif = tif + clock()-to; //time for computing initial form ideal
    5032 #endif
    5033 #ifdef CHECK_IDEAL_MWALK
    5034     idString(Gomega,"Gomega");
    5035 #endif
    5036 #ifndef  BUCHBERGER_ALG
    5037     if(isNolVector(curr_weight) == 0)
    5038     {
    5039       hilb_func = hFirstSeries(Gomega,NULL,NULL,curr_weight,currRing);
    5040     }
    5041     else
    5042     {
    5043       hilb_func = hFirstSeries(Gomega,NULL,NULL,last_omega,currRing);
    5044     }
    5045 #endif
    5046     if(nwalk == 1)
    5047     {
    5048       newRing = VMrDefault(curr_weight); // define a new ring with ordering "(a(curr_weight),lp)
    5049     }
    5050     else
    5051     {
    5052       newRing = VMrRefine(curr_weight,target_weight); //define a new ring with ordering "(a(curr_weight),Wp(target_weight))"
    5053     }
    5054     rChangeCurrRing(newRing);
    5055     Gomega1 = idrMoveR(Gomega, baseRing,currRing);
    5056     idDelete(&Gomega);
    5057     // compute a Groebner basis of <Gomega> w.r.t. "newRing"
    5058 #ifdef TIME_TEST
    5059     to = clock();
    5060 #endif
    5061 #ifndef  BUCHBERGER_ALG
    5062     M=kStd(Gomega1,NULL,isHomog,NULL,hilb_func,0,NULL,curr_weight);
    5063     delete hilb_func;
    5064 #else
    5065     M = kStd(Gomega1,NULL,testHomog,NULL,NULL,0,0,NULL);
    5066 #endif
    5067     idSkipZeroes(M);
    5068 #ifdef TIME_TEST
    5069     tstd = tstd + clock() - to;
    5070 #endif
    5071 #ifdef CHECK_IDEAL_MWALK
    5072     idString(M, "M");
    5073 #endif
    5074     //change the ring to baseRing
    5075     rChangeCurrRing(baseRing);
    5076     M1 =  idrMoveR(M, newRing,currRing);
    5077     idDelete(&M);
    5078     Gomega2 = idrMoveR(Gomega1, newRing,currRing);
    5079     idDelete(&Gomega1);
    5080     to = clock();
    5081     // compute a representation of the generators of submod (M) with respect to those of mod (Gomega), where Gomega is a reduced Groebner basis w.r.t. the current ring
    5082     F = MLifttwoIdeal(Gomega2, M1, G);
    5083     idSkipZeroes(F);
    5084 #ifdef TIME_TEST
    5085     tlift = tlift + clock() - to;
    5086 #endif
    5087 #ifdef CHECK_IDEAL_MWALK
    5088     idString(F,"F");
    5089 #endif
    5090     rChangeCurrRing(newRing); // change the ring to newRing
    5091     G = idrMoveR(F,baseRing,currRing);
    5092     idDelete(&F);
    5093     baseRing = currRing; // set baseRing equal to newRing
    5094 #ifdef CHECK_IDEAL_MWALK
    5095     idString(G,"G");
    5096 #endif
    5097 #ifdef TIME_TEST
    5098     to = clock();
    5099 #endif
    5100     intvec* next_weight = MWalkRandomNextWeight(G, curr_weight, target_weight, weight_rad, op_deg);
    5101 #ifdef TIME_TEST
    5102     tnw = tnw + clock() - to;
    5103 #endif
    5104 #ifdef PRINT_VECTORS
    5105     MivString(curr_weight, target_weight, next_weight);
    5106 #endif
    5107     if(Overflow_Error == TRUE)
    5108     {
    5109       PrintS("\n//**Mprwalk: OVERFLOW: The computed vector does not stay in cone, the result may be wrong.\n");
    5110       delete next_weight;
    5111       break;
    5112     }
    5113 
    5114     if(test_w_in_ConeCC(G,target_weight) == 1 || MivComp(next_weight, ivNull) == 1)
    5115     {
    5116       delete next_weight;
    5117       break;
    5118     }
    5119     //update tmp_weight and curr_weight
    5120     for(i=nV-1; i>=0; i--)
    5121     {
    5122       (*tmp_weight)[i] = (*curr_weight)[i];
    5123       (*curr_weight)[i] = (*next_weight)[i];
    5124     }
    5125     delete next_weight;
    5126   } //end of while-loop
    5127   Print("\n// Mprwalk took %d steps. Ring= %s;\n", nwalk, rString(currRing));
    5128   idSkipZeroes(G);
    5129   si_opt_1 = save1; //set original options, e. g. option(RedSB)
    5130   baseRing = currRing;
    5131   rChangeCurrRing(XXRing);
    5132   ideal Res = idrMoveR(G,baseRing,currRing);
    5133   delete tmp_weight;
    5134   delete ivNull;
    5135   delete exivlp;
    5136 #ifndef BUCHBERGER_ALG
    5137   delete last_omega;
    5138 #endif
    5139 #ifdef TIME_TEST
    5140   TimeString(tinput, tostd, tif, tstd, tlift, tred, tnw, nstep);
    5141 #endif
    5142   return(Res);
    5143 }
    5144 
    5145 poly div_with_remainder( poly f, poly g, ring r)
    5146 {
    5147 return singclap_pdivide ( f, g, r );
     5306  TimeStringFractal(tinput, tostd, tif+xtif, tstd+xtstd,0, tlift+xtlift, tred+xtred,
     5307             tnw+xtnw);
     5308
     5309  Print("\n// pSetm_Error = (%d)", ErrorCheck());
     5310  Print("\n// It took %d steps and Overflow_Error? (%d)\n", nstep,  Overflow_Error);
     5311#endif
     5312  return(Eresult);
    51485313}
    51495314
     
    51675332int Xngleich;
    51685333
    5169 /******************************
    5170  * THE FRACTAL WALK ALGORITHM *
    5171  ******************************/
    5172 ideal Mfwalk(ideal Go, intvec* curr_weight, intvec* orig_target_weight, int pert_deg, ring XXRing)
    5173 {
    5174   Set_Error(FALSE);
    5175   Overflow_Error = FALSE;
    5176 #ifdef TIME_TEST
    5177   clock_t tinput, tostd, tif=0, tstd=0, tlift=0, tred=0, tnw=0;
    5178   xtif=0; xtstd=0; xtlift=0; xtred=0; xtnw=0;
    5179   tinput = clock();
    5180   clock_t tim;
    5181 #endif
    5182   nstep=0;
    5183   ring newRing;
    5184   ring baseRing = currRing;
    5185   rChangeCurrRing(XXRing);
    5186   Go = idrMoveR(Go,baseRing,currRing);
    5187   int i,nwalk,endwalks = 0;
     5334/***********************************************************************
     5335 * Perturb the start weight vector at the top level, i.e. nlev = 1     *
     5336 ***********************************************************************/
     5337static ideal rec_fractal_call(ideal G, int nlev, intvec* omtmp)
     5338{
     5339  Overflow_Error =  FALSE;
     5340  //Print("\n\n// Entering the %d-th recursion:", nlev);
     5341
     5342  int i, nV = currRing->N;
     5343  ring new_ring, testring;
     5344  //ring extoRing;
     5345  ideal Gomega, Gomega1, Gomega2, F, F1, Gresult, Gresult1, G1, Gt;
     5346  int nwalks = 0;
     5347  intvec* Mwlp;
     5348#ifndef BUCHBERGER_ALG
     5349  intvec* hilb_func;
     5350#endif
     5351//  intvec* extXtau;
     5352  intvec* next_vect;
     5353  intvec* omega2 = new intvec(nV);
     5354  intvec* altomega = new intvec(nV);
     5355
     5356  //BOOLEAN isnewtarget = FALSE;
     5357
     5358  // to avoid (1,0,...,0) as the target vector (Hans)
     5359  intvec* last_omega = new intvec(nV);
     5360  for(i=nV-1; i>0; i--)
     5361    (*last_omega)[i] = 1;
     5362  (*last_omega)[0] = 10000;
     5363
     5364  intvec* omega = new intvec(nV);
     5365  for(i=0; i<nV; i++) {
     5366    if(Xsigma->length() == nV)
     5367      (*omega)[i] =  (*Xsigma)[i];
     5368    else
     5369      (*omega)[i] = (*Xsigma)[(nV*(nlev-1))+i];
     5370
     5371    (*omega2)[i] = (*Xtau)[(nlev-1)*nV+i];
     5372  }
     5373
     5374   if(nlev == 1)  Xcall = 1;
     5375   else Xcall = 0;
     5376
     5377  ring oRing = currRing;
     5378
     5379  while(1)
     5380  {
     5381#ifdef FIRST_STEP_FRACTAL
     5382    // perturb the current weight vector only on the top level or
     5383    // after perturbation of the both vectors, nlev = 2 as the top level
     5384    if((nlev == 1 && Xcall == 0) || (nlev == 2 && Xngleich == 1))
     5385      if(islengthpoly2(G) == 1)
     5386      {
     5387        Mwlp = MivWeightOrderlp(omega);
     5388        Xsigma = Mfpertvector(G, Mwlp);
     5389        delete Mwlp;
     5390        Overflow_Error = FALSE;
     5391      }
     5392#endif
     5393    nwalks ++;
     5394    NEXT_VECTOR_FRACTAL:
     5395    to=clock();
     5396    /* determine the next border */
     5397    next_vect = MkInterRedNextWeight(omega,omega2,G);
     5398    xtnw=xtnw+clock()-to;
     5399#ifdef PRINT_VECTORS
     5400    MivString(omega, omega2, next_vect);
     5401#endif
     5402    oRing = currRing;
     5403
     5404    /* We only perturb the current target vector at the recursion  level 1 */
     5405    if(Xngleich == 0 && nlev == 1) //(ngleich == 0) important, e.g. ex2, ex3
     5406      if (MivComp(next_vect, omega2) == 1)
     5407      {
     5408        /* to dispense with taking initial (and lifting/interreducing
     5409           after the call of recursion */
     5410        //Print("\n\n// ** Perturb the both vectors with degree %d with",nlev);
     5411        //idElements(G, "G");
     5412
     5413        Xngleich = 1;
     5414        nlev +=1;
     5415
     5416        if (rParameter(currRing) != NULL)
     5417          DefRingPar(omtmp);
     5418        else
     5419          VMrDefault(omtmp);
     5420
     5421        testring = currRing;
     5422        Gt = idrMoveR(G, oRing,currRing);
     5423
     5424        /* perturb the original target vector w.r.t. the current GB */
     5425        delete Xtau;
     5426        Xtau = NewVectorlp(Gt);
     5427
     5428        rChangeCurrRing(oRing);
     5429        G = idrMoveR(Gt, testring,currRing);
     5430
     5431        /* perturb the current vector w.r.t. the current GB */
     5432        Mwlp = MivWeightOrderlp(omega);
     5433        Xsigma = Mfpertvector(G, Mwlp);
     5434        delete Mwlp;
     5435
     5436        for(i=nV-1; i>=0; i--) {
     5437          (*omega2)[i] = (*Xtau)[nV+i];
     5438          (*omega)[i] = (*Xsigma)[nV+i];
     5439        }
     5440
     5441        delete next_vect;
     5442        to=clock();
     5443
     5444        /* to avoid the value of Overflow_Error that occur in Mfpertvector*/
     5445        Overflow_Error = FALSE;
     5446
     5447        next_vect = MkInterRedNextWeight(omega,omega2,G);
     5448        xtnw=xtnw+clock()-to;
     5449
     5450#ifdef PRINT_VECTORS
     5451        MivString(omega, omega2, next_vect);
     5452#endif
     5453      }
     5454
     5455
     5456    /* check whether the the computed vector is in the correct cone */
     5457    /* If no, the reduced GB of an omega-homogeneous ideal will be
     5458       computed by Buchberger algorithm and stop this recursion step*/
     5459    //if(test_w_in_ConeCC(G, next_vect) != 1) //e.g. Example s7, cyc6
     5460    if(Overflow_Error == TRUE)
     5461    {
     5462      delete next_vect;
     5463      if (rParameter(currRing) != NULL)
     5464      {
     5465        DefRingPar(omtmp);
     5466      }
     5467      else
     5468      {
     5469        VMrDefault(omtmp);
     5470      }
     5471#ifdef TEST_OVERFLOW
     5472      Gt = idrMoveR(G, oRing,currRing);
     5473      Gt = NULL; return(Gt);
     5474#endif
     5475
     5476      //Print("\n\n// apply BB's alg. in ring r = %s;", rString(currRing));
     5477      to=clock();
     5478      Gt = idrMoveR(G, oRing,currRing);
     5479      G1 = MstdCC(Gt);
     5480      xtextra=xtextra+clock()-to;
     5481      Gt = NULL;
     5482
     5483      delete omega2;
     5484      delete altomega;
     5485
     5486      //Print("\n// Leaving the %d-th recursion with %d steps", nlev, nwalks);
     5487      //Print("  ** Overflow_Error? (%d)", Overflow_Error);
     5488      nnflow ++;
     5489
     5490      Overflow_Error = FALSE;
     5491      return (G1);
     5492    }
     5493
     5494
     5495    /* If the perturbed target vector stays in the correct cone,
     5496       return the current GB,
     5497       otherwise, return the computed  GB by the Buchberger-algorithm.
     5498       Then we update the perturbed target vectors w.r.t. this GB. */
     5499
     5500    /* the computed vector is equal to the origin vector, since
     5501       t is not defined */
     5502    if (MivComp(next_vect, XivNull) == 1)
     5503    {
     5504      if (rParameter(currRing) != NULL)
     5505        DefRingPar(omtmp);
     5506      else
     5507        VMrDefault(omtmp);
     5508
     5509      testring = currRing;
     5510      Gt = idrMoveR(G, oRing,currRing);
     5511
     5512      if(test_w_in_ConeCC(Gt, omega2) == 1) {
     5513        delete omega2;
     5514        delete next_vect;
     5515        delete altomega;
     5516        //Print("\n// Leaving the %d-th recursion with %d steps ",nlev, nwalks);
     5517        //Print(" ** Overflow_Error? (%d)", Overflow_Error);
     5518
     5519        return (Gt);
     5520      }
     5521      else
     5522      {
     5523        //ivString(omega2, "tau'");
     5524        //Print("\n//  tau' doesn't stay in the correct cone!!");
     5525
     5526#ifndef  MSTDCC_FRACTAL
     5527        //07.08.03
     5528        //ivString(Xtau, "old Xtau");
     5529        intvec* Xtautmp = Mfpertvector(Gt, MivMatrixOrder(omtmp));
     5530#ifdef TEST_OVERFLOW
     5531      if(Overflow_Error == TRUE)
     5532      Gt = NULL; return(Gt);
     5533#endif
     5534
     5535        if(MivSame(Xtau, Xtautmp) == 1)
     5536        {
     5537          //PrintS("\n// Update vectors are equal to the old vectors!!");
     5538          delete Xtautmp;
     5539          goto FRACTAL_MSTDCC;
     5540        }
     5541
     5542        Xtau = Xtautmp;
     5543        Xtautmp = NULL;
     5544        //ivString(Xtau, "new  Xtau");
     5545
     5546        for(i=nV-1; i>=0; i--)
     5547          (*omega2)[i] = (*Xtau)[(nlev-1)*nV+i];
     5548
     5549        //Print("\n//  ring tau = %s;", rString(currRing));
     5550        rChangeCurrRing(oRing);
     5551        G = idrMoveR(Gt, testring,currRing);
     5552
     5553        goto NEXT_VECTOR_FRACTAL;
     5554#endif
     5555
     5556      FRACTAL_MSTDCC:
     5557        //Print("\n//  apply BB-Alg in ring = %s;", rString(currRing));
     5558        to=clock();
     5559        G = MstdCC(Gt);
     5560        xtextra=xtextra+clock()-to;
     5561
     5562        oRing = currRing;
     5563
     5564        // update the original target vector w.r.t. the current GB
     5565        if(MivSame(Xivinput, Xivlp) == 1)
     5566          if (rParameter(currRing) != NULL)
     5567            DefRingParlp();
     5568          else
     5569            VMrDefaultlp();
     5570        else
     5571          if (rParameter(currRing) != NULL)
     5572            DefRingPar(Xivinput);
     5573          else
     5574            VMrDefault(Xivinput);
     5575
     5576        testring = currRing;
     5577        Gt = idrMoveR(G, oRing,currRing);
     5578
     5579        delete Xtau;
     5580        Xtau = NewVectorlp(Gt);
     5581
     5582        rChangeCurrRing(oRing);
     5583        G = idrMoveR(Gt, testring,currRing);
     5584
     5585        delete omega2;
     5586        delete next_vect;
     5587        delete altomega;
     5588        /*
     5589          Print("\n// Leaving the %d-th recursion with %d steps,", nlev,nwalks);
     5590          Print(" ** Overflow_Error? (%d)", Overflow_Error);
     5591        */
     5592        if(Overflow_Error == TRUE)
     5593          nnflow ++;
     5594
     5595        Overflow_Error = FALSE;
     5596        return(G);
     5597      }
     5598    }
     5599
     5600    for(i=nV-1; i>=0; i--) {
     5601      (*altomega)[i] = (*omega)[i];
     5602      (*omega)[i] = (*next_vect)[i];
     5603    }
     5604    delete next_vect;
     5605
     5606    to=clock();
     5607    /* Take the initial form of <G> w.r.t. omega */
     5608    Gomega = MwalkInitialForm(G, omega);
     5609    xtif=xtif+clock()-to;
     5610
     5611#ifndef  BUCHBERGER_ALG
     5612    if(isNolVector(omega) == 0)
     5613      hilb_func = hFirstSeries(Gomega,NULL,NULL,omega,currRing);
     5614    else
     5615      hilb_func = hFirstSeries(Gomega,NULL,NULL,last_omega,currRing);
     5616#endif // BUCHBERGER_ALG
     5617
     5618    if (rParameter(currRing) != NULL)
     5619      DefRingPar(omega);
     5620    else
     5621      VMrDefault(omega);
     5622
     5623    Gomega1 = idrMoveR(Gomega, oRing,currRing);
     5624
     5625    /* Maximal recursion depth, to compute a red. GB */
     5626    /* Fractal walk with the alternative recursion */
     5627    /* alternative recursion */
     5628    // if(nlev == nV || lengthpoly(Gomega1) == 0)
     5629    if(nlev == Xnlev || lengthpoly(Gomega1) == 0)
     5630      //if(nlev == nV) // blind recursion
     5631    {
     5632      /*
     5633      if(Xnlev != nV)
     5634      {
     5635        Print("\n// ** Xnlev = %d", Xnlev);
     5636        ivString(Xtau, "Xtau");
     5637      }
     5638      */
     5639      to=clock();
     5640#ifdef  BUCHBERGER_ALG
     5641      Gresult = MstdhomCC(Gomega1);
     5642#else
     5643      Gresult =kStd(Gomega1,NULL,isHomog,NULL,hilb_func,0,NULL,omega);
     5644      delete hilb_func;
     5645#endif // BUCHBERGER_ALG
     5646      xtstd=xtstd+clock()-to;
     5647    }
     5648    else {
     5649      rChangeCurrRing(oRing);
     5650      Gomega1 = idrMoveR(Gomega1, oRing,currRing);
     5651      Gresult = rec_fractal_call(idCopy(Gomega1),nlev+1,omega);
     5652    }
     5653
     5654    //convert a Groebner basis from a ring to another ring,
     5655    new_ring = currRing;
     5656
     5657    rChangeCurrRing(oRing);
     5658    Gresult1 = idrMoveR(Gresult, new_ring,currRing);
     5659    Gomega2 = idrMoveR(Gomega1, new_ring,currRing);
     5660
     5661    to=clock();
     5662    /* Lifting process */
     5663    F = MLifttwoIdeal(Gomega2, Gresult1, G);
     5664    xtlift=xtlift+clock()-to;
     5665    idDelete(&Gresult1);
     5666    idDelete(&Gomega2);
     5667    idDelete(&G);
     5668
     5669    rChangeCurrRing(new_ring);
     5670    F1 = idrMoveR(F, oRing,currRing);
     5671
     5672    to=clock();
     5673    /* Interreduce G */
     5674    G = kInterRedCC(F1, NULL);
     5675    xtred=xtred+clock()-to;
     5676    idDelete(&F1);
     5677  }
     5678}
     5679
     5680/************************************************************************
     5681 * Perturb the start weight vector at the top level with random element *
     5682 ************************************************************************/
     5683static ideal rec_r_fractal_call(ideal G, int nlev, intvec* omtmp, int weight_rad)
     5684{
     5685  Overflow_Error =  FALSE;
     5686  //Print("\n\n// Entering the %d-th recursion:", nlev);
     5687
     5688  int i,k,weight_norm;
    51885689  int nV = currRing->N;
    5189   ideal Gomega, M, F, Gomega1, Gomega2, M1, F1;
    5190   intvec* ivNull = new intvec(nV);
    5191   intvec* next_weight = new intvec(nV);
    5192   intvec* tmp_weight = new intvec(nV);
    5193 #ifdef TIME_TEST
    5194   to = clock();
    5195 #endif
    5196   ideal G = MstdCC(Go);
    5197 #ifdef TIME_TEST
    5198   tostd = clock()-to;
    5199 #endif
    5200 #ifdef CHECK_IDEAL_MWALK
    5201     idString(Go,"Go");
    5202 #endif
    5203   intvec* target_weight = MPertVectors(G, MivMatrixOrder(orig_target_weight), pert_deg);
    5204   for(i=0; i<nV; i++)
    5205   {
    5206     (*tmp_weight)[i] = (*target_weight)[i];
    5207   }
    5208   nwalk = 0;
     5690  ring new_ring, testring;
     5691  //ring extoRing;
     5692  ideal Gomega, Gomega1, Gomega2, F, F1, Gresult, Gresult1, G1, Gt;
     5693  int nwalks = 0;
     5694  intvec* Mwlp;
     5695#ifndef BUCHBERGER_ALG
     5696  intvec* hilb_func;
     5697#endif
     5698  //intvec* extXtau;
     5699  intvec* next_vect;
     5700  intvec* omega2 = new intvec(nV);
     5701  intvec* altomega = new intvec(nV);
     5702
     5703  //BOOLEAN isnewtarget = FALSE;
     5704
     5705  // to avoid (1,0,...,0) as the target vector (Hans)
     5706  intvec* last_omega = new intvec(nV);
     5707  for(i=nV-1; i>0; i--)
     5708    (*last_omega)[i] = 1;
     5709  (*last_omega)[0] = 10000;
     5710
     5711  intvec* omega = new intvec(nV);
     5712  for(i=0; i<nV; i++) {
     5713    if(Xsigma->length() == nV)
     5714      (*omega)[i] =  (*Xsigma)[i];
     5715    else
     5716      (*omega)[i] = (*Xsigma)[(nV*(nlev-1))+i];
     5717
     5718    (*omega2)[i] = (*Xtau)[(nlev-1)*nV+i];
     5719  }
     5720
     5721   if(nlev == 1)  Xcall = 1;
     5722   else Xcall = 0;
     5723
     5724  ring oRing = currRing;
     5725
    52095726  while(1)
    52105727  {
    5211     nwalk++;
    5212     baseRing = currRing;
    5213     Gomega = MwalkInitialForm(G, curr_weight);
    5214 #ifdef CHECK_IDEAL_MWALK
    5215     idString(Gomega,"Gomega");
    5216 #endif
    5217     next_weight = MwalkNextWeightCC(curr_weight,target_weight,G);
    5218     MivString(curr_weight, target_weight, tmp_weight);
    5219     if(MivSame(curr_weight,target_weight) ==1)
    5220     {
    5221       break;
    5222     }
    5223     if( MivSame(tmp_weight, curr_weight) == 1 )
    5224     {
    5225       if(test_w_in_ConeCC(G, target_weight) == 1)
    5226       {
    5227         break;
     5728#ifdef FIRST_STEP_FRACTAL
     5729    // perturb the current weight vector only on the top level or
     5730    // after perturbation of the both vectors, nlev = 2 as the top level
     5731    if((nlev == 1 && Xcall == 0) || (nlev == 2 && Xngleich == 1))
     5732      if(islengthpoly2(G) == 1)
     5733      {
     5734        Mwlp = MivWeightOrderlp(omega);
     5735        Xsigma = Mfpertvector(G, Mwlp);
     5736        delete Mwlp;
     5737        Overflow_Error = FALSE;
     5738      }
     5739#endif
     5740    nwalks ++;
     5741    NEXT_VECTOR_FRACTAL:
     5742    to=clock();
     5743    /* determine the next border */
     5744    next_vect = MkInterRedNextWeight(omega,omega2,G);
     5745    xtnw=xtnw+clock()-to;
     5746#ifdef PRINT_VECTORS
     5747    MivString(omega, omega2, next_vect);
     5748#endif
     5749    oRing = currRing;
     5750
     5751    /* We only perturb the current target vector at the recursion  level 1 */
     5752    if(Xngleich == 0 && nlev == 1) //(ngleich == 0) important, e.g. ex2, ex3
     5753    {
     5754      if (MivComp(next_vect, omega2) == 1)
     5755      {
     5756        /* to dispense with taking initial (and lifting/interreducing
     5757           after the call of recursion */
     5758        //Print("\n\n// ** Perturb the both vectors with degree %d with",nlev);
     5759        //idElements(G, "G");
     5760
     5761        Xngleich = 1;
     5762        nlev +=1;
     5763
     5764        if (rParameter(currRing) != NULL)
     5765          DefRingPar(omtmp);
     5766        else
     5767          VMrDefault(omtmp);
     5768
     5769        testring = currRing;
     5770        Gt = idrMoveR(G, oRing,currRing);
     5771
     5772        /* perturb the original target vector w.r.t. the current GB */
     5773        delete Xtau;
     5774        Xtau = NewVectorlp(Gt);
     5775
     5776        rChangeCurrRing(oRing);
     5777        G = idrMoveR(Gt, testring,currRing);
     5778
     5779        /* perturb the current vector w.r.t. the current GB */
     5780        Mwlp = MivWeightOrderlp(omega);
     5781        Xsigma = Mfpertvector(G, Mwlp);
     5782        delete Mwlp;
     5783
     5784        for(i=nV-1; i>=0; i--) {
     5785          (*omega2)[i] = (*Xtau)[nV+i];
     5786          (*omega)[i] = (*Xsigma)[nV+i];
     5787        }
     5788
     5789        delete next_vect;
     5790        to=clock();
     5791
     5792        /* to avoid the value of Overflow_Error that occur in Mfpertvector*/
     5793        Overflow_Error = FALSE;
     5794
     5795        next_vect = MkInterRedNextWeight(omega,omega2,G);
     5796        xtnw=xtnw+clock()-to;
     5797
     5798#ifdef PRINT_VECTORS
     5799        MivString(omega, omega2, next_vect);
     5800#endif
    52285801      }
    52295802      else
    52305803      {
    5231         pert_deg++;
    5232         target_weight = MPertVectors(G, MivMatrixOrder(orig_target_weight), pert_deg);
    5233         next_weight = MwalkNextWeightCC(curr_weight,target_weight,G);
    5234       }
    5235     }
    5236     for(i=0; i<nV; i++)
    5237     {
    5238       (*tmp_weight)[i] = (*curr_weight)[i];
    5239       (*curr_weight)[i] = (*next_weight)[i];
    5240     }
    5241     MivString(curr_weight, target_weight, tmp_weight);
    5242     newRing = VMrRefine(tmp_weight,target_weight); //define a new ring with ordering "(a(curr_weight),Wp(target_weight))"
    5243     rChangeCurrRing(newRing); // change the ring to newRing
    5244     Gomega1 =  idrMoveR(Gomega, baseRing,currRing);
    5245     idDelete(&Gomega);
    5246     if(pert_deg >= nV)
    5247     {
    5248       M = kStd(Gomega1,NULL,testHomog,NULL,NULL,0,0,NULL);
    5249     }
     5804        // compute a perturbed next weight vector "next_vect1"
     5805        intvec* next_vect11 = MPertVectors(G, MivMatrixOrder(omega), nlev);
     5806        intvec* next_vect1 = MkInterRedNextWeight(next_vect11, omega2, G);
     5807        // Print("\n//  size of next_vect  = %d", sizeof((*next_vect)));
     5808        // Print("\n//  size of next_vect1  = %d", sizeof((*next_vect1)));
     5809        // Print("\n//  size of next_vect11  = %d", sizeof((*next_vect11)));
     5810        delete next_vect11;
     5811
     5812        // compare next_vect and next_vect1
     5813        ideal G_test = MwalkInitialForm(G, next_vect);
     5814        ideal G_test1 = MwalkInitialForm(G, next_vect1);
     5815        // Print("\n// G_test, G_test 1 erzeugt");
     5816        if(IDELEMS(G_test1) <= IDELEMS(G_test))
     5817          {
     5818          next_vect = ivCopy(next_vect1);
     5819          }
     5820        delete next_vect1;
     5821        // compute a random next weight vector "next_vect2"
     5822        intvec* next_vect22 = ivCopy(omega2);
     5823        // Print("\n//  size of next_vect22  = %d", sizeof((*next_vect22)));
     5824        k = 0;
     5825        while(test_w_in_ConeCC(G, next_vect22) == 0)
     5826        {
     5827          k = k + 1;
     5828          if(k>10)
     5829          {
     5830            break;
     5831          }
     5832          weight_norm = 0;
     5833          while(weight_norm == 0)
     5834          {
     5835            for(i=nV-1; i>=0; i--)
     5836            {
     5837              (*next_vect22)[i] = rand() % 60000 - 30000;
     5838              weight_norm = weight_norm + (*next_vect22)[i]*(*next_vect22)[i];
     5839            }
     5840            weight_norm = 1 + floor(sqrt(weight_norm));
     5841          }
     5842          for(i=nV-1; i>=0; i--)
     5843          {
     5844            if((*next_vect22)[i] < 0)
     5845            {
     5846              (*next_vect22)[i] = 1 + (*omega)[i] + floor(weight_rad*(*next_vect22)[i]/weight_norm);
     5847            }
     5848            else
     5849            {
     5850              (*next_vect22)[i] = (*omega)[i] + floor(weight_rad*(*next_vect22)[i]/weight_norm);
     5851            }
     5852          }
     5853        }
     5854        if(test_w_in_ConeCC(G, next_vect22) == 1)
     5855        {
     5856          //compare next_weight and next_vect2
     5857          intvec* next_vect2 = MkInterRedNextWeight(next_vect22, omega2, G);
     5858          // Print("\n//  size of next_vect2  = %d", sizeof((*next_vect2)));
     5859          ideal G_test2 = MwalkInitialForm(G, next_vect2);
     5860          if(IDELEMS(G_test2) <= IDELEMS(G_test))
     5861          {
     5862            if(IDELEMS(G_test2) <= IDELEMS(G_test1))
     5863            {
     5864              next_vect = ivCopy(next_vect2);
     5865            }
     5866          }
     5867          idDelete(&G_test2);
     5868          delete next_vect2;
     5869        }
     5870        delete next_vect22;
     5871        idDelete(&G_test);
     5872        idDelete(&G_test1);
     5873#ifdef PRINT_VECTORS
     5874        MivString(omega, omega2, next_vect);
     5875#endif
     5876      }
     5877    }
     5878
     5879
     5880    /* check whether the the computed vector is in the correct cone */
     5881    /* If no, the reduced GB of an omega-homogeneous ideal will be
     5882       computed by Buchberger algorithm and stop this recursion step*/
     5883    //if(test_w_in_ConeCC(G, next_vect) != 1) //e.g. Example s7, cyc6
     5884    if(Overflow_Error == TRUE)
     5885    {
     5886      delete next_vect;
     5887
     5888    //OVERFLOW_IN_NEXT_VECTOR:
     5889      if (rParameter(currRing) != NULL)
     5890        DefRingPar(omtmp);
     5891      else
     5892        VMrDefault(omtmp);
     5893
     5894#ifdef TEST_OVERFLOW
     5895      Gt = idrMoveR(G, oRing,currRing);
     5896      Gt = NULL; return(Gt);
     5897#endif
     5898
     5899      //Print("\n\n// apply BB's alg. in ring r = %s;", rString(currRing));
     5900      to=clock();
     5901      Gt = idrMoveR(G, oRing,currRing);
     5902      G1 = MstdCC(Gt);
     5903      xtextra=xtextra+clock()-to;
     5904      Gt = NULL;
     5905
     5906      delete omega2;
     5907      delete altomega;
     5908
     5909      //Print("\n// Leaving the %d-th recursion with %d steps", nlev, nwalks);
     5910      //Print("  ** Overflow_Error? (%d)", Overflow_Error);
     5911      nnflow ++;
     5912
     5913      Overflow_Error = FALSE;
     5914      return (G1);
     5915    }
     5916
     5917
     5918    /* If the perturbed target vector stays in the correct cone,
     5919       return the current GB,
     5920       otherwise, return the computed  GB by the Buchberger-algorithm.
     5921       Then we update the perturbed target vectors w.r.t. this GB. */
     5922
     5923    /* the computed vector is equal to the origin vector, since
     5924       t is not defined */
     5925    if (MivComp(next_vect, XivNull) == 1)
     5926    {
     5927      if (rParameter(currRing) != NULL)
     5928        DefRingPar(omtmp);
     5929      else
     5930        VMrDefault(omtmp);
     5931
     5932      testring = currRing;
     5933      Gt = idrMoveR(G, oRing,currRing);
     5934
     5935      if(test_w_in_ConeCC(Gt, omega2) == 1) {
     5936        delete omega2;
     5937        delete next_vect;
     5938        delete altomega;
     5939        //Print("\n// Leaving the %d-th recursion with %d steps ",nlev, nwalks);
     5940        //Print(" ** Overflow_Error? (%d)", Overflow_Error);
     5941
     5942        return (Gt);
     5943      }
     5944      else
     5945      {
     5946        //ivString(omega2, "tau'");
     5947        //Print("\n//  tau' doesn't stay in the correct cone!!");
     5948
     5949#ifndef  MSTDCC_FRACTAL
     5950        //07.08.03
     5951        //ivString(Xtau, "old Xtau");
     5952        intvec* Xtautmp = Mfpertvector(Gt, MivMatrixOrder(omtmp));
     5953#ifdef TEST_OVERFLOW
     5954      if(Overflow_Error == TRUE)
     5955      Gt = NULL; return(Gt);
     5956#endif
     5957
     5958        if(MivSame(Xtau, Xtautmp) == 1)
     5959        {
     5960          //PrintS("\n// Update vectors are equal to the old vectors!!");
     5961          delete Xtautmp;
     5962          goto FRACTAL_MSTDCC;
     5963        }
     5964
     5965        Xtau = Xtautmp;
     5966        Xtautmp = NULL;
     5967        //ivString(Xtau, "new  Xtau");
     5968
     5969        for(i=nV-1; i>=0; i--)
     5970          (*omega2)[i] = (*Xtau)[(nlev-1)*nV+i];
     5971
     5972        //Print("\n//  ring tau = %s;", rString(currRing));
     5973        rChangeCurrRing(oRing);
     5974        G = idrMoveR(Gt, testring,currRing);
     5975
     5976        goto NEXT_VECTOR_FRACTAL;
     5977#endif
     5978
     5979      FRACTAL_MSTDCC:
     5980        //Print("\n//  apply BB-Alg in ring = %s;", rString(currRing));
     5981        to=clock();
     5982        G = MstdCC(Gt);
     5983        xtextra=xtextra+clock()-to;
     5984
     5985        oRing = currRing;
     5986
     5987        // update the original target vector w.r.t. the current GB
     5988        if(MivSame(Xivinput, Xivlp) == 1)
     5989          if (rParameter(currRing) != NULL)
     5990            DefRingParlp();
     5991          else
     5992            VMrDefaultlp();
     5993        else
     5994          if (rParameter(currRing) != NULL)
     5995            DefRingPar(Xivinput);
     5996          else
     5997            VMrDefault(Xivinput);
     5998
     5999        testring = currRing;
     6000        Gt = idrMoveR(G, oRing,currRing);
     6001
     6002        delete Xtau;
     6003        Xtau = NewVectorlp(Gt);
     6004
     6005        rChangeCurrRing(oRing);
     6006        G = idrMoveR(Gt, testring,currRing);
     6007
     6008        delete omega2;
     6009        delete next_vect;
     6010        delete altomega;
     6011        /*
     6012          Print("\n// Leaving the %d-th recursion with %d steps,", nlev,nwalks);
     6013          Print(" ** Overflow_Error? (%d)", Overflow_Error);
     6014        */
     6015        if(Overflow_Error == TRUE)
     6016          nnflow ++;
     6017
     6018        Overflow_Error = FALSE;
     6019        return(G);
     6020      }
     6021    }
     6022
     6023    for(i=nV-1; i>=0; i--) {
     6024      (*altomega)[i] = (*omega)[i];
     6025      (*omega)[i] = (*next_vect)[i];
     6026    }
     6027    delete next_vect;
     6028
     6029    to=clock();
     6030    /* Take the initial form of <G> w.r.t. omega */
     6031    Gomega = MwalkInitialForm(G, omega);
     6032    xtif=xtif+clock()-to;
     6033
     6034#ifndef  BUCHBERGER_ALG
     6035    if(isNolVector(omega) == 0)
     6036      hilb_func = hFirstSeries(Gomega,NULL,NULL,omega,currRing);
    52506037    else
    5251     {
    5252       M = idCopy(Mwalk(idCopy(Gomega1), tmp_weight, curr_weight, currRing));
    5253     }
    5254     idSkipZeroes(M);
    5255 #ifdef CHECK_IDEAL_MWALK
    5256     idString(M,"M");
    5257 #endif
    5258     rChangeCurrRing(baseRing); //change ring to baseRing
    5259     M1 =  idrMoveR(M, newRing,currRing);
    5260     Gomega2 = idrMoveR(Gomega1, newRing,currRing);
    5261     idDelete(&M);
    5262     idDelete(&Gomega1);
    5263     F = MLifttwoIdeal(Gomega2, M1, G);
     6038      hilb_func = hFirstSeries(Gomega,NULL,NULL,last_omega,currRing);
     6039#endif // BUCHBERGER_ALG
     6040
     6041    if (rParameter(currRing) != NULL)
     6042      DefRingPar(omega);
     6043    else
     6044      VMrDefault(omega);
     6045
     6046    Gomega1 = idrMoveR(Gomega, oRing,currRing);
     6047
     6048    /* Maximal recursion depth, to compute a red. GB */
     6049    /* Fractal walk with the alternative recursion */
     6050    /* alternative recursion */
     6051    // if(nlev == nV || lengthpoly(Gomega1) == 0)
     6052    if(nlev == Xnlev || lengthpoly(Gomega1) == 0)
     6053      //if(nlev == nV) // blind recursion
     6054    {
     6055      /*
     6056      if(Xnlev != nV)
     6057      {
     6058        Print("\n// ** Xnlev = %d", Xnlev);
     6059        ivString(Xtau, "Xtau");
     6060      }
     6061      */
     6062      to=clock();
     6063#ifdef  BUCHBERGER_ALG
     6064      Gresult = MstdhomCC(Gomega1);
     6065#else
     6066      Gresult =kStd(Gomega1,NULL,isHomog,NULL,hilb_func,0,NULL,omega);
     6067      delete hilb_func;
     6068#endif // BUCHBERGER_ALG
     6069      xtstd=xtstd+clock()-to;
     6070    }
     6071    else {
     6072      rChangeCurrRing(oRing);
     6073      Gomega1 = idrMoveR(Gomega1, oRing,currRing);
     6074      Gresult = rec_fractal_call(idCopy(Gomega1),nlev+1,omega);
     6075    }
     6076
     6077    //convert a Groebner basis from a ring to another ring,
     6078    new_ring = currRing;
     6079
     6080    rChangeCurrRing(oRing);
     6081    Gresult1 = idrMoveR(Gresult, new_ring,currRing);
     6082    Gomega2 = idrMoveR(Gomega1, new_ring,currRing);
     6083
     6084    to=clock();
     6085    /* Lifting process */
     6086    F = MLifttwoIdeal(Gomega2, Gresult1, G);
     6087    xtlift=xtlift+clock()-to;
     6088    idDelete(&Gresult1);
    52646089    idDelete(&Gomega2);
    5265     idSkipZeroes(F);
    5266 #ifdef CHECK_IDEAL_MWALK
    5267     idString(F,"F");
    5268 #endif
    5269     rChangeCurrRing(newRing); // change the ring to newRing
    5270     G = idrMoveR(F,baseRing,currRing);
    5271     idDelete(&F);
    5272     baseRing = currRing;
    5273     if(MivComp(next_weight, ivNull) == 1)
    5274     {
    5275       delete next_weight;
    5276       break;
    5277     }
    5278   }
    5279   delete next_weight;
    5280   if(test_w_in_ConeCC(G, orig_target_weight) != 1)
    5281   {
    5282     newRing = VMrDefault(orig_target_weight);  //define a new ring with ordering "(a(target_weight),lp)
    5283     rChangeCurrRing(newRing);
    5284     G = idrMoveR(G,baseRing,currRing);
    5285     M = idCopy(Mwalk(idCopy(G), tmp_weight, orig_target_weight, currRing));
    5286     idSkipZeroes(G);
    5287     baseRing = currRing;
    5288   }
    5289   ENDE:
    5290   rChangeCurrRing(XXRing);
    5291   G = idrMoveR(G,baseRing,currRing);
    5292   delete tmp_weight;
    5293   delete ivNull;
    5294 #ifdef TIME_TEST
    5295   TimeString(tinput, tostd, tif, tstd, tlift, tred, tnw, nstep);
    5296 #endif
    5297   return(G);
    5298 }
    5299 
    5300 /***********************************************
    5301 THE FRACTAL WALK ALGORITHM WITH RANDOM ELEMENT *
    5302 ************************************************/
    5303 ideal Mfrwalk(ideal Go, intvec* curr_weight, intvec* orig_target_weight, int pert_deg, int weight_rad, ring XXRing)
     6090    idDelete(&G);
     6091
     6092    rChangeCurrRing(new_ring);
     6093    F1 = idrMoveR(F, oRing,currRing);
     6094
     6095    to=clock();
     6096    /* Interreduce G */
     6097    G = kInterRedCC(F1, NULL);
     6098    xtred=xtred+clock()-to;
     6099    idDelete(&F1);
     6100  }
     6101}
     6102
     6103
     6104
     6105/*******************************************************************************
     6106 * The implementation of the fractal walk algorithm                            *
     6107 *                                                                             *
     6108 * The main procedur Mfwalk calls the recursive Subroutine                     *
     6109 * rec_fractal_call to compute the wanted Grï¿œbner basis.                       *
     6110 * At the main procedur we compute the reduced Grï¿œbner basis w.r.t. a "fast"   *
     6111 * order, e.g. "dp" and a sequence of weight vectors which are row vectors     *
     6112 * of a matrix. This matrix defines the given monomial order, e.g. "lp"        *
     6113 *******************************************************************************/
     6114ideal Mfwalk(ideal G, intvec* ivstart, intvec* ivtarget)
    53046115{
    53056116  Set_Error(FALSE);
    53066117  Overflow_Error = FALSE;
    53076118  //Print("// pSetm_Error = (%d)", ErrorCheck());
     6119  //Print("\n// ring ro = %s;", rString(currRing));
     6120
     6121  nnflow = 0;
     6122  Xngleich = 0;
     6123  Xcall = 0;
     6124  xtif=0; xtstd=0; xtlift=0; xtred=0; xtnw=0; xtextra=0;
     6125  xftinput = clock();
     6126
     6127  ring  oldRing = currRing;
     6128  int i, nV = currRing->N;
     6129  XivNull = new intvec(nV);
     6130  Xivinput = ivtarget;
     6131  ngleich = 0;
     6132  to=clock();
     6133  ideal I = MstdCC(G);
     6134  G = NULL;
     6135  xftostd=clock()-to;
     6136  Xsigma = ivstart;
     6137
     6138  Xnlev=nV;
     6139
     6140#ifdef FIRST_STEP_FRACTAL
     6141  ideal Gw = MwalkInitialForm(I, ivstart);
     6142  for(i=IDELEMS(Gw)-1; i>=0; i--)
     6143  {
     6144    if((Gw->m[i]!=NULL) // len >=0
     6145    && (Gw->m[i]->next!=NULL) // len >=1
     6146    && (Gw->m[i]->next->next!=NULL)) // len >=2
     6147    {
     6148      intvec* iv_dp = MivUnit(nV); // define (1,1,...,1)
     6149      intvec* Mdp;
     6150
     6151      if(MivSame(ivstart, iv_dp) != 1)
     6152        Mdp = MivWeightOrderdp(ivstart);
     6153      else
     6154        Mdp = MivMatrixOrderdp(nV);
     6155
     6156      Xsigma = Mfpertvector(I, Mdp);
     6157      Overflow_Error = FALSE;
     6158
     6159      delete Mdp;
     6160      delete iv_dp;
     6161      break;
     6162    }
     6163  }
     6164  idDelete(&Gw);
     6165#endif
     6166
     6167  ideal I1;
     6168  intvec* Mlp;
     6169  Xivlp = Mivlp(nV);
     6170
     6171  if(MivComp(ivtarget, Xivlp)  != 1)
     6172  {
     6173    if (rParameter(currRing) != NULL)
     6174      DefRingPar(ivtarget);
     6175    else
     6176      VMrDefault(ivtarget);
     6177
     6178    I1 = idrMoveR(I, oldRing,currRing);
     6179    Mlp = MivWeightOrderlp(ivtarget);
     6180    Xtau = Mfpertvector(I1, Mlp);
     6181  }
     6182  else
     6183  {
     6184    if (rParameter(currRing) != NULL)
     6185      DefRingParlp();
     6186    else
     6187      VMrDefaultlp();
     6188
     6189    I1 = idrMoveR(I, oldRing,currRing);
     6190    Mlp =  MivMatrixOrderlp(nV);
     6191    Xtau = Mfpertvector(I1, Mlp);
     6192  }
     6193  delete Mlp;
     6194  Overflow_Error = FALSE;
     6195
     6196  //ivString(Xsigma, "Xsigma");
     6197  //ivString(Xtau, "Xtau");
     6198
     6199  id_Delete(&I, oldRing);
     6200  ring tRing = currRing;
     6201
     6202  if (rParameter(currRing) != NULL)
     6203    DefRingPar(ivstart);
     6204  else
     6205    VMrDefault(ivstart);
     6206
     6207  I = idrMoveR(I1,tRing,currRing);
     6208  to=clock();
     6209  ideal J = MstdCC(I);
     6210  idDelete(&I);
     6211  xftostd=xftostd+clock()-to;
     6212
     6213  ideal resF;
     6214  ring helpRing = currRing;
     6215
     6216  J = rec_fractal_call(J, 1, ivtarget);
     6217
     6218  rChangeCurrRing(oldRing);
     6219  resF = idrMoveR(J, helpRing,currRing);
     6220  idSkipZeroes(resF);
     6221
     6222  delete Xivlp;
     6223  delete Xsigma;
     6224  delete Xtau;
     6225  delete XivNull;
     6226
    53086227#ifdef TIME_TEST
    5309   clock_t tinput, tostd, tif=0, tstd=0, tlift=0, tred=0, tnw=0;
    5310   xtif=0; xtstd=0; xtlift=0; xtred=0; xtnw=0;
    5311   tinput = clock();
    5312   clock_t tim;
    5313 #endif
    5314   nstep=0;
    5315   ring newRing;
    5316   ring baseRing = currRing;
    5317   rChangeCurrRing(XXRing);
    5318   Go = idrMoveR(Go,baseRing,currRing);
    5319   int i,nwalk,endwalks = 0;
    5320   int nV = currRing->N;
    5321   ideal Gomega, M, F, Gomega1, Gomega2, M1, F1;
    5322   intvec* ivNull = new intvec(nV);