source: trunk/Cbc/examples/CbcBranchUser.cpp

Last change on this file was 1898, checked in by stefan, 6 years ago

fixup examples

  • Property svn:eol-style set to native
  • Property svn:keywords set to Author Date Id Revision
File size: 12.6 KB
Line 
1// $Id: CbcBranchUser.cpp 1898 2013-04-09 18:06:04Z tkr $
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#include <cassert>
7#include <cmath>
8#include <cfloat>
9
10#include "CoinPragma.hpp"
11#include "OsiSolverInterface.hpp"
12#include "CbcModel.hpp"
13#include "CbcMessage.hpp"
14#include "CbcBranchUser.hpp"
15#include "CoinSort.hpp"
16
17// Default Constructor // Default Constructor
18CbcBranchUserDecision::CbcBranchUserDecision()
19  :CbcBranchDecision()
20{
21}
22
23// Copy constructor
24CbcBranchUserDecision::CbcBranchUserDecision (
25                                    const CbcBranchUserDecision & rhs)
26  :CbcBranchDecision(rhs)
27{
28}
29
30CbcBranchUserDecision::~CbcBranchUserDecision()
31{
32}
33
34// Clone
35CbcBranchDecision * 
36CbcBranchUserDecision::clone() const
37{
38  return new CbcBranchUserDecision(*this);
39}
40
41// Initialize i.e. before start of choosing at a node
42void 
43CbcBranchUserDecision::initialize(CbcModel * model)
44{
45}
46
47
48/* Returns nonzero if branching on first object is "better" than on
49   second (if second NULL first wins). User can play with decision object.
50   This is only used after strong branching.  The initial selection
51   is done by infeasibility() for each CbcObject
52   return code +1 for up branch preferred, -1 for down
53   
54*/
55int
56CbcBranchUserDecision::betterBranch(CbcBranchingObject * thisOne,
57                            CbcBranchingObject * bestSoFar,
58                            double changeUp, int numberInfeasibilitiesUp,
59                            double changeDown, int numberInfeasibilitiesDown)
60{
61  printf("Now obsolete CbcBranchUserDecision::betterBranch\n");
62  abort();
63  return 0;
64}
65/* Compare N branching objects. Return index of best
66   and sets way of branching in chosen object.
67   
68   This routine is used only after strong branching.
69   This is reccommended version as it can be more sophisticated
70*/
71
72int
73CbcBranchUserDecision::bestBranch (CbcBranchingObject ** objects, int numberObjects,
74                                   int numberUnsatisfied,
75                                   double * changeUp, int * numberInfeasibilitiesUp,
76                                   double * changeDown, int * numberInfeasibilitiesDown,
77                                   double objectiveValue) 
78{
79
80  int bestWay=0;
81  int whichObject = -1;
82  if (numberObjects) {
83    CbcModel * model = objects[0]->model();
84    // at continuous
85    //double continuousObjective = model->getContinuousObjective();
86    //int continuousInfeasibilities = model->getContinuousInfeasibilities();
87   
88    // average cost to get rid of infeasibility
89    //double averageCostPerInfeasibility =
90    //(objectiveValue-continuousObjective)/
91    //(double) (abs(continuousInfeasibilities-numberUnsatisfied)+1);
92    /* beforeSolution is :
93       0 - before any solution
94       n - n heuristic solutions but no branched one
95       -1 - branched solution found
96    */
97    int numberSolutions = model->getSolutionCount();
98    double cutoff = model->getCutoff();
99    int method=0;
100    int i;
101    if (numberSolutions) {
102      int numberHeuristic = model->getNumberHeuristicSolutions();
103      if (numberHeuristic<numberSolutions) {
104        method = 1;
105      } else {
106        method = 2;
107        // look further
108        for ( i = 0 ; i < numberObjects ; i++) {
109          int numberNext = numberInfeasibilitiesUp[i];
110         
111          if (numberNext<numberUnsatisfied) {
112            int numberUp = numberUnsatisfied - numberInfeasibilitiesUp[i];
113            double perUnsatisfied = changeUp[i]/(double) numberUp;
114            double estimatedObjective = objectiveValue + numberUnsatisfied * perUnsatisfied;
115            if (estimatedObjective<cutoff) 
116              method=3;
117          }
118          numberNext = numberInfeasibilitiesDown[i];
119          if (numberNext<numberUnsatisfied) {
120            int numberDown = numberUnsatisfied - numberInfeasibilitiesDown[i];
121            double perUnsatisfied = changeDown[i]/(double) numberDown;
122            double estimatedObjective = objectiveValue + numberUnsatisfied * perUnsatisfied;
123            if (estimatedObjective<cutoff) 
124              method=3;
125          }
126        }
127      }
128      method=2;
129    } else {
130      method = 0;
131    }
132    // Uncomment next to force method 4
133    //method=4;
134    /* Methods :
135       0 - fewest infeasibilities
136       1 - largest min change in objective
137       2 - as 1 but use sum of changes if min close
138       3 - predicted best solution
139       4 - take cheapest up branch if infeasibilities same
140    */
141    int bestNumber=COIN_INT_MAX;
142    double bestCriterion=-1.0e50;
143    double alternativeCriterion = -1.0;
144    double bestEstimate = 1.0e100;
145    switch (method) {
146    case 0:
147      // could add in depth as well
148      for ( i = 0 ; i < numberObjects ; i++) {
149        int thisNumber = CoinMin(numberInfeasibilitiesUp[i],numberInfeasibilitiesDown[i]);
150        if (thisNumber<=bestNumber) {
151          int betterWay=0;
152          if (numberInfeasibilitiesUp[i]<numberInfeasibilitiesDown[i]) {
153            if (numberInfeasibilitiesUp[i]<bestNumber) {
154              betterWay = 1;
155            } else {
156              if (changeUp[i]<bestCriterion)
157                betterWay=1;
158            }
159          } else if (numberInfeasibilitiesUp[i]>numberInfeasibilitiesDown[i]) {
160            if (numberInfeasibilitiesDown[i]<bestNumber) {
161              betterWay = -1;
162            } else {
163              if (changeDown[i]<bestCriterion)
164                betterWay=-1;
165            }
166          } else {
167            // up and down have same number
168            bool better=false;
169            if (numberInfeasibilitiesUp[i]<bestNumber) {
170              better=true;
171            } else if (numberInfeasibilitiesUp[i]==bestNumber) {
172              if (CoinMin(changeUp[i],changeDown[i])<bestCriterion)
173                better=true;;
174            }
175            if (better) {
176              // see which way
177              if (changeUp[i]<=changeDown[i])
178                betterWay=1;
179              else
180                betterWay=-1;
181            }
182          }
183          if (betterWay) {
184            bestCriterion = CoinMin(changeUp[i],changeDown[i]);
185            bestNumber = thisNumber;
186            whichObject = i;
187            bestWay = betterWay;
188          }
189        }
190      }
191      break;
192    case 1:
193      for ( i = 0 ; i < numberObjects ; i++) {
194        int betterWay=0;
195        if (changeUp[i]<=changeDown[i]) {
196          if (changeUp[i]>bestCriterion)
197            betterWay=1;
198        } else {
199          if (changeDown[i]>bestCriterion)
200            betterWay=-1;
201        }
202        if (betterWay) {
203          bestCriterion = CoinMin(changeUp[i],changeDown[i]);
204          whichObject = i;
205          bestWay = betterWay;
206        }
207      }
208      break;
209    case 2:
210      for ( i = 0 ; i < numberObjects ; i++) {
211        double change = CoinMin(changeUp[i],changeDown[i]);
212        double sum = changeUp[i] + changeDown[i];
213        bool take=false;
214        if (change>1.1*bestCriterion) 
215          take=true;
216        else if (change>0.9*bestCriterion&&sum+change>bestCriterion+alternativeCriterion) 
217          take=true;
218        if (take) {
219          if (changeUp[i]<=changeDown[i]) {
220            if (changeUp[i]>bestCriterion)
221              bestWay=1;
222          } else {
223            if (changeDown[i]>bestCriterion)
224              bestWay=-1;
225          }
226          bestCriterion = change;
227          alternativeCriterion = sum;
228          whichObject = i;
229        }
230      }
231      break;
232    case 3:
233      for ( i = 0 ; i < numberObjects ; i++) {
234        int numberNext = numberInfeasibilitiesUp[i];
235       
236        if (numberNext<numberUnsatisfied) {
237          int numberUp = numberUnsatisfied - numberInfeasibilitiesUp[i];
238          double perUnsatisfied = changeUp[i]/(double) numberUp;
239          double estimatedObjective = objectiveValue + numberUnsatisfied * perUnsatisfied;
240          if (estimatedObjective<bestEstimate) {
241            bestEstimate = estimatedObjective;
242            bestWay=1;
243            whichObject=i;
244          }
245        }
246        numberNext = numberInfeasibilitiesDown[i];
247        if (numberNext<numberUnsatisfied) {
248          int numberDown = numberUnsatisfied - numberInfeasibilitiesDown[i];
249          double perUnsatisfied = changeDown[i]/(double) numberDown;
250          double estimatedObjective = objectiveValue + numberUnsatisfied * perUnsatisfied;
251          if (estimatedObjective<bestEstimate) {
252            bestEstimate = estimatedObjective;
253            bestWay=-1;
254            whichObject=i;
255          }
256        }
257      }
258      break;
259    case 4:
260      // if number infeas same then cheapest up
261      // first get best number or when going down
262      // now choose smallest change up amongst equal number infeas
263      for ( i = 0 ; i < numberObjects ; i++) {
264        int thisNumber = CoinMin(numberInfeasibilitiesUp[i],numberInfeasibilitiesDown[i]);
265        if (thisNumber<=bestNumber) {
266          int betterWay=0;
267          if (numberInfeasibilitiesUp[i]<numberInfeasibilitiesDown[i]) {
268            if (numberInfeasibilitiesUp[i]<bestNumber) {
269              betterWay = 1;
270            } else {
271              if (changeUp[i]<bestCriterion)
272                betterWay=1;
273            }
274          } else if (numberInfeasibilitiesUp[i]>numberInfeasibilitiesDown[i]) {
275            if (numberInfeasibilitiesDown[i]<bestNumber) {
276              betterWay = -1;
277            } else {
278              if (changeDown[i]<bestCriterion)
279                betterWay=-1;
280            }
281          } else {
282            // up and down have same number
283            bool better=false;
284            if (numberInfeasibilitiesUp[i]<bestNumber) {
285              better=true;
286            } else if (numberInfeasibilitiesUp[i]==bestNumber) {
287              if (CoinMin(changeUp[i],changeDown[i])<bestCriterion)
288                better=true;;
289            }
290            if (better) {
291              // see which way
292              if (changeUp[i]<=changeDown[i])
293                betterWay=1;
294              else
295                betterWay=-1;
296            }
297          }
298          if (betterWay) {
299            bestCriterion = CoinMin(changeUp[i],changeDown[i]);
300            bestNumber = thisNumber;
301            whichObject = i;
302            bestWay = betterWay;
303          }
304        }
305      }
306      bestCriterion=1.0e50;
307      for ( i = 0 ; i < numberObjects ; i++) {
308        int thisNumber = numberInfeasibilitiesUp[i];
309        if (thisNumber==bestNumber&&changeUp) {
310          if (changeUp[i]<bestCriterion) {
311            bestCriterion = changeUp[i];
312            whichObject = i;
313            bestWay = 1;
314          }
315        }
316      }
317      break;
318    }
319    // set way in best
320    if (whichObject>=0)
321      objects[whichObject]->way(bestWay);
322  }
323  return whichObject;
324}
325/** Default Constructor
326
327  Equivalent to an unspecified binary variable.
328*/
329CbcSimpleIntegerFixed::CbcSimpleIntegerFixed ()
330  : CbcSimpleInteger()
331{
332}
333
334/** Useful constructor
335
336  Loads actual upper & lower bounds for the specified variable.
337*/
338CbcSimpleIntegerFixed::CbcSimpleIntegerFixed (CbcModel * model,
339                                    int iColumn, double breakEven)
340  : CbcSimpleInteger(model,iColumn,breakEven)
341{
342}
343// Constructor from simple
344CbcSimpleIntegerFixed::CbcSimpleIntegerFixed (const CbcSimpleInteger & rhs)
345  : CbcSimpleInteger(rhs)
346{
347}
348
349// Copy constructor
350CbcSimpleIntegerFixed::CbcSimpleIntegerFixed ( const CbcSimpleIntegerFixed & rhs)
351  :CbcSimpleInteger(rhs)
352
353{
354}
355
356// Clone
357CbcObject *
358CbcSimpleIntegerFixed::clone() const
359{
360  return new CbcSimpleIntegerFixed(*this);
361}
362
363// Assignment operator
364CbcSimpleIntegerFixed & 
365CbcSimpleIntegerFixed::operator=( const CbcSimpleIntegerFixed& rhs)
366{
367  if (this!=&rhs) {
368    CbcSimpleInteger::operator=(rhs);
369  }
370  return *this;
371}
372
373// Destructor
374CbcSimpleIntegerFixed::~CbcSimpleIntegerFixed ()
375{
376}
377
378// Infeasibility - large is 0.5
379double 
380CbcSimpleIntegerFixed::infeasibility(int & preferredWay) const
381{
382  OsiSolverInterface * solver = model_->solver();
383  const double * solution = model_->testSolution();
384  const double * lower = solver->getColLower();
385  const double * upper = solver->getColUpper();
386  double value = solution[columnNumber_];
387  value = CoinMax(value, lower[columnNumber_]);
388  value = CoinMin(value, upper[columnNumber_]);
389  /*printf("%d %g %g %g %g\n",columnNumber_,value,lower[columnNumber_],
390    solution[columnNumber_],upper[columnNumber_]);*/
391  double nearest = floor(value+(1.0-breakEven_));
392  assert (breakEven_>0.0&&breakEven_<1.0);
393  double integerTolerance = 
394    model_->getDblParam(CbcModel::CbcIntegerTolerance);
395  if (nearest>value) 
396    preferredWay=1;
397  else
398    preferredWay=-1;
399  if (preferredWay_)
400    preferredWay=preferredWay_;
401  double weight = fabs(value-nearest);
402  // normalize so weight is 0.5 at break even
403  if (nearest<value)
404    weight = (0.5/breakEven_)*weight;
405  else
406    weight = (0.5/(1.0-breakEven_))*weight;
407  if (fabs(value-nearest)<=integerTolerance) {
408    if (upper[columnNumber_]==lower[columnNumber_])
409      return 0.0;
410    else
411      return 1.0e-5;
412  } else {
413    return weight;
414  }
415}
416// Creates a branching object
417CbcBranchingObject * 
418CbcSimpleIntegerFixed::createBranch(OsiSolverInterface * solver,
419                                            const OsiBranchingInformation * info, int way) 
420{
421  const double * solution = model_->testSolution();
422  const double * lower = solver->getColLower();
423  const double * upper = solver->getColUpper();
424  double value = solution[columnNumber_];
425  value = CoinMax(value, lower[columnNumber_]);
426  value = CoinMin(value, upper[columnNumber_]);
427  assert (upper[columnNumber_]>lower[columnNumber_]);
428  if (!model_->hotstartSolution()) {
429    double nearest = floor(value+0.5);
430    double integerTolerance = 
431    model_->getDblParam(CbcModel::CbcIntegerTolerance);
432    if (fabs(value-nearest)<integerTolerance) {
433      // adjust value
434      if (nearest!=upper[columnNumber_])
435        value = nearest+2.0*integerTolerance;
436      else
437        value = nearest-2.0*integerTolerance;
438    }
439  } else {
440    const double * hotstartSolution = model_->hotstartSolution();
441    double targetValue = hotstartSolution[columnNumber_];
442    if (way>0)
443      value = targetValue-0.1;
444    else
445      value = targetValue+0.1;
446  }
447  CbcBranchingObject * branch = new CbcIntegerBranchingObject(model_,columnNumber_,way,
448                                             value);
449  branch->setOriginalObject(this);
450  return branch;
451}
Note: See TracBrowser for help on using the repository browser.