Changeset 229d22 in git


Ignore:
Timestamp:
Apr 24, 2009, 5:23:16 PM (15 years ago)
Author:
Martin Monerjan
Branches:
(u'spielwiese', 'fe61d9c35bf7c61f2b6cbf1b56e25e2f08d536cc')
Children:
fa30ac48535b5e4377fcd1c33b14209a59454ae3
Parents:
fb5992617c66ac356cd294471f9005dad76ed84e
Message:
Fixed bug in ring creation (copied from ring.cc::rAssure_SyzComp)
Began method gcone::rCopyAndChangeWeight; does not work yet. Therefore code
is doubled in gcone::flip


git-svn-id: file:///usr/local/Singular/svn/trunk@11752 2c84dea3-7e68-4137-9b89-c4e89433aadc
File:
1 edited

Legend:

Unmodified
Added
Removed

  • kernel/gfan.cc

    rfb5992 r229d22  
    22Compute the Groebner fan of an ideal
    33$Author: monerjan $
    4 $Date: 2009-04-23 09:44:58 $
    5 $Header: /exports/cvsroot-2/cvsroot/kernel/gfan.cc,v 1.38 2009-04-23 09:44:58 monerjan Exp $
    6 $Id: gfan.cc,v 1.38 2009-04-23 09:44:58 monerjan Exp $
     4$Date: 2009-04-24 15:23:16 $
     5$Header: /exports/cvsroot-2/cvsroot/kernel/gfan.cc,v 1.39 2009-04-24 15:23:16 monerjan Exp $
     6$Id: gfan.cc,v 1.39 2009-04-24 15:23:16 monerjan Exp $
    77*/
    88
    … …  
    136136        private:
    137137                int numFacets;          //#of facets of the cone
    138                 ring rootRing;         
     138                ring rootRing;          //good to know this -> generic walk
    139139                ideal inputIdeal;       //the original
    140140        public:
    … …  
    149149                }
    150150               
     151                /** \brief Constructor with ring and ideal
     152                *
     153                * This constructor takes the root ring and the root ideal as parameters and stores
     154                * them in the private members gcone::rootRing and gcone::inputIdeal
     155                */
    151156                gcone(ring r, ideal I)
    152157                {
    … …  
    162167                facet *facetPtr;        //Will hold the adress of the first facet; set by gcone::getConeNormals
    163168                poly gcMarkedTerm;      //marked terms of the cone's Groebner basis
     169                /** # of variables in the ring */
    164170                int numVars;            //#of variables in the ring
    165171               
    … …  
    306312                                        double foo;
    307313                                        foo = mpq_get_d(ddineq->matrix[kk][jj]);
    308 #ifdef gfan_DEBUG
     314/*#ifdef gfan_DEBUG
    309315                                        std::cout << "fAct is " << foo << " at " << fAct << std::endl;
    310 #endif
     316#endif*/
    311317                                        (*load)[jj-1] = (int)foo;               //store typecasted entry at pos jj-1 of load
    312318#endif
    … …  
    314320                                        double *foo;
    315321                                        foo = (double*)ddineq->matrix[kk][jj];  //get entry from actual position#endif
    316 #ifdef gfan_DEBUG
     322/*#ifdef gfan_DEBUG
    317323                                        std::cout << "fAct is " << *foo << " at " << fAct << std::endl;
    318 #endif
     324#endif*/
    319325                                        (*load)[jj-1] = (int)*foo;              //store typecasted entry at pos jj-1 of load
    320326#endif //GMPRATIONAL                                   
    … …  
    383389                * is parallel to \f$ leadexp(g) \f$
    384390                * Parallelity is checked using basic linear algebra. See gcone::isParallel.
    385                 * Other possibilities includes computing the rank of the matrix consisting of the vectors in question and
     391                * Other possibilities include computing the rank of the matrix consisting of the vectors in question and
    386392                * computing an interior point of the facet and taking all terms having the same weight with respect
    387393                * to this interior point.
    … …  
    460466                        */
    461467                        ring srcRing=currRing;
     468                       
     469                        /* copied and modified from ring.cc::rAssureSyzComp */
     470                        //ring tmpRing=rCopyAndChangeWeight(srcRing,fNormal);
    462471                        ring tmpRing=rCopy0(srcRing);
     472                        int i=rBlocks(srcRing);
     473                        int j;
     474                        tmpRing->order=(int *)omAlloc((i+1)*sizeof(int));
     475                        /*NOTE This should probably be set, but as of now crashes Singular*/
     476                        /*tmpRing->block0=(int *)omAlloc0((i+1)*sizeof(int));
     477                        tmpRing->block1=(int *)omAlloc0((i+1)*sizeof(int));*/
     478                        tmpRing->wvhdl[0] =( int *)omAlloc((fNormal->length())*sizeof(int)); //found in Singular/ipshell.cc
     479                        for(j=i;j>0;j--)
     480                        {
     481                                tmpRing->order[j]=srcRing->order[j-1];
     482                                tmpRing->block0[j]=srcRing->block0[j-1];
     483                                tmpRing->block1[j]=srcRing->block1[j-1];
     484                                if (srcRing->wvhdl[j-1] != NULL)
     485                                {
     486                                        tmpRing->wvhdl[j] = (int*) omMemDup(srcRing->wvhdl[j-1]);
     487                                }
     488                        }
    463489                        tmpRing->order[0]=ringorder_a;
    464                         tmpRing->wvhdl[0] =( int *)omAlloc((fNormal->length())*sizeof(int)); //found in Singular/ipshell.cc
     490                        tmpRing->order[1]=ringorder_dp;
     491                        tmpRing->order[2]=ringorder_C;
     492                        //tmpRing->wvhdl[0] =( int *)omAlloc((fNormal->length())*sizeof(int)); //found in Singular/ipshell.cc
    465493                       
    466494                        for (int ii=0;ii<this->numVars;ii++)
    … …  
    492520                        ideal srcRing_HH;                       
    493521                        srcRing_H=idrCopyR(H,tmpRing);
    494                         idShow(srcRing_H);                     
     522#ifdef gfan_DEBUG
     523                        cout << "srcRing_H = ";
     524                        idShow(srcRing_H); cout << endl;
     525#endif
    495526                        srcRing_HH=ffG(srcRing_H,this->gcBasis);               
    496                         idShow(srcRing_HH);
     527#ifdef gfan_DEBUG
     528                        cout << "srcRing_HH = ";
     529                        idShow(srcRing_HH); cout << endl;
     530#endif
    497531                        /*Substep 2.2.1
    498532                        Mark according to G_-\alpha
    … …  
    606640                        */
    607641                        ring dstRing=rCopy0(srcRing);
     642                        i=rBlocks(srcRing);
     643                       
     644                        dstRing->order=(int *)omAlloc((i+1)*sizeof(int));
     645                        for(j=i;j>0;j--)
     646                        {
     647                                dstRing->order[j]=srcRing->order[j-1];
     648                                dstRing->block0[j]=srcRing->block0[j-1];
     649                                dstRing->block1[j]=srcRing->block1[j-1];
     650                                if (srcRing->wvhdl[j-1] != NULL)
     651                                {
     652                                        dstRing->wvhdl[j] = (int*) omMemDup(srcRing->wvhdl[j-1]);
     653                                }
     654                        }
    608655                        dstRing->order[0]=ringorder_a;
    609                         //dstRing->order[1]=ringorder_dp;
     656                        dstRing->order[1]=ringorder_dp;
     657                        dstRing->order[2]=ringorder_C;                 
    610658                        dstRing->wvhdl[0] =( int *)omAlloc((iv_weight->length())*sizeof(int));
     659                       
    611660                        for (int ii=0;ii<this->numVars;ii++)
    612661                        {                               
    … …  
    615664                        rComplete(dstRing);
    616665                        rChangeCurrRing(dstRing);
     666#ifdef gfan_DEBUG
    617667                        rWrite(dstRing); cout << endl;
    618                         ideal dstRing_I;
    619                         //dstRing_I=idrCopyR(gb,srcRing);
     668#endif
     669                        ideal dstRing_I;                       
    620670                        dstRing_I=idrCopyR(srcRing_HH,srcRing);                 
    621671                        //validOpts<1>=TRUE;
     672#ifdef gfan_DEBUG
    622673                        idShow(dstRing_I);
    623                         ideal foo;
     674#endif                 
    624675                        BITSET save=test;
    625676                        test|=Sy_bit(OPT_REDSB);
    626                         test|=Sy_bit(6);        //OPT_DEBUG
    627                         //foo=kStd(idrCopyR(this->inputIdeal,this->rootRing),NULL,testHomog,NULL);
    628                         foo=kStd(dstRing_I,NULL,testHomog,NULL);
    629                         idShow(foo);
    630                         kInterRed(foo);
    631                         idSkipZeroes(foo);
    632                         idShow(foo);
     677                        test|=Sy_bit(6);        //OPT_DEBUG                                     
     678                        dstRing_I=kStd(idrCopyR(this->inputIdeal,this->rootRing),NULL,testHomog,NULL);                                 
     679                        kInterRed(dstRing_I);
     680                        idSkipZeroes(dstRing_I);
    633681                        test=save;
     682                        /*End of step 3 - reduction*/
     683                       
     684                        f->setFlipGB(dstRing_I);//store the flipped GB
     685                        f->printFlipGB();
    634686                       
    635687                }//void flip(ideal gb, facet *f)
    … …  
    784836                }//bool isParallel
    785837               
     838                /** \brief Compute an interior point of a given cone
     839                */
    786840                void interiorPoint(dd_MatrixPtr const &M, intvec &iv) //no const &M here since we want to remove redundant rows
    787841                {
    … …  
    828882                        if (err!=dd_NoError){cout << "Error during dd_CopyLPSolution" << endl;}                 
    829883                        //cout << "Tick 6" << endl;
    830                        
     884#ifdef gfan_DEBUG
    831885                        cout << "Interior point: ";
     886#endif
    832887                        for (int ii=1; ii<(lpSol->d)-1;ii++)
    833888                        {
     889#ifdef gfan_DEBUG
    834890                                dd_WriteNumber(stdout,lpSol->sol[ii]);
     891#endif
     892                                /* NOTE This works only as long as gmp returns fractions with the same denominator*/
    835893                                (iv)[ii-1]=(int)mpz_get_d(mpq_numref(lpSol->sol[ii]));  //looks evil, but does the trick
    836894                        }
    … …  
    839897                        dd_FreeLPData(lp);
    840898                        set_free(ddlinset);
    841                         set_free(ddredrows);
    842                        
    843                         /*At this point we have an interior point of type mpq_t
    844                         Need to convert to an intvec
    845                         NOTE: Resolved
    846                         */
     899                        set_free(ddredrows);                   
    847900                       
    848901                }//void interiorPoint(dd_MatrixPtr const &M)
     902               
     903                ring rCopyAndChangeWeight(ring const r, intvec const *ivw)                             
     904                {
     905                        ring res=rCopy0(r, FALSE, FALSE);
     906                        int i=rBlocks(r);
     907                        int j;
     908
     909                        res->order=(int *)omAlloc((i+1)*sizeof(int));
     910                        /*res->block0=(int *)omAlloc0((i+1)*sizeof(int));
     911                        res->block1=(int *)omAlloc0((i+1)*sizeof(int));
     912                        int ** wvhdl =(int **)omAlloc0((i+1)*sizeof(int**));*/
     913                        for(j=i;j>0;j--)
     914                        {
     915                                res->order[j]=r->order[j-1];
     916                                res->block0[j]=r->block0[j-1];
     917                                res->block1[j]=r->block1[j-1];
     918                                if (r->wvhdl[j-1] != NULL)
     919                                {
     920                                        res->wvhdl[j] = (int*) omMemDup(r->wvhdl[j-1]);
     921                                }
     922                        }
     923                        res->order[0]=ringorder_a;
     924                        res->order[1]=ringorder_dp;
     925                        res->order[2]=ringorder_C;
     926                        //res->wvhdl = wvhdl;
     927                       
     928                        res->wvhdl[0] =( int *)omAlloc((ivw->length())*sizeof(int));
     929                        for (int ii=0;ii<this->numVars;ii++)
     930                        {                               
     931                                res->wvhdl[0][ii]=(*ivw)[ii];                           
     932                        }                       
     933                       
     934                        rComplete(res);
     935                        return res;
     936                }//rCopyAndChange
    849937};//class gcone
    850938
    … …  
    858946        ring inputRing=currRing;        // The ring the user entered
    859947        ring rootRing;                  // The ring associated to the target ordering
    860         ideal res;
    861         //matrix ineq;                  //Matrix containing the boundary inequalities
     948        ideal res;     
    862949        facet *fRoot;
    863950       
    … …  
    874961        rootRing=rCopy0(currRing);
    875962        rootRing->order[0]=ringorder_lp;
     963        /*rootRing->order[0]=ringorder_a;
     964        rootRing->order[1]=ringorder_lp;
     965        rootRing->wvhdl[0] =( int *)omAlloc(numvar*sizeof(int));
     966        rootRing->wvhdl[0][1]=1;
     967        rootRing->wvhdl[0][2]=1;*/
    876968        rComplete(rootRing);
    877969        rChangeCurrRing(rootRing);
Note: See TracChangeset for help on using the changeset viewer.