Changeset 1212


Ignore:
Timestamp:
Aug 21, 2009 12:19:13 PM (10 years ago)
Author:
forrest
Message:

fixes

Location:
trunk/Cbc/src
Files:
15 edited

Legend:

Unmodified
Added
Removed
  • trunk/Cbc/src/CbcBranchDynamic.cpp

    r1205 r1212  
    670670#define WEIGHT_AFTER 0.8
    671671#define WEIGHT_BEFORE 0.1
     672  //Stolen from Constraint Integer Programming book
    672673#define WEIGHT_PRODUCT
    673674  if (fabs(value-nearest)<=integerTolerance) {
  • trunk/Cbc/src/CbcCutGenerator.cpp

    r1211 r1212  
    3838    numberTimes_(0),
    3939    numberCuts_(0),
     40    numberElements_(0),
    4041    numberColumnCuts_(0),
    4142    numberCutsActive_(0),
     
    6061    numberTimes_(0),
    6162    numberCuts_(0),
     63    numberElements_(0),
    6264    numberColumnCuts_(0),
    6365    numberCutsActive_(0),
     
    108110  numberTimes_ = rhs.numberTimes_;
    109111  numberCuts_ = rhs.numberCuts_;
     112  numberElements_ = rhs.numberElements_;
    110113  numberColumnCuts_ = rhs.numberColumnCuts_;
    111114  numberCutsActive_ = rhs.numberCutsActive_;
     
    137140    numberTimes_ = rhs.numberTimes_;
    138141    numberCuts_ = rhs.numberCuts_;
     142    numberElements_ = rhs.numberElements_;
    139143    numberColumnCuts_ = rhs.numberColumnCuts_;
    140144    numberCutsActive_ = rhs.numberCutsActive_;
     
    354358        }
    355359      }
    356       if (doCuts) {
     360      if (doCuts&&generator->tightLower()) {
    357361        // probing may have tightened bounds - check
    358362        const double * tightLower = generator->tightLower();
     
    541545    {
    542546      int numberRowCutsAfter = cs.sizeRowCuts() ;
    543       if (numberRowCutsBefore<numberRowCutsAfter) {
    544 #if 0
    545         printf("generator %s generated %d row cuts\n",
    546                generatorName_,numberRowCutsAfter-numberRowCutsBefore);
    547 #endif
    548         numberCuts_ += numberRowCutsAfter-numberRowCutsBefore;
    549       }
    550       int numberColumnCutsAfter = cs.sizeColCuts() ;
    551       if (numberColumnCutsBefore<numberColumnCutsAfter) {
    552 #if 0
    553         printf("generator %s generated %d column cuts\n",
    554                generatorName_,numberColumnCutsAfter-numberColumnCutsBefore);
    555 #endif
    556         numberColumnCuts_ += numberColumnCutsAfter-numberColumnCutsBefore;
    557       }
    558     }
    559     {
    560       int numberRowCutsAfter = cs.sizeRowCuts() ;
    561547      int k ;
    562548      int nEls=0;
     
    650636        //printf("need to remove cuts\n");
    651637        // just add most effective
     638#if 1
    652639        int nDelete = nEls - nReasonable;
    653640       
     
    657644        // For parallel cuts
    658645        double * element2 = new double [numberColumns];
    659         //#define USE_OBJECTIVE
     646        //#define USE_OBJECTIVE 2
    660647#ifdef USE_OBJECTIVE
    661648        const double *objective = solver->getObjCoefficients() ;
     649#if USE_OBJECTIVE>1
     650        double objNorm=0.0;
     651        for (int i=0;i<numberColumns;i++)
     652          objNorm += objective[i]*objective[i];
     653        if (objNorm)
     654          objNorm = 1.0/sqrt(objNorm);
     655        else
     656          objNorm=1.0;
     657        objNorm *= 0.01; // downgrade
     658#endif
    662659#endif
    663660        CoinZeroN(element2,numberColumns);
     
    700697            double norm=1.0e-3;
    701698#ifdef USE_OBJECTIVE
    702             double normObj=0.0;
     699            double obj=0.0;
    703700#endif
    704701            for (int i=0;i<n;i++) {
     
    708705              norm += value*value;
    709706#ifdef USE_OBJECTIVE
    710               normObj += value*objective[iColumn];
     707              obj += value*objective[iColumn];
    711708#endif
    712709            }
     
    721718#ifdef USE_OBJECTIVE
    722719            if (sum) {
    723               normObj = CoinMax(1.0e-6,fabs(normObj));
    724               norm=sqrt(normObj*norm);
    725               //sum += fabs(normObj)*invObjNorm;
     720#if USE_OBJECTIVE==1
     721              obj = CoinMax(1.0e-6,fabs(obj));
     722              norm=sqrt(obj*norm);
     723              //sum += fabs(obj)*invObjNorm;
    726724              //printf("sum %g norm %g normobj %g invNorm %g mod %g\n",
    727               //     sum,norm,normObj,invObjNorm,normObj*invObjNorm);
    728             }
    729 #endif
     725              //     sum,norm,obj,invObjNorm,obj*invObjNorm);
     726              // normalize
     727              sum /= sqrt(norm);
     728#else
     729              // normalize
     730              norm = 1.0/sqrt(norm);
     731              sum = (sum+objNorm*obj)*norm;
     732#endif
     733            }
     734#else
    730735            // normalize
    731736            sum /= sqrt(norm);
     737#endif
    732738            //sum /= pow(norm,0.3);
    733739            // adjust for length
     
    836842        delete [] which;
    837843        numberRowCutsAfter = cs.sizeRowCuts() ;
     844#else
     845        double * norm = new double [nCuts];
     846        int * which = new int [2*nCuts];
     847        double * score = new double [nCuts];
     848        double * ortho = new double [nCuts];
     849        int nIn=0;
     850        int nOut=nCuts;
     851        // For parallel cuts
     852        double * element2 = new double [numberColumns];
     853        const double *objective = solver->getObjCoefficients() ;
     854        double objNorm=0.0;
     855        for (int i=0;i<numberColumns;i++)
     856          objNorm += objective[i]*objective[i];
     857        if (objNorm)
     858          objNorm = 1.0/sqrt(objNorm);
     859        else
     860          objNorm=1.0;
     861        objNorm *= 0.1; // weight of 0.1
     862        CoinZeroN(element2,numberColumns);
     863        int numberRowCuts = numberRowCutsAfter-numberRowCutsBefore;
     864        int iBest=-1;
     865        double best=0.0;
     866        int nPossible=0;
     867        double testValue = (depth>1) ? 0.7 : 0.5;
     868        for (k = 0;k<numberRowCuts;k++) {
     869          const OsiRowCut * thisCut = cs.rowCutPtr(k+numberRowCutsBefore) ;
     870          double sum=0.0;
     871          if (thisCut->lb()<=thisCut->ub()) {
     872            int n=thisCut->row().getNumElements();
     873            const int * column = thisCut->row().getIndices();
     874            const double * element = thisCut->row().getElements();
     875            assert (n);
     876            double normThis=1.0e-6;
     877            double obj=0.0;
     878            for (int i=0;i<n;i++) {
     879              int iColumn = column[i];
     880              double value = element[i];
     881              sum += value*solution[iColumn];
     882              normThis += value*value;
     883              obj += value*objective[iColumn];
     884            }
     885            if (sum>thisCut->ub()) {
     886              sum= sum-thisCut->ub();
     887            } else if (sum<thisCut->lb()) {
     888              sum= thisCut->lb()-sum;
     889            } else {
     890              sum=0.0;
     891            }
     892            if (sum) {
     893              normThis =1.0/sqrt(normThis);
     894              norm[k]=normThis;
     895              sum *= normThis;
     896              obj *= normThis;
     897              score[k]=sum+obj*objNorm;
     898              ortho[k]=1.0;
     899            }
     900          } else {
     901            // keep and discard others
     902            nIn=1;
     903            which[0]=k;
     904            for (int j = 0;j<numberRowCuts;j++) {
     905              if (j!=k)
     906                which[nOut++]=j;
     907            }
     908            iBest=-1;
     909            break;
     910          }
     911          if (sum) {
     912            if (score[k]>best) {
     913              best=score[k];
     914              iBest=nPossible;
     915            }
     916            which[nPossible++]=k;
     917          } else {
     918            which[nOut++]=k;
     919          }
     920        }
     921        while (iBest>=0) {
     922          int kBest=which[iBest];
     923          int j=which[nIn];
     924          which[iBest]=j;
     925          which[nIn++]=kBest;
     926          const OsiRowCut * thisCut = cs.rowCutPtr(kBest+numberRowCutsBefore) ;
     927          int n=thisCut->row().getNumElements();
     928          nReasonable -= n;
     929          if (nReasonable<=0) {
     930            for (k=nIn;k<nPossible;k++)
     931              which[nOut++]=which[k];
     932            break;
     933          }
     934          // Now see which ones are too similar and choose next
     935          iBest=-1;
     936          best=0.0;
     937          int nOld=nPossible;
     938          nPossible=nIn;
     939          const int * column = thisCut->row().getIndices();
     940          const double * element = thisCut->row().getElements();
     941          assert (n);
     942          double normNew = norm[kBest];
     943          for (int i=0;i<n;i++) {
     944            double value = element[i];
     945            element2[column[i]]=value;
     946          }
     947          for (int j=nIn;j<nOld;j++) {
     948            k = which[j];
     949            const OsiRowCut * thisCut2 = cs.rowCutPtr(k+numberRowCutsBefore) ;
     950            int nB=thisCut2->row().getNumElements();
     951            const int * columnB = thisCut2->row().getIndices();
     952            const double * elementB = thisCut2->row().getElements();
     953            assert (nB);
     954            double normB=norm[k];
     955            double product=0.0;
     956            for (int i=0;i<nB;i++) {
     957              double value = elementB[i];
     958              product += value*element2[columnB[i]];
     959            }
     960            double orthoScore = 1.0-product*normNew*normB;
     961            if (orthoScore>=testValue) {
     962              ortho[k]=CoinMin(orthoScore,ortho[k]);
     963              double test=score[k]+ortho[k];
     964              if (test>best) {
     965                best=score[k];
     966                iBest=nPossible;
     967              }
     968              which[nPossible++]=k;
     969            } else {
     970              which[nOut++]=k;
     971            }
     972          }
     973          for (int i=0;i<n;i++) {
     974            element2[column[i]]=0.0;
     975          }
     976        }
     977        delete [] score;
     978        delete [] ortho;
     979        std::sort(which+nCuts,which+nOut);
     980        k=nOut-1;
     981        for (;k>=nCuts;k--) {
     982          cs.eraseRowCut(which[k]+numberRowCutsBefore);
     983        }
     984        delete [] norm;
     985        delete [] which;
     986        numberRowCutsAfter = cs.sizeRowCuts() ;
     987#endif
    838988      }
    839989    }
     
    8701020    }
    8711021#endif
     1022    {
     1023      int numberRowCutsAfter = cs.sizeRowCuts() ;
     1024      if (numberRowCutsBefore<numberRowCutsAfter) {
     1025        for (int k = numberRowCutsBefore;k<numberRowCutsAfter;k++) {
     1026          OsiRowCut thisCut = cs.rowCut(k) ;
     1027          int n=thisCut.row().getNumElements();
     1028          numberElements_ += n;
     1029        }
     1030#if 0
     1031        printf("generator %s generated %d row cuts\n",
     1032               generatorName_,numberRowCutsAfter-numberRowCutsBefore);
     1033#endif
     1034        numberCuts_ += numberRowCutsAfter-numberRowCutsBefore;
     1035      }
     1036      int numberColumnCutsAfter = cs.sizeColCuts() ;
     1037      if (numberColumnCutsBefore<numberColumnCutsAfter) {
     1038#if 0
     1039        printf("generator %s generated %d column cuts\n",
     1040               generatorName_,numberColumnCutsAfter-numberColumnCutsBefore);
     1041#endif
     1042        numberColumnCuts_ += numberColumnCutsAfter-numberColumnCutsBefore;
     1043      }
     1044    }
    8721045    if (timing())
    8731046      timeInCutGenerator_ += CoinCpuTime()-time1;
  • trunk/Cbc/src/CbcCutGenerator.hpp

    r1211 r1212  
    204204  inline void incrementNumberCutsInTotal(int value=1)
    205205  { numberCuts_ += value;}
     206  /// Total number of elements added
     207  inline int numberElementsInTotal() const
     208  { return numberElements_;}
     209  inline void setNumberElementsInTotal(int value)
     210  { numberElements_ = value;}
     211  inline void incrementNumberElementsInTotal(int value=1)
     212  { numberElements_ += value;}
    206213  /// Total number of column cuts
    207214  inline int numberColumnCuts() const
     
    246253  inline void setIneffectualCuts(bool yesNo)
    247254  { switches_&=~512;switches_ |= yesNo ? 512 : 0;}
     255  /// Whether to use if nay cuts generated
     256  inline bool whetherToUse() const
     257  { return (switches_&1024)!=0;}
     258  /// Set whether to use if any cuts generated
     259  inline void setWhetherToUse(bool yesNo)
     260  { switches_&=~1024;switches_ |= yesNo ? 1024 : 0;}
    248261  /// Number of cuts generated at root
    249262  inline int numberCutsAtRoot() const
     
    322335  /// Total number of cuts added
    323336  int numberCuts_;
     337  /// Total number of elements added
     338  int numberElements_;
    324339  /// Total number of column cuts added
    325340  int numberColumnCuts_;
  • trunk/Cbc/src/CbcHeuristic.cpp

    r1211 r1212  
    11821182  model_->setSpecialOptions(saveModelOptions);
    11831183  model_->setLogLevel(logLevel);
     1184  if (returnCode==1||returnCode==2) {
     1185    OsiSolverInterface * solverC = model_->continuousSolver();
     1186    if (false&&solverC) {
     1187      const double * lower = solver->getColLower();
     1188      const double * upper = solver->getColUpper();
     1189      const double * lowerC = solverC->getColLower();
     1190      const double * upperC = solverC->getColUpper();
     1191      bool good=true;
     1192      for (int iColumn=0;iColumn<numberColumns;iColumn++) {
     1193        if (solverC->isInteger(iColumn)) {
     1194          if (lower[iColumn]>lowerC[iColumn]&&
     1195              upper[iColumn]<upperC[iColumn]) {
     1196            good=false;
     1197            printf("CUT - can't add\n");
     1198            break;
     1199          }
     1200        }
     1201      }
     1202      if (good) {
     1203        double * cut = new double [numberColumns];
     1204        int * which = new int [numberColumns];
     1205        double rhs=-1.0;
     1206        int n=0;
     1207        for (int iColumn=0;iColumn<numberColumns;iColumn++) {
     1208          if (solverC->isInteger(iColumn)) {
     1209            if (lower[iColumn]==upperC[iColumn]) {
     1210              rhs += lower[iColumn];
     1211              cut[n]=1.0;
     1212              which[n++]=iColumn;
     1213            } else if (upper[iColumn]==lowerC[iColumn]) {
     1214              rhs -= upper[iColumn];
     1215              cut[n]=-1.0;
     1216              which[n++]=iColumn;
     1217            }
     1218          }
     1219        }
     1220        printf("CUT has %d entries\n",n);
     1221        OsiRowCut newCut;
     1222        newCut.setLb(-COIN_DBL_MAX);
     1223        newCut.setUb(rhs);
     1224        newCut.setRow(n,which,cut,false);
     1225        model_->makeGlobalCut(newCut);
     1226        delete [] cut;
     1227        delete [] which;
     1228      }
     1229    }
     1230#ifdef COIN_DEVELOP
     1231    if (status==1)
     1232      printf("heuristic could add cut because infeasible (%s)\n",heuristicName_.c_str());
     1233    else if (status==2)
     1234      printf("heuristic could add cut because optimal (%s)\n",heuristicName_.c_str());
     1235#endif
     1236  }
    11841237  if (reset) {
    11851238    for (int iColumn=0;iColumn<numberColumns;iColumn++) {
     
    11911244#ifdef HISTORY_STATISTICS
    11921245  getHistoryStatistics_=true;
    1193 #endif
    1194 #ifdef COIN_DEVELOP
    1195   if (returnCode==1||returnCode==2) {
    1196     if (status==1)
    1197       printf("heuristic could add cut because infeasible (%s)\n",heuristicName_.c_str());
    1198     else if (status==2)
    1199       printf("heuristic could add cut because optimal (%s)\n",heuristicName_.c_str());
    1200   }
    12011246#endif
    12021247  solver->setHintParam(OsiDoReducePrint,takeHint,strength);
     
    14641509  // Get a copy of original matrix (and by row for rounding);
    14651510  assert(model.solver());
    1466   matrix_ = *model.solver()->getMatrixByCol();
    1467   matrixByRow_ = *model.solver()->getMatrixByRow();
    1468   validate();
     1511  if (model.solver()->getNumRows()) {
     1512    matrix_ = *model.solver()->getMatrixByCol();
     1513    matrixByRow_ = *model.solver()->getMatrixByRow();
     1514    validate();
     1515  }
    14691516  seed_=1;
    14701517}
     
    23122359  // Get a copy of original matrix (and by row for rounding);
    23132360  assert(model_->solver());
    2314   matrix_ = *model_->solver()->getMatrixByCol();
    2315   matrixByRow_ = *model_->solver()->getMatrixByRow();
    2316   // make sure model okay for heuristic
    2317   validate();
     2361  if (model_->solver()->getNumRows()) {
     2362    matrix_ = *model_->solver()->getMatrixByCol();
     2363    matrixByRow_ = *model_->solver()->getMatrixByRow();
     2364    // make sure model okay for heuristic
     2365    validate();
     2366  }
    23182367}
    23192368// Validate model i.e. sets when_ to 0 if necessary (may be NULL)
  • trunk/Cbc/src/CbcHeuristicFPump.cpp

    r1211 r1212  
    350350  CoinWarmStartBasis bestBasis;
    351351  bool exitAll=false;
    352   double saveBestObjective = model_->getMinimizationObjValue();
     352  //double saveBestObjective = model_->getMinimizationObjValue();
    353353  int numberSolutions=0;
    354354  OsiSolverInterface * solver = NULL;
     
    595595      int newNumberInfeas=0;
    596596      returnCode=0;
    597       if (model_->getCurrentSeconds()>model_->getMaximumSeconds()) {
     597      if (model_->maximumSecondsReached()) {
    598598        exitAll=true;
    599599        break;
     
    609609      if (numberIterationsPass1>=0) {
    610610        int n = totalNumberIterations - numberIterationsLastPass;
    611         double perPass = totalNumberIterations/(totalNumberPasses+numberPasses);
     611        double perPass = totalNumberIterations/
     612          (totalNumberPasses+numberPasses+1.0e-5);
    612613        perPass /= (solver->getNumRows()+numberColumns);
    613614        double test = moreIterations ? 0.3 : 0.0;
     
    737738                  model_->setBestSolution(saveOldSolution,numberColumns,saveObjectiveValue);
    738739                delete [] saveOldSolution;
    739                 if (!action||model_->getCurrentSeconds()>model_->getMaximumSeconds()) {
     740                if (!action||model_->maximumSecondsReached()) {
    740741                  exitAll=true; // exit
    741742                  break;
     
    979980                    model_->setBestSolution(saveOldSolution,numberColumns,saveObjectiveValue);
    980981                  delete [] saveOldSolution;
    981                   if (!action||model_->getCurrentSeconds()>model_->getMaximumSeconds()) {
     982                  if (!action||model_->maximumSecondsReached()) {
    982983                    exitAll=true; // exit
    983984                    break;
     
    18651866#endif
    18661867  delete solver; // probably NULL but do anyway
    1867   if (!finalReturnCode&&closestSolution&&closestObjectiveValue <= 10.0&&usedColumn) {
     1868  if (!finalReturnCode&&closestSolution&&closestObjectiveValue <= 10.0&&
     1869      usedColumn&&!model_->maximumSecondsReached()) {
    18681870    // try a bit of branch and bound
    18691871    OsiSolverInterface * newSolver = cloneBut(1); // was model_->continuousSolver()->clone();
     
    19291931  if (bestBasis.getNumStructural())
    19301932    model_->setBestSolutionBasis(bestBasis);
    1931   model_->setMinimizationObjValue(saveBestObjective);
     1933  //model_->setMinimizationObjValue(saveBestObjective);
    19321934  if ((accumulate_&1)!=0&&numberSolutions>1&&!model_->getSolutionCount()) {
    19331935    model_->setSolutionCount(1); // for local search
  • trunk/Cbc/src/CbcHeuristicGreedy.cpp

    r1205 r1212  
    5555  // Get a copy of original matrix
    5656  assert(model->solver());
    57   matrix_ = *model->solver()->getMatrixByCol();
     57  if (model->solver()->getNumRows()) {
     58    matrix_ = *model->solver()->getMatrixByCol();
     59  }
    5860  originalNumberRows_=model->solver()->getNumRows();
    5961}
     
    462464  // Get a copy of original matrix
    463465  assert(model->solver());
    464   matrix_ = *model->solver()->getMatrixByCol();
     466  if (model->solver()->getNumRows()) {
     467    matrix_ = *model->solver()->getMatrixByCol();
     468  }
    465469  originalNumberRows_=model->solver()->getNumRows();
    466470}
  • trunk/Cbc/src/CbcHeuristicLocal.cpp

    r1211 r1212  
    3737  // Get a copy of original matrix
    3838  assert(model.solver());
    39   matrix_ = *model.solver()->getMatrixByCol();
     39  if (model.solver()->getNumRows()) {
     40    matrix_ = *model.solver()->getMatrixByCol();
     41  }
    4042  int numberColumns = model.solver()->getNumCols();
    4143  used_ = new char[numberColumns];
     
    599601  // Get a copy of original matrix
    600602  assert(model_->solver());
    601   matrix_ = *model_->solver()->getMatrixByCol();
     603  if (model_->solver()->getNumRows()) {
     604    matrix_ = *model_->solver()->getMatrixByCol();
     605  }
    602606  delete [] used_;
    603607  int numberColumns = model->solver()->getNumCols();
     
    636640{
    637641  CbcHeuristicNaive other;
    638   fprintf(fp,"0#include \"CbcHeuristicNaive.hpp\"\n");
     642  fprintf(fp,"0#include \"CbcHeuristicLocal.hpp\"\n");
    639643  fprintf(fp,"3  CbcHeuristicNaive naive(*cbcModel);\n");
    640644  CbcHeuristic::generateCpp(fp,"naive");
     
    876880  model_ = model;
    877881}
     882// Default Constructor
     883CbcHeuristicCrossover::CbcHeuristicCrossover()
     884  :CbcHeuristic(),
     885   numberSolutions_(0),
     886   useNumber_(3)
     887{
     888  setWhen(1);
     889}
     890
     891// Constructor with model - assumed before cuts
     892
     893CbcHeuristicCrossover::CbcHeuristicCrossover(CbcModel & model)
     894  :CbcHeuristic(model),
     895   numberSolutions_(0),
     896   useNumber_(3)
     897{
     898  setWhen(1);
     899  for (int i=0;i<10;i++)
     900    random_[i]=model.randomNumberGenerator()->randomDouble();
     901}
     902
     903// Destructor
     904CbcHeuristicCrossover::~CbcHeuristicCrossover ()
     905{
     906}
     907
     908// Clone
     909CbcHeuristic *
     910CbcHeuristicCrossover::clone() const
     911{
     912  return new CbcHeuristicCrossover(*this);
     913}
     914// Create C++ lines to get to current state
     915void
     916CbcHeuristicCrossover::generateCpp( FILE * fp)
     917{
     918  CbcHeuristicCrossover other;
     919  fprintf(fp,"0#include \"CbcHeuristicLocal.hpp\"\n");
     920  fprintf(fp,"3  CbcHeuristicCrossover crossover(*cbcModel);\n");
     921  CbcHeuristic::generateCpp(fp,"crossover");
     922  if (useNumber_!=other.useNumber_)
     923    fprintf(fp,"3  crossover.setNumberSolutions(%d);\n",useNumber_);
     924  else
     925    fprintf(fp,"4  crossover.setNumberSolutions(%d);\n",useNumber_);
     926  fprintf(fp,"3  cbcModel->addHeuristic(&crossover);\n");
     927}
     928
     929// Copy constructor
     930CbcHeuristicCrossover::CbcHeuristicCrossover(const CbcHeuristicCrossover & rhs)
     931:
     932  CbcHeuristic(rhs),
     933  attempts_(rhs.attempts_),
     934  numberSolutions_(rhs.numberSolutions_),
     935  useNumber_(rhs.useNumber_)
     936{
     937  memcpy(random_,rhs.random_,10*sizeof(double));
     938}
     939
     940// Assignment operator
     941CbcHeuristicCrossover &
     942CbcHeuristicCrossover::operator=( const CbcHeuristicCrossover& rhs)
     943{
     944  if (this!=&rhs) {
     945    CbcHeuristic::operator=(rhs);
     946    useNumber_ = rhs.useNumber_;
     947    attempts_ = rhs.attempts_;
     948    numberSolutions_ = rhs.numberSolutions_;
     949    memcpy(random_,rhs.random_,10*sizeof(double));
     950  }
     951  return *this;
     952}
     953
     954// Resets stuff if model changes
     955void
     956CbcHeuristicCrossover::resetModel(CbcModel * model)
     957{
     958  CbcHeuristic::resetModel(model);
     959}
     960int
     961CbcHeuristicCrossover::solution(double & solutionValue,
     962                         double * betterSolution)
     963{
     964  if (when_==0)
     965    return 0;
     966  numCouldRun_++;
     967  bool useBest=(numberSolutions_!=model_->getSolutionCount());
     968  if (!useBest&&(when_%10)==1)
     969    return 0;
     970  numberSolutions_=model_->getSolutionCount();
     971  OsiSolverInterface * continuousSolver = model_->continuousSolver();
     972  int useNumber =CoinMin(model_->numberSavedSolutions(),useNumber_);
     973  if (useNumber<2||!continuousSolver)
     974    return 0;
     975  // Fix later
     976  if (!useBest)
     977    abort();
     978  numRuns_++;
     979  double cutoff;
     980  model_->solver()->getDblParam(OsiDualObjectiveLimit,cutoff);
     981  double direction = model_->solver()->getObjSense();
     982  cutoff *= direction;
     983  cutoff = CoinMin(cutoff,solutionValue);
     984  OsiSolverInterface * solver = cloneBut(2);
     985  // But reset bounds
     986  solver->setColLower(continuousSolver->getColLower());
     987  solver->setColUpper(continuousSolver->getColUpper());
     988  int numberColumns = solver->getNumCols();
     989  // Fixed
     990  double * fixed =new double [numberColumns];
     991  for (int i=0;i<numberColumns;i++)
     992    fixed[i]=-COIN_DBL_MAX;
     993  int whichSolution[10];
     994  for (int i=0;i<useNumber;i++)
     995    whichSolution[i]=i;
     996  for (int i=0;i<useNumber;i++) {
     997    int k = whichSolution[i];
     998    const double * solution = model_->savedSolution(k);
     999    for (int j=0;j<numberColumns;j++) {
     1000      if (solver->isInteger(j)) {
     1001        if (fixed[j]==-COIN_DBL_MAX)
     1002          fixed[j]=floor(solution[j]+0.5);
     1003        else if (fabs(fixed[j]-solution[j])>1.0e-7)
     1004          fixed[j]=COIN_DBL_MAX;
     1005      }
     1006    }
     1007  }
     1008  const double * colLower = solver->getColLower();
     1009  for (int i=0;i<numberColumns;i++) {
     1010    if (solver->isInteger(i)) {
     1011      double value=fixed[i];
     1012      if (value!=COIN_DBL_MAX) {
     1013        if (when_<10) {
     1014          solver->setColLower(i,value);
     1015          solver->setColUpper(i,value);
     1016        } else if (value==colLower[i]) {
     1017          solver->setColUpper(i,value);
     1018        }
     1019      }
     1020    }
     1021  }
     1022  int returnCode = smallBranchAndBound(solver,numberNodes_,betterSolution,
     1023                                       solutionValue,
     1024                                       solutionValue,"CbcHeuristicCrossover");
     1025  if (returnCode<0)
     1026    returnCode=0; // returned on size
     1027  if ((returnCode&2)!=0) {
     1028    // could add cut
     1029    returnCode &= ~2;
     1030  }
     1031
     1032  delete solver;
     1033  return returnCode;
     1034}
     1035// update model
     1036void CbcHeuristicCrossover::setModel(CbcModel * model)
     1037{
     1038  model_ = model;
     1039  if (model) {
     1040    for (int i=0;i<10;i++)
     1041      random_[i]=model->randomNumberGenerator()->randomDouble();
     1042  }
     1043}
    8781044
    8791045 
  • trunk/Cbc/src/CbcHeuristicLocal.hpp

    r1173 r1212  
    140140};
    141141
     142/** Crossover Search class
     143 */
     144
     145class CbcHeuristicCrossover : public CbcHeuristic {
     146public:
     147
     148  // Default Constructor
     149  CbcHeuristicCrossover ();
     150
     151  /* Constructor with model - assumed before cuts
     152     Initial version does not do Lps
     153  */
     154  CbcHeuristicCrossover (CbcModel & model);
     155 
     156  // Copy constructor
     157  CbcHeuristicCrossover ( const CbcHeuristicCrossover &);
     158   
     159  // Destructor
     160  ~CbcHeuristicCrossover ();
     161 
     162  /// Clone
     163  virtual CbcHeuristic * clone() const;
     164
     165  /// Assignment operator
     166  CbcHeuristicCrossover & operator=(const CbcHeuristicCrossover& rhs);
     167
     168  /// Create C++ lines to get to current state
     169  virtual void generateCpp( FILE * fp) ;
     170
     171  /// Resets stuff if model changes
     172  virtual void resetModel(CbcModel * model);
     173
     174  /// update model (This is needed if cliques update matrix etc)
     175  virtual void setModel(CbcModel * model);
     176 
     177  using CbcHeuristic::solution ;
     178  /** returns 0 if no solution, 1 if valid solution.
     179      Fix variables if agree in useNumber_ solutions
     180      when_ 0 off, 1 only at new solutions, 2 also every now and then
     181      add 10 to make only if agree at lower bound
     182  */
     183  virtual int solution(double & objectiveValue,
     184                       double * newSolution);
     185
     186  /// Sets number of solutions to use
     187  inline void setNumberSolutions(int value)
     188  {
     189    if (value>0&&value<=10)
     190      useNumber_=value;
     191  }
     192
     193protected:
     194  // Data
     195  /// Attempts
     196  std::vector <double> attempts_;
     197  /// Random numbers to stop same search happening
     198  double random_[10];
     199  /// Number of solutions so we only do after new solution
     200  int numberSolutions_;
     201  /// Number of solutions to use
     202  int useNumber_;
     203};
    142204
    143205#endif
  • trunk/Cbc/src/CbcMessage.cpp

    r1211 r1212  
    2020  {CBC_END_GOOD,1,1,"Search completed - best objective %.16g, took %d iterations and %d nodes (%.2f seconds)"},
    2121  {CBC_MAXNODES,3,1,"Exiting on maximum nodes"},
    22   {CBC_MAXTIME,20,1,"Exiting on maximum time"},
    23   {CBC_MAXSOLS,19,1,"Exiting on maximum solutions"},
    24   {CBC_EVENT,27,1,"Exiting on user event"},
    2522  {CBC_SOLUTION,4,1,"Integer solution of %g found after %d iterations and %d nodes (%.2f seconds)"},
    26   {CBC_END_SOLUTION,34,2,"Final check on integer solution of %g found after %d iterations and %d nodes (%.2f seconds)"},
    27   {CBC_SOLUTION2,33,1,"Integer solution of %g found (by alternate solver) after %d iterations and %d nodes (%.2f seconds)"},
    2823  {CBC_END,5,1,"Partial search - best objective %g (best possible %g), took %d iterations and %d nodes (%.2f seconds)"},
    2924  {CBC_INFEAS,6,1,"The LP relaxation is infeasible or too expensive"},
     
    3429  {CBC_GAP,11,1,"Exiting as integer gap of %g less than %g or %g%%"},
    3530  {CBC_ROUNDING,12,1,"Integer solution of %g found by %s after %d iterations and %d nodes (%.2f seconds)"},
    36   {CBC_TREE_SOL,24,1,"Integer solution of %g found by subtree after %d iterations and %d nodes (%.2f seconds)"},
    3731  {CBC_ROOT,13,1,"At root node, %d cuts changed objective from %g to %g in %d passes"},
    38   {CBC_GENERATOR,14,1,"Cut generator %d (%s) - %d row cuts, %d column cuts (%d active) %? in %.3f seconds - new frequency is %d"},
     32  {CBC_GENERATOR,14,1,"Cut generator %d (%s) - %d row cuts average %.1f elements, %d column cuts (%d active) %? in %.3f seconds - new frequency is %d"},
    3933  {CBC_BRANCH,15,2,"Node %d Obj %g Unsat %d depth %d"},
    4034  {CBC_STRONGSOL,16,1,"Integer solution of %g found by strong branching after %d iterations and %d nodes (%.2f seconds)"},
    41   {CBC_NOINT,3007,1,"No integer variables - nothing to do"},
    4235  {CBC_VUB_PASS,17,1,"%d solved, %d variables fixed, %d tightened"},
    4336  {CBC_VUB_END,18,1,"After tightenVubs, %d variables fixed, %d tightened"},
     37  {CBC_MAXSOLS,19,1,"Exiting on maximum solutions"},
     38  {CBC_MAXTIME,20,1,"Exiting on maximum time"},
    4439  {CBC_NOTFEAS1,21,2,"On closer inspection node is infeasible"},
    4540  {CBC_NOTFEAS2,22,2,"On closer inspection objective value of %g above cutoff of %g"},
    4641  {CBC_NOTFEAS3,23,2,"Allowing solution, even though largest row infeasibility is %g"},
    47   {CBC_CUTOFF_WARNING1,23,1,"Cutoff set to %g - equivalent to best solution of %g"},
     42  {CBC_TREE_SOL,24,1,"Integer solution of %g found by subtree after %d iterations and %d nodes (%.2f seconds)"},
    4843  {CBC_ITERATE_STRONG,25,3,"%d cleanup iterations before strong branching"},
    4944  {CBC_PRIORITY,26,1,"Setting priorities for objects %d to %d inclusive (out of %d)"},
    50   {CBC_WARNING_STRONG,3008,1,"Strong branching is fixing too many variables, too expensively!"},
     45  {CBC_EVENT,27,1,"Exiting on user event"},
    5146  {CBC_START_SUB,28,1,"Starting sub-tree for %s - maximum nodes %d"},
    5247  {CBC_END_SUB,29,1,"Ending sub-tree for %s"},
     
    5449  {CBC_CUTS_STATS,31,1,"%d added rows had average density of %g"},
    5550  {CBC_STRONG_STATS,32,1,"Strong branching done %d times (%d iterations), fathomed %d nodes and fixed %d variables"},
     51  {CBC_SOLUTION2,33,1,"Integer solution of %g found (by alternate solver) after %d iterations and %d nodes (%.2f seconds)"},
    5652  {CBC_UNBOUNDED,34,1,"The LP relaxation is unbounded!"},
    5753  {CBC_OTHER_STATS,35,1,"Maximum depth %d, %g variables fixed on reduced cost"},
     
    6763  {CBC_GENERAL,45,1,"%s"},
    6864  {CBC_ROOT_DETAIL,46,2,"Root node pass %d, %d rows, %d total tight cuts  -  objective %g"},
     65  {CBC_CUTOFF_WARNING1,47,1,"Cutoff set to %g - equivalent to best solution of %g"},
     66  {CBC_END_SOLUTION,48,2,"Final check on integer solution of %g found after %d iterations and %d nodes (%.2f seconds)"},
     67  {CBC_NOINT,3007,1,"No integer variables - nothing to do"},
     68  {CBC_WARNING_STRONG,3008,1,"Strong branching is fixing too many variables, too expensively!"},
    6969  {CBC_DUMMY_END,999999,0,""}
    7070};
  • trunk/Cbc/src/CbcModel.cpp

    r1211 r1212  
    577577#endif
    578578    } else if (largestObj<smallestObj*5.0&&!parentModel_&&
    579                numberIntegerObj*1.2<numberColumns) {
     579               !numberContinuousObj&&
     580               numberIntegerObj*2<numberColumns) {
    580581      // up priorities on costed
    581582      int iPriority=-1;
     
    14671468    feasible=false; // pretend infeasible
    14681469  }
     1470  numberSavedSolutions_=0;
    14691471  int saveNumberStrong = numberStrong_;
    14701472  int saveNumberBeforeTrust = numberBeforeTrust_;
     
    16601662  if(solverCharacteristics_->reducedCostsAccurate())
    16611663    analyzeObjective() ;
     1664  {
     1665    // may be able to change cutoff now
     1666    double cutoff = getCutoff();
     1667    double increment = getDblParam(CbcModel::CbcCutoffIncrement) ;
     1668    if (cutoff > bestObjective_-increment) {
     1669      cutoff = bestObjective_-increment ;
     1670      setCutoff(cutoff) ;
     1671    }
     1672  }
    16621673#ifdef COIN_HAS_CLP
    16631674  // Possible save of pivot method
     
    19781989      OsiClpSolverInterface * clpSolver
    19791990        = dynamic_cast<OsiClpSolverInterface *> (copy1);
    1980       if (clpSolver) {
     1991      if (false&&clpSolver) {
    19811992        numberRows = clpSolver->getNumRows();
    19821993        char * rotate = new char[numberRows];
    19831994        int n = clpSolver->getModelPtr()->findNetwork(rotate,1.0);
     1995        delete [] rotate;
    19841996#ifdef CLP_INVESTIGATE
    19851997        printf("INTA network %d rows out of %d\n",n,numberRows);
     
    22432255    }
    22442256    cpxSolver.branchAndBound();
    2245     static double totalTime=0.0;
    22462257    double timeTaken = CoinCpuTime()-time1;
    2247     totalTime += timeTaken;
    2248     sprintf(printBuffer,"Cplex took %g seconds - total %g",
    2249             timeTaken,totalTime);
     2258    sprintf(printBuffer,"Cplex took %g seconds",
     2259            timeTaken);
    22502260    messageHandler()->message(CBC_GENERAL,messages())
    22512261      << printBuffer << CoinMessageEol ;
     
    22662276  }
    22672277#endif
     2278  if(fastNodeDepth_==1000&&/*!parentModel_*/(specialOptions_&2048)==0) {
     2279    fastNodeDepth_=-1;
     2280    CbcObject * obj =
     2281      new CbcFollowOn(this);
     2282    obj->setPriority(1);
     2283    addObjects(1,&obj);
     2284  }
    22682285  int saveNumberSolves=numberSolves_;
    22692286  int saveNumberIterations=numberIterations_;
     
    27972814  */
    27982815  numberLongStrong_=0;
    2799   double totalTime = 0.0;
    28002816  CbcNode * createdNode=NULL;
    28012817#ifdef CBC_THREAD
     
    29342950      }
    29352951      //abort();
     2952    }
     2953  }
     2954  {
     2955    // may be able to change cutoff now
     2956    double cutoff = getCutoff();
     2957    double increment = getDblParam(CbcModel::CbcCutoffIncrement) ;
     2958    if (cutoff > bestObjective_-increment) {
     2959      cutoff = bestObjective_-increment ;
     2960      setCutoff(cutoff) ;
    29362961    }
    29372962  }
     
    32553280  Check for abort on limits: node count, solution count, time, integrality gap.
    32563281*/
    3257     totalTime = getCurrentSeconds() ;
    3258     double maxSeconds = getMaximumSeconds();
    3259     if (parentModel_)
    3260       maxSeconds=CoinMin(maxSeconds,parentModel_->getMaximumSeconds());
    32613282    if (!(numberNodes_ < intParam_[CbcMaxNumNode] &&
    32623283          numberSolutions_ < intParam_[CbcMaxNumSol] &&
    3263           totalTime < maxSeconds &&
     3284          !maximumSecondsReached() &&
    32643285          !stoppedOnGap_&&!eventHappened_&&(maximumNumberIterations_<0||
    32653286                                            numberIterations_<maximumNumberIterations_))) {
     
    34163437      if (parallelMode()>0)
    34173438        unlockThread();
     3439#ifdef CLP_INVESTIGATE
     3440      if (getCutoff()<1.0e20) {
     3441        if (fabs(getCutoff()-(bestObjective_-getCutoffIncrement()))>1.0e-6&&
     3442            !parentModel_)
     3443          printf("model cutoff in status %g, best %g, increment %g\n",
     3444                 getCutoff(),bestObjective_,getCutoffIncrement());
     3445        assert (getCutoff()<bestObjective_-getCutoffIncrement()+1.0e-6);
     3446      }
     3447#endif
    34183448      if (!intParam_[CbcPrinting]) {
    34193449        messageHandler()->message(CBC_STATUS,messages())
     
    37933823  it'll be deleted in cleanTree. We need to check.
    37943824*/
    3795     double maxSeconds = getMaximumSeconds();
    3796     if (parentModel_)
    3797       maxSeconds=CoinMin(maxSeconds,parentModel_->getMaximumSeconds());
    37983825    if (!(numberNodes_ < intParam_[CbcMaxNumNode] &&
    37993826          numberSolutions_ < intParam_[CbcMaxNumSol] &&
    3800           totalTime < maxSeconds &&
     3827          !maximumSecondsReached() &&
    38013828          !stoppedOnGap_&&!eventHappened_&&(maximumNumberIterations_<0||
    38023829                                            numberIterations_<maximumNumberIterations_))) {
     
    38203847            status_ = 1 ; }
    38213848          else
    3822         if (totalTime >= dblParam_[CbcMaximumSeconds])
     3849            if (maximumSecondsReached())
    38233850          { handler_->message(CBC_MAXTIME,messages_) << CoinMessageEol ;
    38243851          secondaryStatus_ = 4;
     
    42894316  sumChangeObjective2_(0.0),
    42904317  bestSolution_(NULL),
     4318  savedSolutions_(NULL),
    42914319  currentSolution_(NULL),
    42924320  testSolution_(NULL),
    42934321  minimumDrop_(1.0e-4),
    42944322  numberSolutions_(0),
     4323  numberSavedSolutions_(0),
     4324  maximumSavedSolutions_(0),
    42954325  stateOfSearch_(0),
    42964326  whenCuts_(-1),
     
    44564486  minimumDrop_(1.0e-4),
    44574487  numberSolutions_(0),
     4488  numberSavedSolutions_(0),
     4489  maximumSavedSolutions_(0),
    44584490  stateOfSearch_(0),
    44594491  whenCuts_(-1),
     
    45764608  strongInfo_[6]=0;
    45774609  solverCharacteristics_ = NULL;
    4578 
    45794610  nodeCompare_=new CbcCompareDefault();;
    45804611  problemFeasibility_=new CbcFeasibilityBase();
     
    46034634  // Initialize solution and integer variable vectors
    46044635  bestSolution_ = NULL; // to say no solution found
     4636  savedSolutions_ = NULL;
    46054637  numberIntegers_=0;
    46064638  int numberColumns = solver_->getNumCols();
     
    47244756  minimumDrop_(rhs.minimumDrop_),
    47254757  numberSolutions_(rhs.numberSolutions_),
     4758  numberSavedSolutions_(rhs.numberSavedSolutions_),
     4759  maximumSavedSolutions_(rhs.maximumSavedSolutions_),
    47264760  stateOfSearch_(rhs.stateOfSearch_),
    47274761  whenCuts_(rhs.whenCuts_),
     
    49304964  }
    49314965  int numberColumns = solver_->getNumCols();
     4966  if (maximumSavedSolutions_&&rhs.savedSolutions_) {
     4967    savedSolutions_ = new double * [maximumSavedSolutions_];
     4968    for (int i=0;i<maximumSavedSolutions_;i++)
     4969      savedSolutions_[i]=CoinCopyOfArray(rhs.savedSolutions_[i],numberColumns+2);
     4970  } else {
     4971    savedSolutions_=NULL;
     4972  }
    49324973  // Space for current solution
    49334974  currentSolution_ = new double[numberColumns];
     
    50285069      bestSolution_=NULL;
    50295070    }
     5071    for (int i=0;i<maximumSavedSolutions_;i++)
     5072      delete [] savedSolutions_[i];
     5073    delete [] savedSolutions_;
     5074    savedSolutions_=NULL;
    50305075    int numberColumns = rhs.getNumCols();
    50315076    if (numberColumns) {
     
    50405085      usedInSolution_=NULL;
    50415086    }
     5087    if (maximumSavedSolutions_) {
     5088      savedSolutions_ = new double * [maximumSavedSolutions_];
     5089      for (int i=0;i<maximumSavedSolutions_;i++)
     5090        savedSolutions_[i]=CoinCopyOfArray(rhs.savedSolutions_[i],numberColumns+2);
     5091    } else {
     5092      savedSolutions_=NULL;
     5093    }
    50425094    testSolution_=currentSolution_;
    50435095    minimumDrop_ = rhs.minimumDrop_;
    50445096    numberSolutions_=rhs.numberSolutions_;
     5097    numberSavedSolutions_=rhs.numberSavedSolutions_;
     5098    maximumSavedSolutions_=rhs.maximumSavedSolutions_;
    50455099    stateOfSearch_= rhs.stateOfSearch_;
    50465100    whenCuts_ = rhs.whenCuts_;
     
    53645418  delete continuousSolver_;
    53655419  continuousSolver_=NULL;
     5420  numberSavedSolutions_=0;
    53665421  delete [] bestSolution_;
    53675422  bestSolution_=NULL;
     5423  if (savedSolutions_) {
     5424    for (int i=0;i<maximumSavedSolutions_;i++)
     5425      delete [] savedSolutions_[i];
     5426    delete [] savedSolutions_;
     5427    savedSolutions_=NULL;
     5428  }
    53685429  delete [] currentSolution_;
    53695430  currentSolution_=NULL;
     
    54745535  resolveAfterTakeOffCuts_ = rhs.resolveAfterTakeOffCuts_;
    54755536  maximumNumberIterations_ = rhs.maximumNumberIterations_;
     5537  numberSavedSolutions_=rhs.numberSavedSolutions_;
     5538  maximumSavedSolutions_=rhs.maximumSavedSolutions_;
     5539  if (maximumSavedSolutions_) {
     5540    int n=solver_->getNumCols();
     5541    savedSolutions_ = new double * [maximumSavedSolutions_];
     5542    for (int i=0;i<maximumSavedSolutions_;i++)
     5543      savedSolutions_[i]=CoinCopyOfArray(rhs.savedSolutions_[i],n+2);
     5544  }
    54765545  continuousPriority_ = rhs.continuousPriority_;
    54775546  numberThreads_ = rhs.numberThreads_;
     
    65696638  sumChangeObjective1_ += solver_->getObjValue()*solver_->getObjSense()
    65706639    - objectiveValue ;
    6571   if ( getCurrentSeconds() > dblParam_[CbcMaximumSeconds] )
     6640  if ( maximumSecondsReached() )
    65726641    numberTries=0; // exit
    65736642  //if ((numberNodes_%100)==0)
     
    75907659      feasible = ( resolve(node ? node->nodeInfo() : NULL,2) != 0) ;
    75917660      //solver_->setHintParam(OsiDoDualInResolve,true,OsiHintTry);
    7592       if ( getCurrentSeconds() > dblParam_[CbcMaximumSeconds] )
     7661      if ( maximumSecondsReached() ) {
    75937662        numberTries=0; // exit
     7663        break;
     7664      }
    75947665#     ifdef CBC_DEBUG
    75957666      printf("Obj value after cuts %g %d rows\n",solver_->getObjValue(),
     
    78147885    reducedCostFix() ;
    78157886  // If at root node do heuristics
    7816   if (!numberNodes_) {
     7887  if (!numberNodes_&&!maximumSecondsReached()) {
    78177888    // First see if any cuts are slack
    78187889    int numberRows = solver_->getNumRows();
     
    81558226              }
    81568227            }
    8157           } else if (thisCuts+generator_[i]->numberColumnCuts()<smallProblem) {
     8228          } else if ((thisCuts+generator_[i]->numberColumnCuts()<smallProblem)
     8229                     &&!generator_[i] ->whetherToUse()) {
    81588230            if (howOften!=1&&!probingWasOnBut) {
    81598231              if (generator_[i]->whatDepth()<0||howOften!=-1) {
     
    82118283      if (stored&&!generator_[i]->numberCutsInTotal())
    82128284        continue;
     8285      double average=0.0;
     8286      int n=generator_[i]->numberCutsInTotal();
     8287      if (n) {
     8288        average = generator_[i]->numberElementsInTotal();
     8289        average /= n;
     8290      }
    82138291      if (handler_->logLevel()>1||!numberNodes_) {
    82148292        handler_->message(CBC_GENERATOR,messages_)
     
    82168294          <<generator_[i]->cutGeneratorName()
    82178295          //<<generator_[i]->numberCutsInTotal()<<count[i]
    8218           << generator_[i]->numberCutsInTotal()
     8296          << n
     8297          << average
    82198298          <<generator_[i]->numberColumnCuts()
    82208299          <<generator_[i]->numberCutsActive()
     
    1049810577      */
    1049910578      specialOptions_ |= 256; // mark as full cut scan should be done
    10500       bestObjective_ = objectiveValue;
    10501       int numberColumns = solver_->getNumCols();
    10502       if (!bestSolution_)
    10503         bestSolution_ = new double[numberColumns];
    10504       CoinCopyN(solution,numberColumns,bestSolution_);
     10579      saveBestSolution(solution,objectiveValue);
     10580      //bestObjective_ = objectiveValue;
     10581      //int numberColumns = solver_->getNumCols();
     10582      //if (!bestSolution_)
     10583      //bestSolution_ = new double[numberColumns];
     10584      //CoinCopyN(solution,numberColumns,bestSolution_);
    1050510585     
    1050610586      cutoff = bestObjective_-dblParam_[CbcCutoffIncrement];
     
    1067810758        NB - Not all of this if from solve with cuts
    1067910759      */
    10680       bestObjective_ = objectiveValue;
    10681       int numberColumns = solver_->getNumCols();
    10682       if (!bestSolution_)
    10683         bestSolution_ = new double[numberColumns];
    10684       CoinCopyN(candidate,numberColumns,bestSolution_);
     10760      saveBestSolution(solution,objectiveValue);
     10761      //bestObjective_ = objectiveValue;
     10762      //int numberColumns = solver_->getNumCols();
     10763      //if (!bestSolution_)
     10764      //bestSolution_ = new double[numberColumns];
     10765      //CoinCopyN(candidate,numberColumns,bestSolution_);
    1068510766
    1068610767      // don't update if from solveWithCuts
     
    1208112162      messageHandler()->message(CBC_GENERAL,messages())
    1208212163        << printBuffer << CoinMessageEol ;
    12083       // may be able to change cutoff now
    12084       double cutoff = getCutoff();
    12085       double increment = getDblParam(CbcModel::CbcCutoffIncrement) ;
    12086       if (cutoff > objectiveValue-increment) {
    12087         cutoff = objectiveValue-increment ;
    12088         setCutoff(cutoff) ;
    12089       }
    12090     }
    12091   }
     12164    }
     12165  }
     12166  if (bestSolution_)
     12167    saveExtraSolution(bestSolution_,bestObjective_);
    1209212168  bestObjective_ = objectiveValue;
     12169  // may be able to change cutoff now
     12170  double cutoff = getCutoff();
     12171  double increment = getDblParam(CbcModel::CbcCutoffIncrement) ;
     12172  if (cutoff > objectiveValue-increment) {
     12173    cutoff = objectiveValue-increment ;
     12174    setCutoff(cutoff) ;
     12175  }
    1209312176  int n = CoinMax(numberColumns,solver_->getNumCols());
    1209412177  delete [] bestSolution_;
     
    1222212305          if (!heuristic_[i]->shouldHeurRun(whereFrom))
    1222312306            continue;
     12307          if (maximumSecondsReached())
     12308            break;
    1222412309          // see if heuristic will do anything
    1222512310          double saveValue = heuristicValue ;
     
    1274912834            delete [] up;
    1275012835            delete [] numberDown;
     12836            delete [] priority;
    1275112837            delete [] numberUp;
    1275212838            delete [] numberDownInfeasible;
     
    1275612842            solver_->getHintParam(OsiDoReducePrint,takeHint,strength);
    1275712843            ClpSimplex * simplex = clpSolver->getModelPtr();
     12844            //printf("mod cutoff %g solver %g offset %g\n",
     12845            //   getCutoff(),simplex->dualObjectiveLimit(),simplex->objectiveOffset());
    1275812846            int saveLevel = simplex->logLevel();
    1275912847            if (strength!=OsiHintIgnore&&takeHint&&saveLevel==1)
     
    1282612914              double largest=-1.0;
    1282712915              for (int i=0;i<numberIntegers_;i++) {
    12828 #ifndef NDEBUG
    1282912916                CbcSimpleIntegerDynamicPseudoCost * obj =
    1283012917                  dynamic_cast <CbcSimpleIntegerDynamicPseudoCost *>(object_[i]) ;
    12831                 assert (obj&&obj->columnNumber()==integerVariable_[i]);
    12832 #else
    12833                 CbcSimpleIntegerDynamicPseudoCost * obj =
    12834                   static_cast <CbcSimpleIntegerDynamicPseudoCost *>(object_[i]) ;
    12835 #endif
     12918                if (!obj)
     12919                  continue;
     12920                assert (obj->columnNumber()==integerVariable_[i]);
    1283612921                if (info->numberUp_[i]>0) {
    1283712922                  if (info->downPseudo_[i]>largest)
     
    1446114546  }
    1446214547}
     14548// Number of saved solutions (including best)
     14549int
     14550CbcModel::numberSavedSolutions() const
     14551{
     14552  if (!bestSolution_)
     14553    return 0;
     14554  else
     14555    return numberSavedSolutions_+1;
     14556}
     14557// Set maximum number of extra saved solutions
     14558void
     14559CbcModel::setMaximumSavedSolutions(int value)
     14560{
     14561  if (value<maximumSavedSolutions_) {
     14562    for (int i=value;i<maximumSavedSolutions_;i++)
     14563      delete [] savedSolutions_[i];
     14564    maximumSavedSolutions_=value;
     14565    numberSavedSolutions_ = CoinMin(numberSavedSolutions_,
     14566                                    maximumSavedSolutions_);
     14567    if (!maximumSavedSolutions_)
     14568      delete [] savedSolutions_;
     14569  } else if (value>maximumSavedSolutions_) {
     14570    double ** temp = new double * [value];
     14571    int i;
     14572    for ( i=0;i<maximumSavedSolutions_;i++)
     14573      temp[i]=savedSolutions_[i];
     14574    for ( ;i<value;i++)
     14575      temp[i]=NULL;
     14576    delete [] savedSolutions_;
     14577    maximumSavedSolutions_=value;
     14578    savedSolutions_ = temp;
     14579  }
     14580}
     14581// Return a saved solution objective (0==best) - COIN_DBL_MAX if off end
     14582double
     14583CbcModel::savedSolutionObjective(int which) const
     14584{
     14585  if (which==0) {
     14586    return bestObjective_;
     14587  } else if (which<=numberSavedSolutions_) {
     14588    double * sol = savedSolutions_[which-1];
     14589    assert (static_cast<int>(sol[0])==solver_->getNumCols());
     14590    return sol[1];
     14591  } else {
     14592    return COIN_DBL_MAX;
     14593  }
     14594}
     14595// Return a saved solution (0==best) - NULL if off end
     14596const double *
     14597CbcModel::savedSolution(int which) const
     14598{
     14599  if (which==0) {
     14600    return bestSolution_;
     14601  } else if (which<=numberSavedSolutions_) {
     14602    double * sol = savedSolutions_[which-1];
     14603    assert (static_cast<int>(sol[0])==solver_->getNumCols());
     14604    return sol+2;
     14605  } else {
     14606    return NULL;
     14607  }
     14608}
     14609// Save a solution
     14610void
     14611CbcModel::saveExtraSolution(const double * solution, double objectiveValue)
     14612{
     14613  double * save=NULL;
     14614  if (maximumSavedSolutions_) {
     14615    if (!savedSolutions_) {
     14616      savedSolutions_ = new double * [maximumSavedSolutions_];
     14617      for (int i=0;i<maximumSavedSolutions_;i++)
     14618        savedSolutions_[i]=NULL;
     14619    }
     14620    int n=solver_->getNumCols();
     14621    int k;
     14622    for (k=numberSavedSolutions_-1;k>=0;k--) {
     14623      double * sol = savedSolutions_[k];
     14624      assert (static_cast<int>(sol[0])==n);
     14625      if (objectiveValue>sol[1])
     14626        break;
     14627    }
     14628    k++; // where to put
     14629    if (k<maximumSavedSolutions_) {
     14630      if (numberSavedSolutions_==maximumSavedSolutions_) {
     14631        save = savedSolutions_[numberSavedSolutions_-1];
     14632      } else {
     14633        save = new double [n+2];
     14634        numberSavedSolutions_++;
     14635      }
     14636      // move up
     14637      for (int j=maximumSavedSolutions_-1;j>k;j--)
     14638        savedSolutions_[j]=savedSolutions_[j-1];
     14639      savedSolutions_[k]=save;
     14640      save[0]=n;
     14641      save[1]=objectiveValue;
     14642      memcpy(save+2,solution,n*sizeof(double));
     14643    }
     14644  }
     14645}
     14646// Save a solution to best and move current to saved
     14647void
     14648CbcModel::saveBestSolution(const double * solution, double objectiveValue)
     14649{
     14650  int n=solver_->getNumCols();
     14651  if (bestSolution_)
     14652    saveExtraSolution(bestSolution_,bestObjective_);
     14653  else
     14654    bestSolution_ = new double [n];
     14655  bestObjective_=objectiveValue;
     14656  memcpy(bestSolution_,solution,n*sizeof(double));
     14657}
     14658// Delete best and saved solutions
     14659void
     14660CbcModel::deleteSolutions()
     14661{
     14662  delete [] bestSolution_;
     14663  bestSolution_=NULL;
     14664  for (int i=0;i<maximumSavedSolutions_;i++) {
     14665    delete [] savedSolutions_[i];
     14666    savedSolutions_[i]=NULL;
     14667  }
     14668}
    1446314669// Below this is deprecated or at least fairly deprecated
    1446414670/*
     
    1492515131    solver_->setInteger(index);
    1492615132}
     15133// Return true if maximum time reached
     15134bool
     15135CbcModel::maximumSecondsReached() const
     15136{
     15137  double totalTime = getCurrentSeconds() ;
     15138  double maxSeconds = getMaximumSeconds();
     15139  if (totalTime>=maxSeconds)
     15140    return true;
     15141  if (!parentModel_) {
     15142    return false;
     15143  } else {
     15144    totalTime += parentModel_->getCurrentSeconds();
     15145    maxSeconds=parentModel_->getMaximumSeconds();
     15146    return (totalTime>=maxSeconds);
     15147  }
     15148}
    1492715149// Check original model before it gets messed up
    1492815150void
  • trunk/Cbc/src/CbcModel.hpp

    r1211 r1212  
    199199     void branchAndBound(int doStatistics=0);
    200200private:
     201
    201202    /** \brief Evaluate a subproblem using cutting planes and heuristics
    202203
     
    562563  /// Current time since start of branchAndbound
    563564  double getCurrentSeconds() const ;
     565
     566  /// Return true if maximum time reached
     567  bool maximumSecondsReached() const ;
    564568
    565569  /** Set the
     
    11581162  inline void setSolutionCount(int value)
    11591163  { numberSolutions_=value;}
     1164  /// Number of saved solutions (including best)
     1165  int numberSavedSolutions() const;
     1166  /// Maximum number of extra saved solutions
     1167  inline int maximumSavedSolutions() const
     1168  { return maximumSavedSolutions_;}
     1169  /// Set maximum number of extra saved solutions
     1170  void setMaximumSavedSolutions(int value);
     1171  /// Return a saved solution (0==best) - NULL if off end
     1172  const double * savedSolution(int which) const;
     1173  /// Return a saved solution objective (0==best) - COIN_DBL_MAX if off end
     1174  double savedSolutionObjective(int which) const;
     1175
    11601176  /** Current phase (so heuristics etc etc can find out).
    11611177      0 - initial solve
     
    16911707  */
    16921708  void synchronizeHandlers(int makeDefault);
    1693      
     1709  /// Save a solution to saved list
     1710  void saveExtraSolution(const double * solution, double objectiveValue);
     1711  /// Save a solution to best and move current to saved
     1712  void saveBestSolution(const double * solution, double objectiveValue);
     1713  /// Delete best and saved solutions
     1714  void deleteSolutions();
    16941715  /// Encapsulates solver resolve
    16951716  int resolve(OsiSolverInterface * solver);
     
    19511972  /// Array holding the incumbent (best) solution.
    19521973  double * bestSolution_;
     1974  /// Arrays holding other solutions.
     1975  double ** savedSolutions_;
    19531976
    19541977  /** Array holding the current solution.
     
    19751998  /// Number of solutions
    19761999  int numberSolutions_;
     2000  /// Number of saved solutions
     2001  int numberSavedSolutions_;
     2002  /// Maximum number of saved solutions
     2003  int maximumSavedSolutions_;
    19772004  /** State of search
    19782005      0 - no solution
  • trunk/Cbc/src/CbcNode.cpp

    r1205 r1212  
    24072407  // If very odd set of objects then use older chooseBranch
    24082408  bool useOldWay = false;
     2409  // point to useful information
     2410  OsiBranchingInformation usefulInfo = model->usefulInformation();
    24092411  if (numberObjects>model->numberIntegers()) {
    24102412    for (int i=model->numberIntegers();i<numberObjects;i++) {
     
    24122414      CbcObject * obj = dynamic_cast <CbcObject *>(object) ;
    24132415      if (!obj || !obj->optionalObject()) {
    2414         useOldWay=true;
    2415         break;
     2416        int preferredWay;
     2417        double infeasibility = object->infeasibility(&usefulInfo,preferredWay);
     2418        if (infeasibility) {
     2419          useOldWay=true;
     2420          break;
     2421        }
    24162422      }
    24172423    }
     
    24222428  if (useOldWay)
    24232429    return -3;
    2424   // point to useful information
    2425   OsiBranchingInformation usefulInfo = model->usefulInformation();
    2426   // and modify
     2430  // Modify useful info
    24272431  usefulInfo.depth_=depth_;
    24282432  if ((model->specialOptions()&128)!=0) {
  • trunk/Cbc/src/CbcSolver.cpp

    r1211 r1212  
    521521  parameters_[whichParam(GREEDY,numberParameters_,parameters_)].setCurrentOption("on");
    522522  parameters_[whichParam(COMBINE,numberParameters_,parameters_)].setCurrentOption("on");
     523  parameters_[whichParam(CROSSOVER2,numberParameters_,parameters_)].setCurrentOption("off");
    523524  parameters_[whichParam(PIVOTANDFIX,numberParameters_,parameters_)].setCurrentOption("off");
    524525  parameters_[whichParam(RANDROUND,numberParameters_,parameters_)].setCurrentOption("off");
     
    867868int CbcOrClpRead_mode=1;
    868869FILE * CbcOrClpReadCommand=stdin;
     870extern int CbcOrClpEnvironmentIndex;
    869871static bool noPrinting=false;
    870872#ifndef CBC_OTHER_SOLVER
     
    32613263  parameters[whichParam(GREEDY,numberParameters,parameters)].setCurrentOption("on");
    32623264  parameters[whichParam(COMBINE,numberParameters,parameters)].setCurrentOption("on");
     3265  parameters[whichParam(CROSSOVER2,numberParameters,parameters)].setCurrentOption("off");
    32633266  parameters[whichParam(PIVOTANDFIX,numberParameters,parameters)].setCurrentOption("off");
    32643267  parameters[whichParam(RANDROUND,numberParameters,parameters)].setCurrentOption("off");
     
    32973300  int useGreedy = parameters_[whichParam(GREEDY,numberParameters_,parameters_)].currentOptionAsInteger();
    32983301  int useCombine = parameters_[whichParam(COMBINE,numberParameters_,parameters_)].currentOptionAsInteger();
     3302  int useCrossover = parameters_[whichParam(CROSSOVER2,numberParameters_,parameters_)].currentOptionAsInteger();
    32993303  int usePivot = parameters_[whichParam(PIVOTANDFIX,numberParameters_,parameters_)].currentOptionAsInteger();
    33003304  int useRand = parameters_[whichParam(RANDROUND,numberParameters_,parameters_)].currentOptionAsInteger();
     
    34753479    anyToDo=true;
    34763480  }
    3477   if (useCombine>=type&&useCombine>=kType&&useCombine<=kType+1) {
    3478     CbcHeuristicLocal heuristic2(*model);
    3479     heuristic2.setHeuristicName("combine solutions");
    3480     heuristic2.setFractionSmall(0.5);
    3481     heuristic2.setSearchType(1);
    3482     model->addHeuristic(&heuristic2);
    3483     anyToDo=true;
    3484   }
    34853481  if (useGreedy>=type&&useGreedy>=kType&&useGreedy<=kType+1) {
    34863482    CbcHeuristicGreedyCover heuristic3(*model);
     
    36383634    anyToDo=true;
    36393635  }
     3636  if (useCombine>=type&&useCombine>=kType&&useCombine<=kType+1) {
     3637    CbcHeuristicLocal heuristic2(*model);
     3638    heuristic2.setHeuristicName("combine solutions");
     3639    heuristic2.setFractionSmall(0.5);
     3640    heuristic2.setSearchType(1);
     3641    model->addHeuristic(&heuristic2);
     3642    anyToDo=true;
     3643  }
     3644  if (useCrossover>=kType&&useCrossover<=kType+1) {
     3645    CbcHeuristicCrossover heuristic2a(*model);
     3646    heuristic2a.setHeuristicName("crossover");
     3647    heuristic2a.setFractionSmall(0.3);
     3648    // just fix at lower
     3649    heuristic2a.setWhen(11);
     3650    model->addHeuristic(&heuristic2a);
     3651    model->setMaximumSavedSolutions(5);
     3652    anyToDo=true;
     3653  }
    36403654  int heurSwitches=parameters_[whichParam(HOPTIONS,numberParameters_,parameters_)].intValue()%100;
    36413655  if (heurSwitches) {
     
    37033717
    37043718int CbcClpUnitTest (const CbcModel & saveModel,
    3705                     std::string& dirMiplib, bool unitTestOnly,
     3719                    std::string& dirMiplib, int testSwitch,
    37063720                    double * stuff);
    37073721#if NEW_STYLE_SOLVER
     
    68346848                if (howOften==-98||howOften==-99)
    68356849                  generator->setSwitchOffIfLessThan(switches[iGenerator]);
     6850                // Use if any at root as more likely later and fairly cheap
     6851                //if (switches[iGenerator]==-2)
     6852                //generator->setWhetherToUse(true);
    68366853                generator->setInaccuracy(accuracyFlag[iGenerator]);
    68376854                generator->setTiming(true);
     
    82388255                   says -miplib
    82398256                */
    8240                 int extra1 = parameters_[whichParam(EXTRA1,numberParameters_,parameters_)].intValue();
     8257                int extra2 = parameters_[whichParam(EXTRA2,numberParameters_,parameters_)].intValue();
    82418258                double stuff[11];
    82428259                stuff[0]=parameters_[whichParam(FAKEINCREMENT,numberParameters_,parameters_)].doubleValue();
     
    82838300                  }
    82848301                }
    8285                 int returnCode=CbcClpUnitTest(model_, dirMiplib, extra1==1,stuff);
     8302                int returnCode=CbcClpUnitTest(model_, dirMiplib, extra2,stuff);
    82868303                babModel_=NULL;
    82878304                return returnCode;
     
    84078424                  for (int i=0;i<numberColumns2;i++) {
    84088425                    int jColumn = originalColumns[i];
    8409                     solution2[jColumn]=solution[i];
    8410                     lower2[jColumn]=columnLower[i];
    8411                     upper2[jColumn]=columnUpper[i];
     8426                    if (jColumn<n) {
     8427                      solution2[jColumn]=solution[i];
     8428                      lower2[jColumn]=columnLower[i];
     8429                      upper2[jColumn]=columnUpper[i];
     8430                    }
    84128431                  }
    84138432#ifndef CBC_OTHER_SOLVER
     
    84168435#endif
    84178436                  process.postProcess(*babModel_->solver());
     8437#ifdef COIN_DEVELOP
     8438                  if (model_.bestSolution()&&fabs(model_.getMinimizationObjValue()-
     8439                                                  babModel_->getMinimizationObjValue())<1.0e-8) {
     8440                    const double * b1 =model_.bestSolution();
     8441                    const double * b2 = saveSolver->getColSolution();
     8442                    const double * columnLower = saveSolver->getColLower() ;
     8443                    const double * columnUpper = saveSolver->getColUpper() ;
     8444                    for (int i=0;i<n;i++) {
     8445                      if (fabs(b1[i]-b2[i])>1.0e-7) {
     8446                        printf("%d %g %g %g %g\n",i,b1[i],b2[i],
     8447                               columnLower[i],columnUpper[i]);
     8448                      }
     8449                    }
     8450                  }
     8451#endif
    84188452                  // Solution now back in saveSolver
    84198453                  // Double check bounds
     
    84238457                  int numberChanged=0;
    84248458                  for (int i=0;i<n;i++) {
     8459                    if (!saveSolver->isInteger(i))
     8460                      continue;
    84258461                    if (lower2[i]!=COIN_DBL_MAX) {
    84268462                      if (lower2[i]!=columnLower[i]||
     
    84438479                          saveSolver->setColLower(i,lower2[i]);
    84448480                          saveSolver->setColUpper(i,upper2[i]);
     8481                        }
     8482                      }
     8483                    }
     8484                  }
     8485                  // See if sos so we can fix
     8486                  OsiClpSolverInterface * osiclp = dynamic_cast< OsiClpSolverInterface*> (saveSolver);
     8487                  if (osiclp&&osiclp->numberSOS()) {
     8488                    // SOS
     8489                    numberSOS = osiclp->numberSOS();
     8490                    const CoinSet * setInfo = osiclp->setInfo();
     8491                    int i;
     8492                    for ( i=0;i<numberSOS;i++) {
     8493                      int type = setInfo[i].setType();
     8494                      int n=setInfo[i].numberEntries();
     8495                      const int * which = setInfo[i].which();
     8496                      int first=-1;
     8497                      int last=-1;
     8498                      for (int j=0;j<n;j++) {
     8499                        int iColumn = which[j];
     8500                        if (fabs(solution[iColumn])>1.0e-7) {
     8501                          last=j;
     8502                          if (first<0)
     8503                            first=j;
     8504                        }
     8505                      }
     8506                      assert (last-first<type);
     8507                      for (int j=0;j<n;j++) {
     8508                        if (j<first||j>last) {
     8509                          int iColumn = which[j];
     8510                          saveSolver->setColLower(iColumn,0.0);
     8511                          saveSolver->setColUpper(iColumn,0.0);
    84458512                        }
    84468513                      }
     
    84728539#ifndef CBC_OTHER_SOLVER
    84738540                  // and original solver
     8541                  originalSolver->setDblParam(OsiDualObjectiveLimit,COIN_DBL_MAX);
    84748542                  assert (n>=originalSolver->getNumCols());
    84758543                  n=originalSolver->getNumCols();
     
    84998567                }
    85008568#if NEW_STYLE_SOLVER==0
    8501                 if (returnMode==1)
     8569                if (returnMode==1) {
     8570                  model_.deleteSolutions();
    85028571                  model_.setBestSolution(bestSolution,n,babModel_->getMinimizationObjValue());
    8503 #endif
     8572                }
     8573#endif
     8574                babModel_->deleteSolutions();
    85048575                babModel_->setBestSolution(bestSolution,n,babModel_->getMinimizationObjValue());
    85058576#if NEW_STYLE_SOLVER==0
     
    85458616                bestSolution = new double [n];
    85468617                // Put solution now back in saveSolver
     8618                saveSolver->setColSolution(model_.bestSolution());
    85478619                babModel_->assignSolver(saveSolver);
    8548                 saveSolver->setColSolution(model_.bestSolution());
     8620                babModel_->setMinimizationObjValue(model_.getMinimizationObjValue());
    85498621                memcpy(bestSolution,babModel_->solver()->getColSolution(),n*sizeof(double));
    85508622#if NEW_STYLE_SOLVER==0
     
    99179989          case UNITTEST:
    99189990            {
    9919               CbcClpUnitTest(model_, dirSample, true,NULL);
     9991              CbcClpUnitTest(model_, dirSample, -2,NULL);
    99209992            }
    99219993            break;
     
    1056410636          case DUMMY:
    1056510637            break;
     10638          case ENVIRONMENT:
     10639            CbcOrClpEnvironmentIndex=0;
     10640            break;
    1056610641          default:
    1056710642            abort();
  • trunk/Cbc/src/CbcSolver.hpp

    r1200 r1212  
    224224               1 - errors
    225225  */
    226   virtual int importData(CbcSolver * /*model*/, int & /*argc*/, char * /*argv[]*/) {return -1;}
     226  virtual int importData(CbcSolver * /*model*/, int & /*argc*/, char ** /*argv[]*/) {return -1;}
    227227  /// Export 1 OsiClpSolver, 2 CbcModel - add 10 if infeasible from odd situation
    228228  virtual void exportSolution(CbcSolver * /*model*/,
  • trunk/Cbc/src/unitTestClp.cpp

    r1200 r1212  
    7979
    8080//#############################################################################
    81 
     81// testSwitch -2 unitTest, -1 normal (==2)
    8282int CbcClpUnitTest (const CbcModel & saveModel, std::string& dirMiplib,
    83                     bool unitTestOnly,
     83                    int testSwitch,
    8484                    double * stuff)
    8585{
     
    107107  std::vector<double> objValueC ;
    108108  std::vector<double> objValue ;
    109   std::vector<int> strategy ;
     109  std::vector<int> testSet ;
    110110  /*
    111111    And a macro to make the vector creation marginally readable.
     
    113113#define PUSH_MPS(zz_mpsName_zz,\
    114114                 zz_nRows_zz,zz_nCols_zz,zz_objValue_zz,zz_objValueC_zz, \
    115                  zz_strategy_zz) \
     115                 zz_testSet_zz) \
    116116  mpsName.push_back(zz_mpsName_zz) ; \
    117117  nRows.push_back(zz_nRows_zz) ; \
    118118  nCols.push_back(zz_nCols_zz) ; \
    119119  objValueC.push_back(zz_objValueC_zz) ; \
    120   strategy.push_back(zz_strategy_zz) ; \
     120  testSet.push_back(zz_testSet_zz) ; \
    121121  objValue.push_back(zz_objValue_zz) ;
    122 
    123   if (unitTestOnly) {
    124     PUSH_MPS("p0033",16,33,3089,2520.57,7);
    125     PUSH_MPS("p0201",133,201,7615,6875.0,7);
     122  int loSwitch=0;
     123  if (testSwitch==-2) {
     124    PUSH_MPS("p0033",16,33,3089,2520.57,0);
     125    PUSH_MPS("p0201",133,201,7615,6875.0,0);
    126126  } else {
     127    if (testSwitch==-1) {
     128      testSwitch=1;
     129    } else {
     130      loSwitch=static_cast<int>(stuff[6]);
     131      printf("Solving miplib problems in sets >= %d and <=%d\n",
     132             loSwitch,testSwitch);
     133    }
    127134    /*
    128135      Load up the problem vector. Note that the row counts here include the
     
    130137    */
    131138    // 0 for no test, 1 for some, 2 for many, 3 for all
    132     //PUSH_MPS("blend2",274,353,7.598985,6.9156751140,7);
    133     //PUSH_MPS("p2756",755,2756,3124,2688.75,7);
    134     //PUSH_MPS("seymour_1",4944,1372,410.7637014,404.35152,7);
    135     //PUSH_MPS("enigma",21,100,0.0,0.0,7);
    136     //PUSH_MPS("misc03",96,160,3360,1910.,7);
    137     //PUSH_MPS("p0201",133,201,7615,6875.0,7);
    138 #define HOWMANY 2
     139    //PUSH_MPS("blend2",274,353,7.598985,6.9156751140,0);
     140    //PUSH_MPS("p2756",755,2756,3124,2688.75,0);
     141    //PUSH_MPS("seymour_1",4944,1372,410.7637014,404.35152,0);
     142    //PUSH_MPS("enigma",21,100,0.0,0.0,0);
     143    //PUSH_MPS("misc03",96,160,3360,1910.,0);
     144    //PUSH_MPS("p0201",133,201,7615,6875.0,0);
     145#define HOWMANY 6
    139146#if HOWMANY
    140 #if HOWMANY>1
    141     PUSH_MPS("10teams",230,2025,924,917,7);
    142 #endif
    143     PUSH_MPS("air03",124,10757,340160,338864.25,7);
    144 #if HOWMANY>2
     147    PUSH_MPS("10teams",230,2025,924,917,1);
     148    PUSH_MPS("air03",124,10757,340160,338864.25,0);
    145149    PUSH_MPS("air04",823,8904,56137,55535.436,8);
    146150    PUSH_MPS("air05",426,7195,26374,25877.609,8);
    147 #endif
    148   //    PUSH_MPS("arki001",1048,1388,7580813.0459,7579599.80787,7);
    149     PUSH_MPS("bell3a",123,133,878430.32,862578.64,7);
    150 #if HOWMANY>1
    151     PUSH_MPS("bell5",91,104,8966406.49,8608417.95,7);
    152 #endif
    153     PUSH_MPS("blend2",274,353,7.598985,6.9156751140,7);
    154 #if HOWMANY>1
    155     PUSH_MPS("cap6000",2176,6000,-2451377,-2451537.325,7);
    156 #endif
    157     //    PUSH_MPS("dano3mip",3202,13873,728.1111,576.23162474,7);
    158     //PUSH_MPS("danoint",664,521,65.67,62.637280418,7);
    159     PUSH_MPS("dcmulti",290,548,188182,183975.5397,7);
    160     PUSH_MPS("dsbmip",1182,1886,-305.19817501,-305.19817501,7);
    161     PUSH_MPS("egout",98,141,568.101,149.589,7);
    162     PUSH_MPS("enigma",21,100,0.0,0.0,7);
    163 #if HOWMANY>3
    164     PUSH_MPS("fast0507",507,63009,174,172.14556668,7);
    165 #endif
    166     PUSH_MPS("fiber",363,1298,405935.18000,156082.51759,7);
    167 #if HOWMANY>1
    168     PUSH_MPS("fixnet6",478,878,3983,1200.88,7);
    169 #endif
    170     PUSH_MPS("flugpl",18,18,1201500,1167185.7,7);
    171     PUSH_MPS("gen",780,870,112313,112130.0,7);
    172 #if HOWMANY>1
    173     PUSH_MPS("gesa2",1392,1224,25779856.372,25476489.678,7);
    174     PUSH_MPS("gesa2_o",1248,1224,25779856.372,25476489.678,7);
    175 #endif
    176     PUSH_MPS("gesa3",1368,1152,27991042.648,27833632.451,7);
    177     PUSH_MPS("gesa3_o",1224,1152,27991042.648,27833632.451,7);
    178     PUSH_MPS("gt2",29,188,21166.000,13460.233074,7);
    179 #if HOWMANY>2
    180     PUSH_MPS("harp2",112,2993,-73899798.00,-74353341.502,7);
    181 #endif
    182     PUSH_MPS("khb05250",101,1350,106940226,95919464.0,7);
    183 #if HOWMANY>1
    184     PUSH_MPS("l152lav",97,1989,4722,4656.36,7);
    185 #endif
    186     PUSH_MPS("lseu",28,89,1120,834.68,7);
    187     PUSH_MPS("misc03",96,160,3360,1910.,7);
    188     PUSH_MPS("misc06",820,1808,12850.8607,12841.6,7);
    189 #if HOWMANY>1
    190     PUSH_MPS("misc07",212,260,2810,1415.0,7);
    191     PUSH_MPS("mitre",2054,10724,115155,114740.5184,7);
    192 #endif
    193     PUSH_MPS("mod008",6,319,307,290.9,7);
    194     PUSH_MPS("mod010",146,2655,6548,6532.08,7);
    195 #if HOWMANY>2
    196     PUSH_MPS("mod011",4480,10958,-54558535,-62121982.55,7);
    197     PUSH_MPS("modglob",291,422,20740508,20430947.,7);
    198 #endif
    199 #if HOWMANY>3
    200     PUSH_MPS("noswot",182,128,-43,-43.0,7);
    201 #endif
    202 #if HOWMANY>1
    203     PUSH_MPS("nw04",36,87482,16862,16310.66667,7);
    204 #endif
    205     PUSH_MPS("p0033",16,33,3089,2520.57,7);
    206     PUSH_MPS("p0201",133,201,7615,6875.0,7);
    207     PUSH_MPS("p0282",241,282,258411,176867.50,7);
    208     PUSH_MPS("p0548",176,548,8691,315.29,7);
    209     PUSH_MPS("p2756",755,2756,3124,2688.75,7);
    210 #if HOWMANY>2
    211     PUSH_MPS("pk1",45,86,11.0,0.0,7);
    212 #endif
    213 #if HOWMANY>1
    214     PUSH_MPS("pp08a",136,240,7350.0,2748.3452381,7);
    215     PUSH_MPS("pp08aCUTS",246,240,7350.0,5480.6061563,7);
    216 #endif
    217 #if HOWMANY>2
    218     PUSH_MPS("qiu",1192,840,-132.873137,-931.638857,7);
    219 #endif
    220     PUSH_MPS("qnet1",503,1541,16029.692681,14274.102667,7);
    221     PUSH_MPS("qnet1_o",456,1541,16029.692681,12095.571667,7);
    222     PUSH_MPS("rentacar",6803,9557,30356761,28806137.644,7);
    223     PUSH_MPS("rgn",24,180,82.1999,48.7999,7);
    224 #if HOWMANY>3
    225     PUSH_MPS("rout",291,556,1077.56,981.86428571,7);
    226     PUSH_MPS("set1ch",492,712,54537.75,32007.73,7);
    227 #endif
    228     //    PUSH_MPS("seymour",4944,1372,423,403.84647413,7);
    229     PUSH_MPS("stein27",118,27,18,13.0,7);
    230 #if HOWMANY>1
    231     PUSH_MPS("stein45",331,45,30,22.0,7);
    232 #endif
    233     PUSH_MPS("vpm1",234,378,20,15.4167,7);
    234     PUSH_MPS("vpm2",234,378,13.75,9.8892645972,7);
     151    PUSH_MPS("arki001",1048,1388,7580813.0459,7579599.80787,7);
     152    PUSH_MPS("bell3a",123,133,878430.32,862578.64,0);
     153    PUSH_MPS("bell5",91,104,8966406.49,8608417.95,1);
     154    PUSH_MPS("blend2",274,353,7.598985,6.9156751140,0);
     155    PUSH_MPS("cap6000",2176,6000,-2451377,-2451537.325,1);
     156    PUSH_MPS("dano3mip",3202,13873,728.1111,576.23162474,7);
     157    PUSH_MPS("danoint",664,521,65.67,62.637280418,6);
     158    PUSH_MPS("dcmulti",290,548,188182,183975.5397,0);
     159    PUSH_MPS("dsbmip",1182,1886,-305.19817501,-305.19817501,0);
     160    PUSH_MPS("egout",98,141,568.101,149.589,0);
     161    PUSH_MPS("enigma",21,100,0.0,0.0,0);
     162    PUSH_MPS("fast0507",507,63009,174,172.14556668,5);
     163    PUSH_MPS("fiber",363,1298,405935.18000,156082.51759,0);
     164    PUSH_MPS("fixnet6",478,878,3983,1200.88,1);
     165    PUSH_MPS("flugpl",18,18,1201500,1167185.7,0);
     166    PUSH_MPS("gen",780,870,112313,112130.0,0);
     167    PUSH_MPS("gesa2",1392,1224,25779856.372,25476489.678,1);
     168    PUSH_MPS("gesa2_o",1248,1224,25779856.372,25476489.678,1);
     169    PUSH_MPS("gesa3",1368,1152,27991042.648,27833632.451,0);
     170    PUSH_MPS("gesa3_o",1224,1152,27991042.648,27833632.451,0);
     171    PUSH_MPS("gt2",29,188,21166.000,13460.233074,0);
     172    PUSH_MPS("harp2",112,2993,-73899798.00,-74353341.502,6);
     173    PUSH_MPS("khb05250",101,1350,106940226,95919464.0,0);
     174    PUSH_MPS("l152lav",97,1989,4722,4656.36,1);
     175    PUSH_MPS("lseu",28,89,1120,834.68,0);
     176    PUSH_MPS("mas74",13,151,11801.18573,10482.79528,3);
     177    PUSH_MPS("mas76",12,151,40005.05414,38893.9036,2);
     178    PUSH_MPS("misc03",96,160,3360,1910.,0);
     179    PUSH_MPS("misc06",820,1808,12850.8607,12841.6,0);
     180    PUSH_MPS("misc07",212,260,2810,1415.0,1);
     181    PUSH_MPS("mitre",2054,10724,115155,114740.5184,1);
     182    PUSH_MPS("mkc",3411,5325,-553.75,-611.85,7); // this is suboptimal
     183    PUSH_MPS("mod008",6,319,307,290.9,0);
     184    PUSH_MPS("mod010",146,2655,6548,6532.08,0);
     185    PUSH_MPS("mod011",4480,10958,-54558535,-62121982.55,2);
     186    PUSH_MPS("modglob",291,422,20740508,20430947.,2);
     187    PUSH_MPS("noswot",182,128,-43,-43.0,6);
     188    PUSH_MPS("nw04",36,87482,16862,16310.66667,1);
     189    PUSH_MPS("p0033",16,33,3089,2520.57,0);
     190    PUSH_MPS("p0201",133,201,7615,6875.0,0);
     191    PUSH_MPS("p0282",241,282,258411,176867.50,0);
     192    PUSH_MPS("p0548",176,548,8691,315.29,0);
     193    PUSH_MPS("p2756",755,2756,3124,2688.75,0);
     194    PUSH_MPS("pk1",45,86,11.0,0.0,2);
     195    PUSH_MPS("pp08a",136,240,7350.0,2748.3452381,1);
     196    PUSH_MPS("pp08aCUTS",246,240,7350.0,5480.6061563,1);
     197    PUSH_MPS("qiu",1192,840,-132.873137,-931.638857,3);
     198    PUSH_MPS("qnet1",503,1541,16029.692681,14274.102667,0);
     199    PUSH_MPS("qnet1_o",456,1541,16029.692681,12095.571667,0);
     200    PUSH_MPS("rentacar",6803,9557,30356761,28806137.644,0);
     201    PUSH_MPS("rgn",24,180,82.1999,48.7999,0);
     202    PUSH_MPS("rout",291,556,1077.56,981.86428571,3);
     203    PUSH_MPS("set1ch",492,712,54537.75,32007.73,5);
     204    PUSH_MPS("seymour",4944,1372,423,403.84647413,7);
     205    PUSH_MPS("seymour_1",4944,1372,410.76370,403.84647413,5);
     206    PUSH_MPS("stein27",118,27,18,13.0,0);
     207    PUSH_MPS("stein45",331,45,30,22.0,1);
     208    PUSH_MPS("swath",884,6805,497.603,334.4968581,7);
     209    PUSH_MPS("vpm1",234,378,20,15.4167,0);
     210    PUSH_MPS("vpm2",234,378,13.75,9.8892645972,0);
    235211#endif
    236212  }
     
    244220#endif
    245221  int numberFailures=0;
     222  int numberAttempts=0;
    246223 
    247224  /*
     
    249226  */
    250227  for (m = 0 ; m < mpsName.size() ; m++) {
    251     std::cout << "  processing mps file: " << mpsName[m]
    252               << " (" << m+1 << " out of " << mpsName.size() << ")\n";
    253     /*
    254       Stage 1: Read the MPS
    255       and make sure the size of the constraint matrix is correct.
    256     */
    257     CbcModel * model = new CbcModel(saveModel);
     228    int test = testSet[m];
     229    if(testSwitch>=test&&loSwitch<=test) {
     230      numberAttempts++;
     231      std::cout << "  processing mps file: " << mpsName[m]
     232                << " (" << m+1 << " out of " << mpsName.size() << ")\n";
     233      /*
     234        Stage 1: Read the MPS
     235        and make sure the size of the constraint matrix is correct.
     236      */
     237      CbcModel * model = new CbcModel(saveModel);
    258238     
    259     std::string fn = dirMiplib+mpsName[m] ;
    260     if (!CbcTestMpsFile(fn)) {
    261       std::cout << "ERROR: Cannot find MPS file " << fn << "\n";
    262       continue;
    263     }
    264     CoinDrand48(true,1234567);
    265     //printf("RAND1 %g %g\n",CoinDrand48(true,1234567),model->randomNumberGenerator()->randomDouble());
    266     //printf("RAND1 %g\n",CoinDrand48(true,1234567));
    267     model->solver()->readMps(fn.c_str(),"") ;
    268     assert(model->getNumRows() == nRows[m]) ;
    269     assert(model->getNumCols() == nCols[m]) ;
    270 
    271     /*
    272       Stage 2: Call solver to solve the problem.
    273       then check the return code and objective.
    274     */
    275 
     239      std::string fn = dirMiplib+mpsName[m] ;
     240      if (!CbcTestMpsFile(fn)) {
     241        std::cout << "ERROR: Cannot find MPS file " << fn << "\n";
     242        continue;
     243      }
     244      CoinDrand48(true,1234567);
     245      //printf("RAND1 %g %g\n",CoinDrand48(true,1234567),model->randomNumberGenerator()->randomDouble());
     246      //printf("RAND1 %g\n",CoinDrand48(true,1234567));
     247      model->solver()->readMps(fn.c_str(),"") ;
     248      assert(model->getNumRows() == nRows[m]) ;
     249      assert(model->getNumCols() == nCols[m]) ;
     250     
     251      /*
     252        Stage 2: Call solver to solve the problem.
     253        then check the return code and objective.
     254      */
     255     
    276256#ifdef CLP_FACTORIZATION_INSTRUMENT
    277     extern double factorization_instrument(int type);
    278     double facTime1=factorization_instrument(0);
    279     printf("Factorization - initial solve %g seconds\n",
    280            facTime1);
    281     timeTakenFac += facTime1;
    282 #endif
    283     double startTime = CoinCpuTime()+CoinCpuTimeJustChildren();
    284     if (model->getMaximumNodes()>200000) {
    285       model->setMaximumNodes(200000);
    286     }
    287     OsiClpSolverInterface * si =
    288       dynamic_cast<OsiClpSolverInterface *>(model->solver()) ;
    289     ClpSimplex * modelC = NULL;
    290     if (si) {
    291       // get clp itself
    292       modelC = si->getModelPtr();
    293       if (stuff&&stuff[9]) {
    294         // vector matrix!
     257      extern double factorization_instrument(int type);
     258      double facTime1=factorization_instrument(0);
     259      printf("Factorization - initial solve %g seconds\n",
     260             facTime1);
     261      timeTakenFac += facTime1;
     262#endif
     263      double startTime = CoinCpuTime()+CoinCpuTimeJustChildren();
     264      int testMaximumNodes=200000;
     265      if (testSwitch>1)
     266        testMaximumNodes=20000000;
     267      if (model->getMaximumNodes()>testMaximumNodes) {
     268        model->setMaximumNodes(testMaximumNodes);
     269      }
     270      OsiClpSolverInterface * si =
     271        dynamic_cast<OsiClpSolverInterface *>(model->solver()) ;
     272      ClpSimplex * modelC = NULL;
     273      if (si) {
     274        // get clp itself
     275        modelC = si->getModelPtr();
    295276        ClpMatrixBase * matrix = modelC->clpMatrix();
    296         if (dynamic_cast< ClpPackedMatrix*>(matrix)) {
    297           ClpPackedMatrix * clpMatrix = dynamic_cast< ClpPackedMatrix*>(matrix);
     277        ClpPackedMatrix * clpMatrix = dynamic_cast< ClpPackedMatrix*>(matrix);
     278        if (stuff&&stuff[9]&&clpMatrix) {
     279          // vector matrix!
    298280          clpMatrix->makeSpecialColumnCopy();
    299281        }
    300       }
    301       model->checkModel();
    302       modelC->tightenPrimalBounds(0.0,0,true);
    303       model->initialSolve();
    304       if (modelC->dualBound()==1.0e10) {
    305         // user did not set - so modify
    306         // get largest scaled away from bound
    307         ClpSimplex temp=*modelC;
    308         temp.dual(0,7);
    309         double largestScaled=1.0e-12;
    310         double largest=1.0e-12;
    311         int numberRows = temp.numberRows();
    312         const double * rowPrimal = temp.primalRowSolution();
    313         const double * rowLower = temp.rowLower();
    314         const double * rowUpper = temp.rowUpper();
    315         const double * rowScale = temp.rowScale();
    316         int iRow;
    317         for (iRow=0;iRow<numberRows;iRow++) {
    318           double value = rowPrimal[iRow];
    319           double above = value-rowLower[iRow];
    320           double below = rowUpper[iRow]-value;
    321           if (above<1.0e12) {
    322             largest = CoinMax(largest,above);
    323           }
    324           if (below<1.0e12) {
    325             largest = CoinMax(largest,below);
    326           }
    327           if (rowScale) {
    328             double multiplier = rowScale[iRow];
    329             above *= multiplier;
    330             below *= multiplier;
    331           }
    332           if (above<1.0e12) {
    333             largestScaled = CoinMax(largestScaled,above);
    334           }
    335           if (below<1.0e12) {
    336             largestScaled = CoinMax(largestScaled,below);
    337           }
    338         }
    339        
    340         int numberColumns = temp.numberColumns();
    341         const double * columnPrimal = temp.primalColumnSolution();
    342         const double * columnLower = temp.columnLower();
    343         const double * columnUpper = temp.columnUpper();
    344         const double * columnScale = temp.columnScale();
    345         int iColumn;
    346         for (iColumn=0;iColumn<numberColumns;iColumn++) {
    347           double value = columnPrimal[iColumn];
    348           double above = value-columnLower[iColumn];
    349           double below = columnUpper[iColumn]-value;
    350           if (above<1.0e12) {
    351             largest = CoinMax(largest,above);
    352           }
    353           if (below<1.0e12) {
    354             largest = CoinMax(largest,below);
    355           }
    356           if (columnScale) {
    357             double multiplier = 1.0/columnScale[iColumn];
    358             above *= multiplier;
    359             below *= multiplier;
    360           }
    361           if (above<1.0e12) {
    362             largestScaled = CoinMax(largestScaled,above);
    363           }
    364           if (below<1.0e12) {
    365             largestScaled = CoinMax(largestScaled,below);
    366           }
    367         }
    368         std::cout<<"Largest (scaled) away from bound "<<largestScaled
    369                  <<" unscaled "<<largest<<std::endl;
    370282#if 0
    371         modelC->setDualBound(CoinMax(1.0001e8,
    372                                      CoinMin(1000.0*largestScaled,1.00001e10)));
     283        if (clpMatrix) {
     284          int numberRows = clpMatrix->getNumRows();
     285          int numberColumns = clpMatrix->getNumCols();
     286          double * elements = clpMatrix->getMutableElements();
     287          const int * row = clpMatrix->getIndices();
     288          const CoinBigIndex * columnStart = clpMatrix->getVectorStarts();
     289          const int * columnLength = clpMatrix->getVectorLengths();
     290          double * smallest = new double [numberRows];
     291          double * largest = new double [numberRows];
     292          char * flag = new char [numberRows];
     293          CoinZeroN(flag,numberRows);
     294          for (int i=0;i<numberRows;i++) {
     295            smallest[i]=COIN_DBL_MAX;
     296            largest[i]=0.0;
     297          }
     298          for (int iColumn=0;iColumn<numberColumns;iColumn++) {
     299            bool isInteger = modelC->isInteger(iColumn);
     300            CoinBigIndex j;
     301            for (j=columnStart[iColumn];
     302                 j<columnStart[iColumn]+columnLength[iColumn];j++) {
     303              int iRow=row[j];
     304              double value = fabs(elements[j]);
     305              if (!isInteger)
     306                flag[iRow]=1;
     307              smallest[iRow]=CoinMin(smallest[iRow],value);
     308              largest[iRow]=CoinMax(largest[iRow],value);
     309            }
     310          }
     311          double * rowLower = modelC->rowLower();
     312          double * rowUpper = modelC->rowUpper();
     313          bool changed=false;
     314          for (int i=0;i<numberRows;i++) {
     315            if (flag[i]&&smallest[i]>10.0&&false) {
     316              smallest[i]=1.0/smallest[i];
     317              if (rowLower[i]>-1.0e20)
     318                rowLower[i] *= smallest[i];
     319              if (rowUpper[i]<1.0e20)
     320                rowUpper[i] *= smallest[i];
     321              changed = true;
     322            } else {
     323              smallest[i]=0.0;
     324            }
     325          }
     326          if (changed) {
     327            printf("SCALED\n");
     328            for (int iColumn=0;iColumn<numberColumns;iColumn++) {
     329              CoinBigIndex j;
     330              for (j=columnStart[iColumn];
     331                   j<columnStart[iColumn]+columnLength[iColumn];j++) {
     332                int iRow=row[j];
     333                if (smallest[iRow])
     334                  elements[j] *= smallest[iRow];
     335              }
     336            }
     337          }
     338          delete [] smallest;
     339          delete [] largest;
     340          delete [] flag;
     341        }
     342#endif
     343        model->checkModel();
     344        modelC->tightenPrimalBounds(0.0,0,true);
     345        model->initialSolve();
     346        if (modelC->dualBound()==1.0e10) {
     347          // user did not set - so modify
     348          // get largest scaled away from bound
     349          ClpSimplex temp=*modelC;
     350          temp.dual(0,7);
     351          double largestScaled=1.0e-12;
     352          double largest=1.0e-12;
     353          int numberRows = temp.numberRows();
     354          const double * rowPrimal = temp.primalRowSolution();
     355          const double * rowLower = temp.rowLower();
     356          const double * rowUpper = temp.rowUpper();
     357          const double * rowScale = temp.rowScale();
     358          int iRow;
     359          for (iRow=0;iRow<numberRows;iRow++) {
     360            double value = rowPrimal[iRow];
     361            double above = value-rowLower[iRow];
     362            double below = rowUpper[iRow]-value;
     363            if (above<1.0e12) {
     364              largest = CoinMax(largest,above);
     365            }
     366            if (below<1.0e12) {
     367              largest = CoinMax(largest,below);
     368            }
     369            if (rowScale) {
     370              double multiplier = rowScale[iRow];
     371              above *= multiplier;
     372              below *= multiplier;
     373            }
     374            if (above<1.0e12) {
     375              largestScaled = CoinMax(largestScaled,above);
     376            }
     377            if (below<1.0e12) {
     378              largestScaled = CoinMax(largestScaled,below);
     379            }
     380          }
     381         
     382          int numberColumns = temp.numberColumns();
     383          const double * columnPrimal = temp.primalColumnSolution();
     384          const double * columnLower = temp.columnLower();
     385          const double * columnUpper = temp.columnUpper();
     386          const double * columnScale = temp.columnScale();
     387          int iColumn;
     388          for (iColumn=0;iColumn<numberColumns;iColumn++) {
     389            double value = columnPrimal[iColumn];
     390            double above = value-columnLower[iColumn];
     391            double below = columnUpper[iColumn]-value;
     392            if (above<1.0e12) {
     393              largest = CoinMax(largest,above);
     394            }
     395            if (below<1.0e12) {
     396              largest = CoinMax(largest,below);
     397            }
     398            if (columnScale) {
     399              double multiplier = 1.0/columnScale[iColumn];
     400              above *= multiplier;
     401              below *= multiplier;
     402            }
     403            if (above<1.0e12) {
     404              largestScaled = CoinMax(largestScaled,above);
     405            }
     406            if (below<1.0e12) {
     407              largestScaled = CoinMax(largestScaled,below);
     408            }
     409          }
     410          std::cout<<"Largest (scaled) away from bound "<<largestScaled
     411                   <<" unscaled "<<largest<<std::endl;
     412#if 0
     413          modelC->setDualBound(CoinMax(1.0001e8,
     414                                       CoinMin(1000.0*largestScaled,1.00001e10)));
    373415#else
    374         modelC->setDualBound(CoinMax(1.0001e9,
    375                                      CoinMin(1000.0*largestScaled,1.0001e10)));
    376 #endif
    377       }
    378     }
    379     model->setMinimumDrop(CoinMin(5.0e-2,
    380                               fabs(model->getMinimizationObjValue())*1.0e-3+1.0e-4));
    381     if (CoinAbs(model->getMaximumCutPassesAtRoot())<=100) {
    382       if (model->getNumCols()<500) {
    383         model->setMaximumCutPassesAtRoot(-100); // always do 100 if possible
    384       } else if (model->getNumCols()<5000) {
    385         model->setMaximumCutPassesAtRoot(100); // use minimum drop
    386       } else {
    387         model->setMaximumCutPassesAtRoot(20);
    388       }
    389     }
    390     // If defaults then increase trust for small models
    391     if (model->numberStrong()==5&&model->numberBeforeTrust()==10) {
    392       int numberColumns = model->getNumCols();
    393       if (numberColumns<=50) {
    394         model->setNumberBeforeTrust(1000);
    395       } else if (numberColumns<=100) {
    396         model->setNumberBeforeTrust(100);
    397       } else if (numberColumns<=300) {
    398         model->setNumberBeforeTrust(50);
    399       }
    400     }
    401     //if (model->getNumCols()>=500) {
     416          modelC->setDualBound(CoinMax(1.0001e9,
     417                                       CoinMin(1000.0*largestScaled,1.0001e10)));
     418#endif
     419        }
     420      }
     421      model->setMinimumDrop(CoinMin(5.0e-2,
     422                                    fabs(model->getMinimizationObjValue())*1.0e-3+1.0e-4));
     423      if (CoinAbs(model->getMaximumCutPassesAtRoot())<=100) {
     424        if (model->getNumCols()<500) {
     425          model->setMaximumCutPassesAtRoot(-100); // always do 100 if possible
     426        } else if (model->getNumCols()<5000) {
     427          model->setMaximumCutPassesAtRoot(100); // use minimum drop
     428        } else {
     429          model->setMaximumCutPassesAtRoot(20);
     430        }
     431      }
     432      // If defaults then increase trust for small models
     433      if (model->numberStrong()==5&&model->numberBeforeTrust()==10) {
     434        int numberColumns = model->getNumCols();
     435        if (numberColumns<=50) {
     436          model->setNumberBeforeTrust(1000);
     437        } else if (numberColumns<=100) {
     438          model->setNumberBeforeTrust(100);
     439        } else if (numberColumns<=300) {
     440          model->setNumberBeforeTrust(50);
     441        }
     442      }
     443      //if (model->getNumCols()>=500) {
    402444      // switch off Clp stuff
    403445      //model->setFastNodeDepth(-1);
    404     //}
    405     if (model->getNumCols()==-2756) {
    406       // p2756
    407       std::string problemName ;
    408       model->solver()->getStrParam(OsiProbName,problemName) ;
    409       model->solver()->activateRowCutDebugger(problemName.c_str()) ;
     446      //}
     447      if (model->getNumCols()==-2756) {
     448        // p2756
     449        std::string problemName ;
     450        model->solver()->getStrParam(OsiProbName,problemName) ;
     451        model->solver()->activateRowCutDebugger(problemName.c_str()) ;
     452      }
     453      if (model->getNumCols()==-201) {
     454        // p201
     455        std::string problemName ;
     456        model->solver()->getStrParam(OsiProbName,problemName) ;
     457        model->solver()->activateRowCutDebugger(problemName.c_str()) ;
     458      }
     459      if (model->getNumCols()==-104) {
     460        // bell5
     461        std::string problemName ;
     462        model->solver()->getStrParam(OsiProbName,problemName) ;
     463        model->solver()->activateRowCutDebugger(problemName.c_str()) ;
     464      }
     465      if (model->getNumCols()==-548&&model->getNumRows()==176) {
     466        // p0548
     467        std::string problemName ;
     468        model->solver()->getStrParam(OsiProbName,problemName) ;
     469        model->solver()->activateRowCutDebugger(problemName.c_str()) ;
     470      }
     471      if (model->getNumCols()==-160) {
     472        // misc03
     473        std::string problemName ;
     474        model->solver()->getStrParam(OsiProbName,problemName) ;
     475        model->solver()->activateRowCutDebugger(problemName.c_str()) ;
     476      }
     477      if (model->getNumCols()==-353) {
     478        // blend2
     479        std::string problemName ;
     480        model->solver()->getStrParam(OsiProbName,problemName) ;
     481        model->solver()->activateRowCutDebugger(problemName.c_str()) ;
     482      }
     483      if (model->getNumCols()==-100&&model->getNumRows()==21) {
     484        // enigma
     485        std::string problemName ;
     486        model->solver()->getStrParam(OsiProbName,problemName) ;
     487        model->solver()->activateRowCutDebugger(problemName.c_str()) ;
     488      }
     489      if (model->getNumCols()==-1541) {
     490        // qnet1
     491        std::string problemName ;
     492        model->solver()->getStrParam(OsiProbName,problemName) ;
     493        model->solver()->activateRowCutDebugger(problemName.c_str()) ;
     494      }
     495      if (model->getNumCols()==-10724) {
     496        // mitre
     497        std::string problemName ;
     498        model->solver()->getStrParam(OsiProbName,problemName) ;
     499        model->solver()->activateRowCutDebugger(problemName.c_str()) ;
     500      }
     501      if (model->getNumCols()==-1224) {
     502        //PUSH_MPS("gesa2",1392,1224,25779856.372,25476489.678,7);
     503        // gesa2
     504        std::string problemName ;
     505        model->solver()->getStrParam(OsiProbName,problemName) ;
     506        model->solver()->activateRowCutDebugger(problemName.c_str()) ;
     507      }
     508      if (model->getNumCols()==-1152&&model->getNumRows()==1368) {
     509        //PUSH_MPS("gesa3",1368,1152,27991042.648,27833632.451,7);
     510        // gesa3
     511        std::string problemName ;
     512        model->solver()->getStrParam(OsiProbName,problemName) ;
     513        model->solver()->activateRowCutDebugger(problemName.c_str()) ;
     514      }
     515      if (model->getNumCols()==-1152&&model->getNumRows()==1224) {
     516        //PUSH_MPS("gesa3_o",1224,1152,27991042.648,27833632.451,7);
     517        // gesa3
     518        std::string problemName ;
     519        model->solver()->getStrParam(OsiProbName,problemName) ;
     520        model->solver()->activateRowCutDebugger(problemName.c_str()) ;
     521      }
     522      if (model->getNumCols()==-282) {
     523        //PUSH_MPS("p0282",241,282,258411,176867.50,7);
     524        // p0282
     525        std::string problemName ;
     526        model->solver()->getStrParam(OsiProbName,problemName) ;
     527        model->solver()->activateRowCutDebugger(problemName.c_str()) ;
     528      }
     529      if (model->getNumCols()==-141) {
     530        // egout
     531        std::string problemName ;
     532        model->solver()->getStrParam(OsiProbName,problemName) ;
     533        model->solver()->activateRowCutDebugger(problemName.c_str()) ;
     534      }
     535      if (model->getNumCols()==-378) {
     536        // vpm2
     537        std::string problemName ;
     538        model->solver()->getStrParam(OsiProbName,problemName) ;
     539        model->solver()->activateRowCutDebugger(problemName.c_str()) ;
     540      }
     541      if (model->getNumCols()==-240&&model->getNumRows()==246) {
     542        // pp08aCUTS
     543        std::string problemName ;
     544        model->solver()->getStrParam(OsiProbName,problemName) ;
     545        model->solver()->activateRowCutDebugger(problemName.c_str()) ;
     546      }
     547      if (model->getNumCols()==-240&&model->getNumRows()==136) {
     548        // pp08a
     549        std::string problemName ;
     550        model->solver()->getStrParam(OsiProbName,problemName) ;
     551        model->solver()->activateRowCutDebugger(problemName.c_str()) ;
     552      }
     553      if (model->getNumCols()==-1372&&model->getNumRows()==4944) {
     554        // seymour1
     555        std::string problemName ;
     556        model->solver()->getStrParam(OsiProbName,problemName) ;
     557        model->solver()->activateRowCutDebugger(problemName.c_str()) ;
     558      }
     559      setCutAndHeuristicOptions(*model);
     560      if (si) {
     561#ifdef CLP_MULTIPLE_FACTORIZATIONS
     562        if (!modelC->factorization()->isDenseOrSmall()) {
     563          int denseCode = stuff ? static_cast<int> (stuff[4]) : -1;
     564          int smallCode = stuff ? static_cast<int> (stuff[10]) : -1;
     565          if (stuff&&stuff[8]>=1) {
     566            if (denseCode<0)
     567              denseCode=40;
     568            if (smallCode<0)
     569              smallCode=40;
     570          }
     571          if (denseCode>0)
     572            modelC->factorization()->setGoDenseThreshold(denseCode);
     573          if (smallCode>0)
     574            modelC->factorization()->setGoSmallThreshold(smallCode);
     575          if (denseCode>=modelC->numberRows()) {
     576            //printf("problem going dense\n");
     577            //modelC->factorization()->goDenseOrSmall(modelC->numberRows());
     578          }
     579        }
     580#endif
     581        if (stuff&&stuff[8]>=1) {
     582          if (modelC->numberColumns()+modelC->numberRows()<=10000&&
     583              model->fastNodeDepth()==-1)
     584            model->setFastNodeDepth(-9);
     585        }
     586      }
     587      //OsiObject * obj = new CbcBranchToFixLots(model,0.3,0.0,3,3000003);
     588      //model->addObjects(1,&obj);
     589      //delete obj;
     590      model->branchAndBound();
     591#ifdef CLP_FACTORIZATION_INSTRUMENT
     592      double facTime=factorization_instrument(0);
     593      printf("Factorization %g seconds\n",
     594             facTime);
     595      timeTakenFac += facTime;
     596#endif
     597     
     598      double timeOfSolution = CoinCpuTime()+CoinCpuTimeJustChildren()-startTime;
     599      // Print more statistics
     600      std::cout<<"Cuts at root node changed objective from "<<model->getContinuousObjective()
     601               <<" to "<<model->rootObjectiveAfterCuts()<<std::endl;
     602      int numberGenerators = model->numberCutGenerators();
     603      for (int iGenerator=0;iGenerator<numberGenerators;iGenerator++) {
     604        CbcCutGenerator * generator = model->cutGenerator(iGenerator);
     605#ifndef CLP_INVESTIGATE
     606        CglImplication * implication = dynamic_cast<CglImplication*>(generator->generator());
     607        if (implication)
     608          continue;
     609#endif
     610        std::cout<<generator->cutGeneratorName()<<" was tried "
     611                 <<generator->numberTimesEntered()<<" times and created "
     612                 <<generator->numberCutsInTotal()<<" cuts of which "
     613                 <<generator->numberCutsActive()<<" were active after adding rounds of cuts";
     614        if (generator->timing())
     615          std::cout<<" ( "<<generator->timeInCutGenerator()<<" seconds)"<<std::endl;
     616        else
     617          std::cout<<std::endl;
     618      }
     619      if (!model->status()) {
     620        double soln = model->getObjValue(); 
     621        double tolerance = CoinMax(1.0e-5,1.0e-5*CoinMin(fabs(soln),fabs(objValue[m])));
     622        //CoinRelFltEq eq(1.0e-3) ;
     623        if (fabs(soln-objValue[m])<tolerance) {
     624          std::cout
     625            <<"cbc_clp"<<" "
     626            << soln << " = " << objValue[m] << " ; okay";
     627          numProbSolved++;
     628        } else  {
     629          std::cout <<"cbc_clp" <<" " <<soln << " != " <<objValue[m]
     630                    << "; error=" << fabs(objValue[m] - soln);
     631          numberFailures++;
     632          //#ifdef COIN_DEVELOP
     633          //abort();
     634          //#endif
     635        }
     636      } else {
     637        std::cout << "cbc_clp error; too many nodes" ;
     638      }
     639      timeTaken += timeOfSolution;
     640      std::cout<<" - took " <<timeOfSolution<<" seconds.("<<
     641        model->getNodeCount()<<" / "<<model->getIterationCount()<<
     642        " ) subtotal "<<timeTaken
     643               <<" ("<<mpsName[m]<<")"<<std::endl;
     644      delete model;
    410645    }
    411     if (model->getNumCols()==-201) {
    412       // p201
    413       std::string problemName ;
    414       model->solver()->getStrParam(OsiProbName,problemName) ;
    415       model->solver()->activateRowCutDebugger(problemName.c_str()) ;
    416     }
    417     if (model->getNumCols()==-104) {
    418       // bell5
    419       std::string problemName ;
    420       model->solver()->getStrParam(OsiProbName,problemName) ;
    421       model->solver()->activateRowCutDebugger(problemName.c_str()) ;
    422     }
    423     if (model->getNumCols()==-548&&model->getNumRows()==176) {
    424       // p0548
    425       std::string problemName ;
    426       model->solver()->getStrParam(OsiProbName,problemName) ;
    427       model->solver()->activateRowCutDebugger(problemName.c_str()) ;
    428     }
    429     if (model->getNumCols()==-160) {
    430       // misc03
    431       std::string problemName ;
    432       model->solver()->getStrParam(OsiProbName,problemName) ;
    433       model->solver()->activateRowCutDebugger(problemName.c_str()) ;
    434     }
    435     if (model->getNumCols()==-353) {
    436       // blend2
    437       std::string problemName ;
    438       model->solver()->getStrParam(OsiProbName,problemName) ;
    439       model->solver()->activateRowCutDebugger(problemName.c_str()) ;
    440     }
    441     if (model->getNumCols()==-100&&model->getNumRows()==21) {
    442       // enigma
    443       std::string problemName ;
    444       model->solver()->getStrParam(OsiProbName,problemName) ;
    445       model->solver()->activateRowCutDebugger(problemName.c_str()) ;
    446     }
    447     if (model->getNumCols()==-1541) {
    448       // qnet1
    449       std::string problemName ;
    450       model->solver()->getStrParam(OsiProbName,problemName) ;
    451       model->solver()->activateRowCutDebugger(problemName.c_str()) ;
    452     }
    453     if (model->getNumCols()==-10724) {
    454       // mitre
    455       std::string problemName ;
    456       model->solver()->getStrParam(OsiProbName,problemName) ;
    457       model->solver()->activateRowCutDebugger(problemName.c_str()) ;
    458     }
    459     if (model->getNumCols()==-1224) {
    460       //PUSH_MPS("gesa2",1392,1224,25779856.372,25476489.678,7);
    461       // gesa2
    462       std::string problemName ;
    463       model->solver()->getStrParam(OsiProbName,problemName) ;
    464       model->solver()->activateRowCutDebugger(problemName.c_str()) ;
    465     }
    466     if (model->getNumCols()==-1152) {
    467       //PUSH_MPS("gesa3",1368,1152,27991042.648,27833632.451,7);
    468       // gesa3
    469       std::string problemName ;
    470       model->solver()->getStrParam(OsiProbName,problemName) ;
    471       model->solver()->activateRowCutDebugger(problemName.c_str()) ;
    472     }
    473     if (model->getNumCols()==-282) {
    474       //PUSH_MPS("p0282",241,282,258411,176867.50,7);
    475       // p0282
    476       std::string problemName ;
    477       model->solver()->getStrParam(OsiProbName,problemName) ;
    478       model->solver()->activateRowCutDebugger(problemName.c_str()) ;
    479     }
    480     if (model->getNumCols()==-141) {
    481       // egout
    482       std::string problemName ;
    483       model->solver()->getStrParam(OsiProbName,problemName) ;
    484       model->solver()->activateRowCutDebugger(problemName.c_str()) ;
    485     }
    486     if (model->getNumCols()==-378) {
    487       // vpm2
    488       std::string problemName ;
    489       model->solver()->getStrParam(OsiProbName,problemName) ;
    490       model->solver()->activateRowCutDebugger(problemName.c_str()) ;
    491     }
    492     if (model->getNumCols()==-240&&model->getNumRows()==246) {
    493       // pp08aCUTS
    494       std::string problemName ;
    495       model->solver()->getStrParam(OsiProbName,problemName) ;
    496       model->solver()->activateRowCutDebugger(problemName.c_str()) ;
    497     }
    498     if (model->getNumCols()==-240&&model->getNumRows()==136) {
    499       // pp08a
    500       std::string problemName ;
    501       model->solver()->getStrParam(OsiProbName,problemName) ;
    502       model->solver()->activateRowCutDebugger(problemName.c_str()) ;
    503     }
    504     if (model->getNumCols()==-1372&&model->getNumRows()==4944) {
    505       // seymour1
    506       std::string problemName ;
    507       model->solver()->getStrParam(OsiProbName,problemName) ;
    508       model->solver()->activateRowCutDebugger(problemName.c_str()) ;
    509     }
    510     setCutAndHeuristicOptions(*model);
    511     if (si) {
    512 #ifdef CLP_MULTIPLE_FACTORIZATIONS
    513       if (!modelC->factorization()->isDenseOrSmall()) {
    514         int denseCode = stuff ? static_cast<int> (stuff[4]) : -1;
    515         int smallCode = stuff ? static_cast<int> (stuff[10]) : -1;
    516         if (stuff&&stuff[8]>=1) {
    517           if (denseCode<0)
    518             denseCode=40;
    519           if (smallCode<0)
    520             smallCode=40;
    521         }
    522         if (denseCode>0)
    523           modelC->factorization()->setGoDenseThreshold(denseCode);
    524         if (smallCode>0)
    525           modelC->factorization()->setGoSmallThreshold(smallCode);
    526         if (denseCode>=modelC->numberRows()) {
    527           //printf("problem going dense\n");
    528           //modelC->factorization()->goDenseOrSmall(modelC->numberRows());
    529         }
    530       }
    531 #endif
    532       if (stuff&&stuff[8]>=1) {
    533         if (modelC->numberColumns()+modelC->numberRows()<=100000&&
    534             model->fastNodeDepth()==-1)
    535           model->setFastNodeDepth(-9);
    536       }
    537     }
    538     //OsiObject * obj = new CbcBranchToFixLots(model,0.3,0.0,3,3000003);
    539     //model->addObjects(1,&obj);
    540     //delete obj;
    541     model->branchAndBound();
    542 #ifdef CLP_FACTORIZATION_INSTRUMENT
    543     double facTime=factorization_instrument(0);
    544     printf("Factorization %g seconds\n",
    545            facTime);
    546     timeTakenFac += facTime;
    547 #endif
    548      
    549     double timeOfSolution = CoinCpuTime()+CoinCpuTimeJustChildren()-startTime;
    550     // Print more statistics
    551     std::cout<<"Cuts at root node changed objective from "<<model->getContinuousObjective()
    552              <<" to "<<model->rootObjectiveAfterCuts()<<std::endl;
    553     int numberGenerators = model->numberCutGenerators();
    554     for (int iGenerator=0;iGenerator<numberGenerators;iGenerator++) {
    555       CbcCutGenerator * generator = model->cutGenerator(iGenerator);
    556 #ifndef CLP_INVESTIGATE
    557       CglImplication * implication = dynamic_cast<CglImplication*>(generator->generator());
    558       if (implication)
    559         continue;
    560 #endif
    561       std::cout<<generator->cutGeneratorName()<<" was tried "
    562                <<generator->numberTimesEntered()<<" times and created "
    563                <<generator->numberCutsInTotal()<<" cuts of which "
    564                <<generator->numberCutsActive()<<" were active after adding rounds of cuts";
    565       if (generator->timing())
    566         std::cout<<" ( "<<generator->timeInCutGenerator()<<" seconds)"<<std::endl;
    567       else
    568         std::cout<<std::endl;
    569     }
    570     if (!model->status()) {
    571       double soln = model->getObjValue(); 
    572       double tolerance = CoinMax(1.0e-5,1.0e-5*CoinMin(fabs(soln),fabs(objValue[m])));
    573       //CoinRelFltEq eq(1.0e-3) ;
    574       if (fabs(soln-objValue[m])<tolerance) {
    575         std::cout
    576           <<"cbc_clp"<<" "
    577           << soln << " = " << objValue[m] << " ; okay";
    578         numProbSolved++;
    579       } else  {
    580         std::cout <<"cbc_clp" <<" " <<soln << " != " <<objValue[m]
    581                   << "; error=" << fabs(objValue[m] - soln);
    582         numberFailures++;
    583         //#ifdef COIN_DEVELOP
    584         //abort();
    585         //#endif
    586       }
    587     } else {
    588       std::cout << "error; too many nodes" ;
    589     }
    590     timeTaken += timeOfSolution;
    591     std::cout<<" - took " <<timeOfSolution<<" seconds.("<<
    592       model->getNodeCount()<<" / "<<model->getIterationCount()<<
    593       " ) subtotal "<<timeTaken<<std::endl;
    594     delete model;
    595646  }
    596647  int returnCode=0;
     
    600651    <<numProbSolved
    601652    <<" out of "
    602     <<objValue.size();
    603   int numberOnNodes = objValue.size()-numProbSolved-numberFailures;
     653    <<numberAttempts;
     654  int numberOnNodes = numberAttempts-numProbSolved-numberFailures;
    604655  if (numberFailures||numberOnNodes) {
    605656    if (numberOnNodes) {
     
    616667    <<" seconds."
    617668    <<std::endl;
    618   if (unitTestOnly) {
     669  if (testSwitch==-2) {
    619670    if(numberFailures||numberOnNodes) {
    620671      printf("****** Unit Test failed\n");
Note: See TracChangeset for help on using the changeset viewer.