Changeset a7efc97 in git


Ignore:
Timestamp:
Jan 11, 2010, 1:19:18 PM (13 years ago)
Author:
Martin Monerjan
Branches:
(u'spielwiese', '0d6b7fcd9813a1ca1ed4220cfa2b104b97a0a003')
Children:
230d827d0f5f0d0063755ae4baf1c60f8945d447
Parents:
9524c6afef1a341d40fe3a65abdd25dd91db5ff2
Message:
Speedup by removal of several cddlib calls
Simpler check for flippability in getCodim2Normals


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

Legend:

Unmodified
Added
Removed

  • kernel/gfan.cc

    r9524c6 ra7efc97  
    661661                        aktpoly=pNext(aktpoly);
    662662                        int *tailexpv=(int*)omAlloc(((this->numVars)+1)*sizeof(int));
    663                         pGetExpV(aktpoly,tailexpv);
     663                        pGetExpV(aktpoly,tailexpv);                     
    664664                        for(int kk=1;kk<=this->numVars;kk++)
    665                         {
     665                        {                               
    666666                                dd_set_si(ddineq->matrix[(dd_rowrange)aktmatrixrow][kk],leadexpv[kk]-tailexpv[kk]);
    667667                        }
    … …  
    990990                dd_MatrixPtr ddakt;
    991991                ddakt = dd_CopyMatrix(ddineq);
    992                 ddakt->representation=dd_Inequality;
     992//              ddakt->representation=dd_Inequality;    //Not using this makes it faster. But why does the quick check below still work?
    993993                set_addelem(ddakt->linset,ii+1);/*Now set appropriate linearity*/
    994994#ifdef gfanp
    … …  
    996996                gettimeofday(&t_ddMC_start,0);
    997997#endif                         
    998                 /*dd_MatrixCanonicalize(&ddakt, &impl_linset, &redset, &newpos, &err);*/
     998                //dd_MatrixCanonicalize(&ddakt, &impl_linset, &redset, &newpos, &err);
    999999//              set_copy(LL,ddakt->linset);
    10001000                dd_PolyhedraPtr ddpolyh;
    10011001                ddpolyh=dd_DDMatrix2Poly(ddakt, &err);
     1002//              ddpolyh=dd_DDMatrix2Poly2(ddakt, dd_MaxCutoff, &err);
    10021003                dd_MatrixPtr P;
    1003                 P=dd_CopyGenerators(ddpolyh);
     1004                P=dd_CopyGenerators(ddpolyh);           
    10041005                dd_FreePolyhedra(ddpolyh);
     1006//              ddpolyh=dd_DDMatrix2Poly(P,&err);
     1007//              dd_MatrixPtr R;
     1008//              R=dd_CopyInequalities(ddpolyh);
     1009//              dd_FreePolyhedra(ddpolyh);
    10051010#ifdef gfanp
    10061011                gettimeofday(&t_ddMC_end,0);
    … …  
    10471052                * linearity for ddakt above.
    10481053                */
    1049                 intvec *iv_intPoint = new intvec(this->numVars);
    1050                 dd_MatrixPtr shiftMatrix;
    1051                 dd_MatrixPtr intPointMatrix;
    1052                 shiftMatrix = dd_CreateMatrix(this->numVars,this->numVars+1);
    1053                 for(int kk=0;kk<this->numVars;kk++)
    1054                 {
    1055                         dd_set_si(shiftMatrix->matrix[kk][0],1);
    1056                         dd_set_si(shiftMatrix->matrix[kk][kk+1],1);
    1057                 }
    1058                 intPointMatrix=dd_MatrixAppend(ddakt,shiftMatrix);
     1054//              intvec *iv_intPoint = new intvec(this->numVars);
     1055//              dd_MatrixPtr shiftMatrix;
     1056//              dd_MatrixPtr intPointMatrix;
     1057//              shiftMatrix = dd_CreateMatrix(this->numVars,this->numVars+1);
     1058//              for(int kk=0;kk<this->numVars;kk++)
     1059//              {
     1060//                      dd_set_si(shiftMatrix->matrix[kk][0],1);
     1061//                      dd_set_si(shiftMatrix->matrix[kk][kk+1],1);
     1062//              }
     1063//              intPointMatrix=dd_MatrixAppend(ddakt,shiftMatrix);
     1064// #ifdef gfanp
     1065//              timeval t_iP_start, t_iP_end;
     1066//              gettimeofday(&t_iP_start, 0);
     1067// #endif               
     1068//              interiorPoint(intPointMatrix,*iv_intPoint);
     1069// //           dd_rowset impl_linste,lbasis;
     1070// //           dd_LPSolutionPtr lps=NULL;
     1071// //           dd_ErrorType err;
     1072// //           dd_FindRelativeInterior(intPointMatrix, &impl_linset, &lbasis, &lps, &err);
     1073// #ifdef gfanp
     1074//              gettimeofday(&t_iP_end, 0);
     1075//              t_iP += (t_iP_end.tv_sec - t_iP_start.tv_sec + 1e-6*(t_iP_end.tv_usec - t_iP_start.tv_usec));
     1076// #endif
     1077//              for(int ll=0;ll<this->numVars;ll++)
     1078//              {
     1079//                      if( (*iv_intPoint)[ll] < 0 )
     1080//                      {
     1081//                              fAct->isFlippable=FALSE;
     1082//                              break;
     1083//                      }
     1084//              }
     1085                /*End of check*/
     1086                /*This test should be way less time consuming*/
    10591087#ifdef gfanp
    10601088                timeval t_iP_start, t_iP_end;
    10611089                gettimeofday(&t_iP_start, 0);
    1062 #endif         
    1063                 interiorPoint(intPointMatrix,*iv_intPoint);
    1064 //              dd_rowset impl_linste,lbasis;
    1065 //              dd_LPSolutionPtr lps=NULL;
    1066 //              dd_ErrorType err;
    1067 //              dd_FindRelativeInterior(intPointMatrix, &impl_linset, &lbasis, &lps, &err);
     1090#endif
     1091                bool containsStrictlyPosRay=TRUE;
     1092                for(int ii=0;ii<ddakt->rowsize;ii++)
     1093                {
     1094                        containsStrictlyPosRay=TRUE;
     1095                        for(int jj=1;jj<this->numVars;jj++)
     1096                        {
     1097                                if(ddakt->matrix[ii][jj]<=0)
     1098                                {
     1099                                        containsStrictlyPosRay=FALSE;
     1100                                        break;
     1101                                }
     1102                        }
     1103                        if(containsStrictlyPosRay==TRUE)
     1104                                break;
     1105                }
     1106                if(containsStrictlyPosRay==FALSE)
     1107                        fAct->isFlippable=FALSE;
    10681108#ifdef gfanp
    10691109                gettimeofday(&t_iP_end, 0);
    10701110                t_iP += (t_iP_end.tv_sec - t_iP_start.tv_sec + 1e-6*(t_iP_end.tv_usec - t_iP_start.tv_usec));
    10711111#endif
    1072                 for(int ll=0;ll<this->numVars;ll++)
    1073                 {
    1074                         if( (*iv_intPoint)[ll] < 0 )
    1075                         {
    1076                                 fAct->isFlippable=FALSE;
    1077                                 break;
    1078                         }
    1079                 }
    1080                 /*End of check*/                                       
     1112                /**/
    10811113                fAct = fAct->next;     
    10821114                dd_FreeMatrix(ddakt);
    1083                 dd_FreeMatrix(shiftMatrix);
    1084                 dd_FreeMatrix(intPointMatrix);         
    1085                 delete iv_intPoint;
     1115//              dd_FreeMatrix(shiftMatrix);
     1116//              dd_FreeMatrix(intPointMatrix);         
     1117//              delete iv_intPoint;
    10861118                dd_FreeMatrix(P);
    10871119//              set_free(impl_linset);
    … …  
    12881320                                markingsAreCorrect=TRUE;
    12891321                        }
    1290                        
     1322                        int ctr=0;
    12911323                        for (int jj=0;jj<this->numVars;jj++)
    12921324                        {                               
    1293                                 /*Store into ddMatrix*/                                                         
     1325                                /*Store into ddMatrix*/                         
     1326                                if(leadExpV[jj+1]-v[jj+1])
     1327                                        ctr++;
    12941328                                dd_set_si(intPointMatrix->matrix[aktrow][jj+1],leadExpV[jj+1]-v[jj+1]);
    12951329                        }
    1296                         aktrow +=1;
     1330                        /*It ought to be more sensible to avoid 0-rows in the first place*/
     1331//                      if(ctr==this->numVars)//We have a 0-row
     1332//                              dd_MatrixRowRemove(&intPointMatrix,aktrow);
     1333//                      else
     1334                                aktrow +=1;
    12971335                        omFree(v);
    12981336                }
    … …  
    13061344        /*Now it is safe to idDelete(H)*/
    13071345        idDelete(&H);
     1346        /*Preprocessing goes here since otherwise we would delete the constraint
     1347        * for the standard simplex.
     1348        */
     1349        preprocessInequalities(intPointMatrix);
    13081350        /*Now we add the constraint for the standard simplex*/
    1309 // #ifdef gfanp
    1310 //      timeval t_dd_start, t_dd_end;
    1311 //      gettimeofday(&t_dd_start, 0);
    1312 // #endif
    1313         dd_set_si(intPointMatrix->matrix[aktrow][0],-1);
    1314         for (int jj=1;jj<=this->numVars;jj++)
    1315         {
    1316                 dd_set_si(intPointMatrix->matrix[aktrow][jj],1);
    1317         }
     1351//      dd_set_si(intPointMatrix->matrix[aktrow][0],-1);
     1352//      for (int jj=1;jj<=this->numVars;jj++)
     1353//      {
     1354//              dd_set_si(intPointMatrix->matrix[aktrow][jj],1);
     1355//      }
    13181356        //Let's make sure we compute interior points from the positive orthant
    1319         dd_MatrixPtr posRestr=dd_CreateMatrix(this->numVars,this->numVars+1);
    1320        
    1321         int jj=1;
    1322         for (int ii=0;ii<this->numVars;ii++)
    1323         {
    1324                 dd_set_si(posRestr->matrix[ii][jj],1);
    1325                 jj++;                                                   
    1326         }
     1357//      dd_MatrixPtr posRestr=dd_CreateMatrix(this->numVars,this->numVars+1);
     1358//     
     1359//      int jj=1;
     1360//      for (int ii=0;ii<this->numVars;ii++)
     1361//      {
     1362//              dd_set_si(posRestr->matrix[ii][jj],1);
     1363//              jj++;                                                   
     1364//      }
     1365        /*We create a matrix containing the standard simplex
     1366        * and constraints to assure a strictly positive point
     1367        * is computed
     1368        */
     1369        dd_MatrixPtr posRestr = dd_CreateMatrix(this->numVars+1, this->numVars+1);
     1370        for(int ii=0;ii<posRestr->rowsize;ii++)
     1371        {
     1372                if(ii==0)
     1373                {
     1374                        dd_set_si(posRestr->matrix[ii][0],-1);
     1375                        for(int jj=1;jj<=this->numVars;jj++)                   
     1376                                dd_set_si(posRestr->matrix[ii][jj],1);                 
     1377                }
     1378                else
     1379                {
     1380                        dd_set_si(posRestr->matrix[ii][0],-1);
     1381                        dd_set_si(posRestr->matrix[ii][ii],1);
     1382                }
     1383        }
     1384        /**/
    13271385        dd_MatrixAppendTo(&intPointMatrix,posRestr);
    13281386        dd_FreeMatrix(posRestr);
    … …  
    13301388        * Otherwise things will go pear shaped if called after the standard simplex constraints are added
    13311389        */
    1332         preprocessInequalities(intPointMatrix);
     1390//      preprocessInequalities(intPointMatrix);
    13331391        intvec *iv_weight = new intvec(this->numVars);
    13341392#ifdef gfanp
    … …  
    13411399       
    13421400//      dd_MatrixCanonicalize(&intPointMatrix,&implLin,&redrows,&newpos,&err);
    1343         dd_MatrixCanonicalizeLinearity(&intPointMatrix,&implLin, &newpos,&err);
     1401//      dd_MatrixCanonicalizeLinearity(&intPointMatrix,&implLin, &newpos,&err);
    13441402        //dd_MatrixRedundancyRemove is our time sink!
    13451403//      dd_MatrixRedundancyRemove(&intPointMatrix,&redrows,&newpos,&err);
    … …  
    15381596void gcone::preprocessInequalities(dd_MatrixPtr &ddineq)
    15391597{
    1540 //Remove strictly positive rows
    1541         int *posRowsArray=NULL;
     1598/*      int *posRowsArray=NULL;
    15421599        int num_alloc=0;
    15431600        int num_elts=0;
    1544 //      for(int ii=0;ii<ddineq->rowsize;ii++)
    1545 //      {
    1546 //              intvec *ivPos = new intvec(pVariables);
    1547 //              for(int jj=0;jj<pVariables;jj++)
    1548 //                      (*ivPos)[jj]=(int)mpq_get_d(ddineq->matrix[ii][jj+1]);
    1549 //              bool isStrictlyPos=FALSE;
    1550 //              int posCtr=0;           
    1551 //              for(int jj=0;jj<pVariables;jj++)
    1552 //              {
    1553 //                      intvec *ivCanonical = new intvec(pVariables);
    1554 //                      jj==0 ? (*ivCanonical)[ivPos->length()-1]=1 : (*ivCanonical)[jj-1]=1;
    1555 //                      if(dotProduct(*ivCanonical,*ivPos)!=0)
    1556 //                      {
    1557 //                              if ((*ivPos)[jj]>=0)
    1558 //                              {                               
    1559 //                                      posCtr++;                               
    1560 //                              }
    1561 //                      }                       
    1562 //                      delete ivCanonical;
    1563 //              }
    1564 //              if(posCtr==ivPos->length())
    1565 //                      isStrictlyPos=TRUE;
    1566 //              if(isStrictlyPos==TRUE)
    1567 //              {
    1568 //                      if(num_alloc==0)
    1569 //                              num_alloc += 1;
    1570 //                      else
    1571 //                              num_alloc += 1;
    1572 //                      void *tmp = realloc(posRowsArray,(num_alloc*sizeof(int)));
    1573 //                      if(!tmp)
    1574 //                      {
    1575 //                              WerrorS("Woah dude! Couldn't realloc memory\n");
    1576 //                              exit(-1);
    1577 //                      }
    1578 //                      posRowsArray = (int*)tmp;
    1579 //                      posRowsArray[num_elts]=ii;
    1580 //                      num_elts++;     
    1581 //              }
    1582 //              delete ivPos;
    1583 //      }
    1584         int offset=0;
    1585 //      for(int ii=0;ii<num_elts;ii++)
    1586 //      {
    1587 //              dd_MatrixRowRemove(&ddineq,posRowsArray[ii]+1-offset);
    1588 //              offset++;
    1589 //      }
    1590 //      free(posRowsArray);
    1591         //Remove zeroes
    1592 //      int rowsize=ddineq->rowsize;
     1601        int offset=0;*/
     1602        //Remove zeroes (and strictly pos rows?)
    15931603        for(int ii=0;ii<ddineq->rowsize;ii++)
    15941604        {
    … …  
    15981608                {
    15991609                        (*iv)[jj]=(int)mpq_get_d(ddineq->matrix[ii][jj+1]);
    1600 //                      if((*iv)[ii]>0)
    1601 //                              posCtr++;
     1610                        if((*iv)[jj]>0)//check for strictly pos rows
     1611                                posCtr++;
     1612                        //Behold! This will delete the row for the standard simplex!
    16021613                }
    1603                 if( (iv->compare(0)==0) )//|| (posCtr==iv->length()) )
     1614//              if( (iv->compare(0)==0) || (posCtr==iv->length()) )
     1615                if( (posCtr==iv->length()) || (iv->compare(0)==0) )             
    16041616                {
    16051617                        dd_MatrixRowRemove(&ddineq,ii+1);
    1606                         ii--;           
    1607 //                      rowsize=ddineq->rowsize;
     1618                        ii--;//Yes. This is on purpose
    16081619                }
    16091620                delete iv;
Note: See TracChangeset for help on using the changeset viewer.