1 /* $Id$ */
2 // Copyright (C) 2002, International Business Machines
3 // Corporation and others.  All Rights Reserved.
4 // This code is licensed under the terms of the Eclipse Public License (EPL).
5 
6 #if defined(_MSC_VER)
7 // Turn off compiler warning about long names
8 #pragma warning(disable : 4786)
9 #endif
10 #include <cassert>
11 #include <cstdlib>
12 #include <cmath>
13 #include <cfloat>
14 //#define CBC_DEBUG
15 //#define TRACE_ONE 19
16 #include "OsiSolverInterface.hpp"
17 #include "OsiSolverBranch.hpp"
18 #include "CbcModel.hpp"
19 #include "CbcMessage.hpp"
20 #include "CbcBranchDynamic.hpp"
21 #include "CoinSort.hpp"
22 #include "CoinError.hpp"
23 
24 // Removing magic constants.
25 
26 // This is a very small number, added to something to make sure it's non-zero.
27 // Useful, for example in denominators of ratios to avoid any possible division by zero
28 #define nonZeroAmount 1.0e-30
29 
30 // Increasing the size of an array when it grows to the end of its alloted space.
31 // In this file, only used for the history of the outcome of a branch.
32 // New size is   size_scale_numerator* <old value> / size_scale_denominator + additive_size_increase.
33 
34 #define size_scale_numerator 3
35 #define size_scale_denominator 2
36 #define additive_size_increase 100
37 
38 // Explanation of options used in this file
39 
40 // TYPE2 defines a strategy for computing pseudocosts
41 // 0 means to just use the absolute change in objective
42 // 1 means use the relative change in objective
43 // 2 means use a convex combination of absolute and relative objective changes
44 
45 // For option 2 (TYPE2 == 2), the specific combination is controlled by TYPERATIO
46 // Includes a TYPERATIO fraction of the absolute change and (1 - TYPERATIO) fraction of
47 // the relative change.  So should in general have 0 <= TYPERATIO <= 1.  But for the
48 // equality cases, you're better off using the other strategy (TYPE2) options.
49 
50 #ifdef COIN_DEVELOP
51 typedef struct {
52   double sumUp_;
53   double upEst_; // or change in obj in update
54   double sumDown_;
55   double downEst_; // or movement in value in update
56   int sequence_;
57   int numberUp_;
58   int numberUpInf_;
59   int numberDown_;
60   int numberDownInf_;
61   char where_;
62   char status_;
63 } History;
64 static History *history = NULL;
65 static int numberHistory = 0;
66 static int maxHistory = 0;
67 static bool getHistoryStatistics_ = true;
increaseHistory()68 static void increaseHistory()
69 {
70   if (numberHistory == maxHistory) {
71     // This was originally 3 * maxHistory/2 + 100
72     maxHistory = additive_size_increase + (size_scale_numerator * maxHistory) / size_scale_denominator;
73     History *temp = new History[maxHistory];
74     memcpy(temp, history, numberHistory * sizeof(History));
75     delete[] history;
76     history = temp;
77   }
78 }
addRecord(History newOne)79 static bool addRecord(History newOne)
80 {
81   //if (!getHistoryStatistics_)
82   return false;
83   bool fromCompare = false;
84   int i;
85   for (i = numberHistory - 1; i >= 0; i--) {
86     if (newOne.sequence_ != history[i].sequence_)
87       continue;
88     if (newOne.where_ != history[i].where_)
89       continue;
90     if (newOne.numberUp_ != history[i].numberUp_)
91       continue;
92     if (newOne.sumUp_ != history[i].sumUp_)
93       continue;
94     if (newOne.numberUpInf_ != history[i].numberUpInf_)
95       continue;
96     if (newOne.upEst_ != history[i].upEst_)
97       continue;
98     if (newOne.numberDown_ != history[i].numberDown_)
99       continue;
100     if (newOne.sumDown_ != history[i].sumDown_)
101       continue;
102     if (newOne.numberDownInf_ != history[i].numberDownInf_)
103       continue;
104     if (newOne.downEst_ != history[i].downEst_)
105       continue;
106     // If B knock out previous B
107     if (newOne.where_ == 'C') {
108       fromCompare = true;
109       if (newOne.status_ == 'B') {
110         int j;
111         for (j = i - 1; j >= 0; j--) {
112           if (history[j].where_ == 'C') {
113             if (history[j].status_ == 'I') {
114               break;
115             } else if (history[j].status_ == 'B') {
116               history[j].status_ = ' ';
117               break;
118             }
119           }
120         }
121       }
122       break;
123     }
124   }
125   if (i == -1 || fromCompare) {
126     //add
127     increaseHistory();
128     history[numberHistory++] = newOne;
129     return true;
130   } else {
131     return false;
132   }
133 }
134 #endif
135 
136 // Default Constructor
CbcBranchDynamicDecision()137 CbcBranchDynamicDecision::CbcBranchDynamicDecision()
138   : CbcBranchDecision()
139 {
140   bestCriterion_ = 0.0;
141   bestChangeUp_ = 0.0;
142   bestNumberUp_ = 0;
143   bestChangeDown_ = 0.0;
144   bestNumberDown_ = 0;
145   bestObject_ = NULL;
146 }
147 
148 // Copy constructor
CbcBranchDynamicDecision(const CbcBranchDynamicDecision & rhs)149 CbcBranchDynamicDecision::CbcBranchDynamicDecision(
150   const CbcBranchDynamicDecision &rhs)
151   : CbcBranchDecision()
152 {
153   bestCriterion_ = rhs.bestCriterion_;
154   bestChangeUp_ = rhs.bestChangeUp_;
155   bestNumberUp_ = rhs.bestNumberUp_;
156   bestChangeDown_ = rhs.bestChangeDown_;
157   bestNumberDown_ = rhs.bestNumberDown_;
158   bestObject_ = rhs.bestObject_;
159 }
160 
~CbcBranchDynamicDecision()161 CbcBranchDynamicDecision::~CbcBranchDynamicDecision()
162 {
163 }
164 
165 // Clone
166 CbcBranchDecision *
clone() const167 CbcBranchDynamicDecision::clone() const
168 {
169   return new CbcBranchDynamicDecision(*this);
170 }
171 
172 // Initialize i.e. before start of choosing at a node
initialize(CbcModel *)173 void CbcBranchDynamicDecision::initialize(CbcModel * /*model*/)
174 {
175   bestCriterion_ = 0.0;
176   bestChangeUp_ = 0.0;
177   bestNumberUp_ = 0;
178   bestChangeDown_ = 0.0;
179   bestNumberDown_ = 0;
180   bestObject_ = NULL;
181 #ifdef COIN_DEVELOP
182   History hist;
183   hist.where_ = 'C';
184   hist.status_ = 'I';
185   hist.sequence_ = 55555;
186   hist.numberUp_ = 0;
187   hist.numberUpInf_ = 0;
188   hist.sumUp_ = 0.0;
189   hist.upEst_ = 0.0;
190   hist.numberDown_ = 0;
191   hist.numberDownInf_ = 0;
192   hist.sumDown_ = 0.0;
193   hist.downEst_ = 0.0;
194   addRecord(hist);
195 #endif
196 }
197 
198 /* Saves a clone of current branching object.  Can be used to update
199       information on object causing branch - after branch */
saveBranchingObject(OsiBranchingObject * object)200 void CbcBranchDynamicDecision::saveBranchingObject(OsiBranchingObject *object)
201 {
202   OsiBranchingObject *obj = object->clone();
203 #ifndef NDEBUG
204   CbcBranchingObject *obj2 = dynamic_cast< CbcBranchingObject * >(obj);
205   assert(obj2);
206 #if COIN_DEVELOP > 1
207   CbcDynamicPseudoCostBranchingObject *branchingObject = dynamic_cast< CbcDynamicPseudoCostBranchingObject * >(obj);
208   if (!branchingObject)
209     printf("no dynamic branching object Dynamic Decision\n");
210 #endif
211 #else
212   CbcBranchingObject *obj2 = static_cast< CbcBranchingObject * >(obj);
213 #endif
214   //object_=branchingObject;
215   object_ = obj2;
216 }
217 /* Pass in information on branch just done.
218    assumes object can get information from solver */
219 /*
220   The expectation is that this method will be called after the branch has been
221   imposed on the constraint system and resolve() has executed.
222 
223   Note that the CbcBranchDecision is a property of the CbcModel. Note also that
224   this method is reaching right through the CbcBranchingObject to update
225   information in the underlying CbcObject. That's why we delete the
226   branchingObject at the end of the method --- the next time we're called,
227   the CbcObject will be different.
228 */
updateInformation(OsiSolverInterface * solver,const CbcNode * node)229 void CbcBranchDynamicDecision::updateInformation(OsiSolverInterface *solver,
230   const CbcNode *node)
231 {
232   assert(object_);
233   const CbcModel *model = object_->model();
234   double originalValue = node->objectiveValue();
235   int originalUnsatisfied = node->numberUnsatisfied();
236   double objectiveValue = solver->getObjValue() * model->getObjSense();
237   int unsatisfied = 0;
238   int i;
239   int numberIntegers = model->numberIntegers();
240   ;
241   const double *solution = solver->getColSolution();
242   //const double * lower = solver->getColLower();
243   //const double * upper = solver->getColUpper();
244   /*
245 	 Gain access to the associated CbcBranchingObject and its underlying
246 	 CbcObject.
247 
248 	 Seems like we'd want to distinguish between no branching object and a
249 	 branching object of the wrong type. Just deleting an object of the wrong
250 	 type hides many sins.
251 
252 	 Hmmm ... if we're using the OSI side of the hierarchy, is this indicated by a
253 	 null object_? Nah, then we have an assert failure off the top.
254 	*/
255 
256   CbcDynamicPseudoCostBranchingObject *branchingObject = dynamic_cast< CbcDynamicPseudoCostBranchingObject * >(object_);
257   if (!branchingObject) {
258     delete object_;
259     object_ = NULL;
260     return;
261   }
262   CbcSimpleIntegerDynamicPseudoCost *object = branchingObject->object();
263   /*
264 	change is the change in objective due to the branch we've just imposed. It's
265 	possible we may have gone infeasible.
266 	*/
267   double change = CoinMax(0.0, objectiveValue - originalValue);
268   // probably should also ignore if stopped
269   // FIXME. Could use enum to avoid numbers for iStatus (e.g. optimal, unknown, infeasible)
270   int iStatus;
271   if (solver->isProvenOptimal())
272     iStatus = 0; // optimal
273   else if (solver->isIterationLimitReached()
274     && !solver->isDualObjectiveLimitReached())
275     iStatus = 2; // unknown
276   else
277     iStatus = 1; // infeasible
278   /*
279 	  If we're feasible according to the solver, evaluate integer feasibility.
280 	*/
281   bool feasible = iStatus != 1;
282   if (feasible) {
283     double integerTolerance = model->getDblParam(CbcModel::CbcIntegerTolerance);
284     const int *integerVariable = model->integerVariable();
285     for (i = 0; i < numberIntegers; i++) {
286       int j = integerVariable[i];
287       double value = solution[j];
288       double nearest = floor(value + 0.5);
289       if (fabs(value - nearest) > integerTolerance)
290         unsatisfied++;
291     }
292   }
293   /*
294 	  Finally, update the object. Defaults (080104) are TYPE2 = 0, INFEAS = 1.
295 
296 	  Pseudocosts are at heart the average of actual costs for a branch. We just
297 	  need to update the information used to calculate that average.
298 	*/
299   int way = object_->way();
300   double value = object_->value();
301   //#define TYPE2 1
302   //#define TYPERATIO 0.9
303   if (way < 0) {
304     // down
305     if (feasible) {
306       double movement = value - floor(value);
307       movement = CoinMax(movement, MINIMUM_MOVEMENT);
308       //printf("(down change %g value down %g ",change,movement);
309       object->incrementNumberTimesDown();
310       object->addToSumDownChange(nonZeroAmount + movement);
311       object->addToSumDownDecrease(originalUnsatisfied - unsatisfied);
312 #if TYPE2 == 0
313       object->addToSumDownCost(change / (nonZeroAmount + movement));
314       object->setDownDynamicPseudoCost(object->sumDownCost() / static_cast< double >(object->numberTimesDown()));
315 #elif TYPE2 == 1
316       object->addToSumDownCost(change);
317       object->setDownDynamicPseudoCost(object->sumDownCost() / object->sumDownChange());
318 #elif TYPE2 == 2
319       object->addToSumDownCost(change * TYPERATIO + (1.0 - TYPERATIO) * change / (nonZeroAmount + movement));
320       object->setDownDynamicPseudoCost(object->sumDownCost() * (TYPERATIO / object->sumDownChange() + (1.0 - TYPERATIO) / (double)object->numberTimesDown()));
321 #endif
322     } else {
323       //printf("(down infeasible value down %g ",change,movement);
324       object->incrementNumberTimesDown();
325       object->incrementNumberTimesDownInfeasible();
326 #if INFEAS == 2
327       double distanceToCutoff = 0.0;
328       double objectiveValue = model->getCurrentMinimizationObjValue();
329       distanceToCutoff = model->getCutoff() - originalValue;
330       if (distanceToCutoff < 1.0e20)
331         change = distanceToCutoff * 2.0;
332       else
333         change = object->downDynamicPseudoCost() * movement * 10.0;
334       change = CoinMax(1.0e-12 * (1.0 + fabs(originalValue)), change);
335       object->addToSumDownChange(nonZeroAmount + movement);
336       object->addToSumDownDecrease(originalUnsatisfied - unsatisfied);
337 #if TYPE2 == 0
338       object->addToSumDownCost(change / (nonZeroAmount + movement));
339       object->setDownDynamicPseudoCost(object->sumDownCost() / (double)object->numberTimesDown());
340 #elif TYPE2 == 1
341       object->addToSumDownCost(change);
342       object->setDownDynamicPseudoCost(object->sumDownCost() / object->sumDownChange());
343 #elif TYPE2 == 2
344       object->addToSumDownCost(change * TYPERATIO + (1.0 - TYPERATIO) * change / (nonZeroAmount + movement));
345       object->setDownDynamicPseudoCost(object->sumDownCost() * (TYPERATIO / object->sumDownChange() + (1.0 - TYPERATIO) / (double)object->numberTimesDown()));
346 #endif
347 #endif
348     }
349   } else {
350     // up
351     if (feasible) {
352       double movement = ceil(value) - value;
353       movement = CoinMax(movement, MINIMUM_MOVEMENT);
354       //printf("(up change %g value down %g ",change,movement);
355       object->incrementNumberTimesUp();
356       object->addToSumUpChange(nonZeroAmount + movement);
357       object->addToSumUpDecrease(unsatisfied - originalUnsatisfied);
358 #if TYPE2 == 0
359       object->addToSumUpCost(change / (nonZeroAmount + movement));
360       object->setUpDynamicPseudoCost(object->sumUpCost() / static_cast< double >(object->numberTimesUp()));
361 #elif TYPE2 == 1
362       object->addToSumUpCost(change);
363       object->setUpDynamicPseudoCost(object->sumUpCost() / object->sumUpChange());
364 #elif TYPE2 == 2
365       object->addToSumUpCost(change * TYPERATIO + (1.0 - TYPERATIO) * change / (nonZeroAmount + movement));
366       object->setUpDynamicPseudoCost(object->sumUpCost() * (TYPERATIO / object->sumUpChange() + (1.0 - TYPERATIO) / (double)object->numberTimesUp()));
367 #endif
368     } else {
369       //printf("(up infeasible value down %g ",change,movement);
370       object->incrementNumberTimesUp();
371       object->incrementNumberTimesUpInfeasible();
372 #if INFEAS == 2
373       double distanceToCutoff = 0.0;
374       double objectiveValue = model->getCurrentMinimizationObjValue();
375       distanceToCutoff = model->getCutoff() - originalValue;
376       if (distanceToCutoff < 1.0e20)
377         change = distanceToCutoff * 2.0;
378       else
379         change = object->upDynamicPseudoCost() * movement * 10.0;
380       change = CoinMax(1.0e-12 * (1.0 + fabs(originalValue)), change);
381       object->addToSumUpChange(nonZeroAmount + movement);
382       object->addToSumUpDecrease(unsatisfied - originalUnsatisfied);
383 #if TYPE2 == 0
384       object->addToSumUpCost(change / (nonZeroAmount + movement));
385       object->setUpDynamicPseudoCost(object->sumUpCost() / (double)object->numberTimesUp());
386 #elif TYPE2 == 1
387       object->addToSumUpCost(change);
388       object->setUpDynamicPseudoCost(object->sumUpCost() / object->sumUpChange());
389 #elif TYPE2 == 2
390       object->addToSumUpCost(change * TYPERATIO + (1.0 - TYPERATIO) * change / (nonZeroAmount + movement));
391       object->setUpDynamicPseudoCost(object->sumUpCost() * (TYPERATIO / object->sumUpChange() + (1.0 - TYPERATIO) / (double)object->numberTimesUp()));
392 #endif
393 #endif
394     }
395   }
396   //object->print(1,0.5);
397   delete object_;
398   object_ = NULL;
399 }
400 
401 /*
402   Simple dynamic decision algorithm. Compare based on infeasibility (numInfUp,
403   numInfDown) until a solution is found by search, then switch to change in
404   objective (changeUp, changeDown). Note that bestSoFar is remembered in
405   bestObject_, so the parameter bestSoFar is unused.
406 */
407 
betterBranch(CbcBranchingObject * thisOne,CbcBranchingObject *,double changeUp,int numInfUp,double changeDown,int numInfDown)408 int CbcBranchDynamicDecision::betterBranch(CbcBranchingObject *thisOne,
409   CbcBranchingObject * /*bestSoFar*/,
410   double changeUp, int numInfUp,
411   double changeDown, int numInfDown)
412 {
413   CbcModel *model = thisOne->model();
414   int stateOfSearch = model->stateOfSearch() % 10;
415   int betterWay = 0;
416   double value = 0.0;
417   if (!bestObject_) {
418     bestCriterion_ = -1.0e30;
419     bestNumberUp_ = COIN_INT_MAX;
420     bestNumberDown_ = COIN_INT_MAX;
421   }
422   // maybe branch up more if no solution or not many nodes done?
423   if (stateOfSearch <= 2) {
424     //#define TRY_STUFF 1
425 #ifdef TRY_STUFF
426     // before solution - choose smallest number
427     // could add in depth as well
428     int bestNumber = CoinMin(bestNumberUp_, bestNumberDown_);
429     if (numInfUp < numInfDown) {
430       if (numInfUp < bestNumber) {
431         betterWay = 1;
432       } else if (numInfUp == bestNumber) {
433         if (changeUp < bestChangeUp_)
434           betterWay = 1;
435       }
436     } else if (numInfUp > numInfDown) {
437       if (numInfDown < bestNumber) {
438         betterWay = -1;
439       } else if (numInfDown == bestNumber) {
440         if (changeDown < bestChangeDown_)
441           betterWay = -1;
442       }
443     } else {
444       // up and down have same number
445       bool better = false;
446       if (numInfUp < bestNumber) {
447         better = true;
448       } else if (numInfUp == bestNumber) {
449         if (CoinMin(changeUp, changeDown) < CoinMin(bestChangeUp_, bestChangeDown_) - 1.0e-5)
450           better = true;
451         ;
452       }
453       if (better) {
454         // see which way
455         if (changeUp <= changeDown)
456           betterWay = 1;
457         else
458           betterWay = -1;
459       }
460     }
461     if (betterWay) {
462       value = CoinMin(numInfUp, numInfDown);
463     }
464 #else
465     // use pseudo shadow prices modified by locks
466     // test testosi
467 #ifndef JJF_ONE
468     double objectiveValue = model->getCurrentMinimizationObjValue();
469     double distanceToCutoff = model->getCutoff() - objectiveValue;
470     if (distanceToCutoff < 1.0e20)
471       distanceToCutoff *= 10.0;
472     else
473       distanceToCutoff = 1.0e2 + fabs(objectiveValue);
474     distanceToCutoff = CoinMax(distanceToCutoff, 1.0e-12 * (1.0 + fabs(objectiveValue)));
475     double continuousObjective = model->getContinuousObjective();
476     double distanceToCutoffC = model->getCutoff() - continuousObjective;
477     if (distanceToCutoffC > 1.0e20)
478       distanceToCutoffC = 1.0e2 + fabs(objectiveValue);
479     distanceToCutoffC = CoinMax(distanceToCutoffC, 1.0e-12 * (1.0 + fabs(objectiveValue)));
480     int numberInfC = model->getContinuousInfeasibilities();
481     double perInf = distanceToCutoffC / static_cast< double >(numberInfC);
482     assert(perInf > 0.0);
483     //int numberIntegers = model->numberIntegers();
484     changeDown += perInf * numInfDown;
485     changeUp += perInf * numInfUp;
486 #ifdef JJF_ZERO
487     if (numInfDown == 1) {
488       if (numInfUp == 1) {
489         changeUp += 1.0e6;
490         changeDown += 1.0e6;
491       } else if (changeDown <= 1.5 * changeUp) {
492         changeUp += 1.0e6;
493       }
494     } else if (numInfUp == 1 && changeUp <= 1.5 * changeDown) {
495       changeDown += 1.0e6;
496     }
497 #endif
498 #endif
499     double minValue = CoinMin(changeDown, changeUp);
500     double maxValue = CoinMax(changeDown, changeUp);
501     value = WEIGHT_BEFORE * minValue + (1.0 - WEIGHT_BEFORE) * maxValue;
502     if (value > bestCriterion_ + 1.0e-8) {
503       if (changeUp <= 1.5 * changeDown) {
504         betterWay = 1;
505       } else {
506         betterWay = -1;
507       }
508     }
509 #endif
510   } else {
511 #define TRY_STUFF 2
512 #if TRY_STUFF > 1
513     // Get current number of infeasibilities, cutoff and current objective
514     CbcNode *node = model->currentNode();
515     int numberUnsatisfied = node->numberUnsatisfied();
516     double cutoff = model->getCutoff();
517     double objectiveValue = node->objectiveValue();
518 #endif
519     // got a solution
520     double minValue = CoinMin(changeDown, changeUp);
521     double maxValue = CoinMax(changeDown, changeUp);
522     // Reduce
523 #ifdef TRY_STUFF
524     //maxValue = CoinMin(maxValue,minValue*4.0);
525 #else
526     //maxValue = CoinMin(maxValue,minValue*2.0);
527 #endif
528 #ifndef WEIGHT_PRODUCT
529     value = WEIGHT_AFTER * minValue + (1.0 - WEIGHT_AFTER) * maxValue;
530 #else
531     double minProductWeight = model->getDblParam(CbcModel::CbcSmallChange);
532     value = CoinMax(minValue, minProductWeight) * CoinMax(maxValue, minProductWeight);
533     //value += minProductWeight*minValue;
534 #endif
535     double useValue = value;
536     double useBest = bestCriterion_;
537 #if TRY_STUFF > 1
538     int thisNumber = CoinMin(numInfUp, numInfDown);
539     int bestNumber = CoinMin(bestNumberUp_, bestNumberDown_);
540     double distance = cutoff - objectiveValue;
541     assert(distance >= 0.0);
542     if (useValue + 0.1 * distance > useBest && useValue * 1.1 > useBest && useBest + 0.1 * distance > useValue && useBest * 1.1 > useValue) {
543       // not much in it - look at unsatisfied
544       if (thisNumber < numberUnsatisfied || bestNumber < numberUnsatisfied) {
545         double perInteger = distance / (static_cast< double >(numberUnsatisfied));
546         useValue += thisNumber * perInteger;
547         useBest += bestNumber * perInteger;
548       }
549     }
550 #endif
551     if (useValue > useBest + 1.0e-8) {
552       if (changeUp <= 1.5 * changeDown) {
553         betterWay = 1;
554       } else {
555         betterWay = -1;
556       }
557     }
558   }
559 #ifdef COIN_DEVELOP
560   History hist;
561   {
562     CbcDynamicPseudoCostBranchingObject *branchingObject = dynamic_cast< CbcDynamicPseudoCostBranchingObject * >(thisOne);
563     if (branchingObject) {
564       CbcSimpleIntegerDynamicPseudoCost *object = branchingObject->object();
565       assert(object);
566       hist.where_ = 'C';
567       hist.status_ = ' ';
568       hist.sequence_ = object->columnNumber();
569       hist.numberUp_ = object->numberTimesUp();
570       hist.numberUpInf_ = numInfUp;
571       hist.sumUp_ = object->sumUpCost();
572       hist.upEst_ = changeUp;
573       hist.numberDown_ = object->numberTimesDown();
574       hist.numberDownInf_ = numInfDown;
575       hist.sumDown_ = object->sumDownCost();
576       hist.downEst_ = changeDown;
577     }
578   }
579 #endif
580   if (betterWay) {
581 #ifdef COIN_DEVELOP
582     hist.status_ = 'B';
583 #endif
584     // maybe change better way
585     CbcDynamicPseudoCostBranchingObject *branchingObject = dynamic_cast< CbcDynamicPseudoCostBranchingObject * >(thisOne);
586     if (branchingObject) {
587       CbcSimpleIntegerDynamicPseudoCost *object = branchingObject->object();
588       double separator = object->upDownSeparator();
589       if (separator > 0.0) {
590         const double *solution = thisOne->model()->testSolution();
591         double valueVariable = solution[object->columnNumber()];
592         betterWay = (valueVariable - floor(valueVariable) >= separator) ? 1 : -1;
593       }
594     }
595     bestCriterion_ = value;
596     bestChangeUp_ = changeUp;
597     bestNumberUp_ = numInfUp;
598     bestChangeDown_ = changeDown;
599     bestNumberDown_ = numInfDown;
600     bestObject_ = thisOne;
601     // See if user is overriding way
602     if (thisOne->object() && thisOne->object()->preferredWay())
603       betterWay = thisOne->object()->preferredWay();
604   }
605 #ifdef COIN_DEVELOP
606   addRecord(hist);
607 #endif
608   return betterWay;
609 }
610 /* Sets or gets best criterion so far */
setBestCriterion(double value)611 void CbcBranchDynamicDecision::setBestCriterion(double value)
612 {
613   bestCriterion_ = value;
614 }
615 double
getBestCriterion() const616 CbcBranchDynamicDecision::getBestCriterion() const
617 {
618   return bestCriterion_;
619 }
620 #ifdef COIN_DEVELOP
printHistory(const char * file)621 void printHistory(const char *file)
622 {
623   if (!numberHistory)
624     return;
625   FILE *fp = fopen(file, "w");
626   assert(fp);
627   int numberIntegers = 0;
628   int i;
629   for (i = 0; i < numberHistory; i++) {
630     if (history[i].where_ != 'C' || history[i].status_ != 'I')
631       numberIntegers = CoinMax(numberIntegers, history[i].sequence_);
632   }
633   numberIntegers++;
634   for (int iC = 0; iC < numberIntegers; iC++) {
635     int n = 0;
636     for (i = 0; i < numberHistory; i++) {
637       if (history[i].sequence_ == iC) {
638         if (!n)
639           fprintf(fp, "XXX %d\n", iC);
640         n++;
641         fprintf(fp, "%c%c up %8d %8d %12.5f %12.5f down  %8d %8d %12.5f %12.5f\n",
642           history[i].where_,
643           history[i].status_,
644           history[i].numberUp_,
645           history[i].numberUpInf_,
646           history[i].sumUp_,
647           history[i].upEst_,
648           history[i].numberDown_,
649           history[i].numberDownInf_,
650           history[i].sumDown_,
651           history[i].downEst_);
652       }
653     }
654   }
655   fclose(fp);
656 }
657 #endif
658 
659 // Default Constructor
CbcDynamicPseudoCostBranchingObject()660 CbcDynamicPseudoCostBranchingObject::CbcDynamicPseudoCostBranchingObject()
661   : CbcIntegerBranchingObject()
662 {
663   changeInGuessed_ = 1.0e-5;
664   object_ = NULL;
665 }
666 
667 // Useful constructor
CbcDynamicPseudoCostBranchingObject(CbcModel * model,int variable,int way,double value,CbcSimpleIntegerDynamicPseudoCost * object)668 CbcDynamicPseudoCostBranchingObject::CbcDynamicPseudoCostBranchingObject(CbcModel *model,
669   int variable,
670   int way, double value,
671   CbcSimpleIntegerDynamicPseudoCost *object)
672   : CbcIntegerBranchingObject(model, variable, way, value)
673 {
674   changeInGuessed_ = 1.0e-5;
675   object_ = object;
676 }
677 // Does part of work for constructor
fillPart(int variable,int way,double value,CbcSimpleIntegerDynamicPseudoCost * object)678 void CbcDynamicPseudoCostBranchingObject::fillPart(int variable,
679   int way, double value,
680   CbcSimpleIntegerDynamicPseudoCost *object)
681 {
682   CbcIntegerBranchingObject::fillPart(variable, way, value);
683   changeInGuessed_ = 1.0e-5;
684   object_ = object;
685 }
686 // Useful constructor for fixing
CbcDynamicPseudoCostBranchingObject(CbcModel * model,int variable,int way,double lowerValue,double)687 CbcDynamicPseudoCostBranchingObject::CbcDynamicPseudoCostBranchingObject(CbcModel *model,
688   int variable, int way,
689   double lowerValue,
690   double /*upperValue*/)
691   : CbcIntegerBranchingObject(model, variable, way, lowerValue)
692 {
693   changeInGuessed_ = 1.0e100;
694   object_ = NULL;
695 }
696 
697 // Copy constructor
CbcDynamicPseudoCostBranchingObject(const CbcDynamicPseudoCostBranchingObject & rhs)698 CbcDynamicPseudoCostBranchingObject::CbcDynamicPseudoCostBranchingObject(
699   const CbcDynamicPseudoCostBranchingObject &rhs)
700   : CbcIntegerBranchingObject(rhs)
701 {
702   changeInGuessed_ = rhs.changeInGuessed_;
703   object_ = rhs.object_;
704 }
705 
706 // Assignment operator
707 CbcDynamicPseudoCostBranchingObject &
operator =(const CbcDynamicPseudoCostBranchingObject & rhs)708 CbcDynamicPseudoCostBranchingObject::operator=(const CbcDynamicPseudoCostBranchingObject &rhs)
709 {
710   if (this != &rhs) {
711     CbcIntegerBranchingObject::operator=(rhs);
712     changeInGuessed_ = rhs.changeInGuessed_;
713     object_ = rhs.object_;
714   }
715   return *this;
716 }
717 CbcBranchingObject *
clone() const718 CbcDynamicPseudoCostBranchingObject::clone() const
719 {
720   return (new CbcDynamicPseudoCostBranchingObject(*this));
721 }
722 
723 // Destructor
~CbcDynamicPseudoCostBranchingObject()724 CbcDynamicPseudoCostBranchingObject::~CbcDynamicPseudoCostBranchingObject()
725 {
726 }
727 
728 /*
729   Perform a branch by adjusting the bounds of the specified variable. Note
730   that each arm of the branch advances the object to the next arm by
731   advancing the value of way_.
732 
733   Providing new values for the variable's lower and upper bounds for each
734   branching direction gives a little bit of additional flexibility and will
735   be easily extensible to multi-way branching.
736   Returns change in guessed objective on next branch
737 */
738 double
branch()739 CbcDynamicPseudoCostBranchingObject::branch()
740 {
741   CbcIntegerBranchingObject::branch();
742   return changeInGuessed_;
743 }
744 /* Some branchingObjects may claim to be able to skip
745    strong branching.  If so they have to fill in CbcStrongInfo.
746    The object mention in incoming CbcStrongInfo must match.
747    Returns nonzero if skip is wanted */
fillStrongInfo(CbcStrongInfo & info)748 int CbcDynamicPseudoCostBranchingObject::fillStrongInfo(CbcStrongInfo &info)
749 {
750   assert(object_);
751   assert(info.possibleBranch == this);
752   info.upMovement = object_->upDynamicPseudoCost() * (ceil(value_) - value_);
753   info.downMovement = object_->downDynamicPseudoCost() * (value_ - floor(value_));
754   info.numIntInfeasUp -= static_cast< int >(object_->sumUpDecrease() / (1.0e-12 + static_cast< double >(object_->numberTimesUp())));
755   info.numIntInfeasUp = CoinMax(info.numIntInfeasUp, 0);
756   info.numObjInfeasUp = 0;
757   info.finishedUp = false;
758   info.numItersUp = 0;
759   info.numIntInfeasDown -= static_cast< int >(object_->sumDownDecrease() / (1.0e-12 + static_cast< double >(object_->numberTimesDown())));
760   info.numIntInfeasDown = CoinMax(info.numIntInfeasDown, 0);
761   info.numObjInfeasDown = 0;
762   info.finishedDown = false;
763   info.numItersDown = 0;
764   info.fix = 0;
765   if (object_->numberTimesUp() < object_->numberBeforeTrust() + 2 * object_->numberTimesUpInfeasible() || object_->numberTimesDown() < object_->numberBeforeTrust() + 2 * object_->numberTimesDownInfeasible()) {
766     return 0;
767   } else {
768     return 1;
769   }
770 }
771 
772 /* vi: softtabstop=2 shiftwidth=2 expandtab tabstop=2
773 */
774