00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013 #include "CouenneObject.hpp"
00014 #include "BonChooseVariable.hpp"
00015 #include "CouenneChooseStrong.hpp"
00016 #include "CouenneProblem.hpp"
00017 #include "CouenneProblemElem.hpp"
00018 #include "CouenneBranchingObject.hpp"
00019
00020 #include "CouenneRecordBestSol.hpp"
00021
00022
00023
00024
00025
00026
00027
00028
00029
00030
00031 using namespace Ipopt;
00032 using namespace Couenne;
00033
00034 const CouNumber estProdEps = 1e-6;
00035
00036
00038 CouenneChooseStrong::CouenneChooseStrong (Bonmin::BabSetupBase &b, CouenneProblem* p, JnlstPtr jnlst) :
00039
00040 Bonmin::BonChooseVariable (b, b.continuousSolver()),
00041 problem_ (p),
00042 jnlst_ (jnlst),
00043 branchtime_ (0.) {
00044
00045 std::string s;
00046
00047 b.options () -> GetStringValue ("pseudocost_mult_lp", s, "couenne.");
00048 pseudoUpdateLP_ = (s == "yes");
00049
00050 b.options () -> GetStringValue ("estimate_select", s, "couenne.");
00051 estimateProduct_ = (s == "product");
00052
00053 b.options () -> GetStringValue ("trust_strong", s, "couenne.");
00054
00055
00056
00057 setTrustStrongForSolution (s == "yes");
00058 setTrustStrongForBound (s == "yes");
00059 }
00060
00062 CouenneChooseStrong::CouenneChooseStrong (const CouenneChooseStrong& rhs) :
00063 Bonmin::BonChooseVariable (rhs),
00064 problem_ (rhs.problem_),
00065 pseudoUpdateLP_ (rhs.pseudoUpdateLP_),
00066 estimateProduct_ (rhs.estimateProduct_),
00067 jnlst_ (rhs.jnlst_),
00068 branchtime_ (rhs.branchtime_)
00069 {}
00070
00072 CouenneChooseStrong::~CouenneChooseStrong()
00073 {if (branchtime_ > 1e-9) jnlst_ -> Printf (J_ERROR, J_BRANCHING, "Strong branching: total time %g\n", branchtime_);}
00074
00076 OsiChooseVariable *
00077 CouenneChooseStrong::clone() const
00078 {return new CouenneChooseStrong(*this);}
00079
00081 CouenneChooseStrong&
00082 CouenneChooseStrong::operator=(const CouenneChooseStrong & rhs)
00083 {
00084 if (this != &rhs) {
00085 Bonmin::BonChooseVariable::operator=(rhs);
00086 problem_ = rhs.problem_;
00087 pseudoUpdateLP_ = rhs.pseudoUpdateLP_;
00088 estimateProduct_ = rhs.estimateProduct_;
00089 jnlst_ = rhs.jnlst_;
00090 branchtime_ = rhs.branchtime_;
00091 }
00092 return *this;
00093 }
00094
00095
00096 int CouenneChooseStrong::goodCandidate(OsiSolverInterface *solver,
00097 OsiBranchingInformation *info,
00098 OsiObject **object, const int iObject,
00099 const double prec) {
00100 CouenneObject *co = dynamic_cast <CouenneObject *>(object[iObject]);
00101 OsiSimpleInteger *simpl = dynamic_cast <OsiSimpleInteger *>(object[iObject]);
00102
00103 int vInd = -1;
00104 bool varIsInt = false;
00105 if(co) {
00106 vInd = co->Reference()->Index();
00107 if(vInd >= 0) {
00108 varIsInt = co->Reference()->isInteger();
00109 }
00110 }
00111 else {
00112 if(simpl) {
00113 vInd = object[iObject]->columnNumber();
00114 varIsInt = true;
00115 }
00116 else {
00117
00118
00119
00120
00121
00122
00123
00124 return 3;
00125 }
00126 }
00127
00128 int goodCand = 3;
00129
00130
00131
00132
00133 if((vInd >= 0) && varIsInt) {
00134 double vUp = solver->getColUpper()[vInd];
00135 double vLow = solver->getColLower()[vInd];
00136 double infoVal = info->solution_[vInd];
00137 double distToInt = fabs(infoVal - floor(infoVal + 0.5));
00138 if(distToInt > 0.5) {
00139 distToInt = 1 - distToInt;
00140 }
00141 if(simpl) {
00142 goodCand = 0;
00143 if((distToInt > info->integerTolerance_) &&
00144 (vUp > vLow + prec)) {
00145 goodCand = 3;
00146 if((vUp + prec < infoVal) || (infoVal < vLow - prec)) {
00147 goodCand = 1;
00148 }
00149 }
00150 }
00151 if(co) {
00152 goodCand = 2;
00153 if(vUp > vLow + prec) {
00154 goodCand = 3;
00155 }
00156 }
00157 }
00158 return(goodCand);
00159 }
00160
00161
00162 bool CouenneChooseStrong::saveBestCand(OsiObject **object, const int iObject,
00163 const double value,
00164 const double upEstimate,
00165 const double downEstimate,
00166 double &bestVal1,
00167 double &bestVal2, int &bestIndex,
00168 int &bestWay) {
00169 bool retval = false;
00170 if(value > bestVal1) {
00171 retval = true;
00172 bestVal1 = value;
00173 bestIndex = iObject;
00174 bestWay = upEstimate > downEstimate ? 0 : 1;
00175
00176 const OsiObject * obj = object[iObject];
00177 if (obj->preferredWay() >= 0 && obj->infeasibility()) {
00178 bestWay = obj->preferredWay();
00179 }
00180 }
00181 return(retval);
00182 }
00183
00184
00185
00186
00187
00188
00189
00190
00191
00192
00193
00194
00195
00196
00197
00198 int CouenneChooseStrong::chooseVariable (OsiSolverInterface * solver,
00199 OsiBranchingInformation *info,
00200 bool fixVariables) {
00201
00205
00206 problem_ -> domain () -> push
00207 (problem_ -> nVars (),
00208 info -> solution_,
00209 info -> lower_,
00210 info -> upper_);
00211
00212 #ifdef TRACE_STRONG2
00213 int pv = -1;
00214 if(problem_->doPrint_) {
00215 if(pv > -1) {
00216 printf("CCS: beg: x[%d]: %10.4f lb: %10.4f ub: %10.4f\n",
00217 pv, solver->getColSolution()[pv], solver->getColLower()[pv], solver->getColUpper()[pv]);
00218 printf("CCS: info: x[%d]: %10.4f lb: %10.4f ub: %10.4f\n",
00219 pv, info->solution_[pv], info->lower_[pv], info->upper_[pv]);
00220 printf("CCS: problem: lb: %10.4f ub: %10.4f\n",
00221 problem_->Lb(pv), problem_->Ub(pv));
00222 }
00223 }
00224 #endif
00225
00226 int retval;
00227 const double prec = problem_->getFeasTol();
00228
00229
00230
00231
00232
00233
00234
00235
00236
00237 bool isRoot = isRootNode(info);
00238 int numberStrong = numberStrong_;
00239 if (isRoot) {
00240 numberStrong = CoinMax(numberStrong_, numberStrongRoot_);
00241 }
00242 if (numberUnsatisfied_) {
00243 int cardIndForPseudo = 0, *indForPseudo = new int[numberUnsatisfied_];
00244 OsiObject ** object = solver->objects();
00245 const double* upTotalChange = pseudoCosts_.upTotalChange();
00246 const double* downTotalChange = pseudoCosts_.downTotalChange();
00247 const int* upNumber = pseudoCosts_.upNumber();
00248 const int* downNumber = pseudoCosts_.downNumber();
00249 int numberBeforeTrusted = pseudoCosts_.numberBeforeTrusted();
00250 int numberLeft = CoinMin(numberStrong -numberStrongDone_,numberUnsatisfied_);
00251 results_.clear();
00252 int returnCode = 0;
00253 int returnCodeSB = 0;
00254 bestObjectIndex_ = -1;
00255 bestWhichWay_ = -1;
00256 firstForcedObjectIndex_ = -1;
00257 firstForcedWhichWay_ =-1;
00258 double bestTrustedVal1 = -COIN_DBL_MAX;
00259 double bestTrustedVal2 = -COIN_DBL_MAX;
00260
00261 bool smallGap = false;
00262 bool sbObjPosImp = false;
00263
00264
00265
00266
00267 #ifdef USE_SMALL_GAP
00268 int objInd = problem_ -> Obj (0) -> Body () -> Index ();
00269 double lbGap = info -> lower_ [objInd];
00270 double ubGap = problem_ -> getCutOff ();
00271 double currentGap =
00272 (ubGap > COUENNE_INFINITY / 10 ||
00273 lbGap < -Couenne_large_bound / 10) ? 1e3 :
00274 fabs (ubGap - lbGap) / (1.e-3 + CoinMin (fabs (ubGap), fabs (lbGap)));
00275
00276 if(currentGap < 1e-3) {
00277 smallGap = true;
00278 }
00279 #endif
00280
00281 #ifdef TRACE_STRONG
00282 if((problem_->doPrint_) && (number_not_trusted_ > 0)) {
00283 printf("number_not_trusted: %d\n", number_not_trusted_);
00284 }
00285 #endif
00286
00287 for (int i=0;i<numberLeft;i++) {
00288 int iObject = list_[i];
00289 if (numberBeforeTrusted == 0||
00290 i < minNumberStrongBranch_ ||
00291 (
00292 only_pseudo_when_trusted_ && number_not_trusted_>0 ) ||
00293 ( !isRoot && (upNumber[iObject]<numberBeforeTrusted ||
00294 downNumber[iObject]<numberBeforeTrusted ))||
00295 ( isRoot && (!upNumber[iObject] && !downNumber[iObject])) ) {
00296
00297 #ifdef TRACE_STRONG
00298 if(problem_->doPrint_) {
00299 printf("Push object %d for strong branch\n", iObject);
00300 }
00301 #endif
00302 results_.push_back(Bonmin::HotInfo(solver, info,
00303 object, iObject));
00304 }
00305 else {
00306
00307 #ifdef TRACE_STRONG
00308 if(problem_->doPrint_) {
00309 printf("Use pseudo cost for object %d\n", iObject);
00310 }
00311 #endif
00312 indForPseudo[cardIndForPseudo] = iObject;
00313 cardIndForPseudo++;
00314 }
00315 }
00316
00317 int numberFixed=0;
00318
00319 if (results_.size() > 0) {
00320
00321
00322 returnCodeSB = doStrongBranching (solver, info, results_.size(), 1);
00323
00324 if (bb_log_level_>=3) {
00325 const char* stat_msg[] = {"NOTDON", "FEAS", "INFEAS", "NOFINI"};
00326 message(SB_HEADER)<<CoinMessageEol;
00327 for (unsigned int i = 0; i< results_.size(); i++) {
00328 double up_change = results_[i].upChange();
00329 double down_change = results_[i].downChange();
00330 int up_status = results_[i].upStatus();
00331 int down_status = results_[i].downStatus();
00332 message(SB_RES)<<(int) i<<stat_msg[down_status+1]<<down_change
00333 <<stat_msg[up_status+1]<< up_change<< CoinMessageEol;
00334 }
00335 }
00336
00337 if (returnCodeSB >= 0 && returnCodeSB <= 2) {
00338 if(returnCodeSB > 0) {
00339 returnCode = 4;
00340 }
00341 else {
00342 returnCode = 0;
00343 }
00344 for (unsigned int i=0;i < results_.size (); i++) {
00345
00346 if((results_[i].upStatus() < 0) || (results_[i].downStatus() < 0))
00347 continue;
00348
00349
00350 if((results_[i].upStatus() < 0) || (results_[i].downStatus() < 0)) {
00351 continue;
00352 }
00353
00354 int iObject = results_[i].whichObject();
00355 const OsiObject * obj = object[iObject];
00356 int needBranch = goodCandidate(solver, info, object, iObject, prec);
00357
00359
00360 double upEstimate;
00361
00362 if (results_[i].upStatus()!=1) {
00363 assert (results_[i].upStatus()>=0);
00364 upEstimate = results_[i].upChange();
00365 }
00366 else {
00367
00368 if (info->cutoff_<1.0e50)
00369 upEstimate = 2.0*(info->cutoff_-info->objectiveValue_);
00370 else
00371 upEstimate = 2.0*fabs(info->objectiveValue_);
00372 if (firstForcedObjectIndex_ <0) {
00373
00374 firstForcedObjectIndex_ = iObject;
00375 firstForcedWhichWay_ =0;
00376 }
00377
00378 numberFixed++;
00379 if (fixVariables) {
00380 if(needBranch >= 2) {
00381
00382
00383 OsiBranchingObject * branch =
00384 obj->createBranch(solver, info, 0);
00385 branch -> branch (solver);
00386 delete branch;
00387 }
00388 }
00389 }
00390
00392
00393 double downEstimate;
00394
00395 if (results_[i].downStatus()!=1) {
00396 assert (results_[i].downStatus()>=0);
00397 downEstimate = results_[i].downChange();
00398 }
00399 else {
00400
00401 if (info->cutoff_<1.0e50)
00402 downEstimate = 2.0*(info->cutoff_-info->objectiveValue_);
00403 else
00404 downEstimate = 2.0*fabs(info->objectiveValue_);
00405 if (firstForcedObjectIndex_ <0) {
00406 firstForcedObjectIndex_ = iObject;
00407 firstForcedWhichWay_ =1;
00408 }
00409 numberFixed++;
00410 if (fixVariables) {
00411 if(needBranch >= 2) {
00412
00413
00414 OsiBranchingObject * branch =
00415 obj->createBranch(solver, info, 1);
00416 branch -> branch (solver);
00417 delete branch;
00418 }
00419 }
00420 }
00421
00422 double
00423 MAXMIN_CRITERION = maxminCrit (info),
00424 minVal, maxVal, value;
00425
00426 if (downEstimate < upEstimate) {
00427 minVal = downEstimate;
00428 maxVal = upEstimate;
00429 } else {
00430 minVal = upEstimate;
00431 maxVal = downEstimate;
00432 }
00433
00434 if(smallGap) {
00435 value = minVal;
00436 }
00437 else {
00438 value =
00439 estimateProduct_ ?
00440 ((estProdEps + minVal) * maxVal) :
00441 ( MAXMIN_CRITERION * minVal +
00442 (1.0 - MAXMIN_CRITERION) * maxVal);
00443 }
00444
00445 if((needBranch == 3) &&
00446 saveBestCand(object, iObject, value, upEstimate, downEstimate,
00447 bestTrustedVal1,
00448 bestTrustedVal2, bestObjectIndex_, bestWhichWay_)) {
00449 if(returnCodeSB) {
00450 returnCode = 2;
00451 }
00452
00453 #ifdef USE_SMALL_GAP
00454 if(smallGap && (minVal > 1e-3)) {
00455 sbObjPosImp = true;
00456 }
00457 #endif
00458
00459 }
00460 }
00461 }
00462 else {
00463 if (returnCodeSB == 3) {
00464 if(bestObjectIndex_ < 0) {
00465
00466
00467 bestObjectIndex_ = list_[0];
00468 bestWhichWay_ = 0;
00469 returnCode = 0;
00470 }
00471 }
00472 else {
00473 returnCode = -1;
00474 }
00475 }
00476 #ifdef OLD_STYLE
00477 if(bestObjectIndex_ < 0) {
00478 bestObjectIndex_ = list_[0];
00479 }
00480 cardIndForPseudo = 0;
00481 #endif
00482 }
00483
00484 if((returnCodeSB != -1) &&
00485 ((returnCode != 0) || (!sbObjPosImp))) {
00486
00487
00488
00489
00490
00491
00492
00493 int bestObjectIndex2 = -1;
00494 int bestWhichWay2 = 0;
00495 double bestTrusted2Val1 = -COIN_DBL_MAX;
00496 double bestTrusted2Val2 = -COIN_DBL_MAX;
00497
00498 for(int ips=0; ips<cardIndForPseudo; ips++) {
00499 int iObject = indForPseudo[ips];
00500 const OsiObject * obj = object[iObject];
00501 int needBranch = goodCandidate(solver, info, object, iObject, prec);
00502
00503 double
00504 upEstimate = (upTotalChange[iObject]*obj->upEstimate())/upNumber[iObject],
00505 downEstimate = (downTotalChange[iObject]*obj->downEstimate())/downNumber[iObject],
00506 MAXMIN_CRITERION = maxminCrit (info),
00507 minVal, maxVal, value;
00508
00509 if (downEstimate < upEstimate) {
00510 minVal = downEstimate;
00511 maxVal = upEstimate;
00512 } else {
00513 minVal = upEstimate;
00514 maxVal = downEstimate;
00515 }
00516
00517 value =
00518 estimateProduct_ ?
00519 ((estProdEps + minVal) * maxVal) :
00520 ( MAXMIN_CRITERION * minVal +
00521 (1.0 - MAXMIN_CRITERION) * maxVal);
00522
00523
00524
00525 if(needBranch != 3) {
00526 if(saveBestCand(object, iObject, value,
00527 upEstimate, downEstimate,
00528 bestTrusted2Val1,
00529 bestTrusted2Val2, bestObjectIndex2,
00530 bestWhichWay2)) {
00531
00532 }
00533
00534 #ifdef TRACE_STRONG
00535 if(problem_->doPrint_) {
00536 printf("Object %d skip pseudocost\n", iObject);
00537 }
00538 #endif
00539 }
00540 else {
00541 if(saveBestCand(object, iObject, value,
00542 upEstimate, downEstimate,
00543 bestTrustedVal1,
00544 bestTrustedVal2, bestObjectIndex_, bestWhichWay_)) {
00545 returnCode = (returnCode ? 3 : 0);
00546
00547 }
00548
00549 #ifdef TRACE_STRONG
00550 if(problem_->doPrint_) {
00551 printf("Object %d use pseudocost\n", iObject);
00552 }
00553 #endif
00554 }
00555 }
00556
00557 delete[] indForPseudo;
00558
00559 if((bestObjectIndex_ < 0) && (bestObjectIndex2 >= 0)) {
00560 bestObjectIndex_ = bestObjectIndex2;
00561 bestWhichWay_ = bestWhichWay2;
00562 bestTrustedVal1 = bestTrusted2Val1;
00563 returnCode = 4;
00564 }
00565 }
00566
00567 int objectVarInd = -1;
00568 if(bestObjectIndex_ >= 0) {
00569 OsiObject * obj = object[bestObjectIndex_];
00570 obj->setWhichWay(bestWhichWay_);
00571 CouenneObject *co = dynamic_cast <CouenneObject *>(object[bestObjectIndex_]);
00572 if(co) {
00573 objectVarInd = co->Reference()->Index();
00574 }
00575 else {
00576 objectVarInd = obj->columnNumber();
00577 }
00578
00579 message(BRANCH_VAR)<<bestObjectIndex_<< bestWhichWay_ <<CoinMessageEol;
00580
00581 #ifdef TRACE_STRONG
00582 if(problem_->doPrint_) {
00583 if(objectVarInd >= 0) {
00584 double vUb = solver->getColUpper()[objectVarInd];
00585 double vLb = solver->getColLower()[objectVarInd];
00586 double vSolI = info->solution_[objectVarInd];
00587 double vSolS = solver->getColSolution()[objectVarInd];
00588 printf("Branch on object %d (var: %d): solInfo: %10.4f SolSolver: %10.4f low: %10.4f up: %10.4f\n",
00589 bestObjectIndex_, objectVarInd, vSolI, vSolS, vLb, vUb);
00590 }
00591 else {
00592 printf("Branch on object %d (var: -1)\n", bestObjectIndex_);
00593 }
00594 }
00595 #endif
00596 }
00597 message(CHOSEN_VAR)<<bestObjectIndex_<<CoinMessageEol;
00598
00599 if (numberFixed==numberUnsatisfied_&&numberFixed)
00600 returnCode=4;
00601
00602 if((returnCode == 2) || (returnCode == 3)) {
00603 if((objectVarInd > -1) &&
00604 (goodCandidate(solver, info, object, bestObjectIndex_, prec) != 2)) {
00605
00606
00607
00608 returnCode = 4;
00609 }
00610 }
00611 retval = returnCode;
00612 }
00613 else {
00614 retval = 1;
00615 }
00616
00617
00618
00619
00620 #ifdef TRACE_STRONG2
00621 if(problem_->doPrint_) {
00622 if(pv > -1) {
00623 printf("CCS: end: x[%d]: %10.4f lb: %10.4f ub: %10.4f\n",
00624 pv, solver->getColSolution()[pv], solver->getColLower()[pv], solver->getColUpper()[pv]);
00625 printf("CCS: info: x[%d]: %10.4f lb: %10.4f ub: %10.4f\n",
00626 pv, info->solution_[pv], info->lower_[pv], info->upper_[pv]);
00627 printf("CCS: problem: lb: %10.4f ub: %10.4f\n",
00628 problem_->Lb(pv), problem_->Ub(pv));
00629 }
00630 }
00631 #endif
00632
00633 #ifdef TRACE_STRONG
00634 if(problem_->doPrint_) {
00635 printf("CouenneChooseStrong::ChooseVariable(): retval: %d\n", retval);
00636 }
00637 #endif
00638
00639 problem_ -> domain () -> pop ();
00640
00641 return retval;
00642 }
00643
00644
00645
00646
00647 int CouenneChooseStrong::setupList (OsiBranchingInformation *info, bool initialize) {
00648 static bool warned = false;
00649
00650 initialize = true;
00651
00652 problem_ -> domain () -> push
00653 (problem_ -> nVars (),
00654 info -> solution_,
00655 info -> lower_,
00656 info -> upper_);
00657
00658 #ifdef COIN_HAS_NTY
00659
00660 if (problem_ -> orbitalBranching ()) {
00661
00662 problem_ -> ChangeBounds (info -> lower_,
00663 info -> upper_,
00664 problem_ -> nVars ());
00665
00666 problem_ -> Compute_Symmetry();
00667 }
00668
00669 #endif
00670
00671
00672 jnlst_ -> Printf (J_ITERSUMMARY, J_BRANCHING,
00673 "----------------- (strong) setup list\n");
00674
00675 if (jnlst_ -> ProduceOutput (J_DETAILED, J_BRANCHING)) {
00676 for (int i=0; i<problem_ -> domain () -> current () -> Dimension (); i++)
00677 printf ("%4d %20.4g [%20.4g %20.4g]\n", i,
00678 info -> solution_ [i], info -> lower_ [i], info -> upper_ [i]);
00679 }
00680
00681 #ifdef TRACE_STRONG
00682 OsiObject ** object = info->solver_->objects();
00683 #endif
00684
00685 #ifdef TRACE_STRONG
00686 if(problem_->doPrint_) {
00687 printObjViol(info);
00688 }
00689 #endif
00690
00691 int retval = gutsOfSetupList(info, initialize);
00692
00693 if(retval == 0) {
00694
00695 #ifdef FM_CHECKNLP2
00696 if(!(problem_->checkNLP2(info->solution_,
00697 info->objectiveValue_, true,
00698 false,
00699 true,
00700 problem_->getFeasTol()))) {
00701
00702 if (!warned) {
00703 printf("CouenneChooseStrong::setupList(): ### WARNING: checkNLP2() returns infeasible, no branching object selected\n");
00704 warned = true;
00705 }
00706 }
00707 #else
00708 double ckObj = info->objectiveValue_;
00709 if(!(problem_->checkNLP(info->solution_, ckObj, true))) {
00710 if (!warned) {
00711 printf("CouenneChooseStrong::setupList(): ### WARNING: checkNLP() returns infeasible, no branching object selected\n");
00712 warned = true;
00713 }
00714 }
00715 #endif
00716
00717 #ifdef FM_TRACE_OPTSOL
00718 #ifdef FM_CHECKNLP2
00719 problem_->getRecordBestSol()->update();
00720 #else
00721 problem_->getRecordBestSol()->update(info->solution_, problem_->nVars(),
00722 ckObj, problem_->getFeasTol());
00723 #endif
00724 #endif
00725 }
00726
00727 #ifdef TRACE_STRONG
00728 if(problem_->doPrint_) {
00729 printf("Strong list: (obj_ind var_ind priority useful)\n");
00730 printf("numberStrong: %d numberStrongRoot: %d retval: %d\n",
00731 numberStrong_, numberStrongRoot_, retval);
00732 for(int i=0; i<retval; i++) {
00733 CouenneObject *co = dynamic_cast <CouenneObject *>(object[list_[i]]);
00734 int objectInd = -1;
00735 if(co) {
00736 objectInd = co->Reference()->Index();
00737 }
00738 else {
00739 objectInd = object[list_[i]]->columnNumber();
00740 }
00741 printf(" (%d %d %d %6.4f)", list_[i], objectInd,
00742 object[list_[i]]->priority(), useful_[i]);
00743 }
00744 printf("\n");
00745 }
00746 #endif
00747
00748
00749 jnlst_ -> Printf (J_ITERSUMMARY, J_BRANCHING,
00750 "----------------- (strong) setup list done - %d infeasibilities\n", retval);
00751
00752 problem_ -> domain () -> pop ();
00753 return retval;
00754 }
00755
00756
00757
00758
00759
00760
00761
00762
00763
00764
00765
00766
00767
00768
00769
00770
00771
00772
00773
00774
00775
00776
00777
00778
00779
00780
00781
00782
00783
00784
00785
00786
00787
00788
00789
00790 int CouenneChooseStrong::gutsOfSetupList(OsiBranchingInformation *info,
00791 bool initialize)
00792 {
00793 if (numberBeforeTrustedList_ < 0) {
00794 number_not_trusted_ = 1;
00795 printf("CouenneChooseStrong::gutsOfSetupList(): Did not think we were using this; Please double check ...\n");
00796 exit(1);
00797 return OsiChooseVariable::setupList(info, initialize);
00798 }
00799 if (initialize) {
00800 status_=-2;
00801 delete [] goodSolution_;
00802 bestObjectIndex_=-1;
00803 numberStrongDone_=0;
00804 numberStrongIterations_ = 0;
00805 numberStrongFixed_ = 0;
00806 goodSolution_ = NULL;
00807 goodObjectiveValue_ = COIN_DBL_MAX;
00808 number_not_trusted_=0;
00809 }
00810 else {
00811 throw CoinError(CNAME,"setupList","Should not be called with initialize==false");
00812 }
00813 numberOnList_=0;
00814 numberUnsatisfied_=0;
00815 int numberObjects = solver_->numberObjects();
00816 assert (numberObjects);
00817 if (numberObjects>pseudoCosts_.numberObjects()) {
00818
00819
00820
00821
00822
00823 int saveNumberBeforeTrusted = pseudoCosts_.numberBeforeTrusted();
00824 pseudoCosts_.initialize(numberObjects);
00825 pseudoCosts_.setNumberBeforeTrusted(saveNumberBeforeTrusted);
00826 }
00827 int bestPriority = COIN_INT_MAX;
00828 int i;
00829 #ifdef FM_SORT_STRONG
00830 int numStr = numberStrong_;
00831 if(isRootNode(info)) {
00832 numStr = numberStrongRoot_;
00833 }
00834 int maximumStrong = CoinMin(numStr, numberObjects) ;
00835 int lastPrio = problem_->getLastPrioSort();
00836 int card_vPriority = 0;
00837 int posEnd_vPriority = numberObjects;
00838 double *vPriority = new double[numberObjects];
00839 double *infeasVal = new double[numberObjects];
00840 int max_most_fra = setup_pseudo_frac_ > 0. ? (int)floor(setup_pseudo_frac_*(double)maximumStrong): 0;
00841 if (setup_pseudo_frac_ > 0.) {
00842 max_most_fra = CoinMax(1, max_most_fra);
00843 }
00844 #else
00845 int putOther = numberObjects;
00846 double check = -COIN_DBL_MAX;
00847 int checkIndex = 0;
00848 int maximumStrong = CoinMin(CoinMax(numberStrong_,numberStrongRoot_),
00849 numberObjects) ;
00850 for (i=0;i<numberObjects;i++) {
00851 list_[i]=-1;
00852 useful_[i]=0.0;
00853 }
00854
00855 int* list2 = NULL;
00856 double* useful2 = NULL;
00857 double check2 = -COIN_DBL_MAX;
00858 int checkIndex2=0;
00859 int max_most_fra = setup_pseudo_frac_ > 0. ? (int)floor(setup_pseudo_frac_*(double)maximumStrong): 0;
00860 if (setup_pseudo_frac_ > 0.) {
00861 max_most_fra = CoinMax(1, max_most_fra);
00862 }
00863 if (max_most_fra) {
00864 list2 = new int[max_most_fra];
00865 useful2 = new double[max_most_fra];
00866 for (i=0;i<max_most_fra;i++) {
00867 list2[i]=-1;
00868 useful2[i]=0.0;
00869 }
00870 }
00871 #endif
00872
00873 #ifdef FM_CHECK
00874 const double* upTotalChange = pseudoCosts_.upTotalChange();
00875 const double* downTotalChange = pseudoCosts_.downTotalChange();
00876 int pseudoNum = pseudoCosts_.numberObjects();
00877 for(i=0; i<pseudoNum; i++) {
00878 if(isnan(upTotalChange[i]) || isinf(upTotalChange[i])) {
00879 printf("CouenneChooseStrong::gutsOfSetupList(): upTotalChange[%d]: not a number or infinite\n", i);
00880 exit(1);
00881 }
00882 if(isnan(downTotalChange[i]) || isinf(downTotalChange[i])) {
00883 printf("CouenneChooseStrong::gutsOfSetupList(): downTotalChange[%d]: not a number or infinite\n", i);
00884 exit(1);
00885 }
00886 }
00887 #endif
00888
00889 OsiObject ** object = info->solver_->objects();
00890 double upMultiplier, downMultiplier;
00891 computeMultipliers(upMultiplier, downMultiplier);
00892
00893
00894 bool feasible = true;
00895 const double MAXMIN_CRITERION = maxminCrit(info);
00896
00897 bool firstPass = false;
00898
00899
00900 while(numberOnList_ == 0) {
00901 for(i=0;i<numberObjects;i++) {
00902 int way;
00903 double value = object[i]->infeasibility(info, way);
00904
00905 #ifdef FM_SORT_STRONG
00906 infeasVal[i] = value;
00907 #endif
00908
00909 double lbForInfeas = 0.0;
00910 if(value > lbForInfeas) {
00911 numberUnsatisfied_++;
00912 if(value >= 1e50) {
00913
00914 feasible=false;
00915 break;
00916 }
00917 int priorityLevel = object[i]->priority();
00918
00919 #ifdef FM_SORT_STRONG
00920 if(priorityLevel < bestPriority) {
00921 bestPriority = priorityLevel;
00922 }
00923 if(priorityLevel > lastPrio) {
00924 posEnd_vPriority--;
00925 vPriority[posEnd_vPriority] = priorityLevel;
00926 list_[posEnd_vPriority] = i;
00927 }
00928 else {
00929 vPriority[card_vPriority] = priorityLevel;
00930 list_[card_vPriority] = i;
00931 card_vPriority++;
00932 }
00933 #else
00934
00935 if(priorityLevel < bestPriority) {
00936 for (int j=maximumStrong-1; j>=0; j--) {
00937 if(list_[j] >= 0) {
00938 int iObject = list_[j];
00939 list_[j]=-1;
00940 useful_[j]=0.0;
00941 list_[--putOther]=iObject;
00942 }
00943 }
00944 maximumStrong = CoinMin(maximumStrong,putOther);
00945 bestPriority = priorityLevel;
00946 check = -COIN_DBL_MAX;
00947 checkIndex = 0;
00948 check2 = -COIN_DBL_MAX;
00949 checkIndex2 = 0;
00950 number_not_trusted_ = 0;
00951 if(max_most_fra > 0) {
00952 for(int j=0; j<max_most_fra; j++) {
00953 list2[j]=-1;
00954 useful2[j]=0.0;
00955 }
00956 }
00957 }
00958 if(priorityLevel == bestPriority) {
00959
00960 double value2;
00961 value = computeUsefulness(MAXMIN_CRITERION,
00962 upMultiplier, downMultiplier, value,
00963 object[i], i, value2);
00964 if(value > check) {
00965
00966 int iObject = list_[checkIndex];
00967 if(iObject >= 0) {
00968 assert (list_[putOther-1]<0);
00969 list_[--putOther]=iObject;
00970 }
00971 list_[checkIndex]=i;
00972 assert (checkIndex<putOther);
00973 useful_[checkIndex]=value;
00974
00975 check=COIN_DBL_MAX;
00976 maximumStrong = CoinMin(maximumStrong,putOther);
00977 for (int j=0; j<maximumStrong; j++) {
00978 if(list_[j]>=0) {
00979 if (useful_[j]<check) {
00980 check=useful_[j];
00981 checkIndex=j;
00982 }
00983 }
00984 else {
00985 check=0.0;
00986 checkIndex = j;
00987 break;
00988 }
00989 }
00990 }
00991 else {
00992
00993 assert (list_[putOther-1]<0);
00994 list_[--putOther]=i;
00995 maximumStrong = CoinMin(maximumStrong,putOther);
00996 }
00997
00998 if((max_most_fra > 0) && (value2 > check2)) {
00999
01000 number_not_trusted_++;
01001 list2[checkIndex2]=i;
01002 useful2[checkIndex2]=value2;
01003
01004 check2=COIN_DBL_MAX;
01005 for(int j=0; j<max_most_fra; j++) {
01006 if(list2[j] >= 0) {
01007 if(useful2[j] < check2) {
01008 check2=useful2[j];
01009 checkIndex2=j;
01010 }
01011 }
01012 else {
01013 check2=0.0;
01014 checkIndex2 = j;
01015 break;
01016 }
01017 }
01018 }
01019 }
01020 else {
01021
01022
01023 assert (list_[putOther-1]<0);
01024 list_[--putOther]=i;
01025 maximumStrong = CoinMin(maximumStrong,putOther);
01026 }
01027 #endif
01028 }
01029 }
01030
01031 #ifdef FM_SORT_STRONG
01032
01033 #ifdef FM_CHECK
01034 if(card_vPriority - posEnd_vPriority + numberObjects != numberUnsatisfied_) {
01035 printf("CouenneChooseStrong::gutsOfSetupList(): ### ERROR: card_vPriority: %d posEnd_vPriority: %d numberUnsatisfied: %d numberObjects: %d\n",
01036 card_vPriority, posEnd_vPriority, numberUnsatisfied_, numberObjects);
01037 exit(1);
01038 }
01039 #endif
01040
01041 numberOnList_ = 0;
01042 if(feasible) {
01043 int card_smallerThanPrio = card_vPriority;
01044 if(posEnd_vPriority > card_vPriority) {
01045 for(i=posEnd_vPriority; i<numberObjects; i++) {
01046 list_[card_vPriority] = list_[i];
01047 list_[i] = -1;
01048 vPriority[card_vPriority] = vPriority[i];
01049 card_vPriority++;
01050 }
01051 }
01052 else {
01053 card_vPriority = numberUnsatisfied_;
01054 }
01055
01056 int sortFrom = 0;
01057 int sortUpTo = card_smallerThanPrio;
01058 if(card_smallerThanPrio < maximumStrong) {
01059 sortFrom = card_smallerThanPrio;
01060 sortUpTo = card_vPriority;
01061 }
01062 if(card_vPriority > 0) {
01063 numberOnList_ = (card_vPriority < maximumStrong ? card_vPriority : maximumStrong);
01064
01065 #ifdef FM_ALWAYS_SORT
01066 bool alwaysSort = true;
01067 #else
01068 bool alwaysSort = false;
01069 #endif
01070 if(alwaysSort) {
01071 sortFrom = 0;
01072 sortUpTo = card_vPriority;
01073 }
01074 if((sortUpTo > maximumStrong) || alwaysSort){
01075
01076 CoinSort_2(vPriority + sortFrom, vPriority + sortUpTo,
01077 list_ + sortFrom);
01078 }
01079 for(i=0; i<card_vPriority; i++) {
01080 int indObj = list_[i];
01081 double value = 0, value2;
01082 value = computeUsefulness(MAXMIN_CRITERION,
01083 upMultiplier, downMultiplier, value,
01084 object[indObj], indObj, value2);
01085
01086 #ifdef OLD_USEFULLNESS
01087 useful_[i] = -value;
01088 #else
01089 if ((sortCrit_ & 1) == 0) {
01090 useful_[i] = -value;
01091 }
01092 else {
01093 useful_[i] = value;
01094 }
01095 #endif
01096
01097 #ifdef USE_NOT_TRUSTED
01098 if(value2 < -COIN_DBL_MAX / 10) {
01099 infeasVal[indObj] = -COIN_DBL_MAX;
01100 }
01101 #endif
01102 }
01103
01104 #ifdef USE_NOT_TRUSTED
01105
01106 if((card_vPriority > maximumStrong) &&
01107 (vPriority[maximumStrong] < bestPriority + COUENNE_EPS)) {
01108
01109
01110 int cardFrac = 0;
01111 int *fracInd = new int[card_vPriority];
01112
01113 double *fracVal = new double[card_vPriority];
01114
01115
01116 for(i=0; i<card_vPriority; i++) {
01117 int indObj = list_[i];
01118 if(vPriority[i] < bestPriority + COUENNE_EPS) {
01119 if(infeasVal[indObj] > -COIN_DBL_MAX/10) {
01120 fracInd[cardFrac] = i;
01121 fracVal[cardFrac] = -infeasVal[indObj];
01122 cardFrac++;
01123 }
01124 }
01125 }
01126 if(max_most_fra > 0) {
01127
01128
01129 if(cardFrac > max_most_fra) {
01130 CoinSort_2(fracVal, fracVal+cardFrac, fracInd);
01131 }
01132 for(i=0; i<cardFrac; i++) {
01133 useful_[fracInd[i]] =
01134 -1e150*(1. + infeasVal[list_[fracInd[i]]]);
01135 number_not_trusted_++;
01136 if(i == max_most_fra - 1) {
01137 break;
01138 }
01139 }
01140 }
01141 delete[] fracInd;
01142 delete[] fracVal;
01143 }
01144 #endif
01145 }
01146
01147 #ifdef FM_SEC_SORT_USEFUL
01148 CoinSort_2(useful_, useful_+card_vPriority, list_);
01149 #else
01150 #ifdef FM_ALWAYS_SORT
01151 int from = 0, upto = 1;
01152 while(upto < card_vPriority) {
01153 while(vPriority[upto] == vPriority[from]) {
01154 upto++;
01155 if(upto == card_vPriority) {
01156 break;
01157 }
01158 }
01159 CoinSort_2(useful_+from, useful_+upto, list_+from);
01160 from = upto;
01161 upto = from+1;
01162 }
01163 #else
01164 if(sortUpTo > maximumStrong) {
01165
01166
01167 int from = maximumStrong-1, upto = maximumStrong;
01168 int msPrio = vPriority[maximumStrong-1];
01169 problem_->setLastPrioSort(msPrio);
01170 while((from > -1) && (vPriority[from] == msPrio)) {
01171 from--;
01172 }
01173 from++;
01174 while((upto < sortUpTo) && (vPriority[upto] == msPrio)) {
01175 upto++;
01176 }
01177
01178
01179 CoinSort_2(useful_+from, useful_+upto, list_+from);
01180 }
01181
01182 #endif
01183
01184 #ifdef FM_CHECK
01185
01186 double ckPrio = (card_vPriority < numberUnsatisfied_ ?
01187 vPriority[card_vPriority] : 100000);
01188 double ckUse = (card_vPriority < numberUnsatisfied_ ?
01189 useful_[card_vPriority] : 100000);
01190 for(i=0; i<card_vPriority; i++) {
01191 int indObj = list_[i];
01192 if(object[indObj]->priority() > ckPrio + 1e-3) {
01193 printf("CouenneChooseStrong::gutsOfSetupList(): ### ERROR: object[%d]->priority(): %d > ckPrio: %d\n",
01194 indObj, object[indObj]->priority(), ckPrio);
01195 exit(1);
01196 }
01197 if(fabs(object[indObj]->priority() - ckPrio) < 1e-3) {
01198 if(useful_[i] > ckUse + 1e-3) {
01199 printf("CouenneChooseStrong::gutsOfSetupList(): ### ERROR: object[%d]->useful: %f > ckUse: %d\n",
01200 indObj, useful_[i], ckUse);
01201 exit(1);
01202 }
01203 }
01204 }
01205 for(i=card_vPriority; i<numberUnsatisfied_; i++) {
01206 int indObj = list_[i];
01207 if(object[indObj]->priority() < ckPrio - 1e-3) {
01208 printf("CouenneChooseStrong::gutsOfSetupList(): ### ERROR: object[%d]->priority(): %d < ckPrio: %d\n",
01209 indObj, object[indObj]->priority(), ckPrio);
01210 exit(1);
01211 }
01212 if(fabs(object[indObj]->priority() - ckPrio) < 1e-3) {
01213 if(useful_[i] < ckUse - 1e-3) {
01214 printf("CouenneChooseStrong::gutsOfSetupList(): ### ERROR: object[%d]->useful: %f < ckUse: %d\n",
01215 indObj, useful_[i], ckUse);
01216 exit(1);
01217 }
01218 }
01219 }
01220 #endif
01221 #endif
01222 }
01223 else {
01224 numberUnsatisfied_ = -1;
01225 }
01226 #else
01227
01228 numberOnList_=0;
01229 if (feasible) {
01230 maximumStrong = CoinMin(maximumStrong,putOther);
01231 for (i=0;i<maximumStrong;i++) {
01232 if (list_[i]>=0) {
01233 #ifdef OLD_USEFULLNESS
01234 list_[numberOnList_]=list_[i];
01235 useful_[numberOnList_++]=-useful_[i];
01236
01237 #else
01238 list_[numberOnList_]=list_[i];
01239 if ((sortCrit_ & 1) == 0) {
01240 useful_[numberOnList_++]=-useful_[i];
01241 }
01242 else useful_[numberOnList_++] = useful_[i];
01243 #endif
01244 message(CANDIDATE_LIST2)<<numberOnList_-1
01245 <<list_[numberOnList_-1]<<numberOnList_-1<<useful_[numberOnList_-1]
01246 <<CoinMessageEol;
01247 }
01248 }
01249 if (numberOnList_) {
01250 int tmp_on_list = 0;
01251 if (max_most_fra > 0 && numberOnList_ >= maximumStrong) {
01252
01253
01254 number_not_trusted_=0;
01255 for (i=0;i<max_most_fra;i++) {
01256 if (list2[i]>=0) {
01257 list2[number_not_trusted_] = list2[i];
01258 useful2[number_not_trusted_++] = useful2[i];
01259 message(CANDIDATE_LIST3)<<number_not_trusted_-1
01260 <<list2[number_not_trusted_-1]<<number_not_trusted_-1
01261 <<useful2[number_not_trusted_-1]<<CoinMessageEol;
01262 }
01263 }
01264 if (number_not_trusted_) {
01265 CoinSort_2(list_,list_+numberOnList_,useful_);
01266 CoinSort_2(list2,list2+number_not_trusted_,useful2);
01267 int i1=0;
01268 int i2=0;
01269 for (i=0; i<numberObjects; i++) {
01270 bool found1 = (list_[i1]==i);
01271 bool found2 = (list2[i2]==i);
01272 if (found1 && found2) {
01273
01274 #ifdef OLD_USEFULLNESS
01275 useful_[i1] = 1e150*(1.+useful2[i2]);
01276 #else
01277 useful_[i1] = -1e150*(1.+useful2[i2]);
01278 #endif
01279 list2[i2] = -1;
01280 }
01281 if (found1) i1++;
01282 if (found2) i2++;
01283 if (i2==max_most_fra) break;
01284 }
01285 for (i=0; i<number_not_trusted_; i++) {
01286 if (list2[i] >= 0) {
01287 list_[numberOnList_+tmp_on_list] = list2[i];
01288
01289 #ifdef OLD_USEFULLNESS
01290 useful_[numberOnList_+tmp_on_list] = 1e150*(1.+useful2[i]);
01291 #else
01292 useful_[numberOnList_+tmp_on_list] = -1e150*(1.+useful2[i]);
01293 #endif
01294 tmp_on_list++;
01295 }
01296 }
01297 }
01298 }
01299
01300 CoinSort_2(useful_,useful_+numberOnList_+tmp_on_list,list_);
01301
01302 i = numberOnList_;
01303 for (;putOther<numberObjects;putOther++)
01304 list_[i++]=list_[putOther];
01305 assert (i==numberUnsatisfied_);
01306 if (!CoinMax(numberStrong_,numberStrongRoot_))
01307 numberOnList_=0;
01308 }
01309 }
01310 else {
01311
01312 numberUnsatisfied_=-1;
01313 }
01314 #endif
01315
01316 if(!firstPass) {
01317 break;
01318 }
01319 firstPass = false;
01320 }
01321
01322 #ifdef TRACE_STRONG
01323 if(problem_->doPrint_) {
01324 printf("numberStrong_: %d maximumStrong: %d\n",
01325 numberStrong_, maximumStrong);
01326 }
01327 #endif
01328
01329 #ifdef FM_SORT_STRONG
01330 delete [] vPriority;
01331 delete [] infeasVal;
01332 #else
01333 delete [] list2;
01334 delete [] useful2;
01335 #endif
01336
01337
01338 info->defaultDual_ = -1.0;
01339 delete [] info->usefulRegion_;
01340 delete [] info->indexRegion_;
01341
01342 int way;
01343 if (bb_log_level_>3) {
01344
01345 for (int i=0; i<numberOnList_; i++){
01346 message(CANDIDATE_LIST)<<i<< list_[i]<< i<< useful_[i]
01347 <<object[list_[i]]->infeasibility(info,way)
01348 <<CoinMessageEol;
01349 }
01350 }
01351
01352
01353
01354 if(numberUnsatisfied_ == -1) {
01355 return(-1);
01356 }
01357 return numberOnList_;
01358 }
01359
01360
01362 void CouenneChooseStrong::registerOptions (Ipopt::SmartPtr <Bonmin::RegisteredOptions> roptions) {
01363
01364 roptions -> AddStringOption6
01365 ("pseudocost_mult",
01366 "Multipliers of pseudocosts for estimating and update estimation of bound",
01367 "interval_br_rev",
01368
01369 "infeasibility", "infeasibility returned by object",
01370
01371 "projectDist", "distance between current LP point and resulting branches' LP points",
01372
01373 "interval_lp", "width of the interval between bound and current lp point",
01374 "interval_lp_rev", "similar to interval_lp, reversed",
01375
01376 "interval_br", "width of the interval between bound and branching point",
01377 "interval_br_rev", "similar to interval_br, reversed");
01378
01379 roptions -> AddStringOption2
01380 ("pseudocost_mult_lp",
01381 "Use distance between LP points to update multipliers of pseudocosts "
01382 "after simulating branching",
01383 "no",
01384 "yes", "",
01385 "no", "");
01386
01387 roptions -> AddStringOption2
01388 ("estimate_select",
01389 "How the min/max estimates of the subproblems' bounds are used in strong branching",
01390 "normal",
01391 "normal", "as usual in literature",
01392 "product", "use their product");
01393
01394 roptions -> AddStringOption2
01395 ("trust_strong",
01396 "Fathom strong branching LPs when their bound is above the cutoff",
01397 "yes",
01398 "yes", "",
01399 "no", "");
01400 }
01401
01402
01403
01404 bool CouenneChooseStrong::feasibleSolution (const OsiBranchingInformation * info,
01405 const double * solution,
01406 int numberObjects,
01407 const OsiObject ** objects) {
01408
01409 #ifdef FM_CHECKNLP2
01410 bool isFeas = problem_->checkNLP2(solution, 0, false, true, true,
01411 problem_->getFeasTol());
01412 return isFeas;
01413 #else
01414 double obj = solution [problem_ -> Obj (0) -> Body () -> Index ()];
01415 return problem_ -> checkNLP (solution, obj);
01416 #endif
01417 }
01418
01419 void CouenneChooseStrong::printObjViol(OsiBranchingInformation *info) {
01420
01421 OsiObject ** object = info->solver_->objects();
01422 int numberObjects = info->solver_->numberObjects();
01423
01424 printf("CouenneChooseStrong::printObjViol(): Object violations: (obj_ind var_ind violation)");
01425 double maxViol = 0;
01426 double minPosViol = 1e50;
01427 for(int i=0; i<numberObjects; i++) {
01428 CouenneObject *co = dynamic_cast <CouenneObject *>(object[i]);
01429 int indVar = -1;
01430 if(co) {
01431 indVar = co->Reference()->Index();
01432 }
01433 else {
01434 indVar = object[i]->columnNumber();
01435 }
01436 int way;
01437 double value = object[i]->infeasibility(info,way);
01438 maxViol = (value > maxViol ? value : maxViol);
01439 if(value > 0.0) {
01440 printf("(%d %d %f)", i, indVar, value);
01441 minPosViol = (value < minPosViol ? value : minPosViol);
01442 }
01443 }
01444 printf("\nmaxViol: %g minPosViol: %g\n", maxViol, minPosViol);
01445
01446 }
01447
01448
01449