/home/coin/SVN-release/OS-2.4.0/Bonmin/src/Algorithms/BonSubMipSolver.cpp

Go to the documentation of this file.
00001 // (C) Copyright International Business Machines (IBM) 2006
00002 // All Rights Reserved.
00003 // This code is published under the Common Public License.
00004 //
00005 // Authors :
00006 // P. Bonami, International Business Machines
00007 //
00008 // Date :  12/07/2006
00009 
00010 
00011 // Code separated from BonOaDecBase to try to clarify OAs
00012 #include "BonSubMipSolver.hpp"
00013 #include "CbcModel.hpp"
00014 #include "CbcStrategy.hpp"
00015 #include "OsiAuxInfo.hpp"
00016 #include "OsiClpSolverInterface.hpp"
00017 
00018 #include <climits>
00019 #ifdef COIN_HAS_CPX
00020 #include "OsiCpxSolverInterface.hpp"
00021 #include "cplex.h"
00022 void throw_error(const std::string &s, const std::string &f, const std::string &func){
00023 throw CoinError(s,f,func);
00024 }
00025 #define CHECK_CPX_STAT(a,b) if(b) throw_error("Error in CPLEX call",__FILE__,a);
00026 
00027 #endif
00028 
00029 #include "BonRegisteredOptions.hpp"
00030 #include "BonBabSetupBase.hpp"
00031 #include "BonCbcLpStrategy.hpp"
00032 
00033 
00034 namespace Bonmin {
00036   SubMipSolver::SubMipSolver(BabSetupBase &b, const std::string &prefix):
00037       clp_(NULL),
00038       cpx_(NULL),
00039       lowBound_(-DBL_MAX),
00040       optimal_(false),
00041       integerSolution_(NULL),
00042       strategy_(NULL),
00043       ownClp_(false)
00044   {
00045    int ivalue;
00046    b.options()->GetEnumValue("milp_solver",ivalue,prefix);
00047    if (ivalue <= 0) {//uses cbc
00048      strategy_ = new CbcStrategyDefault;
00049      clp_ = new OsiClpSolverInterface;
00050      ownClp_ = true;
00051    }
00052    else if (ivalue == 1) {
00053      CbcStrategyChooseCuts strategy(b, prefix);
00054      strategy_  = new CbcStrategyChooseCuts(b, prefix);
00055      clp_ = new OsiClpSolverInterface;
00056      ownClp_ = true;
00057    }
00058    else if (ivalue == 2) {
00059 #ifdef COIN_HAS_CPX
00060       OsiCpxSolverInterface * cpxSolver = new OsiCpxSolverInterface;
00061 #if 1
00062       b.options()->GetIntegerValue("number_cpx_threads",ivalue,prefix);
00063       CPXsetintparam(cpxSolver->getEnvironmentPtr(), CPX_PARAM_THREADS, ivalue);
00064       b.options()->GetIntegerValue("cpx_parallel_strategy",ivalue,prefix);
00065       CPXsetintparam(cpxSolver->getEnvironmentPtr(), CPX_PARAM_PARALLELMODE, ivalue);
00066 #endif
00067       cpx_ = cpxSolver;
00068 #else
00069       std::cerr << "You have set an option to use CPLEX as the milp\n"
00070       << "subsolver in oa decomposition. However, apparently\n"
00071       << "CPLEX is not configured to be used in bonmin.\n"
00072       << "See the manual for configuring CPLEX\n";
00073       throw -1;
00074 #endif
00075     }
00076 
00077       b.options()->GetEnumValue("milp_strategy",ivalue,prefix);
00078       if(ivalue == 0){
00079         milp_strat_ = FindGoodSolution;
00080       }
00081       else {
00082         milp_strat_ = GetOptimum;
00083       }
00084 
00085       b.options()->GetNumericValue("allowable_fraction_gap", gap_tol_, prefix);
00086 
00087 
00088   }
00089   SubMipSolver::SubMipSolver(const SubMipSolver &copy):
00090       clp_(NULL),
00091       cpx_(NULL),
00092       lowBound_(-DBL_MAX),
00093       optimal_(false),
00094       integerSolution_(NULL),
00095       strategy_(NULL),
00096       milp_strat_(copy.milp_strat_),
00097       gap_tol_(copy.gap_tol_),
00098       ownClp_(copy.ownClp_)
00099   {
00100 #ifdef COIN_HAS_CPX
00101      if(copy.cpx_ != NULL){
00102        cpx_ = new OsiCpxSolverInterface(*copy.cpx_);
00103       int ival;
00104       CPXgetintparam(copy.cpx_->getEnvironmentPtr(), CPX_PARAM_THREADS, &ival);
00105       CPXsetintparam(cpx_->getEnvironmentPtr(), CPX_PARAM_THREADS, ival);
00106       CPXgetintparam(copy.cpx_->getEnvironmentPtr(), CPX_PARAM_PARALLELMODE, &ival);
00107       CPXsetintparam(cpx_->getEnvironmentPtr(), CPX_PARAM_PARALLELMODE, ival);
00108      }
00109 #endif
00110      if(copy.clp_ != NULL){
00111        if(ownClp_) clp_ = new OsiClpSolverInterface(*copy.clp_);
00112        else clp_ = copy.clp_;
00113      }
00114      if(copy.strategy_){
00115         strategy_ = dynamic_cast<CbcStrategyDefault *>(copy.strategy_->clone());
00116         assert(strategy_);
00117      }
00118   }
00119   SubMipSolver::~SubMipSolver()
00120   {
00121     if (strategy_) delete strategy_;
00122     if (integerSolution_) delete [] integerSolution_;
00123     #ifdef COIN_HAS_CPX
00124     if(cpx_) delete cpx_;
00125     #endif
00126     if(ownClp_) delete clp_;
00127   }
00128 
00130   void
00131   SubMipSolver::setLpSolver(OsiSolverInterface * lp)
00132   {
00133 #ifdef COIN_HAS_CPX
00134     if(cpx_){
00135       clp_ = NULL;
00136       cpx_->loadProblem(*lp->getMatrixByCol(), lp->getColLower(), lp->getColUpper(), lp->getObjCoefficients(), lp->getRowLower(), lp->getRowUpper());
00137       int ncols = lp->getNumCols();
00138       for(int i = 0 ; i < ncols ; i++){
00139         if(lp->isInteger(i) || lp->isBinary(i))
00140            cpx_->setInteger(i);
00141         else
00142            cpx_->setContinuous(i);
00143       }
00144     }
00145     else {
00146 #endif
00147       if(ownClp_) delete clp_;
00148       ownClp_ = false;
00149       clp_ = (lp == NULL) ? NULL :
00150               dynamic_cast<OsiClpSolverInterface *>(lp);
00151       assert(clp_);
00152 #ifdef COIN_HAS_CPX
00153     }
00154 #endif
00155     lowBound_ = -COIN_DBL_MAX;
00156     optimal_ = false;
00157     if (integerSolution_) {
00158       delete [] integerSolution_;
00159       integerSolution_ = NULL;
00160     }
00161   }
00162 
00163   OsiSolverInterface * 
00164   SubMipSolver::solver(){
00165          if(clp_ != NULL)
00166            return clp_;
00167          else
00168 #ifdef COIN_HAS_CPX
00169            return cpx_;
00170 #else
00171          return NULL;
00172 #endif
00173       }
00174 
00175  void 
00176  SubMipSolver::find_good_sol(double cutoff, int loglevel, double max_time){
00177 
00178      if(clp_){
00179       CbcStrategyDefault * strat_default = NULL;
00180       if (!strategy_){
00181         strat_default = new CbcStrategyDefault(1,5,5, loglevel);
00182         strat_default->setupPreProcessing();
00183         strategy_ = strat_default;
00184       }
00185       OsiBabSolver empty;
00186       CbcModel cbc(*clp_);
00187       cbc.solver()->setAuxiliaryInfo(&empty);
00188 
00189       //Change Cbc messages prefixes
00190       strcpy(cbc.messagesPointer()->source_,"OCbc");
00191 
00192       cbc.setLogLevel(loglevel);
00193       cbc.solver()->messageHandler()->setLogLevel(0);
00194       clp_->resolve();
00195       cbc.setStrategy(*strategy_);
00196       cbc.setLogLevel(loglevel);
00197       cbc.solver()->messageHandler()->setLogLevel(0);
00198       cbc.setMaximumSeconds(max_time);
00199       cbc.setMaximumSolutions(1);
00200       cbc.setCutoff(cutoff);
00201 
00202       
00203       cbc.branchAndBound();
00204       lowBound_ = cbc.getBestPossibleObjValue();
00205 
00206       if (cbc.isProvenOptimal() || cbc.isProvenInfeasible())
00207         optimal_ = true;
00208       else optimal_ = false;
00209 
00210       if (cbc.getSolutionCount()) {
00211         if (!integerSolution_)
00212           integerSolution_ = new double[clp_->getNumCols()];
00213         CoinCopyN(cbc.bestSolution(), clp_->getNumCols(), integerSolution_);
00214       }
00215       else if (integerSolution_) {
00216         delete [] integerSolution_;
00217         integerSolution_ = NULL;
00218       }
00219       nodeCount_ = cbc.getNodeCount();
00220       iterationCount_ = cbc.getIterationCount();
00221 
00222       if(strat_default != NULL){
00223         delete strat_default;
00224         strategy_ = NULL;
00225       }
00226      }
00227      else if (cpx_){
00228 #ifndef COIN_HAS_CPX
00229         throw CoinError("Unsuported solver, for local searches you should use clp or cplex",
00230             "performLocalSearch",
00231             "OaDecompositionBase::SubMipSolver");
00232 #else
00233         cpx_->messageHandler()->setLogLevel(loglevel);
00234         cpx_->switchToMIP();
00235         CPXENVptr env = cpx_->getEnvironmentPtr();
00236         CPXLPptr cpxlp = cpx_->getLpPtr(OsiCpxSolverInterface::KEEPCACHED_ALL);
00237         CPXsetdblparam(env, CPX_PARAM_TILIM, max_time);
00238         CPXsetdblparam(env, CPX_PARAM_EPINT, 1e-08);
00239         CPXsetdblparam(env, CPX_PARAM_CUTUP, cutoff);
00240       CPXsetdblparam(env, CPX_PARAM_EPGAP, gap_tol_);
00241 
00242 #if 0
00243       CPXsetintparam(env, CPX_PARAM_THREADS, 16);
00244       CPXsetintparam(env, CPX_PARAM_PARALLELMODE, -1);
00245 #endif
00246 
00247         CPXsetintparam(env,CPX_PARAM_INTSOLLIM, 10000);
00248         CPXsetintparam(env,CPX_PARAM_NODELIM, 1000000);
00249 
00250         nodeCount_ = 0;
00251         iterationCount_ = 0;
00252         int status = CPXmipopt(env,cpxlp);
00253         CHECK_CPX_STAT("mipopt",status)
00254        
00255     
00256         status = CPXgetbestobjval(env, cpxlp, &lowBound_);
00257         CHECK_CPX_STAT("bestobjvalue",status)
00258      
00259         int stat = CPXgetstat( env, cpxlp);
00260         optimal_ = (stat == CPXMIP_INFEASIBLE) || (stat == CPXMIP_OPTIMAL) || (stat == CPXMIP_OPTIMAL_TOL); 
00261         nodeCount_ = CPXgetnodecnt(env , cpxlp);
00262         iterationCount_ = CPXgetmipitcnt(env , cpxlp);
00263 
00264         if(stat == CPXMIP_NODE_LIM_INFEAS){
00265           CPXsetintparam(env, CPX_PARAM_INTSOLLIM, 1);
00266           CPXsetintparam(env,CPX_PARAM_NODELIM, 2100000000);
00267            status = CPXmipopt(env,cpxlp);
00268            CHECK_CPX_STAT("mipopt",status)
00269 
00270            status = CPXgetstat( env, cpxlp);
00271            optimal_ = (stat == CPXMIP_INFEASIBLE) || (stat == CPXMIP_OPTIMAL) || (stat == CPXMIP_OPTIMAL_TOL); 
00272            nodeCount_ = CPXgetnodecnt(env , cpxlp);
00273            iterationCount_ = CPXgetmipitcnt(env , cpxlp);
00274        }
00275        bool infeasible = (stat == CPXMIP_INFEASIBLE) || (stat == CPXMIP_ABORT_INFEAS) || (stat == CPXMIP_TIME_LIM_INFEAS) 
00276                           || (stat == CPXMIP_NODE_LIM_INFEAS) || (stat == CPXMIP_FAIL_INFEAS)  
00277                           || (stat == CPXMIP_MEM_LIM_INFEAS) || (stat == CPXMIP_INForUNBD);
00278        
00279 
00280        if(!infeasible){
00281           nodeCount_ += CPXgetnodecnt(env, cpxlp);
00282           //iterationCount_ += CPXgetitcnt(env, cpxlp);
00283           if(!integerSolution_){
00284             integerSolution_ = new double[cpx_->getNumCols()];
00285           }
00286           CPXgetmipx(env, cpxlp, integerSolution_, 0, cpx_->getNumCols() -1);
00287           CHECK_CPX_STAT("getmipx",status)
00288        }
00289        else {
00290          if (integerSolution_) {
00291            delete [] integerSolution_;
00292            integerSolution_ = NULL;
00293          }
00294        }
00295       cpx_->switchToLP();
00296 #endif
00297     } 
00298     else {
00299         throw CoinError("Unsuported solver, for local searches you should use clp or cplex",
00300             "performLocalSearch",
00301             "OaDecompositionBase::SubMipSolver");
00302    }
00303   }
00304 
00305   void
00306   SubMipSolver::optimize(double cutoff, int loglevel, double maxTime)
00307   {
00308     if (clp_) {
00309       assert(strategy_);
00310       CbcStrategyDefault * strat_default = dynamic_cast<CbcStrategyDefault *>(strategy_->clone());
00311       assert(strat_default);
00312       strat_default->setupPreProcessing();
00313 
00314       OsiBabSolver empty;
00315       CbcModel cbc(*clp_);
00316       cbc.solver()->setAuxiliaryInfo(&empty);
00317 
00318       //Change Cbc messages prefixes
00319       strcpy(cbc.messagesPointer()->source_,"OCbc");
00320 
00321       cbc.setLogLevel(loglevel);
00322       cbc.solver()->messageHandler()->setLogLevel(0);
00323       clp_->resolve();
00324       cbc.setStrategy(*strategy_);
00325       cbc.setLogLevel(loglevel);
00326       cbc.solver()->messageHandler()->setLogLevel(0);
00327       cbc.setMaximumSeconds(maxTime);
00328       cbc.setCutoff(cutoff);
00329       cbc.setDblParam( CbcModel::CbcAllowableFractionGap, gap_tol_);
00330 
00331       //cbc.solver()->writeMpsNative("FP.mps", NULL, NULL, 1);
00332       cbc.branchAndBound();
00333       lowBound_ = cbc.getBestPossibleObjValue();
00334 
00335       if (cbc.isProvenOptimal() || cbc.isProvenInfeasible())
00336         optimal_ = true;
00337       else optimal_ = false;
00338 
00339       if (cbc.getSolutionCount()) {
00340         if (!integerSolution_)
00341           integerSolution_ = new double[clp_->getNumCols()];
00342         CoinCopyN(cbc.bestSolution(), clp_->getNumCols(), integerSolution_);
00343       }
00344       else if (integerSolution_) {
00345         delete [] integerSolution_;
00346         integerSolution_ = NULL;
00347       }
00348       nodeCount_ = cbc.getNodeCount();
00349       iterationCount_ = cbc.getIterationCount();
00350       delete strat_default;
00351     }
00352     else 
00353 #ifdef COIN_HAS_CPX
00354     if (cpx_) {
00355       cpx_->switchToMIP();
00356       CPXENVptr env = cpx_->getEnvironmentPtr();
00357       CPXLPptr cpxlp = cpx_->getLpPtr(OsiCpxSolverInterface::KEEPCACHED_ALL);
00358 
00359       CPXsetdblparam(env, CPX_PARAM_TILIM, maxTime);
00360       CPXsetdblparam(env, CPX_PARAM_CUTUP, cutoff);
00361       CPXsetdblparam(env, CPX_PARAM_EPGAP, gap_tol_);
00362       cpx_->messageHandler()->setLogLevel(loglevel);
00363 #if 0
00364       CPXsetintparam(env, CPX_PARAM_THREADS, 16);
00365       CPXsetintparam(env, CPX_PARAM_PARALLELMODE, -1);
00366 #endif
00367       int status = CPXmipopt(env,cpxlp);
00368       CHECK_CPX_STAT("mipopt",status)
00369 
00370       int stat = CPXgetstat( env, cpxlp);
00371       bool infeasible = (stat == CPXMIP_INFEASIBLE) || (stat == CPXMIP_ABORT_INFEAS) || (stat == CPXMIP_TIME_LIM_INFEAS) || (stat == CPXMIP_NODE_LIM_INFEAS) || (stat == CPXMIP_FAIL_INFEAS)
00372                         || (stat == CPXMIP_MEM_LIM_INFEAS) || (stat == CPXMIP_INForUNBD);
00373       optimal_ = (stat == CPXMIP_INFEASIBLE) || (stat == CPXMIP_OPTIMAL) || (stat == CPXMIP_OPTIMAL_TOL); 
00374       nodeCount_ = CPXgetnodecnt(env , cpxlp);
00375       iterationCount_ = CPXgetmipitcnt(env , cpxlp);
00376       status = CPXgetbestobjval(env, cpxlp, &lowBound_);
00377       CHECK_CPX_STAT("getbestobjval",status)
00378        
00379       if(!infeasible){
00380          if(!integerSolution_){
00381            integerSolution_ = new double[cpx_->getNumCols()];
00382          }
00383          CPXgetmipx(env, cpxlp, integerSolution_, 0, cpx_->getNumCols() -1);
00384          CHECK_CPX_STAT("getmipx",status)
00385       }
00386       else {
00387         if (integerSolution_) {
00388           delete [] integerSolution_;
00389           integerSolution_ = NULL;
00390         }
00391       }
00392       cpx_->switchToLP();
00393     }
00394     else {
00395 #else
00396      {
00397 #endif
00398         throw CoinError("Unsuported solver, for local searches you should use clp or cplex",
00399             "performLocalSearch",
00400             "OaDecompositionBase::SubMipSolver");
00401       }
00402 }
00403 
00405    void 
00406    SubMipSolver::setStrategy(CbcStrategyDefault * strategy)
00407    {
00408      if (strategy_) delete strategy_;
00409      strategy_ = dynamic_cast<CbcStrategyDefault *>(strategy->clone());
00410      assert(strategy_);
00411    }
00412 
00414   void
00415   SubMipSolver::registerOptions(Ipopt::SmartPtr<Bonmin::RegisteredOptions> roptions)
00416   {
00417     roptions->SetRegisteringCategory("Options for MILP solver", RegisteredOptions::BonminCategory);
00418     roptions->AddStringOption3("milp_solver",
00419         "Choose the subsolver to solve MILP sub-problems in OA decompositions.",
00420         "Cbc_D",
00421         "Cbc_D","Coin Branch and Cut with its default",
00422         "Cbc_Par", "Coin Branch and Cut with passed parameters",
00423         "Cplex","Ilog Cplex",
00424         " To use Cplex, a valid license is required and you should have compiled OsiCpx in COIN-OR  (see Osi documentation).");
00425     roptions->setOptionExtraInfo("milp_solver",64);
00426 
00427     roptions->AddBoundedIntegerOption("milp_log_level",
00428         "specify MILP solver log level.",
00429         0,4,0,
00430         "Set the level of output of the MILP subsolver in OA : "
00431         "0 - none, 1 - minimal, 2 - normal low, 3 - normal high"
00432                                      );
00433     roptions->setOptionExtraInfo("milp_log_level",64);
00434 
00435     roptions->AddBoundedIntegerOption("cpx_parallel_strategy",
00436                            "Strategy of parallel search mode in CPLEX.",
00437                            -1, 1, 0,
00438                            "-1 = opportunistic, 0 = automatic, 1 = deterministic (refer to CPLEX documentation)"
00439                            );
00440     roptions->setOptionExtraInfo("cpx_parallel_strategy",64);
00441 
00442     roptions->AddLowerBoundedIntegerOption("number_cpx_threads",
00443                            "Set number of threads to use with cplex.",
00444                            0, 0,
00445                            "(refer to CPLEX documentation)"
00446                            );
00447     roptions->setOptionExtraInfo("number_cpx_threads",64);
00448 
00449     
00450     roptions->AddStringOption2("milp_strategy",
00451         "Choose a strategy for MILPs.",
00452         "find_good_sol",
00453         "find_good_sol","Stop sub milps when a solution improving the incumbent is found",
00454         "solve_to_optimality", "Solve MILPs to optimality",
00455         "");
00456     roptions->setOptionExtraInfo("milp_strategy",64);
00457 
00458   }
00459 }/* Ends Bonmin namespace.*/

Generated on Thu Sep 22 03:05:53 2011 by  doxygen 1.4.7