Changeset 1569 for trunk


Ignore:
Timestamp:
Jan 3, 2011 9:00:42 AM (9 years ago)
Author:
forrest
Message:

more sos heuristics

Location:
trunk/Cbc/src
Files:
3 edited

Legend:

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

    r1514 r1569  
    10031003                    if (!gotPump) {
    10041004                        CbcHeuristicFPump heuristic4;
     1005                        // use any cutoff
     1006                        heuristic4.setFakeCutoff(0.5*COIN_DBL_MAX);
    10051007                        if (fractionSmall_<=1.0)
    10061008                          heuristic4.setMaximumPasses(10);
  • trunk/Cbc/src/CbcHeuristicRENS.cpp

    r1564 r1569  
    9090
    9191    const double * currentSolution = solver->getColSolution();
     92    OsiSolverInterface * newSolver = cloneBut(3); // was model_->continuousSolver()->clone();
     93    newSolver->resolve();
     94    double direction = newSolver->getObjSense();
     95    double cutoff=model_->getCutoff();
     96    //newSolver->getDblParam(OsiDualObjectiveLimit, cutoff);
     97    //cutoff *= direction;
     98    double gap = cutoff - newSolver->getObjValue() * direction ;
     99    double tolerance;
     100    newSolver->getDblParam(OsiDualTolerance, tolerance) ;
     101    if (gap > 0.0 || !newSolver->isProvenOptimal()) {
     102      gap += 100.0 * tolerance;
     103      int nFix = newSolver->reducedCostFix(gap);
     104      if (nFix) {
     105        char line [200];
     106        sprintf(line, "Reduced cost fixing fixed %d variables", nFix);
     107        model_->messageHandler()->message(CBC_FPUMP1, model_->messages())
     108          << line
     109          << CoinMessageEol;
     110      }
     111    } else {
     112      return 0; // finished?
     113    }
    92114    int numberColumns = solver->getNumCols();
    93115    double * dj = CoinCopyOfArray(solver->getReducedCost(),numberColumns);
    94     OsiSolverInterface * newSolver = cloneBut(3); // was model_->continuousSolver()->clone();
    95     double direction = newSolver->getObjSense();
    96116    int type = rensType_&15;
    97117    double djTolerance = (type!=1) ? -1.0e30 : 1.0e-4;
     
    99119    const double * colUpper = newSolver->getColUpper();
    100120    int numberFixed = 0;
    101     if ((type&7)==3) {
     121    if (type==3) {
    102122      double total=0.0;
    103123      int n=0;
     
    117137        delete basis;
    118138      }
    119     } else if ((type&7)==4) {
    120       double * sort = new double [numberColumns];
    121       for (int iColumn = 0; iColumn < numberColumns; iColumn++) {
    122         sort[iColumn]=1.0e30;
    123         if (colUpper[iColumn]>colLower[iColumn]) {
    124           sort[iColumn] = fabs(dj[iColumn]);
    125         }
    126       }
    127       std::sort(sort,sort+numberColumns);
    128       int last = static_cast<int>(numberColumns*fractionSmall_);
    129       djTolerance = sort[last];
    130       delete [] sort;
    131     } else if ((type&7)==5||(type&7)==6) {
     139    } else if (type>=5&&type<=10) {
     140      /* 5 fix sets at one
     141         6 fix on dj but leave unfixed SOS slacks
     142         7 fix sets at one but use pi
     143         8 fix all at zero but leave unfixed SOS slacks
     144         9 as 8 but only fix all at zero if just one in set nonzero
     145         10 as 7 but pi other way
     146      */
    132147      // SOS type fixing
    133       bool fixSets = (type&7)==5;
     148      bool fixSets = (type==5)||(type==7)||(type==10);
    134149      CoinWarmStartBasis * basis =
    135150        dynamic_cast<CoinWarmStartBasis *>(solver->getWarmStart()) ;
     
    160175                   j < columnStart[iColumn] + columnLength[iColumn]; j++) {
    161176                int iRow = row[j];
    162                 bestDj[iRow]=1.0e30;
     177                if (bestDj[iRow]<1.0e30) {
     178                  if (element[j] != 1.0)
     179                    bestDj[iRow]=1.0e30;
     180                  else
     181                    bestDj[iRow]=1.0e25;
     182                }
    163183              }
    164184            } else if ( basis->getStructStatus(iColumn) !=
     
    167187                   j < columnStart[iColumn] + columnLength[iColumn]; j++) {
    168188                int iRow = row[j];
    169                 if (bestDj[iRow]<1.0e30) {
     189                if (bestDj[iRow]<1.0e25) {
    170190                  if (element[j] != 1.0)
    171191                    bestDj[iRow]=1.0e30;
     
    179199        int nSOS=0;
    180200        double * sort = new double [numberRows];
     201        const double * pi = newSolver->getRowPrice();
    181202        for (int i=0;i<numberRows;i++) {
    182203          if (bestDj[i]<1.0e30) {
    183             sort[nSOS++]=bestDj[i];
     204            if (type==5)
     205              sort[nSOS++]=bestDj[i];
     206            else if (type==7)
     207              sort[nSOS++]=-fabs(pi[i]);
     208            else
     209              sort[nSOS++]=fabs(pi[i]);
    184210          }
    185211        }
     
    197223                       j < columnStart[iColumn] + columnLength[iColumn]; j++) {
    198224                    int iRow = row[j];
    199                     if (bestDj[iRow]<1.0e30&&bestDj[iRow]>=tolerance) {
     225                    double useDj;
     226                    if (type==5)
     227                      useDj = bestDj[iRow];
     228                    else if (type==7)
     229                      useDj= -fabs(pi[iRow]);
     230                    else
     231                      useDj= fabs(pi[iRow]);
     232                    if (bestDj[iRow]<1.0e30&&useDj>=tolerance) {
    200233                      numberFixed++;
    201234                      if (currentSolution[iColumn]<=1.0e-6)
     
    212245                    dj[iColumn] *= 0.000001;
    213246                  }
     247                } else if (type==8||type==9) {
     248                  if (currentSolution[iColumn]<=1.0e-6) {
     249                    if (type==8) {
     250                      dj[iColumn] *= 1.0e6;
     251                    } else {
     252                      bool fix=false;
     253                      for (j = columnStart[iColumn];
     254                           j < columnStart[iColumn] + columnLength[iColumn]; j++) {
     255                        int iRow = row[j];
     256                        if (bestDj[iRow]<1.0e25) {
     257                          fix=true;
     258                          break;
     259                        }
     260                      }
     261                      if (fix) {
     262                        dj[iColumn] *= 1.0e6;
     263                      }
     264                    }
     265                  } else {
     266                    dj[iColumn] *= 0.000001;
     267                  }
    214268                }
    215269              }
     
    218272          if (fixSets)
    219273            djTolerance = 1.0e30;
    220           else
    221             djTolerance = tolerance;
    222274        }
    223275        delete basis;
     
    225277        delete [] bestDj;
    226278      }
     279    }
     280    // Do dj to get right number
     281    if (type==4||type==6||type>7) {
     282      double * sort = new double [numberColumns];
     283      for (int iColumn = 0; iColumn < numberColumns; iColumn++) {
     284        sort[iColumn]=1.0e30;
     285        if (colUpper[iColumn]>colLower[iColumn]) {
     286          sort[iColumn] = fabs(dj[iColumn]);
     287        }
     288      }
     289      std::sort(sort,sort+numberColumns);
     290      int last = static_cast<int>(numberColumns*fractionSmall_);
     291      djTolerance = CoinMax(sort[last],1.0e-5);
     292      delete [] sort;
    227293    }
    228294   
  • trunk/Cbc/src/CbcModel.cpp

    r1521 r1569  
    21172117    }
    21182118    // Save objective (just so user can access it)
    2119     originalContinuousObjective_ = solver_->getObjValue();
     2119    originalContinuousObjective_ = solver_->getObjValue()* solver_->getObjSense();
    21202120    bestPossibleObjective_ = originalContinuousObjective_;
    21212121    sumChangeObjective1_ = 0.0;
     
    43064306    secondaryStatus_ = -1;
    43074307    originalContinuousObjective_ = solver_->getObjValue() * solver_->getObjSense();
     4308    bestPossibleObjective_ = originalContinuousObjective_;
    43084309    delete [] continuousSolution_;
    43094310    continuousSolution_ = CoinCopyOfArray(solver_->getColSolution(),
     
    1285812859                    // see if heuristic will do anything
    1285912860                    double saveValue = heuristicValue ;
     12861                    double before = getCurrentSeconds();
    1286012862                    int ifSol = heuristic_[i]->solution(heuristicValue,
    1286112863                                                        newSolution);
     12864                    if (handler_->logLevel()>1) {
     12865                      char line[100];
     12866                      sprintf(line,"Heuristic %s took %g seconds",
     12867                              heuristic_[i]->heuristicName(),
     12868                              getCurrentSeconds()-before);
     12869                      handler_->message(CBC_GENERAL, messages_) <<
     12870                        line << CoinMessageEol ;
     12871                    }
    1286212872                    if (ifSol > 0) {
    1286312873                        // better solution found
     
    1288012890                            if (heuristic_[i]->exitNow(bestObjective_))
    1288112891                                break;
     12892                            if (eventHandler) {
     12893                              if (!eventHandler->event(CbcEventHandler::heuristicSolution)) {
     12894                                eventHappened_ = true; // exit
     12895                                break;
     12896                              }
     12897                            }
     12898                            double testGap = CoinMax(dblParam_[CbcAllowableGap],
     12899                                                     CoinMax(fabs(bestObjective_), fabs(bestPossibleObjective_))
     12900                                                     * dblParam_[CbcAllowableFractionGap]);
     12901                            if (bestObjective_ - bestPossibleObjective_ < testGap && getCutoffIncrement() >= 0.0 &&bestPossibleObjective_ < 1.0e30) {
     12902                              if (bestPossibleObjective_ < getCutoff())
     12903                                stoppedOnGap_ = true ;
     12904                              //eventHappened_=true; // stop as fast as possible
     12905                              break;
     12906                            }
    1288212907                        } else {
    1288312908                            // NOT better solution
Note: See TracChangeset for help on using the changeset viewer.