source: trunk/Couenne/src/branch/CouenneChooseStrong.cpp @ 967

Last change on this file since 967 was 967, checked in by fmargot, 7 years ago

fix bug in checkNLP2 and handling of non-var objects in strong branching

  • Property svn:keywords set to Author Date Id Revision
File size: 29.0 KB
Line 
1/* $Id: CouenneChooseStrong.cpp 967 2013-06-24 10:13:05Z fmargot $
2 *
3 * Name:    CouenneChooseStrong.cpp
4 * Authors: Andreas Waechter, IBM Corp.
5 *          Pietro Belotti, Lehigh University
6 *          Francois Margot, Carnegie Mellon University
7 * Purpose: Strong branching objects for Couenne
8 *
9 * (C) Carnegie-Mellon University, 2006-11.
10 * This file is licensed under the Eclipse Public License (EPL)
11 */
12
13#include "BonChooseVariable.hpp"
14
15#include "CouenneObject.hpp"
16#include "CouenneChooseStrong.hpp"
17#include "CouenneProblem.hpp"
18#include "CouenneProblemElem.hpp"
19#include "CouenneBranchingObject.hpp"
20#include "CouenneRecordBestSol.hpp"
21
22// The recommended ones:
23#define FM_SORT_STRONG
24#define FM_SEC_SORT_USEFUL
25#define USE_NOT_TRUSTED
26
27//#define TRACE_STRONG
28//#define TRACE_STRONG2
29//#define FM_ALWAYS_SORT
30//#define USE_SMALL_GAP
31//#define OLD_STYLE
32
33using namespace Ipopt;
34using namespace Couenne;
35
36const CouNumber estProdEps = 1e-6;
37
38#ifdef COIN_HAS_NTY
39#include "Nauty.h"
40#endif
41
42  /// constructor
43  CouenneChooseStrong::CouenneChooseStrong (Bonmin::BabSetupBase &b, CouenneProblem* p, JnlstPtr jnlst) :
44
45  Bonmin::BonChooseVariable (b, b.continuousSolver()),
46  problem_          (p),
47  jnlst_            (jnlst),
48  branchtime_       (0.) {
49
50    std::string s;
51
52    b.options () -> GetStringValue ("pseudocost_mult_lp", s, "couenne.");
53    pseudoUpdateLP_ = (s == "yes");     
54
55    b.options () -> GetStringValue ("estimate_select", s, "couenne.");
56    estimateProduct_ = (s == "product");
57
58    b.options () -> GetStringValue ("trust_strong", s, "couenne.");
59
60    // trust solution from strong branching to provide right LB
61
62    setTrustStrongForSolution (s == "yes");
63    setTrustStrongForBound    (s == "yes");
64  }
65
66  /// copy constructor
67  CouenneChooseStrong::CouenneChooseStrong (const CouenneChooseStrong& rhs) :
68    Bonmin::BonChooseVariable (rhs),
69    problem_          (rhs.problem_),
70    pseudoUpdateLP_   (rhs.pseudoUpdateLP_),
71    estimateProduct_  (rhs.estimateProduct_),
72    jnlst_            (rhs.jnlst_),
73    branchtime_       (rhs.branchtime_)
74  {}
75
76  /// destructor
77  CouenneChooseStrong::~CouenneChooseStrong()
78  {if (branchtime_ > 1e-9) jnlst_ -> Printf (J_ERROR, J_BRANCHING, "Strong branching: total time %g\n", branchtime_);}
79
80  /// cloning method
81  OsiChooseVariable *
82  CouenneChooseStrong::clone() const
83  {return new CouenneChooseStrong(*this);}
84
85  /// assignment operator
86  CouenneChooseStrong&
87  CouenneChooseStrong::operator=(const CouenneChooseStrong & rhs)
88  {
89    if (this != &rhs) {
90      Bonmin::BonChooseVariable::operator=(rhs);
91      problem_         = rhs.problem_;
92      pseudoUpdateLP_  = rhs.pseudoUpdateLP_;
93      estimateProduct_ = rhs.estimateProduct_;
94      jnlst_           = rhs.jnlst_;
95      branchtime_      = rhs.branchtime_;
96    }
97    return *this;
98  }
99
100/***********************************************************************/
101// returns
102// 4 if not a variable object (deemed good)
103// 3 if variable object and is good
104// 2 if variable object and is not good (i.e. fixed variable)
105// 1 if OsiSimpleBranchingObject and not good (deprecated)
106int CouenneChooseStrong::goodCandidate(OsiSolverInterface *solver,
107                                       OsiBranchingInformation *info,
108                                       OsiObject **object, const int iObject,
109                                       const double prec) {
110
111  int vInd = object [iObject] -> columnNumber ();
112
113  if (vInd < 0) return 4;  // not a variable object, so deem it good
114
115  bool varIsInt = solver -> isInteger (vInd);
116
117  // all object now are CouenneObjects or derivates (CouenneVarObject,
118  // CouenneSOSObject, etc.)
119
120  //CouenneObject    *co    = dynamic_cast <CouenneObject    *> (object [iObject]);
121  //OsiSimpleInteger *simpl = dynamic_cast <OsiSimpleInteger *> (object [iObject]);
122 
123  // int vInd = -1;
124  // bool varIsInt = false;
125  // if(co) {
126  //   vInd = co->Reference()->Index();
127  //   if(vInd >= 0) {
128  //     varIsInt = co->Reference()->isInteger();
129  //   }
130  // }
131  // else {
132  //   if(simpl) {
133  //     vInd = object[iObject]->columnNumber();
134  //     varIsInt = true;
135  //   }
136  //   else {
137
138  //     // printf("CouenneChooseStrong::goodCandidate: ### ERROR: unknown object\n");
139  //     // exit(1);
140
141  //     // this is probably a SOS object, and anyhow we don't want to
142  //     // exit on it.
143
144  //     return 3;
145  //   }
146  // }
147 
148  int goodCand = 3; // good candidate
149 
150  // Must not call branch() for integer variable vInd with
151  // upper == lower or for OsiSimpleInteger with
152  // info->solution[vInd] not between lower and upper
153  if((vInd >= 0) && varIsInt) {
154    double vUp = solver->getColUpper()[vInd];
155    double vLow = solver->getColLower()[vInd];
156    double infoVal = info->solution_[vInd];
157    double distToInt = fabs(infoVal - floor(infoVal + 0.5));
158    if(distToInt > 0.5) {
159      distToInt = 1 - distToInt;
160    }
161    // if(simpl) {
162    //   goodCand = 0; // bad candidate
163    //   if((distToInt > info->integerTolerance_) &&
164    //      (vUp > vLow + prec)) {
165    //     goodCand = 3; // good candidate
166    //     if((vUp + prec < infoVal) || (infoVal < vLow - prec)) {
167    //       goodCand = 1; // bad candidate
168    //     }
169    //   }
170    // }
171    //if(co) {
172      goodCand = 2;
173      if(vUp > vLow + prec) {
174        goodCand = 3; // good candidate
175      }
176      //}
177  }
178   
179  return(goodCand);
180} /* goodCandidate */
181
182/***********************************************************************/
183bool CouenneChooseStrong::saveBestCand(OsiObject **object, const int iObject, 
184                                       const double value, 
185                                       const double upEstimate, 
186                                       const double downEstimate,
187                                       double &bestVal1, 
188                                       double &bestVal2, int &bestIndex,
189                                       int &bestWay) {
190  bool retval = false;
191  if(value > bestVal1) {
192    retval = true;
193    bestVal1 = value;
194    bestIndex = iObject;
195    bestWay = upEstimate > downEstimate ? 0 : 1;
196    // but override if there is a preferred way
197    const OsiObject * obj = object[iObject];
198    if (obj->preferredWay() >= 0 && obj->infeasibility()) {
199      bestWay = obj->preferredWay();
200    }
201  }
202  return(retval);
203} /* saveBestCand */
204
205/*******************************************************************/
206  /* Choose a variable. Returns:
207
208    -1  Node is infeasible
209     0  Normal termination - we have a candidate
210     1  All look satisfied - no candidate
211     2  We can change the bound on a variable - but we also have a strong branching candidate
212     3  We can change the bound on a variable - but we have a non-strong branching candidate
213     4  We can change the bound on a variable - no other candidates
214
215     We can pick up branch from whichObject() and whichWay()
216     We can pick up a forced branch (can change bound) from whichForcedObject() and whichForcedWay()
217     If we have a solution then we can pick up from goodObjectiveValue() and goodSolution()
218
219Note: Must not return 4 for a non variable object
220  */
221  int CouenneChooseStrong::chooseVariable (OsiSolverInterface * solver,
222                                           OsiBranchingInformation *info,
223                                           bool fixVariables) {
224
225    /// Note: had to copy code from
226    /// BonChooseVariable::chooseVariable() in order to test product
227    /// thing
228
229    problem_ -> domain () -> push
230      (problem_ -> nVars (),
231       info -> solution_,
232       info -> lower_,
233       info -> upper_);
234
235#ifdef TRACE_STRONG2
236    int pv = -1;
237    if(problem_->doPrint_) {
238      if(pv > -1) {
239        printf("CCS: beg: x[%d]: %10.4f  lb: %10.4f  ub: %10.4f\n",
240               pv, solver->getColSolution()[pv], solver->getColLower()[pv], solver->getColUpper()[pv]);
241        printf("CCS: info: x[%d]: %10.4f  lb: %10.4f  ub: %10.4f\n",
242               pv, info->solution_[pv], info->lower_[pv], info->upper_[pv]);
243        printf("CCS: problem: lb: %10.4f  ub: %10.4f\n",
244               problem_->Lb(pv), problem_->Ub(pv));
245      }
246    }
247 #endif                                                                         
248
249    int retval;
250    const double prec = problem_->getFeasTol();
251
252    //int retval = BonChooseVariable::chooseVariable (solver, info, fixVariables);
253
254    // COPY of Bonmin starts here ////////////////////////////////////////////
255
256    // We assume here that chooseVariable is called once at the very
257    // beginning with fixVariables set to true.  This is then the root
258    // node.
259
260    bool isRoot = isRootNode(info);
261    int numberStrong = numberStrong_;
262
263    if (isRoot) {
264      numberStrong = numberStrongRoot_;
265    }
266
267    if (numberUnsatisfied_) {
268      int cardIndForPseudo = 0, 
269        *indForPseudo = new int[numberUnsatisfied_];
270      OsiObject ** object = solver->objects();
271      const double* upTotalChange = pseudoCosts_.upTotalChange();
272      const double* downTotalChange = pseudoCosts_.downTotalChange();
273      const int* upNumber = pseudoCosts_.upNumber();
274      const int* downNumber = pseudoCosts_.downNumber();
275      int numberBeforeTrusted = pseudoCosts_.numberBeforeTrusted();
276
277      // number of objects to be chosen by doStrongBranching()
278      int numberLeft = CoinMin (numberStrong - numberStrongDone_, numberUnsatisfied_);
279
280      results_.clear();
281
282      // All look satisfied;
283      //useful if all objects considered for strong branching use pseudo costs
284      int returnCode = 1;
285
286      int returnCodeSB = 0;
287      bestObjectIndex_ = -1;
288      bestWhichWay_ = -1;
289      firstForcedObjectIndex_ = -1;
290      firstForcedWhichWay_ =-1;
291      double bestTrustedVal1 = -COIN_DBL_MAX;
292      double bestTrustedVal2 = -COIN_DBL_MAX;
293
294      bool smallGap = false;
295      bool sbObjPosImp = false; // true if an object on which strong branching
296                                // was performed has a positive improvement
297                                // in both branches; used only when gap is
298                                // deemed small
299
300#ifdef USE_SMALL_GAP
301      int objInd = problem_ -> Obj (0) -> Body () -> Index ();
302      double lbGap = objInd >= 0 ? info -> lower_ [objInd] : problem_ -> Obj (0) -> Body () -> Value ();
303      double ubGap = problem_ -> getCutOff ();
304      double currentGap = 
305        (ubGap >  COUENNE_INFINITY    / 10 ||
306         lbGap < -Couenne_large_bound / 10) ? 1e3 : 
307        fabs (ubGap - lbGap) / (1.e-3 + CoinMin (fabs (ubGap), fabs (lbGap)));
308
309      if(currentGap < 1e-3) {
310        smallGap = true;
311      }
312#endif
313
314#ifdef TRACE_STRONG
315      if((problem_->doPrint_) && (number_not_trusted_ > 0)) {
316        printf("number_not_trusted: %d\n", number_not_trusted_);
317      }
318#endif
319
320      for (int i=0;i<numberLeft;i++) {
321        int iObject = list_[i];
322        if (numberBeforeTrusted == 0||
323            i < minNumberStrongBranch_ ||
324            (
325              only_pseudo_when_trusted_ && number_not_trusted_>0 ) ||
326              ( !isRoot && (upNumber[iObject]<numberBeforeTrusted ||
327                          downNumber[iObject]<numberBeforeTrusted ))||
328              ( isRoot && (!upNumber[iObject] && !downNumber[iObject])) ) {
329         
330#ifdef TRACE_STRONG
331          if(problem_->doPrint_) {
332            printf("Push object %d for strong branch\n", iObject);
333          }
334#endif
335          results_.push_back (Bonmin::HotInfo (solver, info, object, iObject));
336        }
337        else {
338
339#ifdef TRACE_STRONG
340          if(problem_->doPrint_) {
341            printf("Use pseudo cost for object %d\n", iObject);
342          }
343#endif
344          indForPseudo[cardIndForPseudo] = iObject;
345          cardIndForPseudo++;
346        }
347      }
348
349      int numberFixed=0;
350
351      if (results_.size() > 0) {
352
353        //
354        // do strong branching
355        //
356
357        returnCodeSB = doStrongBranching (solver, info, results_.size(), 1);
358
359        if (bb_log_level_>=3) {
360          const char* stat_msg[] = {"NOTDON", "FEAS", "INFEAS", "NOFINI"};
361          message(SB_HEADER)<<CoinMessageEol;
362          for (unsigned int i = 0; i< results_.size(); i++) {
363            double up_change = results_[i].upChange();
364            double down_change = results_[i].downChange();
365            int up_status = results_[i].upStatus();
366            int down_status = results_[i].downStatus();
367            message(SB_RES)<<(int) i<<stat_msg[down_status+1]<<down_change
368            <<stat_msg[up_status+1]<< up_change<< CoinMessageEol;
369          }
370        }
371
372        if (returnCodeSB >= 0 && returnCodeSB <= 2) { // 1, 2: some can be fixed
373          if(returnCodeSB > 0) {
374            returnCode = 4; // no bestObject yet
375          }
376          else {
377            returnCode = 0;
378          }
379          for (unsigned int i=0;i < results_.size (); i++) {
380
381            if((results_[i].upStatus() < 0) || (results_[i].downStatus() < 0))
382              continue;
383
384            // if((results_[i].upStatus() < 0) || (results_[i].downStatus() < 0)) {
385            //   continue;
386            // }
387
388            int iObject = results_[i].whichObject();
389            const OsiObject * obj = object[iObject];
390            int needBranch = goodCandidate(solver, info, object, iObject, prec);
391
392            ///////////////////////////////////////////////////////////////////
393
394            double upEstimate;
395
396            if (results_[i].upStatus()!=1) {
397              assert (results_[i].upStatus()>=0);
398              upEstimate = results_[i].upChange();
399            }
400            else {
401              // infeasible - just say expensive
402              if (info->cutoff_<1.0e50)
403                upEstimate = 2.0*(info->cutoff_-info->objectiveValue_);
404              else
405                upEstimate = 2.0*fabs(info->objectiveValue_);
406              if (firstForcedObjectIndex_ <0) {
407                // first fixed variable
408                firstForcedObjectIndex_ = iObject;
409                firstForcedWhichWay_ =0;
410              }
411
412              numberFixed++;
413              if (fixVariables) {
414                if(needBranch >= 2) { // for OsiSimpleInteger: do not branch
415                  // if upper == lower or if info value is outside bounds
416                  // for other objects: branch
417                  OsiBranchingObject * branch = 
418                                        obj->createBranch(solver, info, 0);
419                  branch -> branch (solver);
420                  delete branch;
421                }
422              }
423            }
424
425            ///////////////////////////////////////////////////////////////////
426
427            double downEstimate;
428
429            if (results_[i].downStatus()!=1) {
430              assert (results_[i].downStatus()>=0);
431              downEstimate = results_[i].downChange();
432            }
433            else {
434              // infeasible - just say expensive
435              if (info->cutoff_<1.0e50)
436                downEstimate = 2.0*(info->cutoff_-info->objectiveValue_);
437              else
438                downEstimate = 2.0*fabs(info->objectiveValue_);
439              if (firstForcedObjectIndex_ <0) {
440                firstForcedObjectIndex_ = iObject;
441                firstForcedWhichWay_ =1;
442              }
443              numberFixed++;
444              if (fixVariables) {
445                if(needBranch >= 2) { // for OsiSimpleInteger: do not branch
446                  // if upper == lower or if info value is outside bounds
447                  // for other objects: branch
448                  OsiBranchingObject * branch = 
449                                         obj->createBranch(solver, info, 1);
450                  branch -> branch (solver);
451                  delete branch;
452                }
453              }
454            }
455         
456            double
457              MAXMIN_CRITERION = maxminCrit (info),
458              minVal, maxVal, value;
459
460            if (downEstimate < upEstimate) {
461              minVal = downEstimate;
462              maxVal = upEstimate;
463            } else {
464              minVal = upEstimate;
465              maxVal = downEstimate;
466            }
467
468            if(smallGap) { // use change in objective value
469              value = minVal;
470            }
471            else {
472              value = 
473                estimateProduct_ ? 
474                ((estProdEps + minVal) * maxVal) :
475                (       MAXMIN_CRITERION  * minVal + 
476                        (1.0 - MAXMIN_CRITERION) * maxVal);
477            }
478
479            if((needBranch >= 3) &&
480               saveBestCand(object, iObject, value, upEstimate, downEstimate,
481                            bestTrustedVal1, 
482                            bestTrustedVal2, bestObjectIndex_, bestWhichWay_)) {
483              if(returnCodeSB) { // 1 or 2
484                returnCode = 2; 
485              }
486
487#ifdef USE_SMALL_GAP
488              if(smallGap && (minVal > 1e-3)) {
489                sbObjPosImp = true;
490              }
491#endif
492
493            }
494          }
495        }
496        else { // returnCodeSB is -1 or 3
497          if (returnCodeSB == 3) { // max time - just choose one
498            if(bestObjectIndex_ < 0) {
499              // should not select an integer var with fixed bounds
500              // taken care of below
501              bestObjectIndex_ = list_[0];
502              bestWhichWay_ = 0;
503              returnCode = 0;
504            }
505          }
506          else {
507            returnCode = -1;
508          }
509        }
510#ifdef OLD_STYLE
511        if(bestObjectIndex_ < 0) {
512          bestObjectIndex_ = list_[0];
513        }
514        cardIndForPseudo = 0; // do not scan other vars using pseudo costs
515#endif
516      }
517
518      if((returnCodeSB != -1) && 
519         ((returnCode != 0) || (!sbObjPosImp))) { 
520                  // if returnCodeSB == -1 (i.e. problem is infeasible)
521                  // no need to scan objects with pseudocosts
522                  // if returnCode == 0 and sbObjPOsImp is true
523                  // we want to branch on that object
524                  // to be sure to improve the lower bound
525       
526        // to keep best object skipped on basis of bounds and branching value
527        int bestObjectIndex2 = -1;
528        int bestWhichWay2 = 0;
529        double bestTrusted2Val1 = -COIN_DBL_MAX;
530        double bestTrusted2Val2 = -COIN_DBL_MAX;
531       
532        for(int ips=0; ips<cardIndForPseudo; ips++) {
533          int iObject = indForPseudo[ips];
534          const OsiObject * obj = object[iObject];
535          int needBranch = goodCandidate(solver, info, object, iObject, prec);
536
537          double
538            upEstimate       = (upTotalChange   [iObject] * obj -> upEstimate   ()) / upNumber   [iObject],
539            downEstimate     = (downTotalChange [iObject] * obj -> downEstimate ()) / downNumber [iObject],
540            MAXMIN_CRITERION = maxminCrit (info),
541            minVal, maxVal, value;
542         
543          if (downEstimate < upEstimate) {
544            minVal = downEstimate;
545            maxVal = upEstimate;
546          } else {
547            minVal = upEstimate;
548            maxVal = downEstimate;
549          }
550         
551          value = 
552            estimateProduct_ ? 
553            ((estProdEps + minVal) * maxVal) :
554            (       MAXMIN_CRITERION  * minVal + 
555                    (1.0 - MAXMIN_CRITERION) * maxVal);
556         
557         
558          // store bad candidates in secondary best
559          if(needBranch < 3) {
560            if(saveBestCand(object, iObject, value, 
561                            upEstimate, downEstimate,
562                            bestTrusted2Val1, 
563                            bestTrusted2Val2, bestObjectIndex2, 
564                            bestWhichWay2)) {
565              // no returnCode change
566            }
567
568#ifdef TRACE_STRONG
569            if(problem_->doPrint_) {
570              printf("Object %d skip pseudocost\n", iObject);
571            }
572#endif
573          }
574          else {
575            if(saveBestCand(object, iObject, value, 
576                            upEstimate, downEstimate,
577                            bestTrustedVal1, 
578                            bestTrustedVal2, bestObjectIndex_, bestWhichWay_)) {
579              if(returnCode == 1) { // first saved object
580                                    // is using pseudo-cost
581                returnCode = 0;
582              }
583
584              returnCode = (returnCode ? 3 : 0); // if returnCode was 2 or 3
585                                                 // it becomes 3
586            }
587
588#ifdef TRACE_STRONG
589            if(problem_->doPrint_) {
590              printf("Object %d use pseudocost\n", iObject);
591            }
592#endif
593          }
594        }
595   
596        if((bestObjectIndex_ < 0) && (bestObjectIndex2 >= 0)) {
597          bestObjectIndex_ = bestObjectIndex2;
598          bestWhichWay_ = bestWhichWay2;
599          bestTrustedVal1 = bestTrusted2Val1;
600          if(goodCandidate(solver, info, object, 
601                           bestObjectIndex_, prec) != 4) {
602            returnCode = 4;
603          }
604          else {
605            returnCode = 3;
606          }
607        }
608      }
609
610      int objectVarInd = -1;
611      if(bestObjectIndex_ >= 0) {
612        OsiObject * obj = object[bestObjectIndex_];
613        obj->setWhichWay(bestWhichWay_);
614        // CouenneObject *co =  dynamic_cast <CouenneObject *>(object[bestObjectIndex_]);
615        // if(co) {
616        //   objectVarInd = co->Reference()->Index();
617        // }
618        // else {
619        objectVarInd = obj->columnNumber();
620        //}
621
622        message(BRANCH_VAR)<<bestObjectIndex_<< bestWhichWay_ <<CoinMessageEol;
623
624#ifdef TRACE_STRONG
625        if(problem_->doPrint_) {
626          if(objectVarInd >= 0) {
627            double vUb = solver->getColUpper()[objectVarInd];
628            char strUb[400], strLb[400];
629            if(vUb > 1e50) {
630              sprintf(strUb, "+inf");
631            }
632            else {
633              sprintf(strUb, "%f", vUb);
634            }
635            double vLb = solver->getColLower()[objectVarInd];
636            if(vLb < -1e50) {
637              sprintf(strLb, "-inf");
638            }
639            else {
640              sprintf(strLb, "%f", vLb);
641            }
642            double vSolI = info->solution_[objectVarInd];
643            double vSolS = solver->getColSolution()[objectVarInd];
644            printf("Branch on object %d (var: %d): solInfo: %10.4f  SolSolver: %10.4f low: %s  up: %s\n",
645                   bestObjectIndex_, objectVarInd, vSolI, vSolS, strLb, strU
646
647b);
648          }
649          else {
650            printf("Branch on object %d (var: -1)\n", bestObjectIndex_);
651          }
652        }
653#endif
654      }
655      message(CHOSEN_VAR)<<bestObjectIndex_<<CoinMessageEol;
656   
657      if((numberFixed==numberUnsatisfied_&&numberFixed) &&
658         (goodCandidate(solver, info, object, bestObjectIndex_, prec) != 4)) {
659        returnCode = 4;
660      }
661
662      if((returnCode == 2) || (returnCode == 3)) {
663        if((objectVarInd > -1) && 
664           (goodCandidate(solver, info, object, bestObjectIndex_, prec) != 4)) {
665          // Can occur: two objects for same var, first scanned object
666          // has both branches feasible and is saved as bestObjectIndex_,
667          // second scanned object fixes the variable
668          returnCode = 4;
669        }
670      }
671      retval = returnCode;
672
673      delete[] indForPseudo;
674    }
675    else {
676      retval = 1;
677    }
678
679    // COPY of Bonmin ends here //////////////////////////////////////////////
680
681
682#ifdef TRACE_STRONG2
683    if(problem_->doPrint_) {
684      if(pv > -1) {
685        printf("CCS: end: x[%d]: %10.4f  lb: %10.4f  ub: %10.4f\n",
686               pv, solver->getColSolution()[pv], solver->getColLower()[pv], solver->getColUpper()[pv]);
687        printf("CCS: info: x[%d]: %10.4f  lb: %10.4f  ub: %10.4f\n",
688               pv, info->solution_[pv], info->lower_[pv], info->upper_[pv]);
689        printf("CCS: problem: lb: %10.4f  ub: %10.4f\n",
690               problem_->Lb(pv), problem_->Ub(pv));
691      }
692    }
693#endif
694
695#ifdef TRACE_STRONG
696    if(problem_->doPrint_) {
697      printf("CouenneChooseStrong::ChooseVariable(): retval: %d\n", retval);
698    }
699#endif
700 
701    problem_ -> domain () -> pop ();
702
703    return retval;
704  }
705
706void eliminateIntegerObjects (OsiSolverInterface *model);
707void eliminateIntegerObjects (CbcModel           *model);
708
709/*******************************************************************/
710  // Sets up strong list and clears all if initialize is true.
711  // Returns number of infeasibilities.
712  int CouenneChooseStrong::setupList (OsiBranchingInformation *info, bool initialize) {
713
714    static bool
715      firstCall = true,
716      warned    = false;
717
718    if (firstCall) {
719
720      eliminateIntegerObjects (const_cast <OsiSolverInterface *> (solver_));
721      eliminateIntegerObjects (const_cast <OsiSolverInterface *> (info -> solver_));
722
723      firstCall = false;
724    }
725
726    initialize = true; // to avoid failed assert in BonChooseVariable::setupList()
727
728    problem_ -> domain () -> push
729      (problem_ -> nVars (),
730       info -> solution_, 
731       info -> lower_, 
732       info -> upper_); // have to alloc+copy
733
734    jnlst_ -> Printf (J_ITERSUMMARY, J_BRANCHING, 
735                      "----------------- (strong) setup list\n");
736
737    if (jnlst_ -> ProduceOutput (J_DETAILED, J_BRANCHING)) {
738      for (int i=0; i<problem_ -> domain () -> current () -> Dimension (); i++)
739        printf ("%4d %20.4g [%20.4g %20.4g]\n", i,
740                info -> solution_ [i], info -> lower_ [i], info -> upper_ [i]);
741    }
742
743#ifdef TRACE_STRONG
744    printf("Start CouenneChooseStrong::setupList()\n");
745    OsiObject ** object = info->solver_->objects();
746    if(problem_->doPrint_) {
747      printObjViol(info);
748    }
749#endif
750
751    // int way;
752    // for (int i=0; i<info->solver_->numberObjects(); ++i)
753    //   printf ("[%d:%d,%g] ",
754    //        info -> solver_ -> objects () [i] -> columnNumber (),
755    //        info -> solver_ -> objects () [i] -> priority (),
756    //        info -> solver_ -> objects () [i] -> infeasibility (info,way));
757    // printf ("\n");
758
759    //
760    // Real list setup
761    //
762
763    int retval = gutsOfSetupList (info, initialize);
764
765    if (retval == 0) { // No branching is possible
766
767#ifdef FM_CHECKNLP2
768      if(!(problem_->checkNLP2(info->solution_, 
769                               info->objectiveValue_, true, // care about obj
770                               false, // do not stop at first viol
771                               true, // checkAll
772                               problem_->getFeasTol()))) {
773                                // false for NOT stopping at first violation
774        if (!warned) {
775          printf("CouenneChooseStrong::setupList(): ### WARNING: checkNLP2() returns infeasible, no branching object selected\n");
776          warned = true;
777        }
778      }
779#else /* not FM_CHECKNLP2 */
780      double ckObj = info->objectiveValue_;
781      if(!(problem_->checkNLP(info->solution_, ckObj, true))) {
782        if (!warned) {
783          printf("CouenneChooseStrong::setupList(): ### WARNING: checkNLP() returns infeasible, no branching object selected\n");
784          warned = true;
785        }
786      }
787#endif /* not FM_CHECKNLP2 */
788       
789#ifdef FM_TRACE_OPTSOL
790#ifdef FM_CHECKNLP2
791      problem_->getRecordBestSol()->update();
792#else /* not FM_CHECKNLP2 */
793      problem_->getRecordBestSol()->update(info->solution_, problem_->nVars(),
794                                           ckObj, problem_->getFeasTol());
795#endif /* not FM_CHECKNLP2 */
796#endif
797    }
798
799#ifdef TRACE_STRONG
800    if(problem_->doPrint_) {
801      printf("Strong list: (obj_ind var_ind priority useful viol)\n");
802      printf("numberStrong: %d  numberStrongRoot: %d  retval: %d\n", 
803             numberStrong_, numberStrongRoot_, retval);
804      for(int i=0; i<retval; i++) {
805        // CouenneObject *co =  dynamic_cast <CouenneObject *>(object[list_[i]]);
806        int objectInd = -1;
807        // if(co) {
808        //   objectInd = co->Reference()->Index();
809        // }
810        // else {
811        objectInd = object[list_[i]]->columnNumber();
812        // }
813
814        int wayprint;
815        double violprint = object[list_[i]]->infeasibility(info, wayprint);
816        if(violprint < COIN_DBL_MAX / 100) {
817          printf(" (%d %d %d %6.4f %6.4f)", list_[i], objectInd, 
818                 object[list_[i]]->priority(), useful_[i], 
819                 violprint);
820        }
821        else {
822          printf(" (%d %d %d %6.4f +inf)", list_[i], objectInd, 
823                 object[list_[i]]->priority(), useful_[i]);
824        }
825        printf("\n");
826      }
827#endif
828
829    // for (int i=0; i < (numberStrong_ ? CoinMin (numberStrong_, solver_ -> numberObjects ()) : 1); i++) {
830    //   printf ("list %3d: %3d ", i, list_ [i]);
831    //   if (!((i+1) % 12)) printf ("\n");
832    // }
833
834    jnlst_ -> Printf (J_ITERSUMMARY, J_BRANCHING, 
835                      "----------------- (strong) setup list done - %d infeasibilities\n", retval);
836
837#if (defined TRACE_STRONG)
838    printf("Done CouenneChooseStrong::setupList()\n");
839#endif
840
841    problem_ -> domain () -> pop ();
842    return retval;
843  }
844
845/****************************************************************************/
846  /// Add list of options to be read from file ////////////////////////////////////////
847  void CouenneChooseStrong::registerOptions (Ipopt::SmartPtr <Bonmin::RegisteredOptions> roptions) {
848
849    roptions -> AddStringOption6
850      ("pseudocost_mult",
851       "Multipliers of pseudocosts for estimating and update estimation of bound",
852       "interval_br_rev",
853
854       "infeasibility", "infeasibility returned by object",
855
856       "projectDist",   "distance between current LP point and resulting branches' LP points",
857
858       "interval_lp",   "width of the interval between bound and current lp point",
859       "interval_lp_rev",   "similar to interval_lp, reversed",
860
861       "interval_br",   "width of the interval between bound and branching point",
862       "interval_br_rev",   "similar to interval_br, reversed");
863
864    roptions -> AddStringOption2
865      ("pseudocost_mult_lp",
866       "Use distance between LP points to update multipliers of pseudocosts " 
867       "after simulating branching",
868       "no",
869       "yes", "",
870       "no",  "");
871
872    roptions -> AddStringOption2
873      ("estimate_select",
874       "How the min/max estimates of the subproblems' bounds are used in strong branching",
875       "normal",
876       "normal",   "as usual in literature",
877       "product",  "use their product");
878
879    roptions -> AddStringOption2
880      ("trust_strong",
881       "Fathom strong branching LPs when their bound is above the cutoff",
882       "yes",
883       "yes", "",
884       "no",  "");
885  }
886
887
888  // Returns true if solution looks feasible against given objects
889  bool CouenneChooseStrong::feasibleSolution (const OsiBranchingInformation * info,
890                                              const double * solution,
891                                              int numberObjects,
892                                              const OsiObject ** objects) {
893
894#ifdef FM_CHECKNLP2
895    return problem_ -> checkNLP2 (solution, 0, false, true, true, 
896                                  problem_->getFeasTol());
897#else
898    int indobj = problem_ -> Obj (0) -> Body () -> Index ();
899    return problem_ -> checkNLP (solution, indobj >= 0 ? solution [indobj] : problem_ -> Obj (0) -> Body () -> Value ());
900#endif
901  }
902
903/****************************************************************************/
904  void CouenneChooseStrong::printObjViol(OsiBranchingInformation *info) {
905
906    OsiObject ** object = info->solver_->objects();
907    int numberObjects = info->solver_->numberObjects();
908
909    printf("CouenneChooseStrong::printObjViol(): Object violations: (obj_ind  var_ind  violation)");
910    double maxViol = 0;
911    double minPosViol = 1e50;
912    for(int i=0; i<numberObjects; i++) {
913      //CouenneObject *co =  dynamic_cast <CouenneObject *>(object[i]);
914      int indVar = -1;
915      // if(co) {
916      //        indVar = co->Reference()->Index();
917      // }
918      // else {
919      indVar = object[i]->columnNumber();
920      // }
921      int way;
922      double value = object[i]->infeasibility(info,way);
923      //double value = object[i]->checkInfeasibility(info);
924      maxViol = (value > maxViol ? value : maxViol);
925      if(value > 0.0) {
926        printf("(%d %d %f)", i, indVar, value);
927        minPosViol = (value < minPosViol ? value : minPosViol);
928      }
929    }
930    printf("\nmaxViol: %g  minPosViol: %g\n", maxViol, minPosViol);
931  } /* printObjViol */
932
933//}
Note: See TracBrowser for help on using the repository browser.