source: trunk/CbcBranchDynamic.cpp @ 216

Last change on this file since 216 was 208, checked in by forrest, 15 years ago

large number of changes

  • Property svn:eol-style set to native
  • Property svn:keywords set to Author Date Id Revision
File size: 23.6 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  const double * solution = model_->testSolution();
177  const double * lower = model_->getCbcColLower();
178  const double * upper = model_->getCbcColUpper();
179  double value = solution[columnNumber_];
180  value = CoinMax(value, lower[columnNumber_]);
181  value = CoinMin(value, upper[columnNumber_]);
182#ifndef NDEBUG
183  double nearest = floor(value+0.5);
184  double integerTolerance = 
185    model_->getDblParam(CbcModel::CbcIntegerTolerance);
186  assert (upper[columnNumber_]>lower[columnNumber_]);
187#endif
188  int hotstartStrategy=model_->getHotstartStrategy();
189  if (hotstartStrategy<=0) {
190    assert (fabs(value-nearest)>integerTolerance);
191  } else {
192    const double * bestSolution = model_->bestSolution();
193    double targetValue = bestSolution[columnNumber_];
194    if (way>0)
195      value = targetValue-0.1;
196    else
197      value = targetValue+0.1;
198  }
199  CbcDynamicPseudoCostBranchingObject * newObject = 
200    new CbcDynamicPseudoCostBranchingObject(model_,sequence_,way,
201                                            value,this);
202  double up =  upDynamicPseudoCost_*(ceil(value)-value);
203  double down =  downDynamicPseudoCost_*(value-floor(value));
204  double changeInGuessed=up-down;
205  if (way>0)
206    changeInGuessed = - changeInGuessed;
207  changeInGuessed=CoinMax(0.0,changeInGuessed);
208  //if (way>0)
209  //changeInGuessed += 1.0e8; // bias to stay up
210  newObject->setChangeInGuessed(changeInGuessed);
211  return newObject;
212}
213// Infeasibility - large is 0.5
214double 
215CbcSimpleIntegerDynamicPseudoCost::infeasibility(int & preferredWay) const
216{
217  const double * solution = model_->testSolution();
218  const double * lower = model_->getCbcColLower();
219  const double * upper = model_->getCbcColUpper();
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 = model_->getCurrentMinimizationObjValue();
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  // weight at 1.0 is max min
276#define WEIGHT_AFTER 0.9
277#define WEIGHT_BEFORE 0.3
278  if (fabs(value-nearest)<=integerTolerance) {
279    return 0.0;
280  } else {
281    int stateOfSearch = model_->stateOfSearch();
282    double returnValue=0.0;
283    double minValue = CoinMin(downCost,upCost);
284    double maxValue = CoinMax(downCost,upCost);
285    if (stateOfSearch<=1||model_->currentNode()->depth()<=10/* was ||maxValue>0.2*distanceToCutoff*/) {
286      // no solution
287      returnValue = WEIGHT_BEFORE*minValue + (1.0-WEIGHT_BEFORE)*maxValue;
288    } else {
289      // some solution
290      returnValue = WEIGHT_AFTER*minValue + (1.0-WEIGHT_AFTER)*maxValue;
291    }
292    if (numberTimesUp_<numberBeforeTrust_||
293        numberTimesDown_<numberBeforeTrust_) {
294      //if (returnValue<1.0e10)
295      //returnValue += 1.0e12;
296      //else
297      returnValue *= 1.0e3;
298      if (!numberTimesUp_&&!numberTimesDown_)
299        returnValue=1.0e50;
300    }
301    return CoinMax(returnValue,1.0e-15);
302  }
303}
304
305// Return "up" estimate
306double 
307CbcSimpleIntegerDynamicPseudoCost::upEstimate() const
308{
309  const double * solution = model_->testSolution();
310  const double * lower = model_->getCbcColLower();
311  const double * upper = model_->getCbcColUpper();
312  double value = solution[columnNumber_];
313  value = CoinMax(value, lower[columnNumber_]);
314  value = CoinMin(value, upper[columnNumber_]);
315  if (upper[columnNumber_]==lower[columnNumber_]) {
316    // fixed
317    return 0.0;
318  }
319  double integerTolerance = 
320    model_->getDblParam(CbcModel::CbcIntegerTolerance);
321  double below = floor(value+integerTolerance);
322  double above = below+1.0;
323  if (above>upper[columnNumber_]) {
324    above=below;
325    below = above -1;
326  }
327  double upCost = CoinMax((above-value)*upDynamicPseudoCost_,0.0);
328  return upCost;
329}
330// Return "down" estimate
331double 
332CbcSimpleIntegerDynamicPseudoCost::downEstimate() const
333{
334  const double * solution = model_->testSolution();
335  const double * lower = model_->getCbcColLower();
336  const double * upper = model_->getCbcColUpper();
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    const double * upper = model_->getCbcColUpper();
391    double integerTolerance = 
392      model_->getDblParam(CbcModel::CbcIntegerTolerance);
393    double below = floor(value+integerTolerance);
394    double above = below+1.0;
395    if (above>upper[columnNumber_]) {
396      above=below;
397      below = above -1;
398    }
399    double objectiveValue = model_->getCurrentMinimizationObjValue();
400    double distanceToCutoff =  model_->getCutoff() - objectiveValue;
401    if (distanceToCutoff<1.0e20) 
402      distanceToCutoff *= 10.0;
403    else 
404      distanceToCutoff = 1.0e2 + fabs(objectiveValue);
405    double sum;
406    int number;
407    double downCost = CoinMax(value-below,0.0);
408    double downCost0 = downCost*downDynamicPseudoCost_;
409    sum = sumDownCost();
410    number = numberTimesDown();
411    sum += numberTimesDownInfeasible()*(distanceToCutoff/(downCost+1.0e-12));
412    if (number>0)
413      downCost *= sum / (double) number;
414    else
415      downCost  *=  downDynamicPseudoCost_;
416    double upCost = CoinMax((above-value),0.0);
417    double upCost0 = upCost*upDynamicPseudoCost_;
418    sum = sumUpCost();
419    number = numberTimesUp();
420    sum += numberTimesUpInfeasible()*(distanceToCutoff/(upCost+1.0e-12));
421    if (number>0)
422      upCost *= sum / (double) number;
423    else
424      upCost  *=  upDynamicPseudoCost_;
425    printf("%d down %d times %g (est %g)  up %d times %g (est %g)\n",
426           columnNumber_,
427           numberTimesDown_,downCost,downCost0,
428           numberTimesUp_,upCost,upCost0);
429  }
430}
431
432// Default Constructor
433CbcDynamicPseudoCostBranchingObject::CbcDynamicPseudoCostBranchingObject()
434  :CbcIntegerBranchingObject()
435{
436  changeInGuessed_=1.0e-5;
437  object_=NULL;
438}
439
440// Useful constructor
441CbcDynamicPseudoCostBranchingObject::CbcDynamicPseudoCostBranchingObject (CbcModel * model, 
442                                                                          int variable, 
443                                                                          int way , double value,
444                                       CbcSimpleIntegerDynamicPseudoCost * object) 
445  :CbcIntegerBranchingObject(model,variable,way,value)
446{
447  changeInGuessed_=1.0e-5;
448  object_=object;
449}
450// Useful constructor for fixing
451CbcDynamicPseudoCostBranchingObject::CbcDynamicPseudoCostBranchingObject (CbcModel * model, 
452                                                      int variable, int way,
453                                                      double lowerValue, 
454                                                      double upperValue)
455  :CbcIntegerBranchingObject(model,variable,way,lowerValue)
456{
457  changeInGuessed_=1.0e100;
458  object_=NULL;
459}
460 
461
462// Copy constructor
463CbcDynamicPseudoCostBranchingObject::CbcDynamicPseudoCostBranchingObject ( 
464                                 const CbcDynamicPseudoCostBranchingObject & rhs)
465  :CbcIntegerBranchingObject(rhs)
466{
467  changeInGuessed_ = rhs.changeInGuessed_;
468  object_=rhs.object_;
469}
470
471// Assignment operator
472CbcDynamicPseudoCostBranchingObject & 
473CbcDynamicPseudoCostBranchingObject::operator=( const CbcDynamicPseudoCostBranchingObject& rhs)
474{
475  if (this != &rhs) {
476    CbcIntegerBranchingObject::operator=(rhs);
477    changeInGuessed_ = rhs.changeInGuessed_;
478    object_=rhs.object_;
479  }
480  return *this;
481}
482CbcBranchingObject * 
483CbcDynamicPseudoCostBranchingObject::clone() const
484{ 
485  return (new CbcDynamicPseudoCostBranchingObject(*this));
486}
487
488
489// Destructor
490CbcDynamicPseudoCostBranchingObject::~CbcDynamicPseudoCostBranchingObject ()
491{
492}
493
494/*
495  Perform a branch by adjusting the bounds of the specified variable. Note
496  that each arm of the branch advances the object to the next arm by
497  advancing the value of way_.
498
499  Providing new values for the variable's lower and upper bounds for each
500  branching direction gives a little bit of additional flexibility and will
501  be easily extensible to multi-way branching.
502  Returns change in guessed objective on next branch
503*/
504double
505CbcDynamicPseudoCostBranchingObject::branch(bool normalBranch)
506{
507  CbcIntegerBranchingObject::branch(normalBranch);
508  return changeInGuessed_;
509}
510/* Some branchingObjects may claim to be able to skip
511   strong branching.  If so they have to fill in CbcStrongInfo.
512   The object mention in incoming CbcStrongInfo must match.
513   Returns nonzero if skip is wanted */
514int 
515CbcDynamicPseudoCostBranchingObject::fillStrongInfo( CbcStrongInfo & info)
516{
517  assert (object_);
518  assert (info.possibleBranch==this);
519  if (object_->numberTimesUp()<object_->numberBeforeTrust()||
520      object_->numberTimesDown()<object_->numberBeforeTrust()) {
521    return 0;
522  } else {
523    info.upMovement = object_->upDynamicPseudoCost()*(ceil(value_)-value_);
524    info.downMovement = object_->downDynamicPseudoCost()*(value_-floor(value_));
525    info.numIntInfeasUp  -= (int) (object_->sumUpDecrease()/
526                                   ((double) object_->numberTimesUp()));
527    info.numObjInfeasUp = 0;
528    info.finishedUp = false;
529    info.numItersUp = 0;
530    info.numIntInfeasDown  -= (int) (object_->sumDownDecrease()/
531                                   ((double) object_->numberTimesDown()));
532    info.numObjInfeasDown = 0;
533    info.finishedDown = false;
534    info.numItersDown = 0;
535    info.fix =0;
536    return 1;
537  }
538}
539 
540// Default Constructor
541CbcBranchDynamicDecision::CbcBranchDynamicDecision()
542  :CbcBranchDecision()
543{
544  bestCriterion_ = 0.0;
545  bestChangeUp_ = 0.0;
546  bestNumberUp_ = 0;
547  bestChangeDown_ = 0.0;
548  bestNumberDown_ = 0;
549  bestObject_ = NULL;
550}
551
552// Copy constructor
553CbcBranchDynamicDecision::CbcBranchDynamicDecision (
554                                    const CbcBranchDynamicDecision & rhs)
555  :CbcBranchDecision()
556{
557  bestCriterion_ = rhs.bestCriterion_;
558  bestChangeUp_ = rhs.bestChangeUp_;
559  bestNumberUp_ = rhs.bestNumberUp_;
560  bestChangeDown_ = rhs.bestChangeDown_;
561  bestNumberDown_ = rhs.bestNumberDown_;
562  bestObject_ = rhs.bestObject_;
563}
564
565CbcBranchDynamicDecision::~CbcBranchDynamicDecision()
566{
567}
568
569// Clone
570CbcBranchDecision * 
571CbcBranchDynamicDecision::clone() const
572{
573  return new CbcBranchDynamicDecision(*this);
574}
575
576// Initialize i.e. before start of choosing at a node
577void 
578CbcBranchDynamicDecision::initialize(CbcModel * model)
579{
580  bestCriterion_ = 0.0;
581  bestChangeUp_ = 0.0;
582  bestNumberUp_ = 0;
583  bestChangeDown_ = 0.0;
584  bestNumberDown_ = 0;
585  bestObject_ = NULL;
586}
587
588/* Saves a clone of current branching object.  Can be used to update
589      information on object causing branch - after branch */
590void 
591CbcBranchDynamicDecision::saveBranchingObject(CbcBranchingObject * object) 
592{
593  object_ = object->clone();
594}
595/* Pass in information on branch just done.
596   assumes object can get information from solver */
597void 
598CbcBranchDynamicDecision::updateInformation(OsiSolverInterface * solver,
599                                            const CbcNode * node)
600{
601  assert (object_);
602  const CbcModel * model = object_->model();
603  double originalValue=node->objectiveValue();
604  int originalUnsatisfied = node->numberUnsatisfied();
605  double objectiveValue = solver->getObjValue()*model->getObjSense();
606  int unsatisfied=0;
607  int i;
608  int numberIntegers = model->numberIntegers();;
609  const double * solution = solver->getColSolution();
610  //const double * lower = solver->getColLower();
611  //const double * upper = solver->getColUpper();
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      model->getDblParam(CbcModel::CbcIntegerTolerance);
631    const int * integerVariable = model->integerVariable();
632    for (i=0;i<numberIntegers;i++) {
633      int j=integerVariable[i];
634      double value = solution[j];
635      double nearest = floor(value+0.5);
636      if (fabs(value-nearest)>integerTolerance) 
637        unsatisfied++;
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<=1&&thisOne->model()->currentNode()->depth()>10) {
684#if 0
685    if (!bestObject_) {
686      bestNumberUp_=INT_MAX;
687      bestNumberDown_=INT_MAX;
688    }
689    // before solution - choose smallest number
690    // could add in depth as well
691    int bestNumber = CoinMin(bestNumberUp_,bestNumberDown_);
692    if (numInfUp<numInfDown) {
693      if (numInfUp<bestNumber) {
694        betterWay = 1;
695      } else if (numInfUp==bestNumber) {
696        if (changeUp<bestChangeUp_)
697          betterWay=1;
698      }
699    } else if (numInfUp>numInfDown) {
700      if (numInfDown<bestNumber) {
701        betterWay = -1;
702      } else if (numInfDown==bestNumber) {
703        if (changeDown<bestChangeDown_)
704          betterWay=-1;
705      }
706    } else {
707      // up and down have same number
708      bool better=false;
709      if (numInfUp<bestNumber) {
710        better=true;
711      } else if (numInfUp==bestNumber) {
712        if (CoinMin(changeUp,changeDown)<CoinMin(bestChangeUp_,bestChangeDown_)-1.0e-5)
713          better=true;;
714      }
715      if (better) {
716        // see which way
717        if (changeUp<=changeDown)
718          betterWay=1;
719        else
720          betterWay=-1;
721      }
722    }
723    if (betterWay) {
724      value = CoinMin(numInfUp,numInfDown);
725    }
726#else
727    if (!bestObject_) {
728      bestCriterion_=-1.0;
729    }
730    // got a solution
731    double minValue = CoinMin(changeDown,changeUp);
732    double maxValue = CoinMax(changeDown,changeUp);
733    value = WEIGHT_BEFORE*minValue + (1.0-WEIGHT_BEFORE)*maxValue;
734    if (value>bestCriterion_+1.0e-8) {
735      if (changeUp<=changeDown) {
736        betterWay=1;
737      } else {
738        betterWay=-1;
739      }
740    }
741#endif
742  } else {
743    if (!bestObject_) {
744      bestCriterion_=-1.0;
745    }
746    // got a solution
747    double minValue = CoinMin(changeDown,changeUp);
748    double maxValue = CoinMax(changeDown,changeUp);
749    // Reduce
750    maxValue = CoinMin(maxValue,minValue*2.0);
751    value = WEIGHT_AFTER*minValue + (1.0-WEIGHT_AFTER)*maxValue;
752    if (value>bestCriterion_+1.0e-8) {
753      if (changeUp<=changeDown) {
754        betterWay=1;
755      } else {
756        betterWay=-1;
757      }
758    }
759  }
760  if (betterWay) {
761    bestCriterion_ = value;
762    bestChangeUp_ = changeUp;
763    bestNumberUp_ = numInfUp;
764    bestChangeDown_ = changeDown;
765    bestNumberDown_ = numInfDown;
766    bestObject_=thisOne;
767  }
768  return betterWay;
769}
Note: See TracBrowser for help on using the repository browser.