Changeset 7a6a60 in git


Ignore:
Timestamp:
Apr 23, 2009, 11:44:58 AM (15 years ago)
Author:
Martin Monerjan
Branches:
(u'spielwiese', 'fe61d9c35bf7c61f2b6cbf1b56e25e2f08d536cc')
Children:
228b1c93b3768f0d85df493aa941c2e25f436cef
Parents:
b85587587523c7ec8f89a09d494710648e48569b
Message:
Started work on reduction step after lifting. Doesn't work though.
kStd somehow screws things up.
test|=Sy_bit(OPT_REDSB) doesn't seem to help anything.
Added more private data to class gcone


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

Legend:

Unmodified
Added
Removed

  • kernel/gfan.cc

    rb85587 r7a6a60  
    22Compute the Groebner fan of an ideal
    33$Author: monerjan $
    4 $Date: 2009-04-21 15:23:54 $
    5 $Header: /exports/cvsroot-2/cvsroot/kernel/gfan.cc,v 1.37 2009-04-21 15:23:54 monerjan Exp $
    6 $Id: gfan.cc,v 1.37 2009-04-21 15:23:54 monerjan Exp $
     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 $
    77*/
    88
    … …  
    2222#include "prCopy.h"
    2323#include <iostream.h>   //deprecated
     24#include <bitset>
    2425
    2526/*remove the following at your own risk*/
    … …  
    134135{
    135136        private:
    136                 int numFacets;          //#of facets of the cone               
    137 
     137                int numFacets;          //#of facets of the cone
     138                ring rootRing;         
     139                ideal inputIdeal;       //the original
    138140        public:
    139141                /** \brief Default constructor.
    … …  
    146148                        this->facetPtr=NULL;
    147149                }
     150               
     151                gcone(ring r, ideal I)
     152                {
     153                        this->next=NULL;
     154                        this->facetPtr=NULL;
     155                        this->rootRing=r;
     156                        this->inputIdeal=I;
     157                }
     158               
    148159                ~gcone();               //destructor
    149160               
    … …  
    449460                        */
    450461                        ring srcRing=currRing;
    451                         ring dstRing=rCopy0(srcRing);
    452                         dstRing->order[0]=ringorder_a;
    453                         dstRing->wvhdl[0] =( int *)omAlloc((fNormal->length())*sizeof(int)); //found in Singular/ipshell.cc
     462                        ring tmpRing=rCopy0(srcRing);
     463                        tmpRing->order[0]=ringorder_a;
     464                        tmpRing->wvhdl[0] =( int *)omAlloc((fNormal->length())*sizeof(int)); //found in Singular/ipshell.cc
    454465                       
    455466                        for (int ii=0;ii<this->numVars;ii++)
    456467                        {                               
    457                                 dstRing->wvhdl[0][ii]=-(*fNormal)[ii];  //What exactly am I doing here?
     468                                tmpRing->wvhdl[0][ii]=-(*fNormal)[ii];  //What exactly am I doing here?
    458469                                //cout << tmpring->wvhdl[0][ii] << endl;
    459470                        }
    460                         rComplete(dstRing);
    461                         rChangeCurrRing(dstRing);
     471                        rComplete(tmpRing);
     472                        rChangeCurrRing(tmpRing);
    462473                       
    463474                        rWrite(currRing); cout << endl;
    … …  
    480491                        ideal srcRing_H;
    481492                        ideal srcRing_HH;                       
    482                         srcRing_H=idrCopyR(H,dstRing);
     493                        srcRing_H=idrCopyR(H,tmpRing);
    483494                        idShow(srcRing_H);                     
    484495                        srcRing_HH=ffG(srcRing_H,this->gcBasis);               
    … …  
    518529                                        int *dst_ExpV = (int *)omAlloc((this->numVars+1)*sizeof(int));
    519530                                        pGetExpV(aktpoly,src_ExpV);
    520                                         rChangeCurrRing(dstRing);       //this ring change is crucial!
     531                                        rChangeCurrRing(tmpRing);       //this ring change is crucial!
    521532                                        pGetExpV(pCopy(H->m[ii]),dst_ExpV);
    522533                                        rChangeCurrRing(srcRing);
    … …  
    548559                                else
    549560                                {
    550                                         rChangeCurrRing(dstRing);
     561                                        rChangeCurrRing(tmpRing);
    551562                                        pGetExpV(pHead(H->m[ii]),leadExpV); //We use H->m[ii] as leading monomial
    552563                                        rChangeCurrRing(srcRing);
    … …  
    586597                        }
    587598                        dd_WriteMatrix(stdout,intPointMatrix);
    588                         interiorPoint(intPointMatrix);  //TODO This should finally return an intvec
     599                        intvec *iv_weight = new intvec(this->numVars);
     600                        interiorPoint(intPointMatrix, *iv_weight);      //iv_weight now contains the interior point
    589601                        dd_FreeMatrix(intPointMatrix);
    590                         //dd_free_global_constants();
     602                        dd_free_global_constants();
    591603                       
    592604                        /*Step 3
    593605                        turn the minimal basis into a reduced one
    594606                        */
     607                        ring dstRing=rCopy0(srcRing);
     608                        dstRing->order[0]=ringorder_a;
     609                        //dstRing->order[1]=ringorder_dp;
     610                        dstRing->wvhdl[0] =( int *)omAlloc((iv_weight->length())*sizeof(int));
     611                        for (int ii=0;ii<this->numVars;ii++)
     612                        {                               
     613                                dstRing->wvhdl[0][ii]=(*iv_weight)[ii];                         
     614                        }
     615                        rComplete(dstRing);
     616                        rChangeCurrRing(dstRing);
     617                        rWrite(dstRing); cout << endl;
     618                        ideal dstRing_I;
     619                        //dstRing_I=idrCopyR(gb,srcRing);
     620                        dstRing_I=idrCopyR(srcRing_HH,srcRing);                 
     621                        //validOpts<1>=TRUE;
     622                        idShow(dstRing_I);
     623                        ideal foo;
     624                        BITSET save=test;
     625                        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);
     633                        test=save;
    595634                       
    596635                }//void flip(ideal gb, facet *f)
    … …  
    745784                }//bool isParallel
    746785               
    747                 void interiorPoint(dd_MatrixPtr &M) //no const &M here since we want to remove redundant rows
     786                void interiorPoint(dd_MatrixPtr const &M, intvec &iv) //no const &M here since we want to remove redundant rows
    748787                {
    749788                        dd_LPPtr lp,lpInt;
    … …  
    782821                        }
    783822                       
    784                         dd_LPSolve(lpInt,solver,&err);  //This will not result in a point from the relative interior
     823                        //dd_LPSolve(lpInt,solver,&err);        //This will not result in a point from the relative interior
    785824                        if (err!=dd_NoError){cout << "Error during dd_LPSolve" << endl;}
    786825                        //cout << "Tick 5" << endl;
    … …  
    793832                        for (int ii=1; ii<(lpSol->d)-1;ii++)
    794833                        {
    795                                 dd_WriteNumber(stdout,lpSol->sol[ii]);                         
     834                                dd_WriteNumber(stdout,lpSol->sol[ii]);
     835                                (iv)[ii-1]=(int)mpz_get_d(mpq_numref(lpSol->sol[ii]));  //looks evil, but does the trick
    796836                        }
    797837                        dd_FreeLPSolution(lpSol);
    … …  
    800840                        set_free(ddlinset);
    801841                        set_free(ddredrows);
     842                       
    802843                        /*At this point we have an interior point of type mpq_t
    803844                        Need to convert to an intvec
     845                        NOTE: Resolved
    804846                        */
    805847                       
    … …  
    848890#endif 
    849891       
    850         gcone *gcRoot = new gcone();    //Instantiate the sink
     892        //gcone *gcRoot = new gcone();  //Instantiate the sink
     893        gcone *gcRoot = new gcone(rootRing,rootIdeal);
    851894        gcone *gcAct;
    852895        gcAct = gcRoot;
Note: See TracChangeset for help on using the changeset viewer.