Changeset e92c6a in git


Ignore:
Timestamp:
Apr 19, 2013, 11:52:27 AM (11 years ago)
Author:
Stephan Oberfranz <oberfran@…>
Branches:
(u'spielwiese', '5b153614cbc72bfa198d75b1e9e33dab2645d9fe')
Children:
4d998e87e374f7d2eeb8255db038e35cecbd06f6
Parents:
8659a91de7a2bca7e3533e0b2270e4ecd33c1809
git-author:
Stephan Oberfranz <oberfran@pipo.mathematik.uni-kl.de>2013-04-19 11:52:27+02:00
git-committer:
Hans Schoenemann <hannes@mathematik.uni-kl.de>2013-04-19 11:53:54+02:00
Message:
walk stuff
Location:
Singular
Files:
2 added
3 edited

Legend:

Unmodified
Added
Removed
  • Singular/extra.cc

    r8659a9 re92c6a  
    20212021      }
    20222022      else
     2023      if (strcmp(sys_cmd, "Mrwalk") == 0)
     2024      { // Random Walk
     2025        if (h == NULL || h->Typ() != IDEAL_CMD ||
     2026            h->next == NULL || h->next->Typ() != INTVEC_CMD ||
     2027            h->next->next == NULL || h->next->next->Typ() != INTVEC_CMD ||
     2028            h->next->next->next == NULL || h->next->next->next->Typ() != INT_CMD ||
     2029            h->next->next->next->next == NULL || h->next->next->next->next->Typ() != INT_CMD)
     2030        {
     2031          Werror("system(\"Mrwalk\", ideal, intvec, intvec, int, int) expected");
     2032          return TRUE;
     2033        }
     2034
     2035        if (((intvec*) h->next->Data())->length() != currRing->N &&
     2036            ((intvec*) h->next->next->Data())->length() != currRing->N )
     2037        {
     2038          Werror("system(\"Mrwalk\" ...) intvecs not of length %d\n",
     2039                 currRing->N);
     2040          return TRUE;
     2041        }
     2042        ideal arg1 = (ideal) h->Data();
     2043        intvec* arg2 = (intvec*) h->next->Data();
     2044        intvec* arg3 =  (intvec*) h->next->next->Data();
     2045        int arg4 = (int)(long) h->next->next->next->Data();
     2046        int arg5 = (int)(long) h->next->next->next->next->Data();
     2047
     2048
     2049        ideal result = (ideal) Mrwalk(arg1, arg2, arg3, arg4, arg5);
     2050
     2051        res->rtyp = IDEAL_CMD;
     2052        res->data =  result;
     2053
     2054        return FALSE;
     2055      }
     2056      else
    20232057      if (strcmp(sys_cmd, "MAltwalk1") == 0)
    20242058      {
     
    21192153      }
    21202154      else
     2155      if (strcmp(sys_cmd, "Mfrwalk") == 0)
     2156      {
     2157        if (h == NULL || h->Typ() != IDEAL_CMD ||
     2158            h->next == NULL || h->next->Typ() != INTVEC_CMD ||
     2159            h->next->next == NULL || h->next->next->Typ() != INTVEC_CMD ||
     2160            h->next->next->next == NULL || h->next->next->next->Typ() != INT_CMD)
     2161        {
     2162          Werror("system(\"Mfrwalk\", ideal, intvec, intvec, int) expected");
     2163          return TRUE;
     2164        }
     2165
     2166        if (((intvec*) h->next->Data())->length() != currRing->N &&
     2167            ((intvec*) h->next->next->Data())->length() != currRing->N )
     2168        {
     2169          Werror("system(\"Mfrwalk\" ...) intvecs not of length %d\n",
     2170                 currRing->N);
     2171          return TRUE;
     2172        }
     2173        ideal arg1 = (ideal) h->Data();
     2174        intvec* arg2 = (intvec*) h->next->Data();
     2175        intvec* arg3 = (intvec*) h->next->next->Data();
     2176        int arg4 = (int)(long) h->next->next->next->Data();
     2177       
     2178        ideal result = (ideal) Mfrwalk(arg1, arg2, arg3, arg4);
     2179
     2180        res->rtyp = IDEAL_CMD;
     2181        res->data =  result;
     2182
     2183        return FALSE;
     2184      }
     2185      else
    21212186
    21222187  #ifdef TRAN_Orig
     
    22132278      }
    22142279      else
     2280      if (strcmp(sys_cmd, "TranMrImprovwalk") == 0)
     2281      {
     2282        if (h == NULL || h->Typ() != IDEAL_CMD ||
     2283            h->next == NULL || h->next->Typ() != INTVEC_CMD ||
     2284            h->next->next == NULL || h->next->next->Typ() != INTVEC_CMD ||
     2285            h->next->next->next == NULL || h->next->next->next->Typ() != INT_CMD ||
     2286            h->next->next->next == NULL || h->next->next->next->next->Typ() != INT_CMD ||
     2287            h->next->next->next == NULL || h->next->next->next->next->next->Typ() != INT_CMD)
     2288        {
     2289          Werror("system(\"TranMrImprovwalk\", ideal, intvec, intvec) expected");
     2290          return TRUE;
     2291        }
     2292
     2293        if (((intvec*) h->next->Data())->length() != currRing->N &&
     2294            ((intvec*) h->next->next->Data())->length() != currRing->N )
     2295        {
     2296          Werror("system(\"TranMrImprovwalk\" ...) intvecs not of length %d\n", currRing->N);
     2297          return TRUE;
     2298        }
     2299        ideal arg1 = (ideal) h->Data();
     2300        intvec* arg2 = (intvec*) h->next->Data();
     2301        intvec* arg3 = (intvec*) h->next->next->Data();
     2302        int arg4 = (int)(long) h->next->next->next->Data();
     2303        int arg5 = (int)(long) h->next->next->next->next->Data();
     2304        int arg6 = (int)(long) h->next->next->next->next->next->Data();
     2305
     2306        ideal result = (ideal) TranMrImprovwalk(arg1, arg2, arg3, arg4, arg5, arg6);
     2307
     2308        res->rtyp = IDEAL_CMD;
     2309        res->data =  result;
     2310
     2311        return FALSE;
     2312      }
     2313      else
     2314
    22152315  #endif
    22162316  /*================= Extended system call ========================*/
  • Singular/walk.cc

    r8659a9 re92c6a  
    22*  Computer Algebra System SINGULAR      *
    33*****************************************/
     4/* $Id$ */
    45/*
    56* ABSTRACT: Implementation of the Groebner walk
     
    3031
    3132//#define TIME_TEST // print the used time of each subroutine
    32 //#define ENDWALKS //print the size of the last omega-homogenoues Grï¿œbner basis
     33//#define ENDWALKS //print the size of the last omega-homogenoues Groebner basis
    3334
    3435/* includes */
     
    3839#include <time.h>
    3940#include <sys/time.h>
     41#include <math.h>
    4042#include <sys/stat.h>
    4143#include <unistd.h>
    42 #include <stdio.h>
    4344#include <float.h>
    4445#include <misc/mylimits.h>
     
    106107clock_t xftostd, xtextra, xftinput, to;
    107108
    108 /*2
    109 *utilities for TSet, LSet
    110 */
     109/****************************
     110 * utilities for TSet, LSet *
     111 ****************************/
    111112inline static intset initec (int maxnr)
    112113{
     
    123124}
    124125
    125 #if 0 /*unused*/
    126 /*2
    127 *construct the set s from F u {P}
    128 */
     126/************************************
     127 * construct the set s from F u {P} *
     128 ************************************/
     129// unused
     130#if 0
    129131static void initSSpecialCC (ideal F, ideal Q, ideal P,kStrategy strat)
    130132{
     
    141143  strat->S=strat->Shdl->m;
    142144
    143   /*- put polys into S -*/
     145  // - put polys into S -
    144146  if (Q!=NULL)
    145147  {
     
    183185    }
    184186  }
    185   /*- put polys into S -*/
     187  //- put polys into S -
    186188  for (i=0; i<IDELEMS(F); i++)
    187189  {
     
    276278#endif
    277279
    278 /*2
    279 *interreduces F
    280 */
     280/*****************
     281 *interreduce F  *
     282 *****************/
    281283static ideal kInterRedCC(ideal F, ideal Q)
    282284{
     
    295297  initBuchMoraCrit(strat);
    296298  strat->NotUsedAxis = (BOOLEAN *)omAlloc((currRing->N+1)*sizeof(BOOLEAN));
    297   for (j=currRing->N; j>0; j--) strat->NotUsedAxis[j] = TRUE;
     299  for(j=currRing->N; j>0; j--)
     300  {
     301    strat->NotUsedAxis[j] = TRUE;
     302  }
    298303  strat->enterS      = enterSBba;
    299304  strat->posInT      = posInT0;
     
    305310  strat->R           = initR();
    306311  strat->sevT        = initsevT();
    307   if (rHasLocalOrMixedOrdering_currRing())   strat->honey = TRUE;
    308 
     312  if(rHasLocalOrMixedOrdering_currRing())
     313  {
     314    strat->honey = TRUE;
     315  }
    309316
    310317  //initSCC(F,Q,strat);
     
    316323  tininitS=tininitS+clock()-timetmp;//22.01.02
    317324  */
    318   if (TEST_OPT_REDSB)
     325  if(TEST_OPT_REDSB)
     326  {
    319327    strat->noTailReduction=FALSE;
    320 
     328  }
    321329  updateS(TRUE,strat);
    322330
    323   if (TEST_OPT_REDSB && TEST_OPT_INTSTRATEGY)
     331  if(TEST_OPT_REDSB && TEST_OPT_INTSTRATEGY)
     332  {
    324333    completeReduce(strat);
     334  }
    325335  pDelete(&strat->kHEdge);
    326336  omFreeSize((ADDRESS)strat->T,strat->tmax*sizeof(TObject));
     
    332342  omfree(strat->R);
    333343
    334   if (strat->fromQ)
    335   {
    336     for (j=0;j<IDELEMS(strat->Shdl);j++)
    337     {
    338       if(strat->fromQ[j]) pDelete(&strat->Shdl->m[j]);
     344  if(strat->fromQ)
     345  {
     346    for(j=0; j<IDELEMS(strat->Shdl); j++)
     347    {
     348      if(strat->fromQ[j])
     349      {
     350        pDelete(&strat->Shdl->m[j]);
     351      }
    339352    }
    340353    omFreeSize((ADDRESS)strat->fromQ,IDELEMS(strat->Shdl)*sizeof(int));
    341     strat->fromQ=NULL;
     354    strat->fromQ = NULL;
    342355  }
    343356//  if (TEST_OPT_PROT)
     
    353366}
    354367
    355 #if 0 /*unused*/
     368//unused
     369#if 0
    356370static void TimeString(clock_t tinput, clock_t tostd, clock_t tif,clock_t tstd,
    357371                       clock_t tlf,clock_t tred, clock_t tnw, int step)
    358372{
    359373  double totm = ((double) (clock() - tinput))/1000000;
    360   double ostd,mostd, mif, mstd, mextra, mlf, mred, mnw, mxif,mxstd,mxlf,mxred,mxnw,tot;
    361 
     374  double ostd,mostd, mif, mstd, mlf, mred, mnw, mxif,mxstd,mxlf,mxred,mxnw,tot;
     375  // double mextra
    362376  Print("\n// total time = %.2f sec", totm);
    363377  Print("\n// tostd = %.2f sec = %.2f", ostd=((double) tostd)/1000000,
     
    393407#endif
    394408
    395 #if 0 /*unused*/
     409//unused
     410#if 0
    396411static void TimeStringFractal(clock_t tinput, clock_t tostd, clock_t tif,clock_t tstd,
    397412                       clock_t textra, clock_t tlf,clock_t tred, clock_t tnw)
     
    428443  Print("\n//  ideal %s =  ", st);
    429444  for(i=0; i<nL-1; i++)
     445  {
    430446    Print(" %s, ", pString(L->m[i]));
    431 
     447  }
    432448  Print(" %s;", pString(L->m[nL-1]));
    433449}
    434450
    435 #if 0 /*unused*/
     451//unused
     452//#if 0
    436453static void headidString(ideal L, char* st)
    437454{
     
    440457  Print("\n//  ideal %s =  ", st);
    441458  for(i=0; i<nL-1; i++)
     459  {
    442460    Print(" %s, ", pString(pHead(L->m[i])));
    443 
     461  }
    444462  Print(" %s;", pString(pHead(L->m[nL-1])));
    445463}
    446 #endif
    447 
    448 #if 0 /*unused*/
     464//#endif
     465
     466//unused
     467//#if 0
    449468static void idElements(ideal L, char* st)
    450469{
     
    454473  Print("\n//  #monoms of %s =  ", st);
    455474  for(i=0; i<nL; i++)
     475  {
    456476    K[i] = pLength(L->m[i]);
    457 
    458   int j, nsame, nk=0;
     477  }
     478  int j, nsame;
     479  // int  nk=0;
    459480  for(i=0; i<nL; i++)
    460481  {
     
    462483    {
    463484      nsame = 1;
    464       for(j=i+1; j<nL; j++){
    465         if(K[j]==K[i]){
     485      for(j=i+1; j<nL; j++)
     486      {
     487        if(K[j]==K[i])
     488        {
    466489          nsame ++;
    467490          K[j]=0;
     
    469492      }
    470493      if(nsame == 1)
     494      {
    471495        Print("%d, ",K[i]);
     496      }
    472497      else
     498      {
    473499        Print("%d[%d], ", K[i], nsame);
     500      }
    474501    }
    475502  }
    476503  omFree(K);
    477504}
    478 #endif
     505//#endif
    479506
    480507
     
    482509{
    483510  int nV = iv->length()-1;
    484   //Print("\n// vector %s =  (", ch);
    485511  Print("\n// intvec %s =  ", ch);
    486512
    487513  for(int i=0; i<nV; i++)
     514  {
    488515    Print("%d, ", (*iv)[i]);
     516  }
    489517  Print("%d;", (*iv)[nV]);
    490518}
    491519
    492 #if 0 /*unused*/
     520//unused
     521//#if 0
    493522static void MivString(intvec* iva, intvec* ivb, intvec* ivc)
    494523{
     
    497526  PrintS("\n//  (");
    498527  for(i=0; i<nV; i++)
     528  {
    499529    Print("%d, ", (*iva)[i]);
     530  }
    500531  Print("%d) ==> (", (*iva)[nV]);
    501 
    502532  for(i=0; i<nV; i++)
     533  {
    503534    Print("%d, ", (*ivb)[i]);
     535  }
    504536  Print("%d) := (", (*ivb)[nV]);
    505537
    506538  for(i=0; i<nV; i++)
     539  {
    507540    Print("%d, ", (*ivc)[i]);
     541  }
    508542  Print("%d)", (*ivc)[nV]);
    509543}
    510 #endif
    511 
    512 // returns gcd of integers a and b
     544//#endif
     545
     546/********************************************************************
     547 * returns gcd of integers a and b                                  *
     548 ********************************************************************/
    513549static inline long gcd(const long a, const long b)
    514550{
     
    516552  //assume(p0 >= 0 && p1 >= 0);
    517553  if(p0 < 0)
     554  {
    518555    p0 = -p0;
    519 
     556  }
    520557  if(p1 < 0)
     558  {
    521559    p1 = -p1;
    522 
     560  }
    523561  while(p1 != 0)
    524562  {
     
    530568}
    531569
    532 // cancel gcd of integers zaehler and nenner
     570/*********************************************
     571 * cancel gcd of integers zaehler and nenner *
     572 *********************************************/
    533573static void cancel(mpz_t zaehler, mpz_t nenner)
    534574{
     
    544584}
    545585
    546 #if 0 /*unused*/
    547 /* 23.07.03 */
     586//unused
     587#if 0
    548588static int isVectorNeg(intvec* omega)
    549589{
     
    551591
    552592  for(i=omega->length(); i>=0; i--)
     593  {
    553594    if((*omega)[i]<0)
     595    {
    554596      return 1;
    555 
     597    }
     598  }
    556599  return 0;
    557600}
     
    587630  if(mpz_cmp(zsum, sing_int)>0)
    588631  {
    589     if(Overflow_Error ==  FALSE) {
     632    if(Overflow_Error ==  FALSE)
     633    {
    590634      PrintLn();
    591635      PrintS("\n// ** OVERFLOW in \"MwalkInitialForm\": ");
     
    596640  }
    597641
     642  mpz_clear(zmul);
     643  mpz_clear(zvec);
     644  mpz_clear(zsum);
     645  mpz_clear(sing_int);
     646
    598647  return wgrad;
    599648}
     
    613662
    614663    if (maxtemp > max)
     664    {
    615665      max = maxtemp;
     666    }
    616667  }
    617668  return max;
     
    620671
    621672/********************************************************************
    622  * compute a weight degree of a monomial p w.r.t. a weight_vector *
     673 * compute a weight degree of a monomial p w.r.t. a weight_vector   *
    623674 ********************************************************************/
    624675static  void  MLmWeightedDegree_gmp(mpz_t result, const poly p, intvec* weight)
     
    628679  mpz_init_set_ui(sing_int,  2147483647);
    629680
    630   int i, wgrad;
     681  int i;
    631682
    632683  mpz_t zmul;
     
    657708{
    658709  if(g == NULL)
     710  {
    659711    return NULL;
    660 
     712  }
    661713  mpz_t max; mpz_init(max);
    662714  mpz_t maxtmp; mpz_init(maxtmp);
     
    676728      in_w_g = pHead(hg);
    677729    }
    678     else {
     730    else
     731    {
    679732      if(mpz_cmp(maxtmp, max)==0)
     733      {
    680734        in_w_g = pAdd(in_w_g, pHead(hg));
     735      }
    681736    }
    682737  }
     
    696751
    697752  for(i=nG-1; i>=0; i--)
     753  {
    698754    Gomega->m[i] = MpolyInitialForm(G->m[i], ivw);
    699 
     755  }
    700756  if(Overflow_Error == FALSE)
     757  {
    701758    Overflow_Error = nError;
    702 
     759  }
    703760  return Gomega;
    704761}
     
    706763/************************************************************************
    707764 * test whether the weight vector iv is in the cone of the ideal G      *
    708  *            i.e. are in(in_w(g)) =? in(g), for all g in G             *
     765 *     i.e. test whether in(in_w(g)) = in(g) for all g in G             *
    709766 ************************************************************************/
    710767
     
    731788    {
    732789      pDelete(&mi);
    733 
    734790      if(Overflow_Error == FALSE)
     791      {
    735792        Overflow_Error = nError;
    736 
     793      }
    737794      return 0;
    738795    }
     
    740797    {
    741798      pDelete(&mi);
    742 
    743799      if(Overflow_Error == FALSE)
     800      {
    744801        Overflow_Error = nError;
    745 
     802      }
    746803      return 0;
    747804    }
    748 
    749805    pDelete(&mi);
    750806  }
    751807
    752808  if(Overflow_Error == FALSE)
     809  {
    753810    Overflow_Error = nError;
    754 
     811  }
    755812  return 1;
    756813}
    757814
    758 
    759 //compute a least common multiple of two integers
     815/***************************************************
     816 * compute a least common multiple of two integers *
     817 ***************************************************/
    760818static inline long Mlcm(long &i1, long &i2)
    761819{
     
    775833
    776834  for(i=n-1; i>=0; i--)
     835    {
    777836    result += (*a)[i] * (*b)[i];
    778 
     837    }
    779838  return result;
    780839}
    781840
    782 
     841/*****************************************************
     842 * Substract two given intvecs componentwise         *
     843 *****************************************************/
    783844static intvec* MivSub(intvec* a, intvec* b)
    784845{
     
    788849
    789850  for(i=n-1; i>=0; i--)
     851  {
    790852    (*result)[i] = (*a)[i] - (*b)[i];
    791 
     853  }
    792854  return result;
    793855}
    794856
    795 /**21.10.00*******************************************
     857/*****************************************************
    796858 * return the "intvec" lead exponent of a polynomial *
    797859 *****************************************************/
     
    802864
    803865  for(i=nR-1; i>=0; i--)
     866  {
    804867    (*result)[i] = pGetExp(f,i+1);
    805 
     868  }
    806869  return result;
    807870}
    808871
    809 /* return 1, if two given intvecs are the same, otherwise 0*/
     872/*****************************************************
     873 * Compare two given intvecs and return 1, if they   *
     874 * are the same, otherwise 0                         *
     875 *****************************************************/
    810876int MivSame(intvec* u , intvec* v)
    811877{
     
    815881
    816882  for (i=0; i<niv; i++)
     883  {
    817884    if ((*u)[i] != (*v)[i])
     885    {
    818886      return 0;
    819 
     887    }
     888  }
    820889  return 1;
    821890}
    822891
     892/******************************************************
     893 * Compare 3 given intvecs and return 0, if the first *
     894 * and the second are the same. Return 1, if the      *
     895 * the second and the third are the same, otherwise 2 *
     896 ******************************************************/
    823897int M3ivSame(intvec* temp, intvec* u , intvec* v)
    824898{
     
    826900
    827901  if((MivSame(temp, u)) == 1)
     902  {
    828903    return 0;
    829 
     904  }
    830905  if((MivSame(temp, v)) == 1)
     906  {
    831907    return 1;
    832 
     908  }
    833909  return 2;
    834910}
    835911
    836 
    837 /* compute a Groebner basis of an ideal */
     912/*****************************************************
     913 * compute a Groebner basis of an ideal              *
     914 *****************************************************/
    838915static ideal MstdCC(ideal G)
    839916{
     
    848925}
    849926
    850 
    851 /* compute a Groebner basis of a homogenoues ideal */
     927/*****************************************************
     928 * compute a Groebner basis of an homogeneous ideal  *
     929 *****************************************************/
    852930static ideal MstdhomCC(ideal G)
    853931{
     
    868946intvec* MivMatrixOrder(intvec* iv)
    869947{
    870   int i,j, nR = iv->length();
     948  int i, nR = iv->length();
     949 
    871950  intvec* ivm = new intvec(nR*nR);
    872951
    873952  for(i=0; i<nR; i++)
     953  {
    874954    (*ivm)[i] = (*iv)[i];
    875 
     955  }
    876956  for(i=1; i<nR; i++)
     957  {
    877958    (*ivm)[i*nR+i-1] = 1;
    878 
     959  }
    879960  return ivm;
    880961}
    881962
    882 /* return intvec = (1, ..., 1) */
     963/*******************************
     964 * return intvec = (1, ..., 1) *
     965 *******************************/
    883966intvec* Mivdp(int nR)
    884967{
     
    887970
    888971  for(i=nR-1; i>=0; i--)
     972  {
    889973    (*ivm)[i] = 1;
    890 
     974  }
    891975  return ivm;
    892976}
    893977
    894 /* return intvvec = (1,0, ..., 0) */
     978/**********************************
     979 * return intvvec = (1,0, ..., 0) *
     980 **********************************/
    895981intvec* Mivlp(int nR)
    896982{
    897   int i;
    898983  intvec* ivm = new intvec(nR);
    899984  (*ivm)[0] = 1;
     
    902987}
    903988
    904 /**** 28.10.02  print the max total degree and the max coefficient of G***/
    905 #if 0 /*unused*/
     989//unused
     990/*****************************************************************************
     991 * print the max total degree and the max coefficient of G                   *
     992 *****************************************************************************/
     993#if 0
    906994static void checkComplexity(ideal G, char* cG)
    907995{
    908996  int nV = currRing->N;
    909997  int nG = IDELEMS(G);
    910   intvec* ivUnit = Mivdp(nV);//19.02
    911   int i,j, tmpdeg, maxdeg=0;
     998  intvec* ivUnit = Mivdp(nV);
     999  int i, tmpdeg, maxdeg=0;
    9121000  number tmpcoeff , maxcoeff=currRing->cf->nNULL;
    9131001  poly p;
     
    9151003  {
    9161004    tmpdeg = MwalkWeightDegree(G->m[i], ivUnit);
    917     if (tmpdeg > maxdeg )
     1005    if(tmpdeg > maxdeg )
     1006    {
    9181007      maxdeg = tmpdeg;
     1008    }
    9191009  }
    9201010
     
    9271017      tmpcoeff = pGetCoeff(p);
    9281018      if(nGreater(tmpcoeff,maxcoeff))
     1019      {
    9291020         maxcoeff = nCopy(tmpcoeff);
     1021      }
    9301022      pIter(p);
    9311023    }
     
    9341026  p = pNSet(maxcoeff);
    9351027  char* pStr = pString(p);
     1028  delete ivUnit;
    9361029  Print("// max total degree of %s = %d\n",cG, maxdeg);
    9371030  Print("// max coefficient of %s  = %s", cG, pStr);//ing(p));
     
    9401033}
    9411034#endif
    942 
    9431035
    9441036/*****************************************************************************
     
    9531045* degree which smaller than the numbers of variables                         *
    9541046******************************************************************************/
    955 /* ivtarget is a matrix order of a degree reverse lex. order */
    9561047intvec* MPertVectors(ideal G, intvec* ivtarget, int pdeg)
    9571048{
     1049  // ivtarget is a matrix order of a degree reverse lex. order
    9581050  int nV = currRing->N;
    9591051  //assume(pdeg <= nV && pdeg >= 0);
     
    9631055
    9641056
    965   //Checking that the perturbed degree is valid
     1057  // Check that the perturbed degree is valid
    9661058  if(pdeg > nV || pdeg <= 0)
    9671059  {
     
    9721064
    9731065  if(pdeg == 1)
     1066  {
    9741067    return ivtarget;
    975 
    976   mpz_t *pert_vector=(mpz_t*)omAlloc(nV*sizeof(mpz_t));
     1068  }
     1069  mpz_t *pert_vector = (mpz_t*)omAlloc(nV*sizeof(mpz_t));
     1070  //mpz_t *pert_vector1 = (mpz_t*)omAlloc(nV*sizeof(mpz_t));
    9771071
    9781072  for(i=0; i<nV; i++)
     1073  {
    9791074    mpz_init_set_si(pert_vector[i], (*ivtarget)[i]);
    980 
    981 
     1075   // mpz_init_set_si(pert_vector1[i], (*ivtarget)[i]);
     1076  }
    9821077  // Calculate max1 = Max(A2)+Max(A3)+...+Max(Apdeg),
    9831078  // where the Ai are the i-te rows of the matrix target_ord.
    984 
    9851079  int ntemp, maxAi, maxA=0;
    9861080  for(i=1; i<pdeg; i++)
    9871081  {
    9881082    maxAi = (*ivtarget)[i*nV];
    989     if(maxAi<0) maxAi = -maxAi;
    990 
     1083    if(maxAi<0)
     1084    {
     1085      maxAi = -maxAi;
     1086    }
    9911087    for(j=i*nV+1; j<(i+1)*nV; j++)
    9921088    {
    9931089      ntemp = (*ivtarget)[j];
    994       if(ntemp < 0) ntemp = -ntemp;
    995 
     1090      if(ntemp < 0)
     1091      {
     1092        ntemp = -ntemp;
     1093      }
    9961094      if(ntemp > maxAi)
     1095      {
    9971096        maxAi = ntemp;
     1097      }
    9981098    }
    9991099    maxA += maxAi;
     
    10111111  for(i=nG-1; i>=0; i--)
    10121112  {
    1013     mpz_set_ui(maxdeg, MwalkWeightDegree(G->m[i], ivUnit));
    1014     if (mpz_cmp(maxdeg,  tot_deg) > 0 )
    1015       mpz_set(tot_deg, maxdeg);
     1113     mpz_set_ui(maxdeg, MwalkWeightDegree(G->m[i], ivUnit));
     1114     if (mpz_cmp(maxdeg,  tot_deg) > 0 )
     1115     {
     1116       mpz_set(tot_deg, maxdeg);
     1117     }
    10161118  }
    10171119
     
    10211123
    10221124
    1023   //xx1.06.02 takes  "small" inveps
     1125  // takes  "small" inveps
    10241126#ifdef INVEPS_SMALL_IN_MPERTVECTOR
    10251127  if(mpz_cmp_ui(inveps, pdeg)>0 && pdeg > 3)
    10261128  {
    1027     /*
    1028       Print("\n// choose the\"small\" inverse epsilon := %d / %d = ",
    1029       mpz_get_si(inveps), pdeg);
    1030     */
     1129    //  Print("\n// choose the\"small\" inverse epsilon := %d / %d = ", mpz_get_si(inveps), pdeg);
    10311130    mpz_fdiv_q_ui(inveps, inveps, pdeg);
    1032     //mpz_out_str(stdout, 10, inveps);
     1131    // mpz_out_str(stdout, 10, inveps);
    10331132  }
    10341133#else
    1035   //PrintS("\n// the \"big\" inverse epsilon: ");
     1134  // PrintS("\n// the \"big\" inverse epsilon: ");
    10361135  mpz_out_str(stdout, 10, inveps);
    10371136#endif
     
    10391138  // pert(A1) = inveps^(pdeg-1)*A1 + inveps^(pdeg-2)*A2+...+A_pdeg,
    10401139  // pert_vector := A1
    1041   for ( i=1; i < pdeg; i++ )
     1140  for( i=1; i < pdeg; i++ )
     1141  {
    10421142    for(j=0; j<nV; j++)
    10431143    {
    10441144      mpz_mul(pert_vector[j], pert_vector[j], inveps);
    10451145      if((*ivtarget)[i*nV+j]<0)
     1146      {
    10461147        mpz_sub_ui(pert_vector[j], pert_vector[j],-(*ivtarget)[i*nV+j]);
     1148      }
    10471149      else
     1150      {
    10481151        mpz_add_ui(pert_vector[j], pert_vector[j],(*ivtarget)[i*nV+j]);
    1049     }
    1050 
     1152      }
     1153    }
     1154  }
    10511155  mpz_t ztemp;
    10521156  mpz_init(ztemp);
     
    10561160    mpz_gcd(ztemp, ztemp, pert_vector[i]);
    10571161    if(mpz_cmp_si(ztemp, 1)  == 0)
     1162    {
    10581163      break;
     1164    }
    10591165  }
    10601166  if(mpz_cmp_si(ztemp, 1) != 0)
     1167  {
    10611168    for(i=0; i<nV; i++)
     1169    {
    10621170      mpz_divexact(pert_vector[i], pert_vector[i], ztemp);
    1063 
     1171    }
     1172  }
     1173
     1174  intvec *pert_vector1= new intvec(nV);
     1175  j = 0;
     1176  for(i=0; i<nV; i++)
     1177  {
     1178    (* pert_vector1)[i] = mpz_get_si(pert_vector[i]);
     1179    (* pert_vector1)[i] = 0.1*(* pert_vector1)[i];
     1180    (* pert_vector1)[i] = floor((* pert_vector1)[i] + 0.5);
     1181    if((* pert_vector1)[i] == 0)
     1182    {
     1183      j++;
     1184    }
     1185  }
     1186  if(j > nV - 1)
     1187  {
     1188    // Print("\n//  MPertVectors: geaenderter vector gleich Null! \n");
     1189    delete pert_vector1;
     1190    goto CHECK_OVERFLOW;
     1191  }
     1192
     1193// check that the perturbed weight vector lies in the Groebner cone
     1194  if(test_w_in_ConeCC(G,pert_vector1) != 0)
     1195  {
     1196    // Print("\n//  MPertVectors: geaenderter vector liegt in Groebnerkegel! \n");
     1197    for(i=0; i<nV; i++)
     1198    {
     1199      mpz_set_si(pert_vector[i], (*pert_vector1)[i]);
     1200    }
     1201  }
     1202  else
     1203  {
     1204    //Print("\n// MpertVectors: geaenderter vector liegt nicht in Groebnerkegel! \n");
     1205  }
     1206  delete pert_vector1;
     1207
     1208  CHECK_OVERFLOW:
    10641209  intvec* result = new intvec(nV);
     1210
    10651211  /* 2147483647 is max. integer representation in SINGULAR */
    10661212  mpz_t sing_int;
     
    10711217  {
    10721218    (*result)[i] = mpz_get_si(pert_vector[i]);
    1073 
    10741219    if(mpz_cmp(pert_vector[i], sing_int)>=0)
    10751220    {
     
    10951240  mpz_clear(sing_int);
    10961241  omFree(pert_vector);
     1242  //omFree(pert_vector1);
     1243  mpz_clear(tot_deg);
     1244  mpz_clear(maxdeg);
     1245  mpz_clear(inveps);
     1246
     1247  rComplete(currRing);
     1248  for(j=0; j<IDELEMS(G); j++)
     1249  {
     1250    poly p=G->m[j];
     1251    while(p!=NULL)
     1252    {
     1253      p_Setm(p,currRing); pIter(p);
     1254    }
     1255  }
    10971256  return result;
    10981257}
    10991258
    1100 
    1101 /* ivtarget is a matrix order of the lex. order */
     1259/*****************************************************************************
     1260 * The following procedure returns                                           *
     1261 *     Pert(A1) = 1/eps^(pdeg-1)*A_1 + 1/eps^(pdeg-2)*A_2+...+A_pdeg,        *
     1262 * where the A_i are the i-th rows of the matrix target_ord and              *
     1263 *     1/eps > deg(p)*(max(A_2) + max(A_3)+...+max(A_pdeg))                  *
     1264 *****************************************************************************/
    11021265intvec* MPertVectorslp(ideal G, intvec* ivtarget, int pdeg)
    11031266{
     1267  // ivtarget is a matrix order of the lex. order
    11041268  int nV = currRing->N;
    11051269  //assume(pdeg <= nV && pdeg >= 0);
     
    11151279  }
    11161280  for(i=0; i<nV; i++)
     1281  {
    11171282    (*pert_vector)[i]=(*ivtarget)[i];
    1118 
     1283  }
    11191284  if(pdeg == 1)
     1285  {
    11201286    return pert_vector;
    1121 
     1287  }
    11221288  // Calculate max1 = Max(A2)+Max(A3)+...+Max(Apdeg),
    11231289  // where the Ai are the i-te rows of the matrix target_ord.
     
    11301296      ntemp = (*ivtarget)[j];
    11311297      if(ntemp > maxAi)
     1298      {
    11321299        maxAi = ntemp;
     1300      }
    11331301    }
    11341302    maxA += maxAi;
     
    11411309  for(i=nG-1; i>=0; i--)
    11421310  {
    1143     //maxdeg = pTotaldegree(G->m[i], currRing); //it's wrong for ex1,2,rose
     1311    // maxdeg = pTotaldegree(G->m[i], currRing); //it's wrong for ex1,2,rose
    11441312    maxdeg = MwalkWeightDegree(G->m[i], ivUnit);
    11451313    if (maxdeg > tot_deg )
     1314    {
    11461315      tot_deg = maxdeg;
     1316    }
    11471317  }
    11481318  delete ivUnit;
     
    11501320  inveps = (tot_deg * maxA) + 1;
    11511321
    1152   //9.10.01
    11531322#ifdef INVEPS_SMALL_IN_FRACTAL
    1154   /*
    1155     Print("\n// choose the\"small\" inverse epsilon := %d / %d = ",
    1156     inveps, pdeg);
    1157   */
     1323  //  Print("\n// choose the\"small\" inverse epsilon := %d / %d = ", inveps, pdeg);
    11581324  if(inveps > pdeg && pdeg > 3)
     1325  {
    11591326    inveps = inveps / pdeg;
    1160 
    1161   //Print(" %d", inveps);
     1327  }
     1328  // Print(" %d", inveps);
    11621329#else
    11631330  PrintS("\n// the \"big\" inverse epsilon %d", inveps);
    11641331#endif
    11651332
    1166   // Pert(A1) = inveps^(pdeg-1)*A1 + inveps^(pdeg-2)*A2+...+A_pdeg,
     1333  // Pert(A1) = inveps^(pdeg-1)*A1 + inveps^(pdeg-2)*A2+...+A_pdeg
    11671334  for ( i=1; i < pdeg; i++ )
     1335  {
    11681336    for(j=0; j<nV; j++)
     1337    {
    11691338      (*pert_vector)[j] = inveps*((*pert_vector)[j]) + (*ivtarget)[i*nV+j];
     1339    }
     1340  }
    11701341
    11711342  int temp = (*pert_vector)[0];
     
    11741345    temp = gcd(temp, (*pert_vector)[i]);
    11751346    if(temp == 1)
     1347    {
    11761348      break;
     1349    }
    11771350  }
    11781351  if(temp != 1)
     1352  {
    11791353    for(i=0; i<nV; i++)
     1354    {
    11801355      (*pert_vector)[i] = (*pert_vector)[i] / temp;
     1356    }
     1357  }
    11811358
    11821359  intvec* result = pert_vector;
    1183   pert_vector = NULL;
     1360  delete pert_vector;
    11841361  return result;
    11851362}
    11861363
    1187 
    1188 /* define a lexicographic order matrix as intvec */
     1364/*****************************************************************************
     1365 * define a lexicographic order matrix as intvec                             *
     1366 *****************************************************************************/
    11891367intvec* MivMatrixOrderlp(int nV)
    11901368{
     
    11931371
    11941372  for(i=0; i<nV; i++)
     1373  {
    11951374    (*ivM)[i*nV + i] = 1;
    1196 
     1375  }
    11971376  return(ivM);
    11981377}
    11991378
    1200 /* define a rlex order (dp) matrix as intvec */
     1379
     1380/*****************************************************************************
     1381 * define a reverse lexicographic order (dp) matrix as intvec                *
     1382 *****************************************************************************/
    12011383intvec* MivMatrixOrderdp(int nV)
    12021384{
     
    12051387
    12061388  for(i=0; i<nV; i++)
     1389  {
    12071390    (*ivM)[i] = 1;
    1208 
     1391  }
    12091392  for(i=1; i<nV; i++)
     1393  {
    12101394    (*ivM)[(i+1)*nV - i] = -1;
    1211 
     1395  }
    12121396  return(ivM);
    12131397}
    12141398
    1215 //creates an intvec of the monomial order Wp(ivstart)
     1399/*****************************************************************************
     1400 * creates an intvec of the monomial order Wp(ivstart)                       *
     1401 *****************************************************************************/
    12161402intvec* MivWeightOrderlp(intvec* ivstart)
    12171403{
     
    12211407
    12221408  for(i=0; i<nV; i++)
     1409  {
    12231410    (*ivM)[i] = (*ivstart)[i];
    1224 
     1411  }
    12251412  for(i=1; i<nV; i++)
     1413  {
    12261414    (*ivM)[i*nV + i-1] = 1;
    1227 
     1415  }
    12281416  return(ivM);
    12291417}
    12301418
     1419/*****************************************************************************
     1420 * creates an intvec of the monomial order dp(ivstart)                       *
     1421 *****************************************************************************/
    12311422intvec* MivWeightOrderdp(intvec* ivstart)
    12321423{
     
    12361427
    12371428  for(i=0; i<nV; i++)
     1429  {
    12381430    (*ivM)[i] = (*ivstart)[i];
    1239 
     1431  }
    12401432  for(i=0; i<nV; i++)
     1433  {
    12411434    (*ivM)[nV+i] = 1;
    1242 
     1435  }
    12431436  for(i=2; i<nV; i++)
     1437  {
    12441438    (*ivM)[(i+1)*nV - i] = -1;
    1245 
     1439  }
    12461440  return(ivM);
    12471441}
    12481442
    1249 #if 0 /*unused*/
     1443//unused
     1444#if 0
    12501445static intvec* MatrixOrderdp(int nV)
    12511446{
     
    12541449
    12551450  for(i=0; i<nV; i++)
     1451  {
    12561452    (*ivM)[i] = 1;
    1257 
     1453  }
    12581454  for(i=1; i<nV; i++)
     1455  {
    12591456    (*ivM)[(i+1)*nV - i] = -1;
    1260 
     1457  }
    12611458  return(ivM);
    12621459}
     
    12671464  int i;
    12681465  intvec* ivM = new intvec(nV);
    1269 
    12701466  for(i=nV-1; i>=0; i--)
     1467  {
    12711468    (*ivM)[i] = 1;
    1272 
     1469  }
    12731470  return(ivM);
    12741471}
     
    12791476*************************************************************************/
    12801477int Xnlev;
    1281 
    12821478intvec* Mfpertvector(ideal G, intvec* ivtarget)
    12831479{
     
    12861482  int niv = nV*nV;
    12871483
    1288   // Calculate max1 = Max(A2) + Max(A3) + ... + Max(AnV),
     1484
     1485  // Calculate maxA = Max(A2) + Max(A3) + ... + Max(AnV),
    12891486  // where the Ai are the i-te rows of the matrix 'targer_ord'.
    12901487  int ntemp, maxAi, maxA=0;
     
    12921489  {
    12931490    maxAi = (*ivtarget)[i*nV];
    1294     if(maxAi<0) maxAi = -maxAi;
    1295 
     1491    if(maxAi<0)
     1492    {
     1493      maxAi = -maxAi;
     1494    }
    12961495    for(j=i*nV+1; j<(i+1)*nV; j++)
    12971496    {
    12981497      ntemp = (*ivtarget)[j];
    1299       if(ntemp < 0) ntemp = -ntemp;
    1300 
     1498      if(ntemp < 0)
     1499      {
     1500        ntemp = -ntemp;
     1501      }
    13011502      if(ntemp > maxAi)
     1503      {
    13021504        maxAi = ntemp;
     1505      }
    13031506    }
    13041507    maxA = maxA + maxAi;
     
    13061509  intvec* ivUnit = Mivdp(nV);
    13071510
    1308   // Calculate inveps = 1/eps, where 1/eps > deg(p)*max1 for all p in G.
     1511  // Calculate inveps = 1/eps, where 1/eps > deg(p)*maxA for all p in G.
    13091512  mpz_t tot_deg; mpz_init(tot_deg);
    13101513  mpz_t maxdeg; mpz_init(maxdeg);
     
    13161519    mpz_set_ui(maxdeg, MwalkWeightDegree(G->m[i], ivUnit));
    13171520    if (mpz_cmp(maxdeg,  tot_deg) > 0 )
     1521    {
    13181522      mpz_set(tot_deg, maxdeg);
     1523    }
    13191524  }
    13201525
     
    13241529  mpz_add_ui(inveps, inveps, 1);
    13251530
    1326   //xx1.06.02 takes  "small" inveps
     1531  // takes  "small" inveps
    13271532#ifdef INVEPS_SMALL_IN_FRACTAL
    13281533  if(mpz_cmp_ui(inveps, nV)>0 && nV > 3)
     1534  {
    13291535    mpz_cdiv_q_ui(inveps, inveps, nV);
    1330 
     1536  }
    13311537  //PrintS("\n// choose the \"small\" inverse epsilon!");
    13321538#endif
     
    13381544  mpz_t *pert_vector=(mpz_t *)omAlloc(niv*sizeof(mpz_t));
    13391545
    1340   for(i=0; i<nV; i++)
     1546  for(i=0; i < nV; i++)
    13411547  {
    13421548    mpz_init_set_si(ivtemp[i], (*ivtarget)[i]);
     
    13451551
    13461552  mpz_t ztmp; mpz_init(ztmp);
    1347   BOOLEAN isneg = FALSE;
     1553  // BOOLEAN isneg = FALSE;
    13481554
    13491555  for(i=1; i<nV; i++)
     
    13521558    {
    13531559      mpz_mul(ztmp, inveps, ivtemp[j]);
    1354 
    13551560      if((*ivtarget)[i*nV+j]<0)
     1561      {
    13561562        mpz_sub_ui(ivtemp[j], ztmp, -(*ivtarget)[i*nV+j]);
     1563      }
    13571564      else
     1565      {
    13581566        mpz_add_ui(ivtemp[j], ztmp,(*ivtarget)[i*nV+j]);
     1567      }
    13591568    }
    13601569
    13611570    for(j=0; j<nV; j++)
     1571      {
    13621572      mpz_init_set(pert_vector[i*nV+j],ivtemp[j]);
     1573      }
    13631574  }
    13641575
     
    13681579
    13691580  intvec* result = new intvec(niv);
     1581  intvec* result1 = new intvec(niv);
    13701582  BOOLEAN nflow = FALSE;
    13711583
     
    13761588    mpz_gcd(ztmp, ztmp, pert_vector[i]);
    13771589    if(mpz_cmp_si(ztmp, 1)==0)
     1590    {
    13781591      break;
     1592    }
    13791593  }
    13801594
     
    13831597    mpz_divexact(pert_vector[i], pert_vector[i], ztmp);
    13841598    (* result)[i] = mpz_get_si(pert_vector[i]);
    1385 
     1599  }
     1600
     1601  j = 0;
     1602  for(i=0; i<nV; i++)
     1603  {
     1604    (* result1)[i] = mpz_get_si(pert_vector[i]);
     1605    (* result1)[i] = 0.1*(* result1)[i];
     1606    (* result1)[i] = floor((* result1)[i] + 0.5);
     1607    if((* result1)[i] == 0)
     1608    {
     1609      j++;
     1610    }
     1611  }
     1612  if(j > nV - 1)
     1613  {
     1614    // Print("\n//  MfPertwalk: geaenderter vector gleich Null! \n");
     1615    delete result1;
     1616    goto CHECK_OVERFLOW;
     1617  }
     1618
     1619// check that the perturbed weight vector lies in the Groebner cone
     1620  if(test_w_in_ConeCC(G,result1) != 0)
     1621  {
     1622    // Print("\n//  MfPertwalk: geaenderter vector liegt in Groebnerkegel! \n");
     1623    delete result;
     1624    result = result1;
     1625    for(i=0; i<nV; i++)
     1626    {
     1627      mpz_set_si(pert_vector[i], (*result1)[i]);
     1628    }
     1629  }
     1630  else
     1631  {
     1632    delete result1;
     1633    // Print("\n// Mfpertwalk: geaenderter vector liegt nicht in Groebnerkegel! \n");
     1634  }
     1635
     1636  CHECK_OVERFLOW:
     1637
     1638  for(i=0; i<niv; i++)
     1639  {
    13861640    if(mpz_cmp(pert_vector[i], sing_int)>0)
     1641    {
    13871642      if(nflow == FALSE)
    13881643      {
     
    13901645        nflow = TRUE;
    13911646        Overflow_Error = TRUE;
    1392 
    13931647        Print("\n// Xlev = %d and the %d-th element is", Xnlev,  i+1);
    13941648        PrintS("\n// ** OVERFLOW in \"Mfpertvector\": ");
     
    13971651        Print("\n//  So vector[%d] := %d is wrong!!", i+1, (*result)[i]);
    13981652      }
    1399   }
    1400 
     1653    }
     1654  }
    14011655  if(Overflow_Error == TRUE)
     1656  {
    14021657    ivString(result, "new_vector");
    1403 
     1658  }
    14041659  omFree(pert_vector);
    14051660  omFree(ivtemp);
    14061661  mpz_clear(ztmp);
    1407 
     1662  mpz_clear(tot_deg);
     1663  mpz_clear(maxdeg);
     1664  mpz_clear(inveps);
     1665  mpz_clear(sing_int);
     1666
     1667  rComplete(currRing);
     1668  for(j=0; j<IDELEMS(G); j++)
     1669  {
     1670    poly p=G->m[j];
     1671    while(p!=NULL)
     1672    {
     1673      p_Setm(p,currRing); pIter(p);
     1674    }
     1675  }
    14081676  return result;
    14091677}
    14101678
    14111679/****************************************************************
    1412  * Multiplikation of two ideals by elementwise                  *
     1680 * Multiplication of two ideals element by element              *
    14131681 * i.e. Let be A := (a_i) and B := (b_i), return C := (a_i*b_i) *
    14141682 *  destroy A, keeps B                                          *
     
    14191687
    14201688  if(A==NULL || B==NULL)
     1689  {
    14211690    return NULL;
    1422 
     1691  }
    14231692  if(mB < mA)
     1693  {
    14241694    mA = mB;
    1425 
     1695  }
    14261696  ideal result = idInit(mA, 1);
    14271697
     
    14311701      result->m[k] = pMult(A->m[i], pCopy(B->m[i]));
    14321702      A->m[i]=NULL;
    1433       if (result->m[k]!=NULL) k++;
     1703      if (result->m[k]!=NULL)
     1704      {
     1705        k++;
     1706      }
    14341707    }
    14351708
     
    14511724  ideal Mtmp = idLift(Gw, M, NULL, FALSE, TRUE, TRUE, NULL);
    14521725
    1453   //3.12.02 Note: if Gw is a GB, then isSB = TRUE, otherwise FALSE
    1454   //So, it is better, if one tests whether Gw is a GB
    1455   //in ideals.cc:
    1456   //idLift (ideal mod, ideal submod,ideal * rest, BOOLEAN goodShape,
     1726  // If Gw is a GB, then isSB = TRUE, otherwise FALSE
     1727  // So, it is better, if one tests whether Gw is a GB
     1728  // in ideals.cc:
     1729  // idLift (ideal mod, ideal submod,ideal * rest, BOOLEAN goodShape,
    14571730  //           BOOLEAN isSB,BOOLEAN divide,matrix * unit)
    14581731
    1459   /* Let be Mtmp = {m1,...,ms}, where mi=sum hij.in_gj, for all i=1,...,s
    1460      We compute F = {f1,...,fs}, where fi=sum hij.gj */
     1732  // Let be Mtmp = {m1,...,ms}, where mi=sum hij.in_gj, for all i=1,...,s
     1733  // We compute F = {f1,...,fs}, where fi=sum hij.gj
    14611734  int i, j, nM = IDELEMS(Mtmp);
    14621735  ideal idpol, idLG;
     
    14801753}
    14811754
    1482 
    1483 #if 0 /*unused*/
     1755//unused
     1756#if 0
    14841757static void checkidealCC(ideal G, char* Ch)
    14851758{
     
    14951768
    14961769    if(i != n)
     1770    {
    14971771      Print("%d, ", ntmp);
     1772    }
    14981773    else
     1774    {
    14991775      Print(" bzw. %d ", ntmp);
     1776    }
    15001777  }
    15011778  PrintS(" Monomen.\n");
     
    15051782#endif
    15061783
    1507 #if 0 /*unused*/
     1784//unused
     1785#if 0
    15081786static void HeadidString(ideal L, char* st)
    15091787{
     
    15121790  Print("//  The head terms of the ideal %s = ", st);
    15131791  for(i=0; i<nL; i++)
     1792  {
    15141793    Print(" %s, ", pString(pHead(L->m[i])));
    1515 
     1794  }
    15161795  Print(" %s;\n", pString(pHead(L->m[nL])));
    15171796}
     
    15221801  assume(iva->length() == ivb->length());
    15231802  int i;
    1524 
    15251803  for(i=iva->length()-1; i>=0; i--)
     1804  {
    15261805    if((*iva)[i] - (*ivb)[i] != 0)
     1806    {
    15271807      return 0;
    1528 
     1808    }
     1809  }
    15291810  return 1;
    15301811}
    15311812
    1532 
    1533 /*
    1534   compute a next weight vector between curr_weight and target_weight
    1535   with respect to an ideal G.
    1536 */
     1813/**********************************************
     1814 * Look for the smallest absolut value in vec *
     1815 **********************************************/
     1816static int MivAbsMax(intvec* vec)
     1817{
     1818  int i,k;
     1819  if((*vec)[0] < 0)
     1820  {
     1821    k = -(*vec)[0];
     1822  }
     1823  else
     1824  {
     1825    k = (*vec)[0];
     1826  }
     1827  for(i=1; i < (vec->length()); i++)
     1828  {
     1829    if((*vec)[i] < 0)
     1830    {
     1831      if(-(*vec)[i] > k)
     1832      {
     1833        k = -(*vec)[i];
     1834      }
     1835    }
     1836    else
     1837    {
     1838      if((*vec)[i] > k)
     1839      {
     1840        k = (*vec)[i];
     1841      }
     1842    }
     1843  }
     1844  return k;
     1845}
     1846
     1847/**********************************************************************
     1848 * Compute a next weight vector between curr_weight and target_weight *
     1849 * with respect to an ideal <G>.                                      *
     1850**********************************************************************/
    15371851static intvec* MwalkNextWeightCC(intvec* curr_weight, intvec* target_weight,
    15381852                                 ideal G)
     
    15451859
    15461860  int nRing = currRing->N;
    1547   int j, nG = IDELEMS(G);
     1861  int checkRed, j, kkk, nG = IDELEMS(G);
    15481862  intvec* ivtemp;
    15491863
     
    15581872  mpz_init(MwWd);
    15591873
    1560 
     1874  mpz_t sing_int;
     1875  mpz_init(sing_int);
     1876  mpz_set_si(sing_int,  2147483647);
     1877
     1878  mpz_t sing_int_half;
     1879  mpz_init(sing_int_half);
     1880  mpz_set_si(sing_int_half,  3*(1073741824/2));
     1881 
    15611882  mpz_t deg_w0_p1, deg_d0_p1;
    15621883  mpz_init(deg_w0_p1);
     
    15661887  mpz_init(sztn);
    15671888  mpz_init(sntz);
     1889
    15681890  mpz_t t_null;
    15691891  mpz_init(t_null);
     
    15721894  mpz_init(ggt);
    15731895
    1574   int tn0, tn1, tz1, ncmp, gcd_tmp, ntmp;
     1896  mpz_t dcw;
     1897  mpz_init(dcw);
     1898
     1899  //int tn0, tn1, tz1, ncmp, gcd_tmp, ntmp;
     1900  int gcd_tmp;
    15751901  intvec* diff_weight = MivSub(target_weight, curr_weight);
    15761902
    1577   poly g, gw;
     1903  intvec* diff_weight1 = MivSub(target_weight, curr_weight);
     1904  poly g;
     1905  //poly g, gw;
    15781906  for (j=0; j<nG; j++)
    15791907  {
     
    16111939            }
    16121940
    1613             //compute a simply fraction of s
     1941            //compute a simple fraction of s
    16141942            cancel(s_zaehler, s_nenner);
    16151943
     
    16371965    }
    16381966  }
    1639 
     1967//Print("\n// Alloc Size = %d \n", nRing*sizeof(mpz_t));
    16401968  mpz_t *vec=(mpz_t*)omAlloc(nRing*sizeof(mpz_t));
    1641 
    1642   /* there is no 0<t<1 and define the next weight vector that is equal to
    1643      the current weight vector */
     1969 
     1970
     1971  // there is no 0<t<1 and define the next weight vector that is equal to the current weight vector
    16441972  if(mpz_cmp(t_nenner, t_null) == 0)
    16451973  {
     1974    Print("\n//MwalkNextWeightCC: t_nenner ist Null!");
    16461975    delete diff_weight;
    16471976    diff_weight = ivCopy(curr_weight);//take memory
     
    16491978  }
    16501979
    1651   /* define the target vector as the next weight vector, if t = 1 */
     1980  // define the target vector as the next weight vector, if t = 1
    16521981  if(mpz_cmp_si(t_nenner, 1)==0 && mpz_cmp_si(t_zaehler,1)==0)
    16531982  {
     
    16571986  }
    16581987
    1659 
    1660   //14.08.03 simplify the both vectors  curr_weight and diff_weight (C-int)
     1988   checkRed = 0;
     1989
     1990  SIMPLIFY_GCD:
     1991
     1992  // simplify the vectors curr_weight and diff_weight (C-int)
    16611993  gcd_tmp = (*curr_weight)[0];
    16621994
     
    16651997    gcd_tmp = gcd(gcd_tmp, (*curr_weight)[j]);
    16661998    if(gcd_tmp == 1)
     1999    {
    16672000      break;
    1668   }
    1669 
     2001    }
     2002  }
    16702003  if(gcd_tmp != 1)
     2004  {
    16712005    for (j=0; j<nRing; j++)
    16722006    {
    16732007      gcd_tmp = gcd(gcd_tmp, (*diff_weight)[j]);
    16742008      if(gcd_tmp == 1)
     2009      {
    16752010        break;
    1676     }
    1677 
     2011      }
     2012    }
     2013  }
    16782014  if(gcd_tmp != 1)
     2015  {
    16792016    for (j=0; j<nRing; j++)
    16802017    {
     
    16822019      (*diff_weight)[j] =  (*diff_weight)[j]/gcd_tmp;
    16832020    }
     2021  }
     2022  if(checkRed > 0)
     2023  {
     2024    for (j=0; j<nRing; j++)
     2025    {
     2026      mpz_set_si(vec[j], (*diff_weight)[j]);
     2027    }
     2028    goto TEST_OVERFLOW;
     2029  }
     2030
    16842031#ifdef  NEXT_VECTORS_CC
    16852032  Print("\n// gcd of the weight vectors (current and target) = %d", gcd_tmp);
     
    16912038#endif
    16922039
    1693   mpz_t ddf; mpz_init(ddf);
    1694   mpz_t dcw; mpz_init(dcw);
    1695   BOOLEAN isdwpos;
     2040  //  BOOLEAN isdwpos;
    16962041
    16972042  // construct a new weight vector
     
    17022047
    17032048    if( (*diff_weight)[j]>0)
     2049    {
    17042050      mpz_mul_ui(s_zaehler, t_zaehler, (*diff_weight)[j]);
     2051    }
    17052052    else
    17062053    {
     
    17082055      mpz_neg(s_zaehler, s_zaehler);
    17092056    }
    1710 
    17112057    mpz_add(sntz, s_nenner, s_zaehler);
    1712 
    17132058    mpz_init_set(vec[j], sntz);
    1714 
     2059   
    17152060#ifdef NEXT_VECTORS_CC
    17162061    Print("\n//   j = %d ==> ", j);
     
    17282073
    17292074    if(j==0)
     2075    {
    17302076      mpz_set(ggt, sntz);
     2077    }
    17312078    else
     2079    {
    17322080      if(mpz_cmp_si(ggt,1) != 0)
     2081      {
    17332082        mpz_gcd(ggt, ggt, sntz);
    1734 
     2083      }
     2084    }
    17352085  }
    17362086
     
    17402090#endif
    17412091
    1742   mpz_t omega;
    1743   mpz_t sing_int;
    1744   mpz_init_set_ui(sing_int,  2147483647);
    1745 
    1746   /* construct a new weight vector and check whether vec[j] is overflow!!
    1747      i.e. vec[j] > 2^31.
    1748      If vec[j] doesn't overflow, define a weight vector
    1749      otherwise, report that overflow appears.
    1750      In the second case test whether the defined new vector correct is
    1751      plays an important rolle */
    1752 
     2092/**********************************************************************
     2093* construct a new weight vector and check whether vec[j] is overflow, *
     2094* i.e. vec[j] > 2^31.                                                 *
     2095* If vec[j] doesn't overflow, define a weight vector. Otherwise,      *
     2096* report that overflow appears. In the second case, test whether the  *
     2097* the correctness of the new vector plays an important role           *
     2098**********************************************************************/
     2099  kkk=0;
     2100  for(j=0; j<nRing; j++)
     2101    {
     2102    if(mpz_cmp(vec[j], sing_int_half) >= 0)
     2103      {
     2104      goto REDUCTION;
     2105      }
     2106    }
     2107  checkRed = 1;
    17532108  for (j=0; j<nRing; j++)
    1754   {
    1755     if(mpz_cmp_si(ggt,1)==0)
     2109    {
    17562110      (*diff_weight)[j] = mpz_get_si(vec[j]);
     2111    }
     2112  goto SIMPLIFY_GCD;
     2113
     2114  REDUCTION:
     2115  for (j=0; j<nRing; j++)
     2116  {
     2117    (*diff_weight)[j] = mpz_get_si(vec[j]);
     2118  }
     2119  while(MivAbsMax(diff_weight) >= 5)
     2120  {
     2121    for (j=0; j<nRing; j++)
     2122    {
     2123      if(mpz_cmp_si(ggt,1)==0)
     2124      {     
     2125        (*diff_weight1)[j] = floor(0.1*(*diff_weight)[j] + 0.5);
     2126        // Print("\n//  vector[%d] = %d \n",j+1, (*diff_weight1)[j]);
     2127      }
     2128      else
     2129      {
     2130        mpz_divexact(vec[j], vec[j], ggt);
     2131        (*diff_weight1)[j] = floor(0.1*(*diff_weight)[j] + 0.5);
     2132        // Print("\n//  vector[%d] = %d \n",j+1, (*diff_weight1)[j]);
     2133      }
     2134/*
     2135      if((*diff_weight1)[j] == 0)
     2136      {
     2137        kkk = kkk + 1;
     2138      }
     2139*/
     2140    }
     2141
     2142
     2143/*
     2144    if(kkk > nRing - 1)
     2145    {
     2146      // diff_weight was reduced to zero
     2147      // Print("\n // MwalkNextWeightCC: geaenderter Vector gleich Null! \n");
     2148      goto TEST_OVERFLOW;
     2149    }
     2150*/
     2151 
     2152    if(test_w_in_ConeCC(G,diff_weight1) != 0)
     2153    {
     2154      Print("\n// MwalkNextWeightCC: geaenderter vector liegt in Groebnerkegel! \n");
     2155      for (j=0; j<nRing; j++)
     2156      {
     2157        (*diff_weight)[j] = (*diff_weight1)[j];
     2158      }
     2159      if(MivAbsMax(diff_weight) < 5)
     2160      {
     2161        checkRed = 1;
     2162        goto SIMPLIFY_GCD;
     2163      }
     2164    }
    17572165    else
    17582166    {
    1759       mpz_divexact(vec[j], vec[j], ggt);
    1760       (*diff_weight)[j] = mpz_get_si(vec[j]);
    1761     }
    1762 
     2167      // Print("\n// MwalkNextWeightCC: geaenderter vector liegt nicht in Groebnerkegel! \n");
     2168      break;
     2169    }
     2170  }
     2171
     2172 TEST_OVERFLOW:
     2173
     2174  for (j=0; j<nRing; j++)
     2175  {
    17632176    if(mpz_cmp(vec[j], sing_int)>=0)
     2177    {
    17642178      if(Overflow_Error == FALSE)
    17652179      {
    17662180        Overflow_Error = TRUE;
    1767 
    1768         PrintS("\n// ** OVERFLOW in \"NextVector\": ");
     2181        PrintS("\n// ** OVERFLOW in \"MwalkNextWeightCC\": ");
    17692182        mpz_out_str( stdout, 10, vec[j]);
    17702183        PrintS(" is greater than 2147483647 (max. integer representation)");
    1771         Print("\n//  So vector[%d] := %d is wrong!!\n",j+1, (*diff_weight)[j]);
    1772       }
     2184        Print("\n//  So vector[%d] := %d is wrong!!\n",j+1, vec[j]);
     2185      }
     2186    }
    17732187  }
    17742188
    17752189 FINISH:
    1776 
    1777    mpz_clear(ggt);
     2190   delete diff_weight1;
    17782191   mpz_clear(t_zaehler);
    17792192   mpz_clear(t_nenner);
     2193   mpz_clear(s_zaehler);
     2194   mpz_clear(s_nenner);
    17802195   mpz_clear(sntz);
    17812196   mpz_clear(sztn);
     
    17842199   mpz_clear(deg_w0_p1);
    17852200   mpz_clear(deg_d0_p1);
     2201   mpz_clear(ggt);
    17862202   omFree(vec);
    1787 
     2203   mpz_clear(sing_int_half);
     2204   mpz_clear(sing_int);
     2205   mpz_clear(dcw);
     2206   mpz_clear(t_null);
     2207
     2208
     2209 
    17882210  if(Overflow_Error == FALSE)
     2211  {
    17892212    Overflow_Error = nError;
    1790 
    1791    return diff_weight;
    1792 }
    1793 
    1794 /*
    1795    compute an intermediate weight vector from iva to ivb w.r.t.
    1796    the reduced Groebner basis G.
    1797    Return NULL, if it is equal to iva or iva = avb.
    1798 */
     2213  }
     2214 rComplete(currRing);
     2215  for(kkk=0; kkk<IDELEMS(G);kkk++)
     2216  {
     2217    poly p=G->m[kkk];
     2218    while(p!=NULL)
     2219    {
     2220      p_Setm(p,currRing);
     2221      pIter(p);
     2222    }
     2223  }
     2224return diff_weight;
     2225}
     2226
     2227/**********************************************************************
     2228* Compute an intermediate weight vector from iva to ivb w.r.t.        *
     2229* the reduced Groebner basis G.                                       *
     2230* Return NULL, if it is equal to iva or iva = avb.                    *
     2231**********************************************************************/
    17992232intvec* MkInterRedNextWeight(intvec* iva, intvec* ivb, ideal G)
    18002233{
     
    18032236
    18042237  if(G == NULL)
     2238  {
    18052239    return tmp;
    1806 
     2240  }
    18072241  if(MivComp(iva, ivb) == 1)
     2242  {
    18082243    return tmp;
    1809 
     2244  }
    18102245  result = MwalkNextWeightCC(iva, ivb, G);
    18112246
     
    18202255}
    18212256
    1822 /* 01.11.01 */
    1823 /* define and execute a new ring which order is (a(va),lp,C) */
    1824 static void VMrDefault(intvec* va)
     2257/**************************************************************
     2258 * define and execute a new ring which order is (a(vb),a(va),lp,C)  *
     2259 * ************************************************************/
     2260static void VMrHomogeneous(intvec* va, intvec* vb)
    18252261{
    18262262
    18272263  if ((currRing->ppNoether)!=NULL)
     2264  {
    18282265    pDelete(&(currRing->ppNoether));
    1829 
     2266  }
    18302267  if (((sLastPrinted.rtyp>BEGIN_RING) && (sLastPrinted.rtyp<END_RING)) ||
    18312268      ((sLastPrinted.rtyp==LIST_CMD)&&(lRingDependend((lists)sLastPrinted.data))))
    1832 
    18332269  {
    18342270    sLastPrinted.CleanUp();
     
    18422278  int nb = 4;
    18432279
    1844   /*names*/
    1845   char* Q; //30.10.01 to avoid the corrupted memory, NOT change!!
     2280
     2281  //names
     2282  char* Q; // In order to avoid the corrupted memory, do not change.
     2283  r->names = (char **) omAlloc0(nv * sizeof(char_ptr));
     2284  for(i=0; i<nv; i++)
     2285  {
     2286    Q = currRing->names[i];
     2287    r->names[i]  = omStrDup(Q);
     2288  }
     2289
     2290  //weights: entries for 3 blocks: NULL Made:???
     2291  r->wvhdl = (int **)omAlloc0(nb * sizeof(int_ptr));
     2292  r->wvhdl[0] = (int*) omAlloc(nv*sizeof(int));
     2293  r->wvhdl[1] = (int*) omAlloc((nv-1)*sizeof(int));
     2294
     2295  for(i=0; i<nv-1; i++)
     2296  {
     2297    r->wvhdl[1][i] = (*vb)[i];
     2298    r->wvhdl[0][i] = (*va)[i];
     2299  }
     2300  r->wvhdl[0][nv] = (*va)[nv];
     2301
     2302  // order: (1..1),a,lp,C
     2303  r->order = (int *) omAlloc(nb * sizeof(int *));
     2304  r->block0 = (int *)omAlloc0(nb * sizeof(int *));
     2305  r->block1 = (int *)omAlloc0(nb * sizeof(int *));
     2306
     2307  // ringorder a for the first block: var 1..nv
     2308  r->order[0]  = ringorder_a;
     2309  r->block0[0] = 1;
     2310  r->block1[0] = nv;
     2311
     2312 // ringorder a for the second block: var 2..nv
     2313  r->order[1]  = ringorder_a;
     2314  r->block0[1] = 2;
     2315  r->block1[1] = nv;
     2316
     2317  // ringorder lp for the third block: var 2..nv
     2318  r->order[2]  = ringorder_lp;
     2319  r->block0[2] = 2;
     2320  r->block1[2] = nv;
     2321
     2322  // ringorder C for the 4th block
     2323  // it is very important within "idLift",
     2324  // especially, by ring syz_ring=rCurrRingAssure_SyzComp();
     2325  // therefore, nb  must be (nBlocks(currRing)  + 1)
     2326  r->order[3]  = ringorder_C;
     2327
     2328  // polynomial ring
     2329  r->OrdSgn    = 1;
     2330
     2331  // complete ring intializations
     2332 
     2333  rComplete(r);
     2334
     2335  rChangeCurrRing(r);
     2336}
     2337
     2338
     2339/**************************************************************
     2340 * define and execute a new ring which order is (a(va),lp,C)  *
     2341 * ************************************************************/
     2342static void VMrDefault(intvec* va)
     2343{
     2344
     2345  if ((currRing->ppNoether)!=NULL)
     2346  {
     2347    pDelete(&(currRing->ppNoether));
     2348  }
     2349  if (((sLastPrinted.rtyp>BEGIN_RING) && (sLastPrinted.rtyp<END_RING)) ||
     2350      ((sLastPrinted.rtyp==LIST_CMD)&&(lRingDependend((lists)sLastPrinted.data))))
     2351  {
     2352    sLastPrinted.CleanUp();
     2353  }
     2354
     2355  ring r = (ring) omAlloc0Bin(sip_sring_bin);
     2356  int i, nv = currRing->N;
     2357
     2358  r->cf  = currRing->cf;
     2359  r->N   = currRing->N;
     2360
     2361  int nb = 4;
     2362
     2363  //names
     2364  char* Q; // In order to avoid the corrupted memory, do not change.
    18462365  r->names = (char **) omAlloc0(nv * sizeof(char_ptr));
    18472366  for(i=0; i<nv; i++)
     
    18622381  r->block1 = (int *)omAlloc0(nb * sizeof(int *));
    18632382
    1864   /* ringorder a for the first block: var 1..nv */
     2383  // ringorder a for the first block: var 1..nv
    18652384  r->order[0]  = ringorder_a;
    18662385  r->block0[0] = 1;
    18672386  r->block1[0] = nv;
    18682387
    1869   /* ringorder lp for the second block: var 1..nv */
     2388  // ringorder lp for the second block: var 1..nv
    18702389  r->order[1]  = ringorder_lp;
    18712390  r->block0[1] = 1;
    18722391  r->block1[1] = nv;
    18732392
    1874   /* ringorder C for the third block */
     2393  // ringorder C for the third block
    18752394  // it is very important within "idLift",
    18762395  // especially, by ring syz_ring=rCurrRingAssure_SyzComp();
     
    18782397  r->order[2]  = ringorder_C;
    18792398
    1880   /* the last block: everything is 0 */
     2399  // the last block: everything is 0
    18812400  r->order[3]  = 0;
    18822401
    1883   /*polynomial ring*/
     2402  // polynomial ring
    18842403  r->OrdSgn    = 1;
    18852404
    1886   /* complete ring intializations */
     2405  // complete ring intializations
     2406 
    18872407  rComplete(r);
    18882408
     
    18902410}
    18912411
    1892 /* 03.11.01 */
    1893 /* define and execute a new ring which order is  a lexicographic order */
     2412/**********************************************************************
     2413* define and execute a new ring which order is  a lexicographic order *
     2414***********************************************************************/
    18942415static void VMrDefaultlp(void)
    18952416{
    18962417
    18972418  if ((currRing->ppNoether)!=NULL)
     2419  {
    18982420    pDelete(&(currRing->ppNoether));
    1899 
     2421  }
    19002422  if (((sLastPrinted.rtyp>BEGIN_RING) && (sLastPrinted.rtyp<END_RING)) ||
    19012423      ((sLastPrinted.rtyp==LIST_CMD)&&(lRingDependend((lists)sLastPrinted.data))))
     
    19102432  r->cf  = currRing->cf;
    19112433  r->N   = currRing->N;
    1912   int nb = rBlocks(currRing) + 1;//31.10.01 (+1)
    1913 
    1914   /*names*/
    1915   char* Q; //30.10.01 to avoid the corrupted memory, NOT change!!
     2434  int nb = rBlocks(currRing) + 1;
     2435
     2436  // names
     2437  char* Q; // to avoid the corrupted memory, do not change!!
    19162438  r->names = (char **) omAlloc0(nv * sizeof(char_ptr));
    19172439  for(i=0; i<nv; i++)
     
    19452467
    19462468  /* complete ring intializations */
     2469 
    19472470  rComplete(r);
    19482471
     
    19672490  res->cf = currRing->cf; currRing->cf->ref++;
    19682491
    1969 //   if (currRing->cf->extRing!=NULL)
    1970 //     currRing->cf->extRing->ref++;
    1971 //
    1972 //   if (rParameter (currRing)!=NULL)
    1973 //   {
    1974 //     res->cf->extRing->qideal->m[0]=p_Copy(currRing->cf->extRing->qideal->m[0],currRing->cf->extRing);
    1975 //     int l=rPar(currRing);
    1976 //
    1977 //     res->cf->extRing->names=(char **)omAlloc(l*sizeof(char_ptr));
    1978 //
    1979 //     for(i=l-1;i>=0;i--)
    1980 //       rParameter (res)[i]=omStrDup(rParameter (currRing)[i]);
    1981 //   }
    19822492
    19832493  /*weights: entries for 3 blocks: NULL Made:???*/
     
    19882498
    19892499  /* order: a,lp,C,0 */
     2500
    19902501  res->order = (int *) omAlloc(nb * sizeof(int *));
    19912502  res->block0 = (int *)omAlloc0(nb * sizeof(int *));
    19922503  res->block1 = (int *)omAlloc0(nb * sizeof(int *));
    19932504
    1994   /* ringorder a for the first block: var 1..nv */
     2505  // ringorder a for the first block: var 1..nv
    19952506  res->order[0]  = ringorder_a;
    19962507  res->block0[0] = 1;
    19972508  res->block1[0] = nv;
    19982509
    1999   /* ringorder lp for the second block: var 1..nv */
     2510  // ringorder lp for the second block: var 1..nv
    20002511  res->order[1]  = ringorder_lp;
    20012512  res->block0[1] = 1;
    20022513  res->block1[1] = nv;
    20032514
    2004   /* ringorder C for the third block */
     2515  // ringorder C for the third block
    20052516  // it is very important within "idLift",
    20062517  // especially, by ring syz_ring=rCurrRingAssure_SyzComp();
     
    20082519  res->order[2]  = ringorder_C;
    20092520
    2010   /* the last block: everything is 0 */
     2521  // the last block: everything is 0
    20112522  res->order[3]  = 0;
    20122523
    2013   /*polynomial ring*/
     2524  // polynomial ring
    20142525  res->OrdSgn    = 1;
    20152526
     
    20172528  res->names   = (char **)omAlloc0(nv * sizeof(char_ptr));
    20182529  for (i=nv-1; i>=0; i--)
     2530  {
    20192531    res->names[i] = omStrDup(currRing->names[i]);
    2020 
    2021   /* complete ring intializations */
    2022    rComplete(res);
    2023 
     2532  }
     2533  // complete ring intializations
     2534  rComplete(res);
    20242535
    20252536  // clean up history
     
    20302541
    20312542
    2032   /* execute the created ring */
     2543  // execute the created ring
    20332544  rChangeCurrRing(res);
    20342545}
     
    20482559  r->cf = currRing->cf; currRing->cf->ref++;
    20492560
    2050 //   if (currRing->cf->extRing!=NULL)
    2051 //     currRing->cf->extRing->ref++;
    2052 //
    2053 //   if (rParameter (currRing)!=NULL)
    2054 //   {
    2055 //     r->cf->extRing->qideal->m[0]=p_Copy(currRing->cf->extRing->qideal->m[0], currRing->cf->extRing);
    2056 //     int l=rPar(currRing);
    2057 //     r->cf->extRing->names=(char **)omAlloc(l*sizeof(char_ptr));
    2058 //
    2059 //     for(i=l-1;i>=0;i--)
    2060 //       rParameter(r)[i]=omStrDup(rParameter (currRing)[i]);
    2061 //   }
    2062 
    20632561
    20642562  r->cf  = currRing->cf;
    20652563  r->N   = currRing->N;
    2066   int nb = rBlocks(currRing) + 1;//31.10.01 (+1)
    2067 
    2068   /*names*/
     2564  int nb = rBlocks(currRing) + 1;
     2565
     2566  // names
    20692567  char* Q;
    20702568  r->names = (char **) omAlloc0(nv * sizeof(char_ptr));
     
    21062604//
    21072605//     for(i=l-1;i>=0;i--)
     2606//     {
    21082607//       rParameter(r)[i]=omStrDup(rParameter(currRing)[i]);
     2608//     }
    21092609//   }
    21102610
    2111   /* complete ring intializations */
     2611  // complete ring intializations
     2612 
    21122613  rComplete(r);
    21132614
     
    21182619  }
    21192620
    2120   /* execute the created ring */
     2621  // execute the created ring
    21212622  rChangeCurrRing(r);
    21222623}
    21232624
    2124 #if 0 /*unused*/
    2125 /* check wheather one or more components of a vector are zero */
     2625//unused
     2626/**************************************************************
     2627 * check wheather one or more components of a vector are zero *
     2628 **************************************************************/
     2629//#if 0
    21262630static int isNolVector(intvec* hilb)
    21272631{
    21282632  int i;
    21292633  for(i=hilb->length()-1; i>=0; i--)
     2634  {
    21302635    if((* hilb)[i]==0)
     2636    {
    21312637      return 1;
    2132 
     2638    }
     2639  }
    21332640  return 0;
    21342641}
    2135 #endif
    2136 
     2642//#endif
    21372643
    21382644/******************************  Februar 2002  ****************************
     
    21412647 * its perturbation degree is tp_deg                                      *
    21422648 * We call the following subfunction LastGB, if                           *
    2143  *     the computed intermediate weight vector or                         *
    2144  *     the perturbed target weight vector                                 *
    2145  * does NOT in the correct cone.                                          *
     2649 * the computed intermediate weight vector or                             *
     2650 * if the perturbed target weight vector does NOT lie n the correct cone  *
    21462651 **************************************************************************/
    21472652
     
    21622667  intvec* ivNull = new intvec(nV);
    21632668  intvec* extra_curr_weight = new intvec(nV);
     2669  intvec* next_weight;
     2670
     2671#ifndef  BUCHBERGER_ALG
    21642672  intvec* hilb_func;
    2165   intvec* next_weight;
    2166 
    2167   /* to avoid (1,0,...,0) as the target vector */
     2673#endif
     2674
     2675  // to avoid (1,0,...,0) as the target vector
    21682676  intvec* last_omega = new intvec(nV);
    21692677  for(i=nV-1; i>0; i--)
     2678  {
    21702679    (*last_omega)[i] = 1;
     2680  }
    21712681  (*last_omega)[0] = 10000;
    21722682
    21732683  ring EXXRing = currRing;
    21742684
    2175   /* compute a pertubed weight vector of the target weight vector */
     2685  // compute a pertubed weight vector of the target weight vector
    21762686  if(tp_deg > 1 && tp_deg <= nV)
    21772687  {
    21782688    //..25.03.03 VMrDefaultlp();//    VMrDefault(target_weight);
    21792689    if (rParameter (currRing) != NULL)
     2690    {
    21802691      DefRingParlp();
     2692    }
    21812693    else
     2694    {
    21822695      VMrDefaultlp();
    2183 
     2696    }
    21842697    TargetRing = currRing;
    21852698    ssG = idrMoveR(G,EXXRing,currRing);
     
    21942707  }
    21952708  else
     2709  {
    21962710    target_weight = Mivlp(nV);
    2197 
     2711  }
    21982712  //Print("\n// ring r%d_%d = %s;\n", tp_deg, nwalk, rString(currRing));
    21992713
     
    22032717    nstep++;
    22042718    to=clock();
    2205     /* compute a next weight vector */
     2719   // compute a next weight vector
    22062720    next_weight = MkInterRedNextWeight(curr_weight,target_weight, G);
    22072721    xtnw=xtnw+clock()-to;
     2722
    22082723#ifdef PRINT_VECTORS
    22092724    MivString(curr_weight, target_weight, next_weight);
    22102725#endif
    22112726
    2212     if(Overflow_Error == TRUE){
     2727    if(Overflow_Error == TRUE)
     2728    {
    22132729      newRing = currRing;
    22142730      nnwinC = 0;
    22152731      if(tp_deg == 1)
     2732      {
    22162733        nlast = 1;
     2734      }
    22172735      delete next_weight;
    22182736
     
    22232741    }
    22242742
    2225     if(MivComp(next_weight, ivNull) == 1){
     2743    if(MivComp(next_weight, ivNull) == 1)
     2744    {
    22262745      //Print("\n// ring r%d_%d = %s;\n", tp_deg, nwalk, rString(currRing));
    22272746      newRing = currRing;
     
    22342753
    22352754    for(i=nV-1; i>=0; i--)
    2236         (*extra_curr_weight)[i] = (*curr_weight)[i];
    2237 
     2755    {
     2756      (*extra_curr_weight)[i] = (*curr_weight)[i];
     2757    }
    22382758    /* 06.11.01 NOT Changed */
    22392759    for(i=nV-1; i>=0; i--)
     2760    {
    22402761      (*curr_weight)[i] = (*next_weight)[i];
    2241 
     2762    }
    22422763    oldRing = currRing;
    22432764    to=clock();
    2244     /* compute an initial form ideal of <G> w.r.t. "curr_vector" */
     2765    // compute an initial form ideal of <G> w.r.t. "curr_vector"
    22452766    Gomega = MwalkInitialForm(G, curr_weight);
    22462767    xtif=xtif+clock()-to;
    22472768
    22482769#ifdef ENDWALKS
    2249     if(endwalks == 1){
     2770    if(endwalks == 1)
     2771    {
    22502772      Print("\n// ring r%d_%d = %s;\n", tp_deg, nwalk, rString(currRing));
    22512773      idElements(Gomega, "Gw");
     
    22562778#ifndef  BUCHBERGER_ALG
    22572779    if(isNolVector(curr_weight) == 0)
     2780    {
    22582781      hilb_func = hFirstSeries(Gomega,NULL,NULL,curr_weight,currRing);
     2782    }
    22592783    else
     2784    {
    22602785      hilb_func = hFirstSeries(Gomega,NULL,NULL,last_omega,currRing);
     2786    }
    22612787#endif // BUCHBERGER_ALG
    22622788
     
    22642790    //..25.03.03 VMrDefault(curr_weight);
    22652791    if (rParameter (currRing) != NULL)
     2792    {
    22662793      DefRingPar(curr_weight);
     2794    }
    22672795    else
     2796    {
    22682797      VMrDefault(curr_weight);
    2269 
     2798    }
    22702799    newRing = currRing;
    22712800    Gomega1 = idrMoveR(Gomega, oldRing,currRing);
     
    23032832    idDelete(&F1);
    23042833
    2305     if(endwalks == 1){
     2834    if(endwalks == 1)
     2835    {
    23062836      //Print("\n// ring r%d_%d = %s;\n", tp_deg, nwalk, rString(currRing));
    23072837      break;
     
    23172847    //..25.03.03 VMrDefaultlp();//define and execute the ring "lp"
    23182848    if (rParameter (currRing) != NULL)
     2849    {
    23192850      DefRingParlp();
     2851    }
    23202852    else
     2853    {
    23212854      VMrDefaultlp();
    2322 
     2855    }
    23232856    F1 = idrMoveR(G, newRing,currRing);
    23242857
     
    23582891        //..25.03.03 VMrDefaultlp();//define and execute the ring "lp"
    23592892        if (rParameter (currRing) != NULL)
     2893        {
    23602894          DefRingParlp();
     2895        }
    23612896        else
     2897        {
    23622898          VMrDefaultlp();
    2363 
     2899        }
    23642900
    23652901      F1 = idrMoveR(G, newRing,currRing);
     
    23792915
    23802916  if(Overflow_Error == FALSE)
     2917  {
    23812918    Overflow_Error = nError;
    2382 
     2919  }
    23832920  return(result);
    23842921}
    23852922
    2386 
    2387 /* check whether a polynomial of G has least 3 monomials */
     2923/**********************************************************
     2924 * check whether a polynomial of G has least 3 monomials  *
     2925 **********************************************************/
    23882926static int lengthpoly(ideal G)
    23892927{
    23902928  int i;
    23912929  for(i=IDELEMS(G)-1; i>=0; i--)
     2930  {
    23922931#if 0
    23932932    if(pLength(G->m[i])>2)
     2933    {
    23942934      return 1;
     2935    }
    23952936#else
    23962937    if((G->m[i]!=NULL) /* len >=0 */
     
    23992940       && (G->m[i]->next->next->next!=NULL) /* len >=3 */
    24002941      //&& (G->m[i]->next->next->next->next!=NULL) /* len >=4 */
    2401        ) return 1;
    2402 #endif
     2942       )
     2943    {
     2944    return 1;
     2945    }
     2946#endif
     2947  }
    24032948  return 0;
    24042949}
    24052950
    2406 /* check whether a polynomial of G has least 2 monomials */
     2951/*********************************************************
     2952 * check whether a polynomial of G has least 2 monomials *
     2953**********************************************************/
    24072954static int islengthpoly2(ideal G)
    24082955{
    24092956  int i;
    24102957  for(i=IDELEMS(G)-1; i>=0; i--)
     2958  {
    24112959    if((G->m[i]!=NULL) /* len >=0 */
    24122960       && (G->m[i]->next!=NULL) /* len >=1 */
    24132961       && (G->m[i]->next->next!=NULL)) /* len >=2 */
     2962    {
    24142963      return 1;
    2415 
     2964    }
     2965  }
    24162966  return 0;
    24172967}
     
    24352985 */
    24362986
    2437 /*2
    2438 * return the initial term of an ideal
    2439 */
     2987/***************************************
     2988 * return the initial term of an ideal *
     2989 ***************************************/
    24402990static ideal idHeadCC(ideal h)
    24412991{
     
    24472997  {
    24482998    if (h->m[i]!=NULL)
     2999    {
    24493000      m->m[i]=pHead(h->m[i]);
     3001    }
    24503002  }
    24513003  return m;
    24523004}
    24533005
    2454 /* check whether two head-ideals are the same */
     3006/**********************************************
     3007 * check whether two head-ideals are the same *
     3008 **********************************************/
    24553009static inline int test_G_GB_walk(ideal H0, ideal H1)
    24563010{
     
    24583012
    24593013  if(nG != IDELEMS(H1))
     3014  {
    24603015    return 0;
    2461 
     3016  }
    24623017  for(i=nG-1; i>=0; i--)
    24633018  {
     
    24723027#else
    24733028    if(!pEqualPolys(H0->m[i],H1->m[i]))
     3029    {
    24743030      return 0;
    2475 #endif
    2476   }
    2477 
     3031    }
     3032#endif
     3033  }
    24783034  return 1;
    24793035}
    24803036
    2481 #if 0 /*unused*/
    2482 /* 19.11.01 */
    2483 /* find the maximal total degree of polynomials in G */
     3037//unused
     3038/*****************************************************
     3039 * find the maximal total degree of polynomials in G *
     3040 *****************************************************/
     3041#if 0
    24843042static int Trandegreebound(ideal G)
    24853043{
    24863044  int i, nG = IDELEMS(G);
    2487   int np=1, nV = currRing->N;
     3045  // int np=1;
     3046  int nV = currRing->N;
    24883047  int degtmp, result = 0;
    24893048  intvec* ivUnit = Mivdp(nV);
     
    24913050  for(i=nG-1; i>=0; i--)
    24923051  {
    2493     /* find the maximal total degree of the polynomial G[i] */
     3052    // find the maximal total degree of the polynomial G[i]
    24943053    degtmp = MwalkWeightDegree(G->m[i], ivUnit);
    24953054    if(degtmp > result)
     3055    {
    24963056      result = degtmp;
     3057    }
    24973058  }
    24983059  delete ivUnit;
     
    25013062#endif
    25023063
    2503 /* perturb the weight vector iva w.r.t. the ideal G.
    2504    the monomial order of the current ring is the w_1 weight lex. order.
    2505    define w := d^(n-1)w_1+ d^(n-2)w_2, ...+ dw_(n-1)+ w_n
    2506    where d := 1 + max{totdeg(g):g in G}*m, or
    2507    d := (2*maxdeg*maxdeg + (nV+1)*maxdeg)*m;
    2508 */
    2509 
    2510 #if 0 /*unused*/
    2511 //GMP
     3064//unused
     3065/************************************************************************
     3066 * perturb the weight vector iva w.r.t. the ideal G.                    *
     3067 *  the monomial order of the current ring is the w_1 weight lex. order *
     3068 *  define w := d^(n-1)w_1+ d^(n-2)w_2, ...+ dw_(n-1)+ w_n              *
     3069 *  where d := 1 + max{totdeg(g):g in G}*m, or                          *
     3070 *  d := (2*maxdeg*maxdeg + (nV+1)*maxdeg)*m;                           *
     3071 ************************************************************************/
     3072#if 0
    25123073static intvec* TranPertVector(ideal G, intvec* iva)
    25133074{
     
    25153076  Overflow_Error = FALSE;
    25163077
    2517   int i, j,  nG = IDELEMS(G);
     3078  int i, j;
     3079  // int nG = IDELEMS(G);
    25183080  int nV = currRing->N;
    25193081
     
    25273089  {
    25283090    mtmp = (*ivMat)[i];
    2529 
    2530     if(mtmp <0) mtmp = -mtmp;
    2531 
     3091    if(mtmp <0)
     3092    {
     3093      mtmp = -mtmp;
     3094    }
    25323095    if(mtmp > m)
     3096    {
    25333097      m = mtmp;
    2534   }
    2535 
    2536   /* define the maximal total degree of polynomials of G */
     3098    }
     3099  }
     3100
     3101  // define the maximal total degree of polynomials of G
    25373102  mpz_t ndeg;
    25383103  mpz_init(ndeg);
     
    25553120  mpz_mul_ui(ndeg, ndeg, m);
    25563121
     3122  mpz_clear(ztmp);
     3123
    25573124  //PrintS("\n// with the new upper degree bound (2d^2+(n+1)d)*m ");
    25583125  //Print("\n//         where d = %d, n = %d and bound = %d", maxdeg, nV, ndeg);
    25593126#endif //UPPER_BOUND
    25603127
    2561   /* 29.08.03*/
    25623128#ifdef INVEPS_SMALL_IN_TRAN
    2563  if(mpz_cmp_ui(ndeg, nV)>0 && nV > 3)
     3129  if(mpz_cmp_ui(ndeg, nV)>0 && nV > 3)
     3130  {
    25643131    mpz_cdiv_q_ui(ndeg, ndeg, nV);
    2565 
     3132  }
    25663133 //PrintS("\n// choose the \"small\" inverse epsilon:");
    25673134 //mpz_out_str(stdout, 10, ndeg);
     
    25813148  mpz_t *ivtmp=(mpz_t *)omAlloc(nV*sizeof(mpz_t));
    25823149  for(i=0; i<nV; i++)
     3150  {
    25833151    mpz_init(ivtmp[i]);
    2584 
     3152  }
    25853153  mpz_t sing_int;
    25863154  mpz_init_set_ui(sing_int,  2147483647);
     
    25883156  intvec* repr_vector = new intvec(nV);
    25893157
    2590   /* define ivtmp := ndeg^(n-1).w_1 + ndeg^(n-2).w_2 + ... + w_n */
     3158  // define ivtmp := ndeg^(n-1).w_1 + ndeg^(n-2).w_2 + ... + w_n
    25913159  for(i=0; i<nV; i++)
     3160  {
    25923161    for(j=0; j<nV; j++)
    25933162    {
    25943163      if( (*ivMat)[i*nV+j] >= 0 )
     3164      {
    25953165        mpz_mul_ui(ivres[i], ivres[i], (*ivMat)[i*nV+j]);
     3166      }
    25963167      else
    25973168      {
     
    26013172      mpz_add(ivtmp[j], ivtmp[j], ivres[i]);
    26023173    }
    2603 
     3174  }
    26043175  delete ivMat;
    26053176
     
    26333204  omFree(ivres);
    26343205  omFree(ivtmp);
     3206
     3207  mpz_clear(sing_int);
     3208  mpz_clear(deg_tmp);
     3209  mpz_clear(ndeg);
     3210
    26353211  return repr_vector;
    26363212}
    26373213#endif
    26383214
    2639 
    2640 
    2641 #if 0 /*unused*/
     3215//unused
     3216#if 0
    26423217static intvec* TranPertVector_lp(ideal G)
    26433218{
    26443219  BOOLEAN nError = Overflow_Error;
    26453220  Overflow_Error = FALSE;
    2646 
    2647   int i, j,  nG = IDELEMS(G);
     3221  // int j, nG = IDELEMS(G);
     3222  int i;
    26483223  int nV = currRing->N;
    26493224
    2650   /* define the maximal total degree of polynomials of G */
     3225  // define the maximal total degree of polynomials of G
    26513226  mpz_t ndeg;
    26523227  mpz_init(ndeg);
     
    26673242  mpz_mul_ui(maxdeg, maxdeg, nV+1);
    26683243  mpz_add(ndeg, ztmp, maxdeg);
    2669   /*
    2670     PrintS("\n// with the new upper degree bound (2d^2+(n+1)d)*m ");
    2671     Print("\n//         where d = %d, n = %d and bound = %d",
    2672     mpz_get_si(maxdeg), nV, mpz_get_si(ndeg));
    2673   */
    2674 #endif //UPPER_BOUND
     3244  // PrintS("\n// with the new upper degree bound (2d^2+(n+1)d)*m ");
     3245  // Print("\n//         where d = %d, n = %d and bound = %d",
     3246  // mpz_get_si(maxdeg), nV, mpz_get_si(ndeg));
     3247
     3248 mpz_clear(ztmp);
     3249
     3250#endif
    26753251
    26763252#ifdef INVEPS_SMALL_IN_TRAN
     
    27253301
    27263302  omFree(ivres);
     3303
     3304  mpz_clear(ndeg);
     3305  mpz_clear(sing_int);
     3306
    27273307  return repr_vector;
    27283308}
    27293309#endif
    27303310
    2731 
    2732 #if 0 /*unused*/
    2733 //GMP
     3311//unused
     3312#if 0
    27343313static intvec* RepresentationMatrix_Dp(ideal G, intvec* M)
    27353314{
     
    27453324  for(i=IDELEMS(G)-1; i>=0; i--)
    27463325  {
    2747     /* find the maximal total degree of the polynomial G[i] */
     3326    // find the maximal total degree of the polynomial G[i]
    27483327    degtmp = MwalkWeightDegree(G->m[i], ivUnit);
    27493328    if(degtmp > maxdeg)
     
    27663345    mpz_init(ivtmp[i]);
    27673346
    2768   /* define ivtmp := ndeg^(n-1).w_1 + ndeg^(n-2).w_2 + ... + w_n */
     3347  // define ivtmp := ndeg^(n-1).w_1 + ndeg^(n-2).w_2 + ... + w_n
    27693348  for(i=0; i<nV; i++)
    27703349    for(j=0; j<nV; j++)
     
    27803359      mpz_add(ivtmp[j], ivtmp[j], ztmp);
    27813360    }
    2782 
     3361  delete ivres;
    27833362  mpz_t sing_int;
    27843363  mpz_init_set_ui(sing_int,  2147483647);
     
    28113390    Overflow_Error = nError;
    28123391
     3392  mpz_clear(sing_int);
    28133393  mpz_clear(ztmp);
    28143394  omFree(ivtmp);
     
    28183398#endif
    28193399
    2820 
    2821 
    2822 
    2823 
    2824 /* The following subroutine is the implementation of our first improved
    2825    Groebner walk algorithm, i.e. the first altervative algorithm.
    2826    First we use the Grobner walk algorithm and then we call the changed
    2827    perturbation walk algorithm with decreased degree, if an intermediate
    2828    weight vector is equal to the current target weight vector.
    2829    This call will be only repeated until we get the wanted reduced Groebner
    2830    basis or n times, where n is the numbers of variables.
    2831 */
    2832 
    2833 // 19 Juni 2003
    2834 #if 0 /* unused*/
     3400/*****************************************************************************
     3401 * The following subroutine is the implementation of our first improved      *
     3402 * Groebner walk algorithm, i.e. the first altervative algorithm.            *
     3403 * First we use the Grobner walk algorithm and then we call the changed      *
     3404 * perturbation walk algorithm with decreased degree, if an intermediate     *
     3405 * weight vector is equal to the current target weight vector.               *
     3406 * This call will be only repeated until we get the wanted reduced Groebner  *
     3407 * basis or n times, where n is the numbers of variables.                    *
     3408 *****************************************************************************/
     3409
     3410//unused
     3411#if 0
    28353412static int testnegintvec(intvec* v)
    28363413{
    28373414  int n = v->length();
    28383415  int i;
    2839   for(i=0; i<n; i++){
     3416  for(i=0; i<n; i++)
     3417  {
    28403418    if((*v)[i]<0)
     3419    {
    28413420      return(1);
     3421    }
    28423422  }
    28433423  return(0);
     
    28453425#endif
    28463426
    2847 
    2848 /* 7 Februar 2002 */
    2849 /* npwinc = 0, if curr_weight doesn't stay in the correct Groebner cone */
     3427// npwinc = 0, if curr_weight doesn't stay in the correct Groebner cone
    28503428static ideal Rec_LastGB(ideal G, intvec* curr_weight,
    28513429                        intvec* orig_target_weight, int tp_deg, int npwinc)
     
    28533431  BOOLEAN nError = Overflow_Error;
    28543432  Overflow_Error = FALSE;
    2855   BOOLEAN nOverflow_Error = FALSE;
     3433  // BOOLEAN nOverflow_Error = FALSE;
    28563434
    28573435  clock_t tproc=0;
     
    28673445  intvec* ivNull = new intvec(nV); //define (0,...,0)
    28683446  ring EXXRing = currRing;
    2869   int NEG=0; //19 juni 03
    2870   intvec* extra_curr_weight = new intvec(nV);
     3447  //int NEG=0; //19 juni 03
    28713448  intvec* next_weight;
    2872 
     3449#ifndef  BUCHBERGER_ALG
    28733450  //08 Juli 03
    28743451  intvec* hilb_func;
    2875   /* to avoid (1,0,...,0) as the target vector */
     3452#endif
     3453  // to avoid (1,0,...,0) as the target vector
    28763454  intvec* last_omega = new intvec(nV);
    28773455  for(i=nV-1; i>0; i--)
     
    28813459  BOOLEAN isGB = FALSE;
    28823460
    2883   /* compute a pertubed weight vector of the target weight vector */
    2884   if(tp_deg > 1 && tp_deg <= nV) {
     3461  // compute a pertubed weight vector of the target weight vector
     3462  if(tp_deg > 1 && tp_deg <= nV)
     3463  {
    28853464    ideal H0 = idHeadCC(G);
    28863465
    28873466    if (rParameter (currRing) != NULL)
     3467    {
    28883468      DefRingParlp();
     3469    }
    28893470    else
     3471    {
    28903472      VMrDefaultlp();
    2891 
     3473    }
    28923474    TargetRing = currRing;
    28933475    ssG = idrMoveR(G,EXXRing,currRing);
     
    28963478    ideal H1 = idHeadCC(ssG);
    28973479
    2898     /* Apply Lemma 2.2 in Collart et. al (1997) to check whether
    2899        cone(k-1) is equal to cone(k) */
    2900     if(test_G_GB_walk(H0_tmp,H1)==1) {
     3480    // Apply Lemma 2.2 in Collart et. al (1997) to check whether cone(k-1) is equal to cone(k)
     3481    if(test_G_GB_walk(H0_tmp,H1)==1)
     3482    {
    29013483      idDelete(&H0_tmp);
    29023484      idDelete(&H1);
     
    29073489
    29083490      if(npwinc != 0)
     3491      {
    29093492        goto LastGB_Finish;
    2910       else {
     3493      }
     3494      else
     3495      {
    29113496        isGB = TRUE;
    29123497        goto KSTD_Finish;
     
    29233508    G = idrMoveR(ssG, TargetRing,currRing);
    29243509
    2925     if(Overflow_Error == TRUE)  {
    2926       nOverflow_Error = Overflow_Error;
    2927       NEG = 1;
     3510    if(Overflow_Error == TRUE)
     3511    {
     3512      //nOverflow_Error = Overflow_Error;
     3513      //NEG = 1;
    29283514      newRing = currRing;
    29293515      goto JUNI_STD;
     
    29373523
    29383524    if(nwalk==1)
     3525    {
    29393526      goto FIRST_STEP;
    2940 
     3527    }
    29413528    to=clock();
    2942     /* compute an initial form ideal of <G> w.r.t. "curr_vector" */
     3529    // compute an initial form ideal of <G> w.r.t. "curr_vector"
    29433530    Gomega = MwalkInitialForm(G, curr_weight);
    29443531    xtif=xtif+clock()-to;
     
    29463533#ifndef  BUCHBERGER_ALG
    29473534    if(isNolVector(curr_weight) == 0)
     3535    {
    29483536      hilb_func = hFirstSeries(Gomega,NULL,NULL,curr_weight,currRing);
     3537    }
    29493538    else
     3539    {
    29503540      hilb_func = hFirstSeries(Gomega,NULL,NULL,last_omega,currRing);
     3541    }
    29513542#endif // BUCHBERGER_ALG
    29523543
    29533544    oldRing = currRing;
    29543545
    2955     /* defiNe a new ring that its ordering is "(a(curr_weight),lp) */
     3546    // defiNe a new ring that its ordering is "(a(curr_weight),lp)
    29563547    if (rParameter(currRing) != NULL)
     3548    {
    29573549      DefRingPar(curr_weight);
     3550    }
    29583551    else
     3552    {
    29593553      VMrDefault(curr_weight);
    2960 
     3554    }
    29613555    newRing = currRing;
    29623556    Gomega1 = idrMoveR(Gomega, oldRing,currRing);
    29633557    to=clock();
    2964     /* compute a reduced Groebner basis of <Gomega> w.r.t. "newRing" */
     3558    // compute a reduced Groebner basis of <Gomega> w.r.t. "newRing"
    29653559#ifdef  BUCHBERGER_ALG
    29663560    M = MstdhomCC(Gomega1);
     
    29703564#endif // BUCHBERGER_ALG
    29713565    xtstd=xtstd+clock()-to;
    2972     /* change the ring to oldRing */
     3566    // change the ring to oldRing
    29733567    rChangeCurrRing(oldRing);
    29743568    M1 =  idrMoveR(M, newRing,currRing);
     
    29763570
    29773571     to=clock();
    2978     /* compute a reduced Groebner basis of <G> w.r.t. "newRing" by the
    2979        lifting process*/
     3572    // compute a reduced Groebner basis of <G> w.r.t. "newRing" by the lifting process
    29803573    F = MLifttwoIdeal(Gomega2, M1, G);
    29813574    xtlift=xtlift+clock()-to;
     
    29843577    idDelete(&G);
    29853578
    2986     /* change the ring to newRing */
     3579    // change the ring to newRing
    29873580    rChangeCurrRing(newRing);
    29883581    F1 = idrMoveR(F, oldRing,currRing);
    29893582
    29903583    to=clock();
    2991     /* reduce the Groebner basis <G> w.r.t. new ring */
     3584    // reduce the Groebner basis <G> w.r.t. new ring
    29923585    G = kInterRedCC(F1, NULL);
    29933586    xtred=xtred+clock()-to;
     
    29953588
    29963589    if(endwalks == 1)
     3590    {
    29973591      break;
    2998 
     3592    }
    29993593  FIRST_STEP:
    30003594    to=clock();
    30013595    Overflow_Error = FALSE;
    3002     /* compute a next weight vector */
     3596    // compute a next weight vector
    30033597    next_weight = MkInterRedNextWeight(curr_weight,target_weight, G);
    30043598    xtnw=xtnw+clock()-to;
     
    30073601#endif
    30083602
    3009     if(Overflow_Error == TRUE) {
     3603    if(Overflow_Error == TRUE)
     3604    {
    30103605      //PrintS("\n// ** The next vector does NOT stay in Cone!!\n");
    30113606#ifdef TEST_OVERFLOW
     
    30153610      nnwinC = 0;
    30163611      if(tp_deg == nV)
     3612      {
    30173613        nlast = 1;
    3018 
     3614      }
    30193615      delete next_weight;
    30203616      break;
    30213617    }
    30223618
    3023     if(MivComp(next_weight, ivNull) == 1) {
     3619    if(MivComp(next_weight, ivNull) == 1)
     3620    {
    30243621      //newRing = currRing;
    30253622      delete next_weight;
     
    30273624    }
    30283625
    3029     if(MivComp(next_weight, target_weight) == 1)  {
     3626    if(MivComp(next_weight, target_weight) == 1)
     3627    {
    30303628      if(tp_deg == nV)
     3629      {
    30313630        endwalks = 1;
    3032       else {
    3033       REC_LAST_GB_ALT2:
    3034         nOverflow_Error = Overflow_Error;
     3631      }
     3632      else
     3633      {
     3634        // REC_LAST_GB_ALT2:
     3635        //nOverflow_Error = Overflow_Error;
    30353636        tproc=tproc+clock()-tinput;
    30363637        /*
     
    30453646    }
    30463647
    3047     for(i=nV-1; i>=0; i--) {
    3048       //(*extra_curr_weight)[i] = (*curr_weight)[i];
     3648    for(i=nV-1; i>=0; i--)
     3649    {
    30493650      (*curr_weight)[i] = (*next_weight)[i];
    30503651    }
     
    30543655  delete ivNull;
    30553656
    3056   if(tp_deg != nV) {
     3657  if(tp_deg != nV)
     3658  {
    30573659    newRing = currRing;
    30583660
    30593661    if (rParameter(currRing) != NULL)
     3662    {
    30603663      DefRingParlp();
     3664    }
    30613665    else
     3666    {
    30623667      VMrDefaultlp();
    3063 
     3668    }
    30643669    F1 = idrMoveR(G, newRing,currRing);
    30653670
    3066     if(nnwinC == 0 || test_w_in_ConeCC(F1, target_weight) != 1 ) {
    3067       nOverflow_Error = Overflow_Error;
     3671    if(nnwinC == 0 || test_w_in_ConeCC(F1, target_weight) != 1 )
     3672    {
     3673      // nOverflow_Error = Overflow_Error;
    30683674      //Print("\n//  takes %d steps and calls \"Rec_LastGB (%d):", tp_deg+1);
    30693675      tproc=tproc+clock()-tinput;
     
    30763682    result = idrMoveR(F1, TargetRing,currRing);
    30773683  }
    3078   else  {
    3079     if(nlast == 1) {
     3684  else
     3685  {
     3686    if(nlast == 1)
     3687    {
    30803688      JUNI_STD:
    30813689
    30823690      newRing = currRing;
    30833691      if (rParameter(currRing) != NULL)
     3692      {
    30843693        DefRingParlp();
     3694      }
    30853695      else
     3696      {
    30863697        VMrDefaultlp();
    3087 
     3698      }
    30883699      KSTD_Finish:
    30893700      if(isGB == FALSE)
     3701      {
    30903702        F1 = idrMoveR(G, newRing,currRing);
     3703      }
    30913704      else
     3705      {
    30923706        F1 = G;
     3707      }
    30933708      to=clock();
    30943709      // Print("\n// apply the Buchberger's alg in ring = %s",rString(currRing));
     
    31083723
    31093724  if(Overflow_Error == FALSE)
     3725    {
    31103726    Overflow_Error=nError;
    3111   /*
    3112   Print("\n// \"Rec_LastGB\" (%d) took %d steps and %.2f sec.Overflow_Error (%d)", tp_deg,
    3113         nwalk, ((double) tproc)/1000000, nOverflow_Error);
    3114   */
     3727    }
     3728// Print("\n// \"Rec_LastGB\" (%d) took %d steps and %.2f sec.Overflow_Error (%d)", tp_deg, nwalk, ((double) tproc)/1000000, nOverflow_Error);
    31153729  return(result);
    31163730}
     
    31243738   basis or n times, where n is the numbers of variables.
    31253739*/
    3126 /* walk + recursive LastGB */
     3740
     3741/******************************
     3742 * walk + recursive LastGB    *
     3743 ******************************/
    31273744ideal MAltwalk2(ideal Go, intvec* curr_weight, intvec* target_weight)
    31283745{
    31293746  Set_Error(FALSE);
    31303747  Overflow_Error = FALSE;
    3131   BOOLEAN nOverflow_Error = FALSE;
     3748  //BOOLEAN nOverflow_Error = FALSE;
    31323749  //Print("// pSetm_Error = (%d)", ErrorCheck());
    31333750
     
    31383755  nstep = 0;
    31393756  int i, nV = currRing->N;
    3140   int nwalk=0, endwalks=0, nhilb=1;
    3141 
    3142   ideal Gomega, M, F, Gomega1, Gomega2, M1, F1, G, G1;
    3143   ring endRing, newRing, oldRing;
     3757  int nwalk=0, endwalks=0;
     3758  // int nhilb = 1;
     3759  ideal Gomega, M, F, Gomega1, Gomega2, M1, F1, G;
     3760  //ideal  G1;
     3761  //ring endRing;
     3762  ring newRing, oldRing;
    31443763  intvec* ivNull = new intvec(nV);
    31453764  intvec* next_weight;
     3765#if 0
    31463766  intvec* extra_curr_weight = new intvec(nV);
    3147   intvec* hilb_func;
     3767#endif
     3768  //intvec* hilb_func;
    31483769  intvec* exivlp = Mivlp(nV);
    31493770
     
    31623783  */
    31633784  if(currRing->order[0] == ringorder_a)
     3785  {
    31643786    goto NEXT_VECTOR;
    3165 
     3787  }
    31663788  while(1)
    31673789  {
     
    31853807    /* define a new ring that its ordering is "(a(curr_weight),lp) */
    31863808    if (rParameter(currRing) != NULL)
     3809    {
    31873810      DefRingPar(curr_weight);
     3811    }
    31883812    else
     3813    {
    31893814      VMrDefault(curr_weight);
    3190 
     3815    }
    31913816    newRing = currRing;
    31923817    Gomega1 = idrMoveR(Gomega, oldRing,currRing);
     
    32413866#endif
    32423867
    3243 
    32443868      newRing = currRing;
    32453869      if (rParameter(currRing) != NULL)
     3870      {
    32463871        DefRingPar(target_weight);
     3872      }
    32473873      else
     3874      {
    32483875        VMrDefault(target_weight);
    3249 
     3876      }
    32503877      F1 = idrMoveR(G, newRing,currRing);
    32513878      G = MstdCC(F1);
     
    32663893      if(MivSame(target_weight, exivlp)==1)
    32673894      {
    3268       LAST_GB_ALT2:
    3269         nOverflow_Error = Overflow_Error;
     3895     // LAST_GB_ALT2:
     3896        //nOverflow_Error = Overflow_Error;
    32703897        tproc = clock()-xftinput;
    32713898        //Print("\n// takes %d steps and calls the recursion of level 2:",  nwalk);
     
    32863913    delete next_weight;
    32873914  }
     3915#ifdef TEST_OVERFLOW
    32883916 TEST_OVERFLOW_OI:
     3917#endif
    32893918  rChangeCurrRing(XXRing);
    32903919  G = idrMoveR(G, newRing,currRing);
     
    32933922
    32943923#ifdef TIME_TEST
    3295   Print("\n// \"Main procedure\"  took %d steps dnd %.2f sec. Overflow_Error (%d)", nwalk,
    3296         ((double) tproc)/1000000, nOverflow_Error);
     3924 // Print("\n// \"Main procedure\"  took %d steps dnd %.2f sec. Overflow_Error (%d)", nwalk, ((double) tproc)/1000000, nOverflow_Error);
    32973925
    32983926  TimeStringFractal(xftinput, tostd, xtif, xtstd, xtextra,xtlift, xtred,xtnw);
     
    33073935
    33083936
    3309 /* 5.5.02 */
    3310 /* The implementation of the fractal walk algorithmus */
    3311 
    3312 /* The main procedur Mfwalk calls the recursive Subroutine
    3313    rec_fractal_call to compute the wanted Grï¿œbner basis.
    3314    At the main procedur we compute the reduced Grï¿œbner basis w.r.t. a "fast"
    3315    order, e.g. "dp" and a sequence of weight vectors which are row vectors
    3316    of a matrix. This matrix defines the given monomial order, e.g. "lp"
    3317 */
    3318 
    3319 
    3320 
    3321 
    3322 /* perturb the matrix order of  "lex" */
     3937/**************************************
     3938 * perturb the matrix order of  "lex" *
     3939 **************************************/
    33233940static intvec* NewVectorlp(ideal I)
    33243941{
     
    33373954intvec* Xivlp;
    33383955
    3339 
    3340 
    3341 /***********************************************************************
    3342   The procedur REC_GB_Mwalk computes a GB for <G> w.r.t. the weight order
    3343   otw, where G is a reduced GB w.r.t. the weight order cw.
    3344   The new procedur Mwalk calls REC_GB.
    3345 ************************************************************************/
     3956/********************************
     3957 * compute a next weight vector *
     3958 ********************************/
     3959static intvec* MWalkRandomNextWeight(ideal G, intvec* curr_weight, intvec* target_weight, int weight_rad, int pert_deg)
     3960{
     3961  int i, weight_norm;
     3962  int nV = currRing->N;
     3963  intvec* next_weight2;
     3964  intvec* next_weight22 = new intvec(nV);
     3965  intvec* next_weight = MwalkNextWeightCC(curr_weight,target_weight, G);
     3966  if(MivComp(next_weight, target_weight) == 1)
     3967  {
     3968    return(next_weight);
     3969  }
     3970  else
     3971  {
     3972    //compute a perturbed next weight vector "next_weight1"
     3973    intvec* next_weight1 = MkInterRedNextWeight(MPertVectors(G, MivMatrixOrder(curr_weight), pert_deg), target_weight, G);
     3974    //Print("\n // size of next_weight1 = %d", sizeof((*next_weight1)));
     3975
     3976    //compute a random next weight vector "next_weight2"
     3977    while(1)
     3978    {
     3979      weight_norm = 0;
     3980      while(weight_norm == 0)
     3981      {
     3982        for(i=0; i<nV; i++)
     3983        {
     3984          //Print("\n//  next_weight[%d]  = %d", i, (*next_weight)[i]);
     3985          (*next_weight22)[i] = rand() % 60000 - 30000;
     3986          weight_norm = weight_norm + (*next_weight22)[i]*(*next_weight22)[i];
     3987        }
     3988        weight_norm = 1 + floor(sqrt(weight_norm));
     3989      }
     3990
     3991      for(i=nV-1; i>=0; i--)
     3992      {
     3993        if((*next_weight22)[i] < 0)
     3994        {
     3995          (*next_weight22)[i] = 1 + (*curr_weight)[i] + floor(weight_rad*(*next_weight22)[i]/weight_norm);
     3996        }
     3997        else
     3998        {
     3999          (*next_weight22)[i] = (*curr_weight)[i] + floor(weight_rad*(*next_weight22)[i]/weight_norm);
     4000        }
     4001      //Print("\n//  next_weight22[%d]  = %d", i, (*next_weight22)[i]);
     4002      }
     4003     
     4004      if(test_w_in_ConeCC(G, next_weight22) == 1)
     4005      {
     4006        //Print("\n//MWalkRandomNextWeight: next_weight2 im Kegel\n");
     4007        next_weight2 = MkInterRedNextWeight(next_weight22, target_weight, G);
     4008        delete next_weight22;
     4009        break;
     4010      }
     4011    }
     4012    intvec* result = new intvec(nV);
     4013    ideal G_test = MwalkInitialForm(G, next_weight);
     4014    ideal G_test1 = MwalkInitialForm(G, next_weight1);
     4015    ideal G_test2 = MwalkInitialForm(G, next_weight2);
     4016
     4017    // compare next_weights
     4018    if(IDELEMS(G_test1) < IDELEMS(G_test))
     4019    {
     4020      if(IDELEMS(G_test2) <= IDELEMS(G_test1)) // |G_test2| <= |G_test1| < |G_test|
     4021      {
     4022        for(i=0; i<nV; i++)
     4023        {
     4024          (*result)[i] = (*next_weight2)[i];
     4025        }
     4026      }
     4027      else // |G_test1| < |G_test|, |G_test1| < |G_test2|
     4028      {
     4029        for(i=0; i<nV; i++)
     4030        {
     4031          (*result)[i] = (*next_weight1)[i];
     4032        }
     4033      }   
     4034    }
     4035    else
     4036    {
     4037      if(IDELEMS(G_test2) <= IDELEMS(G_test)) // |G_test2| <= |G_test| <= |G_test1|
     4038      {
     4039        for(i=0; i<nV; i++)
     4040        {
     4041          (*result)[i] = (*next_weight2)[i];
     4042        }
     4043      }
     4044      else // |G_test| <= |G_test1|, |G_test| < |G_test2|
     4045      {
     4046        for(i=0; i<nV; i++)
     4047        {
     4048          (*result)[i] = (*next_weight)[i];
     4049        }
     4050      }
     4051    }
     4052    delete next_weight;
     4053    delete next_weight1;
     4054    idDelete(&G_test);
     4055    idDelete(&G_test1);
     4056    idDelete(&G_test2);
     4057    if(test_w_in_ConeCC(G, result) == 1)
     4058    {
     4059      delete next_weight2;
     4060      return result;
     4061    }
     4062    else
     4063    {
     4064      delete result;
     4065      return next_weight2;
     4066    }
     4067  }
     4068}
     4069
     4070
     4071/***************************************************************************
     4072 * The procedur REC_GB_Mwalk computes a GB for <G> w.r.t. the weight order *
     4073 * otw, where G is a reduced GB w.r.t. the weight order cw.                *
     4074 * The new procedur Mwalk calls REC_GB.                                    *
     4075 ***************************************************************************/
    33464076static ideal REC_GB_Mwalk(ideal G, intvec* curr_weight, intvec* orig_target_weight,
    33474077                          int tp_deg, int npwinc)
     
    33504080  Overflow_Error = FALSE;
    33514081
    3352   int i,  nV = currRing->N, ntwC,  npwinC;
     4082  int i,  nV = currRing->N;
    33534083  int nwalk=0, endwalks=0, nnwinC=1, nlast = 0;
    33544084  ideal Gomega, M, F, Gomega1, Gomega2, M1,F1,result,ssG;
    33554085  ring newRing, oldRing, TargetRing;
    3356   intvec* iv_M_lp;
    33574086  intvec* target_weight;
    33584087  intvec* ivNull = new intvec(nV);
     4088#ifndef BUCHBERGER_ALG
    33594089  intvec* hilb_func;
    3360   BOOLEAN isGB = FALSE;
    3361 
    3362   /* to avoid (1,0,...,0) as the target vector */
     4090  // to avoid (1,0,...,0) as the target vector
    33634091  intvec* last_omega = new intvec(nV);
    33644092  for(i=nV-1; i>0; i--)
     4093  {
    33654094    (*last_omega)[i] = 1;
     4095  }
    33664096  (*last_omega)[0] = 10000;
     4097#endif
     4098  BOOLEAN isGB = FALSE;
    33674099
    33684100  ring EXXRing = currRing;
    33694101
    3370   /* compute a pertubed weight vector of the target weight vector */
     4102  // compute a pertubed weight vector of the target weight vector
    33714103  if(tp_deg > 1 && tp_deg <= nV)
    33724104  {
    33734105    ideal H0 = idHeadCC(G);
    33744106    if (rParameter(currRing) != NULL)
     4107    {
    33754108      DefRingPar(orig_target_weight);
     4109    }
    33764110    else
     4111    {
    33774112      VMrDefault(orig_target_weight);
    3378 
     4113    }
    33794114    TargetRing = currRing;
    33804115    ssG = idrMoveR(G,EXXRing,currRing);
     
    33864121    if(test_G_GB_walk(H0_tmp,H1)==1)
    33874122    {
    3388       //Print("//input in %d-th recursive is a GB",tp_deg);
    3389       idDelete(&H0_tmp);idDelete(&H1);
     4123      Print("\n//REC_GB_Mwalk: input in %d-th recursive is a GB!\n",tp_deg);
     4124      idDelete(&H0_tmp);
     4125      idDelete(&H1);
    33904126      G = ssG;
    33914127      ssG = NULL;
     
    33984134      }
    33994135      else
     4136      {
    34004137        goto LastGB_Finish;
    3401     }
    3402     idDelete(&H0_tmp);   idDelete(&H1);
    3403 
    3404     intvec* ivlp = Mivlp(nV);
    3405     if( MivSame(orig_target_weight, ivlp)==1 )
    3406       iv_M_lp = MivMatrixOrderlp(nV);
    3407     else
    3408       iv_M_lp = MivMatrixOrder(orig_target_weight);
    3409 
    3410     //target_weight  = MPertVectorslp(ssG, iv_M_lp, tp_deg);
    3411     target_weight  = MPertVectors(ssG, iv_M_lp, tp_deg);
    3412 
    3413     delete ivlp;
    3414     delete iv_M_lp;
     4138      }
     4139    }
     4140    idDelete(&H0_tmp);
     4141    idDelete(&H1);
     4142
     4143    target_weight  = MPertVectors(ssG, MivMatrixOrder(orig_target_weight), tp_deg);
    34154144
    34164145    rChangeCurrRing(EXXRing);
     
    34234152    nstep++;
    34244153    if(nwalk == 1)
     4154    {
    34254155      goto NEXT_STEP;
    3426 
    3427     //Print("\n// Entering the %d-th step in the %d-th recursive:",nwalk,tp_deg);
     4156    }
     4157    Print("\n//REC_GB_Mwalk: Entering the %d-th step in the %d-th recursive:\n",nwalk,tp_deg);
    34284158    to = clock();
    3429     /* compute an initial form ideal of <G> w.r.t. "curr_vector" */
     4159    // compute an initial form ideal of <G> w.r.t. "curr_vector"
    34304160    Gomega = MwalkInitialForm(G, curr_weight);
    34314161    xtif = xtif + clock()-to;
     
    34334163#ifndef  BUCHBERGER_ALG
    34344164    if(isNolVector(curr_weight) == 0)
     4165    {
    34354166      hilb_func = hFirstSeries(Gomega,NULL,NULL,curr_weight,currRing);
     4167    }
    34364168    else
     4169    {
    34374170      hilb_func = hFirstSeries(Gomega,NULL,NULL,last_omega,currRing);
    3438 #endif // BUCHBERGER_ALG
     4171    }
     4172#endif
    34394173
    34404174    oldRing = currRing;
    34414175
    3442     /* define a new ring that its ordering is "(a(curr_weight),lp) */
     4176    // define a new ring with ordering "(a(curr_weight),lp)
    34434177    if (rParameter(currRing) != NULL)
     4178    {
    34444179      DefRingPar(curr_weight);
     4180    }
    34454181    else
     4182    {
    34464183      VMrDefault(curr_weight);
    3447 
     4184    }
    34484185    newRing = currRing;
    34494186    Gomega1 = idrMoveR(Gomega, oldRing,currRing);
    34504187
    34514188    to = clock();
    3452     /* compute a reduced Groebner basis of <Gomega> w.r.t. "newRing" */
     4189    // compute a reduced Groebner basis of <Gomega> w.r.t. "newRing"
    34534190#ifdef  BUCHBERGER_ALG
    34544191    M = MstdhomCC(Gomega1);
     
    34564193    M=kStd(Gomega1,NULL,isHomog,NULL,hilb_func,0,NULL,curr_weight);
    34574194    delete hilb_func;
    3458 #endif // BUCHBERGER_ALG
     4195#endif
    34594196    xtstd = xtstd + clock() - to;
    34604197
    3461     /* change the ring to oldRing */
     4198    // change the ring to oldRing
    34624199    rChangeCurrRing(oldRing);
    34634200
     
    34744211
    34754212
    3476     /* change the ring to newRing */
     4213    // change the ring to newRing
    34774214    rChangeCurrRing(newRing);
    34784215    F1 = idrMoveR(F, oldRing,currRing);
    34794216
    34804217    to = clock();
    3481     /* reduce the Groebner basis <G> w.r.t. new ring */
     4218    // reduce the Groebner basis <G> w.r.t. new ring
    34824219    G = kInterRedCC(F1, NULL);
    34834220    xtred = xtred + clock() -to;
     
    34864223
    34874224    if(endwalks == 1)
     4225    {
    34884226      break;
    3489 
     4227    }
    34904228  NEXT_STEP:
    34914229    to = clock();
    3492     /* compute a next weight vector */
     4230    // compute a next weight vector
    34934231    intvec* next_weight = MkInterRedNextWeight(curr_weight,target_weight, G);
     4232 
     4233
    34944234    xtnw = xtnw + clock() - to;
    34954235
     
    34984238#endif
    34994239
    3500     /*check whether the computed vector does in the correct cone */
    3501     //ntwC = test_w_in_ConeCC(G, next_weight);
    3502     //if(ntwC != 1)
    35034240    if(Overflow_Error == TRUE)
    35044241    {
    3505       PrintS("\n// ** The computed vector does NOT stay in Cone!!\n");
     4242      PrintS("\n//REC_GB_Mwalk: The computed vector does NOT stay in the correct cone!!\n");
    35064243      nnwinC = 0;
    35074244      if(tp_deg == nV)
     4245      {
    35084246        nlast = 1;
     4247      }
    35094248      delete next_weight;
    35104249      break;
     
    35204259    {
    35214260      if(tp_deg == nV)
     4261      {
    35224262        endwalks = 1;
    3523       else {
     4263      }
     4264      else
     4265      {
    35244266        G = REC_GB_Mwalk(G,curr_weight, orig_target_weight, tp_deg+1,nnwinC);
    35254267        newRing = currRing;
     
    35294271    }
    35304272
    3531     /* 06.11.01 NOT Changed */
    35324273    for(i=nV-1; i>=0; i--)
     4274    {
    35334275      (*curr_weight)[i] = (*next_weight)[i];
    3534 
     4276    }
    35354277    delete next_weight;
    3536   }//while
     4278  }
    35374279
    35384280  delete ivNull;
     
    35404282  if(tp_deg != nV)
    35414283  {
    3542     //28.07.03
    35434284    newRing = currRing;
    35444285
    35454286    if (rParameter(currRing) != NULL)
    3546       //  DefRingParlp(); //
     4287    {
    35474288      DefRingPar(orig_target_weight);
     4289    }
    35484290    else
     4291    {
    35494292      VMrDefault(orig_target_weight);
    3550 
    3551 
     4293    }
    35524294    F1 = idrMoveR(G, newRing,currRing);
    35534295
    35544296    if(nnwinC == 0)
     4297    {
    35554298      F1 = REC_GB_Mwalk(F1,curr_weight, orig_target_weight, tp_deg+1,nnwinC);
     4299    }
    35564300    else
     4301    {
    35574302      if(test_w_in_ConeCC(F1, target_weight) != 1)
     4303      {
    35584304        F1 = REC_GB_Mwalk(F1,curr_weight, orig_target_weight,tp_deg+1,nnwinC);
    3559 
     4305      }
     4306    }
    35604307    delete target_weight;
    35614308
     
    35694316    {
    35704317      if (rParameter(currRing) != NULL)
     4318      {
    35714319        DefRingPar(orig_target_weight);
     4320      }
    35724321      else
     4322      {
    35734323        VMrDefault(orig_target_weight);
    3574 
    3575 
     4324      }
    35764325    KSTD_Finish:
    35774326      if(isGB == FALSE)
     4327      {
    35784328        F1 = idrMoveR(G, newRing,currRing);
     4329      }
    35794330      else
     4331      {
    35804332        F1 = G;
     4333      }
    35814334      to=clock();
    3582       /* apply Buchberger alg to compute a red. GB of F1 */
     4335      // apply Buchberger alg to compute a red. GB of F1
    35834336      G = MstdCC(F1);
    35844337      xtextra=clock()-to;
     
    35934346
    35944347  if(Overflow_Error == FALSE)
     4348    {
    35954349    Overflow_Error = nError;
    3596 
     4350    }
     4351#ifndef BUCHBERGER_ALG
     4352  delete last_omega;
     4353#endif
    35974354  return(result);
    35984355}
    35994356
    36004357
    3601 /* 08.09.02 */
    3602 /******** THE NEW GRï¿œBNER WALK ALGORITHM **********/
    3603 /* Grï¿œbnerwalk with a recursive "second" alternative GW, REC_GB_Mwalk
    3604    that only computes the last reduced GB */
     4358// THE NEW GROEBNER WALK ALGORITHM
     4359// Groebnerwalk with a recursive "second" alternative GW, called REC_GB_Mwalk that only computes the last reduced GB
    36054360ideal Mwalk(ideal Go, intvec* curr_weight, intvec* target_weight)
    36064361{
     
    36144369  clock_t tim;
    36154370  nstep=0;
    3616   int i, nV = currRing->N;
    3617   int nwalk=0, endwalks=0;
    3618 
    3619   ideal Gomega, M, F, Gomega1, Gomega2, M1, F1, G, G1;
    3620   ring endRing, newRing, oldRing;
     4371  int i;
     4372  int nV = currRing->N;
     4373  int nwalk=0;
     4374  int endwalks=0;
     4375
     4376  ideal Gomega, M, F, Gomega1, Gomega2, M1, F1, G;
     4377  //ideal G1;
     4378  //ring endRing;
     4379  ring newRing, oldRing;
    36214380  intvec* ivNull = new intvec(nV);
    36224381  intvec* exivlp = Mivlp(nV);
     4382#ifndef BUCHBERGER_ALG
    36234383  intvec* hilb_func;
    3624 
     4384#endif
    36254385  intvec* tmp_weight = new intvec(nV);
    36264386  for(i=nV-1; i>=0; i--)
    36274387    (*tmp_weight)[i] = (*curr_weight)[i];
    36284388
    3629    /* to avoid (1,0,...,0) as the target vector */
     4389   // to avoid (1,0,...,0) as the target vector
    36304390  intvec* last_omega = new intvec(nV);
    36314391  for(i=nV-1; i>0; i--)
     
    36364396
    36374397  to = clock();
    3638   /* the monomial ordering of this current ring would be "dp" */
     4398  // the monomial ordering of this current ring would be "dp"
    36394399  G = MstdCC(Go);
    36404400  tostd = clock()-to;
     
    36484408    nstep ++;
    36494409    to = clock();
    3650     /* compute an initial form ideal of <G> w.r.t. "curr_vector" */
     4410    // compute an initial form ideal of <G> w.r.t. "curr_vector"
    36514411    Gomega = MwalkInitialForm(G, curr_weight);
    36524412    tif = tif + clock()-to;
     
    36914451      /* create a new ring newRing */
    36924452       if (rParameter(currRing) != NULL)
     4453       {
    36934454         DefRingPar(curr_weight);
     4455       }
    36944456       else
     4457       {
    36954458         VMrDefault(curr_weight);
    3696 
     4459       }
    36974460      newRing = currRing;
    36984461      F1 = idrMoveR(F, oldRing,currRing);
     
    37034466#ifndef  BUCHBERGER_ALG
    37044467      if(isNolVector(curr_weight) == 0)
     4468      {
    37054469        hilb_func = hFirstSeries(Gomega,NULL,NULL,curr_weight,currRing);
     4470      }
    37064471      else
     4472      {
    37074473        hilb_func = hFirstSeries(Gomega,NULL,NULL,last_omega,currRing);
     4474      }
    37084475#endif // BUCHBERGER_ALG
    37094476
    3710       /* define a new ring that its ordering is "(a(curr_weight),lp) */
     4477      // define a new ring that its ordering is "(a(curr_weight),lp)
    37114478      if (rParameter(currRing) != NULL)
     4479      {
    37124480        DefRingPar(curr_weight);
     4481      }
    37134482      else
     4483      {
    37144484        VMrDefault(curr_weight);
    3715 
     4485      }
    37164486      newRing = currRing;
    37174487      Gomega1 = idrMoveR(Gomega, oldRing,currRing);
    37184488
    37194489      to = clock();
    3720       /* compute a reduced Groebner basis of <Gomega> w.r.t. "newRing" */
     4490      // compute a reduced Groebner basis of <Gomega> w.r.t. "newRing"
    37214491#ifdef  BUCHBERGER_ALG
    37224492      M = MstdhomCC(Gomega1);
     
    37274497      tstd = tstd + clock() - to;
    37284498
    3729       /* change the ring to oldRing */
     4499      // change the ring to oldRing
    37304500      rChangeCurrRing(oldRing);
    37314501      M1 =  idrMoveR(M, newRing,currRing);
     
    37334503
    37344504      to = clock();
    3735       /* compute a representation of the generators of submod (M)
    3736          with respect to those of mod (Gomega).
    3737          Gomega is a reduced Groebner basis w.r.t. the current ring */
     4505      // 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.
    37384506      F = MLifttwoIdeal(Gomega2, M1, G);
    37394507      tlift = tlift + clock() - to;
     
    37434511      idDelete(&G);
    37444512
    3745       /* change the ring to newRing */
     4513      // change the ring to newRing
    37464514      rChangeCurrRing(newRing);
    37474515      F1 = idrMoveR(F, oldRing,currRing);
     
    37494517
    37504518    to = clock();
    3751     /* reduce the Groebner basis <G> w.r.t. new ring */
     4519    // reduce the Groebner basis <G> w.r.t. new ring
    37524520    G = kInterRedCC(F1, NULL);
    37534521    if(endwalks != 1)
     4522    {
    37544523      tred = tred + clock() - to;
     4524    }
    37554525    else
     4526    {
    37564527      xtred = xtred + clock() - to;
    3757 
     4528    }
    37584529    idDelete(&F1);
    37594530    if(endwalks == 1)
     4531    {
    37604532      break;
    3761 
     4533    }
    37624534  NEXT_VECTOR:
    37634535    to = clock();
    3764     /* compute a next weight vector */
     4536    // compute a next weight vector
    37654537    intvec* next_weight = MkInterRedNextWeight(curr_weight,target_weight,G);
    37664538    tnw = tnw + clock() - to;
     
    37764548
    37774549      if (rParameter(currRing) != NULL)
     4550      {
    37784551        DefRingPar(target_weight);
     4552      }
    37794553      else
     4554      {
    37804555        VMrDefault(target_weight);
    3781 
     4556      }
    37824557      F1 = idrMoveR(G, newRing,currRing);
    37834558      G = MstdCC(F1);
     
    37954570    }
    37964571    if(MivComp(next_weight, target_weight) == 1)
     4572    {
    37974573      endwalks = 1;
    3798 
    3799     for(i=nV-1; i>=0; i--) {
     4574    }
     4575    for(i=nV-1; i>=0; i--)
     4576    {
    38004577      (*tmp_weight)[i] = (*curr_weight)[i];
    38014578      (*curr_weight)[i] = (*next_weight)[i];
     
    38194596}
    38204597
    3821   /* 2.12.02*/
     4598// 07.11.2012
     4599// THE RANDOM WALK ALGORITHM
     4600ideal Mrwalk(ideal Go, intvec* curr_weight, intvec* target_weight, int weight_rad, int pert_deg)
     4601{
     4602  Set_Error(FALSE);
     4603  Overflow_Error = FALSE;
     4604  //Print("// pSetm_Error = (%d)", ErrorCheck());
     4605
     4606  clock_t tinput, tostd, tif=0, tstd=0, tlift=0, tred=0, tnw=0;
     4607  xtif=0; xtstd=0; xtlift=0; xtred=0; xtnw=0;
     4608  tinput = clock();
     4609  clock_t tim;
     4610  nstep=0;
     4611  int i,k,nwalk,endwalks = 0;
     4612  int nV = currRing->N;
     4613
     4614  ideal Gomega, M, F, Gomega1, Gomega2, M1, F1, G;
     4615  //ideal G1;
     4616  //ring endRing;
     4617  ring newRing, oldRing;
     4618  intvec* ivNull = new intvec(nV);
     4619  intvec* exivlp = Mivlp(nV);
     4620  intvec* tmp_weight = new intvec(nV);
     4621  for(i=nV-1; i>=0; i--)
     4622  {
     4623    (*tmp_weight)[i] = (*curr_weight)[i];
     4624  }
     4625#ifndef BUCHBERGER_ALG
     4626  intvec* hilb_func;
     4627   // to avoid (1,0,...,0) as the target vector
     4628  intvec* last_omega = new intvec(nV);
     4629  for(i=nV-1; i>0; i--)
     4630  {
     4631    (*last_omega)[i] = 1;
     4632  }
     4633  (*last_omega)[0] = 10000;
     4634#endif
     4635  ring XXRing = currRing;
     4636
     4637  to = clock();
     4638  G = MstdCC(Go);
     4639  tostd = clock()-to;
     4640
     4641  //if(currRing->order[0] == ringorder_a)
     4642  //  goto NEXT_VECTOR;
     4643nwalk = 0;
     4644  while(1)
     4645  {
     4646    nwalk ++;
     4647    nstep ++;
     4648    to = clock();
     4649#ifdef CHECK_IDEAL_MWALK
     4650    idString(G,"G");
     4651#endif
     4652    Gomega = MwalkInitialForm(G, curr_weight);// compute an initial form ideal of <G> w.r.t. "curr_vector"
     4653    tif = tif + clock()-to;
     4654    oldRing = currRing;
     4655#ifdef CHECK_IDEAL_MWALK
     4656    idString(Gomega,"G_omega");
     4657#endif
     4658    if(endwalks == 1)
     4659    {
     4660      // compute a reduced Groebner basis of Gomega w.r.t. >>_cw by the recursive changed perturbation walk alg.
     4661      tim = clock();
     4662      Print("\n//Mrwalk: Groebnerwalk took %d steps.\n", nwalk);
     4663      //PrintS("\n//Mrwalk: call the rec. Pert. Walk to compute a red GB of:");
     4664      //idString(Gomega, "G_omega");
     4665      M = REC_GB_Mwalk(idCopy(Gomega), tmp_weight, curr_weight,pert_deg,1);
     4666      Print("\n//Mrwalk: time for the last std(Gw)  = %.2f sec\n",((double) (clock()-tim)/1000000));
     4667     
     4668#ifdef CHECK_IDEAL_MWALK
     4669      idString(Gomega, "G_omega");
     4670      idString(M, "M");
     4671#endif
     4672      to = clock();
     4673      F = MLifttwoIdeal(Gomega, M, G);
     4674      xtlift = xtlift + clock() - to;
     4675
     4676      idDelete(&Gomega);
     4677      idDelete(&M);
     4678      idDelete(&G);
     4679
     4680      oldRing = currRing;
     4681      VMrDefault(curr_weight); //define a new ring with ordering "(a(curr_weight),lp)
     4682
     4683       /*if (rParameter(currRing) != NULL)
     4684       {
     4685         DefRingPar(curr_weight);
     4686       }
     4687       else
     4688       {
     4689         VMrDefault(curr_weight);
     4690       }*/
     4691      newRing = currRing;
     4692      F1 = idrMoveR(F, oldRing,currRing);
     4693    }
     4694    else
     4695    {
     4696    NORMAL_GW:
     4697#ifndef  BUCHBERGER_ALG
     4698      if(isNolVector(curr_weight) == 0)
     4699      {
     4700        hilb_func = hFirstSeries(Gomega,NULL,NULL,curr_weight,currRing);
     4701      }
     4702      else
     4703      {
     4704        hilb_func = hFirstSeries(Gomega,NULL,NULL,last_omega,currRing);
     4705      }
     4706#endif
     4707/*
     4708      if (rParameter(currRing) != NULL)
     4709      {
     4710        DefRingPar(curr_weight);
     4711      }
     4712      else
     4713      {
     4714        VMrDefault(curr_weight);
     4715      }*/
     4716
     4717      VMrDefault(curr_weight); //define a new ring with ordering "(a(curr_weight),lp)
     4718      newRing = currRing;
     4719      Gomega1 = idrMoveR(Gomega, oldRing,currRing);
     4720#ifdef CHECK_IDEAL_MWALK
     4721      idString(Gomega1, "G_omega1");
     4722#endif
     4723      // compute a reduced Groebner basis of <Gomega> w.r.t. "newRing"
     4724      to = clock();
     4725#ifndef  BUCHBERGER_ALG
     4726      M=kStd(Gomega1,NULL,isHomog,NULL,hilb_func,0,NULL,curr_weight);
     4727      delete hilb_func;
     4728#else
     4729      //M = MstdhomCC(Gomega1);
     4730      //M = MstdCC(Gomega1);
     4731      //M = kStd(Gomega1, NULL, testHomog,NULL, NULL,0,0,curr_weight);
     4732      M = kStd(Gomega1,NULL,testHomog,NULL,NULL,0,0,NULL);
     4733#endif
     4734      tstd = tstd + clock() - to;
     4735#ifdef CHECK_IDEAL_MWALK
     4736      idString(M, "M");
     4737#endif
     4738      //change the ring to oldRing
     4739      rChangeCurrRing(oldRing);
     4740      M1 =  idrMoveR(M, newRing,currRing);
     4741#ifdef CHECK_IDEAL_MWALK
     4742      idString(M1, "M1");
     4743#endif
     4744      Gomega2 =  idrMoveR(Gomega1, newRing,currRing);
     4745
     4746      to = clock();
     4747      // 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
     4748      F = MLifttwoIdeal(Gomega2, M1, G);
     4749#ifdef CHECK_IDEAL_MWALK
     4750      idString(F, "F");
     4751#endif
     4752      tlift = tlift + clock() - to;
     4753
     4754      idDelete(&M1);
     4755      idDelete(&Gomega2);
     4756      idDelete(&G);
     4757
     4758      rChangeCurrRing(newRing); // change the ring to newRing
     4759      F1 = idrMoveR(F,oldRing,currRing);
     4760    }
     4761
     4762    to = clock();
     4763    G = kInterRedCC(F1, NULL); //reduce the Groebner basis <G> w.r.t. new ring
     4764#ifdef CHECK_IDEAL_MWALK
     4765    idString(G, "G");
     4766#endif
     4767    idDelete(&F1);
     4768    if(endwalks == 1)
     4769    {
     4770      xtred = xtred + clock() - to;
     4771      break;
     4772    }
     4773    else
     4774    {
     4775      tred = tred + clock() - to;
     4776    }
     4777
     4778  NEXT_VECTOR:
     4779    to = clock();
     4780    //intvec* next_weight = MkInterRedNextWeight(curr_weight,target_weight,G);
     4781   intvec* next_weight = MWalkRandomNextWeight(G, curr_weight, target_weight, weight_rad, pert_deg);
     4782   
     4783  tnw = tnw + clock() - to;
     4784//#ifdef PRINT_VECTORS
     4785  MivString(curr_weight, target_weight, next_weight);
     4786//#endif
     4787    if(Overflow_Error == TRUE)
     4788    {
     4789      newRing = currRing;
     4790      PrintS("\n//Mrwalk: The computed vector does NOT stay in Cone!!\n");
     4791
     4792      if (rParameter(currRing) != NULL)
     4793      {
     4794        DefRingPar(target_weight);
     4795      }
     4796      else
     4797      {
     4798        VMrDefault(target_weight);  //define a new ring with ordering "(a(curr_weight),lp)
     4799      }
     4800      F1 = idrMoveR(G, newRing,currRing);
     4801      G = MstdCC(F1);
     4802      idDelete(&F1);
     4803
     4804      newRing = currRing;
     4805      break;
     4806    }
     4807
     4808    if(MivComp(next_weight, ivNull) == 1)
     4809    {
     4810      newRing = currRing;
     4811      delete next_weight;
     4812      break;
     4813    }
     4814    if(MivComp(next_weight, target_weight) == 1)
     4815    {
     4816      endwalks = 1;
     4817    }
     4818    for(i=nV-1; i>=0; i--)
     4819    {
     4820      (*tmp_weight)[i] = (*curr_weight)[i];
     4821      (*curr_weight)[i] = (*next_weight)[i];
     4822    }
     4823    delete next_weight;
     4824  }
     4825  rChangeCurrRing(XXRing);
     4826  G = idrMoveR(G, newRing,currRing);
     4827
     4828  delete tmp_weight;
     4829  delete ivNull;
     4830  delete exivlp;
     4831#ifndef BUCHBERGER_ALG
     4832  delete last_omega;
     4833#endif
     4834#ifdef TIME_TEST
     4835  TimeString(tinput, tostd, tif, tstd, tlift, tred, tnw, nstep);
     4836
     4837  Print("\n// pSetm_Error = (%d)", ErrorCheck());
     4838  Print("\n// Overflow_Error? (%d)\n", Overflow_Error);
     4839#endif
     4840  return(G);
     4841}
     4842
     4843//unused
     4844#if 0
    38224845ideal Mwalk_tst(ideal Go, intvec* curr_weight, intvec* target_weight)
    38234846{
    3824   clock_t tinput=clock();
     4847  //clock_t tinput=clock();
    38254848  //idString(Go,"Ginp");
    38264849  int i, nV = currRing->N;
    38274850  int nwalk=0, endwalks=0;
    38284851
    3829   ideal Gomega, M, F, Gomega1, Gomega2, M1, F1, G, G1;
    3830   ring endRing, newRing, oldRing;
     4852  ideal Gomega, M, F, Gomega1, Gomega2, M1, F1, G;
     4853  // ideal G1; ring endRing;
     4854  ring newRing, oldRing;
    38314855  intvec* ivNull = new intvec(nV);
    38324856  ring XXRing = currRing;
     
    38344858  intvec* tmp_weight = new intvec(nV);
    38354859  for(i=nV-1; i>=0; i--)
     4860  {
    38364861    (*tmp_weight)[i] = (*curr_weight)[i];
    3837 
     4862  }
    38384863  /* the monomial ordering of this current ring would be "dp" */
    38394864  G = MstdCC(Go);
    3840 
     4865#ifndef BUCHBERGER_ALG
    38414866  intvec* hilb_func;
    3842 
     4867#endif
    38434868  /* to avoid (1,0,...,0) as the target vector */
    38444869  intvec* last_omega = new intvec(nV);
     
    39404965  return(G);
    39414966}
    3942 
    3943 
     4967#endif
    39444968
    39454969/**************************************************************/
     
    39554979   2) the changed perturbation walk algorithm with a decreased degree.
    39564980*/
    3957 /* use kStd, if nP = 0, else call LastGB */
     4981// use kStd, if nP = 0, else call LastGB
    39584982ideal Mpwalk(ideal Go, int op_deg, int tp_deg,intvec* curr_weight,
    39594983             intvec* target_weight, int nP)
     
    39714995
    39724996  nstep = 0;
    3973   int i, ntwC=1, ntestw=1,  nV = currRing->N, op_tmp = op_deg;
    3974   int endwalks=0, nhilb=0, ntestomega=0;
     4997  int i, ntwC=1, ntestw=1,  nV = currRing->N;
     4998  int endwalks=0;
    39754999
    39765000  ideal Gomega, M, F, G, Gomega1, Gomega2, M1,F1,Eresult,ssG;
     
    39835007  intvec* ivNull = new intvec(nV);
    39845008  intvec* iv_dp = MivUnit(nV);// define (1,1,...,1)
    3985   intvec* cw_tmp = curr_weight;
     5009#ifndef BUCHBERGER_ALG
    39865010  intvec* hilb_func;
     5011#endif
    39875012  intvec* next_weight;
    3988   intvec* extra_curr_weight = new intvec(nV);
    3989 
    3990   /* to avoid (1,0,...,0) as the target vector */
     5013
     5014  // to avoid (1,0,...,0) as the target vector
    39915015  intvec* last_omega = new intvec(nV);
    39925016  for(i=nV-1; i>0; i--)
     
    40115035  else
    40125036  {
    4013     //ring order := (a(curr_weight),lp);
     5037    //define ring order := (a(curr_weight),lp);
    40145038    if (rParameter(currRing) != NULL)
    40155039      DefRingPar(curr_weight);
     
    40545078    }
    40555079    delete iv_M_lp;
    4056     pert_target_vector = target_weight; //vor 19. mai 2003//test 19 Junu 03
     5080    pert_target_vector = target_weight;
    40575081    rChangeCurrRing(HelpRing);
    40585082    G = idrMoveR(ssG, TargetRing,currRing);
     
    40935117    oldRing = currRing;
    40945118
    4095     /* define a new ring that its ordering is "(a(curr_weight),lp) */
     5119    // define a new ring with ordering "(a(curr_weight),lp)
    40965120    if (rParameter(currRing) != NULL)
    40975121      DefRingPar(curr_weight);
     
    41875211    {
    41885212      ntwC = 0;
    4189       ntestomega = 1;
     5213      //ntestomega = 1;
    41905214      //Print("\n// ring r%d = %s;\n", nstep, rString(currRing));
    41915215      //idElements(G, "G");
     
    41955219    if(MivComp(next_weight, ivNull) == 1){
    41965220      newRing = currRing;
    4197       delete next_weight;//16.03.02
     5221      delete next_weight;
    41985222      //Print("\n// ring r%d = %s;\n", nstep, rString(currRing));
    41995223      break;
     
    42105234  if(tp_deg != 1)
    42115235  {
    4212   FINISH_160302://16.03.02
     5236  FINISH_160302:
    42135237    if(MivSame(orig_target, exivlp) == 1)
    42145238      if (rParameter(currRing) != NULL)
     
    42445268      }
    42455269      */
    4246       /* LastGB is "better" than the kStd subroutine */
     5270      // LastGB is "better" than the kStd subroutine
    42475271      to=clock();
    42485272      ideal eF1;
     
    42965320intvec* XivNull;
    42975321
    4298 /* define a matrix (1 ... 1) */
     5322// define a matrix (1 ... 1)
    42995323intvec* MMatrixone(int nV)
    43005324{
     
    43135337int Xngleich;
    43145338
    4315 /* 27.07.03 */
    4316 /* Perturb the start weight vector at the top level, i.e. nlev = 1 */
     5339/***********************************************************************
     5340 * Perturb the start weight vector at the top level, i.e. nlev = 1     *
     5341 ***********************************************************************/
    43175342static ideal rec_fractal_call(ideal G, int nlev, intvec* omtmp)
    43185343{
     
    43215346
    43225347  int i, nV = currRing->N;
    4323   ring new_ring, testring, extoRing;
     5348  ring new_ring, testring;
     5349  //ring extoRing;
    43245350  ideal Gomega, Gomega1, Gomega2, F, F1, Gresult, Gresult1, G1, Gt;
    43255351  int nwalks = 0;
    43265352  intvec* Mwlp;
     5353#ifndef BUCHBERGER_ALG
    43275354  intvec* hilb_func;
    4328   intvec* extXtau;
     5355#endif
     5356//  intvec* extXtau;
    43295357  intvec* next_vect;
    43305358  intvec* omega2 = new intvec(nV);
    43315359  intvec* altomega = new intvec(nV);
    43325360
    4333   BOOLEAN isnewtarget = FALSE;
    4334 
    4335   /* to avoid (1,0,...,0) as the target vector (Hans) */
     5361  //BOOLEAN isnewtarget = FALSE;
     5362
     5363  // to avoid (1,0,...,0) as the target vector (Hans)
    43365364  intvec* last_omega = new intvec(nV);
    43375365  for(i=nV-1; i>0; i--)
     
    44385466    {
    44395467      delete next_vect;
    4440 
    4441     OVERFLOW_IN_NEXT_VECTOR:
    44425468      if (rParameter(currRing) != NULL)
     5469      {
    44435470        DefRingPar(omtmp);
     5471      }
    44445472      else
     5473      {
    44455474        VMrDefault(omtmp);
    4446 
     5475      }
    44475476#ifdef TEST_OVERFLOW
    44485477      Gt = idrMoveR(G, oRing,currRing);
     
    46545683}
    46555684
     5685/************************************************************************
     5686 * Perturb the start weight vector at the top level with random element *
     5687 ************************************************************************/
     5688static ideal rec_r_fractal_call(ideal G, int nlev, intvec* omtmp, int weight_rad)
     5689{
     5690  Overflow_Error =  FALSE;
     5691  //Print("\n\n// Entering the %d-th recursion:", nlev);
     5692
     5693  int i,k,weight_norm;
     5694  int nV = currRing->N;
     5695  ring new_ring, testring;
     5696  //ring extoRing;
     5697  ideal Gomega, Gomega1, Gomega2, F, F1, Gresult, Gresult1, G1, Gt;
     5698  int nwalks = 0;
     5699  intvec* Mwlp;
     5700#ifndef BUCHBERGER_ALG
     5701  intvec* hilb_func;
     5702#endif
     5703  //intvec* extXtau;
     5704  intvec* next_vect;
     5705  intvec* omega2 = new intvec(nV);
     5706  intvec* altomega = new intvec(nV);
     5707
     5708  //BOOLEAN isnewtarget = FALSE;
     5709
     5710  // to avoid (1,0,...,0) as the target vector (Hans)
     5711  intvec* last_omega = new intvec(nV);
     5712  for(i=nV-1; i>0; i--)
     5713    (*last_omega)[i] = 1;
     5714  (*last_omega)[0] = 10000;
     5715
     5716  intvec* omega = new intvec(nV);
     5717  for(i=0; i<nV; i++) {
     5718    if(Xsigma->length() == nV)
     5719      (*omega)[i] =  (*Xsigma)[i];
     5720    else
     5721      (*omega)[i] = (*Xsigma)[(nV*(nlev-1))+i];
     5722
     5723    (*omega2)[i] = (*Xtau)[(nlev-1)*nV+i];
     5724  }
     5725
     5726   if(nlev == 1)  Xcall = 1;
     5727   else Xcall = 0;
     5728
     5729  ring oRing = currRing;
     5730
     5731  while(1)
     5732  {
     5733#ifdef FIRST_STEP_FRACTAL
     5734    // perturb the current weight vector only on the top level or
     5735    // after perturbation of the both vectors, nlev = 2 as the top level
     5736    if((nlev == 1 && Xcall == 0) || (nlev == 2 && Xngleich == 1))
     5737      if(islengthpoly2(G) == 1)
     5738      {
     5739        Mwlp = MivWeightOrderlp(omega);
     5740        Xsigma = Mfpertvector(G, Mwlp);
     5741        delete Mwlp;
     5742        Overflow_Error = FALSE;
     5743      }
     5744#endif
     5745    nwalks ++;
     5746    NEXT_VECTOR_FRACTAL:
     5747    to=clock();
     5748    /* determine the next border */
     5749    next_vect = MkInterRedNextWeight(omega,omega2,G);
     5750    xtnw=xtnw+clock()-to;
     5751#ifdef PRINT_VECTORS
     5752    MivString(omega, omega2, next_vect);
     5753#endif
     5754    oRing = currRing;
     5755
     5756    /* We only perturb the current target vector at the recursion  level 1 */
     5757    if(Xngleich == 0 && nlev == 1) //(ngleich == 0) important, e.g. ex2, ex3
     5758    {
     5759      if (MivComp(next_vect, omega2) == 1)
     5760      {
     5761        /* to dispense with taking initial (and lifting/interreducing
     5762           after the call of recursion */
     5763        //Print("\n\n// ** Perturb the both vectors with degree %d with",nlev);
     5764        //idElements(G, "G");
     5765
     5766        Xngleich = 1;
     5767        nlev +=1;
     5768
     5769        if (rParameter(currRing) != NULL)
     5770          DefRingPar(omtmp);
     5771        else
     5772          VMrDefault(omtmp);
     5773
     5774        testring = currRing;
     5775        Gt = idrMoveR(G, oRing,currRing);
     5776
     5777        /* perturb the original target vector w.r.t. the current GB */
     5778        delete Xtau;
     5779        Xtau = NewVectorlp(Gt);
     5780
     5781        rChangeCurrRing(oRing);
     5782        G = idrMoveR(Gt, testring,currRing);
     5783
     5784        /* perturb the current vector w.r.t. the current GB */
     5785        Mwlp = MivWeightOrderlp(omega);
     5786        Xsigma = Mfpertvector(G, Mwlp);
     5787        delete Mwlp;
     5788
     5789        for(i=nV-1; i>=0; i--) {
     5790          (*omega2)[i] = (*Xtau)[nV+i];
     5791          (*omega)[i] = (*Xsigma)[nV+i];
     5792        }
     5793
     5794        delete next_vect;
     5795        to=clock();
     5796
     5797        /* to avoid the value of Overflow_Error that occur in Mfpertvector*/
     5798        Overflow_Error = FALSE;
     5799
     5800        next_vect = MkInterRedNextWeight(omega,omega2,G);
     5801        xtnw=xtnw+clock()-to;
     5802
     5803#ifdef PRINT_VECTORS
     5804        MivString(omega, omega2, next_vect);
     5805#endif
     5806      }
     5807      else
     5808      {
     5809        // compute a perturbed next weight vector "next_vect1"
     5810        intvec* next_vect11 = MPertVectors(G, MivMatrixOrder(omega), nlev);
     5811        intvec* next_vect1 = MkInterRedNextWeight(next_vect11, omega2, G);
     5812        // Print("\n//  size of next_vect  = %d", sizeof((*next_vect)));
     5813        // Print("\n//  size of next_vect1  = %d", sizeof((*next_vect1)));
     5814        // Print("\n//  size of next_vect11  = %d", sizeof((*next_vect11)));
     5815        delete next_vect11;
     5816
     5817        // compare next_vect and next_vect1
     5818        ideal G_test = MwalkInitialForm(G, next_vect);
     5819        ideal G_test1 = MwalkInitialForm(G, next_vect1);
     5820        // Print("\n// G_test, G_test 1 erzeugt");
     5821        if(IDELEMS(G_test1) <= IDELEMS(G_test))
     5822          {
     5823          next_vect = ivCopy(next_vect1);
     5824          }
     5825        delete next_vect1;
     5826        // compute a random next weight vector "next_vect2"
     5827        intvec* next_vect22 = ivCopy(omega2);
     5828        // Print("\n//  size of next_vect22  = %d", sizeof((*next_vect22)));
     5829        k = 0;
     5830        while(test_w_in_ConeCC(G, next_vect22) == 0)
     5831        {
     5832          k = k + 1;
     5833          if(k>10)
     5834          {
     5835            break;
     5836          }
     5837          weight_norm = 0;
     5838          while(weight_norm == 0)
     5839          {
     5840            for(i=nV-1; i>=0; i--)
     5841            {
     5842              (*next_vect22)[i] = rand() % 60000 - 30000;
     5843              weight_norm = weight_norm + (*next_vect22)[i]*(*next_vect22)[i];
     5844            }
     5845            weight_norm = 1 + floor(sqrt(weight_norm));
     5846          }
     5847          for(i=nV-1; i>=0; i--)
     5848          {
     5849            if((*next_vect22)[i] < 0)
     5850            {
     5851              (*next_vect22)[i] = 1 + (*omega)[i] + floor(weight_rad*(*next_vect22)[i]/weight_norm);
     5852            }
     5853            else
     5854            {
     5855              (*next_vect22)[i] = (*omega)[i] + floor(weight_rad*(*next_vect22)[i]/weight_norm);
     5856            }
     5857          }
     5858        }
     5859        if(test_w_in_ConeCC(G, next_vect22) == 1)
     5860        {
     5861          //compare next_weight and next_vect2
     5862          intvec* next_vect2 = MkInterRedNextWeight(next_vect22, omega2, G);
     5863          // Print("\n//  size of next_vect2  = %d", sizeof((*next_vect2)));
     5864          ideal G_test2 = MwalkInitialForm(G, next_vect2);
     5865          if(IDELEMS(G_test2) <= IDELEMS(G_test))
     5866          {
     5867            if(IDELEMS(G_test2) <= IDELEMS(G_test1))
     5868            {
     5869              next_vect = ivCopy(next_vect2);
     5870            }
     5871          }
     5872          idDelete(&G_test2);
     5873          delete next_vect2;
     5874        }
     5875        delete next_vect22;
     5876        idDelete(&G_test);
     5877        idDelete(&G_test1);
     5878#ifdef PRINT_VECTORS
     5879        MivString(omega, omega2, next_vect);
     5880#endif
     5881      }
     5882    }
     5883
     5884
     5885    /* check whether the the computed vector is in the correct cone */
     5886    /* If no, the reduced GB of an omega-homogeneous ideal will be
     5887       computed by Buchberger algorithm and stop this recursion step*/
     5888    //if(test_w_in_ConeCC(G, next_vect) != 1) //e.g. Example s7, cyc6
     5889    if(Overflow_Error == TRUE)
     5890    {
     5891      delete next_vect;
     5892
     5893    //OVERFLOW_IN_NEXT_VECTOR:
     5894      if (rParameter(currRing) != NULL)
     5895        DefRingPar(omtmp);
     5896      else
     5897        VMrDefault(omtmp);
     5898
     5899#ifdef TEST_OVERFLOW
     5900      Gt = idrMoveR(G, oRing,currRing);
     5901      Gt = NULL; return(Gt);
     5902#endif
     5903
     5904      //Print("\n\n// apply BB's alg. in ring r = %s;", rString(currRing));
     5905      to=clock();
     5906      Gt = idrMoveR(G, oRing,currRing);
     5907      G1 = MstdCC(Gt);
     5908      xtextra=xtextra+clock()-to;
     5909      Gt = NULL;
     5910
     5911      delete omega2;
     5912      delete altomega;
     5913
     5914      //Print("\n// Leaving the %d-th recursion with %d steps", nlev, nwalks);
     5915      //Print("  ** Overflow_Error? (%d)", Overflow_Error);
     5916      nnflow ++;
     5917
     5918      Overflow_Error = FALSE;
     5919      return (G1);
     5920    }
     5921
     5922
     5923    /* If the perturbed target vector stays in the correct cone,
     5924       return the current GB,
     5925       otherwise, return the computed  GB by the Buchberger-algorithm.
     5926       Then we update the perturbed target vectors w.r.t. this GB. */
     5927
     5928    /* the computed vector is equal to the origin vector, since
     5929       t is not defined */
     5930    if (MivComp(next_vect, XivNull) == 1)
     5931    {
     5932      if (rParameter(currRing) != NULL)
     5933        DefRingPar(omtmp);
     5934      else
     5935        VMrDefault(omtmp);
     5936
     5937      testring = currRing;
     5938      Gt = idrMoveR(G, oRing,currRing);
     5939
     5940      if(test_w_in_ConeCC(Gt, omega2) == 1) {
     5941        delete omega2;
     5942        delete next_vect;
     5943        delete altomega;
     5944        //Print("\n// Leaving the %d-th recursion with %d steps ",nlev, nwalks);
     5945        //Print(" ** Overflow_Error? (%d)", Overflow_Error);
     5946
     5947        return (Gt);
     5948      }
     5949      else
     5950      {
     5951        //ivString(omega2, "tau'");
     5952        //Print("\n//  tau' doesn't stay in the correct cone!!");
     5953
     5954#ifndef  MSTDCC_FRACTAL
     5955        //07.08.03
     5956        //ivString(Xtau, "old Xtau");
     5957        intvec* Xtautmp = Mfpertvector(Gt, MivMatrixOrder(omtmp));
     5958#ifdef TEST_OVERFLOW
     5959      if(Overflow_Error == TRUE)
     5960      Gt = NULL; return(Gt);
     5961#endif
     5962
     5963        if(MivSame(Xtau, Xtautmp) == 1)
     5964        {
     5965          //PrintS("\n// Update vectors are equal to the old vectors!!");
     5966          delete Xtautmp;
     5967          goto FRACTAL_MSTDCC;
     5968        }
     5969
     5970        Xtau = Xtautmp;
     5971        Xtautmp = NULL;
     5972        //ivString(Xtau, "new  Xtau");
     5973
     5974        for(i=nV-1; i>=0; i--)
     5975          (*omega2)[i] = (*Xtau)[(nlev-1)*nV+i];
     5976
     5977        //Print("\n//  ring tau = %s;", rString(currRing));
     5978        rChangeCurrRing(oRing);
     5979        G = idrMoveR(Gt, testring,currRing);
     5980
     5981        goto NEXT_VECTOR_FRACTAL;
     5982#endif
     5983
     5984      FRACTAL_MSTDCC:
     5985        //Print("\n//  apply BB-Alg in ring = %s;", rString(currRing));
     5986        to=clock();
     5987        G = MstdCC(Gt);
     5988        xtextra=xtextra+clock()-to;
     5989
     5990        oRing = currRing;
     5991
     5992        // update the original target vector w.r.t. the current GB
     5993        if(MivSame(Xivinput, Xivlp) == 1)
     5994          if (rParameter(currRing) != NULL)
     5995            DefRingParlp();
     5996          else
     5997            VMrDefaultlp();
     5998        else
     5999          if (rParameter(currRing) != NULL)
     6000            DefRingPar(Xivinput);
     6001          else
     6002            VMrDefault(Xivinput);
     6003
     6004        testring = currRing;
     6005        Gt = idrMoveR(G, oRing,currRing);
     6006
     6007        delete Xtau;
     6008        Xtau = NewVectorlp(Gt);
     6009
     6010        rChangeCurrRing(oRing);
     6011        G = idrMoveR(Gt, testring,currRing);
     6012
     6013        delete omega2;
     6014        delete next_vect;
     6015        delete altomega;
     6016        /*
     6017          Print("\n// Leaving the %d-th recursion with %d steps,", nlev,nwalks);
     6018          Print(" ** Overflow_Error? (%d)", Overflow_Error);
     6019        */
     6020        if(Overflow_Error == TRUE)
     6021          nnflow ++;
     6022
     6023        Overflow_Error = FALSE;
     6024        return(G);
     6025      }
     6026    }
     6027
     6028    for(i=nV-1; i>=0; i--) {
     6029      (*altomega)[i] = (*omega)[i];
     6030      (*omega)[i] = (*next_vect)[i];
     6031    }
     6032    delete next_vect;
     6033
     6034    to=clock();
     6035    /* Take the initial form of <G> w.r.t. omega */
     6036    Gomega = MwalkInitialForm(G, omega);
     6037    xtif=xtif+clock()-to;
     6038
     6039#ifndef  BUCHBERGER_ALG
     6040    if(isNolVector(omega) == 0)
     6041      hilb_func = hFirstSeries(Gomega,NULL,NULL,omega,currRing);
     6042    else
     6043      hilb_func = hFirstSeries(Gomega,NULL,NULL,last_omega,currRing);
     6044#endif // BUCHBERGER_ALG
     6045
     6046    if (rParameter(currRing) != NULL)
     6047      DefRingPar(omega);
     6048    else
     6049      VMrDefault(omega);
     6050
     6051    Gomega1 = idrMoveR(Gomega, oRing,currRing);
     6052
     6053    /* Maximal recursion depth, to compute a red. GB */
     6054    /* Fractal walk with the alternative recursion */
     6055    /* alternative recursion */
     6056    // if(nlev == nV || lengthpoly(Gomega1) == 0)
     6057    if(nlev == Xnlev || lengthpoly(Gomega1) == 0)
     6058      //if(nlev == nV) // blind recursion
     6059    {
     6060      /*
     6061      if(Xnlev != nV)
     6062      {
     6063        Print("\n// ** Xnlev = %d", Xnlev);
     6064        ivString(Xtau, "Xtau");
     6065      }
     6066      */
     6067      to=clock();
     6068#ifdef  BUCHBERGER_ALG
     6069      Gresult = MstdhomCC(Gomega1);
     6070#else
     6071      Gresult =kStd(Gomega1,NULL,isHomog,NULL,hilb_func,0,NULL,omega);
     6072      delete hilb_func;
     6073#endif // BUCHBERGER_ALG
     6074      xtstd=xtstd+clock()-to;
     6075    }
     6076    else {
     6077      rChangeCurrRing(oRing);
     6078      Gomega1 = idrMoveR(Gomega1, oRing,currRing);
     6079      Gresult = rec_fractal_call(idCopy(Gomega1),nlev+1,omega);
     6080    }
     6081
     6082    //convert a Groebner basis from a ring to another ring,
     6083    new_ring = currRing;
     6084
     6085    rChangeCurrRing(oRing);
     6086    Gresult1 = idrMoveR(Gresult, new_ring,currRing);
     6087    Gomega2 = idrMoveR(Gomega1, new_ring,currRing);
     6088
     6089    to=clock();
     6090    /* Lifting process */
     6091    F = MLifttwoIdeal(Gomega2, Gresult1, G);
     6092    xtlift=xtlift+clock()-to;
     6093    idDelete(&Gresult1);
     6094    idDelete(&Gomega2);
     6095    idDelete(&G);
     6096
     6097    rChangeCurrRing(new_ring);
     6098    F1 = idrMoveR(F, oRing,currRing);
     6099
     6100    to=clock();
     6101    /* Interreduce G */
     6102    G = kInterRedCC(F1, NULL);
     6103    xtred=xtred+clock()-to;
     6104    idDelete(&F1);
     6105  }
     6106}
     6107
     6108
     6109
     6110/*******************************************************************************
     6111 * The implementation of the fractal walk algorithm                            *
     6112 *                                                                             *
     6113 * The main procedur Mfwalk calls the recursive Subroutine                     *
     6114 * rec_fractal_call to compute the wanted Grï¿œbner basis.                       *
     6115 * At the main procedur we compute the reduced Grï¿œbner basis w.r.t. a "fast"   *
     6116 * order, e.g. "dp" and a sequence of weight vectors which are row vectors     *
     6117 * of a matrix. This matrix defines the given monomial order, e.g. "lp"        *
     6118 *******************************************************************************/
    46566119ideal Mfwalk(ideal G, intvec* ivstart, intvec* ivtarget)
    46576120{
     
    46846147  for(i=IDELEMS(Gw)-1; i>=0; i--)
    46856148  {
    4686     if((Gw->m[i]!=NULL) /* len >=0 */
    4687     && (Gw->m[i]->next!=NULL) /* len >=1 */
    4688     && (Gw->m[i]->next->next!=NULL)) /* len >=2 */
    4689     {
    4690       intvec* iv_dp = MivUnit(nV);// define (1,1,...,1)
     6149    if((Gw->m[i]!=NULL) // len >=0
     6150    && (Gw->m[i]->next!=NULL) // len >=1
     6151    && (Gw->m[i]->next->next!=NULL)) // len >=2
     6152    {
     6153      intvec* iv_dp = MivUnit(nV); // define (1,1,...,1)
    46916154      intvec* Mdp;
    46926155
     
    47806243}
    47816244
    4782 /* Tran algorithm */
    4783 /* use kStd, if nP = 0, else call Ab_Rec_Pert (LastGB) */
     6245ideal Mfrwalk(ideal G, intvec* ivstart, intvec* ivtarget,int weight_rad)
     6246{
     6247  Set_Error(FALSE);
     6248  Overflow_Error = FALSE;
     6249  //Print("// pSetm_Error = (%d)", ErrorCheck());
     6250  //Print("\n// ring ro = %s;", rString(currRing));
     6251
     6252  nnflow = 0;
     6253  Xngleich = 0;
     6254  Xcall = 0;
     6255  xtif=0; xtstd=0; xtlift=0; xtred=0; xtnw=0; xtextra=0;
     6256  xftinput = clock();
     6257
     6258  ring  oldRing = currRing;
     6259  int i, nV = currRing->N;
     6260  XivNull = new intvec(nV);
     6261  Xivinput = ivtarget;
     6262  ngleich = 0;
     6263  to=clock();
     6264  ideal I = MstdCC(G);
     6265  G = NULL;
     6266  xftostd=clock()-to;
     6267  Xsigma = ivstart;
     6268
     6269  Xnlev=nV;
     6270
     6271#ifdef FIRST_STEP_FRACTAL
     6272  ideal Gw = MwalkInitialForm(I, ivstart);
     6273  for(i=IDELEMS(Gw)-1; i>=0; i--)
     6274  {
     6275    if((Gw->m[i]!=NULL) // len >=0
     6276    && (Gw->m[i]->next!=NULL) // len >=1
     6277    && (Gw->m[i]->next->next!=NULL)) // len >=2
     6278    {
     6279      intvec* iv_dp = MivUnit(nV); // define (1,1,...,1)
     6280      intvec* Mdp;
     6281
     6282      if(MivSame(ivstart, iv_dp) != 1)
     6283      {
     6284        Mdp = MivWeightOrderdp(ivstart);
     6285      }
     6286      else
     6287      {
     6288        Mdp = MivMatrixOrderdp(nV);
     6289      }
     6290      Xsigma = Mfpertvector(I, Mdp);
     6291      Overflow_Error = FALSE;
     6292
     6293      delete Mdp;
     6294      delete iv_dp;
     6295      break;
     6296    }
     6297  }
     6298  idDelete(&Gw);
     6299#endif
     6300
     6301  ideal I1;
     6302  intvec* Mlp;
     6303  Xivlp = Mivlp(nV);
     6304
     6305  if(MivComp(ivtarget, Xivlp)  != 1)
     6306  {
     6307    if (rParameter(currRing) != NULL)
     6308      DefRingPar(ivtarget);
     6309    else
     6310      VMrDefault(ivtarget);
     6311
     6312    I1 = idrMoveR(I, oldRing,currRing);
     6313    Mlp = MivWeightOrderlp(ivtarget);
     6314    Xtau = Mfpertvector(I1, Mlp);
     6315  }
     6316  else
     6317  {
     6318    if (rParameter(currRing) != NULL)
     6319      DefRingParlp();
     6320    else
     6321      VMrDefaultlp();
     6322
     6323    I1 = idrMoveR(I, oldRing,currRing);
     6324    Mlp =  MivMatrixOrderlp(nV);
     6325    Xtau = Mfpertvector(I1, Mlp);
     6326  }
     6327  delete Mlp;
     6328  Overflow_Error = FALSE;
     6329
     6330  //ivString(Xsigma, "Xsigma");
     6331  //ivString(Xtau, "Xtau");
     6332
     6333  id_Delete(&I, oldRing);
     6334  ring tRing = currRing;
     6335
     6336  if (rParameter(currRing) != NULL)
     6337    DefRingPar(ivstart);
     6338  else
     6339    VMrDefault(ivstart);
     6340
     6341  I = idrMoveR(I1,tRing,currRing);
     6342  to=clock();
     6343  ideal J = MstdCC(I);
     6344  idDelete(&I);
     6345  xftostd=xftostd+clock()-to;
     6346
     6347  ideal resF;
     6348  ring helpRing = currRing;
     6349  J = rec_r_fractal_call(J,1,ivtarget,weight_rad);
     6350  rChangeCurrRing(oldRing);
     6351  resF = idrMoveR(J, helpRing,currRing);
     6352  idSkipZeroes(resF);
     6353
     6354  delete Xivlp;
     6355  delete Xsigma;
     6356  delete Xtau;
     6357  delete XivNull;
     6358
     6359#ifdef TIME_TEST
     6360  TimeStringFractal(xftinput, xftostd, xtif, xtstd, xtextra,
     6361                    xtlift, xtred, xtnw);
     6362
     6363
     6364  Print("\n// pSetm_Error = (%d)", ErrorCheck());
     6365  Print("\n// Overflow_Error? (%d)\n", Overflow_Error);
     6366  Print("\n// the numbers of Overflow_Error (%d)", nnflow);
     6367#endif
     6368
     6369  return(resF);
     6370}
     6371
     6372/*******************************************************
     6373 * Tran's algorithm                                    *
     6374 *                                                     *
     6375 * use kStd, if nP = 0, else call Ab_Rec_Pert (LastGB) *
     6376 *******************************************************/
    47846377ideal TranMImprovwalk(ideal G,intvec* curr_weight,intvec* target_tmp, int nP)
    47856378{
     6379#ifdef TIME_TEST
    47866380  clock_t mtim = clock();
     6381#endif
    47876382  Set_Error(FALSE  );
    47886383  Overflow_Error =  FALSE;
     
    47916386
    47926387  clock_t tostd, tif=0, tstd=0, tlift=0, tred=0, tnw=0, textra=0;
     6388#ifdef TIME_TEST
    47936389  clock_t tinput = clock();
    4794 
     6390#endif
    47956391  int nsteppert=0, i, nV = currRing->N, nwalk=0, npert_tmp=0;
    47966392  int *npert=(int*)omAlloc(2*nV*sizeof(int));
    47976393  ideal Gomega, M,F,  G1, Gomega1, Gomega2, M1, F1;
    4798   ring endRing, newRing, oldRing, lpRing;
     6394  //ring endRing;
     6395  ring newRing, oldRing, lpRing;
    47996396  intvec* next_weight;
    48006397  intvec* ivNull = new intvec(nV); //define (0,...,0)
    48016398  intvec* iv_dp = MivUnit(nV);// define (1,1,...,1)
    48026399  intvec* iv_lp = Mivlp(nV); //define (1,0,...,0)
    4803   ideal H0, H1, H2, Glp;
     6400  ideal H0;
     6401  //ideal  H1;
     6402  ideal H2, Glp;
    48046403  int nGB, endwalks = 0,  nwalkpert=0,  npertstep=0;
    48056404  intvec* Mlp =  MivMatrixOrderlp(nV);
    48066405  intvec* vector_tmp = new intvec(nV);
     6406#ifndef BUCHBERGER_ALG
    48076407  intvec* hilb_func;
    4808 
     6408#endif
    48096409  /* to avoid (1,0,...,0) as the target vector */
    48106410  intvec* last_omega = new intvec(nV);
     
    49396539    nwalkpert++;
    49406540    to=clock();
    4941     /* compute a next weight vector */
     6541    // compute a next weight vector
    49426542    next_weight = MwalkNextWeightCC(curr_weight,target_weight, G);
    49436543    tnw=tnw+clock()-to;
     
    49476547
    49486548
    4949     /* check whether the computed intermediate weight vector in
    4950        the correct cone is, since sometimes it is very big e.g. s7, cyc7.
    4951        If it is NOT in the cone, then one has directly to compute
     6549    /* check whether the computed intermediate weight vector is in
     6550       the correct cone; sometimes it is very big e.g. s7, cyc7.
     6551       If it is NOT in the correct cone, then compute directly
    49526552       a reduced Groebner basis with respect to the lexicographic ordering
    49536553       for the known Groebner basis that it is computed in the last step.
     
    50056605    }
    50066606
    5007     /* check whether the computed Groebner basis a really Groebner basis is.
    5008        if no, we perturb the target vector with the maximal "perturbation"
     6607    /* check whether the computed Groebner basis is really a Groebner basis.
     6608       If not, we perturb the target vector with the maximal "perturbation"
    50096609       degree.*/
    50106610    if(MivComp(next_weight, target_weight) == 1 ||
     
    51446744    delete next_weight;
    51456745  }//while
    5146 
     6746#ifdef TEST_OVERFLOW
    51476747 BE_FINISH:
     6748#endif
    51486749  rChangeCurrRing(XXRing);
    51496750  G = idrMoveR(G, lpRing,currRing);
     
    51686769}
    51696770
    5170 
    5171 /* compute the reduced Grï¿œbner basis of an ideal <Go> w.r.t. lp */
     6771/*******************************************************
     6772 * Tran's algorithm with random element                *
     6773 *                                                     *
     6774 * use kStd, if nP = 0, else call Ab_Rec_Pert (LastGB) *
     6775 *******************************************************/
     6776ideal TranMrImprovwalk(ideal G,intvec* curr_weight,intvec* target_tmp, int nP, int weight_rad, int pert_deg)
     6777{
     6778#ifdef TIME_TEST
     6779  clock_t mtim = clock();
     6780#endif
     6781  Set_Error(FALSE  );
     6782  Overflow_Error =  FALSE;
     6783  //Print("// pSetm_Error = (%d)", ErrorCheck());
     6784  //Print("\n// ring ro = %s;", rString(currRing));
     6785
     6786  clock_t tostd, tif=0, tstd=0, tlift=0, tred=0, tnw=0, textra=0;
     6787#ifdef TIME_TEST
     6788  clock_t tinput = clock();
     6789#endif
     6790  int nsteppert=0, i, k, nV = currRing->N, nwalk=0, npert_tmp=0;
     6791  int *npert=(int*)omAlloc(2*nV*sizeof(int));
     6792  ideal Gomega, M,F,  G1, Gomega1, Gomega2, M1, F1;
     6793  //ring endRing;
     6794  ring newRing, oldRing, lpRing;
     6795  intvec* next_weight;
     6796  intvec* ivNull = new intvec(nV); //define (0,...,0)
     6797  intvec* iv_dp = MivUnit(nV);// define (1,1,...,1)
     6798  intvec* iv_lp = Mivlp(nV); //define (1,0,...,0)
     6799  ideal H0;
     6800  //ideal H1;
     6801  ideal H2, Glp;
     6802  int weight_norm, nGB, endwalks = 0,  nwalkpert=0,  npertstep=0;
     6803  intvec* Mlp =  MivMatrixOrderlp(nV);
     6804  intvec* vector_tmp = new intvec(nV);
     6805#ifndef BUCHBERGER_ALG
     6806  intvec* hilb_func;
     6807#endif
     6808  // to avoid (1,0,...,0) as the target vector
     6809  intvec* last_omega = new intvec(nV);
     6810  for(i=nV-1; i>0; i--)
     6811  {
     6812    (*last_omega)[i] = 1;
     6813  }
     6814  (*last_omega)[0] = 10000;
     6815
     6816//intvec* extra_curr_weight = new intvec(nV);
     6817  intvec* target_weight = new intvec(nV);
     6818  for(i=nV-1; i>=0; i--)
     6819  {
     6820    (*target_weight)[i] = (*target_tmp)[i];
     6821  }
     6822  ring XXRing = currRing;
     6823  newRing = currRing;
     6824
     6825  to=clock();
     6826  // compute a red. GB w.r.t. the help ring
     6827  if(MivComp(curr_weight, iv_dp) == 1)
     6828  {
     6829    //rOrdStr(currRing) = "dp"
     6830    G = MstdCC(G);
     6831  }
     6832  else
     6833  {
     6834    //rOrdStr(currRing) = (a(.c_w..),lp,C)
     6835    if (rParameter(currRing) != NULL)
     6836    {
     6837      DefRingPar(curr_weight);
     6838    }
     6839    else
     6840    {
     6841      VMrDefault(curr_weight);
     6842    }
     6843    G = idrMoveR(G, XXRing,currRing);
     6844    G = MstdCC(G);
     6845  }
     6846  tostd=clock()-to;
     6847
     6848#ifdef REPRESENTATION_OF_SIGMA
     6849  ideal Gw = MwalkInitialForm(G, curr_weight);
     6850
     6851  if(islengthpoly2(Gw)==1)
     6852  {
     6853    intvec* MDp;
     6854    if(MivComp(curr_weight, iv_dp) == 1)
     6855    {
     6856      MDp = MatrixOrderdp(nV); //MivWeightOrderlp(iv_dp);
     6857    }
     6858    else
     6859    {
     6860      MDp = MivWeightOrderlp(curr_weight);
     6861    }
     6862    curr_weight = RepresentationMatrix_Dp(G, MDp);
     6863
     6864    delete MDp;
     6865
     6866    ring exring = currRing;
     6867
     6868    if (rParameter(currRing) != NULL)
     6869    {
     6870      DefRingPar(curr_weight);
     6871    }
     6872    else
     6873    {
     6874      VMrDefault(curr_weight);
     6875    }
     6876    to=clock();
     6877    Gw = idrMoveR(G, exring,currRing);
     6878    G = MstdCC(Gw);
     6879    Gw = NULL;
     6880    tostd=tostd+clock()-to;
     6881    //ivString(curr_weight,"rep. sigma");
     6882    goto COMPUTE_NEW_VECTOR;
     6883  }
     6884
     6885  idDelete(&Gw);
     6886  delete iv_dp;
     6887#endif
     6888
     6889
     6890  while(1)
     6891  {
     6892    to=clock();
     6893    // compute an initial form ideal of <G> w.r.t. "curr_vector"
     6894    Gomega = MwalkInitialForm(G, curr_weight);
     6895    tif=tif+clock()-to;
     6896
     6897#ifndef  BUCHBERGER_ALG
     6898    if(isNolVector(curr_weight) == 0)
     6899    {
     6900      hilb_func = hFirstSeries(Gomega,NULL,NULL,curr_weight,currRing);
     6901    }
     6902    else
     6903    {
     6904      hilb_func = hFirstSeries(Gomega,NULL,NULL,last_omega,currRing);
     6905    }
     6906#endif // BUCHBERGER_ALG
     6907
     6908    oldRing = currRing;
     6909
     6910    // define a new ring with ordering "(a(curr_weight),lp)
     6911    if (rParameter(currRing) != NULL)
     6912    {
     6913      DefRingPar(curr_weight);
     6914    }
     6915    else
     6916    {
     6917      VMrDefault(curr_weight);
     6918    }
     6919    newRing = currRing;
     6920    Gomega1 = idrMoveR(Gomega, oldRing,currRing);
     6921
     6922    to=clock();
     6923    // compute a reduced Groebner basis of <Gomega> w.r.t. "newRing"
     6924#ifdef  BUCHBERGER_ALG
     6925    M = MstdhomCC(Gomega1);
     6926#else
     6927    M=kStd(Gomega1,NULL,isHomog,NULL,hilb_func,0,NULL,curr_weight);
     6928    delete hilb_func;
     6929#endif
     6930    tstd=tstd+clock()-to;
     6931
     6932    // change the ring to oldRing
     6933    rChangeCurrRing(oldRing);
     6934    M1 =  idrMoveR(M, newRing,currRing);
     6935    Gomega2 =  idrMoveR(Gomega1, newRing,currRing);
     6936
     6937    to=clock();
     6938    // compute a representation of the generators of submod (M) with respect to those of mod (Gomega).
     6939    // Gomega is a reduced Groebner basis w.r.t. the current ring
     6940    F = MLifttwoIdeal(Gomega2, M1, G);
     6941    tlift=tlift+clock()-to;
     6942
     6943    idDelete(&M1);
     6944    idDelete(&Gomega2);
     6945    idDelete(&G);
     6946
     6947    // change the ring to newRing
     6948    rChangeCurrRing(newRing);
     6949    F1 = idrMoveR(F, oldRing,currRing);
     6950
     6951    to=clock();
     6952    // reduce the Groebner basis <G> w.r.t. new ring
     6953    G = kInterRedCC(F1, NULL);
     6954    tred=tred+clock()-to;
     6955    idDelete(&F1);
     6956
     6957  COMPUTE_NEW_VECTOR:
     6958    newRing = currRing;
     6959    nwalk++;
     6960    nwalkpert++;
     6961    to=clock();
     6962    // compute a next weight vector
     6963    //next_weight = MwalkNextWeightCC(curr_weight,target_weight, G);
     6964    next_weight = MWalkRandomNextWeight(G, curr_weight, target_weight, weight_rad, pert_deg);
     6965/*
     6966    next_weight = MkInterRedNextWeight(curr_weight,target_weight,G);
     6967
     6968    if(MivComp(next_weight, target_weight) != 1)
     6969    {
     6970      // compute a perturbed next weight vector "next_weight1"
     6971      intvec* next_weight1 = MkInterRedNextWeight(MPertVectors(G, MivMatrixOrder(curr_weight), pert_deg), target_weight, G);
     6972   
     6973      // compare next_weight and next_weight1
     6974      ideal G_test = MwalkInitialForm(G, next_weight);
     6975      ideal G_test1 = MwalkInitialForm(G, next_weight1);
     6976      if(IDELEMS(G_test1) <= IDELEMS(G_test))
     6977      {
     6978        next_weight = ivCopy(next_weight1);
     6979      }
     6980      delete next_weight1;
     6981      // compute a random next weight vector "next_weight2"
     6982      intvec* next_weight22 = ivCopy(target_weight);
     6983      // Print("\n//  size of target_weight  = %d", sizeof((*target_weight)));
     6984      k = 0;
     6985   
     6986      while(test_w_in_ConeCC(G, next_weight22) == 0 && k < 11)
     6987      {
     6988        k++;
     6989        if(k>10)
     6990        {
     6991          break;
     6992        }
     6993        weight_norm = 0;
     6994        while(weight_norm == 0)
     6995        {
     6996          for(i=nV-1; i>=0; i--)
     6997          {
     6998            // Print("\n//  next_weight[%d]  = %d", i, (*next_weight)[i]);
     6999            (*next_weight22)[i] = rand() % 60000 - 30000;
     7000            weight_norm = weight_norm + (*next_weight22)[i]*(*next_weight22)[i];
     7001          }
     7002          weight_norm = 1 + floor(sqrt(weight_norm));
     7003        }
     7004        for(i=nV-1; i>=0; i--)
     7005        {
     7006          if((*next_weight22)[i] < 0)
     7007          {
     7008            (*next_weight22)[i] = 1 + (*curr_weight)[i] + floor(weight_rad*(*next_weight22)[i]/weight_norm);
     7009          }
     7010          else
     7011          {
     7012            (*next_weight22)[i] = (*curr_weight)[i] + floor(weight_rad*(*next_weight22)[i]/weight_norm);
     7013          }
     7014          // Print("\n//  next_weight22[%d]  = %d", i, (*next_weight22)[i]);
     7015        }
     7016      }
     7017
     7018      if(test_w_in_ConeCC(G, next_weight22) == 1)
     7019      {
     7020        // compare next_weight and next_weight2
     7021        // Print("\n// ZUFALL IM KEGEL");
     7022        intvec* next_weight2 = MkInterRedNextWeight(next_weight22, target_weight, G);
     7023
     7024        ideal G_test2 = MwalkInitialForm(G, next_weight2);
     7025        if(IDELEMS(G_test2) <= IDELEMS(G_test))
     7026        {
     7027          if(IDELEMS(G_test2) <= IDELEMS(G_test1))
     7028          {
     7029             // Print("\n// ZUFALL BENUTZT!\n");
     7030            next_weight = ivCopy(next_weight2);
     7031          }
     7032        }
     7033        idDelete(&G_test2);
     7034        delete next_weight2;
     7035      }
     7036      delete next_weight22;
     7037      idDelete(&G_test);
     7038      idDelete(&G_test1);
     7039    }*/
     7040
     7041    tnw=tnw+clock()-to;
     7042#ifdef PRINT_VECTORS
     7043    MivString(curr_weight, target_weight, next_weight);
     7044#endif
     7045
     7046/*   check whether the computed intermediate weight vector is in
     7047     the correct cone; sometimes it is very big e.g. s7, cyc7.
     7048     If it is NOT in the correct cone, then compute directly
     7049     a reduced Groebner basis with respect to the lexicographic ordering
     7050     for the known Groebner basis that it is computed in the last step.
     7051*/
     7052    //if(test_w_in_ConeCC(G, next_weight) != 1)
     7053    if(Overflow_Error == TRUE)
     7054    {
     7055    OMEGA_OVERFLOW_TRAN_NEW:
     7056      //Print("\n//  takes %d steps!", nwalk-1);
     7057      //Print("\n//ring lastRing = %s;", rString(currRing));
     7058#ifdef TEST_OVERFLOW
     7059      goto  BE_FINISH;
     7060#endif
     7061
     7062#ifdef CHECK_IDEAL_MWALK
     7063      idElements(G, "G");
     7064      //headidString(G, "G");
     7065#endif
     7066
     7067      if(MivSame(target_tmp, iv_lp) == 1)
     7068      {
     7069        if (rParameter(currRing) != NULL)
     7070        {
     7071          DefRingParlp();
     7072        }
     7073        else
     7074        {
     7075          VMrDefaultlp();
     7076        }
     7077      }
     7078      else
     7079      {
     7080        if (rParameter(currRing) != NULL)
     7081        {
     7082          DefRingPar(target_tmp);
     7083        }
     7084        else
     7085        {
     7086          VMrDefault(target_tmp);
     7087        }
     7088      }
     7089      lpRing = currRing;
     7090      G1 = idrMoveR(G, newRing,currRing);
     7091
     7092      to=clock();
     7093      // apply kStd or LastGB to compute  a lex. red. Groebner basis of <G>
     7094      if(nP == 0 || MivSame(target_tmp, iv_lp) == 0)
     7095      {
     7096        //Print("\n\n// calls \"std in ring r_%d = %s;", nwalk, rString(currRing));
     7097        G = MstdCC(G1);//no result for qnt1
     7098      }
     7099      else
     7100      {
     7101        rChangeCurrRing(newRing);
     7102        G1 = idrMoveR(G1, lpRing,currRing);
     7103
     7104        //Print("\n\n// calls \"LastGB\" (%d) to compute a GB", nV-1);
     7105        G = LastGB(G1, curr_weight, nV-1); //no result for kats7
     7106
     7107        rChangeCurrRing(lpRing);
     7108        G = idrMoveR(G, newRing,currRing);
     7109      }
     7110      textra=clock()-to;
     7111      npert[endwalks]=nwalk-npert_tmp;
     7112      npert_tmp = nwalk;
     7113      endwalks ++;
     7114      break;
     7115    }
     7116
     7117    // check whether the computed Groebner basis is really a Groebner basis.
     7118    // If not, we perturb the target vector with the maximal "perturbation" degree.
     7119
     7120    if(MivComp(next_weight, target_weight) == 1 || MivComp(next_weight, curr_weight) == 1 )
     7121    {
     7122      //Print("\n//ring r_%d = %s;", nwalk, rString(currRing));
     7123
     7124
     7125      //compute the number of perturbations and its step
     7126      npert[endwalks]=nwalk-npert_tmp;
     7127      npert_tmp = nwalk;
     7128
     7129      endwalks ++;
     7130
     7131      // it is very important if the walk only uses one step, e.g. Fate, liu
     7132      if(endwalks == 1 && MivComp(next_weight, curr_weight) == 1)
     7133      {
     7134        rChangeCurrRing(XXRing);
     7135        G = idrMoveR(G, newRing,currRing);
     7136        goto FINISH;
     7137      }
     7138      H0 = id_Head(G,currRing);
     7139
     7140      if(MivSame(target_tmp, iv_lp) == 1)
     7141      {
     7142        if (rParameter(currRing) != NULL)
     7143        {
     7144          DefRingParlp();
     7145        }
     7146        else
     7147        {
     7148          VMrDefaultlp();
     7149        }
     7150      }
     7151      else
     7152      {
     7153        if (rParameter(currRing) != NULL)
     7154        {
     7155          DefRingPar(target_tmp);
     7156        }
     7157        else
     7158        {
     7159          VMrDefault(target_tmp);
     7160        }
     7161      }
     7162      lpRing = currRing;
     7163      Glp = idrMoveR(G, newRing,currRing);
     7164      H2 = idrMoveR(H0, newRing,currRing);
     7165
     7166      // Apply Lemma 2.2 in Collart et. al (1997) to check whether cone(k-1) is equal to cone(k)
     7167      nGB = 1;
     7168      for(i=IDELEMS(Glp)-1; i>=0; i--)
     7169      {
     7170        poly t;
     7171        if((t=pSub(pHead(Glp->m[i]), pCopy(H2->m[i]))) != NULL)
     7172        {
     7173          pDelete(&t);
     7174          idDelete(&H2);//5.5.02
     7175          nGB = 0; //i.e. Glp is no reduced Groebner basis
     7176          break;
     7177        }
     7178        pDelete(&t);
     7179      }
     7180
     7181      idDelete(&H2);//5.5.02
     7182
     7183      if(nGB == 1)
     7184      {
     7185        G = Glp;
     7186        Glp = NULL;
     7187        break;
     7188      }
     7189
     7190       // perturb the target weight vector, if the vector target_tmp stays in many cones
     7191      poly p;
     7192      BOOLEAN plength3 = FALSE;
     7193      for(i=IDELEMS(Glp)-1; i>=0; i--)
     7194      {
     7195        p = MpolyInitialForm(Glp->m[i], target_tmp);
     7196        if(p->next != NULL &&
     7197           p->next->next != NULL &&
     7198           p->next->next->next != NULL)
     7199        {
     7200          Overflow_Error = FALSE;
     7201
     7202          for(i=0; i<nV; i++)
     7203          {
     7204            (*vector_tmp)[i] = (*target_weight)[i];
     7205          }
     7206          delete target_weight;
     7207          target_weight = MPertVectors(Glp, Mlp, nV);
     7208
     7209          if(MivComp(vector_tmp, target_weight)==1)
     7210          {
     7211            //PrintS("\n// The old and new representaion vector are the same!!");
     7212            G = Glp;
     7213            newRing = currRing;
     7214            goto OMEGA_OVERFLOW_TRAN_NEW;
     7215           }
     7216
     7217          if(Overflow_Error == TRUE)
     7218          {
     7219            rChangeCurrRing(newRing);
     7220            G = idrMoveR(Glp, lpRing,currRing);
     7221            goto OMEGA_OVERFLOW_TRAN_NEW;
     7222          }
     7223
     7224          plength3 = TRUE;
     7225          pDelete(&p);
     7226          break;
     7227        }
     7228        pDelete(&p);
     7229      }
     7230
     7231      if(plength3 == FALSE)
     7232      {
     7233        rChangeCurrRing(newRing);
     7234        G = idrMoveR(Glp, lpRing,currRing);
     7235        goto TRAN_LIFTING;
     7236      }
     7237
     7238
     7239      npertstep = nwalk;
     7240      nwalkpert = 1;
     7241      nsteppert ++;
     7242
     7243      /*
     7244      Print("\n// Subroutine needs (%d) steps.", nwalk);
     7245      idElements(Glp, "last G in walk:");
     7246      PrintS("\n// ****************************************");
     7247      Print("\n// Perturb the original target vector (%d): ", nsteppert);
     7248      ivString(target_weight, "new target");
     7249      PrintS("\n// ****************************************\n");
     7250      */
     7251      rChangeCurrRing(newRing);
     7252      G = idrMoveR(Glp, lpRing,currRing);
     7253
     7254      delete next_weight;
     7255
     7256      //Print("\n// ring rNEW = %s;", rString(currRing));
     7257      goto COMPUTE_NEW_VECTOR;
     7258    }
     7259
     7260  TRAN_LIFTING:
     7261    for(i=nV-1; i>=0; i--)
     7262    {
     7263      (*curr_weight)[i] = (*next_weight)[i];
     7264    }
     7265    delete next_weight;
     7266  } // end of while
     7267#ifdef TEST_OVERFLOW
     7268 BE_FINISH:
     7269#endif
     7270  rChangeCurrRing(XXRing);
     7271  G = idrMoveR(G, lpRing,currRing);
     7272
     7273 FINISH:
     7274  delete ivNull;
     7275  delete next_weight;
     7276  delete iv_lp;
     7277  omFree(npert);
     7278
     7279#ifdef TIME_TEST
     7280  Print("\n// Computation took %d steps and %.2f sec", nwalk, ((double) (clock()-mtim)/1000000));
     7281
     7282  TimeStringFractal(tinput, tostd, tif, tstd, textra, tlift, tred, tnw);
     7283
     7284  Print("\n// pSetm_Error = (%d)", ErrorCheck());
     7285  Print("\n// Overflow_Error? (%d)\n", Overflow_Error);
     7286#endif
     7287
     7288  return(G);
     7289}
     7290
     7291
     7292/*****************************************************************
     7293 * compute the reduced Groebner basis of an ideal <Go> w.r.t. lp *
     7294 *****************************************************************/
    51727295static ideal Mpwalk_MAltwalk1(ideal Go, intvec* curr_weight, int tp_deg)
    51737296{
    51747297  Overflow_Error = FALSE;
    5175   BOOLEAN nOverflow_Error = FALSE;
     7298 // BOOLEAN nOverflow_Error = FALSE;
    51767299  clock_t tproc=0;
    51777300  clock_t tinput=clock();
     
    51807303  int tp_deg_tmp = tp_deg;
    51817304  ideal Gomega, M, F, G, M1, F1, Gomega1, Gomega2, G1;
    5182   ring endRing, newRing, oldRing, TargetRing;
     7305  ring newRing, oldRing, TargetRing;
    51837306  intvec* next_weight;
    51847307  intvec* ivNull = new intvec(nV);
    5185   intvec* extra_curr_weight = new intvec(nV);
    51867308
    51877309  ring YXXRing = currRing;
     
    51917313  ideal ssG;
    51927314
    5193   /* perturb the target vector */
     7315  // perturb the target vector
    51947316  while(1)
    51957317  {
     
    51977319    {
    51987320      if (rParameter(currRing) != NULL)
     7321      {
    51997322        DefRingParlp();
     7323      }
    52007324      else
     7325      {
    52017326        VMrDefaultlp();
    5202 
     7327      }
    52037328      TargetRing = currRing;
    52047329      ssG = idrMoveR(Go,YXXRing,currRing);
     
    52067331    Overflow_Error = FALSE;
    52077332    if(tp_deg != 1)
     7333    {
    52087334      target_weight = MPertVectors(ssG, iv_M_dpp, tp_deg);
     7335    }
    52097336    else
    52107337    {
     
    52137340    }
    52147341    if(Overflow_Error == FALSE)
     7342    {
    52157343      break;
    5216 
     7344    }
    52177345    Overflow_Error = TRUE;
    52187346    tp_deg --;
     
    52217349  {
    52227350    Overflow_Error = TRUE;
    5223     nOverflow_Error = TRUE;
     7351    //nOverflow_Error = TRUE;
    52247352  }
    52257353
     
    52287356
    52297357  delete iv_M_dpp;
     7358#ifndef  BUCHBERGER_ALG
    52307359  intvec* hilb_func;
    5231 
    5232   /* to avoid (1,0,...,0) as the target vector */
     7360#endif
     7361  // to avoid (1,0,...,0) as the target vector
    52337362  intvec* last_omega = new intvec(nV);
    52347363  for(i=nV-1; i>0; i--)
     7364  {
    52357365    (*last_omega)[i] = 1;
     7366  }
    52367367  (*last_omega)[0] = 10000;
    52377368
     
    52457376
    52467377    if(nwalk==1)
     7378    {
    52477379      goto FIRST_STEP;
    5248 
     7380    }
    52497381    to=clock();
    5250     /* compute an initial form ideal of <G> w.r.t. "curr_vector" */
     7382    // compute an initial form ideal of <G> w.r.t. "curr_vector"
    52517383    Gomega = MwalkInitialForm(G, curr_weight);
    52527384    xtif=xtif+clock()-to;
    5253 #if 0
    5254     if(Overflow_Error == TRUE)
    5255     {
    5256       for(i=nV-1; i>=0; i--)
    5257         (*curr_weight)[i] = (*extra_curr_weight)[i];
    5258       delete extra_curr_weight;
    5259       goto LASTGB_ALT1;
    5260     }
    5261 #endif
     7385
    52627386#ifndef  BUCHBERGER_ALG
    52637387    if(isNolVector(curr_weight) == 0)
     
    52657389    else
    52667390      hilb_func = hFirstSeries(Gomega,NULL,NULL,last_omega,currRing);
    5267 #endif // BUCHBERGER_ALG
     7391#endif
    52687392
    52697393    oldRing = currRing;
    52707394
    5271     /* define a new ring that its ordering is "(a(curr_weight),lp) */
     7395    // define a new ring that its ordering is "(a(curr_weight),lp)
    52727396    if (rParameter(currRing) != NULL)
     7397    {
    52737398      DefRingPar(curr_weight);
     7399    }
    52747400    else
     7401    {
    52757402      VMrDefault(curr_weight);
    5276 
     7403    }
    52777404    newRing = currRing;
    52787405    Gomega1 = idrMoveR(Gomega, oldRing,currRing);
     
    52887415
    52897416    to=clock();
    5290     /* compute a reduced Groebner basis of <Gomega> w.r.t. "newRing" */
     7417    // compute a reduced Groebner basis of <Gomega> w.r.t. "newRing"
    52917418#ifdef  BUCHBERGER_ALG
    52927419    M = MstdhomCC(Gomega1);
     
    52977424    xtstd=xtstd+clock()-to;
    52987425
    5299     /* change the ring to oldRing */
     7426    // change the ring to oldRing
    53007427    rChangeCurrRing(oldRing);
    53017428    M1 =  idrMoveR(M, newRing,currRing);
     
    53037430    to=clock();
    53047431
    5305     // if(endwalks == 1) PrintS("\n//  Lifting is still working:");
    5306 
    5307     /* compute a reduced Groebner basis of <G> w.r.t. "newRing" by the
    5308      lifting process */
     7432    // if(endwalks == 1){PrintS("\n//  Lifting is still working:");}
     7433
     7434    // compute a reduced Groebner basis of <G> w.r.t. "newRing" by the lifting process
    53097435    F = MLifttwoIdeal(Gomega2, M1, G);
    53107436    xtlift=xtlift+clock()-to;
     
    53147440    idDelete(&G);
    53157441
    5316     /* change the ring to newRing */
     7442    // change the ring to newRing
    53177443    rChangeCurrRing(newRing);
    53187444    F1 = idrMoveR(F, oldRing,currRing);
    53197445    to=clock();
    5320     //if(endwalks == 1) PrintS("\n//  InterRed is still working:");
    5321     /* reduce the Groebner basis <G> w.r.t. the new ring */
     7446    //if(endwalks == 1){ PrintS("\n//  InterRed is still working:");}
     7447    // reduce the Groebner basis <G> w.r.t. the new ring
    53227448    G = kInterRedCC(F1, NULL);
    53237449    xtred=xtred+clock()-to;
     
    53307456    Overflow_Error=FALSE;
    53317457    to=clock();
    5332     /* compute a next weight vector */
     7458    // compute a next weight vector
    53337459    next_weight = MkInterRedNextWeight(curr_weight,target_weight, G);
    53347460    xtnw=xtnw+clock()-to;
     
    53407466    {
    53417467      delete next_weight;
    5342 
    5343       LASTGB_ALT1:
    53447468      if(tp_deg > 1){
    5345         nOverflow_Error = Overflow_Error;
     7469        //nOverflow_Error = Overflow_Error;
    53467470        tproc = tproc+clock()-tinput;
    5347         /*
    5348           Print("\n// A subroutine takes %d steps and calls \"Mpwalk\" (1,%d):",
    5349           nwalk, tp_deg-1);
    5350         */
     7471        //Print("\n// A subroutine takes %d steps and calls \"Mpwalk\" (1,%d):", nwalk, tp_deg-1);
    53517472        G1 = Mpwalk_MAltwalk1(G, curr_weight, tp_deg-1);
    53527473        goto MPW_Finish;
     
    53667487    }
    53677488    if(MivComp(next_weight, target_weight) == 1)
     7489    {
    53687490      endwalks = 1;
    5369 
     7491    }
    53707492    for(i=nV-1; i>=0; i--)
    53717493    {
     
    53767498  }//while
    53777499
    5378   /* check wheather the pertubed target vector is correct */
     7500  // check whether the pertubed target vector is correct
    53797501
    53807502  //define and execute ring with lex. order
    53817503  if (rParameter(currRing) != NULL)
     7504  {
    53827505    DefRingParlp();
     7506  }
    53837507  else
     7508  {
    53847509    VMrDefaultlp();
    5385 
     7510  }
    53867511  G1 = idrMoveR(G, newRing,currRing);
    53877512
     
    53897514  {
    53907515    PrintS("\n// The perturbed target vector doesn't STAY in the correct cone!!");
    5391     if(tp_deg == 1){
     7516    if(tp_deg == 1)
     7517    {
    53927518      //Print("\n// subroutine takes %d steps and applys \"std\"", nwalk);
    53937519      to=clock();
     
    53987524      G2 = NULL;
    53997525    }
    5400     else {
    5401       nOverflow_Error = Overflow_Error;
     7526    else
     7527    {
     7528      //nOverflow_Error = Overflow_Error;
    54027529      tproc = tproc+clock()-tinput;
    5403       /*
    5404         Print("\n// B subroutine takes %d steps and calls \"Mpwalk\" (1,%d) :",
    5405             nwalk,  tp_deg-1);
    5406       */
     7530      // Print("\n// B subroutine takes %d steps and calls \"Mpwalk\" (1,%d) :", nwalk,  tp_deg-1);
    54077531      G1 = Mpwalk_MAltwalk1(G1, curr_weight, tp_deg-1);
    54087532    }
     
    54177541  delete target_weight;
    54187542
    5419   /*
    5420   Print("\n// \"Mpwalk\" (1,%d) took %d steps and %.2f sec. Overflow_Error (%d)", tp_deg,
    5421         nwalk, ((double) clock()-tinput)/1000000, nOverflow_Error);
    5422   */
     7543  //Print("\n// \"Mpwalk\" (1,%d) took %d steps and %.2f sec. Overflow_Error (%d)", tp_deg, nwalk, ((double) clock()-tinput)/1000000, nOverflow_Error);
    54237544
    54247545  return(result);
    54257546}
    54267547
    5427 /* August 2003 */
     7548/*******************************************************************
     7549 * Implementation of the first alternative Groebner Walk Algorithm *
     7550 *******************************************************************/
    54287551ideal MAltwalk1(ideal Go, int op_deg, int tp_deg, intvec* curr_weight,
    54297552                intvec* target_weight)
     
    54317554  Set_Error(FALSE  );
    54327555  Overflow_Error = FALSE;
     7556#ifdef TIME_TEST
    54337557  BOOLEAN nOverflow_Error = FALSE;
     7558#endif
    54347559  // Print("// pSetm_Error = (%d)", ErrorCheck());
    54357560
     
    54427567  int nwalk=0, endwalks=0;
    54437568  int op_tmp = op_deg;
    5444   ideal Gomega, M, F, G, G1, Gomega1, Gomega2, M1, F1;
    5445   ring endRing, newRing, oldRing, TargetRing;
     7569  ideal Gomega, M, F, G, Gomega1, Gomega2, M1, F1;
     7570  ring newRing, oldRing;
    54467571  intvec* next_weight;
    54477572  intvec* iv_M_dp;
     
    54497574  intvec* iv_dp = MivUnit(nV);// define (1,1,...,1)
    54507575  intvec* exivlp = Mivlp(nV);
    5451   intvec* extra_curr_weight = new intvec(nV);
     7576  //intvec* extra_curr_weight = new intvec(nV);
     7577#ifndef  BUCHBERGER_ALG
    54527578  intvec* hilb_func;
    5453 
     7579#endif
    54547580  intvec* cw_tmp = curr_weight;
    54557581
    5456   /* to avoid (1,0,...,0) as the target vector */
     7582  // to avoid (1,0,...,0) as the target vector
    54577583  intvec* last_omega = new intvec(nV);
    54587584  for(i=nV-1; i>0; i--)
     7585  {
    54597586    (*last_omega)[i] = 1;
     7587  }
    54607588  (*last_omega)[0] = 10000;
    54617589
     
    54707598    if(Overflow_Error == FALSE)
    54717599    {
    5472       if(MivComp(curr_weight, iv_dp) == 1) //rOrdStr(currRing) = "dp"
    5473         if(op_tmp == op_deg) {
     7600      if(MivComp(curr_weight, iv_dp) == 1)
     7601      {
     7602      //rOrdStr(currRing) = "dp"
     7603        if(op_tmp == op_deg)
     7604        {
    54747605          G = MstdCC(Go);
    54757606          if(op_deg != 1)
     7607          {
    54767608            iv_M_dp = MivMatrixOrderdp(nV);
     7609          }
    54777610        }
     7611      }
    54787612    }
    54797613    else
    54807614    {
    5481       if(op_tmp == op_deg) {
     7615      if(op_tmp == op_deg)
     7616      {
    54827617        //rOrdStr(currRing) = (a(...),lp,C)
    54837618        if (rParameter(currRing) != NULL)
     7619        {
    54847620          DefRingPar(cw_tmp);
     7621        }
    54857622        else
     7623        {
    54867624          VMrDefault(cw_tmp);
    5487 
     7625        }
    54887626        G = idrMoveR(Go, XXRing,currRing);
    54897627        G = MstdCC(G);
     
    54947632    Overflow_Error = FALSE;
    54957633    if(op_deg != 1)
     7634    {
    54967635      curr_weight = MPertVectors(G, iv_M_dp, op_deg);
    5497     else {
     7636    }
     7637    else
     7638    {
    54987639      curr_weight =  cw_tmp;
    54997640      break;
    55007641    }
    55017642    if(Overflow_Error == FALSE)
     7643    {
    55027644      break;
    5503 
     7645    }
    55047646    Overflow_Error = TRUE;
    55057647    op_deg --;
     
    55207662
    55217663    to = clock();
    5522     /* compute an initial form ideal of <G> w.r.t. "curr_vector" */
     7664    // compute an initial form ideal of <G> w.r.t. "curr_vector"
    55237665    Gomega = MwalkInitialForm(G, curr_weight);
    55247666    xtif=xtif+clock()-to;
     
    55367678#ifndef  BUCHBERGER_ALG
    55377679    if(isNolVector(curr_weight) == 0)
     7680    {
    55387681      hilb_func = hFirstSeries(Gomega,NULL,NULL,curr_weight,currRing);
     7682    }
    55397683    else
     7684    {
    55407685      hilb_func = hFirstSeries(Gomega,NULL,NULL,last_omega,currRing);
     7686    }
    55417687#endif // BUCHBERGER_ALG
    55427688
    55437689    oldRing = currRing;
    55447690
    5545     /* define a new ring which ordering is "(a(curr_weight),lp) */
     7691    // define a new ring which ordering is "(a(curr_weight),lp)
    55467692    if (rParameter(currRing) != NULL)
     7693    {
    55477694      DefRingPar(curr_weight);
     7695    }
    55487696    else
     7697    {
    55497698      VMrDefault(curr_weight);
    5550 
     7699    }
    55517700    newRing = currRing;
    55527701    Gomega1 = idrMoveR(Gomega, oldRing,currRing);
    55537702
    55547703    to=clock();
    5555     /* compute a reduced Groebner basis of <Gomega> w.r.t. "newRing" */
     7704    // compute a reduced Groebner basis of <Gomega> w.r.t. "newRing"
    55567705#ifdef  BUCHBERGER_ALG
    55577706    M = MstdhomCC(Gomega1);
     
    55627711    xtstd=xtstd+clock()-to;
    55637712
    5564     /* change the ring to oldRing */
     7713    // change the ring to oldRing
    55657714    rChangeCurrRing(oldRing);
    55667715    M1 =  idrMoveR(M, newRing,currRing);
     
    55687717
    55697718    to=clock();
    5570     /* compute a reduced Groebner basis of <G> w.r.t. "newRing" by the
    5571      lifting process */
     7719    // compute a reduced Groebner basis of <G> w.r.t. "newRing" by the lifting process
    55727720    F = MLifttwoIdeal(Gomega2, M1, G);
    55737721    xtlift=xtlift+clock()-to;
     
    55777725    idDelete(&G);
    55787726
    5579    /* change the ring to newRing */
     7727    // change the ring to newRing
    55807728    rChangeCurrRing(newRing);
    55817729    F1 = idrMoveR(F, oldRing,currRing);
    55827730
    55837731    to=clock();
    5584     /* reduce the Groebner basis <G> w.r.t. new ring */
     7732    // reduce the Groebner basis <G> w.r.t. new ring
    55857733    G = kInterRedCC(F1, NULL);
    55867734    xtred=xtred+clock()-to;
     
    55887736
    55897737    if(endwalks == 1)
     7738    {
    55907739      break;
    5591 
     7740    }
    55927741  NEXT_VECTOR:
    55937742    to=clock();
    5594     /* compute a next weight vector */
     7743    // compute a next weight vector
    55957744    next_weight = MkInterRedNextWeight(curr_weight,target_weight, G);
    55967745    xtnw=xtnw+clock()-to;
     
    56047753
    56057754      if (rParameter(currRing) != NULL)
     7755      {
    56067756        DefRingPar(target_weight);
     7757      }
    56077758      else
     7759      {
    56087760        VMrDefault(target_weight);
    5609 
     7761      }
    56107762      F1 = idrMoveR(G, newRing,currRing);
    56117763      G = MstdCC(F1);
     
    56307782      else
    56317783      {
    5632         MSTD_ALT1:
     7784       // MSTD_ALT1:
     7785#ifdef TIME_TEST
    56337786        nOverflow_Error = Overflow_Error;
     7787#endif
    56347788        tproc = clock()-xftinput;
    5635         /*
    5636         Print("\n//  main routine takes %d steps and calls \"Mpwalk\" (1,%d):",
    5637               nwalk,  tp_deg);
    5638         */
    5639         // compute the red. GB of <G> w.r.t. the lex order by
    5640         // the "recursive-modified" perturbation walk alg (1,tp_deg)
     7789       
     7790        Print("\n//  main routine takes %d steps and calls \"Mpwalk\" (1,%d):", nwalk,  tp_deg);
     7791       
     7792        // compute the red. GB of <G> w.r.t. the lex order by the "recursive-modified" perturbation walk alg (1,tp_deg)
    56417793        G = Mpwalk_MAltwalk1(G, curr_weight, tp_deg);
    56427794        delete next_weight;
     
    56457797    }
    56467798
    5647     /* 06.11.01 NOT Changed, to free memory*/
     7799    //NOT Changed, to free memory
    56487800    for(i=nV-1; i>=0; i--)
    56497801    {
     
    56607812  delete ivNull;
    56617813  if(op_deg != 1 )
     7814  {
    56627815    delete curr_weight;
    5663 
     7816  }
    56647817  delete exivlp;
    56657818#ifdef TIME_TEST
     
    56767829  return(result);
    56777830}
    5678 
  • Singular/walk.h

    r8659a9 re92c6a  
    5555ideal Mwalk(ideal G, intvec* curr_weight, intvec* target_weight);
    5656
     57// random walk algorithm to compute a Groebner basis
     58ideal Mrwalk(ideal Go, intvec* curr_weight, intvec* target_weight, int weight_rad, int pert_deg);
     59
    5760/* the perturbation walk algorithm */
    5861ideal Mpwalk(ideal G,int op,int tp,intvec* curr_weight,intvec* target_weight, int nP);
     
    6164ideal Mfwalk(ideal G, intvec* ivstart, intvec* ivtarget);
    6265
     66/* The fractal walk algorithm with random element */
     67ideal Mfrwalk(ideal G, intvec* ivstart, intvec* ivtarget,int weight_rad);
    6368
    6469/* Implement Tran's idea */
    6570intvec* TranMPertVectorslp(ideal G);
    6671ideal TranMImprovwalk(ideal Go, intvec* curr_weight,intvec* target_weight, int nP);
     72
     73/* Implement Tran's idea with random element*/
     74ideal TranMrImprovwalk(ideal G,intvec* curr_weight,intvec* target_tmp, int nP, int weight_rad, int pert_deg);
    6775
    6876/* the first alternative algorithm */
Note: See TracChangeset for help on using the changeset viewer.