Changeset 20c1a3 in git for Singular/walk.cc


Ignore:
Timestamp:
Jun 3, 2015, 5:16:50 PM (9 years ago)
Author:
Stephan Oberfranz <oberfran@…>
Branches:
(u'spielwiese', '6e5adcba05493683b94648c659a729c189812c77')
Children:
269fc4d27717409c0438dbf31005bd82968fc09d633196803bc6aac699d9f7df3d38f77b6222ea7d
Parents:
2c2d2e89ad6e29cc2a1d2903129b937241868370
Message:
fractal walk fixed
File:
1 edited

Legend:

Unmodified
Added
Removed

  • Singular/walk.cc

    r2c2d2e8 r20c1a3  
    806806  for(i=nG-1; i>=0; i--)
    807807  {
    808     mi = MpolyInitialForm(G->m[i], iv);
    809     gi = G->m[i];
    810 
     808    mi = pHead(MpolyInitialForm(G->m[i], iv));
     809    //Print("\n **// test_w_in_ConeCC: lm(initial)= %s \n",pString(mi));
     810    gi = pHead(G->m[i]);
     811    //Print("\n **// test_w_in_ConeCC: lm(ideal)= %s \n",pString(gi));
    811812    if(mi == NULL)
    812813    {
    … …  
    16461647    (* result)[i] = mpz_get_si(pert_vector[i]);
    16471648  }
    1648 
     1649/*
    16491650  j = 0;
    1650   for(i=0; i<nV; i++)
     1651  for(i=0; i<niv; i++)
    16511652  {
    16521653    (* result1)[i] = mpz_get_si(pert_vector[i]);
    … …  
    16581659    }
    16591660  }
    1660   if(j > nV - 1)
     1661  if(j > niv - 1)
    16611662  {
    16621663    // Print("\n//  MfPertwalk: geaenderter vector gleich Null! \n");
    … …  
    16641665    goto CHECK_OVERFLOW;
    16651666  }
    1666 
     1667/*
    16671668// check that the perturbed weight vector lies in the Groebner cone
     1669  Print("\n========================================\n//** MfPertvector: test in Cone.\n");
    16681670  if(test_w_in_ConeCC(G,result1) != 0)
    16691671  {
    … …  
    16811683    // Print("\n// Mfpertwalk: geaenderter vector liegt nicht in Groebnerkegel! \n");
    16821684  }
    1683 
     1685  Print("\n ========================================\n");*/
    16841686  CHECK_OVERFLOW:
    16851687
    … …  
    50435045  Set_Error(FALSE);
    50445046  Overflow_Error = FALSE;
    5045   BOOLEAN endwalks = FALSE;
     5047  //BOOLEAN endwalks = FALSE;
    50465048#ifdef TIME_TEST
    50475049  clock_t tinput, tostd, tif=0, tstd=0, tlift=0, tred=0, tnw=0;
    … …  
    50575059  ring newRing;
    50585060  ring XXRing = baseRing;
     5061  ring targetRing;
    50595062  intvec* ivNull = new intvec(nV);
    50605063  intvec* curr_weight = new intvec(nV);
    … …  
    50655068  for(i=0; i<nV; i++)
    50665069  {
    5067     (*tmp_weight)[i] = (*target_M)[i];
     5070    (*tmp_weight)[i] = (*orig_M)[i];
    50685071  }
    50695072*/
    … …  
    50905093  }
    50915094//#endif
     5095
     5096  if(target_M->length() == nV)
     5097  {
     5098   // define the target ring
     5099    targetRing = VMrDefault(target_weight);
     5100  }
     5101  else
     5102  {
     5103    targetRing = VMatrDefault(target_M);
     5104  }
     5105  if(orig_M->length() == nV)
     5106  {
     5107    // define a new ring with ordering "(a(curr_weight),lp)
     5108    newRing = VMrDefault(curr_weight);
     5109  }
     5110  else
     5111  {
     5112    newRing = VMatrDefault(orig_M);
     5113  }
     5114  rChangeCurrRing(newRing);
     5115
    50925116#ifdef TIME_TEST
    50935117  to = clock();
    50945118#endif
    5095   if(orig_M->length() == nV)
    5096   {
    5097     // define a new ring with ordering "(a(curr_weight),lp)
    5098     newRing = VMrDefault(curr_weight);
    5099   }
    5100   else
    5101   {
    5102     newRing = VMatrDefault(orig_M);
    5103   }
    5104   rChangeCurrRing(newRing);
     5119
    51055120  ideal G = MstdCC(idrMoveR(Go,baseRing,currRing));
    5106   //ideal G = idrMoveR(Go,baseRing,currRing);
    5107   baseRing = currRing;
     5121
    51085122#ifdef TIME_TEST
    51095123  tostd = clock()-to;
    51105124#endif
    51115125
     5126  baseRing = currRing;
    51125127  nwalk = 0;
     5128
    51135129  while(1)
    51145130  {
    51155131    nwalk ++;
    51165132    nstep ++;
     5133
    51175134#ifdef TIME_TEST
    51185135    to = clock();
    51195136#endif
    5120 //#ifdef CHECK_IDEAL_MWALK
    5121     if(printout > 2)
    5122     {
    5123       idString(G,"//** Mwalk: G");
    5124     }
    5125 //#endif
    5126     // test whether target cone is reached
    5127     if(reduction !=0 && test_w_in_ConeCC(G,target_weight) == 1)
    5128     {
    5129       endwalks = TRUE;
    5130     }
    51315137    // compute an initial form ideal of <G> w.r.t. "curr_vector"
    51325138    Gomega = MwalkInitialForm(G, curr_weight);
    51335139#ifdef TIME_TEST
    5134     //time for computing initial form ideal
    51355140    tif = tif + clock()-to;
    51365141#endif
     5142
    51375143//#ifdef CHECK_IDEAL_MWALK
    51385144    if(printout > 1)
    … …  
    51415147    }
    51425148//#endif
     5149
    51435150    if(reduction == 0)
    51445151    {
    5145       //PrintS("\n//** Mwalk: test middle of cone!\n");
    51465152      FF = middleOfCone(G,Gomega);
    5147       //PrintS("\n//** Mwalk: Test F!\n");
    51485153      if( FF != NULL)
    51495154      {
    … …  
    51515156        G = idCopy(FF);
    51525157        idDelete(&FF);
    5153         //PrintS("\n//** Mwalk: FF nicht NULL! Compue next vector.\n");
    51545158        goto NEXT_VECTOR;
    51555159      }
    51565160    }
     5161
    51575162#ifndef  BUCHBERGER_ALG
    51585163    if(isNolVector(curr_weight) == 0)
    … …  
    51655170    }
    51665171#endif
     5172
    51675173    if(nwalk == 1)
    51685174    {
    … …  
    51815187     if(target_M->length() == nV)
    51825188     {
    5183        //define a new ring with ordering "(a(curr_weight),Wp(target_weight))"
    5184        newRing = VMrRefine(curr_weight,target_weight);
     5189       //define a new ring with ordering "(a(curr_weight),lp)"
     5190       newRing = VMrDefault(curr_weight);
    51855191     }
    51865192     else
    … …  
    52255231    // where Gomega is a reduced Groebner basis w.r.t. the current ring
    52265232    F = MLifttwoIdeal(Gomega2, M1, G);
     5233
    52275234#ifdef TIME_TEST
    52285235    tlift = tlift + clock() - to;
    … …  
    52365243    idDelete(&Gomega2);
    52375244    idDelete(&M1);
     5245
    52385246    rChangeCurrRing(newRing); // change the ring to newRing
    52395247    G = idrMoveR(F,baseRing,currRing);
    52405248    idDelete(&F);
     5249    idSkipZeroes(G);
     5250
     5251//#ifdef CHECK_IDEAL_MWALK
     5252    if(printout > 2)
     5253    {
     5254      idString(G, "//** Mwalk: G");
     5255    }
     5256//#endif
     5257
     5258    rChangeCurrRing(targetRing);
     5259    G = idrMoveR(G,newRing,currRing);
     5260    // test whether target cone is reached
     5261    if(reduction !=0 && test_w_in_ConeCC(G,curr_weight) == 1)
     5262    {
     5263      baseRing = currRing;
     5264      break;
     5265      //endwalks = TRUE;
     5266    }
     5267
     5268    rChangeCurrRing(newRing);
     5269    G = idrMoveR(G,targetRing,currRing);
    52415270    baseRing = currRing;
     5271
     5272/*
    52425273#ifdef TIME_TEST
    52435274    to = clock();
    … …  
    52475278    tstd = tstd + clock() - to;
    52485279#endif
    5249     idSkipZeroes(G);
    5250 //#ifdef CHECK_IDEAL_MWALK
    5251     if(printout > 2)
    5252     {
    5253       idString(G, "//** Mwalk: G");
    5254     }
    5255 //#endif
     5280*/
     5281
     5282
    52565283#ifdef TIME_TEST
    52575284    to = clock();
    … …  
    52685295    }
    52695296//#endif
    5270     if(MivComp(target_weight,curr_weight) == 1 || endwalks == TRUE)
    5271     {
    5272 
     5297    if(MivComp(target_weight,curr_weight) == 1)// || endwalks == TRUE)
     5298    {/*
    52735299//#ifdef CHECK_IDEAL_MWALK
    52745300      if(printout > 0)
    … …  
    53285354      to = clock();
    53295355#endif
    5330       /*PrintS("\n //**Mwalk: Interreduce");
     5356      //PrintS("\n //**Mwalk: Interreduce");
    53315357      //interreduce the Groebner basis <G> w.r.t. currRing
    5332       G = kInterRedCC(G,NULL);*/
     5358      //G = kInterRedCC(G,NULL);
    53335359#ifdef TIME_TEST
    53345360      tred = tred + clock() - to;
    53355361#endif
    53365362      idSkipZeroes(G);
    5337       delete next_weight;
     5363      delete next_weight; */
    53385364      break;
    53395365    }
    … …  
    53485374  rChangeCurrRing(XXRing);
    53495375  ideal result = idrMoveR(G,baseRing,currRing);
     5376  idDelete(&Go);
    53505377  idDelete(&G);
    53515378  //delete tmp_weight;
    … …  
    53765403  Set_Error(FALSE);
    53775404  Overflow_Error = FALSE;
    5378   BOOLEAN endwalks = FALSE;
     5405  //BOOLEAN endwalks = FALSE;
    53795406#ifdef TIME_TEST
    53805407  clock_t tinput, tostd, tif=0, tstd=0, tlift=0, tred=0, tnw=0;
    … …  
    53895416  ideal Gomega, M, F,FF, Gomega1, Gomega2, M1;
    53905417  ring newRing;
     5418  ring targetRing;
    53915419  ring baseRing = currRing;
    53925420  ring XXRing = currRing;
    … …  
    54225450  to = clock();
    54235451#endif
     5452
     5453  if(target_M->length() == nV)
     5454  {
     5455   // define the target ring
     5456    targetRing = VMrDefault(target_weight);
     5457  }
     5458  else
     5459  {
     5460    targetRing = VMatrDefault(target_M);
     5461  }
    54245462  if(orig_M->length() == nV)
    54255463  {
    … …  
    54455483    to = clock();
    54465484#endif
    5447 //#ifdef CHECK_IDEAL_MWALK
    5448     if(printout > 2)
    5449     {
    5450       idString(G,"//** Mrwalk: G");
    5451     }
    5452 //#endif
     5485
    54535486    Gomega = MwalkInitialForm(G, curr_weight); // compute an initial form ideal of <G> w.r.t. "curr_vector"
    54545487    //polylength = 1 if there is a polynomial in Gomega with at least 3 monomials and 0 otherwise
    … …  
    54645497//#endif
    54655498    // test whether target cone is reached
    5466     if(test_w_in_ConeCC(G,target_weight) == 1)
     5499/*    if(test_w_in_ConeCC(G,target_weight) == 1)
    54675500    {
    54685501      endwalks = TRUE;
    5469     }
     5502    }*/
    54705503    if(reduction == 0)
    54715504    {
    5472       //PrintS("\n//** Mrwalk: test middle of cone!\n");
     5505     
    54735506      FF = middleOfCone(G,Gomega);
    5474       //PrintS("\n//** Mrwalk: Test F!\n");
     5507     
    54755508      if(FF != NULL)
    54765509      {
    … …  
    54785511        G = idCopy(FF);
    54795512        idDelete(&FF);
    5480         //PrintS("\n//** Mrwalk: FF nicht NULL! Compue next vector.\n");
     5513       
    54815514        goto NEXT_VECTOR;
    54825515      }
    … …  
    55075540     if(target_M->length() == nV)
    55085541     {
    5509        newRing = VMrRefine(curr_weight,target_weight); //define a new ring with ordering "(a(curr_weight),Wp(target_weight))"
     5542       newRing = VMrDefault(curr_weight); // define a new ring with ordering "(a(curr_weight),lp)
     5543       //newRing = VMrRefine(curr_weight,target_weight); //define a new ring with ordering "(a(curr_weight),Wp(target_weight))"
    55105544     }
    55115545     else
    … …  
    55755609    }
    55765610//#endif
     5611
     5612    rChangeCurrRing(targetRing);
     5613    G = idrMoveR(G,newRing,currRing);
     5614    // test whether target cone is reached
     5615    if(reduction !=0 && test_w_in_ConeCC(G,curr_weight) == 1)
     5616    {
     5617      baseRing = currRing;
     5618      break;
     5619    }
     5620
     5621    rChangeCurrRing(newRing);
     5622    G = idrMoveR(G,targetRing,currRing);
     5623    baseRing = currRing;
     5624
     5625
    55775626    NEXT_VECTOR:
    55785627#ifdef TIME_TEST
    … …  
    55965645    }
    55975646//#endif
    5598     if(MivComp(next_weight, ivNull) == 1 || MivComp(target_weight,curr_weight) == 1 || endwalks == TRUE)
    5599     {
     5647    if(MivComp(next_weight, ivNull) == 1 || MivComp(target_weight,curr_weight) == 1)// || endwalks == TRUE)
     5648    {/*
    56005649//#ifdef CHECK_IDEAL_MWALK
    56015650      if(printout > 0)
    … …  
    56555704      to = clock();
    56565705#endif
    5657       /*PrintS("\n //**Mrwalk: Interreduce");
     5706      //PrintS("\n //**Mrwalk: Interreduce");
    56585707      //interreduce the Groebner basis <G> w.r.t. currRing
    5659       G = kInterRedCC(G,NULL);*/
     5708      //G = kInterRedCC(G,NULL);
    56605709#ifdef TIME_TEST
    56615710      tred = tred + clock() - to;
    56625711#endif
    56635712      idSkipZeroes(G);
    5664       delete next_weight;
     5713      delete next_weight;*/
    56655714      break;
    56665715    }
    … …  
    67246773{
    67256774  Overflow_Error =  FALSE;
    6726   //Print("\n\n// Entering the %d-th recursion:", nlev);
    6727 
     6775  if(printout >0)
     6776  {
     6777    Print("\n\n// Entering the %d-th recursion:", nlev);
     6778  }
    67286779  int i, nV = currRing->N;
    67296780  ring new_ring, testring;
    … …  
    67936844    // We only perturb the current target vector at the recursion level 1
    67946845    if(Xngleich == 0 && nlev == 1) //(ngleich == 0) important, e.g. ex2, ex3
    6795       if (MivComp(next_vect, omega2) != 1)
     6846      if (MivComp(next_vect, omega2) == 1)
    67966847      {
    67976848        // to dispense with taking initial (and lifting/interreducing
    … …  
    68486899        to=clock();
    68496900
    6850         /* to avoid the value of Overflow_Error that occur in Mfpertvector*/
     6901        // to avoid the value of Overflow_Error that occur in Mfpertvector
    68516902        Overflow_Error = FALSE;
    6852 
    68536903        next_vect = MkInterRedNextWeight(omega,omega2,G);
    68546904        xtnw=xtnw+clock()-to;
    … …  
    68626912//#endif
    68636913
    6864     /* check whether the the computed vector is in the correct cone */
    6865     /* If no, the reduced GB of an omega-homogeneous ideal will be
    6866        computed by Buchberger algorithm and stop this recursion step*/
    6867     //if(test_w_in_ConeCC(G, next_vect) != 1) //e.g. Example s7, cyc6
    6868     if(Overflow_Error == TRUE)
     6914    // check whether the the computed vector is in the correct cone.
     6915    // If no, compute the reduced Groebner basis of an omega-homogeneous
     6916    // ideal with Buchberger's algorithm and stop this recursion step
     6917    if(Overflow_Error == TRUE || test_w_in_ConeCC(G, next_vect) != 1)  //e.g. Example s7, cyc6
    68696918    {
    68706919      delete next_vect;
    … …  
    68866935      if(printout > 0)
    68876936      {
    6888         Print("\n//** rec_fractal_call: Overflow. Applying Buchberger's algorithm in ring r = %s;",
     6937        Print("\n//** rec_fractal_call: Applying Buchberger's algorithm in ring r = %s;",
    68896938              rString(currRing));
    68906939      }
    … …  
    69056954
    69066955      nnflow ++;
    6907 
    69086956      Overflow_Error = FALSE;
    69096957      return (G1);
    6910     }//end overflow-check
     6958    }
    69116959
    69126960
    … …  
    69516999        if(printout > 0)
    69527000        {
    6953           Print("\n//** rec_fractal_call: Wrong cone. Tau' doesn't stay in the correct cone.\n");
     7001          Print("\n//** rec_fractal_call: Wrong cone. Tau doesn't stay in the correct cone.\n");
    69547002        }
    6955 /*
     7003
    69567004#ifndef  MSTDCC_FRACTAL
    69577005        intvec* Xtautmp;
    … …  
    69907038        goto NEXT_VECTOR_FRACTAL;
    69917039#endif
    6992 */
     7040
    69937041      FRACTAL_MSTDCC:
    69947042        if(printout > 0)
    … …  
    70637111
    70647112    to=clock();
    7065     /* Take the initial form of <G> w.r.t. omega */
     7113    // Take the initial form of <G> w.r.t. omega
    70667114    Gomega = MwalkInitialForm(G, omega);
    70677115    xtif=xtif+clock()-to;
    … …  
    70737121    if(reduction == 0)
    70747122    {
    7075       /* Check whether the intermediate weight vector lies in the interior of the cone.
    7076        * If so, only perform reductions. Otherwise apply Buchberger's algorithm. */
     7123      // Check whether the intermediate weight vector lies in the interior of the cone.
     7124      // If so, only perform reductions. Otherwise apply Buchberger's algorithm.
    70777125      FF = middleOfCone(G,Gomega);
    70787126      if( FF != NULL)
    … …  
    70817129        G = idCopy(FF);
    70827130        idDelete(&FF);
    7083         /* Compue next vector. */
     7131        // Compue next vector.
    70847132        goto NEXT_VECTOR_FRACTAL;
    70857133      }
    … …  
    71587206    G = idrMoveR(F,oRing,currRing);
    71597207    to=clock();
    7160     /* Interreduce G */
    7161    // G = kInterRedCC(F1, NULL);
     7208    // Interreduce G
     7209    // G = kInterRedCC(F1, NULL);
    71627210    xtred=xtred+clock()-to;
    71637211    //idDelete(&F1);
    … …  
    72397287    if(polylength > 0 && G->m[0] != NULL)
    72407288    {
    7241       /*
    7242       there is a polynomial in Gomega with at least 3 monomials,
    7243       low-dimensional facet of the cone
    7244       */
     7289     
     7290      PrintS("\n**// rec_r_fractal_call: there is a polynomial in Gomega with at least 3 monomials, low-dimensional facet of the cone.\n");
     7291     
    72457292      delete next_vect;
    72467293      next_vect = MWalkRandomNextWeight(G,omega,omega2,weight_rad,nlev);
    … …  
    73247371        if(G->m[0] != NULL && polylength > 0)
    73257372        {
    7326           /*
    7327           there is a polynomial in Gomega with at least 3 monomials,
    7328           low-dimensional facet of the cone
    7329           */
     7373         
     7374          PrintS("//** rec_r_fractal_call: there is a polynomial in Gomega with at least 3 monomials, low-dimensional facet of the cone");
     7375         
    73307376          delete next_vect;
    73317377          next_vect = MWalkRandomNextWeight(G,omega,omega2,weight_rad,nlev);
    … …  
    73497395       computed by Buchberger algorithm and stop this recursion step*/
    73507396    //if(test_w_in_ConeCC(G, next_vect) != 1) //e.g. Example s7, cyc6
    7351     if(Overflow_Error == TRUE)
     7397    if(Overflow_Error == TRUE || test_w_in_ConeCC(G,next_vect) != 1)
    73527398    {
    73537399      delete next_vect;
Note: See TracChangeset for help on using the changeset viewer.