source: branches/devel/Cbc/src/CbcBranchActual.cpp @ 648

Last change on this file since 648 was 619, checked in by forrest, 12 years ago

needed to go with Cbcstrategy

  • Property svn:eol-style set to native
  • Property svn:keywords set to Author Date Id Revision
File size: 87.7 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 "OsiSolverBranch.hpp"
14#include "CbcModel.hpp"
15#include "CbcMessage.hpp"
16#include "CbcBranchActual.hpp"
17#include "CoinSort.hpp"
18#include "CoinError.hpp"
19
20// Default Constructor
21CbcClique::CbcClique ()
22  : CbcObject(),
23    numberMembers_(0),
24    numberNonSOSMembers_(0),
25    members_(NULL),
26    type_(NULL),
27    cliqueType_(-1),
28    slack_(-1)
29{
30}
31
32// Useful constructor (which are integer indices)
33CbcClique::CbcClique (CbcModel * model, int cliqueType, int numberMembers,
34           const int * which, const char * type, int identifier,int slack)
35  : CbcObject(model)
36{
37  id_=identifier;
38  numberMembers_=numberMembers;
39  if (numberMembers_) {
40    members_ = new int[numberMembers_];
41    memcpy(members_,which,numberMembers_*sizeof(int));
42    type_ = new char[numberMembers_];
43    if (type) {
44      memcpy(type_,type,numberMembers_*sizeof(char));
45    } else {
46      for (int i=0;i<numberMembers_;i++)
47        type_[i]=1;
48    }
49  } else {
50    members_ = NULL;
51    type_ = NULL;
52  }
53  // Find out how many non sos
54  int i;
55  numberNonSOSMembers_=0;
56  for (i=0;i<numberMembers_;i++)
57    if (!type_[i])
58      numberNonSOSMembers_++;
59  cliqueType_ = cliqueType;
60  slack_ = slack;
61}
62
63// Copy constructor
64CbcClique::CbcClique ( const CbcClique & rhs)
65  :CbcObject(rhs)
66{
67  numberMembers_ = rhs.numberMembers_;
68  numberNonSOSMembers_ = rhs.numberNonSOSMembers_;
69  if (numberMembers_) {
70    members_ = new int[numberMembers_];
71    memcpy(members_,rhs.members_,numberMembers_*sizeof(int));
72    type_ = new char[numberMembers_];
73    memcpy(type_,rhs.type_,numberMembers_*sizeof(char));
74  } else {
75    members_ = NULL;
76    type_ = NULL;
77  }
78  cliqueType_ = rhs.cliqueType_;
79  slack_ = rhs.slack_;
80}
81
82// Clone
83CbcObject *
84CbcClique::clone() const
85{
86  return new CbcClique(*this);
87}
88
89// Assignment operator
90CbcClique & 
91CbcClique::operator=( const CbcClique& rhs)
92{
93  if (this!=&rhs) {
94    CbcObject::operator=(rhs);
95    delete [] members_;
96    delete [] type_;
97    numberMembers_ = rhs.numberMembers_;
98    numberNonSOSMembers_ = rhs.numberNonSOSMembers_;
99    if (numberMembers_) {
100      members_ = new int[numberMembers_];
101      memcpy(members_,rhs.members_,numberMembers_*sizeof(int));
102      type_ = new char[numberMembers_];
103      memcpy(type_,rhs.type_,numberMembers_*sizeof(char));
104    } else {
105      members_ = NULL;
106      type_ = NULL;
107    }
108    cliqueType_ = rhs.cliqueType_;
109    slack_ = rhs.slack_;
110  }
111  return *this;
112}
113
114// Destructor
115CbcClique::~CbcClique ()
116{
117  delete [] members_;
118  delete [] type_;
119}
120
121// Infeasibility - large is 0.5
122double 
123CbcClique::infeasibility(int & preferredWay) const
124{
125  int numberUnsatis=0, numberFree=0;
126  int j;
127  const int * integer = model_->integerVariable();
128  OsiSolverInterface * solver = model_->solver();
129  const double * solution = model_->testSolution();
130  const double * lower = solver->getColLower();
131  const double * upper = solver->getColUpper();
132  double largestValue=0.0;
133  double integerTolerance = 
134    model_->getDblParam(CbcModel::CbcIntegerTolerance);
135  double * sort = new double[numberMembers_];
136
137  double slackValue=0.0;
138  for (j=0;j<numberMembers_;j++) {
139    int sequence = members_[j];
140    int iColumn = integer[sequence];
141    double value = solution[iColumn];
142    value = CoinMax(value, lower[iColumn]);
143    value = CoinMin(value, upper[iColumn]);
144    double nearest = floor(value+0.5);
145    double distance = fabs(value-nearest);
146    if (distance>integerTolerance) {
147      if (!type_[j])
148        value = 1.0-value; // non SOS
149      // if slack then choose that
150      if (j==slack_&&value>0.05)
151        slackValue = value;
152      largestValue = CoinMax(value,largestValue);
153      sort[numberUnsatis++]=-value;
154    } else if (upper[iColumn]>lower[iColumn]) {
155      numberFree++;
156    }
157  }
158  preferredWay=1;
159  double otherWay = 0.0;
160  if (numberUnsatis) {
161    // sort
162    std::sort(sort,sort+numberUnsatis);
163    for (j=0;j<numberUnsatis;j++) {
164      if ((j&1)!=0)
165        otherWay += -sort[j];
166    }
167    // Need to think more
168    double value = 0.2*numberUnsatis+0.01*(numberMembers_-numberFree);
169    if (fabs(largestValue-0.5)<0.1) {
170      // close to half
171      value +=0.1;
172    }
173    if (slackValue) {
174      // branching on slack
175      value += slackValue;
176    }
177    // scale other way
178    otherWay *= value/(1.0-otherWay);
179    delete [] sort;
180    return value;
181  } else {
182    delete [] sort;
183    return 0.0; // satisfied
184  }
185}
186
187// This looks at solution and sets bounds to contain solution
188void 
189CbcClique::feasibleRegion()
190{
191  int j;
192  const int * integer = model_->integerVariable();
193  OsiSolverInterface * solver = model_->solver();
194  const double * solution = model_->testSolution();
195  const double * lower = solver->getColLower();
196  const double * upper = solver->getColUpper();
197#ifndef NDEBUG
198  double integerTolerance = 
199    model_->getDblParam(CbcModel::CbcIntegerTolerance);
200#endif 
201  for (j=0;j<numberMembers_;j++) {
202    int sequence = members_[j];
203    int iColumn = integer[sequence];
204    double value = solution[iColumn];
205    value = CoinMax(value, lower[iColumn]);
206    value = CoinMin(value, upper[iColumn]);
207    double nearest = floor(value+0.5);
208#ifndef NDEBUG
209    double distance = fabs(value-nearest);
210    assert(distance<=integerTolerance);
211#endif
212    solver->setColLower(iColumn,nearest);
213    solver->setColUpper(iColumn,nearest);
214  }
215}
216// Redoes data when sequence numbers change
217void 
218CbcClique::redoSequenceEtc(CbcModel * model, int numberColumns, const int * originalColumns)
219{
220  model_=model;
221  int n2=0;
222  for (int j=0;j<numberMembers_;j++) {
223    int iColumn = members_[j];
224    int i;
225    for (i=0;i<numberColumns;i++) {
226      if (originalColumns[i]==iColumn)
227        break;
228    }
229    if (i<numberColumns) {
230      members_[n2]=i;
231      type_[n2++]=type_[j];
232    }
233  }
234  if (n2<numberMembers_) {
235    //printf("** SOS number of members reduced from %d to %d!\n",numberMembers_,n2);
236    numberMembers_=n2;
237  }
238  // Find out how many non sos
239  int i;
240  numberNonSOSMembers_=0;
241  for (i=0;i<numberMembers_;i++)
242    if (!type_[i])
243      numberNonSOSMembers_++;
244}
245
246
247// Creates a branching object
248CbcBranchingObject * 
249CbcClique::createBranch(int way) 
250{
251  int numberUnsatis=0;
252  int j;
253  int nUp=0;
254  int nDown=0;
255  int numberFree=numberMembers_;
256  const int * integer = model_->integerVariable();
257  OsiSolverInterface * solver = model_->solver();
258  const double * solution = model_->testSolution();
259  const double * lower = solver->getColLower();
260  const double * upper = solver->getColUpper();
261  int * upList = new int[numberMembers_];
262  int * downList = new int[numberMembers_];
263  double * sort = new double[numberMembers_];
264  double integerTolerance = 
265      model_->getDblParam(CbcModel::CbcIntegerTolerance);
266
267  double slackValue=0.0;
268  for (j=0;j<numberMembers_;j++) {
269    int sequence = members_[j];
270    int iColumn = integer[sequence];
271    double value = solution[iColumn];
272    value = CoinMax(value, lower[iColumn]);
273    value = CoinMin(value, upper[iColumn]);
274    double nearest = floor(value+0.5);
275    double distance = fabs(value-nearest);
276    if (distance>integerTolerance) {
277      if (!type_[j])
278        value = 1.0-value; // non SOS
279      // if slack then choose that
280      if (j==slack_&&value>0.05)
281        slackValue = value;
282      value = -value; // for sort
283      upList[numberUnsatis]=j;
284      sort[numberUnsatis++]=value;
285    } else if (upper[iColumn]>lower[iColumn]) {
286      upList[--numberFree]=j;
287    }
288  }
289  assert (numberUnsatis);
290  if (!slackValue) {
291    // sort
292    CoinSort_2(sort,sort+numberUnsatis,upList);
293    // put first in up etc
294    int kWay=1;
295    for (j=0;j<numberUnsatis;j++) {
296      if (kWay>0) 
297        upList[nUp++]=upList[j];
298      else
299        downList[nDown++]=upList[j];
300      kWay = -kWay;
301    }
302    for (j=numberFree;j<numberMembers_;j++) {
303      if (kWay>0)
304        upList[nUp++]=upList[j];
305      else
306        downList[nDown++]=upList[j];
307      kWay = -kWay;
308    }
309  } else {
310    // put slack to 0 in first way
311    nUp = 1;
312    upList[0]=slack_;
313    for (j=0;j<numberUnsatis;j++) {
314      downList[nDown++]=upList[j];
315    }
316    for (j=numberFree;j<numberMembers_;j++) {
317      downList[nDown++]=upList[j];
318    }
319  }
320  // create object
321  CbcBranchingObject * branch;
322  if (numberMembers_ <=64)
323     branch = new CbcCliqueBranchingObject(model_,this,way,
324                                         nDown,downList,nUp,upList);
325  else
326    branch = new CbcLongCliqueBranchingObject(model_,this,way,
327                                            nDown,downList,nUp,upList);
328  delete [] upList;
329  delete [] downList;
330  delete [] sort;
331  return branch;
332}
333
334// Default Constructor
335CbcSOS::CbcSOS ()
336  : CbcObject(),
337    members_(NULL),
338    weights_(NULL),
339    numberMembers_(0),
340    sosType_(-1),
341    integerValued_(false)
342{
343}
344
345// Useful constructor (which are indices)
346CbcSOS::CbcSOS (CbcModel * model,  int numberMembers,
347           const int * which, const double * weights, int identifier,int type)
348  : CbcObject(model),
349    numberMembers_(numberMembers),
350    sosType_(type)
351{
352  id_=identifier;
353  integerValued_ = type==1;
354  if (numberMembers_) {
355    members_ = new int[numberMembers_];
356    weights_ = new double[numberMembers_];
357    memcpy(members_,which,numberMembers_*sizeof(int));
358    if (weights) {
359      memcpy(weights_,weights,numberMembers_*sizeof(double));
360    } else {
361      for (int i=0;i<numberMembers_;i++)
362        weights_[i]=i;
363    }
364    // sort so weights increasing
365    CoinSort_2(weights_,weights_+numberMembers_,members_);
366    double last = -COIN_DBL_MAX;
367    int i;
368    for (i=0;i<numberMembers_;i++) {
369      double possible = CoinMax(last+1.0e-10,weights_[i]);
370      weights_[i] = possible;
371      last=possible;
372    }
373  } else {
374    members_ = NULL;
375    weights_ = NULL;
376  }
377  assert (sosType_>0&&sosType_<3);
378}
379
380// Copy constructor
381CbcSOS::CbcSOS ( const CbcSOS & rhs)
382  :CbcObject(rhs)
383{
384  numberMembers_ = rhs.numberMembers_;
385  sosType_ = rhs.sosType_;
386  integerValued_ = rhs.integerValued_;
387  if (numberMembers_) {
388    members_ = new int[numberMembers_];
389    weights_ = new double[numberMembers_];
390    memcpy(members_,rhs.members_,numberMembers_*sizeof(int));
391    memcpy(weights_,rhs.weights_,numberMembers_*sizeof(double));
392  } else {
393    members_ = NULL;
394    weights_ = NULL;
395  }
396}
397
398// Clone
399CbcObject *
400CbcSOS::clone() const
401{
402  return new CbcSOS(*this);
403}
404
405// Assignment operator
406CbcSOS & 
407CbcSOS::operator=( const CbcSOS& rhs)
408{
409  if (this!=&rhs) {
410    CbcObject::operator=(rhs);
411    delete [] members_;
412    delete [] weights_;
413    numberMembers_ = rhs.numberMembers_;
414    sosType_ = rhs.sosType_;
415    integerValued_ = rhs.integerValued_;
416    if (numberMembers_) {
417      members_ = new int[numberMembers_];
418      weights_ = new double[numberMembers_];
419      memcpy(members_,rhs.members_,numberMembers_*sizeof(int));
420      memcpy(weights_,rhs.weights_,numberMembers_*sizeof(double));
421    } else {
422      members_ = NULL;
423      weights_ = NULL;
424    }
425  }
426  return *this;
427}
428
429// Destructor
430CbcSOS::~CbcSOS ()
431{
432  delete [] members_;
433  delete [] weights_;
434}
435
436// Infeasibility - large is 0.5
437double 
438CbcSOS::infeasibility(int & preferredWay) const
439{
440  int j;
441  int firstNonZero=-1;
442  int lastNonZero = -1;
443  OsiSolverInterface * solver = model_->solver();
444  const double * solution = model_->testSolution();
445  //const double * lower = solver->getColLower();
446  const double * upper = solver->getColUpper();
447  //double largestValue=0.0;
448  double integerTolerance = 
449    model_->getDblParam(CbcModel::CbcIntegerTolerance);
450  double weight = 0.0;
451  double sum =0.0;
452
453  // check bounds etc
454  double lastWeight=-1.0e100;
455  for (j=0;j<numberMembers_;j++) {
456    int iColumn = members_[j];
457    if (lastWeight>=weights_[j]-1.0e-7)
458      throw CoinError("Weights too close together in SOS","infeasibility","CbcSOS");
459    double value = CoinMax(0.0,solution[iColumn]);
460    sum += value;
461    if (value>integerTolerance&&upper[iColumn]) {
462      // Possibly due to scaling a fixed variable might slip through
463      if (value>upper[iColumn]) {
464        value=upper[iColumn];
465        // Could change to #ifdef CBC_DEBUG
466#ifndef NDEBUG
467        if (model_->messageHandler()->logLevel()>2)
468          printf("** Variable %d (%d) has value %g and upper bound of %g\n",
469                 iColumn,j,value,upper[iColumn]);
470#endif
471      } 
472      weight += weights_[j]*value;
473      if (firstNonZero<0)
474        firstNonZero=j;
475      lastNonZero=j;
476    }
477  }
478  preferredWay=1;
479  if (lastNonZero-firstNonZero>=sosType_) {
480    // find where to branch
481    assert (sum>0.0);
482    weight /= sum;
483    //int iWhere;
484    //for (iWhere=firstNonZero;iWhere<lastNonZero;iWhere++)
485    //if (weight<weights_[iWhere+1])
486    //break;
487    // probably best to use pseudo duals
488    double value = lastNonZero-firstNonZero+1;
489    value *= 0.5/((double) numberMembers_);
490    return value;
491  } else {
492    return 0.0; // satisfied
493  }
494}
495
496// This looks at solution and sets bounds to contain solution
497void 
498CbcSOS::feasibleRegion()
499{
500  int j;
501  int firstNonZero=-1;
502  int lastNonZero = -1;
503  OsiSolverInterface * solver = model_->solver();
504  const double * solution = model_->testSolution();
505  //const double * lower = solver->getColLower();
506  const double * upper = solver->getColUpper();
507  double integerTolerance = 
508    model_->getDblParam(CbcModel::CbcIntegerTolerance);
509  double weight = 0.0;
510  double sum =0.0;
511
512  for (j=0;j<numberMembers_;j++) {
513    int iColumn = members_[j];
514    double value = CoinMax(0.0,solution[iColumn]);
515    sum += value;
516    if (value>integerTolerance&&upper[iColumn]) {
517      weight += weights_[j]*value;
518      if (firstNonZero<0)
519        firstNonZero=j;
520      lastNonZero=j;
521    }
522  }
523  assert (lastNonZero-firstNonZero<sosType_) ;
524  for (j=0;j<firstNonZero;j++) {
525    int iColumn = members_[j];
526    solver->setColUpper(iColumn,0.0);
527  }
528  for (j=lastNonZero+1;j<numberMembers_;j++) {
529    int iColumn = members_[j];
530    solver->setColUpper(iColumn,0.0);
531  }
532}
533// Redoes data when sequence numbers change
534void 
535CbcSOS::redoSequenceEtc(CbcModel * model, int numberColumns, const int * originalColumns)
536{
537  model_=model;
538  int n2=0;
539  for (int j=0;j<numberMembers_;j++) {
540    int iColumn = members_[j];
541    int i;
542    for (i=0;i<numberColumns;i++) {
543      if (originalColumns[i]==iColumn)
544        break;
545    }
546    if (i<numberColumns) {
547      members_[n2]=i;
548      weights_[n2++]=weights_[j];
549    }
550  }
551  if (n2<numberMembers_) {
552    //printf("** SOS number of members reduced from %d to %d!\n",numberMembers_,n2);
553    numberMembers_=n2;
554  }
555}
556
557
558// Creates a branching object
559CbcBranchingObject * 
560CbcSOS::createBranch(int way) 
561{
562  int j;
563  const double * solution = model_->testSolution();
564  double integerTolerance = 
565      model_->getDblParam(CbcModel::CbcIntegerTolerance);
566  OsiSolverInterface * solver = model_->solver();
567  const double * upper = solver->getColUpper();
568  int firstNonFixed=-1;
569  int lastNonFixed=-1;
570  int firstNonZero=-1;
571  int lastNonZero = -1;
572  double weight = 0.0;
573  double sum =0.0;
574  for (j=0;j<numberMembers_;j++) {
575    int iColumn = members_[j];
576    if (upper[iColumn]) {
577      double value = CoinMax(0.0,solution[iColumn]);
578      sum += value;
579      if (firstNonFixed<0)
580        firstNonFixed=j;
581      lastNonFixed=j;
582      if (value>integerTolerance) {
583        weight += weights_[j]*value;
584        if (firstNonZero<0)
585          firstNonZero=j;
586        lastNonZero=j;
587      }
588    }
589  }
590  assert (lastNonZero-firstNonZero>=sosType_) ;
591  // find where to branch
592  assert (sum>0.0);
593  weight /= sum;
594  int iWhere;
595  double separator=0.0;
596  for (iWhere=firstNonZero;iWhere<lastNonZero;iWhere++) 
597    if (weight<weights_[iWhere+1])
598      break;
599  if (sosType_==1) {
600    // SOS 1
601    separator = 0.5 *(weights_[iWhere]+weights_[iWhere+1]);
602  } else {
603    // SOS 2
604    if (iWhere==firstNonFixed)
605      iWhere++;;
606    if (iWhere==lastNonFixed-1)
607      iWhere = lastNonFixed-2;
608    separator = weights_[iWhere+1];
609  }
610  // create object
611  CbcBranchingObject * branch;
612  branch = new CbcSOSBranchingObject(model_,this,way,separator);
613  return branch;
614}
615
616/* Create an OsiSolverBranch object
617   
618This returns NULL if branch not represented by bound changes
619*/
620OsiSolverBranch * 
621CbcSOS::solverBranch() const
622{
623  int j;
624  const double * solution = model_->testSolution();
625  double integerTolerance = 
626      model_->getDblParam(CbcModel::CbcIntegerTolerance);
627  OsiSolverInterface * solver = model_->solver();
628  const double * upper = solver->getColUpper();
629  int firstNonFixed=-1;
630  int lastNonFixed=-1;
631  int firstNonZero=-1;
632  int lastNonZero = -1;
633  double weight = 0.0;
634  double sum =0.0;
635  double * fix = new double[numberMembers_];
636  int * which = new int[numberMembers_];
637  for (j=0;j<numberMembers_;j++) {
638    int iColumn = members_[j];
639    // fix all on one side or other (even if fixed)
640    fix[j]=0.0;
641    which[j]=iColumn;
642    if (upper[iColumn]) {
643      double value = CoinMax(0.0,solution[iColumn]);
644      sum += value;
645      if (firstNonFixed<0)
646        firstNonFixed=j;
647      lastNonFixed=j;
648      if (value>integerTolerance) {
649        weight += weights_[j]*value;
650        if (firstNonZero<0)
651          firstNonZero=j;
652        lastNonZero=j;
653      }
654    }
655  }
656  assert (lastNonZero-firstNonZero>=sosType_) ;
657  // find where to branch
658  assert (sum>0.0);
659  weight /= sum;
660  // down branch fixes ones above weight to 0
661  int iWhere;
662  int iDownStart=0;
663  int iUpEnd=0;
664  for (iWhere=firstNonZero;iWhere<lastNonZero;iWhere++) 
665    if (weight<weights_[iWhere+1])
666      break;
667  if (sosType_==1) {
668    // SOS 1
669    iUpEnd=iWhere+1;
670    iDownStart=iUpEnd;
671  } else {
672    // SOS 2
673    if (iWhere==firstNonFixed)
674      iWhere++;;
675    if (iWhere==lastNonFixed-1)
676      iWhere = lastNonFixed-2;
677    iUpEnd=iWhere+1;
678    iDownStart=iUpEnd+1;
679  }
680  //
681  OsiSolverBranch * branch = new OsiSolverBranch();
682  branch->addBranch(-1,0,NULL,NULL,numberMembers_-iDownStart,which+iDownStart,fix);
683  branch->addBranch(1,0,NULL,NULL,iUpEnd,which,fix);
684  delete [] fix;
685  delete [] which;
686  return branch;
687}
688// Construct an OsiSOS object
689OsiSOS * 
690CbcSOS::osiObject(const OsiSolverInterface * solver) const
691{
692  OsiSOS * obj = new OsiSOS(solver,numberMembers_,members_,weights_,sosType_);
693  obj->setPriority(priority());
694  return obj;
695}
696
697/** Default Constructor
698
699  Equivalent to an unspecified binary variable.
700*/
701CbcSimpleInteger::CbcSimpleInteger ()
702  : CbcObject(),
703    originalLower_(0.0),
704    originalUpper_(1.0),
705    breakEven_(0.5),
706    columnNumber_(-1),
707    preferredWay_(0)
708{
709}
710
711/** Useful constructor
712
713  Loads actual upper & lower bounds for the specified variable.
714*/
715CbcSimpleInteger::CbcSimpleInteger ( CbcModel * model, int iColumn, double breakEven)
716  : CbcObject(model)
717{
718  columnNumber_ = iColumn ;
719  originalLower_ = model->solver()->getColLower()[columnNumber_] ;
720  originalUpper_ = model->solver()->getColUpper()[columnNumber_] ;
721  breakEven_ = breakEven;
722  assert (breakEven_>0.0&&breakEven_<1.0);
723  preferredWay_ = 0;
724}
725
726
727// Copy constructor
728CbcSimpleInteger::CbcSimpleInteger ( const CbcSimpleInteger & rhs)
729  :CbcObject(rhs)
730
731{
732  columnNumber_ = rhs.columnNumber_;
733  originalLower_ = rhs.originalLower_;
734  originalUpper_ = rhs.originalUpper_;
735  breakEven_ = rhs.breakEven_;
736  preferredWay_ = rhs.preferredWay_;
737}
738
739// Clone
740CbcObject *
741CbcSimpleInteger::clone() const
742{
743  return new CbcSimpleInteger(*this);
744}
745
746// Assignment operator
747CbcSimpleInteger & 
748CbcSimpleInteger::operator=( const CbcSimpleInteger& rhs)
749{
750  if (this!=&rhs) {
751    CbcObject::operator=(rhs);
752    columnNumber_ = rhs.columnNumber_;
753    originalLower_ = rhs.originalLower_;
754    originalUpper_ = rhs.originalUpper_;
755    breakEven_ = rhs.breakEven_;
756    preferredWay_ = rhs.preferredWay_;
757  }
758  return *this;
759}
760
761// Destructor
762CbcSimpleInteger::~CbcSimpleInteger ()
763{
764}
765// Construct an OsiSimpleInteger object
766OsiSimpleInteger * 
767CbcSimpleInteger::osiObject() const
768{
769  OsiSimpleInteger * obj = new OsiSimpleInteger(columnNumber_,
770                                                originalLower_,originalUpper_);
771  obj->setPriority(priority());
772  return obj;
773}
774
775double
776CbcSimpleInteger::infeasibility(const OsiSolverInterface * solver, const OsiBranchingInformation * info,
777                         int & preferredWay) const
778{
779  double value = info->solution_[columnNumber_];
780  value = CoinMax(value, info->lower_[columnNumber_]);
781  value = CoinMin(value, info->upper_[columnNumber_]);
782  double nearest = floor(value+(1.0-breakEven_));
783  assert (breakEven_>0.0&&breakEven_<1.0);
784  if (nearest>value) 
785    preferredWay=1;
786  else
787    preferredWay=-1;
788  if (preferredWay_)
789    preferredWay=preferredWay_;
790  double weight = fabs(value-nearest);
791  // normalize so weight is 0.5 at break even
792  if (nearest<value)
793    weight = (0.5/breakEven_)*weight;
794  else
795    weight = (0.5/(1.0-breakEven_))*weight;
796  if (fabs(value-nearest)<=info->integerTolerance_) 
797    return 0.0;
798  else
799    return weight;
800}
801double 
802CbcSimpleInteger::feasibleRegion(OsiSolverInterface * solver, const OsiBranchingInformation * info) const 
803{
804  double value = info->solution_[columnNumber_];
805  double newValue = CoinMax(value, info->lower_[columnNumber_]);
806  newValue = CoinMin(newValue, info->upper_[columnNumber_]);
807  newValue = floor(newValue+0.5);
808  solver->setColLower(columnNumber_,newValue);
809  solver->setColUpper(columnNumber_,newValue);
810  return fabs(value-newValue);
811}
812
813/* Create an OsiSolverBranch object
814   
815This returns NULL if branch not represented by bound changes
816*/
817OsiSolverBranch * 
818CbcSimpleInteger::solverBranch(OsiSolverInterface * solver, const OsiBranchingInformation * info) const
819{
820  double value = info->solution_[columnNumber_];
821  value = CoinMax(value, info->lower_[columnNumber_]);
822  value = CoinMin(value, info->upper_[columnNumber_]);
823  assert (info->upper_[columnNumber_]>info->lower_[columnNumber_]);
824#ifndef NDEBUG
825  double nearest = floor(value+0.5);
826  assert (fabs(value-nearest)>info->integerTolerance_);
827#endif
828  OsiSolverBranch * branch = new OsiSolverBranch();
829  branch->addBranch(columnNumber_,value);
830  return branch;
831}
832// Creates a branching object
833CbcBranchingObject * 
834CbcSimpleInteger::createBranch(OsiSolverInterface * solver, const OsiBranchingInformation * info, int way) 
835{
836  double value = info->solution_[columnNumber_];
837  value = CoinMax(value, info->lower_[columnNumber_]);
838  value = CoinMin(value, info->upper_[columnNumber_]);
839  assert (info->upper_[columnNumber_]>info->lower_[columnNumber_]);
840  if (!info->hotstartSolution_) {
841#ifndef NDEBUG
842    double nearest = floor(value+0.5);
843    assert (fabs(value-nearest)>info->integerTolerance_);
844#endif
845  } else {
846    double targetValue = info->hotstartSolution_[columnNumber_];
847    if (way>0)
848      value = targetValue-0.1;
849    else
850      value = targetValue+0.1;
851  }
852  CbcBranchingObject * branch = new CbcIntegerBranchingObject(model_,columnNumber_,way,
853                                             value);
854  branch->setOriginalObject(this);
855  return branch;
856}
857/* Column number if single column object -1 otherwise,
858   so returns >= 0
859   Used by heuristics
860*/
861int 
862CbcSimpleInteger::columnNumber() const
863{
864  return columnNumber_;
865}
866/* Reset variable bounds to their original values.
867 
868    Bounds may be tightened, so it may be good to be able to set this info in object.
869*/
870void 
871CbcSimpleInteger::resetBounds(const OsiSolverInterface * solver) 
872{
873  originalLower_ = solver->getColLower()[columnNumber_] ;
874  originalUpper_ = solver->getColUpper()[columnNumber_] ;
875}
876
877/*  Change column numbers after preprocessing
878 */
879void 
880CbcSimpleInteger::resetSequenceEtc(int numberColumns, const int * originalColumns) 
881{
882  int iColumn;
883  for (iColumn=0;iColumn<numberColumns;iColumn++) {
884    if (columnNumber_==originalColumns[iColumn])
885      break;
886  }
887  assert (iColumn<numberColumns);
888  columnNumber_ = iColumn;
889}
890
891// Infeasibility - large is 0.5
892double 
893CbcSimpleInteger::infeasibility(int & preferredWay) const
894{
895  OsiBranchingInformation info(model_->solver(),model_->normalSolver(),false);
896  return infeasibility(model_->solver(),&info,preferredWay);
897}
898
899// This looks at solution and sets bounds to contain solution
900/** More precisely: it first forces the variable within the existing
901    bounds, and then tightens the bounds to fix the variable at the
902    nearest integer value.
903*/
904void 
905CbcSimpleInteger::feasibleRegion()
906{
907  abort();
908}
909CbcBranchingObject * 
910CbcSimpleInteger::createBranch( int way) 
911{
912  abort();
913  return NULL;
914}
915// Default Constructor
916CbcIntegerBranchingObject::CbcIntegerBranchingObject()
917  :CbcBranchingObject()
918{
919  down_[0] = 0.0;
920  down_[1] = 0.0;
921  up_[0] = 0.0;
922  up_[1] = 0.0;
923}
924
925// Useful constructor
926CbcIntegerBranchingObject::CbcIntegerBranchingObject (CbcModel * model, 
927                                                      int variable, int way , double value)
928  :CbcBranchingObject(model,variable,way,value)
929{
930  int iColumn = variable;
931  down_[0] = model_->solver()->getColLower()[iColumn];
932  down_[1] = floor(value_);
933  up_[0] = ceil(value_);
934  up_[1] = model->getColUpper()[iColumn];
935}
936// Useful constructor for fixing
937CbcIntegerBranchingObject::CbcIntegerBranchingObject (CbcModel * model, 
938                                                      int variable, int way,
939                                                      double lowerValue, 
940                                                      double upperValue)
941  :CbcBranchingObject(model,variable,way,lowerValue)
942{
943  setNumberBranchesLeft(1);
944  down_[0] = lowerValue;
945  down_[1] = upperValue;
946  up_[0] = lowerValue;
947  up_[1] = upperValue;
948}
949 
950
951// Copy constructor
952CbcIntegerBranchingObject::CbcIntegerBranchingObject ( const CbcIntegerBranchingObject & rhs) :CbcBranchingObject(rhs)
953{
954  down_[0] = rhs.down_[0];
955  down_[1] = rhs.down_[1];
956  up_[0] = rhs.up_[0];
957  up_[1] = rhs.up_[1];
958}
959
960// Assignment operator
961CbcIntegerBranchingObject & 
962CbcIntegerBranchingObject::operator=( const CbcIntegerBranchingObject& rhs)
963{
964  if (this != &rhs) {
965    CbcBranchingObject::operator=(rhs);
966    down_[0] = rhs.down_[0];
967    down_[1] = rhs.down_[1];
968    up_[0] = rhs.up_[0];
969    up_[1] = rhs.up_[1];
970  }
971  return *this;
972}
973CbcBranchingObject * 
974CbcIntegerBranchingObject::clone() const
975{ 
976  return (new CbcIntegerBranchingObject(*this));
977}
978
979
980// Destructor
981CbcIntegerBranchingObject::~CbcIntegerBranchingObject ()
982{
983  // for debugging threads
984  way_=-23456789;
985}
986
987/*
988  Perform a branch by adjusting the bounds of the specified variable. Note
989  that each arm of the branch advances the object to the next arm by
990  advancing the value of way_.
991
992  Providing new values for the variable's lower and upper bounds for each
993  branching direction gives a little bit of additional flexibility and will
994  be easily extensible to multi-way branching.
995  Returns change in guessed objective on next branch
996*/
997double
998CbcIntegerBranchingObject::branch()
999{
1000  // for debugging threads
1001  if (way_<-1||way_>100000) {
1002    printf("way %d, left %d, iCol %d, variable %d\n",
1003           way_,numberBranchesLeft(),
1004           originalCbcObject_->columnNumber(),variable_);
1005    assert (way_!=-23456789);
1006  }
1007  decrementNumberBranchesLeft();
1008  int iColumn = originalCbcObject_->columnNumber();
1009  assert (variable_==iColumn);
1010  double olb,oub ;
1011  olb = model_->solver()->getColLower()[iColumn] ;
1012  oub = model_->solver()->getColUpper()[iColumn] ;
1013#ifdef COIN_DEVELOP
1014  if (olb!=down_[0]||oub!=up_[1]) {
1015    if (way_>0)
1016      printf("branching up on var %d: [%g,%g] => [%g,%g] - other [%g,%g]\n",
1017             iColumn,olb,oub,up_[0],up_[1],down_[0],down_[1]) ; 
1018    else
1019      printf("branching down on var %d: [%g,%g] => [%g,%g] - other [%g,%g]\n",
1020             iColumn,olb,oub,down_[0],down_[1],up_[0],up_[1]) ; 
1021  }
1022#endif
1023  if (way_<0) {
1024#ifdef CBC_DEBUG
1025  { double olb,oub ;
1026    olb = model_->solver()->getColLower()[iColumn] ;
1027    oub = model_->solver()->getColUpper()[iColumn] ;
1028    printf("branching down on var %d: [%g,%g] => [%g,%g]\n",
1029           iColumn,olb,oub,down_[0],down_[1]) ; }
1030#endif
1031    model_->solver()->setColLower(iColumn,down_[0]);
1032    model_->solver()->setColUpper(iColumn,down_[1]);
1033    way_=1;
1034  } else {
1035#ifdef CBC_DEBUG
1036  { double olb,oub ;
1037    olb = model_->solver()->getColLower()[iColumn] ;
1038    oub = model_->solver()->getColUpper()[iColumn] ;
1039    printf("branching up on var %d: [%g,%g] => [%g,%g]\n",
1040           iColumn,olb,oub,up_[0],up_[1]) ; }
1041#endif
1042    model_->solver()->setColLower(iColumn,up_[0]);
1043    model_->solver()->setColUpper(iColumn,up_[1]);
1044    way_=-1;      // Swap direction
1045  }
1046  double nlb = model_->solver()->getColLower()[iColumn];
1047  double nub = model_->solver()->getColUpper()[iColumn];
1048  if (nlb<olb) {
1049#ifndef NDEBUG
1050    printf("bad lb change for column %d from %g to %g\n",iColumn,olb,nlb);
1051#endif
1052    model_->solver()->setColLower(iColumn,CoinMin(olb,nub));
1053    nlb=olb;
1054  }
1055  if (nub>oub) {
1056#ifndef NDEBUG
1057    printf("bad ub change for column %d from %g to %g\n",iColumn,oub,nub);
1058#endif
1059    model_->solver()->setColUpper(iColumn,CoinMax(oub,nlb));
1060  }
1061#ifndef NDEBUG
1062  if (nlb<olb+1.0e-8&&nub>oub-1.0e-8)
1063    printf("bad null change for column %d - bounds %g,%g\n",iColumn,olb,oub);
1064#endif
1065  return 0.0;
1066}
1067// Print what would happen 
1068void
1069CbcIntegerBranchingObject::print()
1070{
1071  int iColumn = originalCbcObject_->columnNumber();
1072  assert (variable_==iColumn);
1073  if (way_<0) {
1074  { double olb,oub ;
1075    olb = model_->solver()->getColLower()[iColumn] ;
1076    oub = model_->solver()->getColUpper()[iColumn] ;
1077    printf("CbcInteger would branch down on var %d (int var %d): [%g,%g] => [%g,%g]\n",
1078           iColumn,variable_,olb,oub,down_[0],down_[1]) ; }
1079  } else {
1080  { double olb,oub ;
1081    olb = model_->solver()->getColLower()[iColumn] ;
1082    oub = model_->solver()->getColUpper()[iColumn] ;
1083    printf("CbcInteger would branch up on var %d (int var %d): [%g,%g] => [%g,%g]\n",
1084           iColumn,variable_,olb,oub,up_[0],up_[1]) ; }
1085  }
1086}
1087
1088
1089/** Default Constructor
1090
1091  Equivalent to an unspecified binary variable.
1092*/
1093CbcSimpleIntegerPseudoCost::CbcSimpleIntegerPseudoCost ()
1094  : CbcSimpleInteger(),
1095    downPseudoCost_(1.0e-5),
1096    upPseudoCost_(1.0e-5),
1097    upDownSeparator_(-1.0),
1098    method_(0)
1099{
1100}
1101
1102/** Useful constructor
1103
1104  Loads actual upper & lower bounds for the specified variable.
1105*/
1106CbcSimpleIntegerPseudoCost::CbcSimpleIntegerPseudoCost (CbcModel * model,
1107                                    int iColumn, double breakEven)
1108  : CbcSimpleInteger(model,iColumn,breakEven)
1109{
1110  const double * cost = model->getObjCoefficients();
1111  double costValue = CoinMax(1.0e-5,fabs(cost[iColumn]));
1112  // treat as if will cost what it says up
1113  upPseudoCost_=costValue;
1114  // and balance at breakeven
1115  downPseudoCost_=((1.0-breakEven_)*upPseudoCost_)/breakEven_;
1116  upDownSeparator_ = -1.0;
1117  method_=0;
1118}
1119
1120/** Useful constructor
1121
1122  Loads actual upper & lower bounds for the specified variable.
1123*/
1124CbcSimpleIntegerPseudoCost::CbcSimpleIntegerPseudoCost (CbcModel * model,
1125                                    int iColumn, double downPseudoCost,
1126                                                        double upPseudoCost)
1127  : CbcSimpleInteger(model,iColumn)
1128{
1129  downPseudoCost_ = CoinMax(1.0e-10,downPseudoCost);
1130  upPseudoCost_ = CoinMax(1.0e-10,upPseudoCost);
1131  breakEven_ = upPseudoCost_/(upPseudoCost_+downPseudoCost_);
1132  upDownSeparator_ = -1.0;
1133  method_=0;
1134}
1135// Useful constructor - passed and model index and pseudo costs
1136CbcSimpleIntegerPseudoCost::CbcSimpleIntegerPseudoCost (CbcModel * model, int dummy,int iColumn, 
1137                                                        double downPseudoCost, double upPseudoCost)
1138{
1139  CbcSimpleIntegerPseudoCost(model,iColumn,downPseudoCost,upPseudoCost);
1140}
1141
1142// Copy constructor
1143CbcSimpleIntegerPseudoCost::CbcSimpleIntegerPseudoCost ( const CbcSimpleIntegerPseudoCost & rhs)
1144  :CbcSimpleInteger(rhs),
1145   downPseudoCost_(rhs.downPseudoCost_),
1146   upPseudoCost_(rhs.upPseudoCost_),
1147   upDownSeparator_(rhs.upDownSeparator_),
1148   method_(rhs.method_)
1149
1150{
1151}
1152
1153// Clone
1154CbcObject *
1155CbcSimpleIntegerPseudoCost::clone() const
1156{
1157  return new CbcSimpleIntegerPseudoCost(*this);
1158}
1159
1160// Assignment operator
1161CbcSimpleIntegerPseudoCost & 
1162CbcSimpleIntegerPseudoCost::operator=( const CbcSimpleIntegerPseudoCost& rhs)
1163{
1164  if (this!=&rhs) {
1165    CbcSimpleInteger::operator=(rhs);
1166    downPseudoCost_=rhs.downPseudoCost_;
1167    upPseudoCost_=rhs.upPseudoCost_;
1168    upDownSeparator_=rhs.upDownSeparator_;
1169    method_=rhs.method_;
1170  }
1171  return *this;
1172}
1173
1174// Destructor
1175CbcSimpleIntegerPseudoCost::~CbcSimpleIntegerPseudoCost ()
1176{
1177}
1178// Creates a branching object
1179CbcBranchingObject * 
1180CbcSimpleIntegerPseudoCost::createBranch(int way) 
1181{
1182  OsiSolverInterface * solver = model_->solver();
1183  const double * solution = model_->testSolution();
1184  const double * lower = solver->getColLower();
1185  const double * upper = solver->getColUpper();
1186  double value = solution[columnNumber_];
1187  value = CoinMax(value, lower[columnNumber_]);
1188  value = CoinMin(value, upper[columnNumber_]);
1189#ifndef NDEBUG
1190  double nearest = floor(value+0.5);
1191  double integerTolerance = 
1192    model_->getDblParam(CbcModel::CbcIntegerTolerance);
1193  assert (upper[columnNumber_]>lower[columnNumber_]);
1194#endif
1195  if (!model_->hotstartSolution()) {
1196    assert (fabs(value-nearest)>integerTolerance);
1197  } else {
1198    const double * hotstartSolution = model_->hotstartSolution();
1199    double targetValue = hotstartSolution[columnNumber_];
1200    if (way>0)
1201      value = targetValue-0.1;
1202    else
1203      value = targetValue+0.1;
1204  }
1205  CbcIntegerPseudoCostBranchingObject * newObject = 
1206    new CbcIntegerPseudoCostBranchingObject(model_,columnNumber_,way,
1207                                            value);
1208  double up =  upPseudoCost_*(ceil(value)-value);
1209  double down =  downPseudoCost_*(value-floor(value));
1210  double changeInGuessed=up-down;
1211  if (way>0)
1212    changeInGuessed = - changeInGuessed;
1213  changeInGuessed=CoinMax(0.0,changeInGuessed);
1214  //if (way>0)
1215  //changeInGuessed += 1.0e8; // bias to stay up
1216  newObject->setChangeInGuessed(changeInGuessed);
1217  newObject->setOriginalObject(this);
1218  return newObject;
1219}
1220// Infeasibility - large is 0.5
1221double 
1222CbcSimpleIntegerPseudoCost::infeasibility(int & preferredWay) const
1223{
1224  OsiSolverInterface * solver = model_->solver();
1225  const double * solution = model_->testSolution();
1226  const double * lower = solver->getColLower();
1227  const double * upper = solver->getColUpper();
1228  if (upper[columnNumber_]==lower[columnNumber_]) {
1229    // fixed
1230    preferredWay=1;
1231    return 0.0;
1232  }
1233  double value = solution[columnNumber_];
1234  value = CoinMax(value, lower[columnNumber_]);
1235  value = CoinMin(value, upper[columnNumber_]);
1236  /*printf("%d %g %g %g %g\n",columnNumber_,value,lower[columnNumber_],
1237    solution[columnNumber_],upper[columnNumber_]);*/
1238  double nearest = floor(value+0.5);
1239  double integerTolerance = 
1240    model_->getDblParam(CbcModel::CbcIntegerTolerance);
1241  double below = floor(value+integerTolerance);
1242  double above = below+1.0;
1243  if (above>upper[columnNumber_]) {
1244    above=below;
1245    below = above -1;
1246  }
1247  double downCost = CoinMax((value-below)*downPseudoCost_,0.0);
1248  double upCost = CoinMax((above-value)*upPseudoCost_,0.0);
1249  if (downCost>=upCost)
1250    preferredWay=1;
1251  else
1252    preferredWay=-1;
1253  // See if up down choice set
1254  if (upDownSeparator_>0.0) {
1255    preferredWay = (value-below>=upDownSeparator_) ? 1 : -1;
1256  }
1257  if (preferredWay_)
1258    preferredWay=preferredWay_;
1259  if (fabs(value-nearest)<=integerTolerance) {
1260    return 0.0;
1261  } else {
1262    // can't get at model so 1,2 don't make sense
1263    assert(method_<1||method_>2);
1264    if (!method_)
1265      return CoinMin(downCost,upCost);
1266    else
1267      return CoinMax(downCost,upCost);
1268  }
1269}
1270
1271// Return "up" estimate
1272double 
1273CbcSimpleIntegerPseudoCost::upEstimate() const
1274{
1275  OsiSolverInterface * solver = model_->solver();
1276  const double * solution = model_->testSolution();
1277  const double * lower = solver->getColLower();
1278  const double * upper = solver->getColUpper();
1279  double value = solution[columnNumber_];
1280  value = CoinMax(value, lower[columnNumber_]);
1281  value = CoinMin(value, upper[columnNumber_]);
1282  if (upper[columnNumber_]==lower[columnNumber_]) {
1283    // fixed
1284    return 0.0;
1285  }
1286  double integerTolerance = 
1287    model_->getDblParam(CbcModel::CbcIntegerTolerance);
1288  double below = floor(value+integerTolerance);
1289  double above = below+1.0;
1290  if (above>upper[columnNumber_]) {
1291    above=below;
1292    below = above -1;
1293  }
1294  double upCost = CoinMax((above-value)*upPseudoCost_,0.0);
1295  return upCost;
1296}
1297// Return "down" estimate
1298double 
1299CbcSimpleIntegerPseudoCost::downEstimate() const
1300{
1301  OsiSolverInterface * solver = model_->solver();
1302  const double * solution = model_->testSolution();
1303  const double * lower = solver->getColLower();
1304  const double * upper = solver->getColUpper();
1305  double value = solution[columnNumber_];
1306  value = CoinMax(value, lower[columnNumber_]);
1307  value = CoinMin(value, upper[columnNumber_]);
1308  if (upper[columnNumber_]==lower[columnNumber_]) {
1309    // fixed
1310    return 0.0;
1311  }
1312  double integerTolerance = 
1313    model_->getDblParam(CbcModel::CbcIntegerTolerance);
1314  double below = floor(value+integerTolerance);
1315  double above = below+1.0;
1316  if (above>upper[columnNumber_]) {
1317    above=below;
1318    below = above -1;
1319  }
1320  double downCost = CoinMax((value-below)*downPseudoCost_,0.0);
1321  return downCost;
1322}
1323
1324// Default Constructor
1325CbcIntegerPseudoCostBranchingObject::CbcIntegerPseudoCostBranchingObject()
1326  :CbcIntegerBranchingObject()
1327{
1328  changeInGuessed_=1.0e-5;
1329}
1330
1331// Useful constructor
1332CbcIntegerPseudoCostBranchingObject::CbcIntegerPseudoCostBranchingObject (CbcModel * model, 
1333                                                      int variable, int way , double value)
1334  :CbcIntegerBranchingObject(model,variable,way,value)
1335{
1336}
1337// Useful constructor for fixing
1338CbcIntegerPseudoCostBranchingObject::CbcIntegerPseudoCostBranchingObject (CbcModel * model, 
1339                                                      int variable, int way,
1340                                                      double lowerValue, 
1341                                                      double upperValue)
1342  :CbcIntegerBranchingObject(model,variable,way,lowerValue)
1343{
1344  changeInGuessed_=1.0e100;
1345}
1346 
1347
1348// Copy constructor
1349CbcIntegerPseudoCostBranchingObject::CbcIntegerPseudoCostBranchingObject ( 
1350                                 const CbcIntegerPseudoCostBranchingObject & rhs)
1351  :CbcIntegerBranchingObject(rhs)
1352{
1353  changeInGuessed_ = rhs.changeInGuessed_;
1354}
1355
1356// Assignment operator
1357CbcIntegerPseudoCostBranchingObject & 
1358CbcIntegerPseudoCostBranchingObject::operator=( const CbcIntegerPseudoCostBranchingObject& rhs)
1359{
1360  if (this != &rhs) {
1361    CbcIntegerBranchingObject::operator=(rhs);
1362    changeInGuessed_ = rhs.changeInGuessed_;
1363  }
1364  return *this;
1365}
1366CbcBranchingObject * 
1367CbcIntegerPseudoCostBranchingObject::clone() const
1368{ 
1369  return (new CbcIntegerPseudoCostBranchingObject(*this));
1370}
1371
1372
1373// Destructor
1374CbcIntegerPseudoCostBranchingObject::~CbcIntegerPseudoCostBranchingObject ()
1375{
1376}
1377
1378/*
1379  Perform a branch by adjusting the bounds of the specified variable. Note
1380  that each arm of the branch advances the object to the next arm by
1381  advancing the value of way_.
1382
1383  Providing new values for the variable's lower and upper bounds for each
1384  branching direction gives a little bit of additional flexibility and will
1385  be easily extensible to multi-way branching.
1386  Returns change in guessed objective on next branch
1387*/
1388double
1389CbcIntegerPseudoCostBranchingObject::branch()
1390{
1391  CbcIntegerBranchingObject::branch();
1392  return changeInGuessed_;
1393}
1394
1395
1396// Default Constructor
1397CbcCliqueBranchingObject::CbcCliqueBranchingObject()
1398  :CbcBranchingObject()
1399{
1400  clique_ = NULL;
1401  downMask_[0]=0;
1402  downMask_[1]=0;
1403  upMask_[0]=0;
1404  upMask_[1]=0;
1405}
1406
1407// Useful constructor
1408CbcCliqueBranchingObject::CbcCliqueBranchingObject (CbcModel * model,
1409                                                    const CbcClique * clique,
1410                                                    int way ,
1411                                                    int numberOnDownSide, const int * down,
1412                                                    int numberOnUpSide, const int * up)
1413  :CbcBranchingObject(model,clique->id(),way,0.5)
1414{
1415  clique_ = clique;
1416  downMask_[0]=0;
1417  downMask_[1]=0;
1418  upMask_[0]=0;
1419  upMask_[1]=0;
1420  int i;
1421  for (i=0;i<numberOnDownSide;i++) {
1422    int sequence = down[i];
1423    int iWord = sequence>>5;
1424    int iBit = sequence - 32*iWord;
1425    unsigned int k = 1<<iBit;
1426    downMask_[iWord] |= k;
1427  }
1428  for (i=0;i<numberOnUpSide;i++) {
1429    int sequence = up[i];
1430    int iWord = sequence>>5;
1431    int iBit = sequence - 32*iWord;
1432    unsigned int k = 1<<iBit;
1433    upMask_[iWord] |= k;
1434  }
1435}
1436
1437// Copy constructor
1438CbcCliqueBranchingObject::CbcCliqueBranchingObject ( const CbcCliqueBranchingObject & rhs) :CbcBranchingObject(rhs)
1439{
1440  clique_=rhs.clique_;
1441  downMask_[0]=rhs.downMask_[0];
1442  downMask_[1]=rhs.downMask_[1];
1443  upMask_[0]=rhs.upMask_[0];
1444  upMask_[1]=rhs.upMask_[1];
1445}
1446
1447// Assignment operator
1448CbcCliqueBranchingObject & 
1449CbcCliqueBranchingObject::operator=( const CbcCliqueBranchingObject& rhs)
1450{
1451  if (this != &rhs) {
1452    CbcBranchingObject::operator=(rhs);
1453    clique_=rhs.clique_;
1454    downMask_[0]=rhs.downMask_[0];
1455    downMask_[1]=rhs.downMask_[1];
1456    upMask_[0]=rhs.upMask_[0];
1457    upMask_[1]=rhs.upMask_[1];
1458  }
1459  return *this;
1460}
1461CbcBranchingObject * 
1462CbcCliqueBranchingObject::clone() const
1463{ 
1464  return (new CbcCliqueBranchingObject(*this));
1465}
1466
1467
1468// Destructor
1469CbcCliqueBranchingObject::~CbcCliqueBranchingObject ()
1470{
1471}
1472double
1473CbcCliqueBranchingObject::branch()
1474{
1475  decrementNumberBranchesLeft();
1476  int iWord;
1477  int numberMembers = clique_->numberMembers();
1478  const int * which = clique_->members();
1479  const int * integerVariables = model_->integerVariable();
1480  int numberWords=(numberMembers+31)>>5;
1481  // *** for way - up means fix all those in down section
1482  if (way_<0) {
1483#ifdef FULL_PRINT
1484    printf("Down Fix ");
1485#endif
1486    for (iWord=0;iWord<numberWords;iWord++) {
1487      int i;
1488      for (i=0;i<32;i++) {
1489        unsigned int k = 1<<i;
1490        if ((upMask_[iWord]&k)!=0) {
1491          int iColumn = which[i+32*iWord];
1492#ifdef FULL_PRINT
1493          printf("%d ",i+32*iWord);
1494#endif
1495          // fix weak way
1496          if (clique_->type(i+32*iWord))
1497            model_->solver()->setColUpper(integerVariables[iColumn],0.0);
1498          else
1499            model_->solver()->setColLower(integerVariables[iColumn],1.0);
1500        }
1501      }
1502    }
1503    way_=1;       // Swap direction
1504  } else {
1505#ifdef FULL_PRINT
1506    printf("Up Fix ");
1507#endif
1508    for (iWord=0;iWord<numberWords;iWord++) {
1509      int i;
1510      for (i=0;i<32;i++) {
1511        unsigned int k = 1<<i;
1512        if ((downMask_[iWord]&k)!=0) {
1513          int iColumn = which[i+32*iWord];
1514#ifdef FULL_PRINT
1515          printf("%d ",i+32*iWord);
1516#endif
1517          // fix weak way
1518          if (clique_->type(i+32*iWord))
1519            model_->solver()->setColUpper(integerVariables[iColumn],0.0);
1520          else
1521            model_->solver()->setColLower(integerVariables[iColumn],1.0);
1522        }
1523      }
1524    }
1525    way_=-1;      // Swap direction
1526  }
1527#ifdef FULL_PRINT
1528  printf("\n");
1529#endif
1530  return 0.0;
1531}
1532// Print what would happen 
1533void
1534CbcCliqueBranchingObject::print()
1535{
1536  int iWord;
1537  int numberMembers = clique_->numberMembers();
1538  const int * which = clique_->members();
1539  const int * integerVariables = model_->integerVariable();
1540  int numberWords=(numberMembers+31)>>5;
1541  // *** for way - up means fix all those in down section
1542  if (way_<0) {
1543    printf("Clique - Down Fix ");
1544    for (iWord=0;iWord<numberWords;iWord++) {
1545      int i;
1546      for (i=0;i<32;i++) {
1547        unsigned int k = 1<<i;
1548        if ((upMask_[iWord]&k)!=0) {
1549          int iColumn = which[i+32*iWord];
1550          printf("%d ",integerVariables[iColumn]);
1551        }
1552      }
1553    }
1554  } else {
1555    printf("Clique - Up Fix ");
1556    for (iWord=0;iWord<numberWords;iWord++) {
1557      int i;
1558      for (i=0;i<32;i++) {
1559        unsigned int k = 1<<i;
1560        if ((downMask_[iWord]&k)!=0) {
1561          int iColumn = which[i+32*iWord];
1562          printf("%d ",integerVariables[iColumn]);
1563        }
1564      }
1565    }
1566  }
1567  printf("\n");
1568}
1569 
1570// Default Constructor
1571CbcLongCliqueBranchingObject::CbcLongCliqueBranchingObject()
1572  :CbcBranchingObject()
1573{
1574  clique_=NULL;
1575  downMask_=NULL;
1576  upMask_=NULL;
1577}
1578
1579// Useful constructor
1580CbcLongCliqueBranchingObject::CbcLongCliqueBranchingObject (CbcModel * model,
1581                                                            const CbcClique * clique, 
1582                                                            int way ,
1583                                                    int numberOnDownSide, const int * down,
1584                                                    int numberOnUpSide, const int * up)
1585  :CbcBranchingObject(model,clique->id(),way,0.5)
1586{
1587  clique_ = clique;
1588  int numberMembers = clique_->numberMembers();
1589  int numberWords=(numberMembers+31)>>5;
1590  downMask_ = new unsigned int [numberWords];
1591  upMask_ = new unsigned int [numberWords];
1592  memset(downMask_,0,numberWords*sizeof(unsigned int));
1593  memset(upMask_,0,numberWords*sizeof(unsigned int));
1594  int i;
1595  for (i=0;i<numberOnDownSide;i++) {
1596    int sequence = down[i];
1597    int iWord = sequence>>5;
1598    int iBit = sequence - 32*iWord;
1599    unsigned int k = 1<<iBit;
1600    downMask_[iWord] |= k;
1601  }
1602  for (i=0;i<numberOnUpSide;i++) {
1603    int sequence = up[i];
1604    int iWord = sequence>>5;
1605    int iBit = sequence - 32*iWord;
1606    unsigned int k = 1<<iBit;
1607    upMask_[iWord] |= k;
1608  }
1609}
1610
1611// Copy constructor
1612CbcLongCliqueBranchingObject::CbcLongCliqueBranchingObject ( const CbcLongCliqueBranchingObject & rhs) :CbcBranchingObject(rhs)
1613{
1614  clique_=rhs.clique_;
1615  if (rhs.downMask_) {
1616    int numberMembers = clique_->numberMembers();
1617    int numberWords=(numberMembers+31)>>5;
1618    downMask_ = new unsigned int [numberWords];
1619    memcpy(downMask_,rhs.downMask_,numberWords*sizeof(unsigned int));
1620    upMask_ = new unsigned int [numberWords];
1621    memcpy(upMask_,rhs.upMask_,numberWords*sizeof(unsigned int));
1622  } else {
1623    downMask_=NULL;
1624    upMask_=NULL;
1625  }   
1626}
1627
1628// Assignment operator
1629CbcLongCliqueBranchingObject & 
1630CbcLongCliqueBranchingObject::operator=( const CbcLongCliqueBranchingObject& rhs)
1631{
1632  if (this != &rhs) {
1633    CbcBranchingObject::operator=(rhs);
1634    clique_=rhs.clique_;
1635    delete [] downMask_;
1636    delete [] upMask_;
1637    if (rhs.downMask_) {
1638      int numberMembers = clique_->numberMembers();
1639      int numberWords=(numberMembers+31)>>5;
1640      downMask_ = new unsigned int [numberWords];
1641      memcpy(downMask_,rhs.downMask_,numberWords*sizeof(unsigned int));
1642      upMask_ = new unsigned int [numberWords];
1643      memcpy(upMask_,rhs.upMask_,numberWords*sizeof(unsigned int));
1644    } else {
1645      downMask_=NULL;
1646      upMask_=NULL;
1647    }   
1648  }
1649  return *this;
1650}
1651CbcBranchingObject * 
1652CbcLongCliqueBranchingObject::clone() const
1653{ 
1654  return (new CbcLongCliqueBranchingObject(*this));
1655}
1656
1657
1658// Destructor
1659CbcLongCliqueBranchingObject::~CbcLongCliqueBranchingObject ()
1660{
1661  delete [] downMask_;
1662  delete [] upMask_;
1663}
1664double
1665CbcLongCliqueBranchingObject::branch()
1666{
1667  decrementNumberBranchesLeft();
1668  int iWord;
1669  int numberMembers = clique_->numberMembers();
1670  const int * which = clique_->members();
1671  const int * integerVariables = model_->integerVariable();
1672  int numberWords=(numberMembers+31)>>5;
1673  // *** for way - up means fix all those in down section
1674  if (way_<0) {
1675#ifdef FULL_PRINT
1676    printf("Down Fix ");
1677#endif
1678    for (iWord=0;iWord<numberWords;iWord++) {
1679      int i;
1680      for (i=0;i<32;i++) {
1681        unsigned int k = 1<<i;
1682        if ((upMask_[iWord]&k)!=0) {
1683          int iColumn = which[i+32*iWord];
1684#ifdef FULL_PRINT
1685          printf("%d ",i+32*iWord);
1686#endif
1687          // fix weak way
1688          if (clique_->type(i+32*iWord))
1689            model_->solver()->setColUpper(integerVariables[iColumn],0.0);
1690          else
1691            model_->solver()->setColLower(integerVariables[iColumn],1.0);
1692        }
1693      }
1694    }
1695    way_=1;       // Swap direction
1696  } else {
1697#ifdef FULL_PRINT
1698    printf("Up Fix ");
1699#endif
1700    for (iWord=0;iWord<numberWords;iWord++) {
1701      int i;
1702      for (i=0;i<32;i++) {
1703        unsigned int k = 1<<i;
1704        if ((downMask_[iWord]&k)!=0) {
1705          int iColumn = which[i+32*iWord];
1706#ifdef FULL_PRINT
1707          printf("%d ",i+32*iWord);
1708#endif
1709          // fix weak way
1710          if (clique_->type(i+32*iWord))
1711            model_->solver()->setColUpper(integerVariables[iColumn],0.0);
1712          else
1713            model_->solver()->setColLower(integerVariables[iColumn],1.0);
1714        }
1715      }
1716    }
1717    way_=-1;      // Swap direction
1718  }
1719#ifdef FULL_PRINT
1720  printf("\n");
1721#endif
1722  return 0.0;
1723}
1724void
1725CbcLongCliqueBranchingObject::print()
1726{
1727  int iWord;
1728  int numberMembers = clique_->numberMembers();
1729  const int * which = clique_->members();
1730  const int * integerVariables = model_->integerVariable();
1731  int numberWords=(numberMembers+31)>>5;
1732  // *** for way - up means fix all those in down section
1733  if (way_<0) {
1734    printf("Clique - Down Fix ");
1735    for (iWord=0;iWord<numberWords;iWord++) {
1736      int i;
1737      for (i=0;i<32;i++) {
1738        unsigned int k = 1<<i;
1739        if ((upMask_[iWord]&k)!=0) {
1740          int iColumn = which[i+32*iWord];
1741          printf("%d ",integerVariables[iColumn]);
1742        }
1743      }
1744    }
1745  } else {
1746    printf("Clique - Up Fix ");
1747    for (iWord=0;iWord<numberWords;iWord++) {
1748      int i;
1749      for (i=0;i<32;i++) {
1750        unsigned int k = 1<<i;
1751        if ((downMask_[iWord]&k)!=0) {
1752          int iColumn = which[i+32*iWord];
1753          printf("%d ",integerVariables[iColumn]);
1754        }
1755      }
1756    }
1757  }
1758  printf("\n");
1759}
1760// Default Constructor
1761CbcSOSBranchingObject::CbcSOSBranchingObject()
1762  :CbcBranchingObject()
1763{
1764  set_ = NULL;
1765  separator_=0.0;
1766}
1767
1768// Useful constructor
1769CbcSOSBranchingObject::CbcSOSBranchingObject (CbcModel * model,
1770                                              const CbcSOS * set,
1771                                              int way ,
1772                                              double separator)
1773  :CbcBranchingObject(model,set->id(),way,0.5)
1774{
1775  set_ = set;
1776  separator_ = separator;
1777}
1778
1779// Copy constructor
1780CbcSOSBranchingObject::CbcSOSBranchingObject ( const CbcSOSBranchingObject & rhs) :CbcBranchingObject(rhs)
1781{
1782  set_=rhs.set_;
1783  separator_ = rhs.separator_;
1784}
1785
1786// Assignment operator
1787CbcSOSBranchingObject & 
1788CbcSOSBranchingObject::operator=( const CbcSOSBranchingObject& rhs)
1789{
1790  if (this != &rhs) {
1791    CbcBranchingObject::operator=(rhs);
1792    set_=rhs.set_;
1793    separator_ = rhs.separator_;
1794  }
1795  return *this;
1796}
1797CbcBranchingObject * 
1798CbcSOSBranchingObject::clone() const
1799{ 
1800  return (new CbcSOSBranchingObject(*this));
1801}
1802
1803
1804// Destructor
1805CbcSOSBranchingObject::~CbcSOSBranchingObject ()
1806{
1807}
1808double
1809CbcSOSBranchingObject::branch()
1810{
1811  decrementNumberBranchesLeft();
1812  int numberMembers = set_->numberMembers();
1813  const int * which = set_->members();
1814  const double * weights = set_->weights();
1815  OsiSolverInterface * solver = model_->solver();
1816  //const double * lower = solver->getColLower();
1817  //const double * upper = solver->getColUpper();
1818  // *** for way - up means fix all those in down section
1819  if (way_<0) {
1820    int i;
1821    for ( i=0;i<numberMembers;i++) {
1822      if (weights[i] > separator_)
1823        break;
1824    }
1825    assert (i<numberMembers);
1826    for (;i<numberMembers;i++) 
1827      solver->setColUpper(which[i],0.0);
1828    way_=1;       // Swap direction
1829  } else {
1830    int i;
1831    for ( i=0;i<numberMembers;i++) {
1832      if (weights[i] >= separator_)
1833        break;
1834      else
1835        solver->setColUpper(which[i],0.0);
1836    }
1837    assert (i<numberMembers);
1838    way_=-1;      // Swap direction
1839  }
1840  return 0.0;
1841}
1842// Print what would happen 
1843void
1844CbcSOSBranchingObject::print()
1845{
1846  int numberMembers = set_->numberMembers();
1847  const int * which = set_->members();
1848  const double * weights = set_->weights();
1849  OsiSolverInterface * solver = model_->solver();
1850  //const double * lower = solver->getColLower();
1851  const double * upper = solver->getColUpper();
1852  int first=numberMembers;
1853  int last=-1;
1854  int numberFixed=0;
1855  int numberOther=0;
1856  int i;
1857  for ( i=0;i<numberMembers;i++) {
1858    double bound = upper[which[i]];
1859    if (bound) {
1860      first = CoinMin(first,i);
1861      last = CoinMax(last,i);
1862    }
1863  }
1864  // *** for way - up means fix all those in down section
1865  if (way_<0) {
1866    printf("SOS Down");
1867    for ( i=0;i<numberMembers;i++) {
1868      double bound = upper[which[i]];
1869      if (weights[i] > separator_)
1870        break;
1871      else if (bound)
1872        numberOther++;
1873    }
1874    assert (i<numberMembers);
1875    for (;i<numberMembers;i++) {
1876      double bound = upper[which[i]];
1877      if (bound)
1878        numberFixed++;
1879    }
1880  } else {
1881    printf("SOS Up");
1882    for ( i=0;i<numberMembers;i++) {
1883      double bound = upper[which[i]];
1884      if (weights[i] >= separator_)
1885        break;
1886      else if (bound)
1887        numberFixed++;
1888    }
1889    assert (i<numberMembers);
1890    for (;i<numberMembers;i++) {
1891      double bound = upper[which[i]];
1892      if (bound)
1893        numberOther++;
1894    }
1895  }
1896  printf(" - at %g, free range %d (%g) => %d (%g), %d would be fixed, %d other way\n",
1897         separator_,which[first],weights[first],which[last],weights[last],numberFixed,numberOther);
1898}
1899 
1900// Default Constructor
1901CbcBranchDefaultDecision::CbcBranchDefaultDecision()
1902  :CbcBranchDecision()
1903{
1904  bestCriterion_ = 0.0;
1905  bestChangeUp_ = 0.0;
1906  bestNumberUp_ = 0;
1907  bestChangeDown_ = 0.0;
1908  bestObject_ = NULL;
1909  bestNumberDown_ = 0;
1910}
1911
1912// Copy constructor
1913CbcBranchDefaultDecision::CbcBranchDefaultDecision (
1914                                    const CbcBranchDefaultDecision & rhs)
1915  :CbcBranchDecision(rhs)
1916{
1917  bestCriterion_ = rhs.bestCriterion_;
1918  bestChangeUp_ = rhs.bestChangeUp_;
1919  bestNumberUp_ = rhs.bestNumberUp_;
1920  bestChangeDown_ = rhs.bestChangeDown_;
1921  bestNumberDown_ = rhs.bestNumberDown_;
1922  bestObject_ = rhs.bestObject_;
1923  model_ = rhs.model_;
1924}
1925
1926CbcBranchDefaultDecision::~CbcBranchDefaultDecision()
1927{
1928}
1929
1930// Clone
1931CbcBranchDecision * 
1932CbcBranchDefaultDecision::clone() const
1933{
1934  return new CbcBranchDefaultDecision(*this);
1935}
1936
1937// Initialize i.e. before start of choosing at a node
1938void 
1939CbcBranchDefaultDecision::initialize(CbcModel * model)
1940{
1941  bestCriterion_ = 0.0;
1942  bestChangeUp_ = 0.0;
1943  bestNumberUp_ = 0;
1944  bestChangeDown_ = 0.0;
1945  bestNumberDown_ = 0;
1946  bestObject_ = NULL;
1947  model_ = model;
1948}
1949
1950
1951/*
1952  Simple default decision algorithm. Compare based on infeasibility (numInfUp,
1953  numInfDn) until a solution is found by search, then switch to change in
1954  objective (changeUp, changeDn). Note that bestSoFar is remembered in
1955  bestObject_, so the parameter bestSoFar is unused.
1956*/
1957
1958int
1959CbcBranchDefaultDecision::betterBranch(CbcBranchingObject * thisOne,
1960                            CbcBranchingObject * bestSoFar,
1961                            double changeUp, int numInfUp,
1962                            double changeDn, int numInfDn)
1963{
1964  bool beforeSolution = cbcModel()->getSolutionCount()==
1965    cbcModel()->getNumberHeuristicSolutions();;
1966  int betterWay=0;
1967  if (beforeSolution) {
1968    if (!bestObject_) {
1969      bestNumberUp_=INT_MAX;
1970      bestNumberDown_=INT_MAX;
1971    }
1972    // before solution - choose smallest number
1973    // could add in depth as well
1974    int bestNumber = CoinMin(bestNumberUp_,bestNumberDown_);
1975    if (numInfUp<numInfDn) {
1976      if (numInfUp<bestNumber) {
1977        betterWay = 1;
1978      } else if (numInfUp==bestNumber) {
1979        if (changeUp<bestCriterion_)
1980          betterWay=1;
1981      }
1982    } else if (numInfUp>numInfDn) {
1983      if (numInfDn<bestNumber) {
1984        betterWay = -1;
1985      } else if (numInfDn==bestNumber) {
1986        if (changeDn<bestCriterion_)
1987          betterWay=-1;
1988      }
1989    } else {
1990      // up and down have same number
1991      bool better=false;
1992      if (numInfUp<bestNumber) {
1993        better=true;
1994      } else if (numInfUp==bestNumber) {
1995        if (min(changeUp,changeDn)<bestCriterion_)
1996          better=true;;
1997      }
1998      if (better) {
1999        // see which way
2000        if (changeUp<=changeDn)
2001          betterWay=1;
2002        else
2003          betterWay=-1;
2004      }
2005    }
2006  } else {
2007    if (!bestObject_) {
2008      bestCriterion_=-1.0;
2009    }
2010    // got a solution
2011    if (changeUp<=changeDn) {
2012      if (changeUp>bestCriterion_)
2013        betterWay=1;
2014    } else {
2015      if (changeDn>bestCriterion_)
2016        betterWay=-1;
2017    }
2018  }
2019  if (betterWay) {
2020    bestCriterion_ = CoinMin(changeUp,changeDn);
2021    bestChangeUp_ = changeUp;
2022    bestNumberUp_ = numInfUp;
2023    bestChangeDown_ = changeDn;
2024    bestNumberDown_ = numInfDn;
2025    bestObject_=thisOne;
2026    // See if user is overriding way
2027    if (thisOne->object()&&thisOne->object()->preferredWay())
2028      betterWay = thisOne->object()->preferredWay();
2029  }
2030  return betterWay;
2031}
2032/* Sets or gets best criterion so far */
2033void 
2034CbcBranchDefaultDecision::setBestCriterion(double value)
2035{ 
2036  bestCriterion_ = value;
2037}
2038double 
2039CbcBranchDefaultDecision::getBestCriterion() const
2040{ 
2041  return bestCriterion_;
2042}
2043
2044/* Compare N branching objects. Return index of best
2045   and sets way of branching in chosen object.
2046   
2047   This routine is used only after strong branching.
2048*/
2049
2050int
2051CbcBranchDefaultDecision::bestBranch (CbcBranchingObject ** objects, int numberObjects,
2052                                   int numberUnsatisfied,
2053                                   double * changeUp, int * numberInfeasibilitiesUp,
2054                                   double * changeDown, int * numberInfeasibilitiesDown,
2055                                   double objectiveValue) 
2056{
2057
2058  int bestWay=0;
2059  int whichObject = -1;
2060  if (numberObjects) {
2061    CbcModel * model = cbcModel();
2062    // at continuous
2063    //double continuousObjective = model->getContinuousObjective();
2064    //int continuousInfeasibilities = model->getContinuousInfeasibilities();
2065   
2066    // average cost to get rid of infeasibility
2067    //double averageCostPerInfeasibility =
2068    //(objectiveValue-continuousObjective)/
2069    //(double) (abs(continuousInfeasibilities-numberUnsatisfied)+1);
2070    /* beforeSolution is :
2071       0 - before any solution
2072       n - n heuristic solutions but no branched one
2073       -1 - branched solution found
2074    */
2075    int numberSolutions = model->getSolutionCount();
2076    double cutoff = model->getCutoff();
2077    int method=0;
2078    int i;
2079    if (numberSolutions) {
2080      int numberHeuristic = model->getNumberHeuristicSolutions();
2081      if (numberHeuristic<numberSolutions) {
2082        method = 1;
2083      } else {
2084        method = 2;
2085        // look further
2086        for ( i = 0 ; i < numberObjects ; i++) {
2087          int numberNext = numberInfeasibilitiesUp[i];
2088         
2089          if (numberNext<numberUnsatisfied) {
2090            int numberUp = numberUnsatisfied - numberInfeasibilitiesUp[i];
2091            double perUnsatisfied = changeUp[i]/(double) numberUp;
2092            double estimatedObjective = objectiveValue + numberUnsatisfied * perUnsatisfied;
2093            if (estimatedObjective<cutoff) 
2094              method=3;
2095          }
2096          numberNext = numberInfeasibilitiesDown[i];
2097          if (numberNext<numberUnsatisfied) {
2098            int numberDown = numberUnsatisfied - numberInfeasibilitiesDown[i];
2099            double perUnsatisfied = changeDown[i]/(double) numberDown;
2100            double estimatedObjective = objectiveValue + numberUnsatisfied * perUnsatisfied;
2101            if (estimatedObjective<cutoff) 
2102              method=3;
2103          }
2104        }
2105      }
2106      method=2;
2107    } else {
2108      method = 0;
2109    }
2110    // Uncomment next to force method 4
2111    //method=4;
2112    /* Methods :
2113       0 - fewest infeasibilities
2114       1 - largest min change in objective
2115       2 - as 1 but use sum of changes if min close
2116       3 - predicted best solution
2117       4 - take cheapest up branch if infeasibilities same
2118    */
2119    int bestNumber=INT_MAX;
2120    double bestCriterion=-1.0e50;
2121    double alternativeCriterion = -1.0;
2122    double bestEstimate = 1.0e100;
2123    switch (method) {
2124    case 0:
2125      // could add in depth as well
2126      for ( i = 0 ; i < numberObjects ; i++) {
2127        int thisNumber = min(numberInfeasibilitiesUp[i],numberInfeasibilitiesDown[i]);
2128        if (thisNumber<=bestNumber) {
2129          int betterWay=0;
2130          if (numberInfeasibilitiesUp[i]<numberInfeasibilitiesDown[i]) {
2131            if (numberInfeasibilitiesUp[i]<bestNumber) {
2132              betterWay = 1;
2133            } else {
2134              if (changeUp[i]<bestCriterion)
2135                betterWay=1;
2136            }
2137          } else if (numberInfeasibilitiesUp[i]>numberInfeasibilitiesDown[i]) {
2138            if (numberInfeasibilitiesDown[i]<bestNumber) {
2139              betterWay = -1;
2140            } else {
2141              if (changeDown[i]<bestCriterion)
2142                betterWay=-1;
2143            }
2144          } else {
2145            // up and down have same number
2146            bool better=false;
2147            if (numberInfeasibilitiesUp[i]<bestNumber) {
2148              better=true;
2149            } else if (numberInfeasibilitiesUp[i]==bestNumber) {
2150              if (min(changeUp[i],changeDown[i])<bestCriterion)
2151                better=true;;
2152            }
2153            if (better) {
2154              // see which way
2155              if (changeUp[i]<=changeDown[i])
2156                betterWay=1;
2157              else
2158                betterWay=-1;
2159            }
2160          }
2161          if (betterWay) {
2162            bestCriterion = min(changeUp[i],changeDown[i]);
2163            bestNumber = thisNumber;
2164            whichObject = i;
2165            bestWay = betterWay;
2166          }
2167        }
2168      }
2169      break;
2170    case 1:
2171      for ( i = 0 ; i < numberObjects ; i++) {
2172        int betterWay=0;
2173        if (changeUp[i]<=changeDown[i]) {
2174          if (changeUp[i]>bestCriterion)
2175            betterWay=1;
2176        } else {
2177          if (changeDown[i]>bestCriterion)
2178            betterWay=-1;
2179        }
2180        if (betterWay) {
2181          bestCriterion = min(changeUp[i],changeDown[i]);
2182          whichObject = i;
2183          bestWay = betterWay;
2184        }
2185      }
2186      break;
2187    case 2:
2188      for ( i = 0 ; i < numberObjects ; i++) {
2189        double change = min(changeUp[i],changeDown[i]);
2190        double sum = changeUp[i] + changeDown[i];
2191        bool take=false;
2192        if (change>1.1*bestCriterion) 
2193          take=true;
2194        else if (change>0.9*bestCriterion&&sum+change>bestCriterion+alternativeCriterion) 
2195          take=true;
2196        if (take) {
2197          if (changeUp[i]<=changeDown[i]) {
2198            if (changeUp[i]>bestCriterion)
2199              bestWay=1;
2200          } else {
2201            if (changeDown[i]>bestCriterion)
2202              bestWay=-1;
2203          }
2204          bestCriterion = change;
2205          alternativeCriterion = sum;
2206          whichObject = i;
2207        }
2208      }
2209      break;
2210    case 3:
2211      for ( i = 0 ; i < numberObjects ; i++) {
2212        int numberNext = numberInfeasibilitiesUp[i];
2213       
2214        if (numberNext<numberUnsatisfied) {
2215          int numberUp = numberUnsatisfied - numberInfeasibilitiesUp[i];
2216          double perUnsatisfied = changeUp[i]/(double) numberUp;
2217          double estimatedObjective = objectiveValue + numberUnsatisfied * perUnsatisfied;
2218          if (estimatedObjective<bestEstimate) {
2219            bestEstimate = estimatedObjective;
2220            bestWay=1;
2221            whichObject=i;
2222          }
2223        }
2224        numberNext = numberInfeasibilitiesDown[i];
2225        if (numberNext<numberUnsatisfied) {
2226          int numberDown = numberUnsatisfied - numberInfeasibilitiesDown[i];
2227          double perUnsatisfied = changeDown[i]/(double) numberDown;
2228          double estimatedObjective = objectiveValue + numberUnsatisfied * perUnsatisfied;
2229          if (estimatedObjective<bestEstimate) {
2230            bestEstimate = estimatedObjective;
2231            bestWay=-1;
2232            whichObject=i;
2233          }
2234        }
2235      }
2236      break;
2237    case 4:
2238      // if number infeas same then cheapest up
2239      // first get best number or when going down
2240      // now choose smallest change up amongst equal number infeas
2241      for ( i = 0 ; i < numberObjects ; i++) {
2242        int thisNumber = min(numberInfeasibilitiesUp[i],numberInfeasibilitiesDown[i]);
2243        if (thisNumber<=bestNumber) {
2244          int betterWay=0;
2245          if (numberInfeasibilitiesUp[i]<numberInfeasibilitiesDown[i]) {
2246            if (numberInfeasibilitiesUp[i]<bestNumber) {
2247              betterWay = 1;
2248            } else {
2249              if (changeUp[i]<bestCriterion)
2250                betterWay=1;
2251            }
2252          } else if (numberInfeasibilitiesUp[i]>numberInfeasibilitiesDown[i]) {
2253            if (numberInfeasibilitiesDown[i]<bestNumber) {
2254              betterWay = -1;
2255            } else {
2256              if (changeDown[i]<bestCriterion)
2257                betterWay=-1;
2258            }
2259          } else {
2260            // up and down have same number
2261            bool better=false;
2262            if (numberInfeasibilitiesUp[i]<bestNumber) {
2263              better=true;
2264            } else if (numberInfeasibilitiesUp[i]==bestNumber) {
2265              if (min(changeUp[i],changeDown[i])<bestCriterion)
2266                better=true;;
2267            }
2268            if (better) {
2269              // see which way
2270              if (changeUp[i]<=changeDown[i])
2271                betterWay=1;
2272              else
2273                betterWay=-1;
2274            }
2275          }
2276          if (betterWay) {
2277            bestCriterion = min(changeUp[i],changeDown[i]);
2278            bestNumber = thisNumber;
2279            whichObject = i;
2280            bestWay = betterWay;
2281          }
2282        }
2283      }
2284      bestCriterion=1.0e50;
2285      for ( i = 0 ; i < numberObjects ; i++) {
2286        int thisNumber = numberInfeasibilitiesUp[i];
2287        if (thisNumber==bestNumber&&changeUp) {
2288          if (changeUp[i]<bestCriterion) {
2289            bestCriterion = changeUp[i];
2290            whichObject = i;
2291            bestWay = 1;
2292          }
2293        }
2294      }
2295      break;
2296    }
2297    // set way in best
2298    if (whichObject>=0) {
2299      CbcBranchingObject * bestObject = objects[whichObject];
2300      if (bestObject->object()&&bestObject->object()->preferredWay()) 
2301        bestWay = bestObject->object()->preferredWay();
2302      bestObject->way(bestWay);
2303    } else {
2304      printf("debug\n");
2305    }
2306  }
2307  return whichObject;
2308}
2309
2310// Default Constructor
2311CbcFollowOn::CbcFollowOn ()
2312  : CbcObject(),
2313    rhs_(NULL)
2314{
2315}
2316
2317// Useful constructor
2318CbcFollowOn::CbcFollowOn (CbcModel * model)
2319  : CbcObject(model)
2320{
2321  assert (model);
2322  OsiSolverInterface * solver = model_->solver();
2323  matrix_ = *solver->getMatrixByCol();
2324  matrix_.removeGaps();
2325  matrix_.setExtraGap(0.0);
2326  matrixByRow_ = *solver->getMatrixByRow();
2327  int numberRows = matrix_.getNumRows();
2328 
2329  rhs_ = new int[numberRows];
2330  int i;
2331  const double * rowLower = solver->getRowLower();
2332  const double * rowUpper = solver->getRowUpper();
2333  // Row copy
2334  const double * elementByRow = matrixByRow_.getElements();
2335  const int * column = matrixByRow_.getIndices();
2336  const CoinBigIndex * rowStart = matrixByRow_.getVectorStarts();
2337  const int * rowLength = matrixByRow_.getVectorLengths();
2338  for (i=0;i<numberRows;i++) {
2339    rhs_[i]=0;
2340    double value = rowLower[i];
2341    if (value==rowUpper[i]) {
2342      if (floor(value)==value&&value>=1.0&&value<10.0) {
2343        // check elements
2344        bool good=true;
2345        for (int j=rowStart[i];j<rowStart[i]+rowLength[i];j++) {
2346          int iColumn = column[j];
2347          if (!solver->isBinary(iColumn))
2348            good=false;
2349          double elValue = elementByRow[j];
2350          if (floor(elValue)!=elValue||value<1.0)
2351            good=false;
2352        }
2353        if (good)
2354          rhs_[i]=(int) value;
2355      }
2356    }
2357  }
2358}
2359
2360// Copy constructor
2361CbcFollowOn::CbcFollowOn ( const CbcFollowOn & rhs)
2362  :CbcObject(rhs),
2363   matrix_(rhs.matrix_),
2364   matrixByRow_(rhs.matrixByRow_)
2365{
2366  int numberRows = matrix_.getNumRows();
2367  rhs_= CoinCopyOfArray(rhs.rhs_,numberRows);
2368}
2369
2370// Clone
2371CbcObject *
2372CbcFollowOn::clone() const
2373{
2374  return new CbcFollowOn(*this);
2375}
2376
2377// Assignment operator
2378CbcFollowOn & 
2379CbcFollowOn::operator=( const CbcFollowOn& rhs)
2380{
2381  if (this!=&rhs) {
2382    CbcObject::operator=(rhs);
2383    delete [] rhs_;
2384    matrix_ = rhs.matrix_;
2385    matrixByRow_ = rhs.matrixByRow_;
2386    int numberRows = matrix_.getNumRows();
2387    rhs_= CoinCopyOfArray(rhs.rhs_,numberRows);
2388  }
2389  return *this;
2390}
2391
2392// Destructor
2393CbcFollowOn::~CbcFollowOn ()
2394{
2395  delete [] rhs_;
2396}
2397// As some computation is needed in more than one place - returns row
2398int 
2399CbcFollowOn::gutsOfFollowOn(int & otherRow, int & preferredWay) const
2400{
2401  int whichRow=-1;
2402  otherRow=-1;
2403  int numberRows = matrix_.getNumRows();
2404 
2405  int i;
2406  // For sorting
2407  int * sort = new int [numberRows];
2408  int * isort = new int [numberRows];
2409  // Column copy
2410  //const double * element = matrix_.getElements();
2411  const int * row = matrix_.getIndices();
2412  const CoinBigIndex * columnStart = matrix_.getVectorStarts();
2413  const int * columnLength = matrix_.getVectorLengths();
2414  // Row copy
2415  const double * elementByRow = matrixByRow_.getElements();
2416  const int * column = matrixByRow_.getIndices();
2417  const CoinBigIndex * rowStart = matrixByRow_.getVectorStarts();
2418  const int * rowLength = matrixByRow_.getVectorLengths();
2419  OsiSolverInterface * solver = model_->solver();
2420  const double * columnLower = solver->getColLower();
2421  const double * columnUpper = solver->getColUpper();
2422  const double * solution = solver->getColSolution();
2423  double integerTolerance = model_->getDblParam(CbcModel::CbcIntegerTolerance);
2424  int nSort=0;
2425  for (i=0;i<numberRows;i++) {
2426    if (rhs_[i]) {
2427      // check elements
2428      double smallest=1.0e10;
2429      double largest=0.0;
2430      int rhsValue=rhs_[i];
2431      int number1=0;
2432      int numberUnsatisfied=0;
2433      for (int j=rowStart[i];j<rowStart[i]+rowLength[i];j++) {
2434        int iColumn = column[j];
2435        double value = elementByRow[j];
2436        double solValue = solution[iColumn];
2437        if (columnLower[iColumn]!=columnUpper[iColumn]) {
2438          smallest = CoinMin(smallest,value);
2439          largest = CoinMax(largest,value);
2440          if (value==1.0)
2441            number1++;
2442          if (solValue<1.0-integerTolerance&&solValue>integerTolerance)
2443            numberUnsatisfied++;
2444        } else {
2445          rhsValue -= (int)(value*floor(solValue+0.5));
2446        }
2447      }
2448      if (numberUnsatisfied>1) {
2449        if (smallest<largest) {
2450          // probably no good but check a few things
2451          assert (largest<=rhsValue);
2452          if (number1==1&&largest==rhsValue)
2453            printf("could fix\n");
2454        } else if (largest==rhsValue) {
2455          sort[nSort]=i;
2456          isort[nSort++]=-numberUnsatisfied;
2457        }
2458      }
2459    }
2460  }
2461  if (nSort>1) {
2462    CoinSort_2(isort,isort+nSort,sort);
2463    CoinZeroN(isort,numberRows);
2464    double * other = new double[numberRows];
2465    CoinZeroN(other,numberRows);
2466    int * which = new int[numberRows];
2467    //#define COUNT
2468#ifndef COUNT
2469    bool beforeSolution = model_->getSolutionCount()==0;
2470#endif
2471    for (int k=0;k<nSort-1;k++) {
2472      i=sort[k];
2473      int numberUnsatisfied = 0;
2474      int n=0;
2475      int j;
2476      for (j=rowStart[i];j<rowStart[i]+rowLength[i];j++) {
2477        int iColumn = column[j];
2478        if (columnLower[iColumn]!=columnUpper[iColumn]) {
2479          double solValue = solution[iColumn]-columnLower[iColumn];
2480          if (solValue<1.0-integerTolerance&&solValue>integerTolerance) {
2481            numberUnsatisfied++;
2482            for (int jj=columnStart[iColumn];jj<columnStart[iColumn]+columnLength[iColumn];jj++) {
2483              int iRow = row[jj];
2484              if (rhs_[iRow]) {
2485                other[iRow]+=solValue;
2486                if (isort[iRow]) {
2487                  isort[iRow]++;
2488                } else {
2489                  isort[iRow]=1;
2490                  which[n++]=iRow;
2491                }
2492              }
2493            }
2494          }
2495        }
2496      }
2497      double total=0.0;
2498      // Take out row
2499      double sumThis=other[i];
2500      other[i]=0.0;
2501      assert (numberUnsatisfied==isort[i]);
2502      // find one nearest half if solution, one if before solution
2503      int iBest=-1;
2504      double dtarget=0.5*total;
2505#ifdef COUNT
2506      int target = (numberUnsatisfied+1)>>1;
2507      int best=numberUnsatisfied;
2508#else
2509      double best;
2510      if (beforeSolution)
2511        best=dtarget;
2512      else
2513        best=1.0e30;
2514#endif
2515      for (j=0;j<n;j++) {
2516        int iRow = which[j];
2517        double dvalue=other[iRow];
2518        other[iRow]=0.0;
2519#ifdef COUNT
2520        int value = isort[iRow];
2521#endif
2522        isort[iRow]=0;
2523        if (fabs(dvalue)<1.0e-8||fabs(sumThis-dvalue)<1.0e-8)
2524          continue;
2525        if (dvalue<integerTolerance||dvalue>1.0-integerTolerance)
2526          continue;
2527#ifdef COUNT
2528        if (abs(value-target)<best&&value!=numberUnsatisfied) {
2529          best=abs(value-target);
2530          iBest=iRow;
2531          if (dvalue<dtarget)
2532            preferredWay=1;
2533          else
2534            preferredWay=-1;
2535        }
2536#else
2537        if (beforeSolution) {
2538          if (fabs(dvalue-dtarget)>best) {
2539            best = fabs(dvalue-dtarget);
2540            iBest=iRow;
2541            if (dvalue<dtarget)
2542              preferredWay=1;
2543            else
2544              preferredWay=-1;
2545          }
2546        } else {
2547          if (fabs(dvalue-dtarget)<best) {
2548            best = fabs(dvalue-dtarget);
2549            iBest=iRow;
2550            if (dvalue<dtarget)
2551              preferredWay=1;
2552            else
2553              preferredWay=-1;
2554          }
2555        }
2556#endif
2557      }
2558      if (iBest>=0) {
2559        whichRow=i;
2560        otherRow=iBest;
2561        break;
2562      }
2563    }
2564    delete [] which;
2565    delete [] other;
2566  }
2567  delete [] sort;
2568  delete [] isort;
2569  return whichRow;
2570}
2571
2572// Infeasibility - large is 0.5
2573double 
2574CbcFollowOn::infeasibility(int & preferredWay) const
2575{
2576  int otherRow=0;
2577  int whichRow = gutsOfFollowOn(otherRow,preferredWay);
2578  if (whichRow<0)
2579    return 0.0;
2580  else
2581  return 2.0* model_->getDblParam(CbcModel::CbcIntegerTolerance);
2582}
2583
2584// This looks at solution and sets bounds to contain solution
2585void 
2586CbcFollowOn::feasibleRegion()
2587{
2588}
2589
2590
2591// Creates a branching object
2592CbcBranchingObject * 
2593CbcFollowOn::createBranch(int way) 
2594{
2595  int otherRow=0;
2596  int preferredWay;
2597  int whichRow = gutsOfFollowOn(otherRow,preferredWay);
2598  assert(way==preferredWay);
2599  assert (whichRow>=0);
2600  int numberColumns = matrix_.getNumCols();
2601 
2602  // Column copy
2603  //const double * element = matrix_.getElements();
2604  const int * row = matrix_.getIndices();
2605  const CoinBigIndex * columnStart = matrix_.getVectorStarts();
2606  const int * columnLength = matrix_.getVectorLengths();
2607  // Row copy
2608  //const double * elementByRow = matrixByRow_.getElements();
2609  const int * column = matrixByRow_.getIndices();
2610  const CoinBigIndex * rowStart = matrixByRow_.getVectorStarts();
2611  const int * rowLength = matrixByRow_.getVectorLengths();
2612  OsiSolverInterface * solver = model_->solver();
2613  const double * columnLower = solver->getColLower();
2614  const double * columnUpper = solver->getColUpper();
2615  //const double * solution = solver->getColSolution();
2616  int nUp=0;
2617  int nDown=0;
2618  int * upList = new int[numberColumns];
2619  int * downList = new int[numberColumns];
2620  int j;
2621  for (j=rowStart[whichRow];j<rowStart[whichRow]+rowLength[whichRow];j++) {
2622    int iColumn = column[j];
2623    if (columnLower[iColumn]!=columnUpper[iColumn]) {
2624      bool up=true;
2625      for (int jj=columnStart[iColumn];jj<columnStart[iColumn]+columnLength[iColumn];jj++) {
2626        int iRow = row[jj];
2627        if (iRow==otherRow) {
2628          up=false;
2629          break;
2630        }
2631      }
2632      if (up)
2633        upList[nUp++]=iColumn;
2634      else
2635        downList[nDown++]=iColumn;
2636    }
2637  }
2638  //printf("way %d\n",way);
2639  // create object
2640  //printf("would fix %d down and %d up\n",nDown,nUp);
2641  CbcBranchingObject * branch
2642     = new CbcFixingBranchingObject(model_,way,
2643                                         nDown,downList,nUp,upList);
2644  delete [] upList;
2645  delete [] downList;
2646  return branch;
2647}
2648// Default Constructor
2649CbcFixingBranchingObject::CbcFixingBranchingObject()
2650  :CbcBranchingObject()
2651{
2652  numberDown_=0;
2653  numberUp_=0;
2654  downList_=NULL;
2655  upList_=NULL;
2656}
2657
2658// Useful constructor
2659CbcFixingBranchingObject::CbcFixingBranchingObject (CbcModel * model,
2660                                                    int way ,
2661                                                    int numberOnDownSide, const int * down,
2662                                                    int numberOnUpSide, const int * up)
2663  :CbcBranchingObject(model,0,way,0.5)
2664{
2665  numberDown_=numberOnDownSide;
2666  numberUp_=numberOnUpSide;
2667  downList_ = CoinCopyOfArray(down,numberDown_);
2668  upList_ = CoinCopyOfArray(up,numberUp_);
2669}
2670
2671// Copy constructor
2672CbcFixingBranchingObject::CbcFixingBranchingObject ( const CbcFixingBranchingObject & rhs) :CbcBranchingObject(rhs)
2673{
2674  numberDown_=rhs.numberDown_;
2675  numberUp_=rhs.numberUp_;
2676  downList_ = CoinCopyOfArray(rhs.downList_,numberDown_);
2677  upList_ = CoinCopyOfArray(rhs.upList_,numberUp_);
2678}
2679
2680// Assignment operator
2681CbcFixingBranchingObject & 
2682CbcFixingBranchingObject::operator=( const CbcFixingBranchingObject& rhs)
2683{
2684  if (this != &rhs) {
2685    CbcBranchingObject::operator=(rhs);
2686    delete [] downList_;
2687    delete [] upList_;
2688    numberDown_=rhs.numberDown_;
2689    numberUp_=rhs.numberUp_;
2690    downList_ = CoinCopyOfArray(rhs.downList_,numberDown_);
2691    upList_ = CoinCopyOfArray(rhs.upList_,numberUp_);
2692  }
2693  return *this;
2694}
2695CbcBranchingObject * 
2696CbcFixingBranchingObject::clone() const
2697{ 
2698  return (new CbcFixingBranchingObject(*this));
2699}
2700
2701
2702// Destructor
2703CbcFixingBranchingObject::~CbcFixingBranchingObject ()
2704{
2705  delete [] downList_;
2706  delete [] upList_;
2707}
2708double
2709CbcFixingBranchingObject::branch()
2710{
2711  decrementNumberBranchesLeft();
2712  OsiSolverInterface * solver = model_->solver();
2713  const double * columnLower = solver->getColLower();
2714  int i;
2715  // *** for way - up means fix all those in up section
2716  if (way_<0) {
2717#ifdef FULL_PRINT
2718    printf("Down Fix ");
2719#endif
2720    //printf("Down Fix %d\n",numberDown_);
2721    for (i=0;i<numberDown_;i++) {
2722      int iColumn = downList_[i];
2723      model_->solver()->setColUpper(iColumn,columnLower[iColumn]);
2724#ifdef FULL_PRINT
2725      printf("Setting bound on %d to lower bound\n",iColumn);
2726#endif
2727    }
2728    way_=1;       // Swap direction
2729  } else {
2730#ifdef FULL_PRINT
2731    printf("Up Fix ");
2732#endif
2733    //printf("Up Fix %d\n",numberUp_);
2734    for (i=0;i<numberUp_;i++) {
2735      int iColumn = upList_[i];
2736      model_->solver()->setColUpper(iColumn,columnLower[iColumn]);
2737#ifdef FULL_PRINT
2738      printf("Setting bound on %d to lower bound\n",iColumn);
2739#endif
2740    }
2741    way_=-1;      // Swap direction
2742  }
2743#ifdef FULL_PRINT
2744  printf("\n");
2745#endif
2746  return 0.0;
2747}
2748void
2749CbcFixingBranchingObject::print()
2750{
2751  int i;
2752  // *** for way - up means fix all those in up section
2753  if (way_<0) {
2754    printf("Down Fix ");
2755    for (i=0;i<numberDown_;i++) {
2756      int iColumn = downList_[i];
2757      printf("%d ",iColumn);
2758    }
2759  } else {
2760    printf("Up Fix ");
2761    for (i=0;i<numberUp_;i++) {
2762      int iColumn = upList_[i];
2763      printf("%d ",iColumn);
2764    }
2765  }
2766  printf("\n");
2767}
2768// Default Constructor
2769CbcNWay::CbcNWay ()
2770  : CbcObject(),
2771    numberMembers_(0),
2772    members_(NULL),
2773    consequence_(NULL)
2774{
2775}
2776
2777// Useful constructor (which are integer indices)
2778CbcNWay::CbcNWay (CbcModel * model, int numberMembers,
2779                  const int * which, int identifier)
2780  : CbcObject(model)
2781{
2782  id_=identifier;
2783  numberMembers_=numberMembers;
2784  if (numberMembers_) {
2785    members_ = new int[numberMembers_];
2786    memcpy(members_,which,numberMembers_*sizeof(int));
2787  } else {
2788    members_ = NULL;
2789  }
2790  consequence_ = NULL;
2791}
2792
2793// Copy constructor
2794CbcNWay::CbcNWay ( const CbcNWay & rhs)
2795  :CbcObject(rhs)
2796{
2797  numberMembers_ = rhs.numberMembers_;
2798  consequence_ = NULL;
2799  if (numberMembers_) {
2800    members_ = new int[numberMembers_];
2801    memcpy(members_,rhs.members_,numberMembers_*sizeof(int));
2802    if (rhs.consequence_) {
2803      consequence_ = new CbcConsequence * [numberMembers_];
2804      for (int i=0;i<numberMembers_;i++) {
2805        if (rhs.consequence_[i])
2806          consequence_[i]= rhs.consequence_[i]->clone();
2807        else
2808          consequence_[i]=NULL;
2809      }
2810    }
2811  } else {
2812    members_ = NULL;
2813  }
2814}
2815
2816// Clone
2817CbcObject *
2818CbcNWay::clone() const
2819{
2820  return new CbcNWay(*this);
2821}
2822
2823// Assignment operator
2824CbcNWay & 
2825CbcNWay::operator=( const CbcNWay& rhs)
2826{
2827  if (this!=&rhs) {
2828    CbcObject::operator=(rhs);
2829    delete [] members_;
2830    numberMembers_ = rhs.numberMembers_;
2831    if (consequence_) {
2832      for (int i=0;i<numberMembers_;i++) 
2833        delete consequence_[i];
2834      delete [] consequence_;
2835      consequence_=NULL;
2836    }
2837    if (numberMembers_) {
2838      members_ = new int[numberMembers_];
2839      memcpy(members_,rhs.members_,numberMembers_*sizeof(int));
2840    } else {
2841      members_ = NULL;
2842    }
2843    if (rhs.consequence_) {
2844      consequence_ = new CbcConsequence * [numberMembers_];
2845      for (int i=0;i<numberMembers_;i++) {
2846        if (rhs.consequence_[i])
2847          consequence_[i]= rhs.consequence_[i]->clone();
2848        else
2849          consequence_[i]=NULL;
2850      }
2851    }
2852  }
2853  return *this;
2854}
2855
2856// Destructor
2857CbcNWay::~CbcNWay ()
2858{
2859  delete [] members_;
2860  if (consequence_) {
2861    for (int i=0;i<numberMembers_;i++) 
2862      delete consequence_[i];
2863    delete [] consequence_;
2864  }
2865}
2866// Set up a consequence for a single member
2867void 
2868CbcNWay::setConsequence(int iColumn, const CbcConsequence & consequence)
2869{
2870  if (!consequence_) {
2871    consequence_ = new CbcConsequence * [numberMembers_];
2872    for (int i=0;i<numberMembers_;i++) 
2873      consequence_[i]=NULL;
2874  }
2875  for (int i=0;i<numberMembers_;i++) {
2876    if (members_[i]==iColumn) {
2877      consequence_[i]=consequence.clone();
2878      break;
2879    }
2880  }
2881}
2882
2883// Applies a consequence for a single member
2884void 
2885CbcNWay::applyConsequence(int iSequence, int state) const
2886{
2887  assert (state==-9999||state==9999);
2888  if (consequence_) {
2889    CbcConsequence * consequence = consequence_[iSequence];
2890    if (consequence) 
2891      consequence->applyToSolver(model_->solver(),state);
2892  }
2893}
2894 
2895// Infeasibility - large is 0.5
2896double 
2897CbcNWay::infeasibility(int & preferredWay) const
2898{
2899  int numberUnsatis=0;
2900  int j;
2901  OsiSolverInterface * solver = model_->solver();
2902  const double * solution = model_->testSolution();
2903  const double * lower = solver->getColLower();
2904  const double * upper = solver->getColUpper();
2905  double largestValue=0.0;
2906 
2907  double integerTolerance = 
2908    model_->getDblParam(CbcModel::CbcIntegerTolerance);
2909
2910  for (j=0;j<numberMembers_;j++) {
2911    int iColumn = members_[j];
2912    double value = solution[iColumn];
2913    value = CoinMax(value, lower[iColumn]);
2914    value = CoinMin(value, upper[iColumn]);
2915    double distance = CoinMin(value-lower[iColumn],upper[iColumn]-value);
2916    if (distance>integerTolerance) {
2917      numberUnsatis++;
2918      largestValue = CoinMax(distance,largestValue);
2919    }
2920  }
2921  preferredWay=1;
2922  if (numberUnsatis) {
2923    return largestValue;
2924  } else {
2925    return 0.0; // satisfied
2926  }
2927}
2928
2929// This looks at solution and sets bounds to contain solution
2930void 
2931CbcNWay::feasibleRegion()
2932{
2933  int j;
2934  OsiSolverInterface * solver = model_->solver();
2935  const double * solution = model_->testSolution();
2936  const double * lower = solver->getColLower();
2937  const double * upper = solver->getColUpper();
2938  double integerTolerance = 
2939    model_->getDblParam(CbcModel::CbcIntegerTolerance);
2940  for (j=0;j<numberMembers_;j++) {
2941    int iColumn = members_[j];
2942    double value = solution[iColumn];
2943    value = CoinMax(value, lower[iColumn]);
2944    value = CoinMin(value, upper[iColumn]);
2945    if (value>=upper[iColumn]-integerTolerance) {
2946      solver->setColLower(iColumn,upper[iColumn]);
2947    } else {
2948      assert (value<=lower[iColumn]+integerTolerance);
2949      solver->setColUpper(iColumn,lower[iColumn]);
2950    }
2951  }
2952}
2953// Redoes data when sequence numbers change
2954void 
2955CbcNWay::redoSequenceEtc(CbcModel * model, int numberColumns, const int * originalColumns)
2956{
2957  model_=model;
2958  int n2=0;
2959  for (int j=0;j<numberMembers_;j++) {
2960    int iColumn = members_[j];
2961    int i;
2962    for (i=0;i<numberColumns;i++) {
2963      if (originalColumns[i]==iColumn)
2964        break;
2965    }
2966    if (i<numberColumns) {
2967      members_[n2]=i;
2968      consequence_[n2++]=consequence_[j];
2969    } else {
2970      delete consequence_[j];
2971    }
2972  }
2973  if (n2<numberMembers_) {
2974    printf("** NWay number of members reduced from %d to %d!\n",numberMembers_,n2);
2975    numberMembers_=n2;
2976  }
2977}
2978
2979
2980// Creates a branching object
2981CbcBranchingObject * 
2982CbcNWay::createBranch(int way) 
2983{
2984  int numberFree=0;
2985  int j;
2986
2987  OsiSolverInterface * solver = model_->solver();
2988  const double * solution = model_->testSolution();
2989  const double * lower = solver->getColLower();
2990  const double * upper = solver->getColUpper();
2991  int * list = new int[numberMembers_];
2992  double * sort = new double[numberMembers_];
2993
2994  for (j=0;j<numberMembers_;j++) {
2995    int iColumn = members_[j];
2996    double value = solution[iColumn];
2997    value = CoinMax(value, lower[iColumn]);
2998    value = CoinMin(value, upper[iColumn]);
2999    if (upper[iColumn]>lower[iColumn]) {
3000      double distance = upper[iColumn]-value;
3001      list[numberFree]=j;
3002      sort[numberFree++]=distance;
3003    }
3004  }
3005  assert (numberFree);
3006  // sort
3007  CoinSort_2(sort,sort+numberFree,list);
3008  // create object
3009  CbcBranchingObject * branch;
3010  branch = new CbcNWayBranchingObject(model_,this,numberFree,list);
3011  branch->setOriginalObject(this);
3012  delete [] list;
3013  delete [] sort;
3014  return branch;
3015}
3016 
3017// Default Constructor
3018CbcNWayBranchingObject::CbcNWayBranchingObject()
3019  :CbcBranchingObject()
3020{
3021  order_=NULL;
3022  object_=NULL;
3023  numberInSet_=0;
3024  way_=0;
3025}
3026
3027// Useful constructor
3028CbcNWayBranchingObject::CbcNWayBranchingObject (CbcModel * model,
3029                                                const CbcNWay * nway, 
3030                                                int number, const int * order)
3031  :CbcBranchingObject(model,nway->id(),-1,0.5)
3032{
3033  numberBranches_ = number;
3034  order_ = new int [number];
3035  object_=nway;
3036  numberInSet_=number;
3037  memcpy(order_,order,number*sizeof(int));
3038}
3039
3040// Copy constructor
3041CbcNWayBranchingObject::CbcNWayBranchingObject ( const CbcNWayBranchingObject & rhs) :CbcBranchingObject(rhs)
3042{
3043  numberInSet_=rhs.numberInSet_;
3044  object_=rhs.object_;
3045  if (numberInSet_) {
3046    order_ = new int [numberInSet_];
3047    memcpy(order_,rhs.order_,numberInSet_*sizeof(int));
3048  } else {
3049    order_=NULL;
3050  }   
3051}
3052
3053// Assignment operator
3054CbcNWayBranchingObject & 
3055CbcNWayBranchingObject::operator=( const CbcNWayBranchingObject& rhs)
3056{
3057  if (this != &rhs) {
3058    CbcBranchingObject::operator=(rhs);
3059    object_=rhs.object_;
3060    delete [] order_;
3061    numberInSet_=rhs.numberInSet_;
3062    if (numberInSet_) {
3063      order_ = new int [numberInSet_];
3064      memcpy(order_,rhs.order_,numberInSet_*sizeof(int));
3065    } else {
3066      order_=NULL;
3067    }   
3068  }
3069  return *this;
3070}
3071CbcBranchingObject * 
3072CbcNWayBranchingObject::clone() const
3073{ 
3074  return (new CbcNWayBranchingObject(*this));
3075}
3076
3077
3078// Destructor
3079CbcNWayBranchingObject::~CbcNWayBranchingObject ()
3080{
3081  delete [] order_;
3082}
3083double
3084CbcNWayBranchingObject::branch()
3085{
3086  int which = branchIndex_;
3087  branchIndex_++;
3088  assert (numberBranchesLeft()>=0);
3089  if (which==0) {
3090    // first branch so way_ may mean something
3091    assert (way_==-1||way_==1);
3092    if (way_==-1)
3093      which++;
3094  } else if (which==1) {
3095    // second branch so way_ may mean something
3096    assert (way_==-1||way_==1);
3097    if (way_==-1)
3098      which--;
3099    // switch way off
3100    way_=0;
3101  }
3102  const double * lower = model_->solver()->getColLower();
3103  const double * upper = model_->solver()->getColUpper();
3104  const int * members = object_->members();
3105  for (int j=0;j<numberInSet_;j++) {
3106    int iSequence = order_[j];
3107    int iColumn = members[iSequence];
3108    if (j!=which) {
3109      model_->solver()->setColUpper(iColumn,lower[iColumn]);
3110      //model_->solver()->setColLower(iColumn,lower[iColumn]);
3111      assert (lower[iColumn]>-1.0e20);
3112      // apply any consequences
3113      object_->applyConsequence(iSequence,-9999);
3114    } else {
3115      model_->solver()->setColLower(iColumn,upper[iColumn]);
3116      //model_->solver()->setColUpper(iColumn,upper[iColumn]);
3117#ifdef FULL_PRINT
3118      printf("Up Fix %d to %g\n",iColumn,upper[iColumn]);
3119#endif
3120      assert (upper[iColumn]<1.0e20);
3121      // apply any consequences
3122      object_->applyConsequence(iSequence,9999);
3123    }
3124  }
3125  return 0.0;
3126}
3127void
3128CbcNWayBranchingObject::print()
3129{
3130  printf("NWay - Up Fix ");
3131  const int * members = object_->members();
3132  for (int j=0;j<way_;j++) {
3133    int iColumn = members[order_[j]];
3134    printf("%d ",iColumn);
3135  }
3136  printf("\n");
3137}
3138
3139// Default Constructor
3140CbcFixVariable::CbcFixVariable ()
3141  : CbcConsequence(),
3142    numberStates_(0),
3143    states_(NULL),
3144    startLower_(NULL),
3145    startUpper_(NULL),
3146    newBound_(NULL),
3147    variable_(NULL)
3148{
3149}
3150
3151// One useful Constructor
3152CbcFixVariable::CbcFixVariable (int numberStates,const int * states, const int * numberNewLower, 
3153                                const int ** newLowerValue,
3154                                const int ** lowerColumn,
3155                                const int * numberNewUpper, const int ** newUpperValue,
3156                                const int ** upperColumn)
3157  : CbcConsequence(),
3158    states_(NULL),
3159    startLower_(NULL),
3160    startUpper_(NULL),
3161    newBound_(NULL),
3162    variable_(NULL)
3163{
3164  // How much space
3165  numberStates_ = numberStates;
3166  if (numberStates_) {
3167    states_ = new int[numberStates_];
3168    memcpy(states_,states,numberStates_*sizeof(int));
3169    int i;
3170    int n=0;
3171    startLower_ = new int[numberStates_+1];
3172    startUpper_ = new int[numberStates_+1];
3173    startLower_[0]=0;
3174    //count
3175    for (i=0;i<numberStates_;i++) {
3176      n += numberNewLower[i];
3177      startUpper_[i]=n;
3178      n += numberNewUpper[i];
3179      startLower_[i+1]=n;
3180    }
3181    newBound_ = new double [n];
3182    variable_ = new int [n];
3183    n=0;
3184    for (i=0;i<numberStates_;i++) {
3185      int j;
3186      int k;
3187      const int * bound;
3188      const int * variable;
3189      k=numberNewLower[i];
3190      bound = newLowerValue[i];
3191      variable = lowerColumn[i];
3192      for (j=0;j<k;j++) {
3193        newBound_[n]=bound[j];
3194        variable_[n++]=variable[j];
3195      }
3196      k=numberNewUpper[i];
3197      bound = newUpperValue[i];
3198      variable = upperColumn[i];
3199      for (j=0;j<k;j++) {
3200        newBound_[n]=bound[j];
3201        variable_[n++]=variable[j];
3202      }
3203    }
3204  }
3205}
3206
3207// Copy constructor
3208CbcFixVariable::CbcFixVariable ( const CbcFixVariable & rhs)
3209  :CbcConsequence(rhs)
3210{
3211  numberStates_ = rhs.numberStates_;
3212  states_ = NULL;
3213  startLower_ = NULL;
3214  startUpper_ = NULL;
3215  newBound_ = NULL;
3216  variable_ = NULL;
3217  if (numberStates_) {
3218    states_ = CoinCopyOfArray(rhs.states_,numberStates_);
3219    startLower_ = CoinCopyOfArray(rhs.startLower_,numberStates_+1);
3220    startUpper_ = CoinCopyOfArray(rhs.startUpper_,numberStates_+1);
3221    int n=startLower_[numberStates_];
3222    newBound_ = CoinCopyOfArray(rhs.newBound_,n);
3223    variable_ = CoinCopyOfArray(rhs.variable_,n);
3224  }
3225}
3226
3227// Clone
3228CbcConsequence *
3229CbcFixVariable::clone() const
3230{
3231  return new CbcFixVariable(*this);
3232}
3233
3234// Assignment operator
3235CbcFixVariable & 
3236CbcFixVariable::operator=( const CbcFixVariable& rhs)
3237{
3238  if (this!=&rhs) {
3239    CbcConsequence::operator=(rhs);
3240    delete [] states_;
3241    delete [] startLower_;
3242    delete [] startUpper_;
3243    delete [] newBound_;
3244    delete [] variable_;
3245    states_ = NULL;
3246    startLower_ = NULL;
3247    startUpper_ = NULL;
3248    newBound_ = NULL;
3249    variable_ = NULL;
3250    numberStates_ = rhs.numberStates_;
3251    if (numberStates_) {
3252      states_ = CoinCopyOfArray(rhs.states_,numberStates_);
3253      startLower_ = CoinCopyOfArray(rhs.startLower_,numberStates_+1);
3254      startUpper_ = CoinCopyOfArray(rhs.startUpper_,numberStates_+1);
3255      int n=startLower_[numberStates_];
3256      newBound_ = CoinCopyOfArray(rhs.newBound_,n);
3257      variable_ = CoinCopyOfArray(rhs.variable_,n);
3258    }
3259  }
3260  return *this;
3261}
3262
3263// Destructor
3264CbcFixVariable::~CbcFixVariable ()
3265{
3266  delete [] states_;
3267  delete [] startLower_;
3268  delete [] startUpper_;
3269  delete [] newBound_;
3270  delete [] variable_;
3271}
3272// Set up a startLower for a single member
3273void 
3274CbcFixVariable::applyToSolver(OsiSolverInterface * solver, int state) const
3275{
3276  assert (state==-9999||state==9999);
3277  // Find state
3278  int find;
3279  for (find=0;find<numberStates_;find++) 
3280    if (states_[find]==state)
3281      break;
3282  if (find==numberStates_)
3283    return;
3284  int i;
3285  // Set new lower bounds
3286  for (i=startLower_[find];i<startUpper_[find];i++) {
3287    int iColumn = variable_[i];
3288    double value = newBound_[i];
3289    double oldValue = solver->getColLower()[iColumn];
3290    //printf("for %d old lower bound %g, new %g",iColumn,oldValue,value);
3291    solver->setColLower(iColumn,CoinMax(value,oldValue));
3292    //printf(" => %g\n",solver->getColLower()[iColumn]);
3293  }
3294  // Set new upper bounds
3295  for (i=startUpper_[find];i<startLower_[find+1];i++) {
3296    int iColumn = variable_[i];
3297    double value = newBound_[i];
3298    double oldValue = solver->getColUpper()[iColumn];
3299    //printf("for %d old upper bound %g, new %g",iColumn,oldValue,value);
3300    solver->setColUpper(iColumn,CoinMin(value,oldValue));
3301    //printf(" => %g\n",solver->getColUpper()[iColumn]);
3302  }
3303}
3304
3305// Default Constructor
3306CbcDummyBranchingObject::CbcDummyBranchingObject(CbcModel * model)
3307  :CbcBranchingObject(model,0,0,0.5)
3308{
3309  setNumberBranchesLeft(1);
3310}
3311
3312
3313// Copy constructor
3314CbcDummyBranchingObject::CbcDummyBranchingObject ( const CbcDummyBranchingObject & rhs) :CbcBranchingObject(rhs)
3315{
3316}
3317
3318// Assignment operator
3319CbcDummyBranchingObject & 
3320CbcDummyBranchingObject::operator=( const CbcDummyBranchingObject& rhs)
3321{
3322  if (this != &rhs) {
3323    CbcBranchingObject::operator=(rhs);
3324  }
3325  return *this;
3326}
3327CbcBranchingObject * 
3328CbcDummyBranchingObject::clone() const
3329{ 
3330  return (new CbcDummyBranchingObject(*this));
3331}
3332
3333
3334// Destructor
3335CbcDummyBranchingObject::~CbcDummyBranchingObject ()
3336{
3337}
3338
3339/*
3340  Perform a dummy branch
3341*/
3342double
3343CbcDummyBranchingObject::branch()
3344{
3345  decrementNumberBranchesLeft();
3346  return 0.0;
3347}
3348// Print what would happen 
3349void
3350CbcDummyBranchingObject::print()
3351{
3352  printf("Dummy branch\n");
3353}
Note: See TracBrowser for help on using the repository browser.