Changeset 173


Ignore:
Timestamp:
Jun 5, 2003 9:20:30 AM (16 years ago)
Author:
forrest
Message:

Improvements to presolve

Location:
trunk
Files:
12 edited

Legend:

Unmodified
Added
Removed
  • trunk/Presolve.cpp

    r145 r173  
    181181                         double feasibilityTolerance,
    182182                         bool keepIntegers,
    183                          int numberPasses)
     183                         int numberPasses,
     184                         bool dropNames)
    184185{
    185186  ncols_ = si.getNumCols();
     
    194195
    195196  // Check matrix
    196   if (!si.clpMatrix()->allElementsInRange(&si,1.0e-20,1.0e20))
     197  if (!si.clpMatrix()->allElementsInRange(&si,si.getSmallElementValue(),
     198                                          1.0e20))
    197199    return NULL;
    198200
     
    208210    delete presolvedModel_;
    209211    presolvedModel_ = new ClpSimplex(si);
    210    
     212    if (dropNames)
     213      presolvedModel_->dropNames();
    211214
    212215    // drop integer information if wanted
     
    576579    const bool dual = ATOI("off")!=0;
    577580#else
     581    // normal
     582#if 1
    578583    const bool slackd = true;
    579584    const bool doubleton = true;
     
    584589    const bool duprow = true;
    585590    const bool dual = doDualStuff;
     591#else
     592    const bool slackd = false;
     593    const bool doubleton = true;
     594    const bool forcing = true;
     595    const bool ifree = false;
     596    const bool zerocost = false;
     597    const bool dupcol = false;
     598    const bool duprow = false;
     599    const bool dual = false;
     600#endif
    586601#endif
    587602   
     
    624639      // look for substitutions with no fill
    625640      int fill_level=2;
     641      //fill_level=10;
     642      //printf("** fill_level == 10 !\n");
    626643      int whichPass=0;
    627644      while (1) {
     
    684701        prob->numberRowsToDo_ = prob->numberNextRowsToDo_;
    685702        int kcheck;
    686         bool found;
     703        bool found=false;
    687704        kcheck=-1;
    688         found=false;
    689705        for (i=0;i<prob->numberNextRowsToDo_;i++) {
    690706          int index = prob->nextRowsToDo_[i];
  • trunk/PresolveDoubleton.cpp

    r151 r173  
    583583         
    584584          s->ncolx      = hincol[icolx];
    585           s->colx       = presolve_duparray(&colels[mcstrt[icolx]], hincol[icolx]);
    586           s->indx       = presolve_duparray(&hrow[mcstrt[icolx]], hincol[icolx]);
    587585         
    588586          s->ncoly      = hincol[icoly];
    589           s->coly       = presolve_duparray(&colels[mcstrt[icoly]], hincol[icoly]);
    590           s->indy       = presolve_duparray(&hrow[mcstrt[icoly]], hincol[icoly]);
     587          if (s->ncoly<s->ncolx) {
     588            s->colel    = presolve_duparray(colels, hrow, hincol[icoly],
     589                                            mcstrt[icoly]);
     590            s->ncolx=0;
     591          } else {
     592            s->colel    = presolve_duparray(colels, hrow, hincol[icolx],
     593                                            mcstrt[icolx]);
     594            s->ncoly=0;
     595          }
    591596        }
    592597       
     
    827832  const double ztoldj   = prob->ztoldj_;
    828833
     834  // Space for accumulating two columns
     835  int nrows = prob->nrows_;
     836  int * index1 = new int[nrows];
     837  double * element1 = new double[nrows];
     838  memset(element1,0,nrows*sizeof(double));
     839
    829840  for (const action *f = &actions[nactions-1]; actions<=f; f--) {
    830841    int irow = f->row;
     
    882893
    883894
    884 
    885     {
     895    if (f->ncoly) {
     896      // y is shorter so was saved - need to reconstruct x
     897      double multiplier = coeffx/coeffy;
     898      //printf("Current colx %d\n",jcolx);
     899      int * indy = (int *) (f->colel+f->ncoly);
    886900      // find the tail
    887901      CoinBigIndex k=mcstrt[jcolx];
    888       for (int i=0; i<hincol[jcolx]-1; ++i) {
     902      int nX=0;
     903      int i,iRow;
     904      for (i=0; i<hincol[jcolx]-1; ++i) {
     905        if (colels[k]) {
     906          iRow = hrow[k];
     907          index1[nX++]=iRow;
     908          element1[iRow]=colels[k];
     909        }
     910        //printf("%d %d %g\n",i,hrow[k],colels[k]);
    889911        k = link[k];
    890912      }
     913      iRow = hrow[k];
     914      index1[nX++]=iRow;
     915      element1[iRow]=colels[k];
     916      //printf("%d %d %g\n",i,hrow[k],colels[k]);
    891917      // append the old x column to the free list, freeing all links
    892918      link[k] = free_list;
    893919      free_list = mcstrt[jcolx];
    894     }
    895 
    896     // use create_row??
    897     {
    898       int xstart = NO_LINK;
    899       for (int i=0; i<f->ncolx; ++i) {
     920      //printf("Coly is %d row %d coeffx %g coeffy %g\n",
     921      //     jcoly,f->row,f->coeffx,f->coeffy);
     922      int ystart = NO_LINK;
     923      for (i=0; i<f->ncoly; ++i) {
     924        int iRow = indy[i];
     925        double yValue = f->colel[i];
     926        //printf("y %d %d %g\n",i,indy[i],f->colel[i]);
    900927        CoinBigIndex k = free_list;
    901928        free_list = link[free_list];
    902 
     929       
    903930        check_free_list(free_list);
    904 
    905         hrow[k] = f->indx[i];
     931       
     932        hrow[k] = iRow;
    906933        PRESOLVEASSERT(rdone[hrow[k]] || hrow[k] == irow);
    907         colels[k] = f->colx[i];
     934        colels[k] = yValue;
     935        link[k] = ystart;
     936        ystart = k;
     937        yValue *= multiplier;
     938        if (!element1[iRow]) {
     939          element1[iRow]=yValue;
     940          index1[nX++]=iRow;
     941        } else {
     942          element1[iRow]+=yValue;
     943        }
     944      }
     945      mcstrt[jcoly] = ystart;
     946      hincol[jcoly] = f->ncoly;
     947      int xstart = NO_LINK;
     948      int n=0;
     949      for (i=0;i<nX;i++) {
     950        int iRow = index1[i];
     951        double xValue = element1[iRow];
     952        element1[iRow]=0.0;
     953        if (fabs(xValue)>=1.0e-12) {
     954          n++;
     955          CoinBigIndex k = free_list;
     956          free_list = link[free_list];
     957         
     958          check_free_list(free_list);
     959         
     960          hrow[k] = iRow;
     961          PRESOLVEASSERT(rdone[hrow[k]] || hrow[k] == irow);
     962          colels[k] = xValue;
     963          link[k] = xstart;
     964          xstart = k;
     965        }
     966      }
     967      mcstrt[jcolx] = xstart;
     968      assert(n);
     969      hincol[jcolx] = n;
     970     
     971    } else {
     972      // will use x
     973      double multiplier = -coeffy/coeffx;
     974      int * indx = (int *) (f->colel+f->ncolx);
     975      //printf("Current colx %d\n",jcolx);
     976      // find the tail
     977      CoinBigIndex k=mcstrt[jcolx];
     978      int nX=0;
     979      int i,iRow;
     980      for (i=0; i<hincol[jcolx]-1; ++i) {
     981        if (colels[k]) {
     982          iRow = hrow[k];
     983          index1[nX++]=iRow;
     984          element1[iRow]=multiplier*colels[k];
     985        }
     986        //printf("%d %d %g\n",i,hrow[k],colels[k]);
     987        k = link[k];
     988      }
     989      iRow = hrow[k];
     990      index1[nX++]=iRow;
     991      element1[iRow]=multiplier*colels[k];
     992      multiplier = - multiplier;
     993      //printf("%d %d %g\n",i,hrow[k],colels[k]);
     994      // append the old x column to the free list, freeing all links
     995      link[k] = free_list;
     996      free_list = mcstrt[jcolx];
     997      //printf("Coly is %d row %d coeffx %g coeffy %g\n",
     998      //     jcoly,f->row,f->coeffx,f->coeffy);
     999      int xstart = NO_LINK;
     1000      for (i=0; i<f->ncolx; ++i) {
     1001        int iRow = indx[i];
     1002        double xValue = f->colel[i];
     1003        //printf("x %d %d %g\n",i,indx[i],f->colel[i]);
     1004        CoinBigIndex k = free_list;
     1005        free_list = link[free_list];
     1006       
     1007        check_free_list(free_list);
     1008       
     1009        hrow[k] = iRow;
     1010        PRESOLVEASSERT(rdone[hrow[k]] || hrow[k] == irow);
     1011        colels[k] = xValue;
    9081012        link[k] = xstart;
    9091013        xstart = k;
     1014        xValue *= multiplier;
     1015        if (!element1[iRow]) {
     1016          element1[iRow]=xValue;
     1017          index1[nX++]=iRow;
     1018        } else {
     1019          element1[iRow]+=xValue;
     1020        }
    9101021      }
    9111022      mcstrt[jcolx] = xstart;
     1023      hincol[jcolx] = f->ncolx;
     1024      int ystart = NO_LINK;
     1025      int n=0;
     1026      for (i=0;i<nX;i++) {
     1027        int iRow = index1[i];
     1028        double yValue = element1[iRow];
     1029        element1[iRow]=0.0;
     1030        if (fabs(yValue)>=1.0e-12) {
     1031          n++;
     1032          CoinBigIndex k = free_list;
     1033          free_list = link[free_list];
     1034         
     1035          check_free_list(free_list);
     1036         
     1037          hrow[k] = iRow;
     1038          PRESOLVEASSERT(rdone[hrow[k]] || hrow[k] == irow);
     1039          colels[k] = yValue;
     1040          link[k] = ystart;
     1041          ystart = k;
     1042        }
     1043      }
     1044      mcstrt[jcoly] = ystart;
     1045      assert(n);
     1046      hincol[jcoly] = n;
     1047       
    9121048    }
    913     {
    914       int ystart = NO_LINK;
    915       for (int i=0; i<f->ncoly; ++i) {
    916         CoinBigIndex k = free_list;
    917         free_list = link[free_list];
    918 
    919         check_free_list(free_list);
    920 
    921         hrow[k] = f->indy[i];
    922         PRESOLVEASSERT(rdone[hrow[k]] || hrow[k] == irow);
    923         colels[k] = f->coly[i];
    924         link[k] = ystart;
    925         ystart = k;
    926       }
    927       mcstrt[jcoly] = ystart;
    928     }
    929 
    930     hincol[jcolx] = f->ncolx;
    931     hincol[jcoly] = f->ncoly;
    9321049
    9331050
     
    9551072    // therefore, we must increase the row bounds by colels[kcoly] * rhs/coeffy,
    9561073    // which is similar to the bias
     1074    double djy = maxmin * dcost[jcoly];
     1075    coeffy=0.0;
     1076    double djx = maxmin * dcost[jcolx];
     1077    coeffx=0.0;
    9571078    {
    9581079      CoinBigIndex k = mcstrt[jcoly];
    9591080      int ny = hincol[jcoly];
    960       double dj = maxmin * dcost[jcoly];
    9611081      double bounds_factor = rhs/coeffy;
    962 
    9631082      // this probably doesn't work (???)
    9641083           
     
    9811100          acts[row] += coeff * bounds_factor;
    9821101
    983           dj -= rowduals[row] * coeff;
    984         }
    985       }
    986 
    987       // this is the coefficient we need to force col y's reduced cost to 0.0;
    988       // for example, this is obviously true if y is a singleton column
    989       rowduals[irow] = dj / coeffy;
    990       rcosts[jcoly] = 0.0;
    991 
    992       // furthermore, this should leave rcosts[jcolx] unchanged.
     1102          djy -= rowduals[row] * coeff;
     1103          //printf("a %d coeff %g dual %g dj %g\n",
     1104          // row,coeff,rowduals[row],dj);
     1105        } else {
     1106          coeffy=coeff;
     1107        }
     1108      }
     1109      k = mcstrt[jcolx];
     1110      int nx = hincol[jcolx];
     1111           
     1112      for (int i=0; i<nx; ++i) {
     1113        int row = hrow[k];
     1114        double coeff = colels[k];
     1115        k = link[k];
     1116
     1117        if (row != irow) {
     1118          djx -= rowduals[row] * coeff;
     1119          //printf("a %d coeff %g dual %g dj %g\n",
     1120          // row,coeff,rowduals[row],dj);
     1121        } else {
     1122          coeffx=coeff;
     1123        }
     1124      }
     1125
    9931126    }
    9941127
     
    9961129    // was that the variable's bound may have moved, requiring it
    9971130    // to become basic.
     1131    //printf("djs x - %g (%g), y - %g (%g)\n",djx,coeffx,djy,coeffy);
    9981132    if (colstat) {
    9991133      if (prob->columnIsBasic(jcolx) ||
     
    10031137       
    10041138        prob->setColumnStatus(jcoly,PrePostsolveMatrix::basic);
    1005         // col y's reduced costs has already been set to 0.0
     1139        // this is the coefficient we need to force col y's reduced cost to 0.0;
     1140        // for example, this is obviously true if y is a singleton column
     1141        rowduals[irow] = djy / coeffy;
     1142        rcosts[jcolx] = djx - rowduals[irow] * coeffx;
     1143        if (prob->columnIsBasic(jcolx))
     1144          assert (fabs(rcosts[jcolx])<1.0e-5);
     1145        rcosts[jcoly] = 0.0;
    10061146      } else {
    10071147        prob->setColumnStatus(jcolx,PrePostsolveMatrix::basic);
    10081148        prob->setColumnStatusUsingValue(jcoly);
    10091149
    1010         if (rcosts[jcolx] != 0.0) {
    1011           // change rowduals[jcolx] enough to cancel out rcosts[jcolx]
    1012           rowduals[irow] += (rcosts[jcolx] / coeffx);
    1013           rcosts[jcoly] -= ((rcosts[jcolx] / coeffx) * coeffy);
    1014           rcosts[jcolx] = 0.0;
    1015         }
    1016       }
     1150        // change rowduals[jcolx] enough to cancel out rcosts[jcolx]
     1151        rowduals[irow] = djx / coeffx;
     1152        rcosts[jcoly] = djy - rowduals[irow] * coeffy;
     1153        rcosts[jcolx] = 0.0;
     1154      }
     1155    } else {
     1156      // No status array
     1157      // this is the coefficient we need to force col y's reduced cost to 0.0;
     1158      // for example, this is obviously true if y is a singleton column
     1159      rowduals[irow] = djy / coeffy;
     1160      rcosts[jcoly] = 0.0;
    10171161    }
    10181162
     
    10341178          printf("BAD DOUBLE X DJ:  %d %d %g %g\n",
    10351179                 irow, jcolx, rcosts[jcolx], dj);
     1180        rcosts[jcolx]=dj;
    10361181      }
    10371182      {
     
    10461191
    10471192          dj -= rowduals[row] * coeff;
     1193          //printf("b %d coeff %g dual %g dj %g\n",
     1194          // row,coeff,rowduals[row],dj);
    10481195        }
    10491196        if (! (fabs(rcosts[jcoly] - dj) < 100*ZTOLDP))
    10501197          printf("BAD DOUBLE Y DJ:  %d %d %g %g\n",
    10511198                 irow, jcoly, rcosts[jcoly], dj);
     1199        rcosts[jcoly]=dj;
     1200        //exit(0);
    10521201      }
    10531202#endif
     
    10561205    rdone[irow] = DOUBLETON;
    10571206  }
     1207  delete [] index1;
     1208  delete [] element1;
    10581209  prob->free_list_ = free_list;
    10591210}
     
    10631214{
    10641215  for (int i=nactions_-1; i>=0; i--) {
    1065     delete[]actions_[i].colx;
    1066     delete[]actions_[i].indx;
    1067     delete[]actions_[i].coly;
    1068     delete[]actions_[i].indy;
     1216    delete[]actions_[i].colel;
    10691217  }
    10701218  deleteAction(actions_,action*);
  • trunk/PresolveDual.cpp

    r144 r173  
    66#include "PresolveDual.hpp"
    77#include "ClpMessage.hpp"
     8//#define PRESOLVE_TIGHTEN_DUALS 1
    89
    910// this looks for "dominated columns"
     
    198199    // I don't know why I stopped doing this.
    199200#if     PRESOLVE_TIGHTEN_DUALS
     201    const double *rowels        = prob->rowels_;
     202    const int *hcol     = prob->hcol_;
     203    const CoinBigIndex *mrstrt  = prob->mrstrt_;
     204    int *hinrow = prob->hinrow_;
    200205    // tighten row dual bounds, as described on p. 229
    201206    for (int i = 0; i < nrows; i++) {
     
    228233#endif
    229234                rdmax[i] = rmax = bnd;
    230                 tightened = 1;
     235                tightened ++;;
    231236              }
    232237            } else if (coeff < -ZTOLDP && djmax0 > -PRESOLVE_INF && rmax0 > -PRESOLVE_INF) {
     
    237242#endif
    238243                rdmin[i] = rmin = bnd;
    239                 tightened = 1;
     244                tightened ++;;
    240245              }         
    241246            }
     
    248253#endif
    249254                rdmin[i] = rmin = bnd;
    250                 tightened = 1;
     255                tightened ++;;
    251256              }
    252257            } else if (coeff < -ZTOLDP && djmin0 < PRESOLVE_INF && rmin0 < PRESOLVE_INF) {
     
    257262#endif
    258263                rdmax[i] = rmax = bnd;
    259                 tightened = 1;
     264                tightened ++;;
    260265              }
    261266            }
     
    270275#if     PRESOLVE_TIGHTEN_DUALS
    271276    else
    272       printf("DUAL TIGHTENED!  %d\n", nactions);
     277      printf("DUAL TIGHTENED!  %d\n", tightened);
    273278#endif
    274279  }
  • trunk/PresolveDupcol.cpp

    r144 r173  
    66#include "PresolveDupcol.hpp"
    77#include "CoinSort.hpp"
     8#include "PresolveUseless.hpp"
     9#include "ClpMessage.hpp"
    810
    911#define DSEED2 2147483647.0
     
    121123  // by the random value associated with the row.
    122124  int nlook = compute_sums(hincol, mcstrt, hrow, colels, workrow, workcol, sort, ncols);
    123 #if 0
    124   int nlook=0;
    125   for (int i = 0; i < ncols; ++i)
    126     if (hincol[i] > 0 /*1?*/) {
    127       CoinBigIndex kcs = mcstrt[i];
    128       CoinBigIndex kce = kcs + hincol[i];
    129 
    130       double value=0.0;
    131 
    132       // sort all columns for easy comparison in next loop
    133       CoinSort_2(hrow+kcs,hrow+kcs+hincol[i],colels+kcs);
    134       //ekk_sort2(hrow+kcs, colels+kcs, hincol[i]);
    135 
    136       for (CoinBigIndex k=kcs;k<kce;k++) {
    137         int irow=hrow[k];
    138         value += workrow[irow]*colels[k];
    139       }
    140       workcol[nlook]=value;
    141       sort[nlook]=i;
    142       ++nlook;
    143     }
    144 #endif
    145125
    146126#if 1
     
    565545                                               const PresolveAction *next)
    566546{
    567    //  double *colels   = prob->colels_;
    568   int *hrow             = prob->hrow_;
    569547  //  CoinBigIndex *mcstrt              = prob->mcstrt_;
    570548  //  int *hincol               = prob->hincol_;
    571549  int ncols             = prob->ncols_;
    572 
    573   //  double *clo       = prob->clo_;
    574   //  double *cup       = prob->cup_;
    575550
    576551  double *rowels        = prob->rowels_;
     
    583558  double *rup   = prob->rup_;
    584559
    585   //  const char *integerType = prob->integerType_;
    586 
    587   //  double maxmin     = prob->maxmin_;
    588 
    589   action *actions       = new action [nrows];
    590   int nactions = 0;
     560  int nuseless_rows = 0;
    591561
    592562  double * workcol = new double[ncols];
     
    597567
    598568  int nlook = compute_sums(hinrow, mrstrt, hcol, rowels, workcol, workrow, sort, nrows);
    599 #if 0
    600   nlook=0;
    601   for (i = 1; i <= nrow; ++i) {
    602     if (mpre[i] == 0) {
    603       double value=0.0;
    604       CoinBigIndex krs=mrstrt[i];
    605       CoinBigIndex kre=mrstrt[i+1];
    606       CoinSort_2(hcol+krs,hcol+kre,dels+krs);
    607       //ekk_sort2(hcol+krs,dels+krs,kre-krs);
    608       for (k=krs;k<kre;k++) {
    609         int icol=hcol[k];
    610         value += workcol[icol]*dels[k];
    611       }
    612       workrow[nlook]=value;
    613       sort[nlook++]=i;
    614     }
    615   }
    616 #endif
    617569  CoinSort_2(workrow,workrow+nlook,sort);
    618   //ekk_sortonDouble(workrow,sort,nlook);
    619570
    620571  double dval = workrow[0];
     
    625576      CoinBigIndex krs = mrstrt[ithis];
    626577      CoinBigIndex kre = krs + hinrow[ithis];
    627       //      int ishift = mrstrt[ilast];
    628578      if (hinrow[ithis] == hinrow[ilast]) {
    629579        int ishift = mrstrt[ilast] - krs;
    630580        CoinBigIndex k;
    631581        for (k=krs;k<kre;k++) {
    632           if (hcol[k] != hrow[k+ishift] ||
     582          if (hcol[k] != hcol[k+ishift] ||
    633583              rowels[k] != rowels[k+ishift]) {
    634584            break;
     
    642592          double rup2=rup[ithis];
    643593
    644           int idelete=0;
    645           if (rlo1<rlo2) {
     594          int idelete=-1;
     595          if (rlo1<=rlo2) {
    646596            if (rup2<=rup1) {
    647597              /* this is strictly tighter than last */
    648598              idelete=ilast;
     599            } else if (fabs(rlo1-rlo2)<1.0e-12) {
     600              /* last is strictly tighter than this */
     601              idelete=ithis;
     602              // swap so can carry on deleting
     603              sort[jj-1]=ithis;
     604              sort[jj]=ilast;
    649605            } else {
    650606              /* overlapping - could merge */
     
    653609                     rlo1,rup1,rlo2,rup2);
    654610#endif
     611              if (rup1<rlo2) {
     612                // infeasible
     613                prob->status_|= 1;
     614                // wrong message - correct if works
     615                prob->messageHandler()->message(CLP_PRESOLVE_ROWINFEAS,
     616                                                prob->messages())
     617                                                  <<ithis
     618                                                  <<rlo[ithis]
     619                                                  <<rup[ithis]
     620                                                  <<CoinMessageEol;
     621                break;
     622              }
    655623            }
    656624          } else {
     
    659627              /* last is strictly tighter than this */
    660628              idelete=ithis;
     629              // swap so can carry on deleting
     630              sort[jj-1]=ithis;
     631              sort[jj]=ilast;
    661632            } else {
    662633              /* overlapping - could merge */
     
    665636                     rlo1,rup1,rlo2,rup2);
    666637#endif
     638              if (rup2<rlo1) {
     639                // infeasible
     640                prob->status_|= 1;
     641                // wrong message - correct if works
     642                prob->messageHandler()->message(CLP_PRESOLVE_ROWINFEAS,
     643                                                prob->messages())
     644                                                  <<ithis
     645                                                  <<rlo[ithis]
     646                                                  <<rup[ithis]
     647                                                  <<CoinMessageEol;
     648                break;
     649              }
    667650            }
    668651          }
    669           if (idelete) {
    670             {
    671               action *s = &actions[nactions];     
    672               nactions++;
    673 
    674               s->row    = idelete;
    675               s->lbound = rlo[idelete];
    676               s->ubound = rup[idelete];
    677             }
    678             // lazy way - let some other transform eliminate the row
    679             rup[idelete] = PRESOLVE_INF;
    680             rlo[idelete] = PRESOLVE_INF;
    681 
    682             nactions++;
    683           }
     652          if (idelete>=0)
     653            sort[nuseless_rows++]=idelete;
    684654        }
    685655      }
     
    690660  delete[]workrow;
    691661  delete[]workcol;
     662
     663
     664  if (nuseless_rows) {
     665#if     PRESOLVE_SUMMARY
     666    printf("DUPLICATE ROWS:  %d\n", nuseless_rows);
     667#endif
     668    next = useless_constraint_action::presolve(prob,
     669                                               sort, nuseless_rows,
     670                                               next);
     671  }
    692672  delete[]sort;
    693673
    694 
    695   if (nactions) {
    696 #if     PRESOLVE_SUMMARY
    697     printf("DUPLICATE ROWS:  %d\n", nactions);
    698 #endif
    699     next = new duprow_action(nactions, copyOfArray(actions,nactions), next);
    700   }
    701   deleteAction(actions,action*);
    702674  return (next);
    703675}
  • trunk/PresolveEmpty.cpp

    r104 r173  
    2828  double *dcost         = prob->cost_;
    2929
     30  const double ztoldj   = prob->ztoldj_;
    3031  //int *mrstrt         = prob->mrstrt_;
    3132  //int *hinrow         = prob->hinrow_;
     
    5758    // there are no more constraints on this variable,
    5859    // so we had better be able to compute the answer now
    59 
     60    if (fabs(dcost[jcol])<ztoldj)
     61      dcost[jcol]=0.0;
    6062    if (dcost[jcol] * maxmin == 0.0) {
    6163      // hopefully, we can make this non-basic
  • trunk/PresolveImpliedFree.cpp

    r144 r173  
    77#include "PresolveImpliedFree.hpp"
    88#include "ClpMessage.hpp"
    9 
    10 
     9#include "CoinHelperFunctions.hpp"
     10#include "CoinSort.hpp"
    1111
    1212// If there is a row with a singleton column such that no matter what
     
    5959  const CoinBigIndex *mrstrt    = prob->mrstrt_;
    6060  int *hinrow   = prob->hinrow_;
    61   //  const int nrows   = prob->nrows_;
     61  int nrows     = prob->nrows_;
    6262
    6363  /*const*/ double *rlo = prob->rlo_;
     
    7878  int nactions = 0;
    7979
    80   char *implied_free = new char[ncols];
    81   memset(implied_free, 0, ncols*sizeof(char));
    82 
    83   double *ilbound = new double[ncols];
    84   double *iubound = new double[ncols];
    85 
    86   double *tclo = new double[ncols];
    87   double *tcup = new double[ncols];
    88 
    89 #if     PRESOLVE_TRY_TIGHTEN
    90   for (int j=0; j<ncols; j++) {
    91     ilbound[j] = -PRESOLVE_INF;
    92     iubound[j] =  PRESOLVE_INF;
    93     tclo[j] = clo[j];
    94     tcup[j] = cup[j];
     80  int *implied_free = new int[ncols];
     81  int i;
     82  for (i=0;i<ncols;i++)
     83    implied_free[i]=-1;
     84
     85  // memory for min max
     86  int * infiniteDown = new int [nrows];
     87  int * infiniteUp = new int [nrows];
     88  double * maxDown = new double[nrows];
     89  double * maxUp = new double[nrows];
     90
     91  // mark as not computed
     92  // -1 => not computed, -2 give up (singleton), -3 give up (other)
     93  for (i=0;i<nrows;i++) {
     94    if (hinrow[i]>1)
     95      infiniteUp[i]=-1;
     96    else
     97      infiniteUp[i]=-2;
    9598  }
    96 
    97   {
    98     int ntightened;
    99     do {
    100       implied_bounds1(rowels, mrstrt, hrow, hinrow, clo, cup, hcol, ncols,
    101                       rlo, rup, integerType, nrows,
    102                       ilbound, iubound);
    103 
    104       ntightened = 0;
    105       for (int j=0; j<ncols; j++) {
    106         if (tclo[j] < ilbound[j]) {
    107           tclo[j] = ilbound[j];
    108           ntightened++;
    109         }
    110         if (tcup[j] > iubound[j]) {
    111           tcup[j] = iubound[j];
    112           ntightened++;
    113         }
    114       }
    115       printf("NTIGHT: %d\n", ntightened);
    116     } while (ntightened);
    117   }
    118 #endif
     99  double large=1.0e20;
    119100
    120101  int numberLook = prob->numberColsToDo_;
     
    126107    look2 = new int[ncols];
    127108    look=look2;
    128     if (!prob->anyProhibited()) {
    129       for (iLook=0;iLook<ncols;iLook++)
    130         look[iLook]=iLook;
    131       numberLook=ncols;
    132     } else {
    133       // some prohibited
    134       numberLook=0;
    135       for (iLook=0;iLook<ncols;iLook++)
    136         if (!prob->colProhibited(iLook))
    137           look[numberLook++]=iLook;
    138     }
     109    for (iLook=0;iLook<ncols;iLook++)
     110      look[iLook]=iLook;
     111    numberLook=ncols;
    139112 }
    140113
    141114  for (iLook=0;iLook<numberLook;iLook++) {
    142115    int j=look[iLook];
    143     if ((hincol[j] >= 1 && hincol[j] <= 3) &&
    144         !integerType[j]) {
    145       CoinBigIndex kcs = mcstrt[j];
    146       CoinBigIndex kce = kcs + hincol[j];
    147 
    148       for (CoinBigIndex k=kcs; k<kce; ++k) {
     116    if ((hincol[j]  <= 3) && !integerType[j]) {
     117      if (hincol[j]>1) {
     118        // extend to > 3 later
     119        CoinBigIndex kcs = mcstrt[j];
     120        CoinBigIndex kce = kcs + hincol[j];
     121        bool possible = false;
     122        bool singleton = false;
     123        CoinBigIndex k;
     124        double largestElement=0.0;
     125        for (k=kcs; k<kce; ++k) {
     126          int row = hrow[k];
     127          double coeffj = colels[k];
     128         
     129          // if its row is an equality constraint...
     130          if (hinrow[row] > 1 ) {
     131            if ( fabs(rlo[row] - rup[row]) < tol &&
     132                 fabs(coeffj) > ZTOLDP) {
     133              possible=true;
     134            }
     135            largestElement = max(largestElement,fabs(coeffj));
     136          } else {
     137            singleton=true;
     138          }
     139        }
     140        if (possible&&!singleton) {
     141          double low=-COIN_DBL_MAX;
     142          double high=COIN_DBL_MAX;
     143          // get bound implied by all rows
     144          for (k=kcs; k<kce; ++k) {
     145            int row = hrow[k];
     146            double coeffj = colels[k];
     147            if (fabs(coeffj) > ZTOLDP) {
     148              if (infiniteUp[row]==-1) {
     149                // compute
     150                CoinBigIndex krs = mrstrt[row];
     151                CoinBigIndex kre = krs + hinrow[row];
     152                int infiniteUpper = 0;
     153                int infiniteLower = 0;
     154                double maximumUp = 0.0;
     155                double maximumDown = 0.0;
     156                CoinBigIndex kk;
     157                // Compute possible lower and upper ranges
     158                for (kk = krs; kk < kre; ++kk) {
     159                  double value=rowels[kk];
     160                  int iColumn = hcol[kk];
     161                  if (value > 0.0) {
     162                    if (cup[iColumn] >= large) {
     163                      ++infiniteUpper;
     164                    } else {
     165                      maximumUp += cup[iColumn] * value;
     166                    }
     167                    if (clo[iColumn] <= -large) {
     168                      ++infiniteLower;
     169                    } else {
     170                      maximumDown += clo[iColumn] * value;
     171                    }
     172                  } else if (value<0.0) {
     173                    if (cup[iColumn] >= large) {
     174                      ++infiniteLower;
     175                    } else {
     176                      maximumDown += cup[iColumn] * value;
     177                    }
     178                    if (clo[iColumn] <= -large) {
     179                      ++infiniteUpper;
     180                    } else {
     181                      maximumUp += clo[iColumn] * value;
     182                    }
     183                  }
     184                }
     185                double maxUpx = maximumUp+infiniteUpper*1.0e31;
     186                double maxDownx = maximumDown-infiniteLower*1.0e31;
     187                if (maxUpx <= rup[row] + tol &&
     188                    maxDownx >= rlo[row] - tol) {
     189         
     190                  // Row is redundant
     191                  infiniteUp[row]=-3;
     192       
     193                } else if (maxUpx < rlo[row] -tol ) {
     194                  /* there is an upper bound and it can't be reached */
     195                  prob->status_|= 1;
     196                  prob->messageHandler()->message(CLP_PRESOLVE_ROWINFEAS,
     197                                                  prob->messages())
     198                                                    <<row
     199                                                    <<rlo[row]
     200                                                    <<rup[row]
     201                                                    <<CoinMessageEol;
     202                  infiniteUp[row]=-3;
     203                  break;
     204                } else if ( maxDownx > rup[row]+tol) {
     205                  /* there is a lower bound and it can't be reached */
     206                  prob->status_|= 1;
     207                  prob->messageHandler()->message(CLP_PRESOLVE_ROWINFEAS,
     208                                                  prob->messages())
     209                                                    <<row
     210                                                    <<rlo[row]
     211                                                    <<rup[row]
     212                                                    <<CoinMessageEol;
     213                  infiniteUp[row]=-3;
     214                  break;
     215                } else {
     216                  infiniteUp[row]=infiniteUpper;
     217                  infiniteDown[row]=infiniteLower;
     218                  maxUp[row]=maximumUp;
     219                  maxDown[row]=maximumDown;
     220                }
     221              }
     222              if (infiniteUp[row]>=0) {
     223                double lower = rlo[row];
     224                double upper = rup[row];
     225                double value=coeffj;
     226                double nowLower = clo[j];
     227                double nowUpper = cup[j];
     228                double newBound;
     229                int infiniteUpper=infiniteUp[row];
     230                int infiniteLower=infiniteDown[row];
     231                double maximumUp = maxUp[row];
     232                double maximumDown = maxDown[row];
     233                if (value > 0.0) {
     234                  // positive value
     235                  if (lower>-large) {
     236                    if (!infiniteUpper) {
     237                      assert(nowUpper < large);
     238                      newBound = nowUpper +
     239                        (lower - maximumUp) / value;
     240                    } else if (infiniteUpper==1&&nowUpper>large) {
     241                      newBound = (lower -maximumUp) / value;
     242                    } else {
     243                      newBound = -COIN_DBL_MAX;
     244                    }
     245                    if (newBound > nowLower + 1.0e-12) {
     246                      // Tighten the lower bound
     247                      // adjust
     248                      double now;
     249                      if (nowLower<-large) {
     250                        now=0.0;
     251                        infiniteLower--;
     252                      } else {
     253                        now = nowLower;
     254                      }
     255                      maximumDown += (newBound-now) * value;
     256                      nowLower = newBound;
     257                    }
     258                    low=max(low,newBound);
     259                  }
     260                  if (upper <large) {
     261                    if (!infiniteLower) {
     262                      assert(nowLower >- large);
     263                      newBound = nowLower +
     264                        (upper - maximumDown) / value;
     265                    } else if (infiniteLower==1&&nowLower<-large) {
     266                      newBound =   (upper - maximumDown) / value;
     267                    } else {
     268                      newBound = COIN_DBL_MAX;
     269                    }
     270                    if (newBound < nowUpper - 1.0e-12) {
     271                      // Tighten the upper bound
     272                      // adjust
     273                      double now;
     274                      if (nowUpper>large) {
     275                        now=0.0;
     276                        infiniteUpper--;
     277                      } else {
     278                        now = nowUpper;
     279                      }
     280                      maximumUp += (newBound-now) * value;
     281                      nowUpper = newBound;
     282                    }
     283                    high=min(high,newBound);
     284                  }
     285                } else {
     286                  // negative value
     287                  if (lower>-large) {
     288                    if (!infiniteUpper) {
     289                      assert(nowLower >- large);
     290                      newBound = nowLower +
     291                        (lower - maximumUp) / value;
     292                    } else if (infiniteUpper==1&&nowLower<-large) {
     293                      newBound = (lower -maximumUp) / value;
     294                    } else {
     295                      newBound = COIN_DBL_MAX;
     296                    }
     297                    if (newBound < nowUpper - 1.0e-12) {
     298                      // Tighten the upper bound
     299                      // adjust
     300                      double now;
     301                      if (nowUpper>large) {
     302                        now=0.0;
     303                        infiniteLower--;
     304                      } else {
     305                        now = nowUpper;
     306                      }
     307                      maximumDown += (newBound-now) * value;
     308                      nowUpper = newBound;
     309                    }
     310                    high=min(high,newBound);
     311                  }
     312                  if (upper <large) {
     313                    if (!infiniteLower) {
     314                      assert(nowUpper < large);
     315                      newBound = nowUpper +
     316                        (upper - maximumDown) / value;
     317                    } else if (infiniteLower==1&&nowUpper>large) {
     318                      newBound =   (upper - maximumDown) / value;
     319                    } else {
     320                      newBound = -COIN_DBL_MAX;
     321                    }
     322                    if (newBound > nowLower + 1.0e-12) {
     323                      // Tighten the lower bound
     324                      // adjust
     325                      double now;
     326                      if (nowLower<-large) {
     327                        now=0.0;
     328                        infiniteUpper--;
     329                      } else {
     330                        now = nowLower;
     331                      }
     332                      maximumUp += (newBound-now) * value;
     333                      nowLower = newBound;
     334                    }
     335                    low = max(low,newBound);
     336                  }
     337                }
     338              }
     339            }
     340          }
     341          if (clo[j] <= low && high <= cup[j]) {
     342             
     343            // both column bounds implied by the constraints of the problem
     344            // get row
     345            largestElement *= 0.1;
     346            int krow=-1;
     347            int ninrow=ncols+1;
     348            for (k=kcs; k<kce; ++k) {
     349              int row = hrow[k];
     350              double coeffj = colels[k];
     351              if ( fabs(rlo[row] - rup[row]) < tol &&
     352                   fabs(coeffj) > largestElement) {
     353                if (hinrow[row]<ninrow) {
     354                  ninrow=hinrow[row];
     355                  krow=row;
     356                }
     357              }
     358            }
     359            if (krow>=0) {
     360              implied_free[j] = krow;
     361              //printf("column %d implied free by row %d hincol %d hinrow %d\n",
     362              //     j,krow,hincol[j],hinrow[krow]);
     363            }
     364          }
     365        }
     366      } else if (hincol[j]) {
     367        // singleton column
     368        CoinBigIndex k = mcstrt[j];
    149369        int row = hrow[k];
    150370        double coeffj = colels[k];
    151 
    152         // if its row is an equality constraint...
    153         if (hinrow[row] > 1 &&  // don't bother with singleton rows
    154 
    155             // either this is a singleton col,
    156             // (for momement - only if no cost)
    157             // or this particular row is an equality
    158             ((hincol[j] == 1 && !cost[j]) || fabs(rlo[row] - rup[row]) < tol) &&
    159 
     371        if ((!cost[j]||rlo[row]==rup[row])&&hinrow[row]>1&&
    160372            fabs(coeffj) > ZTOLDP) {
    161 
     373           
    162374          CoinBigIndex krs = mrstrt[row];
    163375          CoinBigIndex kre = krs + hinrow[row];
    164 
     376         
    165377          double maxup, maxdown, ilow, iup;
    166378          implied_bounds(rowels, clo, cup, hcol,
     
    168380                         &maxup, &maxdown,
    169381                         j, rlo[row], rup[row], &ilow, &iup);
    170 
    171 #if     PRESOLVE_TRY_TIGHTEN
    172           if ((clo[j] <= ilbound[j] && iubound[j] <= cup[j]) &&
    173               ! (clo[j] <= ilow && iup <= cup[j]))
    174             printf("TIGHTER:  %6d %6d\n", row, j);
    175 #endif
    176 
     382         
     383         
    177384          if (maxup < PRESOLVE_INF && maxup + tol < rlo[row]) {
    178385            /* there is an upper bound and it can't be reached */
    179386            prob->status_|= 1;
    180387            prob->messageHandler()->message(CLP_PRESOLVE_ROWINFEAS,
    181                                              prob->messages())
    182                                                <<row
    183                                                <<rlo[row]
    184                                                <<rup[row]
    185                                                <<CoinMessageEol;
     388                                            prob->messages())
     389                                              <<row
     390                                              <<rlo[row]
     391                                              <<rup[row]
     392                                              <<CoinMessageEol;
    186393            break;
    187394          } else if (-PRESOLVE_INF < maxdown && rup[row] < maxdown - tol) {
     
    189396            prob->status_|= 1;
    190397            prob->messageHandler()->message(CLP_PRESOLVE_ROWINFEAS,
    191                                              prob->messages())
    192                                                <<row
    193                                                <<rlo[row]
    194                                                <<rup[row]
    195                                                <<CoinMessageEol;
     398                                            prob->messages())
     399                                              <<row
     400                                              <<rlo[row]
     401                                              <<rup[row]
     402                                              <<CoinMessageEol;
    196403            break;
    197404          } else if (clo[j] <= ilow && iup <= cup[j]) {
    198 
     405           
    199406            // both column bounds implied by the constraints of the problem
    200             implied_free[j] = hincol[j];
    201             break;
     407            implied_free[j] = row;
     408            //printf("column %d implied free by row %d hincol %d hinrow %d\n",
     409            //   j,row,hincol[j],hinrow[row]);
    202410          }
    203411        }
     
    207415  // implied_free[j] == hincol[j] && hincol[j] > 0 ==> j is implied free
    208416
    209 #if 0
    210   // DEBUG
    211   static int nfree = 0;
    212   static int maxfree = atoi(getenv("MAXFREE"));
    213 #endif
     417  delete [] infiniteDown;
     418  delete [] infiniteUp;
     419  delete [] maxDown;
     420  delete [] maxUp;
    214421
    215422  int isolated_row = -1;
     
    221428  for (iLook=0;iLook<numberLook;iLook++) {
    222429    int j=look[iLook];
    223     if (hincol[j] == 1 && implied_free[j] == 1) {
     430    if (hincol[j] == 1 && implied_free[j] >=0) {
    224431      CoinBigIndex kcs = mcstrt[j];
    225432      int row = hrow[kcs];
     
    229436      CoinBigIndex kre = krs + hinrow[row];
    230437
    231 #if 0
    232       if (nfree >= maxfree)
    233         continue;
    234       nfree++;
    235 #endif
    236438
    237439      // isolated rows are weird
     
    321523      hincol[j] = 0;
    322524
    323       implied_free[j] = 0;      // probably unnecessary
     525      implied_free[j] = -1;     // probably unnecessary
    324526    }
    325527  }
     
    334536  deleteAction(actions,action*);
    335537
    336   delete[]ilbound;
    337   delete[]iubound;
    338   delete[]tclo;
    339   delete[]tcup;
    340 
    341538  if (isolated_row != -1) {
    342539    const PresolveAction *nextX = isolated_constraint_action::presolve(prob,
     
    345542      next = nextX; // may fail
    346543  }
    347 
    348544  // try more complex ones
    349   if (fill_level)
     545  if (fill_level) {
    350546    next = subst_constraint_action::presolve(prob, implied_free, next,fill_level);
     547  }
    351548  delete[]implied_free;
    352549
     
    366563  const int nactions = nactions_;
    367564
    368   double *colels        = prob->colels_;
     565  double *elementByColumn       = prob->colels_;
    369566  int *hrow             = prob->hrow_;
    370   CoinBigIndex *mcstrt          = prob->mcstrt_;
    371   int *hincol           = prob->hincol_;
     567  CoinBigIndex *columnStart             = prob->mcstrt_;
     568  int *numberInColumn           = prob->hincol_;
    372569  int *link             = prob->link_;
    373570
     
    421618          check_free_list(free_list);
    422619
    423           link[kk] = mcstrt[jcol];
    424           mcstrt[jcol] = kk;
    425           colels[kk] = coeff;
     620          link[kk] = columnStart[jcol];
     621          columnStart[jcol] = kk;
     622          elementByColumn[kk] = coeff;
    426623          hrow[kk] = irow;
    427624        }
     
    429626        if (jcol == icol) {
    430627          // initialize the singleton column
    431           hincol[jcol] = 1;
     628          numberInColumn[jcol] = 1;
    432629          clo[icol] = f->clo;
    433630          cup[icol] = f->cup;
     
    435632          cdone[icol] = IMPLIED_FREE;
    436633        } else {
    437           hincol[jcol]++;
     634          numberInColumn[jcol]++;
    438635        }
    439636      }
     
    456653        else {
    457654          int jcol = rowcols[k];
    458           PRESOLVE_STMT(CoinBigIndex kk = presolve_find_row2(irow, mcstrt[jcol], hincol[jcol], hrow, link));
     655          PRESOLVE_STMT(CoinBigIndex kk = presolve_find_row2(irow, columnStart[jcol], numberInColumn[jcol], hrow, link));
    459656          act += rowels[k] * sol[jcol];
    460657        }
  • trunk/PresolveMatrix.cpp

    r144 r173  
    7676}
    7777
     78// This one saves in one go to save [] memory
     79double * presolve_duparray(const double * element, const int * index,
     80                           int length, int offset)
     81{
     82  assert (sizeof(double)==2*sizeof(int));
     83  int n = (3*length+1)>>1;
     84  double * dArray = new double [n];
     85  int * iArray = (int *) (dArray+length);
     86  memcpy(dArray,element+offset,length*sizeof(double));
     87  memcpy(iArray,index+offset,length*sizeof(int));
     88  return dArray;
     89}
    7890
    7991
  • trunk/PresolveSubst.cpp

    r172 r173  
    318318}
    319319
    320 
    321 
    322320// add -x/y times row y to row x, thus cancelling out one column of rowx;
    323321// afterwards, that col will be singleton for rowy, so we drop the row.
     
    328326// This implements the functionality of ekkrdc3.
    329327const PresolveAction *subst_constraint_action::presolve(PresolveMatrix *prob,
    330                                          char *implied_free,
     328                                         int *implied_free,
    331329                                        const PresolveAction *next,
    332330                                                        int &try_fill_level)
     
    338336  const int ncols       = prob->ncols_;
    339337
    340   const double *clo     = prob->clo_;
    341   const double *cup     = prob->cup_;
    342 
    343338  double *rowels        = prob->rowels_;
    344339  int *hcol     = prob->hcol_;
     
    391386
    392387  // DEBUGGING
    393   //  int nt = 0;
    394   int ngood = 0;
    395388  int nsubst = 0;
    396389#ifdef  DEBUG_PRESOLVEx
     
    430423  for (iLook=0;iLook<numberLook;iLook++) {
    431424    int jcoly=look[iLook];
     425    bool looksGood=false;
    432426    if (hincol[jcoly] > 1 && hincol[jcoly] <= fill_level &&
    433         implied_free[jcoly] == hincol[jcoly]) {
     427        implied_free[jcoly] >=0) {
     428      looksGood=true;
    434429      CoinBigIndex kcs = mcstrt[jcoly];
    435430      CoinBigIndex kce = kcs + hincol[jcoly];
    436 
     431     
    437432      int bestrowy_size = 0;
    438433      int bestrowy_row=-1;
     
    443438        int row = hrow[k];
    444439        double coeffj = colels[k];
    445 
     440       
    446441        // we don't clean up zeros in the middle of the routine.
    447442        // if there is one, skip this candidate.
     
    450445          break;
    451446        }
    452 
    453         // if its row is an equality constraint...
    454         if (hinrow[row] > 1 &&  // don't bother with singleton rows
    455 
    456             fabs(rlo[row] - rup[row]) < tol) {
    457           CoinBigIndex krs = mrstrt[row];
    458           CoinBigIndex kre = krs + hinrow[row];
    459 
    460           double maxup, maxdown, ilow, iup;
    461 
    462           implied_bounds(rowels, clo, cup, hcol,
    463                          krs, kre,
    464                          &maxup, &maxdown,
    465                          jcoly, rlo[row], rup[row], &ilow, &iup);
    466 
    467           if (maxup < PRESOLVE_INF && maxup + tol < rlo[row]) {
    468             prob->status_|= 1;
    469             prob->messageHandler()->message(CLP_PRESOLVE_ROWINFEAS,
    470                                              prob->messages())
    471                                                <<row
    472                                                <<rlo[row]
    473                                                <<rup[row]
    474                                                <<CoinMessageEol;
    475             break;
    476           } else if (-PRESOLVE_INF < maxdown && rup[row] < maxdown - tol) {
    477             prob->status_|= 1;
    478             prob->messageHandler()->message(CLP_PRESOLVE_ROWINFEAS,
    479                                              prob->messages())
    480                                                <<row
    481                                                <<rlo[row]
    482                                                <<rup[row]
    483                                                <<CoinMessageEol;
    484             break;
    485           } else {
    486             // the row has an implied upper or a lower bound
    487 
    488             if (clo[jcoly] <= ilow && iup <= cup[jcoly]) {
    489               // both column bounds implied by the constraint bounds
    490 
    491               // we want coeffy to be smaller than x, BACKWARDS from in doubleton
    492               if (bestrowy_size == 0 ||
    493                   fabs(coeffj) > fabs(bestrowy_coeff) ||
    494                   (fabs(coeffj) == fabs(bestrowy_coeff) &&
    495                    hinrow[row] < bestrowy_size)) {
    496                 bestrowy_size = hinrow[row];
    497                 bestrowy_row = row;
    498                 bestrowy_coeff = coeffj;
    499                 bestrowy_k = k;
    500               }
    501             }
     447         
     448        if (row==implied_free[jcoly]) {
     449          // if its row is an equality constraint...
     450          if (hinrow[row] > 1 &&        // don't bother with singleton rows
     451             
     452              fabs(rlo[row] - rup[row]) < tol) {
     453            // both column bounds implied by the constraint bounds
     454           
     455            // we want coeffy to be smaller than x, BACKWARDS from in doubleton
     456            bestrowy_size = hinrow[row];
     457            bestrowy_row = row;
     458            bestrowy_coeff = coeffj;
     459            bestrowy_k = k;
    502460          }
    503461        }
     
    510468      for (k=kcs; k<kce; ++k) {
    511469        double coeff_factor = fabs(colels[k] / bestrowy_coeff);
    512         if (fabs(coeff_factor) > 1.3)
     470        if (fabs(coeff_factor) > 10.0)
    513471          all_ok = false;
    514472      }
    515 
     473#if 0
    516474      // check fill-in
    517475      if (all_ok && hincol[jcoly] == 3) {
     
    585543#endif
    586544      }
    587 
     545#endif
    588546      // probably never happens
    589547      if (all_ok && nzerocols + hinrow[bestrowy_row] >= ncols)
     
    911869#endif
    912870      }
     871     
    913872    }
    914873  }
     
    917876  if (nactions < 30&&fill_level==2)
    918877    try_fill_level = -3;
    919 
    920878  if (nactions) {
    921879#if     PRESOLVE_SUMMARY
     
    935893  return (next);
    936894}
    937 
    938895
    939896void subst_constraint_action::postsolve(PostsolveMatrix *prob) const
     
    12021159            }
    12031160#if     DEBUG_PRESOLVE
    1204             PRESOLVEASSERT(rlo[jrowx] - ztolzb <= actx && actx <= rup[jrowx] + ztolzb);
     1161            PRESOLVEASSERT(rlo[jrowx] - prob->ztolzb_ <= actx
     1162                           && actx <= rup[jrowx] + prob->ztolzb_);
    12051163#endif
    12061164            acts[jrowx] = actx;
  • trunk/include/Presolve.hpp

    r144 r173  
    6969      original.  Bounds will be moved by up to feasibilityTolerance
    7070      to try and stay feasible.
     71      Names will be dropped in presolved model if asked
    7172  */
    7273  virtual ClpSimplex * presolvedModel(ClpSimplex & si,
    7374                                      double feasibilityTolerance=0.0,
    7475                                      bool keepIntegers=true,
    75                                       int numberPasses=5);
     76                                      int numberPasses=5,
     77                                      bool dropNames=false);
    7678
    7779  /** Return pointer to presolved model,
  • trunk/include/PresolveDoubleton.hpp

    r50 r173  
    2929
    3030    int ncolx;
    31     double *colx;
    32     int *indx;
     31    double *colel;
    3332
    3433    int ncoly;
    35     double *coly;
    36     int *indy;
    3734  };
    3835
  • trunk/include/PresolveMatrix.hpp

    r144 r173  
    2828char *presolve_duparray(const char *d, int n, int n2);
    2929char *presolve_duparray(const char *d, int n);
     30// This one saves in one go to save [] memory
     31double * presolve_duparray(const double * element, const int * index,
     32                           int length, int offset);
    3033
    3134void presolve_delete_from_row(int row, int col /* thing to delete */,
  • trunk/include/PresolveSubst.hpp

    r50 r173  
    88
    99class subst_constraint_action : public PresolveAction {
     10public:
    1011  struct action {
    1112    int col;
     
    4041
    4142  static const PresolveAction *presolve(PresolveMatrix * prob,
    42                                          char *implied_free,
     43                                         int *implied_free,
    4344                                         const PresolveAction *next,
    4445                                        int & fill_level);
     46  static const PresolveAction *presolveX(PresolveMatrix * prob,
     47                                  const PresolveAction *next,
     48                                  int fillLevel);
    4549
    4650  void postsolve(PostsolveMatrix *prob) const;
Note: See TracChangeset for help on using the changeset viewer.