source: trunk/CbcBranchDynamic.cpp @ 137

Last change on this file since 137 was 137, checked in by forrest, 16 years ago

fiddling with strong branching

  • Property svn:eol-style set to native
  • Property svn:keywords set to Author Date Id Revision
File size: 23.0 KB
Line 
1// Copyright (C) 2002, International Business Machines
2// Corporation and others.  All Rights Reserved.
3#if defined(_MSC_VER)
4// Turn off compiler warning about long names
5#  pragma warning(disable:4786)
6#endif
7#include <cassert>
8#include <cmath>
9#include <cfloat>
10//#define CBC_DEBUG
11
12#include "OsiSolverInterface.hpp"
13#include "CbcModel.hpp"
14#include "CbcMessage.hpp"
15#include "CbcBranchDynamic.hpp"
16#include "CoinSort.hpp"
17#include "CoinError.hpp"
18
19/** Default Constructor
20
21  Equivalent to an unspecified binary variable.
22*/
23CbcSimpleIntegerDynamicPseudoCost::CbcSimpleIntegerDynamicPseudoCost ()
24  : CbcSimpleInteger(),
25    downDynamicPseudoCost_(1.0e-5),
26    upDynamicPseudoCost_(1.0e-5),
27    sumDownCost_(0.0),
28    sumUpCost_(0.0),
29    sumDownCostSquared_(0.0),
30    sumUpCostSquared_(0.0),
31    sumDownDecrease_(0.0),
32    sumUpDecrease_(0.0),
33    lastDownCost_(0.0),
34    lastUpCost_(0.0),
35    lastDownDecrease_(0),
36    lastUpDecrease_(0),
37    numberTimesDown_(0),
38    numberTimesUp_(0),
39    numberTimesDownInfeasible_(0),
40    numberTimesUpInfeasible_(0),
41    numberBeforeTrust_(0),
42    method_(0)
43{
44}
45
46/** Useful constructor
47
48  Loads dynamic upper & lower bounds for the specified variable.
49*/
50CbcSimpleIntegerDynamicPseudoCost::CbcSimpleIntegerDynamicPseudoCost (CbcModel * model, int sequence,
51                                    int iColumn, double breakEven)
52  : CbcSimpleInteger(model,sequence,iColumn,breakEven),
53    sumDownCost_(0.0),
54    sumUpCost_(0.0),
55    sumDownCostSquared_(0.0),
56    sumUpCostSquared_(0.0),
57    sumDownDecrease_(0.0),
58    sumUpDecrease_(0.0),
59    lastDownCost_(0.0),
60    lastUpCost_(0.0),
61    lastDownDecrease_(0),
62    lastUpDecrease_(0),
63    numberTimesDown_(0),
64    numberTimesUp_(0),
65    numberTimesDownInfeasible_(0),
66    numberTimesUpInfeasible_(0),
67    numberBeforeTrust_(0),
68    method_(0)
69{
70  const double * cost = model->getObjCoefficients();
71  double costValue = CoinMax(1.0e-5,fabs(cost[iColumn]));
72  // treat as if will cost what it says up
73  upDynamicPseudoCost_=costValue;
74  // and balance at breakeven
75  downDynamicPseudoCost_=((1.0-breakEven_)*upDynamicPseudoCost_)/breakEven_;
76}
77
78/** Useful constructor
79
80  Loads dynamic upper & lower bounds for the specified variable.
81*/
82CbcSimpleIntegerDynamicPseudoCost::CbcSimpleIntegerDynamicPseudoCost (CbcModel * model, int sequence,
83                                    int iColumn, double downDynamicPseudoCost,
84                                                        double upDynamicPseudoCost)
85  : CbcSimpleInteger(model,sequence,iColumn),
86    sumDownCost_(0.0),
87    sumUpCost_(0.0),
88    sumDownCostSquared_(0.0),
89    sumUpCostSquared_(0.0),
90    sumDownDecrease_(0.0),
91    sumUpDecrease_(0.0),
92    lastDownCost_(0.0),
93    lastUpCost_(0.0),
94    lastDownDecrease_(0),
95    lastUpDecrease_(0),
96    numberTimesDown_(0),
97    numberTimesUp_(0),
98    numberTimesDownInfeasible_(0),
99    numberTimesUpInfeasible_(0),
100    numberBeforeTrust_(0),
101    method_(0)
102{
103  downDynamicPseudoCost_ = downDynamicPseudoCost;
104  upDynamicPseudoCost_ = upDynamicPseudoCost;
105  breakEven_ = upDynamicPseudoCost_/(upDynamicPseudoCost_+downDynamicPseudoCost_);
106}
107
108// Copy constructor
109CbcSimpleIntegerDynamicPseudoCost::CbcSimpleIntegerDynamicPseudoCost ( const CbcSimpleIntegerDynamicPseudoCost & rhs)
110  :CbcSimpleInteger(rhs),
111   downDynamicPseudoCost_(rhs.downDynamicPseudoCost_),
112   upDynamicPseudoCost_(rhs.upDynamicPseudoCost_),
113   sumDownCost_(rhs.sumDownCost_),
114   sumUpCost_(rhs.sumUpCost_),
115   sumDownCostSquared_(rhs.sumDownCostSquared_),
116   sumUpCostSquared_(rhs.sumUpCostSquared_),
117   sumDownDecrease_(rhs.sumDownDecrease_),
118   sumUpDecrease_(rhs.sumUpDecrease_),
119   lastDownCost_(rhs.lastDownCost_),
120   lastUpCost_(rhs.lastUpCost_),
121   lastDownDecrease_(rhs.lastDownDecrease_),
122   lastUpDecrease_(rhs.lastUpDecrease_),
123   numberTimesDown_(rhs.numberTimesDown_),
124   numberTimesUp_(rhs.numberTimesUp_),
125   numberTimesDownInfeasible_(rhs.numberTimesDownInfeasible_),
126   numberTimesUpInfeasible_(rhs.numberTimesUpInfeasible_),
127   numberBeforeTrust_(rhs.numberBeforeTrust_),
128   method_(rhs.method_)
129
130{
131}
132
133// Clone
134CbcObject *
135CbcSimpleIntegerDynamicPseudoCost::clone() const
136{
137  return new CbcSimpleIntegerDynamicPseudoCost(*this);
138}
139
140// Assignment operator
141CbcSimpleIntegerDynamicPseudoCost & 
142CbcSimpleIntegerDynamicPseudoCost::operator=( const CbcSimpleIntegerDynamicPseudoCost& rhs)
143{
144  if (this!=&rhs) {
145    CbcSimpleInteger::operator=(rhs);
146    downDynamicPseudoCost_=rhs.downDynamicPseudoCost_;
147    upDynamicPseudoCost_=rhs.upDynamicPseudoCost_;
148    sumDownCost_ = rhs.sumDownCost_;
149    sumUpCost_ = rhs.sumUpCost_;
150    sumDownCostSquared_ = rhs.sumDownCostSquared_;
151    sumUpCostSquared_ = rhs.sumUpCostSquared_;
152    sumDownDecrease_ = rhs.sumDownDecrease_;
153    sumUpDecrease_ = rhs.sumUpDecrease_;
154    lastDownCost_ = rhs.lastDownCost_;
155    lastUpCost_ = rhs.lastUpCost_;
156    lastDownDecrease_ = rhs.lastDownDecrease_;
157    lastUpDecrease_ = rhs.lastUpDecrease_;
158    numberTimesDown_ = rhs.numberTimesDown_;
159    numberTimesUp_ = rhs.numberTimesUp_;
160    numberTimesDownInfeasible_ = rhs.numberTimesDownInfeasible_;
161    numberTimesUpInfeasible_ = rhs.numberTimesUpInfeasible_;
162    numberBeforeTrust_ = rhs.numberBeforeTrust_;
163    method_=rhs.method_;
164  }
165  return *this;
166}
167
168// Destructor
169CbcSimpleIntegerDynamicPseudoCost::~CbcSimpleIntegerDynamicPseudoCost ()
170{
171}
172// Creates a branching object
173CbcBranchingObject * 
174CbcSimpleIntegerDynamicPseudoCost::createBranch(int way) 
175{
176  OsiSolverInterface * solver = model_->solver();
177  const double * solution = model_->testSolution();
178  const double * lower = solver->getColLower();
179  const double * upper = solver->getColUpper();
180  double value = solution[columnNumber_];
181  value = CoinMax(value, lower[columnNumber_]);
182  value = CoinMin(value, upper[columnNumber_]);
183  double nearest = floor(value+0.5);
184  double integerTolerance = 
185    model_->getDblParam(CbcModel::CbcIntegerTolerance);
186  assert (upper[columnNumber_]>lower[columnNumber_]);
187  int hotstartStrategy=model_->getHotstartStrategy();
188  if (hotstartStrategy<=0) {
189    assert (fabs(value-nearest)>integerTolerance);
190  } else {
191    const double * bestSolution = model_->bestSolution();
192    double targetValue = bestSolution[columnNumber_];
193    if (way>0)
194      value = targetValue-0.1;
195    else
196      value = targetValue+0.1;
197  }
198  CbcDynamicPseudoCostBranchingObject * newObject = 
199    new CbcDynamicPseudoCostBranchingObject(model_,sequence_,way,
200                                            value,this);
201  double up =  upDynamicPseudoCost_*(ceil(value)-value);
202  double down =  downDynamicPseudoCost_*(value-floor(value));
203  double changeInGuessed=up-down;
204  if (way>0)
205    changeInGuessed = - changeInGuessed;
206  changeInGuessed=CoinMax(0.0,changeInGuessed);
207  //if (way>0)
208  //changeInGuessed += 1.0e8; // bias to stay up
209  newObject->setChangeInGuessed(changeInGuessed);
210  return newObject;
211}
212// Infeasibility - large is 0.5
213double 
214CbcSimpleIntegerDynamicPseudoCost::infeasibility(int & preferredWay) const
215{
216  OsiSolverInterface * solver = model_->solver();
217  const double * solution = model_->testSolution();
218  const double * lower = solver->getColLower();
219  const double * upper = solver->getColUpper();
220  if (upper[columnNumber_]==lower[columnNumber_]) {
221    // fixed
222    preferredWay=1;
223    return 0.0;
224  }
225  double value = solution[columnNumber_];
226  value = CoinMax(value, lower[columnNumber_]);
227  value = CoinMin(value, upper[columnNumber_]);
228  /*printf("%d %g %g %g %g\n",columnNumber_,value,lower[columnNumber_],
229    solution[columnNumber_],upper[columnNumber_]);*/
230  double nearest = floor(value+0.5);
231  double integerTolerance = 
232    model_->getDblParam(CbcModel::CbcIntegerTolerance);
233  double below = floor(value+integerTolerance);
234  double above = below+1.0;
235  if (above>upper[columnNumber_]) {
236    above=below;
237    below = above -1;
238  }
239#define INFEAS
240#ifdef INFEAS
241  double distanceToCutoff=0.0;
242  double objectiveValue = solver->getObjSense()*solver->getObjValue();
243  distanceToCutoff =  model_->getCutoff()  - objectiveValue;
244  if (distanceToCutoff<1.0e20) 
245    distanceToCutoff *= 10.0;
246  else 
247    distanceToCutoff = 1.0e2 + fabs(objectiveValue);
248#endif
249  double sum;
250  int number;
251  double downCost = CoinMax(value-below,0.0);
252  sum = sumDownCost();
253  number = numberTimesDown();
254#ifdef INFEAS
255  sum += numberTimesDownInfeasible()*(distanceToCutoff/(downCost+1.0e-12));
256#endif
257  if (number>0)
258    downCost *= sum / (double) number;
259  else
260    downCost  *=  downDynamicPseudoCost_;
261  double upCost = CoinMax((above-value),0.0);
262  sum = sumUpCost();
263  number = numberTimesUp();
264#ifdef INFEAS
265  sum += numberTimesUpInfeasible()*(distanceToCutoff/(upCost+1.0e-12));
266#endif
267  if (number>0)
268    upCost *= sum / (double) number;
269  else
270    upCost  *=  upDynamicPseudoCost_;
271  if (downCost>=upCost)
272    preferredWay=1;
273  else
274    preferredWay=-1;
275#define WEIGHT_AFTER 0.8
276  if (fabs(value-nearest)<=integerTolerance) {
277    return 0.0;
278  } else {
279    int stateOfSearch = model_->stateOfSearch();
280    double returnValue=0.0;
281    double minValue = CoinMin(downCost,upCost);
282    double maxValue = CoinMax(downCost,upCost);
283    if (stateOfSearch==0) {
284      // no solution
285      returnValue = 0.5*minValue + 0.5*maxValue;
286    } else {
287      // some solution
288      returnValue = WEIGHT_AFTER*minValue + (1.0-WEIGHT_AFTER)*maxValue;
289    }
290    if (numberTimesUp_<numberBeforeTrust_||
291        numberTimesDown_<numberBeforeTrust_) {
292      //if (returnValue<1.0e10)
293      //returnValue += 1.0e12;
294      //else
295      returnValue *= 1.0e3;
296      if (!numberTimesUp_&&!numberTimesDown_)
297        returnValue=1.0e50;
298    }
299    return CoinMax(returnValue,1.0e-15);
300  }
301}
302
303// Return "up" estimate
304double 
305CbcSimpleIntegerDynamicPseudoCost::upEstimate() const
306{
307  OsiSolverInterface * solver = model_->solver();
308  const double * solution = model_->testSolution();
309  const double * lower = solver->getColLower();
310  const double * upper = solver->getColUpper();
311  double value = solution[columnNumber_];
312  value = CoinMax(value, lower[columnNumber_]);
313  value = CoinMin(value, upper[columnNumber_]);
314  if (upper[columnNumber_]==lower[columnNumber_]) {
315    // fixed
316    return 0.0;
317  }
318  double integerTolerance = 
319    model_->getDblParam(CbcModel::CbcIntegerTolerance);
320  double below = floor(value+integerTolerance);
321  double above = below+1.0;
322  if (above>upper[columnNumber_]) {
323    above=below;
324    below = above -1;
325  }
326  double upCost = CoinMax((above-value)*upDynamicPseudoCost_,0.0);
327  return upCost;
328}
329// Return "down" estimate
330double 
331CbcSimpleIntegerDynamicPseudoCost::downEstimate() const
332{
333  OsiSolverInterface * solver = model_->solver();
334  const double * solution = model_->testSolution();
335  const double * lower = solver->getColLower();
336  const double * upper = solver->getColUpper();
337  double value = solution[columnNumber_];
338  value = CoinMax(value, lower[columnNumber_]);
339  value = CoinMin(value, upper[columnNumber_]);
340  if (upper[columnNumber_]==lower[columnNumber_]) {
341    // fixed
342    return 0.0;
343  }
344  double integerTolerance = 
345    model_->getDblParam(CbcModel::CbcIntegerTolerance);
346  double below = floor(value+integerTolerance);
347  double above = below+1.0;
348  if (above>upper[columnNumber_]) {
349    above=below;
350    below = above -1;
351  }
352  double downCost = CoinMax((value-below)*downDynamicPseudoCost_,0.0);
353  return downCost;
354}
355// Print
356void 
357CbcSimpleIntegerDynamicPseudoCost::print(int type,double value) const
358{
359  if (!type) {
360    double meanDown =0.0;
361    double devDown =0.0;
362    if (numberTimesDown_) {
363      meanDown = sumDownCost_/(double) numberTimesDown_;
364      devDown = meanDown*meanDown + sumDownCostSquared_ - 
365        2.0*meanDown*sumDownCost_;
366      if (devDown>=0.0)
367        devDown = sqrt(devDown);
368    }
369    double meanUp =0.0;
370    double devUp =0.0;
371    if (numberTimesUp_) {
372      meanUp = sumUpCost_/(double) numberTimesUp_;
373      devUp = meanUp*meanUp + sumUpCostSquared_ - 
374        2.0*meanUp*sumUpCost_;
375      if (devUp>=0.0)
376        devUp = sqrt(devUp);
377    }
378#if 0
379    printf("%d down %d times (%d inf) mean %g (dev %g) up %d times (%d inf) mean %g (dev %g)\n",
380           columnNumber_,
381           numberTimesDown_,numberTimesDownInfeasible_,meanDown,devDown,
382           numberTimesUp_,numberTimesUpInfeasible_,meanUp,devUp);
383#else
384    printf("%d down %d times (%d inf) mean %g  up %d times (%d inf) mean %g\n",
385           columnNumber_,
386           numberTimesDown_,numberTimesDownInfeasible_,meanDown,
387           numberTimesUp_,numberTimesUpInfeasible_,meanUp);
388#endif
389  } else {
390    OsiSolverInterface * solver = model_->solver();
391    const double * upper = solver->getColUpper();
392    double integerTolerance = 
393      model_->getDblParam(CbcModel::CbcIntegerTolerance);
394    double below = floor(value+integerTolerance);
395    double above = below+1.0;
396    if (above>upper[columnNumber_]) {
397      above=below;
398      below = above -1;
399    }
400    double objectiveValue = solver->getObjSense()*solver->getObjValue();
401    double distanceToCutoff =  model_->getCutoff() - objectiveValue;
402    if (distanceToCutoff<1.0e20) 
403      distanceToCutoff *= 10.0;
404    else 
405      distanceToCutoff = 1.0e2 + fabs(objectiveValue);
406    double sum;
407    int number;
408    double downCost = CoinMax(value-below,0.0);
409    double downCost0 = downCost*downDynamicPseudoCost_;
410    sum = sumDownCost();
411    number = numberTimesDown();
412    sum += numberTimesDownInfeasible()*(distanceToCutoff/(downCost+1.0e-12));
413    if (number>0)
414      downCost *= sum / (double) number;
415    else
416      downCost  *=  downDynamicPseudoCost_;
417    double upCost = CoinMax((above-value),0.0);
418    double upCost0 = upCost*upDynamicPseudoCost_;
419    sum = sumUpCost();
420    number = numberTimesUp();
421    sum += numberTimesUpInfeasible()*(distanceToCutoff/(upCost+1.0e-12));
422    if (number>0)
423      upCost *= sum / (double) number;
424    else
425      upCost  *=  upDynamicPseudoCost_;
426    printf("%d down %d times %g (est %g)  up %d times %g (est %g)\n",
427           columnNumber_,
428           numberTimesDown_,downCost,downCost0,
429           numberTimesUp_,upCost,upCost0);
430  }
431}
432
433// Default Constructor
434CbcDynamicPseudoCostBranchingObject::CbcDynamicPseudoCostBranchingObject()
435  :CbcIntegerBranchingObject()
436{
437  changeInGuessed_=1.0e-5;
438  object_=NULL;
439}
440
441// Useful constructor
442CbcDynamicPseudoCostBranchingObject::CbcDynamicPseudoCostBranchingObject (CbcModel * model, 
443                                                                          int variable, 
444                                                                          int way , double value,
445                                       CbcSimpleIntegerDynamicPseudoCost * object) 
446  :CbcIntegerBranchingObject(model,variable,way,value)
447{
448  changeInGuessed_=1.0e-5;
449  object_=object;
450}
451// Useful constructor for fixing
452CbcDynamicPseudoCostBranchingObject::CbcDynamicPseudoCostBranchingObject (CbcModel * model, 
453                                                      int variable, int way,
454                                                      double lowerValue, 
455                                                      double upperValue)
456  :CbcIntegerBranchingObject(model,variable,way,lowerValue)
457{
458  changeInGuessed_=1.0e100;
459  object_=NULL;
460}
461 
462
463// Copy constructor
464CbcDynamicPseudoCostBranchingObject::CbcDynamicPseudoCostBranchingObject ( 
465                                 const CbcDynamicPseudoCostBranchingObject & rhs)
466  :CbcIntegerBranchingObject(rhs)
467{
468  changeInGuessed_ = rhs.changeInGuessed_;
469  object_=rhs.object_;
470}
471
472// Assignment operator
473CbcDynamicPseudoCostBranchingObject & 
474CbcDynamicPseudoCostBranchingObject::operator=( const CbcDynamicPseudoCostBranchingObject& rhs)
475{
476  if (this != &rhs) {
477    CbcIntegerBranchingObject::operator=(rhs);
478    changeInGuessed_ = rhs.changeInGuessed_;
479    object_=rhs.object_;
480  }
481  return *this;
482}
483CbcBranchingObject * 
484CbcDynamicPseudoCostBranchingObject::clone() const
485{ 
486  return (new CbcDynamicPseudoCostBranchingObject(*this));
487}
488
489
490// Destructor
491CbcDynamicPseudoCostBranchingObject::~CbcDynamicPseudoCostBranchingObject ()
492{
493}
494
495/*
496  Perform a branch by adjusting the bounds of the specified variable. Note
497  that each arm of the branch advances the object to the next arm by
498  advancing the value of way_.
499
500  Providing new values for the variable's lower and upper bounds for each
501  branching direction gives a little bit of additional flexibility and will
502  be easily extensible to multi-way branching.
503  Returns change in guessed objective on next branch
504*/
505double
506CbcDynamicPseudoCostBranchingObject::branch(bool normalBranch)
507{
508  CbcIntegerBranchingObject::branch(normalBranch);
509  return changeInGuessed_;
510}
511/* Some branchingObjects may claim to be able to skip
512   strong branching.  If so they have to fill in CbcStrongInfo.
513   The object mention in incoming CbcStrongInfo must match.
514   Returns nonzero if skip is wanted */
515int 
516CbcDynamicPseudoCostBranchingObject::fillStrongInfo( CbcStrongInfo & info)
517{
518  assert (object_);
519  assert (info.possibleBranch==this);
520  if (object_->numberTimesUp()<object_->numberBeforeTrust()||
521      object_->numberTimesDown()<object_->numberBeforeTrust()) {
522    return 0;
523  } else {
524    info.upMovement = object_->upDynamicPseudoCost()*(ceil(value_)-value_);
525    info.downMovement = object_->downDynamicPseudoCost()*(value_-floor(value_));
526    info.numIntInfeasUp  -= (int) (object_->sumUpDecrease()/
527                                   ((double) object_->numberTimesUp()));
528    info.numObjInfeasUp = 0;
529    info.finishedUp = false;
530    info.numItersUp = 0;
531    info.numIntInfeasDown  -= (int) (object_->sumDownDecrease()/
532                                   ((double) object_->numberTimesDown()));
533    info.numObjInfeasDown = 0;
534    info.finishedDown = false;
535    info.numItersDown = 0;
536    info.fix =0;
537    return 1;
538  }
539}
540 
541// Default Constructor
542CbcBranchDynamicDecision::CbcBranchDynamicDecision()
543  :CbcBranchDecision()
544{
545  bestCriterion_ = 0.0;
546  bestChangeUp_ = 0.0;
547  bestNumberUp_ = 0;
548  bestChangeDown_ = 0.0;
549  bestNumberDown_ = 0;
550  bestObject_ = NULL;
551}
552
553// Copy constructor
554CbcBranchDynamicDecision::CbcBranchDynamicDecision (
555                                    const CbcBranchDynamicDecision & rhs)
556  :CbcBranchDecision()
557{
558  bestCriterion_ = rhs.bestCriterion_;
559  bestChangeUp_ = rhs.bestChangeUp_;
560  bestNumberUp_ = rhs.bestNumberUp_;
561  bestChangeDown_ = rhs.bestChangeDown_;
562  bestNumberDown_ = rhs.bestNumberDown_;
563  bestObject_ = rhs.bestObject_;
564}
565
566CbcBranchDynamicDecision::~CbcBranchDynamicDecision()
567{
568}
569
570// Clone
571CbcBranchDecision * 
572CbcBranchDynamicDecision::clone() const
573{
574  return new CbcBranchDynamicDecision(*this);
575}
576
577// Initialize i.e. before start of choosing at a node
578void 
579CbcBranchDynamicDecision::initialize(CbcModel * model)
580{
581  bestCriterion_ = 0.0;
582  bestChangeUp_ = 0.0;
583  bestNumberUp_ = 0;
584  bestChangeDown_ = 0.0;
585  bestNumberDown_ = 0;
586  bestObject_ = NULL;
587}
588
589/* Saves a clone of current branching object.  Can be used to update
590      information on object causing branch - after branch */
591void 
592CbcBranchDynamicDecision::saveBranchingObject(CbcBranchingObject * object) 
593{
594  object_ = object->clone();
595}
596/* Pass in information on branch just done.
597   assumes object can get information from solver */
598void 
599CbcBranchDynamicDecision::updateInformation(OsiSolverInterface * solver,
600                                            const CbcNode * node)
601{
602  double originalValue=node->objectiveValue();
603  int originalUnsatisfied = node->numberUnsatisfied();
604  double objectiveValue = solver->getObjSense()*solver->getObjValue();
605  int unsatisfied=0;
606  int i;
607  int numberColumns = solver->getNumCols();
608  const double * solution = solver->getColSolution();
609  //const double * lower = solver->getColLower();
610  //const double * upper = solver->getColUpper();
611  assert (object_);
612  CbcDynamicPseudoCostBranchingObject * branchingObject =
613    dynamic_cast<CbcDynamicPseudoCostBranchingObject *>(object_);
614  assert (branchingObject);
615  CbcSimpleIntegerDynamicPseudoCost *  object = branchingObject->object();
616  double change = CoinMax(0.0,objectiveValue-originalValue);
617  // probably should also ignore if stopped
618  int iStatus;
619  if (solver->isProvenOptimal())
620    iStatus=0; // optimal
621  else if (solver->isIterationLimitReached()
622           &&!solver->isDualObjectiveLimitReached())
623    iStatus=2; // unknown
624  else
625    iStatus=1; // infeasible
626
627  bool feasible = iStatus!=1;
628  if (feasible) {
629    double integerTolerance = 
630      object->model()->getDblParam(CbcModel::CbcIntegerTolerance);
631    for (i=0;i<numberColumns;i++) {
632      if (solver->isInteger(i)) {
633        double value = solution[i];
634        double nearest = floor(value+0.5);
635        if (fabs(value-nearest)>integerTolerance) 
636          unsatisfied++;
637      }
638    }
639  }
640  int way = object_->way();
641  double value = object_->value();
642  if (way<0) {
643    // down
644    object->incrementNumberTimesDown();
645    if (feasible) {
646      object->addToSumDownCost(change/(1.0e-30+(value-floor(value))));
647      object->addToSumDownDecrease(originalUnsatisfied-unsatisfied);
648      object->setDownDynamicPseudoCost(object->sumDownCost()/(double) object->numberTimesDown());
649    } else {
650      object->incrementNumberTimesDownInfeasible();
651    }
652  } else {
653    // up
654    if (feasible) {
655      object->incrementNumberTimesUp();
656      object->addToSumUpCost(change/(1.0e-30+(ceil(value)-value)));
657      object->addToSumUpDecrease(unsatisfied-originalUnsatisfied);
658      object->setUpDynamicPseudoCost(object->sumUpCost()/(double) object->numberTimesUp());
659    } else {
660      object->incrementNumberTimesUpInfeasible();
661    }
662  }
663  delete object_;
664  object_=NULL;
665}
666
667/*
668  Simple dynamic decision algorithm. Compare based on infeasibility (numInfUp,
669  numInfDown) until a solution is found by search, then switch to change in
670  objective (changeUp, changeDown). Note that bestSoFar is remembered in
671  bestObject_, so the parameter bestSoFar is unused.
672*/
673
674int
675CbcBranchDynamicDecision::betterBranch(CbcBranchingObject * thisOne,
676                            CbcBranchingObject * bestSoFar,
677                            double changeUp, int numInfUp,
678                            double changeDown, int numInfDown)
679{
680  int stateOfSearch = thisOne->model()->stateOfSearch();
681  int betterWay=0;
682  double value=0.0;
683  if (!stateOfSearch) {
684    if (!bestObject_) {
685      bestNumberUp_=INT_MAX;
686      bestNumberDown_=INT_MAX;
687    }
688    // before solution - choose smallest number
689    // could add in depth as well
690    int bestNumber = CoinMin(bestNumberUp_,bestNumberDown_);
691    if (numInfUp<numInfDown) {
692      if (numInfUp<bestNumber) {
693        betterWay = 1;
694      } else if (numInfUp==bestNumber) {
695        if (changeUp<bestChangeUp_)
696          betterWay=1;
697      }
698    } else if (numInfUp>numInfDown) {
699      if (numInfDown<bestNumber) {
700        betterWay = -1;
701      } else if (numInfDown==bestNumber) {
702        if (changeDown<bestChangeDown_)
703          betterWay=-1;
704      }
705    } else {
706      // up and down have same number
707      bool better=false;
708      if (numInfUp<bestNumber) {
709        better=true;
710      } else if (numInfUp==bestNumber) {
711        if (CoinMin(changeUp,changeDown)<CoinMin(bestChangeUp_,bestChangeDown_)-1.0e-5)
712          better=true;;
713      }
714      if (better) {
715        // see which way
716        if (changeUp<=changeDown)
717          betterWay=1;
718        else
719          betterWay=-1;
720      }
721    }
722    if (betterWay) {
723      value = CoinMin(numInfUp,numInfDown);
724    }
725  } else {
726    if (!bestObject_) {
727      bestCriterion_=-1.0;
728    }
729    // got a solution
730    double minValue = CoinMin(changeDown,changeUp);
731    double maxValue = CoinMax(changeDown,changeUp);
732    value = WEIGHT_AFTER*minValue + (1.0-WEIGHT_AFTER)*maxValue;
733    if (value>bestCriterion_+1.0e-8) {
734      if (changeUp<=changeDown) {
735        betterWay=1;
736      } else {
737        betterWay=-1;
738      }
739    }
740  }
741  if (betterWay) {
742    bestCriterion_ = value;
743    bestChangeUp_ = changeUp;
744    bestNumberUp_ = numInfUp;
745    bestChangeDown_ = changeDown;
746    bestNumberDown_ = numInfDown;
747    bestObject_=thisOne;
748  }
749  return betterWay;
750}
Note: See TracBrowser for help on using the repository browser.