Changeset 3b45cf in git


Ignore:
Timestamp:
May 5, 2009, 5:15:02 PM (14 years ago)
Author:
Martin Monerjan
Branches:
(u'spielwiese', '0d6b7fcd9813a1ca1ed4220cfa2b104b97a0a003')
Children:
8ca3c766834164238a40705d1b906a5f5e003d3f
Parents:
def568e9169e67d58842356d8b9718f9a255b706
Message:
Fixed bug in flippability
More work on gcone::isSearchFacet
Crashes after output of polynomials  with:

error: no more memory
System 675k:675k Appl 216k/458k Malloc 163k/0k Valloc 512k/458k Pages 31/97 Regions 1:1

halt 14


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

Legend:

Unmodified
Added
Removed
  • kernel/gfan.cc

    rdef568 r3b45cf  
    22Compute the Groebner fan of an ideal
    33$Author: monerjan $
    4 $Date: 2009-05-04 15:14:38 $
    5 $Header: /exports/cvsroot-2/cvsroot/kernel/gfan.cc,v 1.44 2009-05-04 15:14:38 monerjan Exp $
    6 $Id: gfan.cc,v 1.44 2009-05-04 15:14:38 monerjan Exp $
     4$Date: 2009-05-05 15:15:02 $
     5$Header: /exports/cvsroot-2/cvsroot/kernel/gfan.cc,v 1.45 2009-05-05 15:15:02 monerjan Exp $
     6$Id: gfan.cc,v 1.45 2009-05-05 15:15:02 monerjan Exp $
    77*/
    88
     
    7575
    7676                       
    77         public:
     77        public:         
     78                //bool isFlippable;     //flippable facet? Want to have cone->isflippable.facet[i]
     79                bool isIncoming;        //Is the facet incoming or outgoing?
     80                facet *next;            //Pointer to next facet
     81               
    7882                /** The default constructor. Do I need a constructor of type facet(intvec)? */
    7983                facet()                 
     
    125129                }
    126130               
    127                 bool isFlippable;       //flippable facet? Want to have cone->isflippable.facet[i]
    128                 bool isIncoming;        //Is the facet incoming or outgoing?
    129                 facet *next;            //Pointer to next facet
     131                /*bool isFlippable(intvec &load)
     132                {
     133                        bool res=TRUE;                 
     134                        int jj;
     135                        for (int jj = 0; jj<load.length(); jj++)
     136                        {
     137                                intvec *ivCanonical = new intvec(load.length());
     138                                (*ivCanonical)[jj]=1;                           
     139                                if (ivMult(&load,ivCanonical)<0)
     140                                {
     141                                        res=FALSE;
     142                                        break;
     143                                }
     144                        }
     145                        return res;
     146                       
     147                        /*while (dotProduct(load,ivCanonical)>=0)
     148                        {
     149                                if (jj!=this->numVars)
     150                                {
     151                                        intvec *ivCanonical = new intvec(this->numVars);
     152                                        (*ivCanonical)[jj]=1;                   
     153                                        res=TRUE;
     154                                        jj += 1;
     155                                }
     156                        }
     157                        if (jj==this->numVars)
     158                        {                       
     159                                delete ivCanonical;
     160                                return FALSE;
     161                        }
     162                        else
     163                        {
     164                                delete ivCanonical;
     165                                return TRUE;
     166                }*/                                             
     167                }//bool isFlippable(facet &f)
     168               
    130169               
    131170                friend class gcone;     //Bad style
     
    403442                                }//for (int jj = 1; jj <ddcols; jj++)
    404443                                /*Quick'n'dirty hack for flippability*/
    405                                 bool isFlippable;
     444                                bool isFlippable=FALSE;
    406445                                //NOTE BUG HERE
    407446                                /* The flippability check is wrong!
     
    409448                                REWRITE THIS
    410449                                */                     
    411                                 for (int jj = 0; jj<this->numVars; jj++)
     450                                /*for (int jj = 0; jj<this->numVars; jj++)
    412451                                {                                       
    413452                                        intvec *ivCanonical = new intvec(this->numVars);
     
    423462                                        delete ivCanonical;
    424463                                }//for (int jj = 0; jj<this->numVars; jj++)
     464                                */     
     465                                for (int jj = 0; jj<load->length(); jj++)
     466                                {
     467                                        intvec *ivCanonical = new intvec(load->length());
     468                                        (*ivCanonical)[jj]=1;                           
     469                                        if (dotProduct(load,ivCanonical)<0)
     470                                        {
     471                                                isFlippable=TRUE;
     472                                                break;  //URGHS
     473                                        }
     474                                }       
    425475                                if (isFlippable==FALSE)
    426476                                {
     
    10341084                * If this is the case, then our facet is indeed a search facet and TRUE is retuned.
    10351085                */
    1036                 bool isSearchFacet(gcone &tmpcone, facet &testfacet)
     1086                bool isSearchFacet(gcone &gcTmp, facet &testfacet)
    10371087                {                               
    10381088                        ring actRing=currRing;
    1039                         facet *facetPtr=(facet*)tmpcone.facetPtr;                       
     1089                        facet *facetPtr=(facet*)gcTmp.facetPtr;                 
    10401090                        facet *fMin=new facet(*facetPtr);       //Pointer to the "minimal" facet
    10411091                        //facet *fMin = new facet(tmpcone.facetPtr);
    10421092                        //fMin=tmpcone.facetPtr;                //Initialise to first facet of tmpcone
    1043                         facet *fAct;
    1044                         //fMin = testfacet;
    1045                         //fMin->next=testfacet->next;
     1093                        facet *fAct;    //Ptr to alpha_i
     1094                        facet *fCmp;    //Ptr to alpha_j
    10461095                        fAct = fMin;
    1047                        
    1048                         cout << endl<< fMin->next << endl;
    1049                         cout << tmpcone.facetPtr->next << endl;
    1050                         tmpcone.facetPtr->printNormal();
    1051                         fAct=tmpcone.facetPtr->next;
    1052                         fAct->printNormal();
    1053                         /*fMin->printNormal();
    1054                         fMin->next->printNormal();*/
     1096                        fCmp = fMin->next;
     1097                       
     1098                        //cout << endl<< fMin->next << endl;
     1099                        //cout << gcTmp.facetPtr->next << endl;
     1100                        //gcTmp.facetPtr->printNormal();
     1101                        //fAct=gcTmp.facetPtr->next;
     1102                        //fAct->printNormal();                 
    10551103                       
    10561104                        rChangeCurrRing(this->rootRing);        //because we compare the monomials in the rootring                     
    1057                         poly p=NULL;
    1058                         poly q=NULL;
    1059                         intvec *p_weight = new intvec(this->numVars);                   
    1060                         intvec *q_weight = new intvec(this->numVars);
     1105                        poly p=pInit();
     1106                        poly q=pInit();
     1107                        intvec *alpha_i = new intvec(this->numVars);                   
     1108                        intvec *alpha_j = new intvec(this->numVars);
    10611109                        intvec *sigma = new intvec(this->numVars);
    1062                         sigma=tmpcone.getIntPoint();
    1063                         cout << "pling" << endl;
     1110                        sigma=gcTmp.getIntPoint();
     1111                       
    10641112                        int *u=(int *)omAlloc((this->numVars+1)*sizeof(int));
    10651113                        int *v=(int *)omAlloc((this->numVars+1)*sizeof(int));
     1114                        u[0]=0; v[0]=0;
    10661115                        int weight1,weight2;
    1067                         while(fAct->next!=NULL) //ersetzen durch fAct
     1116                        while(fAct->next->next!=NULL)   //NOTE this is ugly. Can it be done without fCmp?
    10681117                        {
    10691118                                /* Get alpha_i and alpha_{i+1} */
    1070                                 p_weight=fMin->getFacetNormal();
    1071                                 q_weight=fMin->next->getFacetNormal();
    1072                                 p_weight->show();
    1073                                 q_weight->show();
     1119                                alpha_i=fAct->getFacetNormal();
     1120                                alpha_j=fCmp->getFacetNormal();
     1121#ifdef gfan_DEBUG
     1122                                alpha_i->show();
     1123                                alpha_j->show();
     1124#endif
    10741125                                /*Compute the dot product of sigma and alpha_{i,j}*/
    1075                                 weight1=dotProduct(sigma,p_weight);
    1076                                 weight2=dotProduct(sigma,q_weight);
     1126                                weight1=dotProduct(sigma,alpha_i);
     1127                                weight2=dotProduct(sigma,alpha_j);
     1128#ifdef gfan_DEBUG
    10771129                                cout << "weight1=" << weight1 << " " << "weight2=" << weight2 << endl;
     1130#endif
    10781131                                /*Adjust alpha_i and alpha_i+1 accordingly*/
    10791132                                for(int ii=1;ii<=this->numVars;ii++)
    1080                                 {
    1081                                         /*p_weight[ii]=weight1*(*p_weight)[ii];
    1082                                         q_weight[ii]=weight2*(*q_weight)[ii];*/
    1083                                         u[ii]=weight1*(*p_weight)[ii];
    1084                                         v[ii]=weight2*(*q_weight)[ii];
    1085                                 }
    1086                                 cout << "PLONK" << endl;
     1133                                {                                       
     1134                                        u[ii]=weight1*(*alpha_i)[ii-1];
     1135                                        v[ii]=weight2*(*alpha_j)[ii-1];
     1136                                }                               
    10871137                               
    10881138                                /*Now p_weight and q_weight need to be compared as exponent vectors*/
    1089                                 pSetExpV(p,u);
    1090                                 pSetExpV(q,v);
    1091                                 cout << "AARGH" << endl;
     1139                                pSetCoeff0(p,nInit(1));
     1140                                pSetCoeff0(q,nInit(1));
     1141                                pSetExpV(p,u);
     1142                                pSetm(p);                       
     1143                                pSetExpV(q,v);
     1144                                pSetm(q);
     1145#ifdef gfan_DEBUG                               
     1146                                pWrite(p);pWrite(q);
     1147#endif 
    10921148                                /*We want to check whether x^p < x^q
    10931149                                => want to check for return value 1 */
    10941150                                if (pLmCmp(p,q)==1)     //i.e. x^q is smaller
    10951151                                {
    1096                                         fMin=fMin->next;
     1152                                        fMin=fCmp;
    10971153                                        fAct=fMin;
    10981154                                }
    10991155                                else
    11001156                                {
    1101                                         fAct=fAct->next;
     1157                                        //fAct=fAct->next;
     1158                                        if(fCmp->next!=NULL)
     1159                                        {
     1160                                                fCmp=fCmp->next;
     1161                                        }
     1162                                        else
     1163                                        {
     1164                                                fAct=fAct->next;
     1165                                        }
    11021166                                }
    1103                                 fAct=fAct->next;
    1104                         }//while(tmpcone.facetPtr->next!=NULL)
     1167                                //fAct=fAct->next;
     1168                        }//while(fAct.facetPtr->next!=NULL)
    11051169                       
    11061170                        /*If testfacet was minimal then fMin should still point there */
    1107                         //if (isParallel(fMin->getFacetNormal,testfacet.getFacetNormal))
    1108                         if (fMin==tmpcone.facetPtr)
     1171                        intvec *alpha_min = new intvec(this->numVars);
     1172                        alpha_min=fMin->getFacetNormal();
     1173                        intvec *test = new intvec(this->numVars);
     1174                        test=testfacet.getFacetNormal();
     1175                        if (isParallel(alpha_min,test))
     1176                        //if (fMin==gcTmp.facetPtr)
    11091177                        {
    11101178                                rChangeCurrRing(actRing);
     
    11561224                        }//while(fAct->next!=NULL)
    11571225                }//reverseSearch
     1226                friend class facet;
    11581227};//class gcone
    11591228
Note: See TracChangeset for help on using the changeset viewer.