Changeset 1127


Ignore:
Timestamp:
Jan 18, 2009 12:09:55 PM (11 years ago)
Author:
forrest
Message:

changes to fpump

Location:
trunk/Cbc/src
Files:
3 edited

Legend:

Unmodified
Added
Removed
  • trunk/Cbc/src/CbcHeuristic.hpp

    r1121 r1127  
    138138      1 bit - stop once allowable gap on objective reached
    139139      2 bit - always do given number of passes
    140       4 bit - increase cutoff by 5% every 500 passes?
     140      for other possibilities see switches_
    141141  */
    142142  inline void setSwitches(int value)
     
    145145      1 bit - stop once allowable gap on objective reached
    146146      2 bit - always do given number of passes
    147       4 bit - increase cutoff by 5% every 500 passes?
     147      for other possibilities see switches_
    148148  */
    149149  inline int switches() const
     
    240240  /// How much to increase how often
    241241  double decayFactor_;
    242   /** Switches on input (does not apply equally to all heuristics)
     242  /** Switches (does not apply equally to all heuristics)
    243243      1 bit - stop once allowable gap on objective reached
    244244      2 bit - always do given number of passes
    245       4 bit - increase cutoff by 5% every 500 passes?
     245      4 bit - weaken cutoff by 5% every 50 passes?
     246      8 bit - if has cutoff and suminf bobbling for 20 passes then
     247              first try halving distance to best possible then
     248              try keep halving distance to known cutoff
    246249  */
    247250  int switches_;
  • trunk/Cbc/src/CbcHeuristicFPump.cpp

    r1121 r1127  
    369369    }
    370370  }
     371  // guess exact multiple of objective
     372  double exactMultiple = model_->getCutoffIncrement();
     373  exactMultiple *= 2520;
     374  if (fabs(exactMultiple/0.999-floor(exactMultiple/0.999+0.5))<1.0e-9)
     375    exactMultiple /= 2520.0*0.999;
     376  else if (fabs(exactMultiple-floor(exactMultiple+0.5))<1.0e-9)
     377    exactMultiple /= 2520.0;
     378  else
     379    exactMultiple = 0.0;
     380  // check for rounding errors (only for integral case)
     381  if (fabs(exactMultiple-floor(exactMultiple+0.5))<1.0e-8)
     382    exactMultiple = floor(exactMultiple+0.5);
     383  //printf("exact multiple %g\n",exactMultiple);
    371384  while (!exitAll) {
     385    // Cutoff rhs
     386    double useRhs=COIN_DBL_MAX;
     387    double useOffset=0.0;
    372388    int numberPasses=0;
    373389    artificialFactor *= 10.0;
     
    411427    bool useCutoff = (fabs(cutoff)<1.0e20&&(fakeCutoff_!=COIN_DBL_MAX||numberTries>1));
    412428    // but there may be a close one
    413     if (firstCutoff<2.0*solutionValue&&numberTries==1/*&&CoinMin(cutoff,fakeCutoff_)<1.0e20*/)
     429    if (firstCutoff<2.0*solutionValue&&numberTries==1&&CoinMin(cutoff,fakeCutoff_)<1.0e20)
    414430      useCutoff=true;
    415431    if (useCutoff) {
     
    427443        }
    428444      }
    429       double offset;
    430       solver->getDblParam(OsiObjOffset,offset);
     445      solver->getDblParam(OsiObjOffset,useOffset);
    431446#ifdef COIN_DEVELOP
    432       if (offset)
    433         printf("CbcHeuristicFPump obj offset %g\n",offset);
     447      if (useOffset)
     448        printf("CbcHeuristicFPump obj offset %g\n",useOffset);
    434449#endif
    435       // Tweak rhs
    436       double useRhs = rhs+offset*direction;
     450      useOffset *= direction;
     451      // Tweak rhs and save
     452      useRhs = rhs;
    437453#if 0
    438454      double tempValue = 60.0*useRhs;
     
    442458      }
    443459#endif
    444       solver->addRow(nel,which,els,-COIN_DBL_MAX,useRhs);
     460      solver->addRow(nel,which,els,-COIN_DBL_MAX,useRhs+useOffset);
    445461      delete [] which;
    446462      delete [] els;
     
    521537             largestCost,nArtificial);
    522538#endif
     539    // This is an array of sums of infeasibilities so can see if "bobbling"
     540#define SIZE_BOBBLE 20
     541    double saveSumInf[SIZE_BOBBLE];
     542    CoinFillN(saveSumInf,SIZE_BOBBLE,COIN_DBL_MAX);
     543    // 0 before doing anything
     544    int bobbleMode = 0;
    523545    // 5. MAIN WHILE LOOP
    524546    //bool newLineNeeded=false;
     
    10391061            }
    10401062            newTrueSolutionValue *= direction;
    1041             if (newNumberInfeas&&newNumberInfeas<-20) {
    1042 #if 1
     1063            if (newNumberInfeas&&newNumberInfeas<15) {
     1064#if 0
    10431065              roundingObjective=solutionValue;
    10441066              OsiSolverInterface * saveSolver = model_->swapSolver(solver);
     
    10461068                CoinCopyOfArray(solver->getObjCoefficients(),numberColumns);
    10471069              solver->setObjective(saveObjective);
     1070              double saveOffset2;
     1071              solver->getDblParam(OsiObjOffset,saveOffset2);
     1072              solver->setDblParam(OsiObjOffset,saveOffset);
    10481073              int ifSol = roundingHeuristic.solution(roundingObjective,roundingSolution);
    10491074              solver->setObjective(currentObjective);
     1075              solver->setDblParam(OsiObjOffset,saveOffset2);
    10501076              delete [] currentObjective;
    10511077              model_->swapSolver(saveSolver);
     
    11311157                    newObj += contrib2-contrib;
    11321158                    newObj *= direction;
     1159#ifdef COIN_DEVELOP
    11331160                    printf("FFFeasible! - obj %g\n",newObj);
    1134                     if (newObj<roundingObjective) {
     1161#endif
     1162                    if (newObj<roundingObjective-1.0e-6) {
     1163#ifdef COIN_DEVELOP
    11351164                      printf("FBetter\n");
     1165#endif
    11361166                      roundingObjective = newObj;
    11371167                      memcpy(roundingSolution,newSolution,numberColumns*sizeof(double));
     
    11621192            }
    11631193          }
    1164           if (false) {
     1194          if (true) {
    11651195            OsiSolverInterface * saveSolver = model_->swapSolver(solver);
    11661196            double * currentObjective =
    11671197              CoinCopyOfArray(solver->getObjCoefficients(),numberColumns);
     1198            //double * saveSolution = CoinCopyOfArray(solution,numberColumns);
    11681199            solver->setObjective(saveObjective);
    11691200            double saveOffset2;
     
    11811212            solver->setObjective(currentObjective);
    11821213            solver->setDblParam(OsiObjOffset,saveOffset2);
     1214            //if ((switches_&16)==0||numberPasses<-10)
     1215            //solver->setColSolution(saveSolution);
     1216            //delete [] saveSolution;
    11831217            delete [] currentObjective;
    11841218            if (returnCode==1) {
     1219#ifdef COIN_DEVELOP
    11851220              printf("rounding obj of %g?\n",roundingObjective);
     1221#endif
    11861222              //roundingObjective = newSolutionValue;
    11871223            } else {
     
    13331369          closestObjectiveValue = solver->getObjValue();
    13341370        }
    1335         //newLineNeeded=true;
    1336        
     1371        // See if we need to think about changing rhs
     1372        if ((switches_&12)!=0&&useRhs<1.0e50) {
     1373          double oldRhs = useRhs;
     1374          bool trying=false;
     1375          if ((switches_&4)!=0&&numberPasses&&(numberPasses%50)==0) {
     1376            if (solutionValue>1.0e20) {
     1377              // only if no genuine solution
     1378              double gap = useRhs-continuousObjectiveValue;
     1379              useRhs += 0.1*gap;
     1380              if (exactMultiple) {
     1381                useRhs = exactMultiple*ceil(useRhs/exactMultiple);
     1382                useRhs = CoinMax(useRhs,oldRhs+exactMultiple);
     1383              }
     1384              trying=true;
     1385            }
     1386          }
     1387          if ((switches_&8)!=0) {
     1388            // Put in new suminf and check
     1389            double largest=newSumInfeas;
     1390            double smallest=newSumInfeas;
     1391            for (int i=0;i<SIZE_BOBBLE-1;i++) {
     1392              double value = saveSumInf[i+1];
     1393              saveSumInf[i]=value;
     1394              largest = CoinMax(largest,value);
     1395              smallest = CoinMin(smallest,value);
     1396            }
     1397            saveSumInf[SIZE_BOBBLE-1]=newSumInfeas;
     1398            if (smallest*1.5>largest&&smallest>2.0) {
     1399              if (bobbleMode==0) {
     1400                // go closer
     1401                double gap = oldRhs-continuousObjectiveValue;
     1402                useRhs -= 0.4*gap;
     1403                if (exactMultiple) {
     1404                  double value = floor(useRhs/exactMultiple);
     1405                  useRhs = CoinMin(value*exactMultiple,oldRhs-exactMultiple);
     1406                }
     1407                if (useRhs<continuousObjectiveValue) {
     1408                  // skip decrease
     1409                  bobbleMode=1;
     1410                  useRhs=oldRhs;
     1411                }
     1412              }
     1413              if (bobbleMode) {
     1414                trying=true;
     1415                // weaken
     1416                if (solutionValue<1.0e20) {
     1417                  double gap = solutionValue-oldRhs;
     1418                  useRhs += 0.3*gap;
     1419                } else {
     1420                  double gap = oldRhs-continuousObjectiveValue;
     1421                  useRhs += 0.05*gap;
     1422                }
     1423                if (exactMultiple) {
     1424                  double value = ceil(useRhs/exactMultiple);
     1425                  useRhs = CoinMin(value*exactMultiple,
     1426                                   solutionValue-exactMultiple);
     1427                }
     1428              }
     1429              bobbleMode++;
     1430              // reset
     1431              CoinFillN(saveSumInf,SIZE_BOBBLE,COIN_DBL_MAX);
     1432            }
     1433          }
     1434          if (useRhs!=oldRhs) {
     1435            // tidy up
     1436            if (exactMultiple) {
     1437              double value = floor(useRhs/exactMultiple);
     1438              double bestPossible = ceil(continuousObjectiveValue/exactMultiple);
     1439              useRhs = CoinMax(value,bestPossible)*exactMultiple;
     1440            } else {
     1441              useRhs = CoinMax(useRhs,continuousObjectiveValue);
     1442            }
     1443            int k = solver->getNumRows()-1;
     1444            solver->setRowUpper(k,useRhs+useOffset);
     1445            bool takeHint;
     1446            OsiHintStrength strength;
     1447            solver->getHintParam(OsiDoDualInResolve,takeHint,strength);
     1448            if (useRhs<oldRhs) {
     1449              solver->setHintParam(OsiDoDualInResolve,true);
     1450              solver->resolve();
     1451            } else if (useRhs>oldRhs) {
     1452              solver->setHintParam(OsiDoDualInResolve,false);
     1453              solver->resolve();
     1454            }
     1455            solver->setHintParam(OsiDoDualInResolve,takeHint);
     1456            if (!solver->isProvenOptimal()) {
     1457              // presumably max time or some such
     1458              exitAll=true;
     1459              break;
     1460            }
     1461          } else if (trying) {
     1462            // doesn't look good
     1463            break;
     1464          }
     1465        }
    13371466      }
    13381467      // reduce scale factor
     
    16551784      }
    16561785      // But round down
    1657       if (model_->getCutoffIncrement()==0.999||
    1658           model_->getCutoffIncrement()==1.0)
    1659         cutoff = floor(cutoff);
     1786      if (exactMultiple)
     1787        cutoff = exactMultiple*floor(cutoff/exactMultiple);
    16601788      if (cutoff<continuousObjectiveValue)
    16611789        break;
     
    18341962    }
    18351963  }
     1964  int numberColumns = solver->getNumCols();
    18361965#if 0
    18371966  // Do set covering variables
     
    18411970  const CoinBigIndex * rowStart = matrixByRow->getVectorStarts();
    18421971  const int * rowLength = matrixByRow->getVectorLengths();
    1843   int numberColumns = solver->getNumCols();
    18441972  double * sortTemp = new double[numberColumns];
    18451973  int * whichTemp = new int [numberColumns];
     
    19422070  if (i==numberIntegers) {
    19432071    // may be able to use solution even if 0.99999's
    1944     // get row activities
    1945     double * rowActivity = new double[numberRows];
    1946     memset(rowActivity,0,numberRows*sizeof(double));
    1947     solver->getMatrixByCol()->times(solution,rowActivity) ;
    1948     double largestInfeasibility =0.0;
    1949     for (i=0 ; i < numberRows ; i++) {
    1950       largestInfeasibility = CoinMax(largestInfeasibility,
    1951                                  rowLower[i]-rowActivity[i]);
    1952       largestInfeasibility = CoinMax(largestInfeasibility,
    1953                                  rowActivity[i]-rowUpper[i]);
    1954     }
    1955     delete [] rowActivity;
    1956     if (largestInfeasibility<=primalTolerance) {
     2072    double * saveLower = CoinCopyOfArray(columnLower,numberColumns);
     2073    double * saveUpper = CoinCopyOfArray(columnUpper,numberColumns);
     2074    double * saveSolution = CoinCopyOfArray(solution,numberColumns);
     2075    double * tempSolution = CoinCopyOfArray(solution,numberColumns);
     2076    CoinWarmStartBasis  * saveBasis =
     2077      dynamic_cast<CoinWarmStartBasis *>(solver->getWarmStart()) ;
     2078    for (i=0;i<numberIntegers;i++) {
     2079      int iColumn = integerVariable[i];
     2080      double value=solution[iColumn];
     2081      double round = floor(value+0.5);
     2082      solver->setColLower(iColumn,round);
     2083      solver->setColUpper(iColumn,round);
     2084      tempSolution[iColumn] = round;
     2085    }
     2086    solver->setColSolution(tempSolution);
     2087    delete [] tempSolution;
     2088    solver->resolve();
     2089    solver->setColLower(saveLower);
     2090    solver->setColUpper(saveUpper);
     2091    solver->setWarmStart(saveBasis);
     2092    delete [] saveLower;
     2093    delete [] saveUpper;
     2094    delete saveBasis;
     2095    if (!solver->isProvenOptimal()) {
     2096      solver->setColSolution(saveSolution);
     2097    }
     2098    delete [] saveSolution;
     2099    if (solver->isProvenOptimal()) {
    19572100      // feasible
    19582101      delete [] list;
  • trunk/Cbc/src/CbcSolver.cpp

    r1126 r1127  
    35303530    anyToDo=true;
    35313531  }
    3532   int heurSwitches=parameters_[whichParam(HOPTIONS,numberParameters_,parameters_)].intValue()/10;
     3532  int heurSwitches=parameters_[whichParam(HOPTIONS,numberParameters_,parameters_)].intValue()%100;
    35333533  if (heurSwitches) {
    35343534    for (int iHeur=0;iHeur<model->numberHeuristics();iHeur++) {
     
    35543554    bool cleanModel = !model2.numberIntegers()&&!model2.numberObjects();
    35553555    model2.findIntegers(false);
    3556     int heurOptions=parameters_[whichParam(HOPTIONS,numberParameters_,parameters_)].intValue()%10;
     3556    int heurOptions=parameters_[whichParam(HOPTIONS,numberParameters_,parameters_)].intValue()/100;
    35573557    if (heurOptions==0||heurOptions==2) {
    35583558      model2.doHeuristicsAtRoot(1);
     
    76837683                }
    76847684                int heurOptions=parameters_[whichParam(HOPTIONS,numberParameters_,parameters_)].intValue();
    7685                 if (heurOptions>1)
     7685                if (heurOptions>100)
    76867686                  babModel_->setSpecialOptions(babModel_->specialOptions()|8192);
    76877687               
Note: See TracChangeset for help on using the changeset viewer.