/home/coin/SVN-release/OS-2.1.1/Couenne/src/branch/CouenneChooseStrong.cpp

Go to the documentation of this file.
00001 /* $Id: CouenneChooseStrong.cpp 272 2009-11-20 22:21:50Z pbelotti $
00002  *
00003  * Name:    CouenneChooseStrong.cpp
00004  * Authors: Andreas Waechter, IBM Corp.
00005  *          Pietro Belotti, Lehigh University
00006  * Purpose: Strong branching objects for Couenne
00007  *
00008  * (C) Carnegie-Mellon University, 2006-09.
00009  * This file is licensed under the Common Public License (CPL)
00010  */
00011 
00012 #include "BonChooseVariable.hpp"
00013 #include "CouenneChooseStrong.hpp"
00014 #include "CouenneProblem.hpp"
00015 #include "CouenneBranchingObject.hpp"
00016 
00017 const CouNumber estProdEps = 1e-6;
00018 
00019 
00021   CouenneChooseStrong::CouenneChooseStrong (Bonmin::BabSetupBase &b, CouenneProblem* p, JnlstPtr jnlst) :
00022 
00023   Bonmin::BonChooseVariable (b, b.continuousSolver()),
00024   problem_          (p),
00025   jnlst_            (jnlst),
00026   branchtime_       (0.) {
00027 
00028     std::string s;
00029 
00030     b.options () -> GetStringValue ("pseudocost_mult_lp", s, "couenne.");
00031     pseudoUpdateLP_ = (s == "yes");      
00032 
00033     b.options () -> GetStringValue ("estimate_select", s, "couenne.");
00034     estimateProduct_ = (s == "product");
00035   }
00036 
00038   CouenneChooseStrong::CouenneChooseStrong (const CouenneChooseStrong& rhs) :
00039     Bonmin::BonChooseVariable (rhs),
00040     problem_          (rhs.problem_),
00041     pseudoUpdateLP_   (rhs.pseudoUpdateLP_),
00042     estimateProduct_  (rhs.estimateProduct_),
00043     jnlst_            (rhs.jnlst_),
00044     branchtime_       (rhs.branchtime_)
00045   {}
00046 
00048   CouenneChooseStrong::~CouenneChooseStrong()
00049   {if (branchtime_ > 1e-9) jnlst_ -> Printf (J_ERROR, J_BRANCHING, "Strong branching: total time %g\n", branchtime_);}
00050 
00052   OsiChooseVariable *
00053   CouenneChooseStrong::clone() const
00054   {return new CouenneChooseStrong(*this);}
00055 
00057   CouenneChooseStrong&
00058   CouenneChooseStrong::operator=(const CouenneChooseStrong & rhs)
00059   {
00060     if (this != &rhs) {
00061       Bonmin::BonChooseVariable::operator=(rhs);
00062       problem_         = rhs.problem_;
00063       pseudoUpdateLP_  = rhs.pseudoUpdateLP_;
00064       estimateProduct_ = rhs.estimateProduct_;
00065       jnlst_           = rhs.jnlst_;
00066       branchtime_      = rhs.branchtime_;
00067     }
00068     return *this;
00069   }
00070 
00071 
00072   /* Choose a variable. Returns:
00073 
00074     -1  Node is infeasible
00075      0  Normal termination - we have a candidate
00076      1  All look satisfied - no candidate
00077      2  We can change the bound on a variable - but we also have a strong branching candidate
00078      3  We can change the bound on a variable - but we have a non-strong branching candidate
00079      4  We can change the bound on a variable - no other candidates
00080 
00081      We can pick up branch from whichObject() and whichWay()
00082      We can pick up a forced branch (can change bound) from whichForcedObject() and whichForcedWay()
00083      If we have a solution then we can pick up from goodObjectiveValue() and goodSolution()
00084   */
00085   int CouenneChooseStrong::chooseVariable (OsiSolverInterface * solver,
00086                                            OsiBranchingInformation *info,
00087                                            bool fixVariables) {
00088 
00092 
00093     problem_ -> domain () -> push
00094       (problem_ -> nVars (),
00095        info -> solution_,
00096        info -> lower_,
00097        info -> upper_);
00098 
00099     int retval;
00100 
00101     //int retval = BonChooseVariable::chooseVariable (solver, info, fixVariables);
00102 
00103     // COPY of Bonmin starts here ////////////////////////////////////////////
00104 
00105     // We assume here that chooseVariable is called once at the very
00106     // beginning with fixVariables set to true.  This is then the root
00107     // node.
00108 
00109     bool isRoot = isRootNode(info);
00110     int numberStrong = numberStrong_;
00111     if (isRoot) {
00112       numberStrong = CoinMax(numberStrong_, numberStrongRoot_);
00113     }
00114     if (numberUnsatisfied_) {
00115       const double* upTotalChange = pseudoCosts_.upTotalChange();
00116       const double* downTotalChange = pseudoCosts_.downTotalChange();
00117       const int* upNumber = pseudoCosts_.upNumber();
00118       const int* downNumber = pseudoCosts_.downNumber();
00119       int numberBeforeTrusted = pseudoCosts_.numberBeforeTrusted();
00120       int numberLeft = CoinMin(numberStrong -numberStrongDone_,numberUnsatisfied_);
00121       results_.clear();
00122       int returnCode=0;
00123       bestObjectIndex_ = -1;
00124       bestWhichWay_ = -1;
00125       firstForcedObjectIndex_ = -1;
00126       firstForcedWhichWay_ =-1;
00127       double bestTrusted=-COIN_DBL_MAX;
00128       for (int i=0;i<numberLeft;i++) {
00129         int iObject = list_[i];
00130         if (numberBeforeTrusted == 0||
00131             i < minNumberStrongBranch_ ||
00132             (
00133               only_pseudo_when_trusted_ && number_not_trusted_>0 ) ||
00134               ( !isRoot && (upNumber[iObject]<numberBeforeTrusted ||
00135                           downNumber[iObject]<numberBeforeTrusted ))||
00136               ( isRoot && (!upNumber[iObject] && !downNumber[iObject])) ) {
00137          
00138           results_.push_back(Bonmin::HotInfo(solver, info,
00139                                 solver->objects(), iObject));
00140         }
00141         else {
00142           const OsiObject * obj = solver->object(iObject);
00143           double
00144             upEstimate       = (upTotalChange[iObject]*obj->upEstimate())/upNumber[iObject],
00145             downEstimate     = (downTotalChange[iObject]*obj->downEstimate())/downNumber[iObject],
00146             MAXMIN_CRITERION = maxminCrit (info),
00147             minVal, maxVal, value;
00148 
00149           if (downEstimate < upEstimate) {
00150             minVal = downEstimate;
00151             maxVal = upEstimate;
00152           } else {
00153             minVal = upEstimate;
00154             maxVal = downEstimate;
00155           }
00156 
00157           value = 
00158             estimateProduct_ ? 
00159             ((estProdEps + minVal) * maxVal) :
00160             (       MAXMIN_CRITERION  * minVal + 
00161              (1.0 - MAXMIN_CRITERION) * maxVal);
00162 
00163           if (value > bestTrusted) {
00164             bestObjectIndex_=iObject;
00165             bestWhichWay_ = upEstimate>downEstimate ? 0 : 1;
00166             bestTrusted = value;
00167           }
00168         }
00169       }
00170       int numberFixed=0;
00171       if (results_.size() > 0) {
00172 
00173         // do strong branching //////////////////////////////////////////////////
00174         returnCode = doStrongBranching (solver, info, results_.size(), 1);
00175 
00176         if (bb_log_level_>=3) {
00177           const char* stat_msg[] = {"NOTDON", "FEAS", "INFEAS", "NOFINI"};
00178           message(SB_HEADER)<<CoinMessageEol;
00179           for (unsigned int i = 0; i< results_.size(); i++) {
00180             double up_change = results_[i].upChange();
00181             double down_change = results_[i].downChange();
00182             int up_status = results_[i].upStatus();
00183             int down_status = results_[i].downStatus();
00184             message(SB_RES)<<(int) i<<stat_msg[down_status+1]<<down_change
00185             <<stat_msg[up_status+1]<< up_change<< CoinMessageEol;
00186           }
00187         }
00188         if (returnCode >= 0 && 
00189             returnCode <= 2) {
00190           if (returnCode)
00191             returnCode = (bestObjectIndex_>=0) ? 3 : 4;
00192           for (unsigned int i=0;i < results_.size();i++) {
00193             int iObject = results_[i].whichObject();
00194             double upEstimate;
00195             if (results_[i].upStatus()!=1) {
00196               assert (results_[i].upStatus()>=0);
00197               upEstimate = results_[i].upChange();
00198             }
00199             else {
00200               // infeasible - just say expensive
00201               if (info->cutoff_<1.0e50)
00202                 upEstimate = 2.0*(info->cutoff_-info->objectiveValue_);
00203               else
00204                 upEstimate = 2.0*fabs(info->objectiveValue_);
00205               if (firstForcedObjectIndex_ <0) {
00206                 // first fixed variable
00207                 firstForcedObjectIndex_ = iObject;
00208                 firstForcedWhichWay_ =0;
00209               }
00210               numberFixed++;
00211               if (fixVariables) {
00212                 const OsiObject * obj = solver->object(iObject);
00213                 OsiBranchingObject * branch = obj->createBranch(solver,info,0);
00214                 branch->branch(solver);
00215                 delete branch;
00216               }
00217             }
00218             double downEstimate;
00219             if (results_[i].downStatus()!=1) {
00220               assert (results_[i].downStatus()>=0);
00221               downEstimate = results_[i].downChange();
00222             }
00223             else {
00224               // infeasible - just say expensive
00225               if (info->cutoff_<1.0e50)
00226                 downEstimate = 2.0*(info->cutoff_-info->objectiveValue_);
00227               else
00228                 downEstimate = 2.0*fabs(info->objectiveValue_);
00229               if (firstForcedObjectIndex_ <0) {
00230                 firstForcedObjectIndex_ = iObject;
00231                 firstForcedWhichWay_ =1;
00232               }
00233               numberFixed++;
00234               if (fixVariables) {
00235                 const OsiObject * obj = solver->object(iObject);
00236                 OsiBranchingObject * branch = obj->createBranch(solver,info,1);
00237                 branch->branch(solver);
00238                 delete branch;
00239               }
00240             }
00241 
00242             double
00243               MAXMIN_CRITERION = maxminCrit(info),            
00244               minVal, maxVal, value;
00245 
00246             if (downEstimate < upEstimate) {
00247               minVal = downEstimate;
00248               maxVal = upEstimate;
00249             } else {
00250               minVal = upEstimate;
00251               maxVal = downEstimate;
00252             }
00253 
00254             value = 
00255               estimateProduct_ ? 
00256               ((estProdEps + minVal) * maxVal) :
00257               (       MAXMIN_CRITERION  * minVal + 
00258                (1.0 - MAXMIN_CRITERION) * maxVal);
00259 
00260             if (value>bestTrusted) {
00261               bestTrusted = value;
00262               bestObjectIndex_ = iObject;
00263               bestWhichWay_ = upEstimate>downEstimate ? 0 : 1;
00264               // but override if there is a preferred way
00265               const OsiObject * obj = solver->object(iObject);
00266               if (obj->preferredWay()>=0&&obj->infeasibility())
00267                 bestWhichWay_ = obj->preferredWay();
00268               if (returnCode)
00269                 returnCode = 2;
00270             }
00271           }
00272         }
00273         else if (returnCode==3) {
00274           // max time - just choose one
00275           bestObjectIndex_ = list_[0];
00276           bestWhichWay_ = 0;
00277           returnCode=0;
00278         }
00279       }
00280       else {
00281         bestObjectIndex_=list_[0];
00282       }
00283 
00284       if ( bestObjectIndex_ >=0 ) {
00285         OsiObject * obj = solver->objects()[bestObjectIndex_];
00286         obj->setWhichWay(bestWhichWay_);
00287         message(BRANCH_VAR)<<obj->columnNumber()<< bestWhichWay_ <<CoinMessageEol;
00288       }
00289 
00290       message(CHOSEN_VAR)<<bestObjectIndex_<<CoinMessageEol;
00291 
00292       if (numberFixed==numberUnsatisfied_&&numberFixed)
00293         returnCode=4;
00294       retval = returnCode;
00295     }
00296     else {
00297       retval = 1;
00298     }
00299 
00300     // COPY of Bonmin ends here //////////////////////////////////////////////
00301 
00302     problem_ -> domain () -> pop ();
00303 
00304     return retval;
00305   }
00306 
00307 
00308   // Sets up strong list and clears all if initialize is true.
00309   // Returns number of infeasibilities.
00310   int CouenneChooseStrong::setupList (OsiBranchingInformation *info, bool initialize) {
00311 
00312     initialize = true; // to avoid failed assert in BonChooseVariable::setupList()
00313 
00314     problem_ -> domain () -> push 
00315       (problem_ -> nVars (),
00316        info -> solution_, 
00317        info -> lower_, 
00318        info -> upper_); // have to alloc+copy
00319 
00320     jnlst_ -> Printf (J_ITERSUMMARY, J_BRANCHING, 
00321                       "----------------- (strong) setup list\n");
00322 
00323     if (jnlst_ -> ProduceOutput (J_DETAILED, J_BRANCHING)) {
00324       for (int i=0; i<problem_ -> domain () -> current () -> Dimension (); i++)
00325         printf ("%4d %20.4g [%20.4g %20.4g]\n", i,
00326                 info -> solution_ [i], info -> lower_ [i], info -> upper_ [i]);
00327     }
00328 
00329     // call Bonmin's setuplist
00330     int retval = Bonmin::BonChooseVariable::setupList (info, initialize);
00331 
00332     jnlst_ -> Printf (J_ITERSUMMARY, J_BRANCHING, 
00333                       "----------------- (strong) setup list done - %d infeasibilities\n", retval);
00334 
00335     problem_ -> domain () -> pop ();
00336     return retval;
00337   }
00338 
00339 
00341   void CouenneChooseStrong::registerOptions (Ipopt::SmartPtr <Bonmin::RegisteredOptions> roptions) {
00342 
00343     roptions -> AddStringOption6
00344       ("pseudocost_mult",
00345        "Multipliers of pseudocosts for estimating and update estimation of bound",
00346        "interval_br_rev",
00347 
00348        "infeasibility", "infeasibility returned by object",
00349 
00350        "projectDist",   "distance between current LP point and resulting branches' LP points",
00351 
00352        "interval_lp",   "width of the interval between bound and current lp point",
00353        "interval_lp_rev",   "similar to interval_lp, reversed",
00354 
00355        "interval_br",   "width of the interval between bound and branching point",
00356        "interval_br_rev",   "similar to interval_br, reversed");
00357 
00358     roptions -> AddStringOption2
00359       ("pseudocost_mult_lp",
00360        "Use distance between LP points to update multipliers of pseudocosts "  
00361        "after simulating branching",
00362        "no",
00363        "yes", "",
00364        "no",  "");
00365 
00366     roptions -> AddStringOption2
00367       ("estimate_select",
00368        "How the min/max estimates of the subproblems' bounds are used in strong branching",
00369        "normal",
00370        "normal",   "as usual in literature",
00371        "product",  "use their product");
00372   }
00373 
00374 
00375   // Returns true if solution looks feasible against given objects
00376   bool CouenneChooseStrong::feasibleSolution (const OsiBranchingInformation * info,
00377                                               const double * solution,
00378                                               int numberObjects,
00379                                               const OsiObject ** objects) {
00380 
00381     double obj = solution [problem_ -> Obj (0) -> Body () -> Index ()];
00382     return problem_ -> checkNLP (solution, obj);
00383   }
00384 //}

Generated on Mon May 3 03:05:18 2010 by  doxygen 1.4.7