Changeset 8d7ca0 in git for Singular/walk.cc


Ignore:
Timestamp:
Jul 2, 2015, 7:33:23 PM (9 years ago)
Author:
Stephan Oberfranz <oberfran@…>
Branches:
(u'spielwiese', '2a584933abf2a2d3082034c7586d38bb6de1a30a')
Children:
654a230847a8318a2a297f1520550f19f34a4401
Parents:
bbfc4a825e46a6bf0bdd30896d5189f16c4e4d99
Message:
computation of intermediate weight adapted
File:
1 edited

Legend:

Unmodified
Added
Removed
  • Singular/walk.cc

    rbbfc4a8 r8d7ca0  
    1616
    1717//#define TEST_OVERFLOW
    18 //#define CHECK_IDEAL
    19 //#define CHECK_IDEAL_MWALK
     18
     19#define CHECK_IDEAL_MWALK //to print intermediate results
    2020
    2121//#define NEXT_VECTORS_CC
    22 //#define PRINT_VECTORS //to print vectors (sigma, tau, omega)
     22#define PRINT_VECTORS //to print weight vectors
    2323
    2424#define INVEPS_SMALL_IN_FRACTAL  //to choose the small invers of epsilon
     
    2727
    2828#define FIRST_STEP_FRACTAL // to define the first step of the fractal
    29 #define MSTDCC_FRACTAL // apply Buchberger alg to compute a red GB, if
    30 //                          tau doesn't stay in the correct cone
     29#define MSTDCC_FRACTAL // apply Buchberger alg to compute a red GB, if tau doesn't stay in the correct cone
    3130
    3231//#define TIME_TEST // print the used time of each subroutine
     
    122121 ************************************/
    123122// unused
    124 #if 0
     123/*
    125124static void initSSpecialCC (ideal F, ideal Q, ideal P,kStrategy strat)
    126125{
     
    270269#endif
    271270}
    272 #endif
     271*/
    273272
    274273/*****************
     
    279278  int j;
    280279  kStrategy strat = new skStrategy;
    281 
    282 //  if (TEST_OPT_PROT)
    283 //  {
    284 //    writeTime("start InterRed:");
    285 //    mflush();
    286 //  }
    287   //strat->syzComp     = 0;
     280/*
     281  if (TEST_OPT_PROT)
     282  {
     283    writeTime("start InterRed:");
     284    mflush();
     285  }
     286  strat->syzComp     = 0;
     287*/
    288288  strat->kHEdgeFound = (currRing->ppNoether) != NULL;
    289289  strat->kNoether=pCopy((currRing->ppNoether));
     
    348348    strat->fromQ = NULL;
    349349  }
    350 //  if (TEST_OPT_PROT)
    351 //  {
    352 //    writeTime("end Interred:");
    353 //    mflush();
    354 //  }
     350/*
     351  if (TEST_OPT_PROT)
     352  {
     353    writeTime("end Interred:");
     354    mflush();
     355  }
     356*/
    355357  ideal shdl=strat->Shdl;
    356358  idSkipZeroes(shdl);
     
    360362}
    361363
    362 //unused
    363 #if 0
     364#ifdef TIME_TEST
    364365static void TimeString(clock_t tinput, clock_t tostd, clock_t tif,clock_t tstd,
    365366                       clock_t tlf,clock_t tred, clock_t tnw, int step)
     
    399400        ((((double) xtextra)/1000000)/totm)*100);
    400401}
    401 #endif
    402 
    403 //unused
    404 #if 0
     402
    405403static void TimeStringFractal(clock_t tinput, clock_t tostd, clock_t tif,clock_t tstd,
    406404                       clock_t textra, clock_t tlf,clock_t tred, clock_t tnw)
     
    431429#endif
    432430
    433 //#ifdef CHECK_IDEAL_MWALK
     431#ifdef CHECK_IDEAL_MWALK
    434432static void idString(ideal L, const char* st)
    435433{
     
    443441  Print(" %s;", pString(L->m[nL-1]));
    444442}
    445 //#endif
    446 
     443#endif
     444/*
    447445#if defined(CHECK_IDEAL_MWALK) || defined(ENDWALKS)
    448446static void headidString(ideal L, char* st)
     
    498496}
    499497#endif
    500 
     498*/
    501499
    502500static void ivString(intvec* iv, const char* ch)
     
    512510}
    513511
    514 //unused
    515 //#if 0
     512#ifdef PRINT_VECTORS
    516513static void MivString(intvec* iva, intvec* ivb, intvec* ivc)
    517514{
     
    536533  Print("%d)", (*ivc)[nV]);
    537534}
    538 //#endif
     535#endif
    539536
    540537/********************************************************************
     
    10401037 * print the max total degree and the max coefficient of G                   *
    10411038 *****************************************************************************/
    1042 #if 0
     1039/*
    10431040static void checkComplexity(ideal G, char* cG)
    10441041{
     
    10811078  PrintLn();
    10821079}
    1083 #endif
     1080*/
    10841081
    10851082/*****************************************************************************
     
    11031100  intvec* v_null =  new intvec(nV);
    11041101
    1105 
    11061102  // Check that the perturbed degree is valid
    11071103  if(pdeg > nV || pdeg <= 0)
     
    11171113  }
    11181114  mpz_t *pert_vector = (mpz_t*)omAlloc(nV*sizeof(mpz_t));
    1119   //mpz_t *pert_vector1 = (mpz_t*)omAlloc(nV*sizeof(mpz_t));
     1115  mpz_t *pert_vector1 = (mpz_t*)omAlloc(nV*sizeof(mpz_t));
    11201116
    11211117  for(i=0; i<nV; i++)
    11221118  {
    11231119    mpz_init_set_si(pert_vector[i], (*ivtarget)[i]);
    1124    // mpz_init_set_si(pert_vector1[i], (*ivtarget)[i]);
     1120    mpz_init_set_si(pert_vector1[i], (*ivtarget)[i]);
    11251121  }
    11261122  // Calculate max1 = Max(A2)+Max(A3)+...+Max(Apdeg),
     
    12021198    }
    12031199  }
     1200
     1201  // 2147483647 is max. integer representation in SINGULAR
     1202  mpz_t sing_int;
     1203  mpz_init_set_ui(sing_int,  2147483647);
     1204
     1205  mpz_t check_int;
     1206  mpz_init_set_ui(check_int,  100000);
     1207
    12041208  mpz_t ztemp;
    12051209  mpz_init(ztemp);
     
    12211225  }
    12221226
    1223   intvec *pert_vector1= new intvec(nV);
    1224   j = 0;
    12251227  for(i=0; i<nV; i++)
    12261228  {
    1227     (* pert_vector1)[i] = mpz_get_si(pert_vector[i]);
    1228     (* pert_vector1)[i] = 0.1*(* pert_vector1)[i];
    1229     (* pert_vector1)[i] = floor((* pert_vector1)[i] + 0.5);
    1230     if((* pert_vector1)[i] == 0)
    1231     {
    1232       j++;
    1233     }
    1234   }
    1235   if(j > nV - 1)
    1236   {
    1237     // Print("\n//  MPertVectors: geaenderter vector gleich Null! \n");
    1238     delete pert_vector1;
    1239     goto CHECK_OVERFLOW;
    1240   }
    1241 
    1242 // check that the perturbed weight vector lies in the Groebner cone
    1243   if(test_w_in_ConeCC(G,pert_vector1) != 0)
    1244   {
    1245     // Print("\n//  MPertVectors: geaenderter vector liegt in Groebnerkegel! \n");
     1229    if(mpz_cmp(pert_vector[i], check_int)>=0)
     1230    {
     1231      for(j=0; j<nV; j++)
     1232      {
     1233        mpz_fdiv_q_ui(pert_vector1[j], pert_vector[j], 100);
     1234      }
     1235    }
     1236  }
     1237
     1238  intvec* result = new intvec(nV);
     1239
     1240  int ntrue=0;
     1241
     1242  for(i=0; i<nV; i++)
     1243  {
     1244    (*result)[i] = mpz_get_si(pert_vector1[i]);
     1245    if(mpz_cmp(pert_vector1[i], sing_int)>=0)
     1246    {
     1247      ntrue++;
     1248    }
     1249  }
     1250  if(ntrue > 0 || test_w_in_ConeCC(G,result)==0)
     1251  {
     1252    ntrue=0;
    12461253    for(i=0; i<nV; i++)
    12471254    {
    1248       mpz_set_si(pert_vector[i], (*pert_vector1)[i]);
    1249     }
    1250   }
    1251   else
    1252   {
    1253     //Print("\n// MpertVectors: geaenderter vector liegt nicht in Groebnerkegel! \n");
    1254   }
    1255   delete pert_vector1;
    1256 
    1257   CHECK_OVERFLOW:
    1258   intvec* result = new intvec(nV);
    1259 
    1260   /* 2147483647 is max. integer representation in SINGULAR */
    1261   mpz_t sing_int;
    1262   mpz_init_set_ui(sing_int,  2147483647);
    1263 
    1264   int ntrue=0;
    1265   for(i=0; i<nV; i++)
    1266   {
    1267     (*result)[i] = mpz_get_si(pert_vector[i]);
    1268     if(mpz_cmp(pert_vector[i], sing_int)>=0)
    1269     {
    1270       ntrue++;
    1271       if(Overflow_Error == FALSE)
    1272       {
    1273         Overflow_Error = TRUE;
    1274         PrintS("\n// ** OVERFLOW in \"MPertvectors\": ");
    1275         mpz_out_str( stdout, 10, pert_vector[i]);
    1276         PrintS(" is greater than 2147483647 (max. integer representation)");
    1277         Print("\n//  So vector[%d] := %d is wrong!!", i+1, (*result)[i]);
    1278       }
    1279     }
    1280   }
    1281 
    1282   if(Overflow_Error == TRUE)
    1283   {
    1284     ivString(result, "pert_vector");
    1285     Print("\n// %d element(s) of it is overflow!!", ntrue);
     1255      (*result)[i] = mpz_get_si(pert_vector[i]);
     1256      if(mpz_cmp(pert_vector[i], sing_int)>=0)
     1257      {
     1258        ntrue++;
     1259        if(Overflow_Error == FALSE)
     1260        {
     1261          Overflow_Error = TRUE;
     1262          PrintS("\n// ** OVERFLOW in \"MPertvectors\": ");
     1263          mpz_out_str( stdout, 10, pert_vector[i]);
     1264          PrintS(" is greater than 2147483647 (max. integer representation)");
     1265          Print("\n//  So vector[%d] := %d is wrong!!", i+1, (*result)[i]);
     1266        }
     1267      }
     1268    }
     1269
     1270    if(Overflow_Error == TRUE)
     1271    {
     1272      ivString(result, "pert_vector");
     1273      Print("\n// %d element(s) of it is overflow!!", ntrue);
     1274    }
    12861275  }
    12871276
    12881277  mpz_clear(ztemp);
    12891278  mpz_clear(sing_int);
     1279  mpz_clear(check_int);
    12901280  omFree(pert_vector);
    1291   //omFree(pert_vector1);
     1281  omFree(pert_vector1);
    12921282  mpz_clear(tot_deg);
    12931283  mpz_clear(maxdeg);
     
    14911481
    14921482//unused
    1493 #if 0
     1483/*
    14941484static intvec* MatrixOrderdp(int nV)
    14951485{
     
    15071497  return(ivM);
    15081498}
    1509 #endif
     1499*/
    15101500
    15111501intvec* MivUnit(int nV)
     
    15841574    mpz_cdiv_q_ui(inveps, inveps, nV);
    15851575  }
    1586   //PrintS("\n// choose the \"small\" inverse epsilon!");
     1576  // choose the small inverse epsilon
    15871577#endif
    15881578
     
    16181608
    16191609    for(j=0; j<nV; j++)
    1620       {
     1610    {
    16211611      mpz_init_set(pert_vector[i*nV+j],ivtemp[j]);
    1622       }
    1623   }
    1624 
    1625   /* 2147483647 is max. integer representation in SINGULAR */
     1612    }
     1613  }
     1614
     1615  // 2147483647 is max. integer representation in SINGULAR
    16261616  mpz_t sing_int;
    16271617  mpz_init_set_ui(sing_int,  2147483647);
     
    16471637    (* result)[i] = mpz_get_si(pert_vector[i]);
    16481638  }
    1649 /*
    1650   j = 0;
    1651   for(i=0; i<niv; i++)
    1652   {
    1653     (* result1)[i] = mpz_get_si(pert_vector[i]);
    1654     (* result1)[i] = 0.1*(* result1)[i];
    1655     (* result1)[i] = floor((* result1)[i] + 0.5);
    1656     if((* result1)[i] == 0)
    1657     {
    1658       j++;
    1659     }
    1660   }
    1661   if(j > niv - 1)
    1662   {
    1663     // Print("\n//  MfPertwalk: geaenderter vector gleich Null! \n");
    1664     delete result1;
    1665     goto CHECK_OVERFLOW;
    1666   }
    1667 /*
    1668 // check that the perturbed weight vector lies in the Groebner cone
    1669   Print("\n========================================\n//** MfPertvector: test in Cone.\n");
    1670   if(test_w_in_ConeCC(G,result1) != 0)
    1671   {
    1672     // Print("\n//  MfPertwalk: geaenderter vector liegt in Groebnerkegel! \n");
    1673     delete result;
    1674     result = result1;
    1675     for(i=0; i<nV; i++)
    1676     {
    1677       mpz_set_si(pert_vector[i], (*result1)[i]);
    1678     }
    1679   }
    1680   else
    1681   {
    1682     delete result1;
    1683     // Print("\n// Mfpertwalk: geaenderter vector liegt nicht in Groebnerkegel! \n");
    1684   }
    1685   Print("\n ========================================\n");*/
     1639
    16861640  CHECK_OVERFLOW:
    16871641
     
    17211675    while(p!=NULL)
    17221676    {
    1723       p_Setm(p,currRing); pIter(p);
     1677      p_Setm(p,currRing);
     1678      pIter(p);
    17241679    }
    17251680  }
     
    18041759
    18051760//unused
    1806 #if 0
     1761/*
    18071762static void checkidealCC(ideal G, char* Ch)
    18081763{
     
    18301785  PrintLn();
    18311786}
    1832 #endif
     1787*/
    18331788
    18341789//unused
    1835 #if 0
     1790/*
    18361791static void HeadidString(ideal L, char* st)
    18371792{
     
    18451800  Print(" %s;\n", pString(pHead(L->m[nL])));
    18461801}
    1847 #endif
    1848 
     1802
     1803*/
    18491804static inline int MivComp(intvec* iva, intvec* ivb)
    18501805{
     
    19191874 * with respect to an ideal <G>.                                      *
    19201875**********************************************************************/
     1876/*
    19211877static intvec* MwalkNextWeightCC(intvec* curr_weight, intvec* target_weight,
    19221878                                 ideal G)
     
    20401996  if(mpz_cmp(t_nenner, t_null) == 0)
    20411997  {
    2042     #ifndef SING_NDEBUG
    2043     Print("\n//MwalkNextWeightCC: t_nenner ist Null!");
    2044     #endif
     1998#ifndef SING_NDEBUG
     1999    Print("\n//MwalkNextWeightCC: t_nenner=0\n");
     2000#endif
    20452001    delete diff_weight;
    20462002    diff_weight = ivCopy(curr_weight);//take memory
     
    21712127
    21722128  for(j=0; j<nRing; j++)
    2173     {
     2129  {
    21742130    if(mpz_cmp(vec[j], sing_int_half) >= 0)
    2175       {
     2131    {
    21762132      goto REDUCTION;
    2177       }
    2178     }
     2133    }
     2134  }
    21792135  checkRed = 1;
    21802136  for (j=0; j<nRing; j++)
     
    22712227return diff_weight;
    22722228}
     2229*/
     2230/**********************************************************************
     2231 * Compute a next weight vector between curr_weight and target_weight *
     2232 * with respect to an ideal <G>.                                      *
     2233**********************************************************************/
     2234static intvec* MwalkNextWeightCC(intvec* curr_weight, intvec* target_weight,
     2235                                 ideal G)
     2236{
     2237  BOOLEAN nError = Overflow_Error;
     2238  Overflow_Error = FALSE;
     2239
     2240  assume(currRing != NULL && curr_weight != NULL &&
     2241         target_weight != NULL && G != NULL);
     2242
     2243  int nRing = currRing->N;
     2244  int checkRed, j, nG = IDELEMS(G);
     2245  intvec* ivtemp;
     2246
     2247  mpz_t t_zaehler, t_nenner;
     2248  mpz_init(t_zaehler);
     2249  mpz_init(t_nenner);
     2250
     2251  mpz_t s_zaehler, s_nenner, temp, MwWd;
     2252  mpz_init(s_zaehler);
     2253  mpz_init(s_nenner);
     2254  mpz_init(temp);
     2255  mpz_init(MwWd);
     2256
     2257  mpz_t sing_int;
     2258  mpz_init(sing_int);
     2259  mpz_set_si(sing_int,  2147483647);
     2260
     2261  mpz_t sing_int_half;
     2262  mpz_init(sing_int_half);
     2263  mpz_set_si(sing_int_half,  3*(1073741824/2));
     2264
     2265  mpz_t deg_w0_p1, deg_d0_p1;
     2266  mpz_init(deg_w0_p1);
     2267  mpz_init(deg_d0_p1);
     2268
     2269  mpz_t sztn, sntz;
     2270  mpz_init(sztn);
     2271  mpz_init(sntz);
     2272
     2273  mpz_t t_null;
     2274  mpz_init(t_null);
     2275
     2276  mpz_t ggt;
     2277  mpz_init(ggt);
     2278
     2279  mpz_t dcw;
     2280  mpz_init(dcw);
     2281
     2282  int gcd_tmp;
     2283  //intvec* diff_weight = MivSub(target_weight, curr_weight);
     2284
     2285  intvec* diff_weight1 = new intvec(nRing); //MivSub(target_weight, curr_weight);
     2286  poly g;
     2287
     2288  // reduce the size of the entries of the current weight vector
     2289  if(TEST_OPT_REDSB)
     2290  {
     2291    for (j=0; j<nRing; j++)
     2292    {
     2293      (*diff_weight1)[j] = (*curr_weight)[j];
     2294    }
     2295    while(MivAbsMax(diff_weight1)>10000 && test_w_in_ConeCC(G,diff_weight1)==1)
     2296    {
     2297      for(j=0; j<nRing; j++)
     2298      {
     2299        (*curr_weight)[j] = (*diff_weight1)[j]; 
     2300      }
     2301      for(j=0; j<nRing; j++)
     2302      {
     2303        (*diff_weight1)[j] = floor(0.1*(*diff_weight1)[j] + 0.5);
     2304      }
     2305    }
     2306
     2307    if(MivAbsMax(curr_weight)>100000)
     2308    {
     2309      for(j=0; j<nRing; j++)
     2310      {
     2311        (*diff_weight1)[j] = (*curr_weight)[j];
     2312      }
     2313      j = 0;
     2314      while(test_w_in_ConeCC(G,diff_weight1)==1 && MivAbsMax(diff_weight1)>1000)
     2315      {
     2316        (*curr_weight)[j] = (*diff_weight1)[j];
     2317        j = MivAbsMaxArg(diff_weight1);
     2318        (*diff_weight1)[j] = floor(0.1*(*diff_weight1)[j] + 0.5);
     2319      }
     2320    }
     2321
     2322  }
     2323  intvec* diff_weight = MivSub(target_weight, curr_weight);
     2324
     2325  // compute a suitable next weight vector
     2326  for (j=0; j<nG; j++)
     2327  {
     2328    g = G->m[j];
     2329    if (g != NULL)
     2330    {
     2331      ivtemp = MExpPol(g);
     2332      mpz_set_si(deg_w0_p1, MivDotProduct(ivtemp, curr_weight));
     2333      mpz_set_si(deg_d0_p1, MivDotProduct(ivtemp, diff_weight));
     2334      delete ivtemp;
     2335
     2336      pIter(g);
     2337      while (g != NULL)
     2338      {
     2339        ivtemp = MExpPol(g);
     2340        mpz_set_si(MwWd, MivDotProduct(ivtemp, curr_weight));
     2341        mpz_sub(s_zaehler, deg_w0_p1, MwWd);
     2342        if(mpz_cmp(s_zaehler, t_null) != 0)
     2343        {
     2344          mpz_set_si(MwWd, MivDotProduct(ivtemp, diff_weight));
     2345          mpz_sub(s_nenner, MwWd, deg_d0_p1);
     2346          // check for 0 < s <= 1
     2347          if( (mpz_cmp(s_zaehler,t_null) > 0 &&
     2348               mpz_cmp(s_nenner, s_zaehler)>=0) ||
     2349              (mpz_cmp(s_zaehler, t_null) < 0 &&
     2350               mpz_cmp(s_nenner, s_zaehler)<=0))
     2351          {
     2352            // make both positive
     2353            if (mpz_cmp(s_zaehler, t_null) < 0)
     2354            {
     2355              mpz_neg(s_zaehler, s_zaehler);
     2356              mpz_neg(s_nenner, s_nenner);
     2357            }
     2358
     2359            //compute a simple fraction of s
     2360            cancel(s_zaehler, s_nenner);
     2361
     2362            if(mpz_cmp(t_nenner, t_null) != 0)
     2363            {
     2364              mpz_mul(sztn, s_zaehler, t_nenner);
     2365              mpz_mul(sntz, s_nenner, t_zaehler);
     2366
     2367              if(mpz_cmp(sztn,sntz) < 0)
     2368              {
     2369                mpz_add(t_nenner, t_null, s_nenner);
     2370                mpz_add(t_zaehler,t_null, s_zaehler);
     2371              }
     2372            }
     2373            else
     2374            {
     2375              mpz_add(t_nenner, t_null, s_nenner);
     2376              mpz_add(t_zaehler,t_null, s_zaehler);
     2377            }
     2378          }
     2379        }
     2380        pIter(g);
     2381        delete ivtemp;
     2382      }
     2383    }
     2384  }
     2385  //Print("\n// Alloc Size = %d \n", nRing*sizeof(mpz_t));
     2386  mpz_t *vec=(mpz_t*)omAlloc(nRing*sizeof(mpz_t));
     2387
     2388
     2389  // there is no 0<t<1 and define the next weight vector that is equal
     2390  // to the current weight vector
     2391  if(mpz_cmp(t_nenner, t_null) == 0)
     2392  {
     2393#ifndef SING_NDEBUG
     2394    Print("\n//MwalkNextWeightCC: t_nenner=0\n");
     2395#endif
     2396    delete diff_weight;
     2397    diff_weight = ivCopy(curr_weight);//take memory
     2398    goto FINISH;
     2399  }
     2400
     2401  // define the target vector as the next weight vector, if t = 1
     2402  if(mpz_cmp_si(t_nenner, 1)==0 && mpz_cmp_si(t_zaehler,1)==0)
     2403  {
     2404    delete diff_weight;
     2405    diff_weight = ivCopy(target_weight); //this takes memory
     2406    goto FINISH;
     2407  }
     2408
     2409   //checkRed = 0;
     2410
     2411  SIMPLIFY_GCD:
     2412
     2413  // simplify the vectors curr_weight and diff_weight (C-int)
     2414  gcd_tmp = (*curr_weight)[0];
     2415
     2416  for (j=1; j<nRing; j++)
     2417  {
     2418    gcd_tmp = gcd(gcd_tmp, (*curr_weight)[j]);
     2419    if(gcd_tmp == 1)
     2420    {
     2421      break;
     2422    }
     2423  }
     2424  if(gcd_tmp != 1)
     2425  {
     2426    for (j=0; j<nRing; j++)
     2427    {
     2428      gcd_tmp = gcd(gcd_tmp, (*diff_weight)[j]);
     2429      if(gcd_tmp == 1)
     2430      {
     2431        break;
     2432      }
     2433    }
     2434  }
     2435  if(gcd_tmp != 1)
     2436  {
     2437    for (j=0; j<nRing; j++)
     2438    {
     2439      (*curr_weight)[j] =  (*curr_weight)[j]/gcd_tmp;
     2440      (*diff_weight)[j] =  (*diff_weight)[j]/gcd_tmp;
     2441    }
     2442  }
     2443  if(checkRed > 0)
     2444  {
     2445    for (j=0; j<nRing; j++)
     2446    {
     2447      mpz_set_si(vec[j], (*diff_weight)[j]);
     2448    }
     2449    goto TEST_OVERFLOW;
     2450  }
     2451
     2452#ifdef  NEXT_VECTORS_CC
     2453  Print("\n// gcd of the weight vectors (current and target) = %d", gcd_tmp);
     2454  ivString(curr_weight, "new cw");
     2455  ivString(diff_weight, "new dw");
     2456
     2457  PrintS("\n// t_zaehler: ");  mpz_out_str( stdout, 10, t_zaehler);
     2458  PrintS(", t_nenner: ");  mpz_out_str( stdout, 10, t_nenner);
     2459#endif
     2460
     2461// construct a new weight vector and check whether vec[j] is overflow, i.e. vec[j] > 2^31.
     2462// If vec[j] doesn't overflow, define a weight vector. Otherwise, report that overflow
     2463// appears. In the second case, test whether the the correctness of the new vector plays
     2464// an important role
     2465
     2466  for (j=0; j<nRing; j++)
     2467  {
     2468    mpz_set_si(dcw, (*curr_weight)[j]);
     2469    mpz_mul(s_nenner, t_nenner, dcw);
     2470
     2471    if( (*diff_weight)[j]>0)
     2472    {
     2473      mpz_mul_ui(s_zaehler, t_zaehler, (*diff_weight)[j]);
     2474    }
     2475    else
     2476    {
     2477      mpz_mul_ui(s_zaehler, t_zaehler, -(*diff_weight)[j]);
     2478      mpz_neg(s_zaehler, s_zaehler);
     2479    }
     2480    mpz_add(sntz, s_nenner, s_zaehler);
     2481    mpz_init_set(vec[j], sntz);
     2482
     2483#ifdef NEXT_VECTORS_CC
     2484    Print("\n//   j = %d ==> ", j);
     2485    PrintS("(");
     2486    mpz_out_str( stdout, 10, t_nenner);
     2487    Print(" * %d)", (*curr_weight)[j]);
     2488    Print(" + ("); mpz_out_str( stdout, 10, t_zaehler);
     2489    Print(" * %d) =  ",  (*diff_weight)[j]);
     2490    mpz_out_str( stdout, 10, s_nenner);
     2491    PrintS(" + ");
     2492    mpz_out_str( stdout, 10, s_zaehler);
     2493    PrintS(" = "); mpz_out_str( stdout, 10, sntz);
     2494    Print(" ==> vector[%d]: ", j); mpz_out_str(stdout, 10, vec[j]);
     2495#endif
     2496
     2497    if(j==0)
     2498    {
     2499      mpz_set(ggt, sntz);
     2500    }
     2501    else
     2502    {
     2503      if(mpz_cmp_si(ggt,1) != 0)
     2504      {
     2505        mpz_gcd(ggt, ggt, sntz);
     2506      }
     2507    }
     2508  }
     2509  // reduce the vector with the gcd
     2510  if(mpz_cmp_si(ggt,1) != 0)
     2511  {
     2512    for (j=0; j<nRing; j++)
     2513    {
     2514      mpz_divexact(vec[j], vec[j], ggt);
     2515    }
     2516  }
     2517#ifdef  NEXT_VECTORS_CC
     2518  PrintS("\n// gcd of elements of the vector: ");
     2519  mpz_out_str( stdout, 10, ggt);
     2520#endif
     2521
     2522  for (j=0; j<nRing; j++)
     2523  {
     2524    (*diff_weight)[j] = mpz_get_si(vec[j]);
     2525  }
     2526
     2527 TEST_OVERFLOW:
     2528
     2529  for (j=0; j<nRing; j++)
     2530  {
     2531    if(mpz_cmp(vec[j], sing_int)>=0)
     2532    {
     2533      if(Overflow_Error == FALSE)
     2534      {
     2535        Overflow_Error = TRUE;
     2536        PrintS("\n// ** OVERFLOW in \"MwalkNextWeightCC\": ");
     2537        mpz_out_str( stdout, 10, vec[j]);
     2538        PrintS(" is greater than 2147483647 (max. integer representation)\n");
     2539        Print("//  So vector[%d] := %d is wrong!!\n",j+1, vec[j]);// vec[j] is mpz_t
     2540      }
     2541    }
     2542  }
     2543
     2544 FINISH:
     2545   delete diff_weight1;
     2546   mpz_clear(t_zaehler);
     2547   mpz_clear(t_nenner);
     2548   mpz_clear(s_zaehler);
     2549   mpz_clear(s_nenner);
     2550   mpz_clear(sntz);
     2551   mpz_clear(sztn);
     2552   mpz_clear(temp);
     2553   mpz_clear(MwWd);
     2554   mpz_clear(deg_w0_p1);
     2555   mpz_clear(deg_d0_p1);
     2556   mpz_clear(ggt);
     2557   omFree(vec);
     2558   mpz_clear(sing_int_half);
     2559   mpz_clear(sing_int);
     2560   mpz_clear(dcw);
     2561   mpz_clear(t_null);
     2562
     2563  if(Overflow_Error == FALSE)
     2564  {
     2565    Overflow_Error = nError;
     2566  }
     2567  rComplete(currRing);
     2568  for(j=0; j<IDELEMS(G); j++)
     2569  {
     2570    poly p=G->m[j];
     2571    while(p!=NULL)
     2572    {
     2573      p_Setm(p,currRing);
     2574      pIter(p);
     2575    }
     2576  }
     2577return diff_weight;
     2578}
     2579
    22732580
    22742581/**********************************************************************
     
    29563263  int i,j,N = IDELEMS(Gomega);
    29573264  poly p,lm,factor1,factor2;
    2958   //PrintS("\n//** idCopy\n");
     3265
    29593266  ideal Go = idCopy(G);
    29603267 
    2961   //PrintS("\n//** jetzt for-Loop!\n");
    2962 
    29633268  // check whether leading monomials of G and Gomega coincide
    29643269  // and return NULL if not
    29653270  for(i=0; i<N; i++)
    29663271  {
    2967     p = pCopy(Gomega->m[i]);
    2968     lm = pCopy(pHead(G->m[i]));
    2969     if(!pIsConstant(pSub(p,lm)))
    2970     {
    2971       //pDelete(&p);
    2972       //pDelete(&lm);
     3272    if(!pIsConstant(pSub(pCopy(Gomega->m[i]),pCopy(pHead(G->m[i])))))
     3273    {
    29733274      idDelete(&Go);
    29743275      return NULL;
    29753276    }
    2976     //pDelete(&p);
    2977     //pDelete(&lm);
    29783277  }
    29793278  for(i=0; i<N; i++)
     
    30073306    }
    30083307  }
    3009  
    3010   //PrintS("\n//** jetzt Delete!\n");
    3011   //pDelete(&p);
    3012   //pDelete(&factor);
    3013   //pDelete(&lm);
     3308
    30143309  if(middle == TRUE)
    30153310  {
    3016     //PrintS("\n//** middle TRUE!\n");
    30173311    return Go;
    30183312  }
    3019   //PrintS("\n//** middle FALSE!\n");
    30203313  idDelete(&Go);
    30213314  return NULL;
     
    31513444    {
    31523445      Print("\n// ring r%d_%d = %s;\n", tp_deg, nwalk, rString(currRing));
     3446/*
    31533447      idElements(Gomega, "Gw");
    31543448      headidString(Gomega, "Gw");
     3449*/
    31553450    }
    31563451#endif
     
    33023597
    33033598/**********************************************************
    3304  * check whether a polynomial of G has least 3 monomials  *
     3599 * check whether a polynomial of G has least 4 monomials  *
    33053600 **********************************************************/
    33063601static int lengthpoly(ideal G)
     
    33093604  for(i=IDELEMS(G)-1; i>=0; i--)
    33103605  {
    3311 #if 0
    3312     if(pLength(G->m[i])>2)
    3313     {
    3314       return 1;
    3315     }
    3316 #else
    33173606    if((G->m[i]!=NULL) /* len >=0 */
    33183607       && (G->m[i]->next!=NULL) /* len >=1 */
    33193608       && (G->m[i]->next->next!=NULL) /* len >=2 */
    33203609       && (G->m[i]->next->next->next!=NULL) /* len >=3 */
    3321       //&& (G->m[i]->next->next->next->next!=NULL) /* len >=4 */
    3322        )
    3323     {
    3324     return 1;
    3325     }
    3326 #endif
     3610       && (G->m[i]->next->next->next->next!=NULL) /* len >=4*/ )
     3611    {
     3612      return 1;
     3613    }
    33273614  }
    33283615  return 0;
     
    34143701  for(i=nG-1; i>=0; i--)
    34153702  {
    3416 #if 0
     3703/*
    34173704    poly t;
    34183705    if((t=pSub(pCopy(H0->m[i]), pCopy(H1->m[i]))) != NULL)
     
    34223709    }
    34233710    pDelete(&t);
    3424 #else
     3711*/
    34253712    if(!pEqualPolys(H0->m[i],H1->m[i]))
    34263713    {
    34273714      return 0;
    34283715    }
    3429 #endif
    34303716  }
    34313717  return 1;
     
    34363722 * find the maximal total degree of polynomials in G *
    34373723 *****************************************************/
    3438 #if 0
     3724/*
    34393725static int Trandegreebound(ideal G)
    34403726{
     
    34573743  return result;
    34583744}
    3459 #endif
     3745*/
    34603746
    34613747//unused
     
    38044090 * basis or n times, where n is the numbers of variables.                    *
    38054091 *****************************************************************************/
    3806 
    3807 //unused
    3808 #if 0
    3809 static int testnegintvec(intvec* v)
    3810 {
    3811   int n = v->length();
    3812   int i;
    3813   for(i=0; i<n; i++)
    3814   {
    3815     if((*v)[i]<0)
    3816     {
    3817       return(1);
    3818     }
    3819   }
    3820   return(0);
    3821 }
    3822 #endif
    38234092
    38244093// npwinc = 0, if curr_weight doesn't stay in the correct Groebner cone
     
    40324301        //nOverflow_Error = Overflow_Error;
    40334302        tproc=tproc+clock()-tinput;
    4034         /*
    4035           Print("\n// takes %d steps and calls \"Rec_LastGB\" (%d):",
    4036           nwalk, tp_deg+1);
    4037         */
     4303       
     4304        Print("\n// takes %d steps and calls \"Rec_LastGB\" (%d):",
     4305        nwalk, tp_deg+1);
     4306       
    40384307        G = Rec_LastGB(G,curr_weight, orig_target_weight, tp_deg+1,nnwinC);
    40394308        newRing = currRing;
     
    40694338    {
    40704339      // nOverflow_Error = Overflow_Error;
    4071       //Print("\n//  takes %d steps and calls \"Rec_LastGB (%d):", tp_deg+1);
     4340      Print("\n//  takes %d steps and calls \"Rec_LastGB (%d):", tp_deg+1);
    40724341      tproc=tproc+clock()-tinput;
    40734342      F1 = Rec_LastGB(F1,curr_weight, orig_target_weight, tp_deg+1,nnwinC);
     
    41234392    Overflow_Error=nError;
    41244393    }
    4125 // Print("\n// \"Rec_LastGB\" (%d) took %d steps and %.2f sec.Overflow_Error (%d)", tp_deg, nwalk, ((double) tproc)/1000000, nOverflow_Error);
     4394#ifdef TIME_TEST
     4395   Print("\n// \"Rec_LastGB\" (%d) took %d steps and %.2f sec.Overflow_Error (%d)", tp_deg, nwalk, ((double) tproc)/1000000, nOverflow_Error);
     4396#endif
    41264397  return(result);
    41274398}
     
    41454416  //BOOLEAN nOverflow_Error = FALSE;
    41464417  //Print("// pSetm_Error = (%d)", ErrorCheck());
    4147 
     4418#ifdef TIME_TEST
    41484419  xtif=0; xtstd=0; xtlift=0; xtred=0; xtnw=0; xtextra=0;
    41494420  xftinput = clock();
    41504421  clock_t tostd, tproc;
    4151 
     4422#endif
    41524423  nstep = 0;
    41534424  int i, nV = currRing->N;
     
    41604431  intvec* ivNull = new intvec(nV);
    41614432  intvec* next_weight;
    4162 #if 0
    4163   intvec* extra_curr_weight = new intvec(nV);
    4164 #endif
     4433  //intvec* extra_curr_weight = new intvec(nV);
    41654434  //intvec* hilb_func;
    41664435  intvec* exivlp = Mivlp(nV);
    4167 
    41684436  ring XXRing = currRing;
    41694437
    41704438  //Print("\n// ring r_input = %s;", rString(currRing));
     4439#ifdef TIME_TEST
    41714440  to = clock();
     4441#endif
    41724442  /* compute the reduced Groebner basis of the given ideal w.r.t.
    41734443     a "fast" monomial order, e.g. degree reverse lex. order (dp) */
    41744444  G = MstdCC(Go);
     4445#ifdef TIME_TEST
    41754446  tostd=clock()-to;
    41764447
    4177   /*
    41784448  Print("\n// Computation of the first std took = %.2f sec",
    41794449        ((double) tostd)/1000000);
    4180   */
     4450#endif
    41814451  if(currRing->order[0] == ringorder_a)
    41824452  {
     
    41874457    nwalk ++;
    41884458    nstep ++;
     4459#ifdef TIME_TEST
    41894460    to = clock();
     4461#endif
    41904462    /* compute an initial form ideal of <G> w.r.t. "curr_vector" */
    41914463    Gomega = MwalkInitialForm(G, curr_weight);
     4464#ifdef TIME_TEST
    41924465    xtif=xtif+clock()-to;
    4193 #if 0
     4466#endif
     4467/*
    41944468    if(Overflow_Error == TRUE)
    41954469    {
     
    41994473      goto LAST_GB_ALT2;
    42004474    }
    4201 #endif
     4475*/
    42024476    oldRing = currRing;
    42034477
     
    42134487    newRing = currRing;
    42144488    Gomega1 = idrMoveR(Gomega, oldRing,currRing);
     4489#ifdef TIME_TEST
    42154490    to = clock();
     4491#endif
    42164492    /* compute a reduced Groebner basis of <Gomega> w.r.t. "newRing" */
    42174493    M = MstdhomCC(Gomega1);
     4494#ifdef TIME_TEST
    42184495    xtstd=xtstd+clock()-to;
     4496#endif
    42194497    /* change the ring to oldRing */
    42204498    rChangeCurrRing(oldRing);
    42214499    M1 =  idrMoveR(M, newRing,currRing);
    42224500    Gomega2 =  idrMoveR(Gomega1, newRing,currRing);
    4223 
     4501#ifdef TIME_TEST
    42244502    to = clock();
     4503#endif
    42254504    /* compute the reduced Groebner basis of <G> w.r.t. "newRing"
    42264505       by the liftig process */
    42274506    F = MLifttwoIdeal(Gomega2, M1, G);
     4507#ifdef TIME_TEST
    42284508    xtlift=xtlift+clock()-to;
     4509#endif
    42294510    idDelete(&M1);
    42304511    idDelete(&Gomega2);
     
    42344515    rChangeCurrRing(newRing);
    42354516    F1 = idrMoveR(F, oldRing,currRing);
    4236 
     4517#ifdef TIME_TEST
    42374518    to = clock();
     4519#endif
    42384520    /* reduce the Groebner basis <G> w.r.t. newRing */
    42394521    G = kInterRedCC(F1, NULL);
     4522#ifdef TIME_TEST
    42404523    xtred=xtred+clock()-to;
     4524#endif
    42414525    idDelete(&F1);
    42424526
     
    42454529
    42464530  NEXT_VECTOR:
     4531#ifdef TIME_TEST
    42474532    to = clock();
     4533#endif
    42484534    /* compute a next weight vector */
    42494535    next_weight = MkInterRedNextWeight(curr_weight,target_weight, G);
     4536#ifdef TIME_TEST
    42504537    xtnw=xtnw+clock()-to;
     4538#endif
    42514539#ifdef PRINT_VECTORS
    42524540    MivString(curr_weight, target_weight, next_weight);
     
    42924580     // LAST_GB_ALT2:
    42934581        //nOverflow_Error = Overflow_Error;
     4582#ifdef TIME_TEST
    42944583        tproc = clock()-xftinput;
     4584#endif
    42954585        //Print("\n// takes %d steps and calls the recursion of level 2:",  nwalk);
    42964586        /* call the changed perturbation walk algorithm with degree 2 */
     
    43194609
    43204610#ifdef TIME_TEST
    4321  // Print("\n// \"Main procedure\"  took %d steps dnd %.2f sec. Overflow_Error (%d)", nwalk, ((double) tproc)/1000000, nOverflow_Error);
     4611  Print("\n// \"Main procedure\"  took %d steps dnd %.2f sec. Overflow_Error (%d)",
     4612        nwalk, ((double) tproc)/1000000, nOverflow_Error);
    43224613
    43234614  TimeStringFractal(xftinput, tostd, xtif, xtstd, xtextra,xtlift, xtred,xtnw);
     
    43554646 * compute a next weight vector *
    43564647 ********************************/
    4357 static intvec* MWalkRandomNextWeight(ideal G, intvec* curr_weight,
    4358                intvec* target_weight, int weight_rad, int pert_deg)
     4648static intvec* MWalkRandomNextWeight(ideal G, intvec* orig_M, intvec* target_weight,
     4649       int weight_rad, int pert_deg)
    43594650{
    4360   assume(currRing != NULL && curr_weight != NULL &&
     4651  assume(currRing != NULL && orig_M != NULL &&
    43614652         target_weight != NULL && G->m[0] != NULL);
    43624653
     4654  //BOOLEAN nError = Overflow_Error;
     4655  Overflow_Error = FALSE;
     4656
     4657  BOOLEAN found_random_weight = FALSE;
     4658  int i,nV = currRing->N;
     4659  intvec* curr_weight = new intvec(nV);
     4660
     4661  for(i=0; i<nV; i++)
     4662  {
     4663    (*curr_weight)[i] = (*orig_M)[i];
     4664  }
     4665
     4666  int k=0,weight_norm;
     4667  intvec* next_weight;
    43634668  intvec* next_weight1 = MkInterRedNextWeight(curr_weight,target_weight,G);
    4364   if(weight_rad == 0)
    4365   {
    4366     return(next_weight1);
    4367   }
    4368 
    4369   int i,weight_norm,nV = currRing->N;
    43704669  intvec* next_weight2 = new intvec(nV);
    43714670  intvec* next_weight22 = new intvec(nV);
    43724671  intvec* result = new intvec(nV);
    4373 
    4374   ideal G_test, G_test2;
    4375 
    4376   //compute a random next weight vector "next_weight2"
    4377   while(1)
     4672  intvec* curr_weight1;
     4673  ideal G_test, G_test1, G_test2;
     4674
     4675  //try to find a random next weight vector "next_weight2"
     4676  if(weight_rad > 0){ while(k<10)
    43784677  {
    43794678    weight_norm = 0;
     
    44014700    }
    44024701   
    4403     if(test_w_in_ConeCC(G, next_weight2) == 1)
    4404     {
    4405       PrintS("\n Hier 1\n");
    4406       G_test2 = MwalkInitialForm(G, next_weight2);
    4407       PrintS("\n Hier 2\n");
    4408       if(maxlengthpoly(G_test2) <= 1)
     4702    if(test_w_in_ConeCC(G,next_weight2) == 1)
     4703    {
     4704      if(maxlengthpoly(MwalkInitialForm(G,next_weight2))<2)
    44094705      {
    44104706        next_weight2 = MkInterRedNextWeight(next_weight2,target_weight,G);
    44114707      }
    4412       idDelete(&G_test2);
    4413 
     4708/*
    44144709      if(MivAbsMax(next_weight2)>1147483647)
    44154710      {
     
    44284723        delete next_weight22;
    44294724      }
     4725*/
     4726      G_test2 = MwalkInitialForm(G, next_weight2);
     4727      found_random_weight = TRUE;
    44304728      break;
    44314729    }
    4432     weight_rad++;
    4433   }
    4434  
    4435   // compute "usual" next weight vector
    4436   intvec* next_weight = MwalkNextWeightCC(curr_weight,target_weight, G);
    4437   G_test = MwalkInitialForm(G, next_weight);
    4438   G_test2 = MwalkInitialForm(G, next_weight2);
    4439 
     4730    k++;
     4731  }}
     4732Print("\n MWalkRandomNextWeight: compute perurbation...\n");
     4733  // compute "perturbed" next weight vector
     4734  if(pert_deg > 1)
     4735  {
     4736    curr_weight1 = MPertVectors(G,orig_M,pert_deg);
     4737    next_weight = MkInterRedNextWeight(curr_weight1,target_weight,G);
     4738    delete curr_weight1;
     4739  }
     4740  else
     4741  {
     4742    next_weight = MkInterRedNextWeight(curr_weight,target_weight,G);
     4743  }
     4744  if(MivSame(curr_weight,next_weight)==1 || Overflow_Error == TRUE)
     4745  {
     4746    Overflow_Error = FALSE;
     4747    delete next_weight;
     4748    next_weight = MkInterRedNextWeight(curr_weight,target_weight,G);
     4749  }
     4750  G_test=MwalkInitialForm(G,next_weight);
     4751  G_test1=MwalkInitialForm(G,next_weight1);
     4752Print("\n MWalkRandomNextWeight: finished...\n");
    44404753  // compare next weights
    44414754  if(Overflow_Error == FALSE)
    44424755  {
    4443     ideal G_test1 = MwalkInitialForm(G, next_weight1);
    4444     if(G_test1->m[0] != NULL && maxlengthpoly(G_test1) < maxlengthpoly(G_test))
    4445     {
    4446       if(G_test2->m[0] != NULL && maxlengthpoly(G_test2) < maxlengthpoly(G_test1))
    4447       {
    4448         for(i=0; i<nV; i++)
     4756    if(found_random_weight == TRUE)
     4757    {
     4758    // random next weight vector found
     4759      if(G_test1->m[0] != NULL && maxlengthpoly(G_test1) < maxlengthpoly(G_test))
     4760      {
     4761        if(G_test2->m[0] != NULL && maxlengthpoly(G_test2) < maxlengthpoly(G_test1))
    44494762        {
    4450           (*result)[i] = (*next_weight2)[i];
     4763          for(i=0; i<nV; i++)
     4764          {
     4765            (*result)[i] = (*next_weight2)[i];
     4766          }
    44514767        }
     4768        else
     4769        {
     4770          for(i=0; i<nV; i++)
     4771          {
     4772            (*result)[i] = (*next_weight1)[i];
     4773          }
     4774        }   
    44524775      }
    44534776      else
    44544777      {
    4455         for(i=0; i<nV; i++)
     4778        if(G_test2->m[0] != NULL && maxlengthpoly(G_test2) < maxlengthpoly(G_test))
     4779        {
     4780          for(i=0; i<nV; i++)
     4781          {
     4782            (*result)[i] = (*next_weight2)[i];
     4783          }
     4784        }
     4785        else
     4786        {
     4787          for(i=0; i<nV; i++)
     4788          {
     4789            (*result)[i] = (*next_weight)[i];
     4790          }
     4791        }
     4792      }
     4793    }
     4794    else
     4795    {
     4796      // no random next weight vector found
     4797      if(G_test1->m[0] != NULL && maxlengthpoly(G_test1) < maxlengthpoly(G_test))
     4798      {
     4799       for(i=0; i<nV; i++)
    44564800        {
    44574801          (*result)[i] = (*next_weight1)[i];
    4458         }
    4459       }   
    4460     }
    4461     else
    4462     {
    4463       if(G_test2->m[0] != NULL && maxlengthpoly(G_test2) < maxlengthpoly(G_test))
    4464       {
    4465         for(i=0; i<nV; i++)
    4466         {
    4467           (*result)[i] = (*next_weight2)[i];
    44684802        }
    44694803      }
     
    44764810      }
    44774811    }
    4478     idDelete(&G_test1);
    44794812  }
    44804813  else
    44814814  {
    44824815    Overflow_Error = FALSE;
    4483     if(G_test2->m[0] != NULL && maxlengthpoly(G_test2) < maxlengthpoly(G_test))
    4484     {
    4485       for(i=1; i<nV; i++)
    4486       {
    4487         (*result)[i] = (*next_weight2)[i];
     4816    if(found_random_weight == TRUE)
     4817    {
     4818      if(G_test2->m[0] != NULL && maxlengthpoly(G_test2) < maxlengthpoly(G_test))
     4819      {
     4820        for(i=1; i<nV; i++)
     4821        {
     4822          (*result)[i] = (*next_weight2)[i];
     4823        }
     4824      }
     4825      else
     4826      {
     4827        for(i=0; i<nV; i++)
     4828        {
     4829          (*result)[i] = (*next_weight)[i];
     4830        }
    44884831      }
    44894832    }
     
    44964839    }
    44974840  }
    4498 
     4841 // delete curr_weight1;
     4842  delete next_weight;
     4843  delete next_weight2;
    44994844  idDelete(&G_test);
    4500   idDelete(&G_test2);
    4501   if(test_w_in_ConeCC(G, result) == 1)
    4502   {
    4503     delete next_weight2;
    4504     delete next_weight;
     4845  idDelete(&G_test1);
     4846  if(found_random_weight == TRUE)
     4847  {
     4848    idDelete(&G_test2);
     4849  }
     4850  if(test_w_in_ConeCC(G, result) == 1 && MivSame(curr_weight,result)==0)
     4851  {
     4852    delete curr_weight;
    45054853    delete next_weight1;
    45064854    return result;
     
    45084856  else
    45094857  {
     4858    delete curr_weight;
    45104859    delete result;
    4511     delete next_weight2;
    4512     delete next_weight1;
    4513     return next_weight;
     4860    return next_weight1;
    45144861  }
    45154862}
     
    45194866 * The procedur REC_GB_Mwalk computes a GB for <G> w.r.t. the weight order *
    45204867 * otw, where G is a reduced GB w.r.t. the weight order cw.                *
    4521  * The new procedur Mwalk calls REC_GB.                                    *
     4868 * The new procedure Mwalk calls REC_GB.                                   *
    45224869 ***************************************************************************/
    45234870static ideal REC_GB_Mwalk(ideal G, intvec* curr_weight, intvec* orig_target_weight,
     
    48225169
    48235170  ideal Gomega, M, F, Gomega1, Gomega2, M1, F1, G;
    4824   //ideal G1;
    4825   //ring endRing;
     5171
    48265172  ring newRing, oldRing;
    48275173  intvec* ivNull = new intvec(nV);
     
    48655211         the recursive changed perturbation walk alg. */
    48665212      tim = clock();
    4867       /*
     5213#ifdef CHECK_IDEAL_MWALK
    48685214        Print("\n// **** Groebnerwalk took %d steps and ", nwalk);
    48695215        PrintS("\n// **** call the rec. Pert. Walk to compute a red GB of:");
    4870         idElements(Gomega, "G_omega");
    4871       */
     5216        idString(Gomega, "Gomega");
     5217#endif
    48725218
    48735219      if(MivSame(exivlp, target_weight)==1)
     
    48755221      else
    48765222        goto NORMAL_GW;
    4877       /*
     5223#ifdef TIME_TEST
    48785224        Print("\n//  time for the last std(Gw)  = %.2f sec",
    48795225        ((double) (clock()-tim)/1000000));
    4880         PrintS("\n// ***************************************************\n");
    4881       */
     5226#endif
     5227/*
    48825228#ifdef CHECK_IDEAL_MWALK
    48835229      idElements(Gomega, "G_omega");
     
    48865232      //headidString(M, "M");
    48875233#endif
     5234*/
    48885235      to = clock();
    48895236      F = MLifttwoIdeal(Gomega, M, G);
     
    50455392}
    50465393
    5047 
    50485394/*******************************
    50495395 * THE GROEBNER WALK ALGORITHM *
     
    51035449#endif
    51045450  rComplete(currRing);
    5105 //#ifdef CHECK_IDEAL_MWALK
     5451#ifdef CHECK_IDEAL_MWALK
    51065452  if(printout > 2)
    51075453  {
    51085454    idString(Go,"//** Mwalk: Go");
    51095455  }
    5110 //#endif
     5456#endif
    51115457
    51125458  if(target_M->length() == nV)
     
    51335479  to = clock();
    51345480#endif
    5135 
    51365481  ideal G = MstdCC(idrMoveR(Go,baseRing,currRing));
    5137 
    51385482#ifdef TIME_TEST
    51395483  tostd = clock()-to;
     
    51475491    nwalk ++;
    51485492    nstep ++;
    5149 
     5493    //compute an initial form ideal of <G> w.r.t. "curr_vector"
    51505494#ifdef TIME_TEST
    51515495    to = clock();
    51525496#endif
    5153     // compute an initial form ideal of <G> w.r.t. "curr_vector"
    51545497    Gomega = MwalkInitialForm(G, curr_weight);
    51555498#ifdef TIME_TEST
     
    51575500#endif
    51585501
    5159 //#ifdef CHECK_IDEAL_MWALK
     5502#ifdef CHECK_IDEAL_MWALK
    51605503    if(printout > 1)
    51615504    {
    51625505      idString(Gomega,"//** Mwalk: Gomega");
    51635506    }
    5164 //#endif
     5507#endif
    51655508
    51665509    if(reduction == 0)
    51675510    {
    51685511      FF = middleOfCone(G,Gomega);
    5169       if( FF != NULL)
     5512      if(FF != NULL)
    51705513      {
    51715514        idDelete(&G);
     
    52295572#endif
    52305573    idSkipZeroes(M);
    5231 //#ifdef CHECK_IDEAL_MWALK
     5574#ifdef CHECK_IDEAL_MWALK
    52325575    if(printout > 2)
    52335576    {
    52345577      idString(M, "//** Mwalk: M");
    52355578    }
    5236 //#endif
     5579#endif
    52375580    //change the ring to baseRing
    52385581    rChangeCurrRing(baseRing);
     
    52475590    // where Gomega is a reduced Groebner basis w.r.t. the current ring
    52485591    F = MLifttwoIdeal(Gomega2, M1, G);
    5249 
    52505592#ifdef TIME_TEST
    52515593    tlift = tlift + clock() - to;
    52525594#endif
    5253 //#ifdef CHECK_IDEAL_MWALK
     5595#ifdef CHECK_IDEAL_MWALK
    52545596    if(printout > 2)
    52555597    {
    52565598      idString(F, "//** Mwalk: F");
    52575599    }
    5258 //#endif
     5600#endif
    52595601    idDelete(&Gomega2);
    52605602    idDelete(&M1);
     
    52655607    idSkipZeroes(G);
    52665608
    5267 //#ifdef CHECK_IDEAL_MWALK
     5609#ifdef CHECK_IDEAL_MWALK
    52685610    if(printout > 2)
    52695611    {
    52705612      idString(G, "//** Mwalk: G");
    52715613    }
    5272 //#endif
     5614#endif
    52735615
    52745616    rChangeCurrRing(targetRing);
     
    52865628    baseRing = currRing;
    52875629
     5630    NEXT_VECTOR:
     5631#ifdef TIME_TEST
     5632    to = clock();
     5633#endif
     5634    intvec* next_weight = MwalkNextWeightCC(curr_weight,target_weight,G);
     5635#ifdef TIME_TEST
     5636    tnw = tnw + clock() - to;
     5637#endif
     5638#ifdef PRINT_VECTORS
     5639    if(printout > 0)
     5640    {
     5641      MivString(curr_weight, target_weight, next_weight);
     5642    }
     5643#endif
     5644    if(MivComp(target_weight,curr_weight) == 1)// || endwalks == TRUE)
     5645    {
    52885646/*
    5289 #ifdef TIME_TEST
    5290     to = clock();
    5291 #endif
    5292 
    5293 #ifdef TIME_TEST
    5294     tstd = tstd + clock() - to;
    5295 #endif
    5296 */
    5297 
    5298 
    5299 #ifdef TIME_TEST
    5300     to = clock();
    5301 #endif
    5302     NEXT_VECTOR:
    5303     intvec* next_weight = MwalkNextWeightCC(curr_weight,target_weight,G);
    5304 #ifdef TIME_TEST
    5305     tnw = tnw + clock() - to;
    5306 #endif
    5307 //#ifdef PRINT_VECTORS
    5308     if(printout > 0)
    5309     {
    5310       MivString(curr_weight, target_weight, next_weight);
    5311     }
    5312 //#endif
    5313     if(MivComp(target_weight,curr_weight) == 1)// || endwalks == TRUE)
    5314     {/*
    5315 //#ifdef CHECK_IDEAL_MWALK
     5647#ifdef CHECK_IDEAL_MWALK
    53165648      if(printout > 0)
    53175649      {
    53185650        PrintS("\n//** Mwalk: entering last cone.\n");
    53195651      }
    5320 //#endif
     5652#endif
    53215653
    53225654      Gomega = MwalkInitialForm(G, curr_weight); // compute an initial form ideal of <G> w.r.t. "curr_vector"
     
    53325664      Gomega1 = idrMoveR(Gomega, baseRing,currRing);
    53335665      idDelete(&Gomega);
    5334 //#ifdef CHECK_IDEAL_MWALK
     5666#ifdef CHECK_IDEAL_MWALK
    53355667      if(printout > 1)
    53365668      {
     
    53385670      }
    53395671      PrintS("\n //** Mwalk: kStd(Gomega)");
    5340 //#endif
     5672#endif
    53415673      M = kStd(Gomega1,NULL,testHomog,NULL,NULL,0,0,NULL);
    5342 //#ifdef CHECK_IDEAL_MWALK
     5674#ifdef CHECK_IDEAL_MWALK
    53435675      if(printout > 1)
    53445676      {
    53455677        idString(M,"//** Mwalk: M");
    53465678      }
    5347 //#endif
     5679#endif
    53485680      rChangeCurrRing(baseRing);
    53495681      M1 =  idrMoveR(M, newRing,currRing);
     
    53535685      //PrintS("\n //** Mwalk: MLifttwoIdeal");
    53545686      F = MLifttwoIdeal(Gomega2, M1, G);
    5355 //#ifdef CHECK_IDEAL_MWALK
     5687#ifdef CHECK_IDEAL_MWALK
    53565688      if(printout > 2)
    53575689      {
    53585690        idString(F,"//** Mwalk: F");
    53595691      }
    5360 //#endif
     5692#endif
    53615693      idDelete(&Gomega2);
    53625694      idDelete(&M1);
     
    53655697      idDelete(&F);
    53665698      baseRing = currRing;
    5367       si_opt_1 = save1; //set original options, e. g. option(RedSB)
    53685699      idSkipZeroes(G);
    53695700#ifdef TIME_TEST
     
    53775708#endif
    53785709      idSkipZeroes(G);
    5379       delete next_weight; */
     5710      delete next_weight;
     5711*/
    53805712      break;
    53815713    }
     
    54045736#endif
    54055737  Print("\n//** Mwalk: Groebner Walk took %d steps.\n", nstep);
     5738  si_opt_1 = save1; //set original options
    54065739  return(result);
    54075740}
     
    54175750    si_opt_1 &= (~Sy_bit(OPT_REDTAIL)); // not tail reductions
    54185751  }
     5752
    54195753  Set_Error(FALSE);
    54205754  Overflow_Error = FALSE;
     
    54275761#endif
    54285762  nstep=0;
    5429   int i,polylength,nwalk;
     5763  int i,nwalk;//polylength;
    54305764  int nV = currRing->N;
     5765
     5766  //check that weight radius is valid
     5767  if(weight_rad < 0)
     5768  {
     5769    Werror("Invalid radius.\n");
     5770    return NULL;
     5771  }
     5772
     5773  //check that perturbation degree is valid
     5774  if(pert_deg > nV || pert_deg < 1)
     5775  {
     5776    Werror("Invalid perturbation degree.\n");
     5777    return NULL;
     5778  }
    54315779
    54325780  ideal Gomega, M, F,FF, Gomega1, Gomega2, M1;
     
    54355783  ring baseRing = currRing;
    54365784  ring XXRing = currRing;
     5785  intvec* iv_M;
    54375786  intvec* ivNull = new intvec(nV);
    54385787  intvec* curr_weight = new intvec(nV);
     
    54455794    (*target_weight)[i] = (*target_M)[i];
    54465795  }
     5796
    54475797#ifndef BUCHBERGER_ALG
    54485798  intvec* hilb_func;
     
    54565806#endif
    54575807  rComplete(currRing);
     5808
     5809  if(target_M->length() == nV)
     5810  {
     5811    targetRing = VMrDefault(target_weight); // define the target ring
     5812  }
     5813  else
     5814  {
     5815    targetRing = VMatrDefault(target_M);
     5816  }
     5817  if(orig_M->length() == nV)
     5818  {
     5819    newRing = VMrDefault(curr_weight); // define a new ring with ordering "(a(curr_weight),lp)
     5820  }
     5821  else
     5822  {
     5823    newRing = VMatrDefault(orig_M);
     5824  }
     5825  rChangeCurrRing(newRing);
    54585826#ifdef TIME_TEST
    54595827  to = clock();
    54605828#endif
    5461 
    5462   if(target_M->length() == nV)
    5463   {
    5464    // define the target ring
    5465     targetRing = VMrDefault(target_weight);
    5466   }
    5467   else
    5468   {
    5469     targetRing = VMatrDefault(target_M);
    5470   }
    5471   if(orig_M->length() == nV)
    5472   {
    5473     newRing = VMrDefault(curr_weight); // define a new ring with ordering "(a(curr_weight),lp)
    5474   }
    5475   else
    5476   {
    5477     newRing = VMatrDefault(orig_M);
    5478   }
    5479   rChangeCurrRing(newRing);
    54805829  ideal G = MstdCC(idrMoveR(Go,baseRing,currRing));
     5830#ifdef TIME_TEST
     5831  tostd = clock()-to;
     5832#endif
    54815833  baseRing = currRing;
    5482 #ifdef TIME_TEST
    5483   tostd = clock()-to;
    5484 #endif
    5485 
    54865834  nwalk = 0;
    54875835
     
    54985846    nwalk ++;
    54995847    nstep ++;
     5848#ifdef CHECK_IDEAL_MWALK
    55005849    if(printout > 1)
    55015850    {
    55025851      idString(Gomega,"//** Mrwalk: Gomega");
    55035852    }
     5853#endif
    55045854    if(reduction == 0)
    55055855    {
     
    55105860        G = idCopy(FF);
    55115861        idDelete(&FF);
    5512        
    55135862        goto NEXT_VECTOR;
    55145863      }
     
    55535902    to = clock();
    55545903#endif
    5555 #ifndef  BUCHBERGER_ALG
     5904#ifndef BUCHBERGER_ALG
    55565905    M=kStd(Gomega1,NULL,isHomog,NULL,hilb_func,0,NULL,curr_weight);
    55575906    delete hilb_func;
     
    55635912#endif
    55645913    idSkipZeroes(M);
    5565 //#ifdef CHECK_IDEAL_MWALK
     5914#ifdef CHECK_IDEAL_MWALK
    55665915    if(printout > 2)
    55675916    {
    55685917      idString(M, "//** Mrwalk: M");
    55695918    }
    5570 //#endif
     5919#endif
    55715920    //change the ring to baseRing
    55725921    rChangeCurrRing(baseRing);
     
    55845933    tlift = tlift + clock() - to;
    55855934#endif
    5586 //#ifdef CHECK_IDEAL_MWALK
     5935#ifdef CHECK_IDEAL_MWALK
    55875936    if(printout > 2)
    55885937    {
    5589       idString(F, "//** Mrwalk: F");
    5590     }
    5591 //#endif
     5938      idString(F,"//** Mrwalk: F");
     5939    }
     5940#endif
    55925941    idDelete(&Gomega2);
    55935942    idDelete(&M1);
     
    56015950#endif
    56025951    idSkipZeroes(G);
    5603 //#ifdef CHECK_IDEAL_MWALK
     5952#ifdef CHECK_IDEAL_MWALK
    56045953    if(printout > 2)
    56055954    {
    5606       idString(G, "//** Mrwalk: G");
    5607     }
    5608 //#endif
     5955      idString(G,"//** Mrwalk: G");
     5956    }
     5957#endif
    56095958
    56105959    rChangeCurrRing(targetRing);
     
    56395988#endif
    56405989
    5641     //polylength = 1 if there is a polynomial in Gomega with at least 3 monomials and 0 otherwise
    5642     polylength = lengthpoly(Gomega);
    5643     if(polylength > 0)
     5990    //lengthpoly(Gomega) = 1 if there is a polynomial in Gomega with at least 3 monomials and 0 otherwise
     5991    //polylength = lengthpoly(Gomega);
     5992    if(lengthpoly(Gomega) > 0)
    56445993    {
    56455994      //there is a polynomial in Gomega with at least 3 monomials,
    56465995      //low-dimensional facet of the cone
    56475996      delete next_weight;
    5648 #ifdef TIME_TEST
    5649     to = clock();
    5650 #endif
    5651       next_weight = MWalkRandomNextWeight(G, curr_weight, target_weight, weight_rad, pert_deg);
    5652 #ifdef TIME_TEST
    5653     tnw = tnw + clock() - to;
     5997      if(target_M->length() == nV)
     5998      {
     5999        iv_M = MivMatrixOrder(curr_weight);
     6000      }
     6001      else
     6002      {
     6003        iv_M = MivMatrixOrderRefine(curr_weight,target_M);
     6004      }
     6005#ifdef TIME_TEST
     6006      to = clock();
     6007#endif
     6008      next_weight = MWalkRandomNextWeight(G, iv_M, target_weight, weight_rad, pert_deg);
     6009#ifdef TIME_TEST
     6010      tnw = tnw + clock() - to;
    56546011#endif
    56556012      idDelete(&Gomega);
    56566013#ifdef TIME_TEST
    5657     to = clock();
     6014      to = clock();
    56586015#endif
    56596016      Gomega = MwalkInitialForm(G, next_weight);
    56606017#ifdef TIME_TEST
    5661     tif = tif + clock()-to; //time for computing initial form ideal
    5662 #endif
     6018      tif = tif + clock()-to; //time for computing initial form ideal
     6019#endif
     6020      delete iv_M;
    56636021    }
    56646022
     
    56716029    }
    56726030
    5673 //#ifdef PRINT_VECTORS
     6031#ifdef PRINT_VECTORS
    56746032    if(printout > 0)
    56756033    {
    56766034      MivString(curr_weight, target_weight, next_weight);
    56776035    }
    5678 //#endif
     6036#endif
    56796037
    56806038    for(i=nV-1; i>=0; i--)
     
    56886046  ideal result = idrMoveR(G,baseRing,currRing);
    56896047  idDelete(&G);
    5690   si_opt_1 = save1; //set original options, e. g. option(RedSB)
    56916048  delete ivNull;
    56926049#ifndef BUCHBERGER_ALG
     
    56996056  //Print("\n// Overflow_Error? (%d)\n", Overflow_Error);
    57006057#endif
     6058  si_opt_1 = save1; //set original options
    57016059  return(result);
    57026060}
    5703 
    5704 //unused
    5705 #if 0
    5706 ideal Mwalk_tst(ideal Go, intvec* curr_weight, intvec* target_weight)
    5707 {
    5708   //clock_t tinput=clock();
    5709   //idString(Go,"Ginp");
    5710   int i, nV = currRing->N;
    5711   int nwalk=0, endwalks=0;
    5712 
    5713   ideal Gomega, M, F, Gomega1, Gomega2, M1, F1, G;
    5714   // ideal G1; ring endRing;
    5715   ring newRing, oldRing;
    5716   intvec* ivNull = new intvec(nV);
    5717   ring XXRing = currRing;
    5718 
    5719   intvec* tmp_weight = new intvec(nV);
    5720   for(i=nV-1; i>=0; i--)
    5721   {
    5722     (*tmp_weight)[i] = (*curr_weight)[i];
    5723   }
    5724   /* the monomial ordering of this current ring would be "dp" */
    5725   G = MstdCC(Go);
    5726 #ifndef BUCHBERGER_ALG
    5727   intvec* hilb_func;
    5728 #endif
    5729   /* to avoid (1,0,...,0) as the target vector */
    5730   intvec* last_omega = new intvec(nV);
    5731   for(i=nV-1; i>0; i--)
    5732     (*last_omega)[i] = 1;
    5733   (*last_omega)[0] = 10000;
    5734 
    5735   while(1)
    5736   {
    5737     nwalk ++;
    5738     //Print("\n// Entering the %d-th step:", nwalk);
    5739     //Print("\n// ring r[%d] = %s;", nwalk, rString(currRing));
    5740     idString(G,"G");
    5741     /* compute an initial form ideal of <G> w.r.t. "curr_vector" */
    5742     Gomega = MwalkInitialForm(G, curr_weight);
    5743     //ivString(curr_weight, "omega");
    5744     idString(Gomega,"Gw");
    5745 
    5746 #ifndef  BUCHBERGER_ALG
    5747     if(isNolVector(curr_weight) == 0)
    5748       hilb_func = hFirstSeries(Gomega,NULL,NULL,curr_weight,currRing);
    5749     else
    5750       hilb_func = hFirstSeries(Gomega,NULL,NULL,last_omega,currRing);
    5751 #endif // BUCHBERGER_ALG
    5752 
    5753 
    5754     oldRing = currRing;
    5755 
    5756     /* define a new ring that its ordering is "(a(curr_weight),lp) */
    5757     VMrDefault(curr_weight);
    5758     newRing = currRing;
    5759 
    5760     Gomega1 = idrMoveR(Gomega, oldRing,currRing);
    5761 
    5762     /* compute a reduced Groebner basis of <Gomega> w.r.t. "newRing" */
    5763 #ifdef  BUCHBERGER_ALG
    5764     M = MstdhomCC(Gomega1);
    5765 #else
    5766     M=kStd(Gomega1,NULL,isHomog,NULL,hilb_func,0,NULL,curr_weight);
    5767     delete hilb_func;
    5768 #endif // BUCHBERGER_ALG
    5769 
    5770     idString(M,"M");
    5771 
    5772       /* change the ring to oldRing */
    5773     rChangeCurrRing(oldRing);
    5774     M1 =  idrMoveR(M, newRing,currRing);
    5775     Gomega2 =  idrMoveR(Gomega1, newRing,currRing);
    5776 
    5777       /* compute a representation of the generators of submod (M)
    5778          with respect to those of mod (Gomega).
    5779          Gomega is a reduced Groebner basis w.r.t. the current ring */
    5780     F = MLifttwoIdeal(Gomega2, M1, G);
    5781     idDelete(&M1);
    5782     idDelete(&Gomega2);
    5783     idDelete(&G);
    5784     idString(F,"F");
    5785 
    5786     /* change the ring to newRing */
    5787     rChangeCurrRing(newRing);
    5788     F1 = idrMoveR(F, oldRing,currRing);
    5789 
    5790     /* reduce the Groebner basis <G> w.r.t. new ring */
    5791     G = kInterRedCC(F1, NULL);
    5792     //idSkipZeroes(G);//done by kInterRed
    5793     idDelete(&F1);
    5794     idString(G,"G");
    5795     if(endwalks == 1)
    5796       break;
    5797 
    5798     /* compute a next weight vector */
    5799     intvec* next_weight = MkInterRedNextWeight(curr_weight,target_weight,G);
    5800 #ifdef PRINT_VECTORS
    5801     MivString(curr_weight, target_weight, next_weight);
    5802 #endif
    5803 
    5804     if(MivComp(next_weight, ivNull) == 1)
    5805     {
    5806       delete next_weight;
    5807       break;
    5808     }
    5809     if(MivComp(next_weight, target_weight) == 1)
    5810       endwalks = 1;
    5811 
    5812     for(i=nV-1; i>=0; i--)
    5813       (*tmp_weight)[i] = (*curr_weight)[i];
    5814 
    5815     /* 06.11.01 to free the memory: NOT Changed!!*/
    5816     for(i=nV-1; i>=0; i--)
    5817       (*curr_weight)[i] = (*next_weight)[i];
    5818     delete next_weight;
    5819   }
    5820   rChangeCurrRing(XXRing);
    5821   G = idrMoveR(G, newRing,currRing);
    5822 
    5823   delete tmp_weight;
    5824   delete ivNull;
    5825   PrintLn();
    5826   return(G);
    5827 }
    5828 #endif
    58296061
    58306062/**************************************************************/
     
    58406072   2) the changed perturbation walk algorithm with a decreased degree.
    58416073*/
    5842 // use kStd, if nP = 0, else call LastGB
     6074// if nP = 0 use kStd, else call LastGB
    58436075ideal Mpwalk(ideal Go, int op_deg, int tp_deg,intvec* curr_weight,
    58446076             intvec* target_weight, int nP, int reduction, int printout)
     
    58486080  {
    58496081    si_opt_1 &= (~Sy_bit(OPT_REDSB)); // no reduced Groebner basis
    5850     //si_opt_1 &= (~Sy_bit(OPT_REDTAIL)); // not tail reductions
    5851     //si_opt_1|=(Sy_bit(OPT_REDTAIL)|Sy_bit(OPT_REDSB));
     6082    si_opt_1 &= (~Sy_bit(OPT_REDTAIL)); // not tail reductions
    58526083  }
    58536084  Set_Error(FALSE  );
    58546085  Overflow_Error = FALSE;
    58556086  //Print("// pSetm_Error = (%d)", ErrorCheck());
    5856 
     6087#ifdef TIME_TEST
    58576088  clock_t  tinput, tostd, tif=0, tstd=0, tlift=0, tred=0, tnw=0;
    58586089  xtextra=0;
     
    58616092
    58626093  clock_t tim;
    5863 
     6094#endif
    58646095  nstep = 0;
    58656096  int i, ntwC=1, ntestw=1,  nV = currRing->N;
     6097
     6098  //check that perturbation degree is valid
     6099  if(op_deg < 1 || tp_deg < 1 || op_deg > nV || tp_deg > nV)
     6100  {
     6101    Werror("Invalid perturbation degree.\n");
     6102    return NULL;
     6103  }
     6104
    58666105  BOOLEAN endwalks = FALSE;
    5867 
    58686106  ideal Gomega, M, F, FF, G, Gomega1, Gomega2, M1,F1,Eresult,ssG;
    58696107  ring newRing, oldRing, TargetRing;
     
    58876125
    58886126  ring XXRing = currRing;
    5889 
     6127#ifdef TIME_TEST
    58906128  to = clock();
     6129#endif
    58916130  // perturbs the original vector
    58926131  if(MivComp(curr_weight, iv_dp) == 1) //rOrdStr(currRing) := "dp"
    58936132  {
    58946133    G = MstdCC(Go);
     6134#ifdef TIME_TEST
    58956135    tostd = clock()-to;
     6136#endif
    58966137    if(op_deg != 1){
    58976138      iv_M_dp = MivMatrixOrderdp(nV);
     
    59106151    G = idrMoveR(Go, XXRing,currRing);
    59116152    G = MstdCC(G);
     6153#ifdef TIME_TEST
    59126154    tostd = clock()-to;
     6155#endif
    59136156    if(op_deg != 1){
    59146157      iv_M_dp = MivMatrixOrder(curr_weight);
     
    59466189    G = idrMoveR(ssG, TargetRing,currRing);
    59476190  }
    5948     if(printout > 0)
    5949     {
    5950       Print("\n//** Mpwalk: Perturbation Walk of degree (%d,%d):",op_deg,tp_deg);
    5951       ivString(curr_weight, "//** Mpwalk: new current weight");
    5952       ivString(target_weight, "//** Mpwalk: new target weight");
    5953     }
     6191  if(printout > 0)
     6192  {
     6193    Print("\n//** Mpwalk: Perturbation Walk of degree (%d,%d):",op_deg,tp_deg);
     6194#ifdef PRINT_VECTORS
     6195    ivString(curr_weight, "//** Mpwalk: new current weight");
     6196    ivString(target_weight, "//** Mpwalk: new target weight");
     6197#endif
     6198  }
    59546199  while(1)
    59556200  {
    59566201    nstep ++;
     6202#ifdef TIME_TEST
    59576203    to = clock();
     6204#endif
    59586205    // compute an initial form ideal of <G> w.r.t. the weight vector
    59596206    // "curr_weight"
    59606207    Gomega = MwalkInitialForm(G, curr_weight);
    5961 //#ifdef CHECK_IDEAL_MWALK
     6208#ifdef TIME_TEST
     6209    tif = tif + clock()-to;
     6210#endif
     6211#ifdef CHECK_IDEAL_MWALK
    59626212    if(printout > 1)
    59636213    {
    59646214      idString(Gomega,"//** Mpwalk: Gomega");
    59656215    }
    5966 //#endif
     6216#endif
    59676217    if(reduction == 0 && nstep > 1)
    59686218    {
     
    59786228          idDelete(&FF);
    59796229          next_weight = MwalkNextWeightCC(curr_weight,target_weight,G);
     6230#ifdef PRINT_VECTORS
    59806231          if(printout > 0)
    59816232          {
    59826233            MivString(curr_weight, target_weight, next_weight);
    59836234          }
     6235#endif
    59846236        }
    59856237        else
     
    59926244        }
    59936245        Gomega = MwalkInitialForm(G, curr_weight);
     6246#ifdef CHECK_IDEAL_MWALK
    59946247        if(printout > 1)
    59956248        {
    59966249          idString(Gomega,"//** Mpwalk: Gomega");
    59976250        }
     6251#endif
    59986252      }
    59996253*/
     
    60126266    {
    60136267      Print("\n// ring r%d = %s;\n", nstep, rString(currRing));
     6268/*
    60146269      idElements(G, "G");
    60156270      headidString(G, "G");
    6016     }
    6017 #endif
    6018 
    6019     tif = tif + clock()-to;
     6271*/
     6272    }
     6273#endif
    60206274
    60216275#ifndef  BUCHBERGER_ALG
     
    60416295    {
    60426296      Print("\n// ring r%d = %s;\n", nstep, rString(currRing));
     6297/*
    60436298      idElements(Gomega1, "Gw");
    60446299      headidString(Gomega1, "headGw");
     6300*/
    60456301      PrintS("\n// compute a rGB of Gw:\n");
    6046 
    60476302#ifndef  BUCHBERGER_ALG
    60486303      ivString(hilb_func, "w");
     
    60506305    }
    60516306#endif
    6052 
     6307#ifdef TIME_TEST
    60536308    tim = clock();
    60546309    to = clock();
     6310#endif
    60556311    // compute a reduced Groebner basis of <Gomega> w.r.t. "newRing"
    60566312#ifdef  BUCHBERGER_ALG
     
    60606316    delete hilb_func;
    60616317#endif
    6062 //#ifdef CHECK_IDEAL_MWALK
    6063       if(printout > 2)
    6064       {
    6065         idString(M,"//** Mpwalk: M");
    6066       }
    6067 //#endif
    6068 
    6069     if(endwalks == TRUE){
     6318
     6319    if(endwalks == TRUE)
     6320    {
     6321#ifdef TIME_TEST
    60706322      xtstd = xtstd+clock()-to;
     6323#endif
    60716324#ifdef ENDWALKS
    60726325      Print("\n// time for the last std(Gw)  = %.2f sec\n",
     
    60756328    }
    60766329    else
     6330    {
     6331#ifdef TIME_TEST
    60776332      tstd=tstd+clock()-to;
    6078 
     6333#endif
     6334    }
     6335#ifdef CHECK_IDEAL_MWALK
     6336    if(printout > 2)
     6337    {
     6338      idString(M,"//** Mpwalk: M");
     6339    }
     6340#endif
    60796341    // change the ring to oldRing
    60806342    rChangeCurrRing(oldRing);
    60816343    M1 =  idrMoveR(M, newRing,currRing);
    60826344    Gomega2 =  idrMoveR(Gomega1, newRing,currRing);
    6083 
     6345#ifdef TIME_TEST
    60846346    to=clock();
     6347#endif
    60856348    /* compute a representation of the generators of submod (M)
    60866349       with respect to those of mod (Gomega).
    60876350       Gomega is a reduced Groebner basis w.r.t. the current ring */
    60886351    F = MLifttwoIdeal(Gomega2, M1, G);
     6352#ifdef TIME_TEST
    60896353    if(endwalks == FALSE)
    60906354      tlift = tlift+clock()-to;
    60916355    else
    60926356      xtlift=clock()-to;
    6093 
    6094 //#ifdef CHECK_IDEAL_MWALK
     6357#endif
     6358#ifdef CHECK_IDEAL_MWALK
    60956359    if(printout > 2)
    60966360    {
    60976361      idString(F,"//** Mpwalk: F");
    60986362    }
    6099 //#endif
     6363#endif
    61006364
    61016365    idDelete(&M1);
     
    61166380        PrintS("\n //** Mpwalk: reduce the Groebner basis.\n");
    61176381      }
     6382#ifdef TIME_TEST
    61186383      to=clock();
     6384#endif
    61196385      G = kInterRedCC(F1, NULL);
     6386#ifdef TIME_TEST
    61206387      if(endwalks == FALSE)
    61216388        tred = tred+clock()-to;
    61226389      else
    61236390        xtred=clock()-to;
     6391#endif
    61246392      idDelete(&F1);
    61256393    }
     
    61286396
    61296397    NEXT_VECTOR:
     6398#ifdef TIME_TEST
    61306399    to=clock();
     6400#endif
    61316401    // compute a next weight vector
    61326402    next_weight = MkInterRedNextWeight(curr_weight,target_weight, G);
     6403#ifdef TIME_TEST
    61336404    tnw=tnw+clock()-to;
    6134 //#ifdef PRINT_VECTORS
     6405#endif
     6406#ifdef PRINT_VECTORS
    61356407    if(printout > 0)
    61366408    {
    61376409      MivString(curr_weight, target_weight, next_weight);
    61386410    }
    6139 //#endif
     6411#endif
    61406412
    61416413    if(Overflow_Error == TRUE)
     
    61796451    TargetRing=currRing;
    61806452    F1 = idrMoveR(G, newRing,currRing);
    6181 #ifdef CHECK_IDEAL
     6453/*
     6454#ifdef CHECK_IDEAL_MWALK
    61826455      headidString(G, "G");
    61836456#endif
    6184 
     6457*/
    61856458
    61866459    // check whether the pertubed target vector stays in the correct cone
     
    62656538  Overflow_Error = FALSE;
    62666539  //Print("// pSetm_Error = (%d)", ErrorCheck());
    6267 
     6540#ifdef TIME_TEST
    62686541  clock_t  tinput, tostd, tif=0, tstd=0, tlift=0, tred=0, tnw=0;
    62696542  xtextra=0;
     
    62726545
    62736546  clock_t tim;
    6274 
     6547#endif
    62756548  nstep = 0;
    6276   int i, ntwC=1, ntestw=1, polylength, nV = currRing->N;
     6549  int i, ntwC=1, ntestw=1, nV = currRing->N; //polylength
     6550
     6551  //check that weight radius is valid
     6552  if(weight_rad < 0)
     6553  {
     6554    Werror("Invalid radius.\n");
     6555    return NULL;
     6556  }
     6557
     6558  //check that perturbation degree is valid
     6559  if(op_deg < 1 || tp_deg < 1 || op_deg > nV || tp_deg > nV)
     6560  {
     6561    Werror("Invalid perturbation degree.\n");
     6562    return NULL;
     6563  }
     6564
    62776565  BOOLEAN endwalks = FALSE;
    62786566
    62796567  ideal Gomega, M, F, FF, G, Gomega1, Gomega2, M1,F1,Eresult,ssG;
    62806568  ring newRing, oldRing, TargetRing;
     6569  intvec* iv_M;
    62816570  intvec* iv_M_dp;
    62826571  intvec* iv_M_lp;
     
    63066595  ring XXRing = currRing;
    63076596
    6308   to = clock();
    63096597  // perturbs the original vector
    63106598  if(orig_M->length() == nV)
     
    63126600    if(MivComp(curr_weight, iv_dp) == 1) //rOrdStr(currRing) := "dp"
    63136601    {
     6602#ifdef TIME_TEST
     6603  to = clock();
     6604#endif
    63146605      G = MstdCC(Go);
     6606#ifdef TIME_TEST
    63156607      tostd = clock()-to;
     6608#endif
    63166609      if(op_deg != 1)
    63176610      {
     
    63296622
    63306623      G = idrMoveR(Go, XXRing,currRing);
     6624#ifdef TIME_TEST
     6625  to = clock();
     6626#endif
    63316627      G = MstdCC(G);
     6628#ifdef TIME_TEST
    63326629      tostd = clock()-to;
     6630#endif
    63336631      if(op_deg != 1)
    63346632      {
     
    63426640    rChangeCurrRing(VMatrDefault(orig_M));
    63436641    G = idrMoveR(Go, XXRing,currRing);
     6642#ifdef TIME_TEST
     6643    to = clock();
     6644#endif
    63446645    G = MstdCC(G);
     6646#ifdef TIME_TEST
    63456647    tostd = clock()-to;
     6648#endif
    63466649    if(op_deg != 1)
    63476650    {
     
    64116714  {
    64126715    nstep ++;
    6413     to = clock();
    6414     // compute an initial form ideal of G w.r.t. the weight vector "curr_weight"
    6415 /*    Gomega = MwalkInitialForm(G, curr_weight);
    6416     polylength = lengthpoly(Gomega);
    6417 */
    6418 //#ifdef CHECK_IDEAL_MWALK
     6716#ifdef CHECK_IDEAL_MWALK
    64196717    if(printout > 1)
    64206718    {
    64216719      idString(Gomega,"//** Mprwalk: Gomega");
    64226720    }
    6423 //#endif
     6721#endif
    64246722
    64256723    if(reduction == 0 && nstep > 1)
     
    64366734          idDelete(&FF);
    64376735          next_weight = MwalkNextWeightCC(curr_weight,target_weight,G);
     6736#ifdef PRINT_VECTORS
    64386737          if(printout > 0)
    64396738          {
    64406739            MivString(curr_weight, target_weight, next_weight);
    64416740          }
     6741#endif
    64426742        }
    64436743        else
     
    64506750        }
    64516751        Gomega = MwalkInitialForm(G, curr_weight);
     6752#ifdef CHECK_IDEAL_MWALK
    64526753        if(printout > 1)
    64536754        {
    64546755          idString(Gomega,"//** Mprwalk: Gomega");
    64556756        }
     6757#endif
    64566758      }
    64576759*/
     
    64706772    {
    64716773      Print("\n// ring r%d = %s;\n", nstep, rString(currRing));
     6774/*
    64726775      idElements(G, "G");
    64736776      headidString(G, "G");
     6777*/
    64746778    }
    64756779#endif
     
    65026806    {
    65036807      Print("\n// ring r%d = %s;\n", nstep, rString(currRing));
     6808/*
    65046809      idElements(Gomega1, "Gw");
    65056810      headidString(Gomega1, "headGw");
     6811*/
    65066812      PrintS("\n// compute a rGB of Gw:\n");
    65076813
     
    65116817    }
    65126818#endif
    6513 
     6819#ifdef TIME_TEST
    65146820    tim = clock();
    65156821    to = clock();
     6822#endif
    65166823    // compute a reduced Groebner basis of <Gomega> w.r.t. "newRing"
    65176824#ifdef  BUCHBERGER_ALG
     
    65216828    delete hilb_func;
    65226829#endif
    6523 //#ifdef CHECK_IDEAL_MWALK
     6830#ifdef CHECK_IDEAL_MWALK
    65246831    if(printout > 2)
    65256832    {
    65266833      idString(M,"//** Mprwalk: M");
    65276834    }
    6528 //#endif
    6529 
     6835#endif
     6836#ifdef TIME_TEST
    65306837    if(endwalks == TRUE)
    65316838    {
     
    65386845    else
    65396846      tstd=tstd+clock()-to;
    6540 
     6847#endif
    65416848    /* change the ring to oldRing */
    65426849    rChangeCurrRing(oldRing);
    65436850    M1 =  idrMoveR(M, newRing,currRing);
    65446851    Gomega2 =  idrMoveR(Gomega1, newRing,currRing);
    6545 
     6852#ifdef TIME_TEST
    65466853    to=clock();
     6854#endif
    65476855    /* compute a representation of the generators of submod (M)
    65486856       with respect to those of mod (Gomega).
    65496857       Gomega is a reduced Groebner basis w.r.t. the current ring */
    65506858    F = MLifttwoIdeal(Gomega2, M1, G);
     6859#ifdef TIME_TEST
    65516860    if(endwalks == FALSE)
    65526861      tlift = tlift+clock()-to;
    65536862    else
    65546863      xtlift=clock()-to;
    6555 
    6556 //#ifdef CHECK_IDEAL_MWALK
     6864#endif
     6865#ifdef CHECK_IDEAL_MWALK
    65576866    if(printout > 2)
    65586867    {
    65596868      idString(F,"//** Mprwalk: F");
    65606869    }
    6561 //#endif
     6870#endif
    65626871
    65636872    idDelete(&M1);
     
    65786887        PrintS("\n //** Mprwalk: reduce the Groebner basis.\n");
    65796888      }
     6889#ifdef TIME_TEST
    65806890      to=clock();
     6891#endif
    65816892      G = kInterRedCC(F1, NULL);
     6893#ifdef TIME_TEST
    65826894      if(endwalks == FALSE)
    65836895        tred = tred+clock()-to;
    65846896      else
    65856897        xtred=clock()-to;
     6898#endif
    65866899      idDelete(&F1);
    65876900    }
     
    65906903      break;
    65916904
    6592 Print("\n Next weight");
    65936905    NEXT_VECTOR:
    65946906#ifdef TIME_TEST
     
    66036915    to = clock();
    66046916#endif
    6605     Gomega = MwalkInitialForm(G, next_weight); // compute an initial form ideal of <G> w.r.t. "curr_vector"
     6917    // compute an initial form ideal of <G> w.r.t. "next_vector"
     6918    Gomega = MwalkInitialForm(G, next_weight);
    66066919#ifdef TIME_TEST
    66076920    tif = tif + clock()-to; //time for computing initial form ideal
    66086921#endif
    66096922
    6610     //polylength = 1 if there is a polynomial in Gomega with at least 3 monomials and 0 otherwise
    6611     polylength = lengthpoly(Gomega);
    6612     if(polylength > 0)
    6613     {
    6614       Print("\n there is a polynomial in Gomega with at least 3 monomials");
    6615       //low-dimensional facet of the cone
     6923    //lengthpoly(Gomega) = 1 if there is a polynomial in Gomega with at least 3 monomials and 0 otherwise
     6924    if(lengthpoly(Gomega) > 0)
     6925    {
     6926      Print("\n there is a polynomial in Gomega with at least 3 monomials,\n");
     6927      // low-dimensional facet of the cone
    66166928      delete next_weight;
     6929      if(target_M->length() == nV)
     6930      {
     6931        iv_M = MivMatrixOrder(curr_weight);
     6932      }
     6933      else
     6934      {
     6935        iv_M = MivMatrixOrderRefine(curr_weight,target_M);
     6936      }
    66176937#ifdef TIME_TEST
    66186938      to = clock();
    66196939#endif
    6620       next_weight = MWalkRandomNextWeight(G, curr_weight, target_weight, weight_rad, op_deg);
     6940      next_weight = MWalkRandomNextWeight(G, iv_M, target_weight, weight_rad, op_deg);
    66216941#ifdef TIME_TEST
    66226942      tnw = tnw + clock() - to;
     
    66306950      tif = tif + clock()-to; //time for computing initial form ideal
    66316951#endif
     6952  Print("delete\n");
     6953      delete iv_M;
    66326954    }
    66336955
     
    66466968    tnw=tnw+clock()-to;
    66476969*/
    6648 //#ifdef PRINT_VECTORS
     6970#ifdef PRINT_VECTORS
    66496971    if(printout > 0)
    66506972    {
    66516973      MivString(curr_weight, target_weight, next_weight);
    66526974    }
    6653 //#endif
     6975#endif
    66546976
    66556977    if(Overflow_Error == TRUE)
     
    66746996
    66756997    delete next_weight;
    6676   }//while
     6998  }// end of while-loop
    66776999
    66787000  if(tp_deg != 1)
     
    66987020    TargetRing=currRing;
    66997021    F1 = idrMoveR(G, newRing,currRing);
    6700 #ifdef CHECK_IDEAL
    6701       headidString(G, "G");
    6702 #endif
    67037022
    67047023    // check whether the pertubed target vector stays in the correct cone
     
    67117030      if(ntestw != 1 && printout > 2)
    67127031      {
     7032#ifdef PRINT_VECTORS
    67137033        ivString(pert_target_vector, "tau");
     7034#endif
    67147035        PrintS("\n// **Mprwalk: perturbed target vector doesn't stay in cone.");
    67157036        Print("\n// ring r%d = %s;\n", nstep, rString(currRing));
     
    67177038      }
    67187039      // LastGB is "better" than the kStd subroutine
     7040#ifdef TIME_TEST
    67197041      to=clock();
     7042#endif
    67207043      ideal eF1;
    67217044      if(nP == 0 || tp_deg == 1 || MivSame(orig_target, exivlp) != 1 || target_M->length() != nV)
     
    67397062        F2=NULL;
    67407063      }
     7064#ifdef TIME_TEST
    67417065      xtextra=clock()-to;
     7066#endif
    67427067      ring exTargetRing = currRing;
    67437068
     
    68177142  intvec* hilb_func;
    68187143#endif
    6819 //  intvec* extXtau;
     7144  //intvec* extXtau;
    68207145  intvec* next_vect;
    68217146  intvec* omega2 = new intvec(nV);
     
    68667191    nwalks ++;
    68677192    NEXT_VECTOR_FRACTAL:
     7193#ifdef TIME_TEST
    68687194    to=clock();
     7195#endif
    68697196    // determine the next border
    68707197    next_vect = MkInterRedNextWeight(omega,omega2,G);
     7198#ifdef TIME_TEST
    68717199    xtnw=xtnw+clock()-to;
    6872 
     7200#endif
    68737201    oRing = currRing;
    68747202
     
    69287256
    69297257        delete next_vect;
     7258#ifdef TIME_TEST
    69307259        to=clock();
    6931 
     7260#endif
    69327261        // to avoid the value of Overflow_Error that occur in Mfpertvector
    69337262        Overflow_Error = FALSE;
    69347263        next_vect = MkInterRedNextWeight(omega,omega2,G);
     7264#ifdef TIME_TEST
    69357265        xtnw=xtnw+clock()-to;
     7266#endif
    69367267      }// end of (if MivComp(next_vect, omega2) == 1)
    69377268
    6938 //#ifdef PRINT_VECTORS
     7269#ifdef PRINT_VECTORS
    69397270      if(printout > 0)
    69407271      {
    69417272        MivString(omega, omega2, next_vect);
    69427273      }
    6943 //#endif
     7274#endif
    69447275
    69457276    // check whether the the computed vector is in the correct cone.
     
    69697300              rString(currRing));
    69707301      }
     7302#ifdef TIME_TEST
    69717303      to=clock();
     7304#endif
    69727305      Gt = idrMoveR(G, oRing,currRing);
    69737306      G1 = MstdCC(Gt);
     7307#ifdef TIME_TEST
    69747308      xtextra=xtextra+clock()-to;
     7309#endif
    69757310      Gt = NULL;
    69767311
     
    69887323      return (G1);
    69897324    }
    6990 
    69917325
    69927326    /* If the perturbed target vector stays in the correct cone,
     
    70767410                rString(currRing));
    70777411        }
     7412#ifdef TIME_TEST
    70787413        to=clock();
     7414#endif
    70797415        G = MstdCC(Gt);
     7416#ifdef TIME_TEST
    70807417        xtextra=xtextra+clock()-to;
    7081 
     7418#endif
    70827419        oRing = currRing;
    70837420
     
    71407477    }
    71417478    delete next_vect;
    7142 
     7479#ifdef TIME_TEST
    71437480    to=clock();
     7481#endif
    71447482    // Take the initial form of <G> w.r.t. omega
    71457483    Gomega = MwalkInitialForm(G, omega);
     7484#ifdef TIME_TEST
    71467485    xtif=xtif+clock()-to;
     7486#endif
     7487#ifdef CHECK_IDEAL_MWALK
    71477488    if(printout > 1)
    71487489    {
    71497490      idString(Gomega,"//** rec_fractal_call: Gomega");
    71507491    }
    7151 
     7492#endif
    71527493    if(reduction == 0)
    71537494    {
     
    71947535        Print("\n//** rec_fractal_call: Maximal recursion depth.\n");
    71957536      }
    7196 
     7537#ifdef TIME_TEST
    71977538      to=clock();
     7539#endif
    71987540#ifdef  BUCHBERGER_ALG
    71997541      Gresult = MstdhomCC(Gomega1);
     
    72027544      delete hilb_func;
    72037545#endif
     7546#ifdef TIME_TEST
    72047547      xtstd=xtstd+clock()-to;
     7548#endif
    72057549    }
    72067550    else
     
    72107554      Gresult = rec_fractal_call(idCopy(Gomega1),nlev+1,omega,reduction,printout);
    72117555    }
     7556#ifdef CHECK_IDEAL_MWALK
    72127557    if(printout > 2)
    72137558    {
    72147559      idString(Gresult,"//** rec_fractal_call: M");
    72157560    }
     7561#endif
    72167562    //convert a Groebner basis from a ring to another ring
    72177563    new_ring = currRing;
     
    72207566    Gresult1 = idrMoveR(Gresult, new_ring,currRing);
    72217567    Gomega2 = idrMoveR(Gomega1, new_ring,currRing);
    7222 
     7568#ifdef TIME_TEST
    72237569    to=clock();
     7570#endif
    72247571    // Lifting process
    72257572    F = MLifttwoIdeal(Gomega2, Gresult1, G);
     7573#ifdef TIME_TEST
    72267574    xtlift=xtlift+clock()-to;
     7575#endif
     7576#ifdef CHECK_IDEAL_MWALK
    72277577    if(printout > 2)
    72287578    {
    72297579      idString(F,"//** rec_fractal_call: F");
    72307580    }
     7581#endif
    72317582    idDelete(&Gresult1);
    72327583    idDelete(&Gomega2);
     
    72347585
    72357586    rChangeCurrRing(new_ring);
    7236     //F1 = idrMoveR(F, oRing,currRing);
    72377587    G = idrMoveR(F,oRing,currRing);
     7588/*
     7589    F1 = idrMoveR(F, oRing,currRing);
     7590#ifdef TIME_TEST
    72387591    to=clock();
     7592#endif
    72397593    // Interreduce G
    7240     // G = kInterRedCC(F1, NULL);
     7594    G = kInterRedCC(F1, NULL);
     7595#ifdef TIME_TEST
    72417596    xtred=xtred+clock()-to;
    7242     //idDelete(&F1);
     7597#endif
     7598    idDelete(&F1);
     7599*/
    72437600  }
    72447601}
     
    72537610  //Print("\n\n// Entering the %d-th recursion:", nlev);
    72547611
    7255   int i, polylength, nV = currRing->N;
     7612  int nwalks = 0,i,nV=currRing->N;//polylength
    72567613  ring new_ring, testring;
    72577614  //ring extoRing;
    72587615  ideal Gomega, Gomega1, Gomega2, F, FF, F1, Gresult, Gresult1, G1, Gt;
    7259   int nwalks = 0;
    72607616  intvec* Mwlp;
    72617617#ifndef BUCHBERGER_ALG
     
    72647620//  intvec* extXtau;
    72657621  intvec* next_vect;
     7622  intvec* iv_M;
    72667623  intvec* omega2 = new intvec(nV);
    72677624  intvec* omtmp = new intvec(nV);
     
    73137670    nwalks ++;
    73147671    NEXT_VECTOR_FRACTAL:
     7672#ifdef TIME_TEST
    73157673    to=clock();
     7674#endif
    73167675    /* determine the next border */
    73177676    next_vect = MkInterRedNextWeight(omega,omega2,G);
    7318     if(polylength > 0 && G->m[0] != NULL)
     7677#ifdef TIME_TEST
     7678    xtnw=xtnw+clock()-to;
     7679#endif
     7680    if(lengthpoly(MwalkInitialForm(G, next_vect)) > 0 && G->m[0] != NULL)
    73197681    {
    73207682      if(printout > 0)
    73217683      {
    7322         PrintS("\n**// rec_r_fractal_call: there is a polynomial in Gomega with at least 3 monomials, low-dimensional facet of the cone.\n");
     7684        PrintS("\n**// rec_r_fractal_call: there is a polynomial in Gomega with at least 3 monomials.\n");
    73237685      }
    73247686      delete next_vect;
    7325       next_vect = MWalkRandomNextWeight(G,omega,omega2,weight_rad,nlev);
     7687      iv_M = MivMatrixOrder(omega);
     7688#ifdef TIME_TEST
     7689      to=clock();
     7690#endif
     7691      next_vect = MWalkRandomNextWeight(G,iv_M,omega2,weight_rad,nlev);
     7692#ifdef TIME_TEST
     7693      xtnw=xtnw+clock()-to;
     7694#endif
    73267695      if(isNegNolVector(next_vect) == 1)
    73277696      {
    73287697        delete next_vect;
     7698#ifdef TIME_TEST
     7699        to=clock();
     7700#endif
    73297701        next_vect = MkInterRedNextWeight(omega,omega2,G);
    7330       }
    7331     }
    7332     xtnw=xtnw+clock()-to;
    7333 
     7702#ifdef TIME_TEST
     7703        xtnw=xtnw+clock()-to;
     7704#endif
     7705      }
     7706    }
    73347707    oRing = currRing;
    73357708
     
    73867759        delete Mwlp;
    73877760
    7388         for(i=nV-1; i>=0; i--) {
     7761        for(i=nV-1; i>=0; i--)
     7762        {
    73897763          (*omega2)[i] = (*Xtau)[nV+i];
    73907764          (*omega)[i] = (*Xsigma)[nV+i];
     
    73927766
    73937767        delete next_vect;
     7768
     7769   //to avoid the value of Overflow_Error that occur in Mfpertvector
     7770        Overflow_Error = FALSE;
     7771#ifdef TIME_TEST
    73947772        to=clock();
    7395 
    7396         /*
    7397         to avoid the value of Overflow_Error that occur in
    7398         Mfpertvector
    7399         */
    7400         Overflow_Error = FALSE;
    7401 
     7773#endif
    74027774        next_vect = MkInterRedNextWeight(omega,omega2,G);
    7403         if(G->m[0] != NULL && polylength > 0)
     7775#ifdef TIME_TEST
     7776        xtnw=xtnw+clock()-to;
     7777#endif
     7778        if(lengthpoly(MwalkInitialForm(G, next_vect)) > 0 && G->m[0] != NULL)
    74047779        {
    7405          
    7406           PrintS("//** rec_r_fractal_call: there is a polynomial in Gomega with at least 3 monomials, low-dimensional facet of the cone");
    7407          
     7780          // there is a polynomial in Gomega with at least 3 monomials
     7781          iv_M = MivMatrixOrder(omega);
    74087782          delete next_vect;
    7409           next_vect = MWalkRandomNextWeight(G,omega,omega2,weight_rad,nlev);
     7783#ifdef TIME_TEST
     7784          to=clock();
     7785#endif
     7786          next_vect = MWalkRandomNextWeight(G,iv_M,omega2,weight_rad,nlev);
     7787#ifdef TIME_TEST
     7788          xtnw=xtnw+clock()-to;
     7789#endif
     7790          delete iv_M;
    74107791          if(isNegNolVector(next_vect) == 1)
    74117792          {
    74127793            delete next_vect;
     7794#ifdef TIME_TEST
     7795            to=clock();
     7796#endif
    74137797            next_vect = MkInterRedNextWeight(omega,omega2,G);
     7798#ifdef TIME_TEST
     7799        xtnw=xtnw+clock()-to;
     7800#endif
    74147801          }
    74157802        }
    7416         xtnw=xtnw+clock()-to;
    7417       }
    7418 //#ifdef PRINT_VECTORS
     7803      }
     7804#ifdef PRINT_VECTORS
    74197805      if(printout > 0)
    74207806      {
    74217807        MivString(omega, omega2, next_vect);
    74227808      }
    7423 //#endif
    7424 
    7425     /* check whether the the computed vector is in the correct cone
     7809#endif
     7810
     7811/*    check whether the the computed vector is in the correct cone
    74267812       If no, the reduced GB of an omega-homogeneous ideal will be
    7427        computed by Buchberger algorithm and stop this recursion step*/
    7428     //if(test_w_in_ConeCC(G, next_vect) != 1) //e.g. Example s7, cyc6
    7429     if(Overflow_Error == TRUE || test_w_in_ConeCC(G,next_vect) != 1)
     7813       computed by Buchberger algorithm and stop this recursion step
     7814*/
     7815    if(Overflow_Error == TRUE || test_w_in_ConeCC(G,next_vect) != 1)//e.g. Example s7, cyc6
    74307816    {
    74317817      delete next_vect;
     
    74557841              rString(currRing));
    74567842      }
     7843      Gt = idrMoveR(G, oRing,currRing);
     7844#ifdef TIME_TEST
    74577845      to=clock();
    7458       Gt = idrMoveR(G, oRing,currRing);
     7846#endif
    74597847      G1 = MstdCC(Gt);
     7848#ifdef TIME_TEST
    74607849      xtextra=xtextra+clock()-to;
     7850#endif
    74617851      Gt = NULL;
    74627852
     
    75197909
    75207910#ifndef  MSTDCC_FRACTAL
    7521         //ivString(Xtau, "old Xtau");
     7911#ifdef PRINT_VECTORS
     7912        if(printout > 0)
     7913        {
     7914          ivString(Xtau, "old Xtau");
     7915        }
     7916#endif
    75227917        intvec* Xtautmp;
    75237918        if(ivtarget->length() == nV)
     
    75437938        Xtau = Xtautmp;
    75447939        Xtautmp = NULL;
    7545         //ivString(Xtau, "new  Xtau");
     7940#ifdef PRINT_VECTORS
     7941        if(printout > 0)
     7942        {
     7943          ivString(Xtau, "new  Xtau");
     7944        }
     7945#endif
    75467946
    75477947        for(i=nV-1; i>=0; i--)
     
    75617961                rString(currRing));
    75627962        }
     7963#ifdef TIME_TEST
    75637964        to=clock();
     7965#endif
    75647966        G = MstdCC(Gt);
     7967#ifdef TIME_TEST
    75657968        xtextra=xtextra+clock()-to;
    7566 
     7969#endif
    75677970        oRing = currRing;
    75687971
     
    76268029    }
    76278030    delete next_vect;
    7628 
     8031#ifdef TIME_TEST
    76298032    to=clock();
     8033#endif
    76308034    // Take the initial form of <G> w.r.t. omega
    76318035    Gomega = MwalkInitialForm(G, omega);
     8036#ifdef TIME_TEST
    76328037    xtif=xtif+clock()-to;
     8038#endif
    76338039    //polylength = 1 if there is a polynomial in Gomega with at least 3 monomials and 0 otherwise
    7634     polylength = lengthpoly(Gomega);
     8040    //polylength = lengthpoly(Gomega);
     8041#ifdef CHECK_IDEAL_MWALK
    76358042    if(printout > 1)
    76368043    {
    76378044      idString(Gomega,"//** rec_r_fractal_call: Gomega");
    76388045    }
    7639 
     8046#endif
    76408047    if(reduction == 0)
    76418048    {
     
    76778084    if(nlev == Xnlev || lengthpoly(Gomega1) == 0)
    76788085    {
     8086#ifdef TIME_TEST
    76798087      to=clock();
     8088#endif
    76808089#ifdef  BUCHBERGER_ALG
    76818090      Gresult = MstdhomCC(Gomega1);
     
    76848093      delete hilb_func;
    76858094#endif
     8095#ifdef TIME_TEST
    76868096      xtstd=xtstd+clock()-to;
     8097#endif
    76878098    }
    76888099    else
     
    76928103      Gresult = rec_r_fractal_call(idCopy(Gomega1),nlev+1,omega,weight_rad,reduction,printout);
    76938104    }
     8105#ifdef CHECK_IDEAL_MWALK
    76948106    if(printout > 2)
    76958107    {
    76968108      idString(Gresult,"//** rec_r_fractal_call: M");
    76978109    }
     8110#endif
    76988111    //convert a Groebner basis from a ring to another ring
    76998112    new_ring = currRing;
     
    77028115    Gresult1 = idrMoveR(Gresult, new_ring,currRing);
    77038116    Gomega2 = idrMoveR(Gomega1, new_ring,currRing);
    7704 
     8117#ifdef TIME_TEST
    77058118    to=clock();
     8119#endif
    77068120    // Lifting process
    77078121    F = MLifttwoIdeal(Gomega2, Gresult1, G);
     8122#ifdef TIME_TEST
    77088123    xtlift=xtlift+clock()-to;
    7709 
     8124#endif
     8125#ifdef CHECK_IDEAL_MWALK
    77108126    if(printout > 2)
    77118127    {
    77128128      idString(F,"//** rec_r_fractal_call: F");
    77138129    }
    7714 
     8130#endif
    77158131    idDelete(&Gresult1);
    77168132    idDelete(&Gomega2);
     
    77208136    //F1 = idrMoveR(F, oRing,currRing);
    77218137    G = idrMoveR(F,oRing,currRing);
     8138/*
     8139#ifdef TIME_TEST
    77228140    to=clock();
     8141#endif
    77238142    // Interreduce G
    7724     //G = kInterRedCC(F1, NULL);
     8143    G = kInterRedCC(F1, NULL);
     8144#ifdef TIME_TEST
    77258145    xtred=xtred+clock()-to;
    7726     //idDelete(&F1);
     8146#endif
     8147    idDelete(&F1);
     8148*/
    77278149  }
    77288150}
     
    77558177  Xngleich = 0;
    77568178  Xcall = 0;
     8179#ifdef TIME_TEST
    77578180  xtif=0; xtstd=0; xtlift=0; xtred=0; xtnw=0; xtextra=0;
    77588181  xftinput = clock();
    7759 
     8182#endif
    77608183  ring  oldRing = currRing;
    77618184  int i, nV = currRing->N;
     
    77638186  Xivinput = ivtarget;
    77648187  ngleich = 0;
     8188#ifdef TIME_TEST
    77658189  to=clock();
     8190#endif
    77668191  ideal I = MstdCC(G);
    77678192  G = NULL;
     8193#ifdef TIME_TEST
    77688194  xftostd=clock()-to;
     8195#endif
    77698196  Xsigma = ivstart;
    77708197
     
    78618288
    78628289  I = idrMoveR(I1,tRing,currRing);
     8290#ifdef TIME_TEST
    78638291  to=clock();
     8292#endif
    78648293  ideal J = MstdCC(I);
    78658294  idDelete(&I);
     8295#ifdef TIME_TEST
    78668296  xftostd=xftostd+clock()-to;
    7867 
     8297#endif
    78688298  ideal resF;
    78698299  ring helpRing = currRing;
     
    79078337{
    79088338  BITSET save1 = si_opt_1; // save current options
     8339  //check that weight radius is valid
     8340  if(weight_rad < 0)
     8341  {
     8342    Werror("Invalid radius.\n");
     8343    return NULL;
     8344  }
    79098345  if(reduction == 0)
    79108346  {
    79118347    si_opt_1 &= (~Sy_bit(OPT_REDSB)); // no reduced Groebner basis
    7912     //si_opt_1 &= (~Sy_bit(OPT_REDTAIL)); // not tail reductions
     8348    si_opt_1 &= (~Sy_bit(OPT_REDTAIL)); // not tail reductions
    79138349  }
    79148350  Set_Error(FALSE);
     
    79208356  Xngleich = 0;
    79218357  Xcall = 0;
     8358#ifdef TIME_TEST
    79228359  xtif=0; xtstd=0; xtlift=0; xtred=0; xtnw=0; xtextra=0;
    79238360  xftinput = clock();
    7924 
     8361#endif
    79258362  ring  oldRing = currRing;
    79268363  int i, nV = currRing->N;
     
    79288365  Xivinput = ivtarget;
    79298366  ngleich = 0;
     8367#ifdef TIME_TEST
    79308368  to=clock();
     8369#endif
    79318370  ideal I = MstdCC(G);
    79328371  G = NULL;
     8372#ifdef TIME_TEST
    79338373  xftostd=clock()-to;
     8374#endif
    79348375  Xsigma = ivstart;
    79358376
     
    80268467
    80278468  I = idrMoveR(I1,tRing,currRing);
     8469#ifdef TIME_TEST
    80288470  to=clock();
     8471#endif
    80298472  ideal J = MstdCC(I);
    80308473  idDelete(&I);
     8474#ifdef TIME_TEST
    80318475  xftostd=xftostd+clock()-to;
    8032 
     8476#endif
    80338477  ideal resF;
    80348478  ring helpRing = currRing;
     
    82518695      goto  BE_FINISH;
    82528696#endif
    8253 
     8697/*
    82548698#ifdef CHECK_IDEAL_MWALK
    82558699      idElements(G, "G");
    82568700      //headidString(G, "G");
    82578701#endif
    8258 
     8702*/
    82598703      if(MivSame(target_tmp, iv_lp) == 1)
    82608704        if (rParameter(currRing) != NULL)
     
    89909434  clock_t tinput=clock();
    89919435  int i, nV = currRing->N;
     9436
     9437  //check that perturbation degree is valid
     9438  if(tp_deg < 1 || tp_deg > nV)
     9439  {
     9440    Werror("Invalid perturbation degree.\n");
     9441    return NULL;
     9442  }
     9443
    89929444  int nwalk=0, endwalks=0, ntestwinC=1;
    89939445  int tp_deg_tmp = tp_deg;
     
    90999551    {
    91009552      Print("\n//  it is  %d-th step!!", nwalk);
    9101       idElements(Gomega1, "Gw");
     9553      idString(Gomega1, "Gw");
    91029554      PrintS("\n//  compute a rGB of Gw:");
    91039555    }
Note: See TracChangeset for help on using the changeset viewer.