source: trunk/Cbc/src/CbcNode.cpp @ 936

Last change on this file since 936 was 936, checked in by forrest, 11 years ago

maybe I was a bit too cavalier with numberStrong

  • Property svn:eol-style set to native
  • Property svn:keywords set to Author Date Id Revision
File size: 176.0 KB
Line 
1// Copyright (C) 2002, International Business Machines
2// Corporation and others.  All Rights Reserved.
3#if defined(_MSC_VER)
4// Turn off compiler warning about long names
5#  pragma warning(disable:4786)
6#endif
7
8#include "CbcConfig.h"
9
10#include <string>
11//#define CBC_DEBUG 1
12//#define CHECK_CUT_COUNTS
13//#define CHECK_NODE
14//#define CBC_CHECK_BASIS
15#define CBC_WEAK_STRONG
16#include <cassert>
17#include <cfloat>
18#define CUTS
19#include "OsiSolverInterface.hpp"
20#include "OsiChooseVariable.hpp"
21#include "OsiAuxInfo.hpp"
22#include "OsiSolverBranch.hpp"
23#include "CoinWarmStartBasis.hpp"
24#include "CoinTime.hpp"
25#include "CbcModel.hpp"
26#include "CbcNode.hpp"
27#include "CbcStatistics.hpp"
28#include "CbcStrategy.hpp"
29#include "CbcBranchActual.hpp"
30#include "CbcBranchDynamic.hpp"
31#include "OsiRowCut.hpp"
32#include "OsiRowCutDebugger.hpp"
33#include "OsiCuts.hpp"
34#include "CbcCountRowCut.hpp"
35#include "CbcFeasibilityBase.hpp"
36#include "CbcMessage.hpp"
37#ifdef COIN_HAS_CLP
38#include "OsiClpSolverInterface.hpp"
39#include "ClpSimplexOther.hpp"
40#endif
41using namespace std;
42#include "CglCutGenerator.hpp"
43// Default Constructor
44CbcNodeInfo::CbcNodeInfo ()
45  :
46  numberPointingToThis_(0),
47  parent_(NULL),
48  parentBranch_(NULL),
49  owner_(NULL),
50  numberCuts_(0),
51  nodeNumber_(0),
52  cuts_(NULL),
53  numberRows_(0),
54  numberBranchesLeft_(0),
55  active_(7)
56{
57#ifdef CHECK_NODE
58  printf("CbcNodeInfo %x Constructor\n",this);
59#endif
60}
61
62void
63CbcNodeInfo::setParentBasedData()
64{
65  if (parent_) {
66    numberRows_ = parent_->numberRows_+parent_->numberCuts_;
67    //parent_->increment();
68    if (parent_->owner()) {
69      const OsiBranchingObject* br = parent_->owner()->branchingObject();
70      assert(br);
71      parentBranch_ = br->clone();
72    }
73  }
74}
75
76#if 0
77// Constructor given parent
78CbcNodeInfo::CbcNodeInfo (CbcNodeInfo * parent)
79  :
80  numberPointingToThis_(2),
81  parent_(parent),
82  parentBranch_(NULL),
83  owner_(NULL),
84  numberCuts_(0),
85  nodeNumber_(0),
86  cuts_(NULL),
87  numberRows_(0),
88  numberBranchesLeft_(2),
89  active_(7)
90{
91#ifdef CHECK_NODE
92  printf("CbcNodeInfo %x Constructor from parent %x\n",this,parent_);
93#endif
94  setParentBasedData();
95}
96#endif
97
98// Copy Constructor
99CbcNodeInfo::CbcNodeInfo (const CbcNodeInfo & rhs)
100  :
101  numberPointingToThis_(rhs.numberPointingToThis_),
102  parent_(rhs.parent_),
103  parentBranch_(NULL),
104  owner_(rhs.owner_),
105  numberCuts_(rhs.numberCuts_),
106  nodeNumber_(rhs.nodeNumber_),
107  cuts_(NULL),
108  numberRows_(rhs.numberRows_),
109  numberBranchesLeft_(rhs.numberBranchesLeft_),
110  active_(rhs.active_)
111{
112#ifdef CHECK_NODE
113  printf("CbcNodeInfo %x Copy constructor\n",this);
114#endif
115  if (numberCuts_) {
116    cuts_ = new CbcCountRowCut * [numberCuts_];
117    int n=0;
118    for (int i=0;i<numberCuts_;i++) {
119      CbcCountRowCut * thisCut = rhs.cuts_[i];
120      if (thisCut) {
121        // I think this is correct - new one should take priority
122        thisCut->setInfo(this,n);
123        thisCut->increment(numberBranchesLeft_); 
124        cuts_[n++] = thisCut;
125      }
126    }
127    numberCuts_=n;
128  }
129  if (rhs.parentBranch_) {
130    parentBranch_ = rhs.parentBranch_->clone();
131  }
132}
133// Constructor given parent and owner
134CbcNodeInfo::CbcNodeInfo (CbcNodeInfo * parent, CbcNode * owner)
135  :
136  numberPointingToThis_(2),
137  parent_(parent),
138  parentBranch_(NULL),
139  owner_(owner),
140  numberCuts_(0),
141  nodeNumber_(0),
142  cuts_(NULL),
143  numberRows_(0),
144  numberBranchesLeft_(2),
145  active_(7)
146{
147#ifdef CHECK_NODE
148  printf("CbcNodeInfo %x Constructor from parent %x\n",this,parent_);
149#endif
150  setParentBasedData();
151}
152
153/**
154  Take care to detach from the owning CbcNode and decrement the reference
155  count in the parent.  If this is the last nodeInfo object pointing to the
156  parent, make a recursive call to delete the parent.
157*/
158CbcNodeInfo::~CbcNodeInfo()
159{
160#ifdef CHECK_NODE
161  printf("CbcNodeInfo %x Destructor parent %x\n",this,parent_);
162#endif
163
164  assert(!numberPointingToThis_);
165  // But there may be some left (max nodes?)
166  for (int i=0;i<numberCuts_;i++) {
167    if (cuts_[i]) {
168#ifndef GLOBAL_CUTS_JUST_POINTERS
169      delete cuts_[i];
170#else
171      if (cuts_[i]->globallyValidAsInteger()!=2)
172        delete cuts_[i];
173#endif
174    }
175  }
176  delete [] cuts_;
177  if (owner_) 
178    owner_->nullNodeInfo();
179  if (parent_) {
180    int numberLinks = parent_->decrement();
181    if (!numberLinks) delete parent_;
182  }
183  delete parentBranch_;
184}
185
186
187//#define ALLCUTS
188void
189CbcNodeInfo::decrementCuts(int change)
190{
191  int i;
192  // get rid of all remaining if negative
193  int changeThis;
194  if (change<0)
195    changeThis = numberBranchesLeft_;
196  else
197    changeThis = change;
198 // decrement cut counts
199  for (i=0;i<numberCuts_;i++) {
200    if (cuts_[i]) {
201      int number = cuts_[i]->decrement(changeThis);
202      if (!number) {
203        //printf("info %x del cut %d %x\n",this,i,cuts_[i]);
204#ifndef GLOBAL_CUTS_JUST_POINTERS
205        delete cuts_[i];
206#else
207        if (cuts_[i]->globallyValidAsInteger()!=2)
208          delete cuts_[i];
209#endif
210        cuts_[i]=NULL;
211      }
212    }
213  }
214}
215void
216CbcNodeInfo::incrementCuts(int change)
217{
218  int i;
219  assert (change>0);
220  // increment cut counts
221  for (i=0;i<numberCuts_;i++) {
222    if (cuts_[i]) 
223      cuts_[i]->increment(change);
224  }
225}
226void
227CbcNodeInfo::decrementParentCuts(CbcModel * model,int change)
228{
229  if (parent_) {
230    // get rid of all remaining if negative
231    int changeThis;
232    if (change<0)
233      changeThis = numberBranchesLeft_;
234    else
235      changeThis = change;
236    int i;
237    // Get over-estimate of space needed for basis
238    CoinWarmStartBasis & dummy = model->workingBasis();
239    dummy.setSize(0,numberRows_+numberCuts_);
240    buildRowBasis(dummy);
241    /* everything is zero (i.e. free) so we can use to see
242       if latest basis */
243    CbcNodeInfo * thisInfo = parent_;
244    while (thisInfo) 
245      thisInfo = thisInfo->buildRowBasis(dummy);
246    // decrement cut counts
247    thisInfo = parent_;
248    int numberRows=numberRows_;
249    while (thisInfo) {
250      for (i=thisInfo->numberCuts_-1;i>=0;i--) {
251        CoinWarmStartBasis::Status status = dummy.getArtifStatus(--numberRows);
252#ifdef ALLCUTS
253        status = CoinWarmStartBasis::isFree;
254#endif
255        if (thisInfo->cuts_[i]) {
256          int number=1;
257          if (status!=CoinWarmStartBasis::basic) {
258            // tight - drop 1 or 2
259            if (change<0)
260              number = thisInfo->cuts_[i]->decrement(changeThis);
261            else
262              number = thisInfo->cuts_[i]->decrement(change);
263          }
264          if (!number) {
265#ifndef GLOBAL_CUTS_JUST_POINTERS
266            delete thisInfo->cuts_[i];
267#else
268            if (thisInfo->cuts_[i]->globallyValidAsInteger()!=2)
269              delete thisInfo->cuts_[i];
270#endif
271            thisInfo->cuts_[i]=NULL;
272          }
273        }
274      }
275      thisInfo = thisInfo->parent_;
276    }
277  }
278}
279#if 0
280void
281CbcNodeInfo::incrementParentCuts(CbcModel * model, int change)
282{
283  if (parent_) {
284    int i;
285    // Get over-estimate of space needed for basis
286    CoinWarmStartBasis & dummy = model->workingBasis();
287    dummy.setSize(0,numberRows_+numberCuts_);
288    /* everything is zero (i.e. free) so we can use to see
289       if latest basis */
290    buildRowBasis(dummy);
291    CbcNodeInfo * thisInfo = parent_;
292    while (thisInfo)
293      thisInfo = thisInfo->buildRowBasis(dummy);
294    // increment cut counts
295    thisInfo = parent_;
296    int numberRows=numberRows_;
297    while (thisInfo) {
298      for (i=thisInfo->numberCuts_-1;i>=0;i--) {
299        CoinWarmStartBasis::Status status = dummy.getArtifStatus(--numberRows);
300#ifdef ALLCUTS
301        status = CoinWarmStartBasis::isFree;
302#endif
303        if (thisInfo->cuts_[i]&&status!=CoinWarmStartBasis::basic) {
304          thisInfo->cuts_[i]->increment(change);
305        }
306      }
307      thisInfo = thisInfo->parent_;
308    }
309  }
310}
311#endif
312/*
313  Append cuts to the cuts_ array in a nodeInfo. The initial reference count
314  is set to numberToBranchOn, which will normally be the number of arms
315  defined for the CbcBranchingObject attached to the CbcNode that owns this
316  CbcNodeInfo.
317*/
318void
319CbcNodeInfo::addCuts (OsiCuts & cuts, int numberToBranchOn,
320                      int * whichGenerator)
321{
322  int numberCuts = cuts.sizeRowCuts();
323  if (numberCuts) {
324    int i;
325    if (!numberCuts_) {
326      cuts_ = new CbcCountRowCut * [numberCuts];
327    } else {
328      CbcCountRowCut ** temp = new CbcCountRowCut * [numberCuts+numberCuts_];
329      memcpy(temp,cuts_,numberCuts_*sizeof(CbcCountRowCut *));
330      delete [] cuts_;
331      cuts_ = temp;
332    }
333    for (i=0;i<numberCuts;i++) {
334      CbcCountRowCut * thisCut = new CbcCountRowCut(*cuts.rowCutPtr(i),
335                                                    this,numberCuts_);
336      thisCut->increment(numberToBranchOn); 
337      cuts_[numberCuts_++] = thisCut;
338#ifdef CBC_DEBUG
339#if CBC_DEBUG>1
340      int n=thisCut->row().getNumElements();
341      printf("Cut %d has %d entries, rhs %g %g =>",i,n,thisCut->lb(),
342             thisCut->ub());
343      int j;
344      const int * index = thisCut->row().getIndices();
345      const double * element = thisCut->row().getElements();
346      for (j=0;j<n;j++) {
347        printf(" (%d,%g)",index[j],element[j]);
348        assert(fabs(element[j])>1.00e-12);
349      }
350      printf("\n");
351#else
352      int n=thisCut->row().getNumElements();
353      int j;
354      const double * element = thisCut->row().getElements();
355      for (j=0;j<n;j++) {
356        assert(fabs(element[j])>1.00e-12);
357      }
358#endif
359#endif
360    }
361  }
362}
363
364void
365CbcNodeInfo::addCuts(int numberCuts, CbcCountRowCut ** cut, 
366                     int numberToBranchOn)
367{
368  if (numberCuts) {
369    int i;
370    if (!numberCuts_) {
371      cuts_ = new CbcCountRowCut * [numberCuts];
372    } else {
373      CbcCountRowCut ** temp = new CbcCountRowCut * [numberCuts+numberCuts_];
374      memcpy(temp,cuts_,numberCuts_*sizeof(CbcCountRowCut *));
375      delete [] cuts_;
376      cuts_ = temp;
377    }
378    for (i=0;i<numberCuts;i++) {
379      CbcCountRowCut * thisCut = cut[i];
380      thisCut->setInfo(this,numberCuts_);
381      //printf("info %x cut %d %x\n",this,i,thisCut);
382      thisCut->increment(numberToBranchOn); 
383      cuts_[numberCuts_++] = thisCut;
384#ifdef CBC_DEBUG
385      int n=thisCut->row().getNumElements();
386#if CBC_DEBUG>1
387      printf("Cut %d has %d entries, rhs %g %g =>",i,n,thisCut->lb(),
388             thisCut->ub());
389#endif
390      int j;
391#if CBC_DEBUG>1
392      const int * index = thisCut->row().getIndices();
393#endif
394      const double * element = thisCut->row().getElements();
395      for (j=0;j<n;j++) {
396#if CBC_DEBUG>1
397        printf(" (%d,%g)",index[j],element[j]);
398#endif
399        assert(fabs(element[j])>1.00e-12);
400      }
401      printf("\n");
402#endif
403    }
404  }
405}
406
407// delete cuts
408void
409CbcNodeInfo::deleteCuts(int numberToDelete, CbcCountRowCut ** cuts)
410{
411  int i;
412  int j;
413  int last=-1;
414  for (i=0;i<numberToDelete;i++) {
415    CbcCountRowCut * next = cuts[i];
416    for (j=last+1;j<numberCuts_;j++) {
417      if (next==cuts_[j])
418        break;
419    }
420    if (j==numberCuts_) {
421      // start from beginning
422      for (j=0;j<last;j++) {
423        if (next==cuts_[j])
424          break;
425      }
426      assert(j<last);
427    }
428    last=j;
429    int number = cuts_[j]->decrement();
430    if (!number) {
431#ifndef GLOBAL_CUTS_JUST_POINTERS
432      delete cuts_[j];
433#else
434      if (cuts_[j]->globallyValidAsInteger()!=2)
435        delete cuts_[j];
436#endif
437    }
438    cuts_[j]=NULL;
439  }
440  j=0;
441  for (i=0;i<numberCuts_;i++) {
442    if (cuts_[i])
443      cuts_[j++]=cuts_[i];
444  }
445  numberCuts_ = j;
446}
447
448// delete cuts
449void
450CbcNodeInfo::deleteCuts(int numberToDelete, int * which)
451{
452  int i;
453  for (i=0;i<numberToDelete;i++) {
454    int iCut=which[i];
455    int number = cuts_[iCut]->decrement();
456    if (!number) {
457#ifndef GLOBAL_CUTS_JUST_POINTERS
458      delete cuts_[iCut];
459#else
460      if (cuts_[iCut]->globallyValidAsInteger()!=2)
461        delete cuts_[iCut];
462#endif
463    }
464    cuts_[iCut]=NULL;
465  }
466  int j=0;
467  for (i=0;i<numberCuts_;i++) {
468    if (cuts_[i])
469      cuts_[j++]=cuts_[i];
470  }
471  numberCuts_ = j;
472}
473
474// Really delete a cut
475void 
476CbcNodeInfo::deleteCut(int whichOne)
477{
478  assert(whichOne<numberCuts_);
479  cuts_[whichOne]=NULL;
480}
481/* Deactivate node information.
482   1 - bounds
483   2 - cuts
484   4 - basis!
485*/
486void 
487CbcNodeInfo::deactivate(int mode)
488{
489  active_ &= (~mode);
490}
491
492CbcFullNodeInfo::CbcFullNodeInfo() :
493  CbcNodeInfo(),
494  basis_(),
495  numberIntegers_(0),
496  lower_(NULL),
497  upper_(NULL)
498{
499}
500CbcFullNodeInfo::CbcFullNodeInfo(CbcModel * model,
501                                 int numberRowsAtContinuous) :
502  CbcNodeInfo(NULL, model->currentNode())
503{
504  OsiSolverInterface * solver = model->solver();
505  numberRows_ = numberRowsAtContinuous;
506  numberIntegers_ = model->numberIntegers();
507  int numberColumns = model->getNumCols();
508  lower_ = new double [numberColumns];
509  upper_ = new double [numberColumns];
510  const double * lower = solver->getColLower();
511  const double * upper = solver->getColUpper();
512  int i;
513
514  for (i=0;i<numberColumns;i++) {
515    lower_[i]=lower[i];
516    upper_[i]=upper[i];
517  }
518
519  basis_ =  dynamic_cast<CoinWarmStartBasis*>(solver->getWarmStart());
520}
521
522CbcFullNodeInfo::CbcFullNodeInfo(const CbcFullNodeInfo & rhs) :
523  CbcNodeInfo(rhs)
524{ 
525  basis_= dynamic_cast<CoinWarmStartBasis *>(rhs.basis_->clone()) ;
526  numberIntegers_=rhs.numberIntegers_;
527  lower_=NULL;
528  upper_=NULL;
529  if (rhs.lower_!=NULL) {
530    int numberColumns = basis_->getNumStructural();
531    lower_ = new double [numberColumns];
532    upper_ = new double [numberColumns];
533    assert (upper_!=NULL);
534    memcpy(lower_,rhs.lower_,numberColumns*sizeof(double));
535    memcpy(upper_,rhs.upper_,numberColumns*sizeof(double));
536  }
537}
538
539CbcNodeInfo * 
540CbcFullNodeInfo::clone() const
541{ 
542  return (new CbcFullNodeInfo(*this));
543}
544
545CbcFullNodeInfo::~CbcFullNodeInfo ()
546{
547  delete basis_ ;
548  delete [] lower_;
549  delete [] upper_;
550}
551
552/*
553   The basis supplied as a parameter is deleted and replaced with a new basis
554   appropriate for the node, and lower and upper bounds on variables are
555   reset according to the stored bounds arrays. Any cuts associated with this
556   node are added to the list in addCuts, but not actually added to the
557   constraint system in the model.
558
559   Why pass in a basis at all? The short answer is ``We need the parameter to
560   pass out a basis, so might as well use it to pass in the size.''
561   
562   A longer answer is that in practice we take a memory allocation hit up in
563   addCuts1 (the only place applyToModel is called) when we setSize() the
564   basis that's passed in. It's immediately tossed here in favour of a clone
565   of the basis attached to this nodeInfo. This can probably be fixed, given
566   a bit of thought.
567*/
568
569void CbcFullNodeInfo::applyToModel (CbcModel *model,
570                                    CoinWarmStartBasis *&basis,
571                                    CbcCountRowCut **addCuts,
572                                    int &currentNumberCuts) const 
573
574{ OsiSolverInterface *solver = model->solver() ;
575
576  // branch - do bounds
577  assert (active_==7||active_==15);
578  int i;
579  solver->setColLower(lower_);
580  solver->setColUpper(upper_);
581  int numberColumns = model->getNumCols();
582  // move basis - but make sure size stays
583  // for bon-min - should not be needed int numberRows = model->getNumRows();
584  int numberRows=basis->getNumArtificial();
585  delete basis ;
586  if (basis_) {
587    basis = dynamic_cast<CoinWarmStartBasis *>(basis_->clone()) ;
588    basis->resize(numberRows,numberColumns);
589#ifdef CBC_CHECK_BASIS
590    std::cout << "Basis (after applying root " <<this<<") "<< std::endl ;
591    basis->print() ;
592#endif   
593  } else {
594    // We have a solver without a basis
595    basis=NULL;
596  }
597  for (i=0;i<numberCuts_;i++) 
598    addCuts[currentNumberCuts+i]= cuts_[i];
599  currentNumberCuts += numberCuts_;
600  assert(!parent_);
601  return ;
602}
603// Just apply bounds to one variable (1=>infeasible)
604int 
605CbcFullNodeInfo::applyBounds(int iColumn, double & lower, double & upper,int force) 
606{
607  if ((force&&1)==0) {
608    if (lower>lower_[iColumn])
609      printf("%d odd lower going from %g to %g\n",iColumn,lower,lower_[iColumn]);
610    lower = lower_[iColumn];
611  } else {
612    lower_[iColumn]=lower;
613  }
614  if ((force&&2)==0) {
615    if (upper<upper_[iColumn])
616      printf("%d odd upper going from %g to %g\n",iColumn,upper,upper_[iColumn]);
617    upper = upper_[iColumn];
618  } else {
619    upper_[iColumn]=upper;
620  }
621  return (upper_[iColumn]>=lower_[iColumn]) ? 0 : 1;
622}
623
624/* Builds up row basis backwards (until original model).
625   Returns NULL or previous one to apply .
626   Depends on Free being 0 and impossible for cuts
627*/
628CbcNodeInfo * 
629CbcFullNodeInfo::buildRowBasis(CoinWarmStartBasis & basis ) const 
630{
631  const unsigned int * saved = 
632    (const unsigned int *) basis_->getArtificialStatus();
633  unsigned int * now = 
634    (unsigned int *) basis.getArtificialStatus();
635  int number=basis_->getNumArtificial()>>4;;
636  int i;
637  for (i=0;i<number;i++) { 
638    if (!now[i])
639      now[i] = saved[i];
640  }
641  return NULL;
642}
643
644
645// Default constructor
646CbcPartialNodeInfo::CbcPartialNodeInfo()
647
648  : CbcNodeInfo(),
649    basisDiff_(NULL),
650    variables_(NULL),
651    newBounds_(NULL),
652    numberChangedBounds_(0)
653
654{ /* this space intentionally left blank */ }
655
656// Constructor from current state
657CbcPartialNodeInfo::CbcPartialNodeInfo (CbcNodeInfo *parent, CbcNode *owner,
658                                        int numberChangedBounds,
659                                        const int *variables,
660                                        const double *boundChanges,
661                                        const CoinWarmStartDiff *basisDiff)
662 : CbcNodeInfo(parent,owner)
663{
664  basisDiff_ = basisDiff->clone() ;
665
666  numberChangedBounds_ = numberChangedBounds;
667  int size = numberChangedBounds_*(sizeof(double)+sizeof(int));
668  char * temp = new char [size];
669  newBounds_ = (double *) temp;
670  variables_ = (int *) (newBounds_+numberChangedBounds_);
671
672  int i ;
673  for (i=0;i<numberChangedBounds_;i++) {
674    variables_[i]=variables[i];
675    newBounds_[i]=boundChanges[i];
676  }
677}
678
679CbcPartialNodeInfo::CbcPartialNodeInfo (const CbcPartialNodeInfo & rhs)
680
681  : CbcNodeInfo(rhs)
682
683{ basisDiff_ = rhs.basisDiff_->clone() ;
684
685  numberChangedBounds_ = rhs.numberChangedBounds_;
686  int size = numberChangedBounds_*(sizeof(double)+sizeof(int));
687  char * temp = new char [size];
688  newBounds_ = (double *) temp;
689  variables_ = (int *) (newBounds_+numberChangedBounds_);
690
691  int i ;
692  for (i=0;i<numberChangedBounds_;i++) {
693    variables_[i]=rhs.variables_[i];
694    newBounds_[i]=rhs.newBounds_[i];
695  }
696}
697
698CbcNodeInfo * 
699CbcPartialNodeInfo::clone() const
700{ 
701  return (new CbcPartialNodeInfo(*this));
702}
703
704
705CbcPartialNodeInfo::~CbcPartialNodeInfo ()
706{
707  delete basisDiff_ ;
708  delete [] newBounds_;
709}
710
711
712/**
713   The basis supplied as a parameter is incrementally modified, and lower and
714   upper bounds on variables in the model are incrementally modified. Any
715   cuts associated with this node are added to the list in addCuts.
716*/
717
718void CbcPartialNodeInfo::applyToModel (CbcModel *model,
719                                       CoinWarmStartBasis *&basis,
720                                       CbcCountRowCut **addCuts,
721                                       int &currentNumberCuts) const 
722
723{ OsiSolverInterface *solver = model->solver();
724  if ((active_&4)!=0) {
725    basis->applyDiff(basisDiff_) ;
726#ifdef CBC_CHECK_BASIS
727    std::cout << "Basis (after applying " <<this<<") "<< std::endl ;
728    basis->print() ;
729#endif   
730  }
731
732  // branch - do bounds
733  int i;
734  if ((active_&1)!=0) {
735    for (i=0;i<numberChangedBounds_;i++) {
736      int variable = variables_[i];
737      int k = variable&0x3fffffff;
738      if ((variable&0x80000000)==0) {
739        // lower bound changing
740        //#define CBC_PRINT2
741#ifdef CBC_PRINT2
742        if(solver->getColLower()[k]!=newBounds_[i])
743          printf("lower change for column %d - from %g to %g\n",k,solver->getColLower()[k],newBounds_[i]);
744#endif
745#ifndef NDEBUG
746        if ((variable&0x40000000)==0&&false) {
747          double oldValue = solver->getColLower()[k];
748          assert (newBounds_[i]>oldValue-1.0e-8);
749          if (newBounds_[i]<oldValue+1.0e-8)
750            printf("bad null lower change for column %d - bound %g\n",k,oldValue);
751        }
752#endif
753        solver->setColLower(k,newBounds_[i]);
754      } else {
755        // upper bound changing
756#ifdef CBC_PRINT2
757        if(solver->getColUpper()[k]!=newBounds_[i])
758          printf("upper change for column %d - from %g to %g\n",k,solver->getColUpper()[k],newBounds_[i]);
759#endif
760#ifndef NDEBUG
761        if ((variable&0x40000000)==0&&false) {
762          double oldValue = solver->getColUpper()[k];
763          assert (newBounds_[i]<oldValue+1.0e-8);
764          if (newBounds_[i]>oldValue-1.0e-8)
765            printf("bad null upper change for column %d - bound %g\n",k,oldValue);
766        }
767#endif
768        solver->setColUpper(k,newBounds_[i]);
769      }
770    }
771  }
772  if ((active_&2)!=0) {
773    for (i=0;i<numberCuts_;i++) {
774      addCuts[currentNumberCuts+i]= cuts_[i];
775      if (cuts_[i]&&model->messageHandler()->logLevel()>4) {
776        cuts_[i]->print();
777      }
778    }
779   
780    currentNumberCuts += numberCuts_;
781  }
782  return ;
783}
784// Just apply bounds to one variable (1=>infeasible)
785int
786CbcPartialNodeInfo::applyBounds(int iColumn, double & lower, double & upper,int force) 
787{
788  // branch - do bounds
789  int i;
790  int found=0;
791  double newLower = -COIN_DBL_MAX;
792  double newUpper = COIN_DBL_MAX;
793  for (i=0;i<numberChangedBounds_;i++) {
794    int variable = variables_[i];
795    int k = variable&0x3fffffff;
796    if (k==iColumn) {
797      if ((variable&0x80000000)==0) {
798        // lower bound changing
799        found |= 1;
800        newLower = CoinMax(newLower,newBounds_[i]);
801        if ((force&1)==0) {
802          if (lower>newBounds_[i])
803            printf("%d odd lower going from %g to %g\n",iColumn,lower,newBounds_[i]);
804          lower = newBounds_[i];
805        } else {
806          newBounds_[i]=lower;
807          variables_[i] |= 0x40000000; // say can go odd way
808        }
809      } else {
810        // upper bound changing
811        found |= 2;
812        newUpper = CoinMin(newUpper,newBounds_[i]);
813        if ((force&2)==0) {
814          if (upper<newBounds_[i])
815            printf("%d odd upper going from %g to %g\n",iColumn,upper,newBounds_[i]);
816          upper = newBounds_[i];
817        } else {
818          newBounds_[i]=upper;
819          variables_[i] |= 0x40000000; // say can go odd way
820        }
821      }
822    }
823  }
824  newLower = CoinMax(newLower,lower);
825  newUpper = CoinMin(newUpper,upper);
826  int nAdd=0;
827  if ((force&2)!=0&&(found&2)==0) {
828    // need to add new upper
829    nAdd++;
830  }
831  if ((force&1)!=0&&(found&1)==0) {
832    // need to add new lower
833    nAdd++;
834  }
835  if (nAdd) { 
836    int size = (numberChangedBounds_+nAdd)*(sizeof(double)+sizeof(int));
837    char * temp = new char [size];
838    double * newBounds = (double *) temp;
839    int * variables = (int *) (newBounds+numberChangedBounds_+nAdd);
840
841    int i ;
842    for (i=0;i<numberChangedBounds_;i++) {
843      variables[i]=variables_[i];
844      newBounds[i]=newBounds_[i];
845    }
846    delete [] newBounds_;
847    newBounds_ = newBounds;
848    variables_ = variables;
849    if ((force&2)!=0&&(found&2)==0) {
850      // need to add new upper
851      int variable = iColumn | 0x80000000;
852      variables_[numberChangedBounds_]=variable;
853      newBounds_[numberChangedBounds_++]=newUpper;
854    }
855    if ((force&1)!=0&&(found&1)==0) {
856      // need to add new lower
857      int variable = iColumn;
858      variables_[numberChangedBounds_]=variable;
859      newBounds_[numberChangedBounds_++]=newLower;
860    }
861  }
862 
863  return (newUpper>=newLower) ? 0 : 1;
864}
865
866/* Builds up row basis backwards (until original model).
867   Returns NULL or previous one to apply .
868   Depends on Free being 0 and impossible for cuts
869*/
870
871CbcNodeInfo * 
872CbcPartialNodeInfo::buildRowBasis(CoinWarmStartBasis & basis ) const 
873
874{ basis.applyDiff(basisDiff_) ;
875
876  return parent_ ; }
877
878CbcNode::CbcNode() :
879  nodeInfo_(NULL),
880  objectiveValue_(1.0e100),
881  guessedObjectiveValue_(1.0e100),
882  sumInfeasibilities_(0.0),
883  branch_(NULL),
884  depth_(-1),
885  numberUnsatisfied_(0),
886  nodeNumber_(-1),
887  state_(0)
888{
889#ifdef CHECK_NODE
890  printf("CbcNode %x Constructor\n",this);
891#endif
892}
893// Print
894void 
895CbcNode::print() const
896{
897  printf("number %d obj %g depth %d sumun %g nunsat %d state %d\n",
898         nodeNumber_,objectiveValue_,depth_,sumInfeasibilities_,numberUnsatisfied_,state_);
899}
900CbcNode::CbcNode(CbcModel * model,
901                 CbcNode * lastNode) :
902  nodeInfo_(NULL),
903  objectiveValue_(1.0e100),
904  guessedObjectiveValue_(1.0e100),
905  sumInfeasibilities_(0.0),
906  branch_(NULL),
907  depth_(-1),
908  numberUnsatisfied_(0),
909  nodeNumber_(-1),
910  state_(0)
911{
912#ifdef CHECK_NODE
913  printf("CbcNode %x Constructor from model\n",this);
914#endif
915  model->setObjectiveValue(this,lastNode);
916
917  if (lastNode) {
918    if (lastNode->nodeInfo_) {
919       lastNode->nodeInfo_->increment();
920    }
921  }
922  nodeNumber_= model->getNodeCount();
923}
924
925#define CBC_NEW_CREATEINFO
926#ifdef CBC_NEW_CREATEINFO
927
928/*
929  New createInfo, with basis manipulation hidden inside mergeBasis. Allows
930  solvers to override and carry over all information from one basis to
931  another.
932*/
933
934void
935CbcNode::createInfo (CbcModel *model,
936                     CbcNode *lastNode,
937                     const CoinWarmStartBasis *lastws,
938                     const double *lastLower, const double *lastUpper,
939                     int numberOldActiveCuts, int numberNewCuts)
940
941{ OsiSolverInterface *solver = model->solver();
942  CbcStrategy *strategy = model->strategy();
943/*
944  The root --- no parent. Create full basis and bounds information.
945*/
946  if (!lastNode)
947  { 
948    if (!strategy)
949      nodeInfo_=new CbcFullNodeInfo(model,solver->getNumRows());
950    else
951      nodeInfo_ = strategy->fullNodeInfo(model,solver->getNumRows());
952  } else {
953/*
954  Not the root. Create an edit from the parent's basis & bound information.
955  This is not quite as straightforward as it seems. We need to reintroduce
956  cuts we may have dropped out of the basis, in the correct position, because
957  this whole process is strictly positional. Start by grabbing the current
958  basis.
959*/
960    bool mustDeleteBasis;
961    const CoinWarmStartBasis *ws =
962      dynamic_cast<const CoinWarmStartBasis*>(solver->getPointerToWarmStart(mustDeleteBasis));
963    assert(ws!=NULL); // make sure not volume
964    //int numberArtificials = lastws->getNumArtificial();
965    int numberColumns = solver->getNumCols();
966    int numberRowsAtContinuous = model->numberRowsAtContinuous();
967    int currentNumberCuts = model->currentNumberCuts();
968#   ifdef CBC_CHECK_BASIS
969    std::cout
970      << "Before expansion: orig " << numberRowsAtContinuous
971      << ", old " << numberOldActiveCuts
972      << ", new " << numberNewCuts
973      << ", current " << currentNumberCuts << "." << std::endl ;
974    ws->print();
975#   endif
976/*
977  Clone the basis and resize it to hold the structural constraints, plus
978  all the cuts: old cuts, both active and inactive (currentNumberCuts),
979  and new cuts (numberNewCuts). This will become the expanded basis.
980*/
981    CoinWarmStartBasis *expanded = 
982      dynamic_cast<CoinWarmStartBasis *>(ws->clone()) ;
983    int iCompact = numberRowsAtContinuous+numberOldActiveCuts+numberNewCuts ;
984    // int nPartial = numberRowsAtContinuous+currentNumberCuts;
985    int iFull = numberRowsAtContinuous+currentNumberCuts+numberNewCuts;
986    // int maxBasisLength = ((iFull+15)>>4)+((numberColumns+15)>>4);
987    // printf("l %d full %d\n",maxBasisLength,iFull);
988    expanded->resize(iFull,numberColumns);
989#   ifdef CBC_CHECK_BASIS
990    std::cout
991      << "\tFull basis " << iFull << " rows, "
992      << numberColumns << " columns; compact "
993      << iCompact << " rows." << std::endl ;
994#   endif
995/*
996  Now flesh out the expanded basis. The clone already has the
997  correct status information for the variables and for the structural
998  (numberRowsAtContinuous) constraints. Any indices beyond nPartial must be
999  cuts created while processing this node --- they can be copied en bloc
1000  into the correct position in the expanded basis. The space reserved for
1001  xferRows is a gross overestimate.
1002*/
1003    CoinWarmStartBasis::XferVec xferRows ;
1004    xferRows.reserve(iFull-numberRowsAtContinuous+1) ;
1005    if (numberNewCuts) {
1006      xferRows.push_back(
1007          CoinWarmStartBasis::XferEntry(iCompact-numberNewCuts,
1008                                        iFull-numberNewCuts,numberNewCuts)) ;
1009    }
1010/*
1011  From nPartial down, record the entries we want to copy from the current
1012  basis (the entries for the active cuts; non-zero in the list returned
1013  by addedCuts). Fill the expanded basis with entries showing a status of
1014  basic for the deactivated (loose) cuts.
1015*/
1016    CbcCountRowCut **cut = model->addedCuts();
1017    iFull -= (numberNewCuts+1) ;
1018    iCompact -= (numberNewCuts+1) ;
1019    int runLen = 0 ;
1020    CoinWarmStartBasis::XferEntry entry(-1,-1,-1) ;
1021    while (iFull >= numberRowsAtContinuous) { 
1022      for ( ; iFull >= numberRowsAtContinuous &&
1023              cut[iFull-numberRowsAtContinuous] ; iFull--)
1024        runLen++ ;
1025      if (runLen) {
1026        iCompact -= runLen ;
1027        entry.first = iCompact+1 ;
1028        entry.second = iFull+1 ;
1029        entry.third = runLen ;
1030        runLen = 0 ;
1031        xferRows.push_back(entry) ;
1032      }
1033      for ( ; iFull >= numberRowsAtContinuous &&
1034              !cut[iFull-numberRowsAtContinuous] ; iFull--)
1035        expanded->setArtifStatus(iFull,CoinWarmStartBasis::basic);
1036    }
1037/*
1038  Finally, call mergeBasis to copy over entries from the current basis to
1039  the expanded basis. Since we cloned the expanded basis from the active basis
1040  and haven't changed the number of variables, only row status entries need
1041  to be copied.
1042*/
1043    expanded->mergeBasis(ws,&xferRows,0) ;
1044
1045#ifdef CBC_CHECK_BASIS
1046    std::cout << "Expanded basis:" << std::endl ;
1047    expanded->print() ;
1048    std::cout << "Diffing against:" << std::endl ;
1049    lastws->print() ;
1050#endif   
1051    assert (expanded->getNumArtificial()>=lastws->getNumArtificial());
1052    if (!expanded->fullBasis()) {
1053      int iFull = numberRowsAtContinuous+currentNumberCuts+numberNewCuts;
1054      printf("cont %d old %d new %d current %d full inc %d full %d\n",
1055             numberRowsAtContinuous,numberOldActiveCuts,numberNewCuts,
1056             currentNumberCuts,iFull,iFull-numberNewCuts);
1057    }
1058    assert (lastws->fullBasis());
1059
1060/*
1061  Now that we have two bases in proper positional correspondence, creating
1062  the actual diff is dead easy.
1063
1064  Note that we're going to compare the expanded basis here to the stripped
1065  basis (lastws) produced by addCuts. It doesn't affect the correctness (the
1066  diff process has no knowledge of the meaning of an entry) but it does
1067  mean that we'll always generate a whack of diff entries because the expanded
1068  basis is considerably larger than the stripped basis.
1069*/
1070    CoinWarmStartDiff *basisDiff = expanded->generateDiff(lastws) ;
1071/*
1072  Diff the bound vectors. It's assumed the number of structural variables
1073  is not changing. For branching objects that change bounds on integer
1074  variables, we should see at least one bound change as a consequence
1075  of applying the branch that generated this subproblem from its parent.
1076  This need not hold for other types of branching objects (hyperplane
1077  branches, for example).
1078*/
1079    const double * lower = solver->getColLower();
1080    const double * upper = solver->getColUpper();
1081
1082    double *boundChanges = new double [2*numberColumns] ;
1083    int *variables = new int [2*numberColumns] ;
1084    int numberChangedBounds=0;
1085   
1086    int i;
1087    for (i=0;i<numberColumns;i++) {
1088      if (lower[i]!=lastLower[i]) {
1089        variables[numberChangedBounds]=i;
1090        boundChanges[numberChangedBounds++]=lower[i];
1091      }
1092      if (upper[i]!=lastUpper[i]) {
1093        variables[numberChangedBounds]=i|0x80000000;
1094        boundChanges[numberChangedBounds++]=upper[i];
1095      }
1096#ifdef CBC_DEBUG
1097      if (lower[i] != lastLower[i]) {
1098        std::cout
1099          << "lower on " << i << " changed from "
1100          << lastLower[i] << " to " << lower[i] << std::endl ;
1101      }
1102      if (upper[i] != lastUpper[i]) {
1103        std::cout
1104          << "upper on " << i << " changed from "
1105          << lastUpper[i] << " to " << upper[i] << std::endl ;
1106      }
1107#endif
1108    }
1109#ifdef CBC_DEBUG
1110    std::cout << numberChangedBounds << " changed bounds." << std::endl ;
1111#endif
1112    //if (lastNode->branchingObject()->boundBranch())
1113    //assert (numberChangedBounds);
1114/*
1115  Hand the lot over to the CbcPartialNodeInfo constructor, then clean up and
1116  return.
1117*/
1118    if (!strategy)
1119      nodeInfo_ =
1120        new CbcPartialNodeInfo(lastNode->nodeInfo_,this,numberChangedBounds,
1121                               variables,boundChanges,basisDiff) ;
1122    else
1123      nodeInfo_ =
1124        strategy->partialNodeInfo(model,lastNode->nodeInfo_,this,
1125                                  numberChangedBounds,variables,boundChanges,
1126                                  basisDiff) ;
1127    delete basisDiff ;
1128    delete [] boundChanges;
1129    delete [] variables;
1130    delete expanded ;
1131    if  (mustDeleteBasis)
1132      delete ws;
1133  }
1134  // Set node number
1135  nodeInfo_->setNodeNumber(model->getNodeCount2());
1136  state_ |= 2; // say active
1137}
1138
1139#else   // CBC_NEW_CREATEINFO
1140
1141/*
1142  Original createInfo, with bare manipulation of basis vectors. Fails if solver
1143  maintains additional information in basis.
1144*/
1145
1146void
1147CbcNode::createInfo (CbcModel *model,
1148                     CbcNode *lastNode,
1149                     const CoinWarmStartBasis *lastws,
1150                     const double *lastLower, const double *lastUpper,
1151                     int numberOldActiveCuts,int numberNewCuts)
1152{ OsiSolverInterface * solver = model->solver();
1153 CbcStrategy * strategy = model->strategy();
1154/*
1155  The root --- no parent. Create full basis and bounds information.
1156*/
1157  if (!lastNode)
1158  { 
1159    if (!strategy)
1160      nodeInfo_=new CbcFullNodeInfo(model,solver->getNumRows());
1161    else
1162      nodeInfo_ = strategy->fullNodeInfo(model,solver->getNumRows());
1163  }
1164/*
1165  Not the root. Create an edit from the parent's basis & bound information.
1166  This is not quite as straightforward as it seems. We need to reintroduce
1167  cuts we may have dropped out of the basis, in the correct position, because
1168  this whole process is strictly positional. Start by grabbing the current
1169  basis.
1170*/
1171  else
1172  { 
1173    bool mustDeleteBasis;
1174    const CoinWarmStartBasis* ws =
1175      dynamic_cast<const CoinWarmStartBasis*>(solver->getPointerToWarmStart(mustDeleteBasis));
1176    assert(ws!=NULL); // make sure not volume
1177    //int numberArtificials = lastws->getNumArtificial();
1178    int numberColumns = solver->getNumCols();
1179   
1180    const double * lower = solver->getColLower();
1181    const double * upper = solver->getColUpper();
1182
1183    int i;
1184/*
1185  Create a clone and resize it to hold all the structural constraints, plus
1186  all the cuts: old cuts, both active and inactive (currentNumberCuts), and
1187  new cuts (numberNewCuts).
1188
1189  TODO: You'd think that the set of constraints (logicals) in the expanded
1190        basis should match the set represented in lastws. At least, that's
1191        what I thought. But at the point I first looked hard at this bit of
1192        code, it turned out that lastws was the stripped basis produced at
1193        the end of addCuts(), rather than the raw basis handed back by
1194        addCuts1(). The expanded basis here is equivalent to the raw basis of
1195        addCuts1(). I said ``whoa, that's not good, I must have introduced a
1196        bug'' and went back to John's code to see where I'd gone wrong.
1197        And discovered the same `error' in his code.
1198
1199        After a bit of thought, my conclusion is that correctness is not
1200        affected by whether lastws is the stripped or raw basis. The diffs
1201        have no semantics --- just a set of changes that need to be made
1202        to convert lastws into expanded. I think the only effect is that we
1203        store a lot more diffs (everything in expanded that's not covered by
1204        the stripped basis). But I need to give this more thought. There
1205        may well be some subtle error cases.
1206
1207        In the mean time, I've twiddled addCuts() to set lastws to the raw
1208        basis. Makes me (Lou) less nervous to compare apples to apples.
1209*/
1210    CoinWarmStartBasis *expanded = 
1211      dynamic_cast<CoinWarmStartBasis *>(ws->clone()) ;
1212    int numberRowsAtContinuous = model->numberRowsAtContinuous();
1213    int iFull = numberRowsAtContinuous+model->currentNumberCuts()+
1214      numberNewCuts;
1215    //int numberArtificialsNow = iFull;
1216    //int maxBasisLength = ((iFull+15)>>4)+((numberColumns+15)>>4);
1217    //printf("l %d full %d\n",maxBasisLength,iFull);
1218    if (expanded) 
1219      expanded->resize(iFull,numberColumns);
1220#ifdef CBC_CHECK_BASIS
1221    printf("Before expansion: orig %d, old %d, new %d, current %d\n",
1222           numberRowsAtContinuous,numberOldActiveCuts,numberNewCuts,
1223           model->currentNumberCuts()) ;
1224    ws->print();
1225#endif
1226/*
1227  Now fill in the expanded basis. Any indices beyond nPartial must
1228  be cuts created while processing this node --- they can be copied directly
1229  into the expanded basis. From nPartial down, pull the status of active cuts
1230  from ws, interleaving with a B entry for the deactivated (loose) cuts.
1231*/
1232    int numberDropped = model->currentNumberCuts()-numberOldActiveCuts;
1233    int iCompact=iFull-numberDropped;
1234    CbcCountRowCut ** cut = model->addedCuts();
1235    int nPartial = model->currentNumberCuts()+numberRowsAtContinuous;
1236    iFull--;
1237    for (;iFull>=nPartial;iFull--) {
1238      CoinWarmStartBasis::Status status = ws->getArtifStatus(--iCompact);
1239      //assert (status != CoinWarmStartBasis::basic); // may be permanent cut
1240      expanded->setArtifStatus(iFull,status);
1241    }
1242    for (;iFull>=numberRowsAtContinuous;iFull--) {
1243      if (cut[iFull-numberRowsAtContinuous]) {
1244        CoinWarmStartBasis::Status status = ws->getArtifStatus(--iCompact);
1245        // If no cut generator being used then we may have basic variables
1246        //if (model->getMaximumCutPasses()&&
1247        //  status == CoinWarmStartBasis::basic)
1248        //printf("cut basic\n");
1249        expanded->setArtifStatus(iFull,status);
1250      } else {
1251        expanded->setArtifStatus(iFull,CoinWarmStartBasis::basic);
1252      }
1253    }
1254#ifdef CBC_CHECK_BASIS
1255    printf("Expanded basis\n");
1256    expanded->print() ;
1257    printf("Diffing against\n") ;
1258    lastws->print() ;
1259#endif   
1260/*
1261  Now that we have two bases in proper positional correspondence, creating
1262  the actual diff is dead easy.
1263*/
1264
1265    CoinWarmStartDiff *basisDiff = expanded->generateDiff(lastws) ;
1266/*
1267  Diff the bound vectors. It's assumed the number of structural variables is
1268  not changing. Assuming that branching objects all involve integer variables,
1269  we should see at least one bound change as a consequence of processing this
1270  subproblem. Different types of branching objects could break this assertion.
1271  Not true at all - we have not applied current branch - JJF.
1272*/
1273    double *boundChanges = new double [2*numberColumns] ;
1274    int *variables = new int [2*numberColumns] ;
1275    int numberChangedBounds=0;
1276    for (i=0;i<numberColumns;i++) {
1277      if (lower[i]!=lastLower[i]) {
1278        variables[numberChangedBounds]=i;
1279        boundChanges[numberChangedBounds++]=lower[i];
1280      }
1281      if (upper[i]!=lastUpper[i]) {
1282        variables[numberChangedBounds]=i|0x80000000;
1283        boundChanges[numberChangedBounds++]=upper[i];
1284      }
1285#ifdef CBC_DEBUG
1286      if (lower[i]!=lastLower[i])
1287        printf("lower on %d changed from %g to %g\n",
1288               i,lastLower[i],lower[i]);
1289      if (upper[i]!=lastUpper[i])
1290        printf("upper on %d changed from %g to %g\n",
1291               i,lastUpper[i],upper[i]);
1292#endif
1293    }
1294#ifdef CBC_DEBUG
1295    printf("%d changed bounds\n",numberChangedBounds) ;
1296#endif
1297    //if (lastNode->branchingObject()->boundBranch())
1298    //assert (numberChangedBounds);
1299/*
1300  Hand the lot over to the CbcPartialNodeInfo constructor, then clean up and
1301  return.
1302*/
1303    if (!strategy)
1304      nodeInfo_ =
1305        new CbcPartialNodeInfo(lastNode->nodeInfo_,this,numberChangedBounds,
1306                               variables,boundChanges,basisDiff) ;
1307    else
1308      nodeInfo_ = strategy->partialNodeInfo(model, lastNode->nodeInfo_,this,numberChangedBounds,
1309                               variables,boundChanges,basisDiff) ;
1310    delete basisDiff ;
1311    delete [] boundChanges;
1312    delete [] variables;
1313    delete expanded ;
1314    if  (mustDeleteBasis)
1315      delete ws;
1316  }
1317  // Set node number
1318  nodeInfo_->setNodeNumber(model->getNodeCount2());
1319  state_ |= 2; // say active
1320}
1321
1322#endif  // CBC_NEW_CREATEINFO
1323
1324/*
1325  The routine scans through the object list of the model looking for objects
1326  that indicate infeasibility. It tests each object using strong branching
1327  and selects the one with the least objective degradation.  A corresponding
1328  branching object is left attached to lastNode.
1329
1330  If strong branching is disabled, a candidate object is chosen essentially
1331  at random (whatever object ends up in pos'n 0 of the candidate array).
1332
1333  If a branching candidate is found to be monotone, bounds are set to fix the
1334  variable and the routine immediately returns (the caller is expected to
1335  reoptimize).
1336
1337  If a branching candidate is found to result in infeasibility in both
1338  directions, the routine immediately returns an indication of infeasibility.
1339
1340  Returns:  0   both branch directions are feasible
1341           -1   branching variable is monotone
1342           -2   infeasible
1343
1344  Original comments:
1345    Here could go cuts etc etc
1346    For now just fix on objective from strong branching.
1347*/
1348
1349int CbcNode::chooseBranch (CbcModel *model, CbcNode *lastNode,int numberPassesLeft)
1350
1351{ if (lastNode)
1352    depth_ = lastNode->depth_+1;
1353  else
1354    depth_ = 0;
1355  delete branch_;
1356  branch_=NULL;
1357  OsiSolverInterface * solver = model->solver();
1358  double saveObjectiveValue = solver->getObjValue();
1359  double objectiveValue = CoinMax(solver->getObjSense()*saveObjectiveValue,objectiveValue_);
1360  const double * lower = solver->getColLower();
1361  const double * upper = solver->getColUpper();
1362  // See what user thinks
1363  int anyAction=model->problemFeasibility()->feasible(model,0);
1364  if (anyAction) {
1365    // will return -2 if infeasible , 0 if treat as integer
1366    return anyAction-1;
1367  }
1368  double integerTolerance = 
1369    model->getDblParam(CbcModel::CbcIntegerTolerance);
1370  // point to useful information
1371  OsiBranchingInformation usefulInfo = model->usefulInformation();
1372  // and modify
1373  usefulInfo.depth_=depth_;
1374  int i;
1375  bool beforeSolution = model->getSolutionCount()==0;
1376  int numberStrong=model->numberStrong();
1377  // switch off strong if hotstart
1378  if (model->hotstartSolution())
1379    numberStrong=0;
1380  int numberStrongDone=0;
1381  int numberUnfinished=0;
1382  int numberStrongInfeasible=0;
1383  int numberStrongIterations=0;
1384  int saveNumberStrong=numberStrong;
1385  int numberObjects = model->numberObjects();
1386  bool checkFeasibility = numberObjects>model->numberIntegers();
1387  int maximumStrong = CoinMax(CoinMin(numberStrong,numberObjects),1);
1388  int numberColumns = model->getNumCols();
1389  double * saveUpper = new double[numberColumns];
1390  double * saveLower = new double[numberColumns];
1391
1392  // Save solution in case heuristics need good solution later
1393 
1394  double * saveSolution = new double[numberColumns];
1395  memcpy(saveSolution,solver->getColSolution(),numberColumns*sizeof(double));
1396  model->reserveCurrentSolution(saveSolution);
1397  /*
1398    Get a branching decision object. Use the default decision criteria unless
1399    the user has loaded a decision method into the model.
1400  */
1401  CbcBranchDecision *decision = model->branchingMethod();
1402  CbcDynamicPseudoCostBranchingObject * dynamicBranchingObject =
1403    dynamic_cast<CbcDynamicPseudoCostBranchingObject *>(decision);
1404  if (!decision||dynamicBranchingObject)
1405    decision = new CbcBranchDefaultDecision();
1406  decision->initialize(model);
1407  CbcStrongInfo * choice = new CbcStrongInfo[maximumStrong];
1408  for (i=0;i<numberColumns;i++) {
1409    saveLower[i] = lower[i];
1410    saveUpper[i] = upper[i];
1411  }
1412  // May go round twice if strong branching fixes all local candidates
1413  bool finished=false;
1414  double estimatedDegradation=0.0; 
1415  while(!finished) {
1416    finished=true;
1417    // Some objects may compute an estimate of best solution from here
1418    estimatedDegradation=0.0; 
1419    //int numberIntegerInfeasibilities=0; // without odd ones
1420    numberStrongDone=0;
1421    numberUnfinished=0;
1422    numberStrongInfeasible=0;
1423    numberStrongIterations=0;
1424   
1425    // We may go round this loop twice (only if we think we have solution)
1426    for (int iPass=0;iPass<2;iPass++) {
1427     
1428      // compute current state
1429      //int numberObjectInfeasibilities; // just odd ones
1430      //model->feasibleSolution(
1431      //                      numberIntegerInfeasibilities,
1432      //                      numberObjectInfeasibilities);
1433      const double * hotstartSolution = model->hotstartSolution();
1434      const int * hotstartPriorities = model->hotstartPriorities();
1435     
1436      // Some objects may compute an estimate of best solution from here
1437      estimatedDegradation=0.0; 
1438      numberUnsatisfied_ = 0;
1439      // initialize sum of "infeasibilities"
1440      sumInfeasibilities_ = 0.0;
1441      int bestPriority=COIN_INT_MAX;
1442      /*
1443        Scan for branching objects that indicate infeasibility. Choose the best
1444        maximumStrong candidates, using priority as the first criteria, then
1445        integer infeasibility.
1446       
1447        The algorithm is to fill the choice array with a set of good candidates (by
1448        infeasibility) with priority bestPriority.  Finding a candidate with
1449        priority better (less) than bestPriority flushes the choice array. (This
1450        serves as initialization when the first candidate is found.)
1451       
1452        A new candidate is added to choices only if its infeasibility exceeds the
1453        current max infeasibility (mostAway). When a candidate is added, it
1454        replaces the candidate with the smallest infeasibility (tracked by
1455        iSmallest).
1456      */
1457      int iSmallest = 0;
1458      double mostAway = 1.0e-100;
1459      for (i = 0 ; i < maximumStrong ; i++)
1460        choice[i].possibleBranch = NULL ;
1461      numberStrong=0;
1462      bool canDoOneHot=false;
1463      for (i=0;i<numberObjects;i++) {
1464        OsiObject * object = model->modifiableObject(i);
1465        int preferredWay;
1466        double infeasibility = object->infeasibility(&usefulInfo,preferredWay);
1467        int priorityLevel = object->priority();
1468        if (hotstartSolution) {
1469          // we are doing hot start
1470          const CbcSimpleInteger * thisOne = dynamic_cast <const CbcSimpleInteger *> (object);
1471          if (thisOne) {
1472            int iColumn = thisOne->columnNumber();
1473            bool canDoThisHot=true;
1474            double targetValue = hotstartSolution[iColumn];
1475            if (saveUpper[iColumn]>saveLower[iColumn]) {
1476              double value = saveSolution[iColumn];
1477              if (hotstartPriorities)
1478                priorityLevel=hotstartPriorities[iColumn]; 
1479              //double originalLower = thisOne->originalLower();
1480              //double originalUpper = thisOne->originalUpper();
1481              // switch off if not possible
1482              if (targetValue>=saveLower[iColumn]&&targetValue<=saveUpper[iColumn]) {
1483                /* priority outranks rest always if negative
1484                   otherwise can be downgraded if at correct level.
1485                   Infeasibility may be increased to choose 1.0 values first.
1486                   choose one near wanted value
1487                */
1488                if (fabs(value-targetValue)>integerTolerance) {
1489                  infeasibility = 1.0-fabs(value-targetValue);
1490                  if (targetValue==1.0)
1491                    infeasibility += 1.0;
1492                  if (value>targetValue) {
1493                    preferredWay=-1;
1494                  } else {
1495                    preferredWay=1;
1496                  }
1497                  priorityLevel = CoinAbs(priorityLevel);
1498                } else if (priorityLevel<0) {
1499                  priorityLevel = CoinAbs(priorityLevel);
1500                  if (targetValue==saveLower[iColumn]) {
1501                    infeasibility = integerTolerance+1.0e-12;
1502                    preferredWay=-1;
1503                  } else if (targetValue==saveUpper[iColumn]) {
1504                    infeasibility = integerTolerance+1.0e-12;
1505                    preferredWay=1;
1506                  } else {
1507                    // can't
1508                    priorityLevel += 10000000;
1509                    canDoThisHot=false;
1510                  }
1511                } else {
1512                  priorityLevel += 10000000;
1513                  canDoThisHot=false;
1514                }
1515              } else {
1516                // switch off if not possible
1517                canDoThisHot=false;
1518              }
1519              if (canDoThisHot)
1520                canDoOneHot=true;
1521            } else if (targetValue<saveLower[iColumn]||targetValue>saveUpper[iColumn]) {
1522            }
1523          } else {
1524            priorityLevel += 10000000;
1525          }
1526        }
1527        if (infeasibility) {
1528          // Increase estimated degradation to solution
1529          estimatedDegradation += CoinMin(object->upEstimate(),object->downEstimate());
1530          numberUnsatisfied_++;
1531          sumInfeasibilities_ += infeasibility;
1532          // Better priority? Flush choices.
1533          if (priorityLevel<bestPriority) {
1534            int j;
1535            iSmallest=0;
1536            for (j=0;j<maximumStrong;j++) {
1537              choice[j].upMovement=0.0;
1538              delete choice[j].possibleBranch;
1539              choice[j].possibleBranch=NULL;
1540            }
1541            bestPriority = priorityLevel;
1542            mostAway=1.0e-100;
1543            numberStrong=0;
1544          } else if (priorityLevel>bestPriority) {
1545            continue;
1546          }
1547          // Check for suitability based on infeasibility.
1548          if (infeasibility>mostAway) {
1549            //add to list
1550            choice[iSmallest].upMovement=infeasibility;
1551            delete choice[iSmallest].possibleBranch;
1552            CbcSimpleInteger * obj =
1553              dynamic_cast <CbcSimpleInteger *>(object) ;
1554            if (obj) {
1555              choice[iSmallest].possibleBranch=obj->createBranch(solver,&usefulInfo,preferredWay);
1556            } else {
1557              CbcObject * obj =
1558                dynamic_cast <CbcObject *>(object) ;
1559              assert (obj);
1560              choice[iSmallest].possibleBranch=obj->createBranch(preferredWay);
1561            }
1562            numberStrong = CoinMax(numberStrong,iSmallest+1);
1563            // Save which object it was
1564            choice[iSmallest].objectNumber=i;
1565            int j;
1566            iSmallest=-1;
1567            mostAway = 1.0e50;
1568            for (j=0;j<maximumStrong;j++) {
1569              if (choice[j].upMovement<mostAway) {
1570                mostAway=choice[j].upMovement;
1571                iSmallest=j;
1572              }
1573            }
1574          }
1575        }
1576      }
1577      if (!canDoOneHot&&hotstartSolution) {
1578        // switch off as not possible
1579        hotstartSolution=NULL;
1580        model->setHotstartSolution(NULL,NULL);
1581      }
1582      if (numberUnsatisfied_) {
1583        // some infeasibilities - go to next steps
1584        break;
1585      } else if (!iPass) {
1586        // looks like a solution - get paranoid
1587        bool roundAgain=false;
1588        // get basis
1589        CoinWarmStartBasis * ws = dynamic_cast<CoinWarmStartBasis*>(solver->getWarmStart());
1590        if (!ws)
1591          break;
1592        for (i=0;i<numberColumns;i++) {
1593          double value = saveSolution[i];
1594          if (value<lower[i]) {
1595            saveSolution[i]=lower[i];
1596            roundAgain=true;
1597            ws->setStructStatus(i,CoinWarmStartBasis::atLowerBound);
1598          } else if (value>upper[i]) {
1599            saveSolution[i]=upper[i];
1600            roundAgain=true;
1601            ws->setStructStatus(i,CoinWarmStartBasis::atUpperBound);
1602          } 
1603        }
1604        if (roundAgain&&saveNumberStrong) {
1605          // restore basis
1606          solver->setWarmStart(ws);
1607          delete ws;
1608          solver->resolve();
1609          memcpy(saveSolution,solver->getColSolution(),numberColumns*sizeof(double));
1610          model->reserveCurrentSolution(saveSolution);
1611          if (!solver->isProvenOptimal()) {
1612            // infeasible
1613            anyAction=-2;
1614            break;
1615          }
1616        } else {
1617          delete ws;
1618          break;
1619        }
1620      }
1621    }
1622    /* Some solvers can do the strong branching calculations faster if
1623       they do them all at once.  At present only Clp does for ordinary
1624       integers but I think this coding would be easy to modify
1625    */
1626    bool allNormal=true; // to say if we can do fast strong branching
1627    // Say which one will be best
1628    int bestChoice=0;
1629    double worstInfeasibility=0.0;
1630    for (i=0;i<numberStrong;i++) {
1631      choice[i].numIntInfeasUp = numberUnsatisfied_;
1632      choice[i].numIntInfeasDown = numberUnsatisfied_;
1633      choice[i].fix=0; // say not fixed
1634      if (!dynamic_cast <const CbcSimpleInteger *> (model->object(choice[i].objectNumber)))
1635        allNormal=false; // Something odd so lets skip clever fast branching
1636      if ( !model->object(choice[i].objectNumber)->boundBranch())
1637        numberStrong=0; // switch off
1638      if ( choice[i].possibleBranch->numberBranches()>2)
1639        numberStrong=0; // switch off
1640      // Do best choice in case switched off
1641      if (choice[i].upMovement>worstInfeasibility) {
1642        worstInfeasibility=choice[i].upMovement;
1643        bestChoice=i;
1644      }
1645    }
1646    // If we have hit max time don't do strong branching
1647    bool hitMaxTime = ( CoinCpuTime()-model->getDblParam(CbcModel::CbcStartSeconds) > 
1648                        model->getDblParam(CbcModel::CbcMaximumSeconds));
1649    // also give up if we are looping round too much
1650    if (hitMaxTime||numberPassesLeft<=0)
1651      numberStrong=0;
1652    /*
1653      Is strong branching enabled? If so, set up and do it. Otherwise, we'll
1654      fall through to simple branching.
1655     
1656      Setup for strong branching involves saving the current basis (for restoration
1657      afterwards) and setting up for hot starts.
1658    */
1659    if (numberStrong&&saveNumberStrong) {
1660     
1661      bool solveAll=false; // set true to say look at all even if some fixed (experiment)
1662      solveAll=true;
1663      // worth trying if too many times
1664      // Save basis
1665      CoinWarmStart * ws = solver->getWarmStart();
1666      // save limit
1667      int saveLimit;
1668      solver->getIntParam(OsiMaxNumIterationHotStart,saveLimit);
1669      if (beforeSolution&&saveLimit<100)
1670        solver->setIntParam(OsiMaxNumIterationHotStart,100); // go to end
1671#     ifdef COIN_HAS_CLP     
1672      /* If we are doing all strong branching in one go then we create new arrays
1673         to store information.  If clp NULL then doing old way.
1674         Going down -
1675         outputSolution[2*i] is final solution.
1676         outputStuff[2*i] is status (0 - finished, 1 infeas, other unknown
1677         outputStuff[2*i+numberStrong] is number iterations
1678         On entry newUpper[i] is new upper bound, on exit obj change
1679         Going up -
1680         outputSolution[2*i+1] is final solution.
1681         outputStuff[2*i+1] is status (0 - finished, 1 infeas, other unknown
1682         outputStuff[2*i+1+numberStrong] is number iterations
1683       On entry newLower[i] is new lower bound, on exit obj change
1684      */
1685      OsiClpSolverInterface * osiclp = dynamic_cast< OsiClpSolverInterface*> (solver);
1686      ClpSimplex * clp=NULL;
1687      double * newLower = NULL;
1688      double * newUpper = NULL;
1689      double ** outputSolution=NULL;
1690      int * outputStuff=NULL;
1691      // Go back to normal way if user wants it
1692      if (osiclp&&(osiclp->specialOptions()&16)!=0&&osiclp->specialOptions()>0)
1693        allNormal=false;
1694      if (osiclp&&!allNormal) {
1695        // say do fast
1696        int easy=1;
1697        osiclp->setHintParam(OsiDoInBranchAndCut,true,OsiHintDo,&easy) ;
1698      }
1699      if (osiclp&& allNormal) {
1700        clp = osiclp->getModelPtr();
1701        // Clp - do a different way
1702        newLower = new double[numberStrong];
1703        newUpper = new double[numberStrong];
1704        outputSolution = new double * [2*numberStrong];
1705        outputStuff = new int [4*numberStrong];
1706        int * which = new int[numberStrong];
1707        int startFinishOptions;
1708        int specialOptions = osiclp->specialOptions();
1709        int clpOptions = clp->specialOptions();
1710        int returnCode=0;
1711#define CRUNCH
1712#ifdef CRUNCH
1713        // Crunch down problem
1714        int numberRows = clp->numberRows();
1715        // Use dual region
1716        double * rhs = clp->dualRowSolution();
1717        int * whichRow = new int[3*numberRows];
1718        int * whichColumn = new int[2*numberColumns];
1719        int nBound;
1720        ClpSimplex * small = ((ClpSimplexOther *) clp)->crunch(rhs,whichRow,whichColumn,nBound,true);
1721        if (!small) {
1722          anyAction=-2;
1723          //printf("XXXX Inf by inspection\n");
1724          delete [] whichColumn;
1725          whichColumn=NULL;
1726          delete [] whichRow;
1727          whichRow=NULL;
1728          break;
1729        } else {
1730          clp = small;
1731        }
1732#else
1733        int saveLogLevel = clp->logLevel();
1734        int saveMaxIts = clp->maximumIterations();
1735#endif
1736        clp->setLogLevel(0);
1737        if((specialOptions&1)==0) {
1738          startFinishOptions=0;
1739          clp->setSpecialOptions(clpOptions|(64|1024));
1740        } else {
1741          startFinishOptions=1+2+4;
1742          //startFinishOptions=1+4; // for moment re-factorize
1743          if((specialOptions&4)==0) 
1744            clp->setSpecialOptions(clpOptions|(64|128|512|1024|4096));
1745          else
1746            clp->setSpecialOptions(clpOptions|(64|128|512|1024|2048|4096));
1747        }
1748        // User may want to clean up before strong branching
1749        if ((clp->specialOptions()&32)!=0) {
1750          clp->primal(1);
1751          if (clp->numberIterations())
1752            model->messageHandler()->message(CBC_ITERATE_STRONG,*model->messagesPointer())
1753              << clp->numberIterations()
1754              <<CoinMessageEol;
1755        }
1756        clp->setMaximumIterations(saveLimit);
1757#ifdef CRUNCH
1758        int * backColumn = whichColumn+numberColumns;
1759#endif
1760        for (i=0;i<numberStrong;i++) {
1761          int iObject = choice[i].objectNumber;
1762          const OsiObject * object = model->object(iObject);
1763          const CbcSimpleInteger * simple = dynamic_cast <const CbcSimpleInteger *> (object);
1764          int iSequence = simple->columnNumber();
1765          newLower[i]= ceil(saveSolution[iSequence]);
1766          newUpper[i]= floor(saveSolution[iSequence]);
1767#ifdef CRUNCH
1768          iSequence = backColumn[iSequence];
1769          assert (iSequence>=0);
1770#endif
1771          which[i]=iSequence;
1772          outputSolution[2*i]= new double [numberColumns];
1773          outputSolution[2*i+1]= new double [numberColumns];
1774        }
1775        //clp->writeMps("bad");
1776        returnCode=clp->strongBranching(numberStrong,which,
1777                                            newLower, newUpper,outputSolution,
1778                                            outputStuff,outputStuff+2*numberStrong,!solveAll,false,
1779                                            startFinishOptions);
1780#ifndef CRUNCH
1781        clp->setSpecialOptions(clpOptions); // restore
1782        clp->setMaximumIterations(saveMaxIts);
1783        clp->setLogLevel(saveLogLevel);
1784#endif
1785        if (returnCode==-2) {
1786          // bad factorization!!!
1787          // Doing normal way
1788          // Mark hot start
1789          solver->markHotStart();
1790          clp = NULL;
1791        } else {
1792#ifdef CRUNCH
1793          // extract solution
1794          //bool checkSol=true;
1795          for (i=0;i<numberStrong;i++) {
1796            int iObject = choice[i].objectNumber;
1797            const OsiObject * object = model->object(iObject);
1798            const CbcSimpleInteger * simple = dynamic_cast <const CbcSimpleInteger *> (object);
1799            int iSequence = simple->columnNumber();
1800            which[i]=iSequence;
1801            double * sol = outputSolution[2*i];
1802            double * sol2 = outputSolution[2*i+1];
1803            //bool x=true;
1804            //bool x2=true;
1805            for (int iColumn=numberColumns-1;iColumn>=0;iColumn--) {
1806              int jColumn = backColumn[iColumn];
1807              if (jColumn>=0) {
1808                sol[iColumn]=sol[jColumn];
1809                sol2[iColumn]=sol2[jColumn];
1810              } else {
1811                sol[iColumn]=saveSolution[iColumn];
1812                sol2[iColumn]=saveSolution[iColumn];
1813              }
1814            }
1815          }
1816#endif
1817        }
1818#ifdef CRUNCH
1819        delete [] whichColumn;
1820        delete [] whichRow;
1821        delete small;
1822#endif
1823        delete [] which;
1824      } else {
1825        // Doing normal way
1826        // Mark hot start
1827        solver->markHotStart();
1828      }
1829#     else      /* COIN_HAS_CLP */
1830
1831      OsiSolverInterface *clp = NULL ;
1832      double **outputSolution = NULL ;
1833      int *outputStuff = NULL ;
1834      double * newLower = NULL ;
1835      double * newUpper = NULL ;
1836
1837      solver->markHotStart();
1838
1839#     endif     /* COIN_HAS_CLP */
1840      /*
1841        Open a loop to do the strong branching LPs. For each candidate variable,
1842        solve an LP with the variable forced down, then up. If a direction turns
1843        out to be infeasible or monotonic (i.e., over the dual objective cutoff),
1844        force the objective change to be big (1.0e100). If we determine the problem
1845        is infeasible, or find a monotone variable, escape the loop.
1846       
1847        TODO: The `restore bounds' part might be better encapsulated as an
1848        unbranch() method. Branching objects more exotic than simple integers
1849        or cliques might not restrict themselves to variable bounds.
1850
1851        TODO: Virtuous solvers invalidate the current solution (or give bogus
1852        results :-) when the bounds are changed out from under them. So we
1853        need to do all the work associated with finding a new solution before
1854        restoring the bounds.
1855      */
1856      for (i = 0 ; i < numberStrong ; i++)
1857        { double objectiveChange ;
1858        double newObjectiveValue=1.0e100;
1859        // status is 0 finished, 1 infeasible and other
1860        int iStatus;
1861        /*
1862          Try the down direction first. (Specify the initial branching alternative as
1863          down with a call to way(-1). Each subsequent call to branch() performs the
1864          specified branch and advances the branch object state to the next branch
1865          alternative.)
1866        */
1867        if (!clp) {
1868          choice[i].possibleBranch->way(-1) ;
1869          choice[i].possibleBranch->branch() ;
1870          bool feasible=true;
1871          if (checkFeasibility) {
1872            // check branching did not make infeasible
1873            int iColumn;
1874            int numberColumns = solver->getNumCols();
1875            const double * columnLower = solver->getColLower();
1876            const double * columnUpper = solver->getColUpper();
1877            for (iColumn= 0;iColumn<numberColumns;iColumn++) {
1878              if (columnLower[iColumn]>columnUpper[iColumn]+1.0e-5)
1879                feasible=false;
1880            }
1881          }
1882          if (feasible) {
1883            solver->solveFromHotStart() ;
1884            numberStrongDone++;
1885            numberStrongIterations += solver->getIterationCount();
1886            /*
1887              We now have an estimate of objective degradation that we can use for strong
1888              branching. If we're over the cutoff, the variable is monotone up.
1889              If we actually made it to optimality, check for a solution, and if we have
1890              a good one, call setBestSolution to process it. Note that this may reduce the
1891              cutoff, so we check again to see if we can declare this variable monotone.
1892            */
1893            if (solver->isProvenOptimal())
1894              iStatus=0; // optimal
1895            else if (solver->isIterationLimitReached()
1896                     &&!solver->isDualObjectiveLimitReached())
1897              iStatus=2; // unknown
1898            else
1899              iStatus=1; // infeasible
1900            newObjectiveValue = solver->getObjSense()*solver->getObjValue();
1901            choice[i].numItersDown = solver->getIterationCount();
1902          } else {
1903            iStatus=1; // infeasible
1904            newObjectiveValue = 1.0e100;
1905            choice[i].numItersDown = 0;
1906          }
1907        } else {
1908          iStatus = outputStuff[2*i];
1909          choice[i].numItersDown = outputStuff[2*numberStrong+2*i];
1910          numberStrongDone++;
1911          numberStrongIterations += choice[i].numItersDown;
1912          newObjectiveValue = objectiveValue+newUpper[i];
1913          solver->setColSolution(outputSolution[2*i]);
1914        }
1915        objectiveChange = CoinMax(newObjectiveValue  - objectiveValue_,0.0);
1916        if (!iStatus) {
1917          choice[i].finishedDown = true ;
1918          if (newObjectiveValue>=model->getCutoff()) {
1919            objectiveChange = 1.0e100; // say infeasible
1920            numberStrongInfeasible++;
1921          } else {
1922            // See if integer solution
1923            if (model->feasibleSolution(choice[i].numIntInfeasDown,
1924                                        choice[i].numObjInfeasDown)
1925                &&model->problemFeasibility()->feasible(model,-1)>=0) {
1926              model->setBestSolution(CBC_STRONGSOL,
1927                                     newObjectiveValue,
1928                                     solver->getColSolution()) ;
1929              // only needed for odd solvers
1930              newObjectiveValue = solver->getObjSense()*solver->getObjValue();
1931              objectiveChange = CoinMax(newObjectiveValue-objectiveValue_,0.0) ;
1932              model->setLastHeuristic(NULL);
1933              model->incrementUsed(solver->getColSolution());
1934              if (newObjectiveValue >= model->getCutoff()) {    //  *new* cutoff
1935                objectiveChange = 1.0e100 ;
1936                numberStrongInfeasible++;
1937              }
1938            }
1939          }
1940        } else if (iStatus==1) {
1941          objectiveChange = 1.0e100 ;
1942          numberStrongInfeasible++;
1943        } else {
1944          // Can't say much as we did not finish
1945          choice[i].finishedDown = false ;
1946          numberUnfinished++;
1947        }
1948        choice[i].downMovement = objectiveChange ;
1949       
1950        // restore bounds
1951        if (!clp)
1952          { for (int j=0;j<numberColumns;j++) {
1953            if (saveLower[j] != lower[j])
1954              solver->setColLower(j,saveLower[j]);
1955            if (saveUpper[j] != upper[j])
1956              solver->setColUpper(j,saveUpper[j]);
1957          }
1958          }
1959        //printf("Down on %d, status is %d, obj %g its %d cost %g finished %d inf %d infobj %d\n",
1960        //     choice[i].objectNumber,iStatus,newObjectiveValue,choice[i].numItersDown,
1961        //     choice[i].downMovement,choice[i].finishedDown,choice[i].numIntInfeasDown,
1962        //     choice[i].numObjInfeasDown);
1963       
1964        // repeat the whole exercise, forcing the variable up
1965        if (!clp) {
1966          bool feasible=true;
1967          // If odd branching then maybe just one possibility
1968          if(choice[i].possibleBranch->numberBranchesLeft()>0) {
1969            choice[i].possibleBranch->branch();
1970            if (checkFeasibility) {
1971              // check branching did not make infeasible
1972              int iColumn;
1973              int numberColumns = solver->getNumCols();
1974              const double * columnLower = solver->getColLower();
1975              const double * columnUpper = solver->getColUpper();
1976              for (iColumn= 0;iColumn<numberColumns;iColumn++) {
1977                if (columnLower[iColumn]>columnUpper[iColumn]+1.0e-5)
1978                  feasible=false;
1979              }
1980            }
1981          } else {
1982            // second branch infeasible
1983            feasible=false;
1984          }
1985          if (feasible) {
1986            solver->solveFromHotStart() ;
1987            numberStrongDone++;
1988            numberStrongIterations += solver->getIterationCount();
1989            /*
1990              We now have an estimate of objective degradation that we can use for strong
1991              branching. If we're over the cutoff, the variable is monotone up.
1992              If we actually made it to optimality, check for a solution, and if we have
1993              a good one, call setBestSolution to process it. Note that this may reduce the
1994              cutoff, so we check again to see if we can declare this variable monotone.
1995            */
1996            if (solver->isProvenOptimal())
1997              iStatus=0; // optimal
1998            else if (solver->isIterationLimitReached()
1999                     &&!solver->isDualObjectiveLimitReached())
2000              iStatus=2; // unknown
2001            else
2002              iStatus=1; // infeasible
2003            newObjectiveValue = solver->getObjSense()*solver->getObjValue();
2004            choice[i].numItersUp = solver->getIterationCount();
2005          } else {
2006            iStatus=1; // infeasible
2007            newObjectiveValue = 1.0e100;
2008            choice[i].numItersDown = 0;
2009          }
2010        } else {
2011          iStatus = outputStuff[2*i+1];
2012          choice[i].numItersUp = outputStuff[2*numberStrong+2*i+1];
2013          numberStrongDone++;
2014          numberStrongIterations += choice[i].numItersUp;
2015          newObjectiveValue = objectiveValue+newLower[i];
2016          solver->setColSolution(outputSolution[2*i+1]);
2017        }
2018        objectiveChange = CoinMax(newObjectiveValue  - objectiveValue_,0.0);
2019        if (!iStatus) {
2020          choice[i].finishedUp = true ;
2021          if (newObjectiveValue>=model->getCutoff()) {
2022            objectiveChange = 1.0e100; // say infeasible
2023            numberStrongInfeasible++;
2024          } else {
2025            // See if integer solution
2026            if (model->feasibleSolution(choice[i].numIntInfeasUp,
2027                                        choice[i].numObjInfeasUp)
2028                &&model->problemFeasibility()->feasible(model,-1)>=0) {
2029              model->setBestSolution(CBC_STRONGSOL,
2030                                     newObjectiveValue,
2031                                     solver->getColSolution()) ;
2032              // only needed for odd solvers
2033              newObjectiveValue = solver->getObjSense()*solver->getObjValue();
2034              objectiveChange = CoinMax(newObjectiveValue-objectiveValue_,0.0) ;
2035              model->setLastHeuristic(NULL);
2036              model->incrementUsed(solver->getColSolution());
2037              if (newObjectiveValue >= model->getCutoff()) {    //  *new* cutoff
2038                objectiveChange = 1.0e100 ;
2039                numberStrongInfeasible++;
2040              }
2041            }
2042          }
2043        } else if (iStatus==1) {
2044          objectiveChange = 1.0e100 ;
2045          numberStrongInfeasible++;
2046        } else {
2047          // Can't say much as we did not finish
2048          choice[i].finishedUp = false ;
2049          numberUnfinished++;
2050        }
2051        choice[i].upMovement = objectiveChange ;
2052       
2053        // restore bounds
2054        if (!clp)
2055          { for (int j=0;j<numberColumns;j++) {
2056            if (saveLower[j] != lower[j])
2057              solver->setColLower(j,saveLower[j]);
2058            if (saveUpper[j] != upper[j])
2059              solver->setColUpper(j,saveUpper[j]);
2060          }
2061          }
2062       
2063        //printf("Up on %d, status is %d, obj %g its %d cost %g finished %d inf %d infobj %d\n",
2064        //     choice[i].objectNumber,iStatus,newObjectiveValue,choice[i].numItersUp,
2065        //     choice[i].upMovement,choice[i].finishedUp,choice[i].numIntInfeasUp,
2066        //     choice[i].numObjInfeasUp);
2067       
2068        /*
2069          End of evaluation for this candidate variable. Possibilities are:
2070          * Both sides below cutoff; this variable is a candidate for branching.
2071          * Both sides infeasible or above the objective cutoff: no further action
2072          here. Break from the evaluation loop and assume the node will be purged
2073          by the caller.
2074          * One side below cutoff: Install the branch (i.e., fix the variable). Break
2075          from the evaluation loop and assume the node will be reoptimised by the
2076          caller.
2077        */
2078        // reset
2079        choice[i].possibleBranch->resetNumberBranchesLeft();
2080        if (choice[i].upMovement<1.0e100) {
2081          if(choice[i].downMovement<1.0e100) {
2082            // feasible - no action
2083          } else {
2084            // up feasible, down infeasible
2085            anyAction=-1;
2086            //printf("Down infeasible for choice %d sequence %d\n",i,
2087            // model->object(choice[i].objectNumber)->columnNumber());
2088            if (!solveAll) {
2089              choice[i].possibleBranch->way(1);
2090              choice[i].possibleBranch->branch();
2091              break;
2092            } else {
2093              choice[i].fix=1;
2094            }
2095          }
2096        } else {
2097          if(choice[i].downMovement<1.0e100) {
2098            // down feasible, up infeasible
2099            anyAction=-1;
2100            //printf("Up infeasible for choice %d sequence %d\n",i,
2101            // model->object(choice[i].objectNumber)->columnNumber());
2102            if (!solveAll) {
2103              choice[i].possibleBranch->way(-1);
2104              choice[i].possibleBranch->branch();
2105              break;
2106            } else {
2107              choice[i].fix=-1;
2108            }
2109          } else {
2110            // neither side feasible
2111            anyAction=-2;
2112            //printf("Both infeasible for choice %d sequence %d\n",i,
2113            // model->object(choice[i].objectNumber)->columnNumber());
2114            break;
2115          }
2116        }
2117        bool hitMaxTime = ( CoinCpuTime()-model->getDblParam(CbcModel::CbcStartSeconds) > 
2118                            model->getDblParam(CbcModel::CbcMaximumSeconds));
2119        if (hitMaxTime) {
2120          numberStrong=i+1;
2121          break;
2122        }
2123        }
2124      if (!clp) {
2125        // Delete the snapshot
2126        solver->unmarkHotStart();
2127      } else {
2128        delete [] newLower;
2129        delete [] newUpper;
2130        delete [] outputStuff;
2131        int i;
2132        for (i=0;i<2*numberStrong;i++)
2133          delete [] outputSolution[i];
2134        delete [] outputSolution;
2135      }
2136      solver->setIntParam(OsiMaxNumIterationHotStart,saveLimit);
2137      // restore basis
2138      solver->setWarmStart(ws);
2139      // Unless infeasible we will carry on
2140      // But we could fix anyway
2141      if (anyAction==-1&&solveAll) {
2142        // apply and take off
2143        for (i = 0 ; i < numberStrong ; i++) {
2144          if (choice[i].fix) {
2145            choice[i].possibleBranch->way(choice[i].fix) ;
2146            choice[i].possibleBranch->branch() ;
2147          }
2148        }
2149        bool feasible=true;
2150        if (checkFeasibility) {
2151          // check branching did not make infeasible
2152          int iColumn;
2153          int numberColumns = solver->getNumCols();
2154          const double * columnLower = solver->getColLower();
2155          const double * columnUpper = solver->getColUpper();
2156          for (iColumn= 0;iColumn<numberColumns;iColumn++) {
2157            if (columnLower[iColumn]>columnUpper[iColumn]+1.0e-5)
2158              feasible=false;
2159          }
2160        }
2161        if (feasible) {
2162          // can do quick optimality check
2163          int easy=2;
2164          solver->setHintParam(OsiDoInBranchAndCut,true,OsiHintDo,&easy) ;
2165          solver->resolve() ;
2166          solver->setHintParam(OsiDoInBranchAndCut,true,OsiHintDo,NULL) ;
2167          feasible = solver->isProvenOptimal();
2168        }
2169        if (feasible) {
2170          memcpy(saveSolution,solver->getColSolution(),numberColumns*sizeof(double));
2171          model->reserveCurrentSolution(saveSolution);
2172          memcpy(saveLower,solver->getColLower(),numberColumns*sizeof(double));
2173          memcpy(saveUpper,solver->getColUpper(),numberColumns*sizeof(double));
2174          // Clean up all candidates whih are fixed
2175          int numberLeft=0;
2176          for (i = 0 ; i < numberStrong ; i++) {
2177            CbcStrongInfo thisChoice = choice[i];
2178            choice[i].possibleBranch=NULL;
2179            const OsiObject * object = model->object(thisChoice.objectNumber);
2180            int preferredWay;
2181            double infeasibility = object->infeasibility(&usefulInfo,preferredWay);
2182            if (!infeasibility) {
2183              // take out
2184              delete thisChoice.possibleBranch;
2185            } else {
2186              choice[numberLeft++]=thisChoice;
2187            }
2188          }
2189          numberStrong=numberLeft;
2190          for (;i<maximumStrong;i++) {
2191            delete choice[i].possibleBranch;
2192            choice[i].possibleBranch=NULL;
2193          }
2194          // If all fixed then round again
2195          if (!numberLeft) {
2196            finished=false;
2197            numberStrong=0;
2198            saveNumberStrong=0;
2199            maximumStrong=1;
2200          } else {
2201            anyAction=0;
2202          }
2203          // If these two uncommented then different action
2204          anyAction=-1;
2205          finished=true;
2206          //printf("some fixed but continuing %d left\n",numberLeft);
2207        } else {
2208          anyAction=-2; // say infeasible
2209        }
2210      }
2211      delete ws;
2212      int numberNodes = model->getNodeCount();
2213      // update number of strong iterations etc
2214      model->incrementStrongInfo(numberStrongDone,numberStrongIterations,
2215                                 anyAction==-2 ? 0:numberStrongInfeasible,anyAction==-2);
2216     
2217      /*
2218        anyAction >= 0 indicates that strong branching didn't produce any monotone
2219        variables. Sift through the candidates for the best one.
2220       
2221        QUERY: Setting numberNodes looks to be a distributed noop. numberNodes is
2222        local to this code block. Perhaps should be numberNodes_ from model?
2223        Unclear what this calculation is doing.
2224      */
2225      if (anyAction>=0) {
2226       
2227        // get average cost per iteration and assume stopped ones
2228        // would stop after 50% more iterations at average cost??? !!! ???
2229        double averageCostPerIteration=0.0;
2230        double totalNumberIterations=1.0;
2231        int smallestNumberInfeasibilities=COIN_INT_MAX;
2232        for (i=0;i<numberStrong;i++) {
2233          totalNumberIterations += choice[i].numItersDown +
2234            choice[i].numItersUp ;
2235          averageCostPerIteration += choice[i].downMovement +
2236            choice[i].upMovement;
2237          smallestNumberInfeasibilities= 
2238            CoinMin(CoinMin(choice[i].numIntInfeasDown ,
2239                            choice[i].numIntInfeasUp ),
2240                    smallestNumberInfeasibilities);
2241        }
2242        //if (smallestNumberInfeasibilities>=numberIntegerInfeasibilities)
2243        //numberNodes=1000000; // switch off search for better solution
2244        numberNodes=1000000; // switch off anyway
2245        averageCostPerIteration /= totalNumberIterations;
2246        // all feasible - choose best bet
2247       
2248        // New method does all at once so it can be more sophisticated
2249        // in deciding how to balance actions.
2250        // But it does need arrays
2251        double * changeUp = new double [numberStrong];
2252        int * numberInfeasibilitiesUp = new int [numberStrong];
2253        double * changeDown = new double [numberStrong];
2254        int * numberInfeasibilitiesDown = new int [numberStrong];
2255        CbcBranchingObject ** objects = new CbcBranchingObject * [ numberStrong];
2256        for (i = 0 ; i < numberStrong ; i++) {
2257          int iColumn = choice[i].possibleBranch->variable() ;
2258          model->messageHandler()->message(CBC_STRONG,*model->messagesPointer())
2259            << i << iColumn
2260            <<choice[i].downMovement<<choice[i].numIntInfeasDown
2261            <<choice[i].upMovement<<choice[i].numIntInfeasUp
2262            <<choice[i].possibleBranch->value()
2263            <<CoinMessageEol;
2264          changeUp[i]=choice[i].upMovement;
2265          numberInfeasibilitiesUp[i] = choice[i].numIntInfeasUp;
2266          changeDown[i]=choice[i].downMovement;
2267          numberInfeasibilitiesDown[i] = choice[i].numIntInfeasDown;
2268          objects[i] = choice[i].possibleBranch;
2269        }
2270        int whichObject = decision->bestBranch(objects,numberStrong,numberUnsatisfied_,
2271                                               changeUp,numberInfeasibilitiesUp,
2272                                               changeDown,numberInfeasibilitiesDown,
2273                                               objectiveValue_);
2274        // move branching object and make sure it will not be deleted
2275        if (whichObject>=0) {
2276          branch_ = objects[whichObject];
2277          if (model->messageHandler()->logLevel()>3) 
2278            printf("Choosing column %d\n",choice[whichObject].possibleBranch->variable()) ;
2279          choice[whichObject].possibleBranch=NULL;
2280        }
2281        delete [] changeUp;
2282        delete [] numberInfeasibilitiesUp;
2283        delete [] changeDown;
2284        delete [] numberInfeasibilitiesDown;
2285        delete [] objects;
2286      }
2287#     ifdef COIN_HAS_CLP
2288      if (osiclp&&!allNormal) {
2289        // back to normal
2290        osiclp->setHintParam(OsiDoInBranchAndCut,true,OsiHintDo,NULL) ;
2291      }
2292#     endif 
2293    }
2294    /*
2295      Simple branching. Probably just one, but we may have got here
2296      because of an odd branch e.g. a cut
2297    */
2298    else {
2299      // not strong
2300      // C) create branching object
2301      branch_ = choice[bestChoice].possibleBranch;
2302      choice[bestChoice].possibleBranch=NULL;
2303    }
2304  }
2305  // Set guessed solution value
2306  guessedObjectiveValue_ = objectiveValue_+estimatedDegradation;
2307/*
2308  Cleanup, then we're outta here.
2309*/
2310  if (!model->branchingMethod()||dynamicBranchingObject)
2311    delete decision;
2312   
2313  for (i=0;i<maximumStrong;i++)
2314    delete choice[i].possibleBranch;
2315  delete [] choice;
2316  delete [] saveLower;
2317  delete [] saveUpper;
2318 
2319  // restore solution
2320  solver->setColSolution(saveSolution);
2321  delete [] saveSolution;
2322  return anyAction;
2323}
2324
2325/*
2326  Version for dynamic pseudo costs.
2327
2328  **** For now just return if anything odd
2329  later allow even if odd
2330
2331  The routine scans through the object list of the model looking for objects
2332  that indicate infeasibility. It tests each object using strong branching
2333  and selects the one with the least objective degradation.  A corresponding
2334  branching object is left attached to lastNode.
2335  This version gives preference in evaluation to variables which
2336  have not been evaluated many times.  It also uses numberStrong
2337  to say give up if last few tries have not changed incumbent.
2338  See Achterberg, Koch and Martin.
2339
2340  If strong branching is disabled, a candidate object is chosen essentially
2341  at random (whatever object ends up in pos'n 0 of the candidate array).
2342
2343  If a branching candidate is found to be monotone, bounds are set to fix the
2344  variable and the routine immediately returns (the caller is expected to
2345  reoptimize).
2346
2347  If a branching candidate is found to result in infeasibility in both
2348  directions, the routine immediately returns an indication of infeasibility.
2349
2350  Returns:  0   both branch directions are feasible
2351           -1   branching variable is monotone
2352           -2   infeasible
2353           -3   Use another method
2354
2355           For now just fix on objective from strong branching.
2356*/
2357
2358int CbcNode::chooseDynamicBranch (CbcModel *model, CbcNode *lastNode,
2359                                  OsiSolverBranch * & branches,int numberPassesLeft)
2360
2361{ if (lastNode)
2362    depth_ = lastNode->depth_+1;
2363  else
2364    depth_ = 0;
2365  // Go to other choose if hot start
2366  if (model->hotstartSolution()) 
2367    return -3;
2368  delete branch_;
2369  branch_=NULL;
2370  OsiSolverInterface * solver = model->solver();
2371  // get information on solver type
2372  const OsiAuxInfo * auxInfo = solver->getAuxiliaryInfo();
2373  const OsiBabSolver * auxiliaryInfo = dynamic_cast<const OsiBabSolver *> (auxInfo);
2374  if (!auxiliaryInfo) {
2375    // use one from CbcModel
2376    auxiliaryInfo = model->solverCharacteristics();
2377  }
2378  int numberObjects = model->numberObjects();
2379  bool checkFeasibility = numberObjects>model->numberIntegers();
2380  if ((model->specialOptions()&128)!=0)
2381    checkFeasibility=false; // allow
2382  // For now return if not simple
2383  if (checkFeasibility)
2384    return -3;
2385  // point to useful information
2386  OsiBranchingInformation usefulInfo = model->usefulInformation();
2387  // and modify
2388  usefulInfo.depth_=depth_;
2389  if ((model->specialOptions()&128)!=0) {
2390    // SOS - shadow prices
2391    int numberRows = solver->getNumRows();
2392    const double * pi = usefulInfo.pi_;
2393    double sumPi=0.0;
2394    for (int i=0;i<numberRows;i++) 
2395      sumPi += fabs(pi[i]);
2396    sumPi /= ((double) numberRows);
2397    // and scale back
2398    sumPi *= 0.01;
2399    usefulInfo.defaultDual_ = sumPi; // switch on
2400    int numberColumns = solver->getNumCols();
2401    int size = CoinMax(numberColumns,2*numberRows);
2402    usefulInfo.usefulRegion_ = new double [size];
2403    CoinZeroN(usefulInfo.usefulRegion_,size);
2404    usefulInfo.indexRegion_ = new int [size];
2405    // pi may change
2406    usefulInfo.pi_=CoinCopyOfArray(usefulInfo.pi_,numberRows);
2407  }
2408  assert (auxiliaryInfo);
2409  //assert(objectiveValue_ == solver->getObjSense()*solver->getObjValue());
2410  double cutoff =model->getCutoff();
2411  double distanceToCutoff=cutoff-objectiveValue_;
2412  const double * lower = solver->getColLower();
2413  const double * upper = solver->getColUpper();
2414  // See what user thinks
2415  int anyAction=model->problemFeasibility()->feasible(model,0);
2416  if (anyAction) {
2417    // will return -2 if infeasible , 0 if treat as integer
2418    return anyAction-1;
2419  }
2420  int i;
2421  int saveStateOfSearch = model->stateOfSearch();
2422  int numberStrong=model->numberStrong();
2423  if (!auxiliaryInfo->warmStart())
2424    numberStrong=0;
2425  // But make more likely to get out after some times
2426  int changeStrategy=numberStrong;
2427  double changeFactor=1.0;
2428  // Use minimum of this and one stored in objects
2429  //int numberBeforeTrust = model->numberBeforeTrust();
2430  // Return if doing hot start (in BAB sense)
2431  if (model->hotstartSolution()) 
2432    return -3;
2433  //#define RANGING
2434#ifdef RANGING
2435  // Pass number
2436  int kPass=0;
2437  int numberRows = solver->getNumRows();
2438#endif
2439  int numberColumns = model->getNumCols();
2440  double * saveUpper = new double[numberColumns];
2441  double * saveLower = new double[numberColumns];
2442
2443  // Save solution in case heuristics need good solution later
2444 
2445  double * saveSolution = new double[numberColumns];
2446  memcpy(saveSolution,solver->getColSolution(),numberColumns*sizeof(double));
2447  model->reserveCurrentSolution(saveSolution);
2448  /*
2449    Get a branching decision object. Use the default dynamic decision criteria unless
2450    the user has loaded a decision method into the model.
2451  */
2452  CbcBranchDecision *decision = model->branchingMethod();
2453  if (!decision)
2454    decision = new CbcBranchDynamicDecision();
2455  int numberMini=0;
2456  int xPen=0;
2457  int xMark=0;
2458  for (i=0;i<numberColumns;i++) {
2459    saveLower[i] = lower[i];
2460    saveUpper[i] = upper[i];
2461  }
2462  // Get arrays to sort
2463  double * sort = new double[numberObjects];
2464  int * whichObject = new int[numberObjects];
2465  int * objectMark = new int[2*numberObjects+1];
2466  // Arrays with movements
2467  double * upEstimate = new double[numberObjects];
2468  double * downEstimate = new double[numberObjects];
2469  CbcStrongInfo * fixObject = new CbcStrongInfo[numberObjects];
2470  double estimatedDegradation=0.0; 
2471  int numberNodes=model->getNodeCount();
2472  int saveLogLevel = model->logLevel();
2473  if ((numberNodes%500)==0&&false) {
2474    model->setLogLevel(6);
2475    // Get average up and down costs
2476    double averageUp=0.0;
2477    double averageDown=0.0;
2478    int numberUp=0;
2479    int numberDown=0;
2480    int i;
2481    for ( i=0;i<numberObjects;i++) {
2482      OsiObject * object = model->modifiableObject(i);
2483      CbcSimpleIntegerDynamicPseudoCost * dynamicObject =
2484        dynamic_cast <CbcSimpleIntegerDynamicPseudoCost *>(object) ;
2485      assert(dynamicObject);
2486      int  numberUp2=0;
2487      int numberDown2=0;
2488      double up=0.0;
2489      double down=0.0;
2490      if (dynamicObject->numberTimesUp()) {
2491        numberUp++;
2492        averageUp += dynamicObject->upDynamicPseudoCost();
2493        numberUp2 += dynamicObject->numberTimesUp();
2494        up = dynamicObject->upDynamicPseudoCost();
2495      }
2496      if (dynamicObject->numberTimesDown()) {
2497        numberDown++;
2498        averageDown += dynamicObject->downDynamicPseudoCost();
2499        numberDown2 += dynamicObject->numberTimesDown();
2500        down = dynamicObject->downDynamicPseudoCost();
2501      }
2502      if (numberUp2||numberDown2)
2503        printf("col %d - up %d times cost %g, - down %d times cost %g\n",
2504               dynamicObject->columnNumber(),numberUp2,up,numberDown2,down);
2505    }
2506    if (numberUp) 
2507      averageUp /= (double) numberUp;
2508    else
2509      averageUp=1.0;
2510    if (numberDown) 
2511      averageDown /= (double) numberDown;
2512    else
2513      averageDown=1.0;
2514    printf("total - up %d vars average %g, - down %d vars average %g\n",
2515           numberUp,averageUp,numberDown,averageDown);
2516  }
2517  int numberBeforeTrust = model->numberBeforeTrust();
2518  int numberPenalties = model->numberPenalties();
2519  if (numberBeforeTrust>=1000000) {
2520    numberBeforeTrust = numberBeforeTrust % 1000000;
2521    numberPenalties=0;
2522  } else if (numberBeforeTrust<0) {
2523    if (numberBeforeTrust==-1)
2524      numberPenalties=numberColumns;
2525    else if (numberBeforeTrust==-2)
2526      numberPenalties=0;
2527    numberBeforeTrust=0;
2528  }
2529  // May go round twice if strong branching fixes all local candidates
2530  bool finished=false;
2531  int numberToFix=0;
2532# ifdef COIN_HAS_CLP
2533  OsiClpSolverInterface * osiclp = dynamic_cast< OsiClpSolverInterface*> (solver);
2534  int saveClpOptions=0;
2535  if (osiclp) {
2536    // for faster hot start
2537    saveClpOptions = osiclp->specialOptions();
2538    osiclp->setSpecialOptions(saveClpOptions|8192);
2539  }
2540# else
2541  OsiSolverInterface *osiclp = 0 ;
2542# endif
2543  const CglTreeProbingInfo * probingInfo = model->probingInfo();
2544  int saveSearchStrategy2 = model->searchStrategy();
2545#define NODE_NEW  2
2546#ifdef RANGING
2547  bool useRanging;
2548#if NODE_NEW !=3
2549  useRanging = false;
2550#else
2551  useRanging = true;
2552#endif
2553  if (saveSearchStrategy2!=2)
2554    useRanging=false;
2555#endif
2556  if (saveSearchStrategy2<999) {
2557#if NODE_NEW ==4
2558    if (saveSearchStrategy2!=2) {
2559#endif
2560    // Get average up and down costs
2561    double averageUp=0.0;
2562    double averageDown=0.0;
2563    {
2564      int numberUp=0;
2565      int numberDown=0;
2566      int i;
2567      for ( i=0;i<numberObjects;i++) {
2568        OsiObject * object = model->modifiableObject(i);
2569        CbcSimpleIntegerDynamicPseudoCost * dynamicObject =
2570          dynamic_cast <CbcSimpleIntegerDynamicPseudoCost *>(object) ;
2571        if (dynamicObject) {
2572          if (dynamicObject->numberTimesUp()) {
2573            numberUp++;
2574            averageUp += dynamicObject->upDynamicPseudoCost();
2575          }
2576          if (dynamicObject->numberTimesDown()) {
2577            numberDown++;
2578            averageDown += dynamicObject->downDynamicPseudoCost();
2579          }
2580        }
2581      }
2582      if (numberUp) 
2583        averageUp /= (double) numberUp;
2584      else
2585        averageUp=1.0;
2586      if (numberDown) 
2587        averageDown /= (double) numberDown;
2588      else
2589        averageDown=1.0;
2590      for ( i=0;i<numberObjects;i++) {
2591        OsiObject * object = model->modifiableObject(i);
2592        CbcSimpleIntegerDynamicPseudoCost * dynamicObject =
2593          dynamic_cast <CbcSimpleIntegerDynamicPseudoCost *>(object) ;
2594        if (dynamicObject) {
2595          if (!dynamicObject->numberTimesUp()) 
2596            dynamicObject->setUpDynamicPseudoCost(averageUp);
2597          if (!dynamicObject->numberTimesDown()) 
2598            dynamicObject->setDownDynamicPseudoCost(averageDown);
2599        }
2600      }
2601    }
2602#if NODE_NEW ==4
2603    } else {
2604    // pseudo shadow prices
2605    model->pseudoShadow(NULL,NULL);
2606    }
2607#endif
2608  } else if (saveSearchStrategy2<1999) {
2609    // pseudo shadow prices
2610    model->pseudoShadow(NULL,NULL);
2611  } else if (saveSearchStrategy2<2999) {
2612    // leave old ones
2613  } else if (saveSearchStrategy2<3999) {
2614    // pseudo shadow prices at root
2615    if (!numberNodes)
2616      model->pseudoShadow(NULL,NULL);
2617  } else {
2618    abort();
2619  }
2620  if (saveSearchStrategy2>=0)
2621    saveSearchStrategy2 = saveSearchStrategy2 % 1000;
2622  if (saveSearchStrategy2==999)
2623    saveSearchStrategy2=-1;
2624  int px[4]={-1,-1,-1,-1};
2625  int saveSearchStrategy = saveSearchStrategy2<99 ? saveSearchStrategy2 : saveSearchStrategy2-100;
2626  bool newWay = saveSearchStrategy2>98;
2627  int numberNotTrusted=0;
2628  int numberStrongDone=0;
2629  int numberUnfinished=0;
2630  int numberStrongInfeasible=0;
2631  int numberStrongIterations=0;
2632  // so we can save lots of news
2633  CbcStrongInfo choice;
2634  CbcDynamicPseudoCostBranchingObject * choiceObject = NULL;
2635  if (model->allDynamic()) {
2636    CbcSimpleIntegerDynamicPseudoCost * object = NULL;
2637    choiceObject=new CbcDynamicPseudoCostBranchingObject(model,0,-1,0.5,object);
2638  }
2639  choice.possibleBranch=choiceObject;
2640  while(!finished) {
2641    finished=true;
2642    decision->initialize(model);
2643    // Some objects may compute an estimate of best solution from here
2644    estimatedDegradation=0.0; 
2645    numberToFix=0;
2646    int numberIntegerInfeasibilities=0; // without odd ones
2647    int numberToDo=0;
2648    int iBestNot=-1;
2649    int iBestGot=-1;
2650    double best=0.0;
2651    numberNotTrusted=0;
2652    numberStrongDone=0;
2653    numberUnfinished=0;
2654    numberStrongInfeasible=0;
2655    numberStrongIterations=0;
2656    int * which = objectMark+numberObjects+1;
2657    int neededPenalties;
2658    int branchingMethod=-1;
2659    // We may go round this loop three times (only if we think we have solution)
2660    for (int iPass=0;iPass<3;iPass++) {
2661     
2662      // compute current state
2663      int numberObjectInfeasibilities; // just odd ones
2664      model->feasibleSolution(
2665                              numberIntegerInfeasibilities,
2666                              numberObjectInfeasibilities);
2667     
2668      // Some objects may compute an estimate of best solution from here
2669      estimatedDegradation=0.0; 
2670      numberUnsatisfied_ = 0;
2671      // initialize sum of "infeasibilities"
2672      sumInfeasibilities_ = 0.0;
2673      int bestPriority=COIN_INT_MAX;
2674      int number01 = 0;
2675      const fixEntry * entry = NULL;
2676      const int * toZero = NULL;
2677      const int * toOne = NULL;
2678      const int * backward = NULL;
2679      int numberUnsatisProbed=0;
2680      int numberUnsatisNotProbed=0; // 0-1
2681      if (probingInfo) {
2682        number01 = probingInfo->numberIntegers();
2683        entry = probingInfo->fixEntries();
2684        toZero = probingInfo->toZero();
2685        toOne = probingInfo->toOne();
2686        backward = probingInfo->backward();
2687        if (!toZero[number01]||number01<numberObjects||true) {
2688          // no info
2689          probingInfo=NULL;
2690        }
2691      }
2692      /*
2693        Scan for branching objects that indicate infeasibility. Choose candidates
2694        using priority as the first criteria, then integer infeasibility.
2695       
2696        The algorithm is to fill the array with a set of good candidates (by
2697        infeasibility) with priority bestPriority.  Finding a candidate with
2698        priority better (less) than bestPriority flushes the choice array. (This
2699        serves as initialization when the first candidate is found.)
2700       
2701      */
2702      numberToDo=0;
2703      neededPenalties=0;
2704      iBestNot=-1;
2705      double bestNot=0.0;
2706      iBestGot=-1;
2707      best=0.0;
2708      /* Problem type as set by user or found by analysis.  This will be extended
2709         0 - not known
2710         1 - Set partitioning <=
2711         2 - Set partitioning ==
2712         3 - Set covering
2713         4 - all +- 1 or all +1 and odd
2714      */
2715      int problemType = model->problemType();
2716#define PRINT_STUFF -1
2717      for (i=0;i<numberObjects;i++) {
2718        OsiObject * object = model->modifiableObject(i);
2719        CbcSimpleIntegerDynamicPseudoCost * dynamicObject =
2720          dynamic_cast <CbcSimpleIntegerDynamicPseudoCost *>(object) ;
2721        int preferredWay;
2722        double infeasibility = object->infeasibility(&usefulInfo,preferredWay);
2723        int priorityLevel = object->priority();
2724#define ZERO_ONE 0
2725#define ZERO_FAKE 1.0e20;
2726#if ZERO_ONE==1
2727        // branch on 0-1 first (temp)
2728        if (fabs(saveSolution[dynamicObject->columnNumber()])<1.0)
2729          priorityLevel--;
2730#endif
2731#if ZERO_ONE==2
2732        if (fabs(saveSolution[dynamicObject->columnNumber()])<1.0)
2733          infeasibility *= ZERO_FAKE;
2734#endif
2735        if (infeasibility) {
2736          int iColumn = numberColumns+i;
2737          bool gotDown=false;
2738          int numberThisDown = 0;
2739          bool gotUp=false;
2740          int numberThisUp = 0;
2741          // check branching method
2742          if (dynamicObject) {
2743            if (branchingMethod!=dynamicObject->method()) {
2744              if (branchingMethod==-1)
2745                branchingMethod = dynamicObject->method();
2746              else
2747                branchingMethod = 100;
2748            }
2749            iColumn = dynamicObject->columnNumber();
2750            //double gap = saveUpper[iColumn]-saveLower[iColumn];
2751            // Give precedence to ones with gap of 1.0
2752            //assert(gap>0.0);
2753            //infeasibility /= CoinMin(gap,100.0);
2754            if (!depth_&&false) {
2755              // try closest to 0.5
2756              double part =saveSolution[iColumn]-floor(saveSolution[iColumn]);
2757              infeasibility = fabs(0.5-part);
2758            }
2759            if (problemType>0&&problemType<4&&false) {
2760              // try closest to 0.5
2761              double part =saveSolution[iColumn]-floor(saveSolution[iColumn]);
2762              infeasibility = 0.5-fabs(0.5-part);
2763            }
2764            if (probingInfo) {
2765              int iSeq = backward[iColumn];
2766              assert (iSeq>=0);
2767              infeasibility = 1.0 + (toZero[iSeq+1]-toZero[iSeq])+ 
2768                5.0*CoinMin(toOne[iSeq]-toZero[iSeq],toZero[iSeq+1]-toOne[iSeq]);
2769              if (toZero[iSeq+1]>toZero[iSeq]) {
2770                numberUnsatisProbed++;
2771              } else {
2772                numberUnsatisNotProbed++;
2773              }
2774            }
2775            gotDown=false;
2776            numberThisDown = dynamicObject->numberTimesDown();
2777            if (numberThisDown>=numberBeforeTrust)
2778              gotDown=true;
2779            gotUp=false;
2780            numberThisUp = dynamicObject->numberTimesUp();
2781            if (numberThisUp>=numberBeforeTrust)
2782              gotUp=true;
2783            if ((numberNodes%PRINT_STUFF)==0&&PRINT_STUFF>0)
2784              printf("%d down %d %g up %d %g - infeas %g\n",
2785                     i,numberThisDown,object->downEstimate(),numberThisUp,object->upEstimate(),
2786                     infeasibility);
2787          } else {
2788            // see if SOS
2789            CbcSOS * sosObject =
2790              dynamic_cast <CbcSOS *>(object) ;
2791            assert (sosObject);
2792            gotDown=false;
2793            numberThisDown = sosObject->numberTimesDown();
2794            if (numberThisDown>=numberBeforeTrust)
2795              gotDown=true;
2796            gotUp=false;
2797            numberThisUp = sosObject->numberTimesUp();
2798            if (numberThisUp>=numberBeforeTrust)
2799              gotUp=true;
2800          }
2801          // Increase estimated degradation to solution
2802          estimatedDegradation += CoinMin(object->upEstimate(),object->downEstimate());
2803          downEstimate[i]=object->downEstimate();
2804          upEstimate[i]=object->upEstimate();
2805          numberUnsatisfied_++;
2806          sumInfeasibilities_ += infeasibility;
2807          // Better priority? Flush choices.
2808          if (priorityLevel<bestPriority) {
2809            numberToDo=0;
2810            bestPriority = priorityLevel;
2811            iBestGot=-1;
2812            best=0.0;
2813            numberNotTrusted=0;
2814          } else if (priorityLevel>bestPriority) {
2815            continue;
2816          }
2817          if (!gotUp||!gotDown)
2818            numberNotTrusted++;
2819          // Check for suitability based on infeasibility.
2820          if ((gotDown&&gotUp)&&numberStrong>0) {
2821            sort[numberToDo]=-infeasibility;
2822            if (infeasibility>best) {
2823              best=infeasibility;
2824              iBestGot=numberToDo;
2825            }
2826          } else {
2827            objectMark[neededPenalties]=numberToDo;
2828            which[neededPenalties++]=iColumn;
2829            sort[numberToDo]=-10.0*infeasibility;
2830            if (!(numberThisUp+numberThisDown))
2831              sort[numberToDo] *= 100.0; // make even more likely
2832            if (iColumn<numberColumns) {
2833              double part =saveSolution[iColumn]-floor(saveSolution[iColumn]);
2834              if (1.0-fabs(part-0.5)>bestNot) {
2835                iBestNot=numberToDo;
2836                bestNot = 1.0-fabs(part-0.5);
2837              }
2838            } else {
2839              // SOS
2840              if (-sort[numberToDo]>bestNot) {
2841                iBestNot=numberToDo;
2842                bestNot = -sort[numberToDo];
2843              }
2844            }
2845          }
2846          if (model->messageHandler()->logLevel()>3) { 
2847            printf("%d (%d) down %d %g up %d %g - infeas %g - sort %g solution %g\n",
2848                   i,iColumn,numberThisDown,object->downEstimate(),numberThisUp,object->upEstimate(),
2849                   infeasibility,sort[numberToDo],saveSolution[iColumn]);
2850          }
2851          whichObject[numberToDo++]=i;
2852        } else {
2853          // for debug
2854          downEstimate[i]=-1.0;
2855          upEstimate[i]=-1.0;
2856        }
2857      }
2858      if (numberUnsatisfied_) {
2859        if (probingInfo&&false)
2860          printf("nunsat %d, %d probed, %d other 0-1\n",numberUnsatisfied_,
2861                 numberUnsatisProbed,numberUnsatisNotProbed);
2862        // some infeasibilities - go to next steps
2863        break;
2864      } else if (!iPass) {
2865        // may just need resolve
2866        solver->resolve();
2867        memcpy(saveSolution,solver->getColSolution(),numberColumns*sizeof(double));
2868        model->reserveCurrentSolution(saveSolution);
2869        if (!solver->isProvenOptimal()) {
2870          // infeasible
2871          anyAction=-2;
2872          break;
2873        }
2874      } else if (iPass==1) {
2875        // looks like a solution - get paranoid
2876        bool roundAgain=false;
2877        // get basis
2878        CoinWarmStartBasis * ws = dynamic_cast<CoinWarmStartBasis*>(solver->getWarmStart());
2879        if (!ws)
2880          break;
2881        double tolerance;
2882        solver->getDblParam(OsiPrimalTolerance,tolerance);
2883        for (i=0;i<numberColumns;i++) {
2884          double value = saveSolution[i];
2885          if (value<lower[i]-tolerance) {
2886            saveSolution[i]=lower[i];
2887            roundAgain=true;
2888            ws->setStructStatus(i,CoinWarmStartBasis::atLowerBound);
2889          } else if (value>upper[i]+tolerance) {
2890            saveSolution[i]=upper[i];
2891            roundAgain=true;
2892            ws->setStructStatus(i,CoinWarmStartBasis::atUpperBound);
2893          } 
2894        }
2895        if (roundAgain) {
2896          // restore basis
2897          solver->setWarmStart(ws);
2898          solver->setColSolution(saveSolution);
2899          delete ws;
2900          bool takeHint;
2901          OsiHintStrength strength;
2902          solver->getHintParam(OsiDoDualInResolve,takeHint,strength);
2903          solver->setHintParam(OsiDoDualInResolve,false,OsiHintDo) ;
2904          solver->resolve();
2905          solver->setHintParam(OsiDoDualInResolve,takeHint,strength) ;
2906          memcpy(saveSolution,solver->getColSolution(),numberColumns*sizeof(double));
2907          model->reserveCurrentSolution(saveSolution);
2908          if (!solver->isProvenOptimal()) {
2909            // infeasible
2910            anyAction=-2;
2911            break;
2912          }
2913        } else {
2914          delete ws;
2915          break;
2916        }
2917      }
2918    }
2919    if (anyAction==-2) {
2920      break;
2921    }
2922    bool solveAll=false; // set true to say look at all even if some fixed (experiment)
2923    solveAll=true;
2924    // skip if solution
2925    if (!numberUnsatisfied_)
2926      break;
2927    //bool skipAll = (numberBeforeTrust>20&&numberNodes>20000&&numberNotTrusted==0);
2928    bool skipAll = numberNotTrusted==0||numberToDo==1;
2929    bool doneHotStart=false;
2930    int searchStrategy = saveSearchStrategy>=0 ? (saveSearchStrategy%10) : -1;
2931    if (0) {
2932      int numberIterations = model->getIterationCount();
2933      int numberStrongIterations = model->numberStrongIterations();
2934      printf("node %d saveSearch %d search %d - its %d strong %d\n",
2935             numberNodes,saveSearchStrategy,searchStrategy,
2936             numberIterations,numberStrongIterations);
2937    }
2938#ifndef CBC_WEAK_STRONG
2939    if (((numberNodes%20)==0&&searchStrategy!=2)||(model->specialOptions()&8)!=0)
2940      skipAll=false;
2941#endif
2942    if (!newWay) {
2943    // 10 up always use %10, 20 up as 10 and allow penalties
2944    // But adjust depending on ratio of iterations
2945    if (searchStrategy>0&&saveSearchStrategy<10) {
2946      if (numberBeforeTrust>=5&&numberBeforeTrust<=10) {
2947        if (searchStrategy!=2) {
2948          if (depth_>5) {
2949            int numberIterations = model->getIterationCount();
2950            int numberStrongIterations = model->numberStrongIterations();
2951            if (numberStrongIterations>numberIterations+10000) {
2952              searchStrategy=2;
2953              //skipAll=true;
2954            } else if (numberStrongIterations*4+1000<numberIterations||depth_<5) {
2955              searchStrategy=3;
2956              skipAll=false;
2957            }
2958          } else {
2959            searchStrategy=3;
2960            skipAll=false;
2961          }
2962        } else {
2963          //skipAll=true;
2964        }
2965      }
2966    }
2967    } else {
2968    // But adjust depending on ratio of iterations
2969    if (saveSearchStrategy<0) {
2970      // unset
2971      if ((numberNodes%20)==0||(model->specialOptions()&8)!=0) {
2972        // Do numberStrong
2973        searchStrategy=3;
2974      } else if (depth_<5) {
2975        // Do numberStrong
2976        searchStrategy=2;
2977      } else {
2978        int numberIterations = model->getIterationCount();
2979        int numberStrongIterations = model->numberStrongIterations();
2980        int numberRows = solver->getNumRows();
2981        if (numberStrongIterations>numberIterations+CoinMin(10000,10*numberRows)) {
2982          // off
2983          searchStrategy=0;
2984        } else if (numberStrongIterations*4+1000<numberIterations) {
2985          // Do numberStrong if not trusted
2986          searchStrategy=2;
2987        } else {
2988          searchStrategy=1;
2989        }
2990      }
2991    }
2992    if (searchStrategy<3&&(!numberNotTrusted||!searchStrategy))
2993      skipAll=true;
2994    else
2995      skipAll=false;
2996    }
2997    // worth trying if too many times
2998    // Save basis
2999    CoinWarmStart * ws = NULL;
3000    // save limit
3001    int saveLimit=0;
3002    solver->getIntParam(OsiMaxNumIterationHotStart,saveLimit);
3003    if (!skipAll) {
3004      ws = solver->getWarmStart();
3005      int limit=100;
3006#if 0
3007      int averageBranchIterations = model->getIterationCount()/(model->getNodeCount()+1);
3008      if (numberNodes)
3009        limit = CoinMin(CoinMax(limit,2*averageBranchIterations),500);
3010      else
3011        limit = 500;
3012#endif
3013      if ((!saveStateOfSearch||searchStrategy>3)&&saveLimit<limit&&saveLimit==100)
3014        solver->setIntParam(OsiMaxNumIterationHotStart,limit); 
3015    }
3016    // Say which one will be best
3017    int whichChoice=0;
3018    int bestChoice;
3019    if (iBestGot>=0)
3020      bestChoice=iBestGot;
3021    else
3022      bestChoice=iBestNot;
3023    assert (bestChoice>=0);
3024    // If we have hit max time don't do strong branching
3025    bool hitMaxTime = ( CoinCpuTime()-model->getDblParam(CbcModel::CbcStartSeconds) > 
3026                        model->getDblParam(CbcModel::CbcMaximumSeconds));
3027    // also give up if we are looping round too much
3028    if (hitMaxTime||numberPassesLeft<=0||(!numberNotTrusted&&false)||branchingMethod==11) {
3029      int iObject = whichObject[bestChoice];
3030      OsiObject * object = model->modifiableObject(iObject);
3031      int preferredWay;
3032      object->infeasibility(&usefulInfo,preferredWay);
3033      CbcSimpleInteger * obj =
3034        dynamic_cast <CbcSimpleInteger *>(object) ;
3035      if (obj) {
3036        branch_=obj->createBranch(solver,&usefulInfo,preferredWay);
3037      } else {
3038        CbcObject * obj =
3039          dynamic_cast <CbcObject *>(object) ;
3040        assert (obj);
3041        branch_=obj->createBranch(preferredWay);
3042      }
3043      {
3044        CbcBranchingObject * branchObj =
3045          dynamic_cast <CbcBranchingObject *>(branch_) ;
3046        assert (branchObj);
3047        branchObj->way(preferredWay);
3048      }
3049      delete ws;
3050      ws=NULL;
3051      break;
3052    } else {
3053      // say do fast
3054      int easy=1;
3055      if (!skipAll)
3056        solver->setHintParam(OsiDoInBranchAndCut,true,OsiHintDo,&easy) ;
3057      int iDo;
3058#ifdef RANGING
3059      if (useRanging) {
3060      if ((skipAll&&numberBeforeTrust&&saveSearchStrategy<20)||saveSearchStrategy<10)
3061        numberPenalties=0;
3062      {
3063        // off penalties if too much
3064        double needed = neededPenalties;
3065        needed *= numberRows;
3066        if (needed>1.0e6&&numberNodes&&saveSearchStrategy<20) {
3067          numberPenalties=0;
3068          neededPenalties=0;
3069        }
3070      }
3071#     ifdef COIN_HAS_CLP
3072      if (osiclp&&numberPenalties&&neededPenalties) {
3073        xPen += neededPenalties;
3074        which--;
3075        which[0]=neededPenalties;
3076        osiclp->passInRanges(which);
3077        // Mark hot start and get ranges
3078        if (kPass) {
3079          // until can work out why solution can go funny
3080          int save = osiclp->specialOptions();
3081          osiclp->setSpecialOptions(save|256);
3082          solver->markHotStart();
3083          osiclp->setSpecialOptions(save);
3084        } else {
3085          solver->markHotStart();
3086        }
3087        assert (auxiliaryInfo->warmStart());
3088        doneHotStart=true;
3089        xMark++;
3090        kPass++;
3091        osiclp->passInRanges(NULL);
3092        const double * downCost=osiclp->upRange();
3093        const double * upCost=osiclp->downRange();
3094        //printf("numberTodo %d needed %d numberpenalties %d\n",numberToDo,neededPenalties,numberPenalties);
3095        double invTrust = 1.0/((double) numberBeforeTrust);
3096        for (int i=0;i<neededPenalties;i++) {
3097          int j = objectMark[i];
3098          int iObject = whichObject[j];
3099          OsiObject * object = model->modifiableObject(iObject);
3100          CbcSimpleIntegerDynamicPseudoCost * dynamicObject =
3101            dynamic_cast <CbcSimpleIntegerDynamicPseudoCost *>(object) ;
3102          int iSequence=dynamicObject->columnNumber();
3103          double value = saveSolution[iSequence];
3104          value -= floor(value);
3105          double upPenalty = CoinMin(upCost[i],1.0e110)*(1.0-value);
3106          double downPenalty = CoinMin(downCost[i],1.0e110)*value;
3107          if (!numberBeforeTrust) {
3108            // override
3109            downEstimate[iObject]=downPenalty;
3110            upEstimate[iObject]=upPenalty;
3111          } else {
3112            int numberThisDown = dynamicObject->numberTimesDown();
3113            if (numberThisDown<numberBeforeTrust) {
3114              double fraction = ((double) numberThisDown)*invTrust;
3115              downEstimate[iObject] = fraction*downEstimate[iObject]+(1.0-fraction)*downPenalty;
3116            }
3117            int numberThisUp = dynamicObject->numberTimesUp();
3118            if (numberThisUp<numberBeforeTrust) {
3119              double fraction = ((double) numberThisUp)*invTrust;
3120              upEstimate[iObject] = fraction*upEstimate[iObject]+(1.0-fraction)*upPenalty;
3121            }
3122          }
3123          sort[j] = - CoinMin(downEstimate[iObject],upEstimate[iObject]);
3124#ifdef CBC_WEAK_STRONG
3125          sort[j] -= 1.0e10; // make more likely to be chosen
3126#endif
3127          //if ((numberNodes%PRINT_STUFF)==0&&PRINT_STUFF>0)
3128          if (!numberNodes)
3129            printf("%d pen down ps %g -> %g up ps %g -> %g\n",
3130                   iObject,downCost[i],downPenalty,upCost[i],upPenalty);
3131        }
3132      } else
3133#     endif     /* COIN_HAS_CLP */
3134      {
3135        if (!skipAll) {
3136          // Mark hot start
3137          solver->markHotStart();
3138          doneHotStart=true;
3139          assert (auxiliaryInfo->warmStart());
3140          xMark++;
3141          //if (solver->isProvenPrimalInfeasible())
3142          //printf("**** Hot start says node infeasible\n");
3143        }
3144        // make sure best will be first
3145        if (iBestGot>=0)
3146          sort[iBestGot]=-1.0e120;
3147      }
3148      } else {
3149#endif          /* RANGING */
3150      if (!skipAll) {
3151        // Mark hot start
3152        doneHotStart=true;
3153        assert (auxiliaryInfo->warmStart());
3154        solver->markHotStart();
3155        xMark++;
3156      }
3157      // make sure best will be first
3158      if (iBestGot>=0)
3159        sort[iBestGot]=-COIN_DBL_MAX;
3160#ifdef RANGING
3161      }
3162#endif          /* RANGING */
3163      // Actions 0 - exit for repeat, 1 resolve and try old choice,2 exit for continue
3164#define ACTION 0
3165#if ACTION<2
3166      if (anyAction)
3167        numberToDo=0; // skip as we will be trying again
3168#endif
3169      // Sort
3170      CoinSort_2(sort,sort+numberToDo,whichObject);
3171      // Change in objective opposite infeasible
3172      double worstFeasible=0.0;
3173      // Just first if strong off
3174      if (!numberStrong)
3175        numberToDo=CoinMin(numberToDo,1);
3176#if NODE_NEW
3177      if (searchStrategy==2)
3178        numberToDo=CoinMin(numberToDo,20);
3179#endif
3180      iDo=0;
3181      int saveLimit2;
3182      solver->getIntParam(OsiMaxNumIterationHotStart,saveLimit2);
3183      bool doQuickly = false; // numberToDo>2*numberStrong;
3184      if (searchStrategy==2)
3185        doQuickly=true;
3186      //printf("todo %d, strong %d\n",numberToDo,numberStrong);
3187      int numberTest=numberNotTrusted>0 ? numberStrong : (numberStrong+1)/2;
3188      int numberTest2 = 2*numberStrong;
3189      //double distanceToCutoff2 = model->getCutoff()-objectiveValue_;
3190      if (!newWay) {
3191      if (searchStrategy==3) {
3192        // Previously decided we need strong
3193        doQuickly=false;
3194        numberTest = numberStrong;
3195        //numberTest2 = 1000000;
3196      }
3197      //if (searchStrategy<0||searchStrategy==1)
3198        //numberTest2 = 1000000;
3199#if 0
3200      if (numberBeforeTrust>20&&(numberNodes>20000||(numberNodes>200&&numberNotTrusted==0))) {
3201        if ((numberNodes%20)!=0) {
3202          numberTest=0;
3203          doQuickly=true;
3204        }
3205      }
3206#else
3207      // Try nearly always off
3208      if (searchStrategy<2) {
3209        if ((numberNodes%20)!=0) {
3210          if ((model->specialOptions()&8)==0) {
3211            numberTest=0;
3212            doQuickly=true;
3213          }
3214        } else {
3215          doQuickly=false;
3216          numberTest=2*numberStrong;
3217          skipAll=false;
3218        }
3219      } else if (searchStrategy!=3) {
3220        doQuickly=true;
3221        numberTest=numberStrong;
3222      }
3223#endif
3224      if (depth_<8&&numberStrong) {
3225        if (searchStrategy!=2) {
3226          doQuickly=false;
3227          int numberRows = solver->getNumRows();
3228          // whether to do this or not is important - think
3229          if (numberRows<300||numberRows+numberColumns<2500) {
3230            if (depth_<7)
3231              numberStrong = CoinMin(3*numberStrong,numberToDo);
3232            if (!depth_) 
3233              numberStrong=CoinMin(6*numberStrong,numberToDo);
3234          }
3235          numberTest=numberStrong;
3236          skipAll=false;
3237        }
3238        //model->setStateOfSearch(2); // use min min
3239      }
3240      // could adjust using average iterations per branch
3241      // double average = ((double)model->getIterationCount())/
3242      //((double) model->getNodeCount()+1.0);
3243      // if too many and big then just do 10 its
3244      if (!skipAll&&saveStateOfSearch) {
3245        //if (numberNotTrusted>3*numberStrong&&numberRows>250&&numberColumns>1000&&saveLimit==100)
3246          // off solver->setIntParam(OsiMaxNumIterationHotStart,10);
3247      }
3248      // make negative for test
3249      distanceToCutoff = - distanceToCutoff;
3250      if (numberObjects>-100) {
3251        // larger
3252        distanceToCutoff *= 100.0;
3253      }
3254        distanceToCutoff = -COIN_DBL_MAX;
3255      // Do at least 5 strong
3256      if (numberColumns<1000&&(depth_<15||numberNodes<1000000))
3257        numberTest = CoinMax(numberTest,5);
3258      if ((model->specialOptions()&8)==0) {
3259        if (skipAll) {
3260          numberTest=0;
3261          doQuickly=true;
3262        }
3263      } else {
3264        // do 5 as strong is fixing
3265        numberTest = CoinMax(numberTest,5);
3266      }
3267      } else {
3268        if (!depth_&&numberToDo<200&&solver->getNumRows()<2000&&
3269            numberColumns<2000&&false) {
3270          numberStrong = numberToDo;
3271          doQuickly = false;
3272        }
3273      int numberTest=numberNotTrusted>0 ? numberStrong : (numberStrong+1)/2;
3274      int numberTest2 = 2*numberStrong;
3275      if (searchStrategy>=3) {
3276        // Previously decided we need strong
3277        doQuickly=false;
3278        if (depth_<7)
3279          numberStrong *=3;
3280        if (!depth_) 
3281          numberStrong=CoinMin(6*numberStrong,numberToDo);
3282        numberTest = numberStrong;
3283        numberTest2 *= 2;
3284      } else if (searchStrategy==2||(searchStrategy==1&&depth_<6)) {
3285        numberStrong *=2;
3286        if (!depth_) 
3287          numberStrong=CoinMin(2*numberStrong,numberToDo);
3288        numberTest = numberStrong;
3289      } else if (searchStrategy==1&&numberNotTrusted) {
3290        numberTest = numberStrong;
3291      } else {
3292        numberTest=0;
3293        skipAll=true;
3294      }
3295      distanceToCutoff=model->getCutoff()-objectiveValue_;
3296      // make negative for test
3297      distanceToCutoff = - distanceToCutoff;
3298      if (numberObjects>-100) {
3299        // larger
3300        distanceToCutoff *= 100.0;
3301      }
3302      distanceToCutoff = -COIN_DBL_MAX;
3303      if (skipAll) {
3304        numberTest=0;
3305        doQuickly=true;
3306      }
3307      }
3308      // temp - always switch off
3309      if (0) {
3310        int numberIterations = model->getIterationCount();
3311        int numberStrongIterations = model->numberStrongIterations();
3312        int numberRows = solver->getNumRows();
3313        if (numberStrongIterations>numberIterations+CoinMin(100,10*numberRows)&&depth_>=5) {
3314          skipAll=true;
3315          newWay=false;
3316          numberTest=0;
3317          doQuickly=true;
3318        }
3319      }
3320#if 0
3321      // temp - always switch on
3322      if (0) {
3323        int numberIterations = model->getIterationCount();
3324        int numberStrongIterations = model->numberStrongIterations();
3325        if (2*numberStrongIterations<numberIterations||depth_<=5) {
3326          skipAll=false;
3327          newWay=false;
3328          numberTest=CoinMax(numberTest,numberStrong);
3329          doQuickly=false;
3330        }
3331      }
3332#endif
3333      px[0]=numberTest;
3334      px[1]=numberTest2;
3335      px[2]= doQuickly ? 1 : -1;
3336      px[3]=numberStrong;
3337      if (!newWay) {
3338        if (numberColumns>8*solver->getNumRows()&&false) {
3339          printf("skipAll %c doQuickly %c numberTest %d numberTest2 %d numberNot %d\n",
3340                 skipAll ? 'Y' : 'N',doQuickly ? 'Y' : 'N',numberTest,numberTest2,numberNotTrusted);
3341          numberTest = CoinMin(numberTest,model->numberStrong());
3342          numberTest2 = CoinMin(numberTest2,model->numberStrong());
3343          printf("new test,test2 %d %d\n",numberTest,numberTest2);
3344        }
3345      }
3346      //printf("skipAll %c doQuickly %c numberTest %d numberTest2 %d numberNot %d\n",
3347      //   skipAll ? 'Y' : 'N',doQuickly ? 'Y' : 'N',numberTest,numberTest2,numberNotTrusted);
3348      // See if we want mini tree
3349      bool wantMiniTree=false;
3350      if (model->sizeMiniTree()&&depth_>7&&saveStateOfSearch>0)
3351        wantMiniTree=true;
3352      numberMini=0;
3353      //if (skipAll&&numberTest==0&&doQuickly)
3354      //numberToDo = 1; // trust previous stuff
3355      bool couldChooseFirst = false ; //(skipAll&&numberTest==0&&doQuickly);
3356      //skipAll=false;
3357      int realMaxHotIterations=999999;
3358#if 0
3359      if (saveSearchStrategy==-1&&!model->parentModel()&&
3360          !depth_&&saveLimit==100) {
3361        realMaxHotIterations=saveLimit;
3362        saveLimit2=200;
3363        solver->setIntParam(OsiMaxNumIterationHotStart,saveLimit2);
3364      }
3365#endif
3366      for ( iDo=0;iDo<numberToDo;iDo++) {
3367        int iObject = whichObject[iDo];
3368        OsiObject * object = model->modifiableObject(iObject);
3369        CbcSimpleIntegerDynamicPseudoCost * dynamicObject =
3370          dynamic_cast <CbcSimpleIntegerDynamicPseudoCost *>(object) ;
3371        int iColumn = dynamicObject ? dynamicObject->columnNumber() : numberColumns+iObject;
3372        int preferredWay;
3373        double infeasibility = object->infeasibility(&usefulInfo,preferredWay);
3374        // may have become feasible
3375        if (!infeasibility)
3376          continue;
3377        CbcSimpleInteger * obj =
3378          dynamic_cast <CbcSimpleInteger *>(object) ;
3379        if (obj) {
3380          if (choiceObject) {
3381            obj->fillCreateBranch(choiceObject,&usefulInfo,preferredWay);
3382            choiceObject->setObject(dynamicObject);
3383          } else {
3384            choice.possibleBranch=obj->createBranch(solver,&usefulInfo,preferredWay);
3385          }
3386        } else {
3387          CbcObject * obj =
3388            dynamic_cast <CbcObject *>(object) ;
3389          assert (obj);
3390          choice.possibleBranch=obj->createBranch(preferredWay);
3391        }
3392        // Save which object it was
3393        choice.objectNumber=iObject;
3394        choice.numIntInfeasUp = numberUnsatisfied_;
3395        choice.numIntInfeasDown = numberUnsatisfied_;
3396        choice.upMovement = upEstimate[iObject];
3397        choice.downMovement = downEstimate[iObject];
3398        assert (choice.upMovement>=0.0);
3399        assert (choice.downMovement>=0.0);
3400        choice.fix=0; // say not fixed
3401        double maxChange = 0.5*(choice.upMovement+choice.downMovement);
3402        maxChange = CoinMin(choice.upMovement,choice.downMovement);
3403        maxChange = CoinMax(choice.upMovement,choice.downMovement);
3404        if (searchStrategy==2)
3405          maxChange = COIN_DBL_MAX;
3406        //maxChange *= 5.0;
3407        if (dynamicObject&&dynamicObject->method()==1)
3408          maxChange *= 0.1; // probing
3409        // see if can skip strong branching
3410        int canSkip = choice.possibleBranch->fillStrongInfo(choice);
3411#if 0
3412        if (!newWay) {
3413          if ((maxChange>distanceToCutoff2)&&(!doQuickly||(numberTest>0&&searchStrategy!=2)))
3414          canSkip=0;
3415        } else {
3416        if (skipAll)
3417          canSkip=1;
3418        else if (numberTest>0&&searchStrategy>=3)
3419          canSkip=0;
3420        }
3421        if (!numberBeforeTrust) {
3422          canSkip=1;
3423        }
3424        if (sort[iDo]<distanceToCutoff)
3425          canSkip=0;
3426        if (((numberTest2<=0&&numberTest<=0)||skipAll)&&sort[iDo]>distanceToCutoff) {
3427          canSkip=1; // always skip
3428          if (iDo>20) {
3429            if (!choiceObject) {
3430              delete choice.possibleBranch;
3431              choice.possibleBranch=NULL;
3432            }
3433            break; // give up anyway
3434          }
3435        }
3436#else
3437        if (((numberTest2<=0&&numberTest<=0)||skipAll)&&sort[iDo]>distanceToCutoff) {
3438          //canSkip=1; // always skip
3439          if (iDo>20) {
3440            if (!choiceObject) {
3441              delete choice.possibleBranch;
3442              choice.possibleBranch=NULL;
3443            }
3444            break; // give up anyway
3445          }
3446        }
3447#endif
3448        if (model->messageHandler()->logLevel()>3&&numberBeforeTrust&&dynamicObject) 
3449          dynamicObject->print(1,choice.possibleBranch->value());
3450        // was if (!canSkip)
3451        if (newWay)
3452        numberTest2--;
3453        if (!canSkip) {
3454          //#ifndef RANGING
3455          if (!doneHotStart) {
3456            // Mark hot start
3457            doneHotStart=true;
3458            assert (auxiliaryInfo->warmStart());
3459            solver->markHotStart();
3460            xMark++;
3461          }
3462          //#endif
3463          assert (!couldChooseFirst);
3464          numberTest--;
3465          if (!newWay)
3466          numberTest2--;
3467          // just do a few
3468#if NODE_NEW == 5  || NODE_NEW == 2
3469          //if (canSkip)
3470          if (searchStrategy==2)
3471            solver->setIntParam(OsiMaxNumIterationHotStart,10); 
3472#endif
3473          double objectiveChange ;
3474          double newObjectiveValue=1.0e100;
3475          int j;
3476          // status is 0 finished, 1 infeasible and other
3477          int iStatus;
3478          if (0) {
3479            CbcDynamicPseudoCostBranchingObject * cbcobj = dynamic_cast<CbcDynamicPseudoCostBranchingObject *> (choice.possibleBranch);
3480            if (cbcobj) {
3481              CbcSimpleIntegerDynamicPseudoCost * object = cbcobj->object();
3482              printf("strong %d ",iDo);
3483              object->print(1,0.5);
3484            }
3485          }
3486          /*
3487            Try the down direction first. (Specify the initial branching alternative as
3488            down with a call to way(-1). Each subsequent call to branch() performs the
3489            specified branch and advances the branch object state to the next branch
3490            alternative.)
3491          */
3492          choice.possibleBranch->way(-1) ;
3493#if NEW_UPDATE_OBJECT==0
3494          decision->saveBranchingObject( choice.possibleBranch);
3495#endif
3496          choice.possibleBranch->branch() ;
3497          solver->solveFromHotStart() ;
3498          bool needHotStartUpdate=false;
3499          numberStrongDone++;
3500          numberStrongIterations += solver->getIterationCount();
3501          /*
3502            We now have an estimate of objective degradation that we can use for strong
3503            branching. If we're over the cutoff, the variable is monotone up.
3504            If we actually made it to optimality, check for a solution, and if we have
3505            a good one, call setBestSolution to process it. Note that this may reduce the
3506            cutoff, so we check again to see if we can declare this variable monotone.
3507          */
3508          if (solver->isProvenOptimal())
3509            iStatus=0; // optimal
3510          else if (solver->isIterationLimitReached()
3511                   &&!solver->isDualObjectiveLimitReached())
3512            iStatus=2; // unknown
3513          else
3514            iStatus=1; // infeasible
3515          if (iStatus!=2&&solver->getIterationCount()>
3516              realMaxHotIterations)
3517            numberUnfinished++;
3518          newObjectiveValue = solver->getObjSense()*solver->getObjValue();
3519          choice.numItersDown = solver->getIterationCount();
3520          objectiveChange = CoinMax(newObjectiveValue  - objectiveValue_,0.0);
3521          // Update branching information if wanted
3522#if NEW_UPDATE_OBJECT==0
3523          decision->updateInformation( solver,this);
3524#elif NEW_UPDATE_OBJECT<2
3525          CbcBranchingObject * cbcobj = dynamic_cast<CbcBranchingObject *> (choice.possibleBranch);
3526          if (cbcobj) {
3527            CbcObject * object = cbcobj->object();
3528            assert (object);
3529            CbcObjectUpdateData update = object->createUpdateInformation(solver,this,cbcobj);
3530            object->updateInformation(update);
3531          } else {
3532            decision->updateInformation( solver,this);
3533          }
3534#else
3535          CbcBranchingObject * cbcobj = dynamic_cast<CbcBranchingObject *> (choice.possibleBranch);
3536          if (cbcobj) {
3537            CbcObject * object = cbcobj->object();
3538            assert (object) ;
3539            CbcObjectUpdateData update = object->createUpdateInformation(solver,this,cbcobj);
3540            update.objectNumber_ = choice.objectNumber;
3541            model->addUpdateInformation(update);
3542          } else {
3543            decision->updateInformation( solver,this);
3544          }
3545#endif
3546          if (!iStatus) {
3547            choice.finishedDown = true ;
3548            if (newObjectiveValue>=cutoff) {
3549              objectiveChange = 1.0e100; // say infeasible
3550              numberStrongInfeasible++;
3551            } else {
3552              // See if integer solution
3553              if (model->feasibleSolution(choice.numIntInfeasDown,
3554                                          choice.numObjInfeasDown)
3555                  &&model->problemFeasibility()->feasible(model,-1)>=0) {
3556                if (auxiliaryInfo->solutionAddsCuts()) {
3557                  needHotStartUpdate=true;
3558                  solver->unmarkHotStart();
3559                }
3560                model->setBestSolution(CBC_STRONGSOL,
3561                                       newObjectiveValue,
3562                                       solver->getColSolution()) ;
3563                if (needHotStartUpdate) {
3564                  solver->resolve();
3565                  newObjectiveValue = solver->getObjSense()*solver->getObjValue();
3566                  objectiveChange = CoinMax(newObjectiveValue  - objectiveValue_,0.0);
3567                  model->feasibleSolution(choice.numIntInfeasDown,
3568                                          choice.numObjInfeasDown);
3569                }
3570                model->setLastHeuristic(NULL);
3571                model->incrementUsed(solver->getColSolution());
3572                cutoff =model->getCutoff();
3573                if (newObjectiveValue >= cutoff) {      //  *new* cutoff
3574                  objectiveChange = 1.0e100 ;
3575                  numberStrongInfeasible++;
3576                }
3577              }
3578            }
3579          } else if (iStatus==1) {
3580            objectiveChange = 1.0e100 ;
3581            numberStrongInfeasible++;
3582          } else {
3583            // Can't say much as we did not finish
3584            choice.finishedDown = false ;
3585            numberUnfinished++;
3586          }
3587          choice.downMovement = objectiveChange ;
3588         
3589          // restore bounds
3590          for ( j=0;j<numberColumns;j++) {
3591            if (saveLower[j] != lower[j])
3592              solver->setColLower(j,saveLower[j]);
3593            if (saveUpper[j] != upper[j])
3594              solver->setColUpper(j,saveUpper[j]);
3595          }
3596          if(needHotStartUpdate) {
3597            needHotStartUpdate = false;
3598            solver->resolve();
3599            //we may again have an integer feasible solution
3600            int numberIntegerInfeasibilities;
3601            int numberObjectInfeasibilities;
3602            if (model->feasibleSolution(
3603                                        numberIntegerInfeasibilities,
3604                                        numberObjectInfeasibilities)) {
3605#ifdef BONMIN
3606              //In this case node has become integer feasible, let us exit the loop
3607              std::cout<<"Node has become integer feasible"<<std::endl;
3608              numberUnsatisfied_ = 0;
3609              break;
3610#endif
3611              double objValue = solver->getObjValue();
3612              model->setBestSolution(CBC_STRONGSOL,
3613                                     objValue,
3614                                     solver->getColSolution()) ;
3615              solver->resolve();
3616              cutoff =model->getCutoff();
3617            }
3618            solver->markHotStart();
3619          }
3620          //printf("Down on %d, status is %d, obj %g its %d cost %g finished %d inf %d infobj %d\n",
3621          //printf("Down on %d, status is %d, obj %g its %d cost %g finished %d inf %d infobj %d\n",
3622          //     choice.objectNumber,iStatus,newObjectiveValue,choice.numItersDown,
3623          //     choice.downMovement,choice.finishedDown,choice.numIntInfeasDown,
3624          //     choice.numObjInfeasDown);
3625         
3626          // repeat the whole exercise, forcing the variable up
3627#if NEW_UPDATE_OBJECT==0
3628          decision->saveBranchingObject( choice.possibleBranch);
3629#endif
3630          choice.possibleBranch->branch();
3631          solver->solveFromHotStart() ;
3632          numberStrongDone++;
3633          numberStrongIterations += solver->getIterationCount();
3634          /*
3635            We now have an estimate of objective degradation that we can use for strong
3636            branching. If we're over the cutoff, the variable is monotone up.
3637            If we actually made it to optimality, check for a solution, and if we have
3638            a good one, call setBestSolution to process it. Note that this may reduce the
3639            cutoff, so we check again to see if we can declare this variable monotone.
3640          */
3641          if (solver->isProvenOptimal())
3642            iStatus=0; // optimal
3643          else if (solver->isIterationLimitReached()
3644                   &&!solver->isDualObjectiveLimitReached())
3645            iStatus=2; // unknown
3646          else
3647            iStatus=1; // infeasible
3648          if (iStatus!=2&&solver->getIterationCount()>
3649              realMaxHotIterations)
3650            numberUnfinished++;
3651          newObjectiveValue = solver->getObjSense()*solver->getObjValue();
3652          choice.numItersUp = solver->getIterationCount();
3653          objectiveChange = CoinMax(newObjectiveValue  - objectiveValue_,0.0);
3654          // Update branching information if wanted
3655#if NEW_UPDATE_OBJECT==0
3656          decision->updateInformation( solver,this);
3657#elif NEW_UPDATE_OBJECT<2
3658          cbcobj = dynamic_cast<CbcBranchingObject *> (choice.possibleBranch);
3659          if (cbcobj) {
3660            CbcObject * object = cbcobj->object();
3661            CbcObjectUpdateData update = object->createUpdateInformation(solver,this,cbcobj);
3662            object->updateInformation(update);
3663          } else {
3664            decision->updateInformation( solver,this);
3665          }
3666#else
3667          cbcobj = dynamic_cast<CbcBranchingObject *> (choice.possibleBranch);
3668          if (cbcobj) {
3669            CbcObject * object = cbcobj->object();
3670            assert (object) ;
3671            CbcObjectUpdateData update = object->createUpdateInformation(solver,this,cbcobj);
3672            update.objectNumber_ = choice.objectNumber;
3673            model->addUpdateInformation(update);
3674          } else {
3675            decision->updateInformation( solver,this);
3676          }
3677#endif
3678          if (!iStatus) {
3679            choice.finishedUp = true ;
3680            if (newObjectiveValue>=cutoff) {
3681              objectiveChange = 1.0e100; // say infeasible
3682              numberStrongInfeasible++;
3683            } else {
3684              // See if integer solution
3685              if (model->feasibleSolution(choice.numIntInfeasUp,
3686                                          choice.numObjInfeasUp)
3687                  &&model->problemFeasibility()->feasible(model,-1)>=0) {
3688#ifdef BONMIN
3689                std::cout<<"Node has become integer feasible"<<std::endl;
3690                numberUnsatisfied_ = 0;
3691                break;
3692#endif
3693                if (auxiliaryInfo->solutionAddsCuts()) {
3694                  needHotStartUpdate=true;
3695                  solver->unmarkHotStart();
3696                }
3697                model->setBestSolution(CBC_STRONGSOL,
3698                                       newObjectiveValue,
3699                                       solver->getColSolution()) ;
3700                if (needHotStartUpdate) {
3701                  solver->resolve();
3702                  newObjectiveValue = solver->getObjSense()*solver->getObjValue();
3703                  objectiveChange = CoinMax(newObjectiveValue  - objectiveValue_,0.0);
3704                  model->feasibleSolution(choice.numIntInfeasDown,
3705                                          choice.numObjInfeasDown);
3706                }
3707                model->setLastHeuristic(NULL);
3708                model->incrementUsed(solver->getColSolution());
3709                cutoff =model->getCutoff();
3710                if (newObjectiveValue >= cutoff) {      //  *new* cutoff
3711                  objectiveChange = 1.0e100 ;
3712                  numberStrongInfeasible++;
3713                }
3714              }
3715            }
3716          } else if (iStatus==1) {
3717            objectiveChange = 1.0e100 ;
3718            numberStrongInfeasible++;
3719          } else {
3720            // Can't say much as we did not finish
3721            choice.finishedUp = false ;
3722            numberUnfinished++;
3723          }
3724          choice.upMovement = objectiveChange ;
3725         
3726          // restore bounds
3727          for ( j=0;j<numberColumns;j++) {
3728            if (saveLower[j] != lower[j])
3729              solver->setColLower(j,saveLower[j]);
3730            if (saveUpper[j] != upper[j])
3731              solver->setColUpper(j,saveUpper[j]);
3732          }
3733          if(needHotStartUpdate) {
3734            needHotStartUpdate = false;
3735            solver->resolve();
3736            //we may again have an integer feasible solution
3737            int numberIntegerInfeasibilities;
3738            int numberObjectInfeasibilities;
3739            if (model->feasibleSolution(
3740                                        numberIntegerInfeasibilities,
3741                                        numberObjectInfeasibilities)) {
3742              double objValue = solver->getObjValue();
3743              model->setBestSolution(CBC_STRONGSOL,
3744                                     objValue,
3745                                     solver->getColSolution()) ;
3746              solver->resolve();
3747              cutoff =model->getCutoff();
3748            }
3749            solver->markHotStart();
3750          }
3751         
3752          //printf("Up on %d, status is %d, obj %g its %d cost %g finished %d inf %d infobj %d\n",
3753          //     choice.objectNumber,iStatus,newObjectiveValue,choice.numItersUp,
3754          //     choice.upMovement,choice.finishedUp,choice.numIntInfeasUp,
3755          //     choice.numObjInfeasUp);
3756        }
3757   
3758        solver->setIntParam(OsiMaxNumIterationHotStart,saveLimit2); 
3759        /*
3760          End of evaluation for this candidate variable. Possibilities are:
3761          * Both sides below cutoff; this variable is a candidate for branching.
3762          * Both sides infeasible or above the objective cutoff: no further action
3763          here. Break from the evaluation loop and assume the node will be purged
3764          by the caller.
3765          * One side below cutoff: Install the branch (i.e., fix the variable). Break
3766          from the evaluation loop and assume the node will be reoptimised by the
3767          caller.
3768        */
3769        // reset
3770        choice.possibleBranch->resetNumberBranchesLeft();
3771        if (choice.upMovement<1.0e100) {
3772          if(choice.downMovement<1.0e100) {
3773            // In case solution coming in was odd
3774            choice.upMovement = CoinMax(0.0,choice.upMovement);
3775            choice.downMovement = CoinMax(0.0,choice.downMovement);
3776            if (couldChooseFirst)
3777              printf("candidate %d up %g down %g sort %g\n",iDo,choice.upMovement,choice.downMovement,sort[iDo]);
3778#if ZERO_ONE==2
3779            // branch on 0-1 first (temp)
3780            if (fabs(choice.possibleBranch->value())<1.0) {
3781              choice.upMovement *= ZERO_FAKE;
3782              choice.downMovement *= ZERO_FAKE;
3783            }
3784#endif
3785            // feasible - see which best
3786            if (!canSkip) {
3787              if (iColumn==-46) {
3788                printf("sort %g downest %g upest %g ",sort[iDo],downEstimate[iObject],
3789                     upEstimate[iObject]);
3790                printf("downMove %g upMove %g value %g current pseudo %g %g\n",
3791                       choice.downMovement,choice.upMovement,choice.possibleBranch->value(),
3792                       dynamicObject->downDynamicPseudoCost(),dynamicObject->upDynamicPseudoCost());
3793              }
3794              if (model->messageHandler()->logLevel()>3) 
3795                printf("sort %g downest %g upest %g ",sort[iDo],downEstimate[iObject],
3796                     upEstimate[iObject]);
3797              model->messageHandler()->message(CBC_STRONG,*model->messagesPointer())
3798                << iObject << iColumn
3799                <<choice.downMovement<<choice.numIntInfeasDown
3800                <<choice.upMovement<<choice.numIntInfeasUp
3801                <<choice.possibleBranch->value()
3802                <<CoinMessageEol;
3803            }
3804            //if (!stateOfSearch)
3805            //choice.numIntInfeasDown=99999; // temp fudge
3806            if (wantMiniTree)
3807              decision->setBestCriterion(-1.0);
3808            double bestCriterion = -1.0;
3809            //double gap = saveUpper[iColumn]-saveLower[iColumn];
3810            // Give precedence to ones with gap of 1.0
3811            //assert(gap>0.0);
3812            double factor = 1.0; //changeFactor/CoinMin(gap,100.0);
3813            int betterWay;
3814            {
3815              CbcBranchingObject * branchObj =
3816                dynamic_cast <CbcBranchingObject *>(branch_) ;
3817              if (branch_)
3818                assert (branchObj);
3819              betterWay = decision->betterBranch(choice.possibleBranch,
3820                                                     branchObj,
3821                                                     choice.upMovement*factor,
3822                                                     choice.numIntInfeasUp ,
3823                                                     choice.downMovement*factor,
3824                                                     choice.numIntInfeasDown );
3825            }
3826            if (wantMiniTree) {
3827              double criterion = decision->getBestCriterion();
3828              sort[numberMini]=-criterion;
3829              whichObject[numberMini++]=whichObject[iDo];
3830              assert (betterWay);
3831              if (criterion>bestCriterion) 
3832                bestCriterion=criterion;
3833              else
3834                betterWay=0;
3835            }
3836            if (iDo>=changeStrategy) {
3837              // make less likely
3838              changeStrategy+=numberStrong;
3839              changeFactor *= 0.9;
3840            }
3841            if (betterWay) {
3842              // C) create branching object
3843              if (choiceObject) {
3844                delete branch_;
3845                branch_ = choice.possibleBranch->clone();
3846              } else {
3847                delete branch_;
3848                branch_ = choice.possibleBranch;
3849                choice.possibleBranch=NULL;
3850              }
3851              {
3852                CbcBranchingObject * branchObj =
3853                  dynamic_cast <CbcBranchingObject *>(branch_) ;
3854                assert (branchObj);
3855                //branchObj->way(preferredWay);
3856                branchObj->way(betterWay);
3857              }
3858              if (couldChooseFirst)
3859                printf("choosing %d way %d\n",iDo,betterWay);
3860              bestChoice = choice.objectNumber;
3861              whichChoice = iDo;
3862              if (numberStrong<=1) {
3863                delete ws;
3864                ws=NULL;
3865                break;
3866              }
3867            } else {
3868              if (!choiceObject) {
3869                delete choice.possibleBranch;
3870                choice.possibleBranch=NULL;
3871              }
3872              if (iDo>=2*numberStrong) {
3873                delete ws;
3874                ws=NULL;
3875                break;
3876              }
3877              if (!dynamicObject||dynamicObject->numberTimesUp()>1) {
3878                if (iDo-whichChoice>=numberStrong) {
3879                  if (!choiceObject) {
3880                    delete choice.possibleBranch;
3881                    choice.possibleBranch=NULL;
3882                  }
3883                  break; // give up
3884                }
3885              } else {
3886                if (iDo-whichChoice>=2*numberStrong) {
3887                  delete ws;
3888                  ws=NULL;
3889                  if (!choiceObject) {
3890                    delete choice.possibleBranch;
3891                    choice.possibleBranch=NULL;
3892                  }
3893                  break; // give up
3894                }
3895              }
3896            }
3897          } else {
3898            // up feasible, down infeasible
3899            anyAction=-1;
3900            worstFeasible = CoinMax(worstFeasible,choice.upMovement);
3901            model->messageHandler()->message(CBC_STRONG,*model->messagesPointer())
3902              << iObject << iColumn
3903              <<choice.downMovement<<choice.numIntInfeasDown
3904              <<choice.upMovement<<choice.numIntInfeasUp
3905              <<choice.possibleBranch->value()
3906              <<CoinMessageEol;
3907            //printf("Down infeasible for choice %d sequence %d\n",i,
3908            // model->object(choice.objectNumber)->columnNumber());
3909            if (!solveAll) {
3910              choice.possibleBranch->way(1);
3911              choice.possibleBranch->branch();
3912              if (!choiceObject) {
3913                delete choice.possibleBranch;
3914                choice.possibleBranch=NULL;
3915              }
3916              delete ws;
3917              ws=NULL;
3918              break;
3919            } else {
3920              choice.fix=1;
3921              fixObject[numberToFix++]=choice;
3922#define FIXNOW
3923#ifdef FIXNOW
3924#if 0
3925              double value = ceil(saveSolution[iColumn]);
3926              saveLower[iColumn]=value;
3927              solver->setColLower(iColumn,value);
3928#else
3929              choice.possibleBranch->fix(solver,saveLower,saveUpper,1);
3930#endif
3931#endif
3932              if (!choiceObject) {
3933                choice.possibleBranch=NULL;
3934              } else {
3935                choiceObject = new CbcDynamicPseudoCostBranchingObject(*choiceObject);
3936                choice.possibleBranch=choiceObject;
3937              }
3938#ifdef FIXNOW
3939              assert(doneHotStart);
3940              solver->unmarkHotStart();
3941              solver->resolve();
3942              solver->markHotStart();
3943              // may be infeasible (if other way stopped on iterations)
3944              if (!solver->isProvenOptimal()) {
3945                // neither side feasible
3946                anyAction=-2;
3947                if (!choiceObject) {
3948                  delete choice.possibleBranch;
3949                  choice.possibleBranch=NULL;
3950                }
3951                //printf("Both infeasible for choice %d sequence %d\n",i,
3952                // model->object(choice.objectNumber)->columnNumber());
3953                delete ws;
3954                ws=NULL;
3955                break;
3956              }
3957#endif
3958            }
3959          }
3960        } else {
3961          if(choice.downMovement<1.0e100) {
3962            // down feasible, up infeasible
3963            anyAction=-1;
3964            worstFeasible = CoinMax(worstFeasible,choice.downMovement);
3965            model->messageHandler()->message(CBC_STRONG,*model->messagesPointer())
3966              << iObject << iColumn
3967              <<choice.downMovement<<choice.numIntInfeasDown
3968              <<choice.upMovement<<choice.numIntInfeasUp
3969              <<choice.possibleBranch->value()
3970              <<CoinMessageEol;
3971            //printf("Up infeasible for choice %d sequence %d\n",i,
3972            // model->object(choice.objectNumber)->columnNumber());
3973            if (!solveAll) {
3974              choice.possibleBranch->way(-1);
3975              choice.possibleBranch->branch();
3976              if (!choiceObject) {
3977                delete choice.possibleBranch;
3978                choice.possibleBranch=NULL;
3979              }
3980              delete ws;
3981              ws=NULL;
3982              break;
3983            } else {
3984              choice.fix=-1;
3985              fixObject[numberToFix++]=choice;
3986#ifdef FIXNOW
3987#if 0
3988              double value = floor(saveSolution[iColumn]);
3989              saveUpper[iColumn]=value;
3990              solver->setColUpper(iColumn,value);
3991#else
3992              choice.possibleBranch->fix(solver,saveLower,saveUpper,-1);
3993#endif
3994#endif
3995              if (!choiceObject) {
3996                choice.possibleBranch=NULL;
3997              } else {
3998                choiceObject = new CbcDynamicPseudoCostBranchingObject(*choiceObject);
3999                choice.possibleBranch=choiceObject;
4000              }
4001#ifdef FIXNOW
4002              assert(doneHotStart);
4003              solver->unmarkHotStart();
4004              solver->resolve();
4005              solver->markHotStart();
4006              // may be infeasible (if other way stopped on iterations)
4007              if (!solver->isProvenOptimal()) {
4008                // neither side feasible
4009                anyAction=-2;
4010                if (!choiceObject) {
4011                  delete choice.possibleBranch;
4012                  choice.possibleBranch=NULL;
4013                }
4014                //printf("Both infeasible for choice %d sequence %d\n",i,
4015                // model->object(choice.objectNumber)->columnNumber());
4016                delete ws;
4017                ws=NULL;
4018                break;
4019              }
4020#endif
4021            }
4022          } else {
4023            // neither side feasible
4024            anyAction=-2;
4025            if (!choiceObject) {
4026              delete choice.possibleBranch;
4027              choice.possibleBranch=NULL;
4028            }
4029            //printf("Both infeasible for choice %d sequence %d\n",i,
4030            // model->object(choice.objectNumber)->columnNumber());
4031            delete ws;
4032            ws=NULL;
4033            break;
4034          }
4035        }
4036        // Check max time
4037        hitMaxTime = ( CoinCpuTime()-model->getDblParam(CbcModel::CbcStartSeconds) > 
4038                       model->getDblParam(CbcModel::CbcMaximumSeconds));
4039        if (hitMaxTime) {
4040          // make sure rest are fast
4041          doQuickly=true;
4042          for ( int jDo=iDo+1;jDo<numberToDo;jDo++) {
4043            int iObject = whichObject[iDo];
4044            OsiObject * object = model->modifiableObject(iObject);
4045            CbcSimpleIntegerDynamicPseudoCost * dynamicObject =
4046              dynamic_cast <CbcSimpleIntegerDynamicPseudoCost *>(object) ;
4047            if (dynamicObject) 
4048              dynamicObject->setNumberBeforeTrust(0);
4049          }
4050          numberTest=0;
4051          distanceToCutoff=-COIN_DBL_MAX;
4052        }
4053        if (!choiceObject) {
4054          delete choice.possibleBranch;
4055        }
4056      }
4057      double averageChange = model->sumChangeObjective()/((double) model->getNodeCount());
4058      if (depth_<10||worstFeasible>0.2*averageChange) 
4059        solveAll=false;
4060      if (model->messageHandler()->logLevel()>3||false) { 
4061        if (anyAction==-2) {
4062          printf("infeasible\n");
4063        } else if(anyAction==-1) {
4064          if (!solveAll)
4065            printf("%d fixed\n",numberToFix);
4066          else
4067            printf("%d fixed AND choosing %d iDo %d iChosenWhen %d numberToDo %d\n",numberToFix,bestChoice,
4068                   iDo,whichChoice,numberToDo);
4069        } else {
4070          int iObject = whichObject[whichChoice];
4071          OsiObject * object = model->modifiableObject(iObject);
4072          CbcSimpleIntegerDynamicPseudoCost * dynamicObject =
4073            dynamic_cast <CbcSimpleIntegerDynamicPseudoCost *>(object) ;
4074          if (dynamicObject) {
4075            int iColumn = dynamicObject->columnNumber();
4076            printf("choosing %d (column %d) iChosenWhen %d numberToDo %d\n",bestChoice,
4077                   iColumn,whichChoice,numberToDo);
4078          }
4079        }
4080      }
4081      if (doneHotStart) {
4082        // Delete the snapshot
4083        solver->unmarkHotStart();
4084        // back to normal
4085        solver->setHintParam(OsiDoInBranchAndCut,true,OsiHintDo,NULL) ;
4086        // restore basis
4087        solver->setWarmStart(ws);
4088      }
4089      solver->setIntParam(OsiMaxNumIterationHotStart,saveLimit);
4090      // Unless infeasible we will carry on
4091      // But we could fix anyway
4092      if (numberToFix&&!hitMaxTime) {
4093        if (anyAction==-2) {
4094          // take off
4095          for (i = 0 ; i < numberToFix ; i++) {
4096            delete fixObject[i].possibleBranch;
4097          }
4098        } else {
4099          // apply and take off
4100          for (i = 0 ; i < numberToFix ; i++) {
4101#ifndef FIXNOW
4102            fixObject[i].possibleBranch->way(fixObject[i].fix) ;
4103            fixObject[i].possibleBranch->branch() ;
4104#endif
4105            delete fixObject[i].possibleBranch;
4106          }
4107          bool feasible=true;
4108#if ACTION <2
4109          if (solveAll) {
4110            // can do quick optimality check
4111            int easy=2;
4112            solver->setHintParam(OsiDoInBranchAndCut,true,OsiHintDo,&easy) ;
4113            solver->resolve() ;
4114            solver->setHintParam(OsiDoInBranchAndCut,true,OsiHintDo,NULL) ;
4115            feasible = solver->isProvenOptimal();
4116            if (feasible) {
4117              anyAction=0;
4118              numberMini=0;
4119              memcpy(saveSolution,solver->getColSolution(),numberColumns*sizeof(double));
4120              model->reserveCurrentSolution(saveSolution);
4121              memcpy(saveLower,solver->getColLower(),numberColumns*sizeof(double));
4122              memcpy(saveUpper,solver->getColUpper(),numberColumns*sizeof(double));
4123              model->setPointers(solver);
4124              // See if candidate still possible
4125              if (branch_) {
4126                const OsiObject * object = model->object(bestChoice);
4127                int preferredWay;
4128                double infeasibility = object->infeasibility(&usefulInfo,preferredWay);
4129                if (!infeasibility) {
4130                  // take out
4131                  delete branch_;
4132                  branch_=NULL;
4133                } else {
4134                  CbcBranchingObject * branchObj =
4135                    dynamic_cast <CbcBranchingObject *>(branch_) ;
4136                  assert (branchObj);
4137                  branchObj->way(preferredWay);
4138                }
4139              }
4140            } else {
4141              anyAction=-2;
4142              finished=true;
4143            }
4144          }
4145#endif
4146          // If  fixed then round again
4147          if (!branch_&&anyAction!=-2) {
4148            finished=false;
4149          }
4150          // If these in then different action
4151#if ACTION == 1
4152          if (!anyAction)
4153            anyAction=-1;
4154          finished=true;
4155#endif
4156        }
4157      }
4158      delete ws;
4159    }
4160  }
4161  if (model->messageHandler()->logLevel()>2) 
4162    printf("%d strong, %d iters, %d pen, %d mark, %d fixed, action %d nnott %d nt %d, %d dq %s ns %d\n",
4163         numberStrongDone,numberStrongIterations,xPen,xMark,
4164           numberToFix,anyAction,numberNotTrusted,px[0],px[1],px[2]>0 ? "y" : "n",px[3]);
4165  // update number of strong iterations etc
4166  model->incrementStrongInfo(numberStrongDone,numberStrongIterations,
4167                             anyAction==-2 ? 0:numberToFix,anyAction==-2);
4168  if (!newWay) {
4169  if (((model->searchStrategy()+1)%1000)==0) {
4170#ifndef COIN_DEVELOP
4171    if (solver->messageHandler()->logLevel()>1)
4172#endif
4173      printf("%d strong, %d iters, %d inf, %d not finished, %d not trusted\n",
4174             numberStrongDone,numberStrongIterations,numberStrongInfeasible,numberUnfinished,
4175             numberNotTrusted);
4176    // decide what to do
4177    int strategy=1;
4178    if (((numberUnfinished*4>numberStrongDone&&
4179         numberStrongInfeasible*40<numberStrongDone)||
4180         numberStrongInfeasible<0)&&model->numberStrong()<10&&model->numberBeforeTrust()<20) {
4181      strategy=2;
4182#ifdef COIN_DEVELOP
4183      //if (model->logLevel()>1)
4184        printf("going to strategy 2\n");
4185#endif
4186#if NODE_NEW==1
4187      // Basically off
4188      model->setNumberStrong(0);
4189      model->setNumberBeforeTrust(-2);
4190#elif NODE_NEW==2
4191      // Weaken
4192      model->setNumberStrong(2);
4193      model->setNumberBeforeTrust(1);
4194#elif NODE_NEW==3
4195      // Off and ranging
4196      model->setNumberStrong(0);
4197      model->setNumberBeforeTrust(-1);
4198#elif NODE_NEW==4
4199      // Off and pseudo shadow prices
4200      model->setNumberStrong(0);
4201      model->setNumberBeforeTrust(-2);
4202#elif NODE_NEW==5
4203      // Just fewer iterations
4204#endif
4205      model->synchronizeNumberBeforeTrust();
4206    }
4207    if (numberNodes)
4208      strategy=1;  // should only happen after hot start
4209    if (model->searchStrategy()<999)
4210      model->setSearchStrategy(strategy);
4211  }
4212  if (model->searchStrategy()==1&&numberNodes>500&&numberNodes<-510) {
4213#ifndef COIN_DEVELOP
4214    if (solver->messageHandler()->logLevel()>1)
4215#endif
4216      printf("after %d nodes - %d strong, %d iters, %d inf, %d not finished, %d not trusted\n",
4217             numberNodes,numberStrongDone,numberStrongIterations,numberStrongInfeasible,numberUnfinished,
4218             numberNotTrusted);
4219    // decide what to do
4220    if (numberUnfinished*10>numberStrongDone+1||
4221        !numberStrongInfeasible) {
4222      //if (model->logLevel()>1)
4223        printf("going to strategy 2\n");
4224#if NODE_NEW==1
4225      // Basically off
4226      model->setNumberStrong(0);
4227      model->setNumberBeforeTrust(-2);
4228#elif NODE_NEW==2
4229      // Weaken
4230      model->setNumberStrong(2);
4231      model->setNumberBeforeTrust(1);
4232#elif NODE_NEW==3
4233      // Off and ranging
4234      model->setNumberStrong(0);
4235      model->setNumberBeforeTrust(-1);
4236#elif NODE_NEW==4
4237      // Off and pseudo shadow prices
4238      model->setNumberStrong(0);
4239      model->setNumberBeforeTrust(-2);
4240#elif NODE_NEW==5
4241      // Just fewer iterations
4242#endif
4243      model->synchronizeNumberBeforeTrust();
4244      model->setSearchStrategy(2);
4245    }
4246  }
4247  }
4248  //if (numberToFix&&depth_<5)
4249  //printf("%d fixed by strong at depth %d\n",numberToFix,depth_);
4250  // Set guessed solution value
4251  guessedObjectiveValue_ = objectiveValue_+estimatedDegradation;
4252 
4253  // Get collection of branches if mini tree wanted
4254  if (anyAction==0&&numberMini&&numberMini>1) {
4255    // Sort
4256    CoinSort_2(sort,sort+numberMini,whichObject);
4257    delete branch_;
4258    branch_=NULL;
4259    numberMini = CoinMin(numberMini,model->sizeMiniTree());
4260    anyAction=numberMini;
4261    branches = new OsiSolverBranch[numberMini];
4262    for (int iDo=0;iDo<numberMini;iDo++) {
4263      int iObject = whichObject[iDo];
4264      OsiObject * object = model->modifiableObject(iObject);
4265      CbcSimpleInteger * obj =
4266        dynamic_cast <CbcSimpleInteger *>(object) ;
4267      OsiSolverBranch * oneBranch;
4268      if (obj) {
4269        oneBranch = obj->solverBranch(solver,&usefulInfo);
4270      } else {
4271        CbcObject * obj =
4272          dynamic_cast <CbcObject *>(object) ;
4273        assert (obj);
4274        oneBranch = obj->solverBranch();
4275      }
4276      branches[iDo]=*oneBranch;
4277      delete oneBranch;
4278    }
4279  }
4280/*
4281  Cleanup, then we're finished
4282*/
4283  if (!model->branchingMethod())
4284    delete decision;
4285
4286  delete choiceObject;
4287  delete [] fixObject;
4288  delete [] sort;
4289  delete [] whichObject;
4290  delete [] objectMark;
4291  delete [] saveLower;
4292  delete [] saveUpper;
4293  delete [] upEstimate;
4294  delete [] downEstimate;
4295# ifdef COIN_HAS_CLP
4296  if (osiclp) 
4297    osiclp->setSpecialOptions(saveClpOptions);
4298# endif
4299  // restore solution
4300  solver->setColSolution(saveSolution);
4301  model->reserveCurrentSolution(saveSolution);
4302  delete [] saveSolution;
4303  model->setStateOfSearch(saveStateOfSearch);
4304  model->setLogLevel(saveLogLevel);
4305  // delete extra regions
4306  if (usefulInfo.usefulRegion_) {
4307    delete [] usefulInfo.usefulRegion_;
4308    delete [] usefulInfo.indexRegion_;
4309    delete [] usefulInfo.pi_;
4310    usefulInfo.usefulRegion_ = NULL;
4311    usefulInfo.indexRegion_ = NULL;
4312    usefulInfo.pi_ = NULL;
4313  }
4314  return anyAction;
4315}
4316int CbcNode::analyze (CbcModel *model, double * results)
4317{
4318  int i;
4319  int numberIterationsAllowed = model->numberAnalyzeIterations();
4320  OsiSolverInterface * solver = model->solver();
4321  objectiveValue_ = solver->getObjSense()*solver->getObjValue();
4322  double cutoff =model->getCutoff();
4323  const double * lower = solver->getColLower();
4324  const double * upper = solver->getColUpper();
4325  const double * dj = solver->getReducedCost();
4326  int numberObjects = model->numberObjects();
4327  int numberColumns = model->getNumCols();
4328  // Initialize arrays
4329  int numberIntegers = model->numberIntegers();
4330  int * back = new int[numberColumns];
4331  const int * integerVariable = model->integerVariable();
4332  for (i=0;i<numberColumns;i++) 
4333    back[i]=-1;
4334  // What results is
4335  double * newLower = results;
4336  double * objLower = newLower+numberIntegers;
4337  double * newUpper = objLower+numberIntegers;
4338  double * objUpper = newUpper+numberIntegers;
4339  for (i=0;i<numberIntegers;i++) {
4340    int iColumn = integerVariable[i];
4341    back[iColumn]=i;
4342    newLower[i]=0.0;
4343    objLower[i]=-COIN_DBL_MAX;
4344    newUpper[i]=0.0;
4345    objUpper[i]=-COIN_DBL_MAX;
4346  }
4347  double * saveUpper = new double[numberColumns];
4348  double * saveLower = new double[numberColumns];
4349  int anyAction=0;
4350  // Save solution in case heuristics need good solution later
4351 
4352  double * saveSolution = new double[numberColumns];
4353  memcpy(saveSolution,solver->getColSolution(),numberColumns*sizeof(double));
4354  model->reserveCurrentSolution(saveSolution);
4355  for (i=0;i<numberColumns;i++) {
4356    saveLower[i] = lower[i];
4357    saveUpper[i] = upper[i];
4358  }
4359  // Get arrays to sort
4360  double * sort = new double[numberObjects];
4361  int * whichObject = new int[numberObjects];
4362  int numberToFix=0;
4363  int numberToDo=0;
4364  double integerTolerance = 
4365    model->getDblParam(CbcModel::CbcIntegerTolerance);
4366  // point to useful information
4367  OsiBranchingInformation usefulInfo = model->usefulInformation();
4368  // and modify
4369  usefulInfo.depth_=depth_;
4370     
4371  // compute current state
4372  int numberObjectInfeasibilities; // just odd ones
4373  int numberIntegerInfeasibilities;
4374  model->feasibleSolution(
4375                          numberIntegerInfeasibilities,
4376                          numberObjectInfeasibilities);
4377# ifdef COIN_HAS_CLP
4378  OsiClpSolverInterface * osiclp = dynamic_cast< OsiClpSolverInterface*> (solver);
4379  int saveClpOptions=0;
4380  bool fastIterations = (model->specialOptions()&8)!=0;
4381  if (osiclp&&fastIterations) {
4382    // for faster hot start
4383    saveClpOptions = osiclp->specialOptions();
4384    osiclp->setSpecialOptions(saveClpOptions|8192);
4385  }
4386# else
4387  bool fastIterations = false ;
4388# endif
4389  /*
4390    Scan for branching objects that indicate infeasibility. Choose candidates
4391    using priority as the first criteria, then integer infeasibility.
4392   
4393    The algorithm is to fill the array with a set of good candidates (by
4394    infeasibility) with priority bestPriority.  Finding a candidate with
4395    priority better (less) than bestPriority flushes the choice array. (This
4396    serves as initialization when the first candidate is found.)
4397   
4398  */
4399  numberToDo=0;
4400  for (i=0;i<numberObjects;i++) {
4401    OsiObject * object = model->modifiableObject(i);
4402    CbcSimpleIntegerDynamicPseudoCost * dynamicObject =
4403      dynamic_cast <CbcSimpleIntegerDynamicPseudoCost *>(object) ;
4404    if(!dynamicObject)
4405      continue;
4406    int preferredWay;
4407    double infeasibility = object->infeasibility(&usefulInfo,preferredWay);
4408    int iColumn = dynamicObject->columnNumber();
4409    if (saveUpper[iColumn]==saveLower[iColumn])
4410      continue;
4411    if (infeasibility)
4412      sort[numberToDo]=-1.0e10-infeasibility;
4413    else
4414      sort[numberToDo]=-fabs(dj[iColumn]);
4415    whichObject[numberToDo++]=i;
4416  }
4417  // Save basis
4418  CoinWarmStart * ws = solver->getWarmStart();
4419  int saveLimit;
4420  solver->getIntParam(OsiMaxNumIterationHotStart,saveLimit);
4421  int targetIterations = CoinMax(500,numberIterationsAllowed/numberObjects);
4422  if (saveLimit<targetIterations)
4423    solver->setIntParam(OsiMaxNumIterationHotStart,targetIterations); 
4424  // Mark hot start
4425  solver->markHotStart();
4426  // Sort
4427  CoinSort_2(sort,sort+numberToDo,whichObject);
4428  //double distanceToCutoff=model->getCutoff()-objectiveValue_;
4429  double * currentSolution = model->currentSolution();
4430  double objMin = 1.0e50;
4431  double objMax = -1.0e50;
4432  bool needResolve=false;
4433  int iDo;
4434  for (iDo=0;iDo<numberToDo;iDo++) {
4435    CbcStrongInfo choice;
4436    int iObject = whichObject[iDo];
4437    OsiObject * object = model->modifiableObject(iObject);
4438    CbcSimpleIntegerDynamicPseudoCost * dynamicObject =
4439      dynamic_cast <CbcSimpleIntegerDynamicPseudoCost *>(object) ;
4440    int iColumn = dynamicObject->columnNumber();
4441    int preferredWay;
4442    object->infeasibility(&usefulInfo,preferredWay);
4443    double value = currentSolution[iColumn];
4444    double nearest = floor(value+0.5);
4445    double lowerValue = floor(value);
4446    bool satisfied=false;
4447    if (fabs(value-nearest)<=integerTolerance||value<saveLower[iColumn]||value>saveUpper[iColumn]) {
4448      satisfied=true;
4449      double newValue;
4450      if (nearest<saveUpper[iColumn]) {
4451        newValue = nearest + 1.0001*integerTolerance;
4452        lowerValue = nearest;
4453      } else {
4454        newValue = nearest - 1.0001*integerTolerance;
4455        lowerValue = nearest-1;
4456      }
4457      currentSolution[iColumn]=newValue;
4458    }
4459    double upperValue = lowerValue+1.0;
4460    CbcSimpleInteger * obj =
4461      dynamic_cast <CbcSimpleInteger *>(object) ;
4462    if (obj) {
4463      choice.possibleBranch=obj->createBranch(solver,&usefulInfo,preferredWay);
4464    } else {
4465      CbcObject * obj =
4466        dynamic_cast <CbcObject *>(object) ;
4467      assert (obj);
4468      choice.possibleBranch=obj->createBranch(preferredWay);
4469    }
4470    currentSolution[iColumn]=value;
4471    // Save which object it was
4472    choice.objectNumber=iObject;
4473    choice.numIntInfeasUp = numberUnsatisfied_;
4474    choice.numIntInfeasDown = numberUnsatisfied_;
4475    choice.downMovement = 0.0;
4476    choice.upMovement = 0.0;
4477    choice.numItersDown = 0;
4478    choice.numItersUp = 0;
4479    choice.fix=0; // say not fixed
4480    double objectiveChange ;
4481    double newObjectiveValue=1.0e100;
4482    int j;
4483    // status is 0 finished, 1 infeasible and other
4484    int iStatus;
4485    /*
4486      Try the down direction first. (Specify the initial branching alternative as
4487      down with a call to way(-1). Each subsequent call to branch() performs the
4488      specified branch and advances the branch object state to the next branch
4489      alternative.)
4490    */
4491    choice.possibleBranch->way(-1) ;
4492    choice.possibleBranch->branch() ;
4493    if (fabs(value-lowerValue)>integerTolerance) {
4494      solver->solveFromHotStart() ;
4495      /*
4496        We now have an estimate of objective degradation that we can use for strong
4497        branching. If we're over the cutoff, the variable is monotone up.
4498        If we actually made it to optimality, check for a solution, and if we have
4499        a good one, call setBestSolution to process it. Note that this may reduce the
4500        cutoff, so we check again to see if we can declare this variable monotone.
4501      */
4502      if (solver->isProvenOptimal())
4503        iStatus=0; // optimal
4504      else if (solver->isIterationLimitReached()
4505               &&!solver->isDualObjectiveLimitReached())
4506        iStatus=2; // unknown
4507      else
4508        iStatus=1; // infeasible
4509      newObjectiveValue = solver->getObjSense()*solver->getObjValue();
4510      choice.numItersDown = solver->getIterationCount();
4511      numberIterationsAllowed -= choice.numItersDown;
4512      objectiveChange = newObjectiveValue  - objectiveValue_;
4513      if (!iStatus) {
4514        choice.finishedDown = true ;
4515        if (newObjectiveValue>=cutoff) {
4516          objectiveChange = 1.0e100; // say infeasible
4517        } else {
4518          // See if integer solution
4519          if (model->feasibleSolution(choice.numIntInfeasDown,
4520                                      choice.numObjInfeasDown)
4521              &&model->problemFeasibility()->feasible(model,-1)>=0) {
4522            model->setBestSolution(CBC_STRONGSOL,
4523                                   newObjectiveValue,
4524                                   solver->getColSolution()) ;
4525            model->setLastHeuristic(NULL);
4526            model->incrementUsed(solver->getColSolution());
4527            cutoff =model->getCutoff();
4528            if (newObjectiveValue >= cutoff)    //  *new* cutoff
4529              objectiveChange = 1.0e100 ;
4530          }
4531        }
4532      } else if (iStatus==1) {
4533        objectiveChange = 1.0e100 ;
4534      } else {
4535        // Can't say much as we did not finish
4536        choice.finishedDown = false ;
4537      }
4538      choice.downMovement = objectiveChange ;
4539    }
4540    // restore bounds
4541    for ( j=0;j<numberColumns;j++) {
4542      if (saveLower[j] != lower[j])
4543        solver->setColLower(j,saveLower[j]);
4544      if (saveUpper[j] != upper[j])
4545        solver->setColUpper(j,saveUpper[j]);
4546    }
4547    // repeat the whole exercise, forcing the variable up
4548    choice.possibleBranch->branch();
4549    if (fabs(value-upperValue)>integerTolerance) {
4550      solver->solveFromHotStart() ;
4551      /*
4552        We now have an estimate of objective degradation that we can use for strong
4553        branching. If we're over the cutoff, the variable is monotone up.
4554        If we actually made it to optimality, check for a solution, and if we have
4555        a good one, call setBestSolution to process it. Note that this may reduce the
4556        cutoff, so we check again to see if we can declare this variable monotone.
4557      */
4558      if (solver->isProvenOptimal())
4559        iStatus=0; // optimal
4560      else if (solver->isIterationLimitReached()
4561               &&!solver->isDualObjectiveLimitReached())
4562        iStatus=2; // unknown
4563      else
4564        iStatus=1; // infeasible
4565      newObjectiveValue = solver->getObjSense()*solver->getObjValue();
4566      choice.numItersUp = solver->getIterationCount();
4567      numberIterationsAllowed -= choice.numItersUp;
4568      objectiveChange = newObjectiveValue  - objectiveValue_;
4569      if (!iStatus) {
4570        choice.finishedUp = true ;
4571        if (newObjectiveValue>=cutoff) {
4572          objectiveChange = 1.0e100; // say infeasible
4573        } else {
4574          // See if integer solution
4575          if (model->feasibleSolution(choice.numIntInfeasUp,
4576                                      choice.numObjInfeasUp)
4577              &&model->problemFeasibility()->feasible(model,-1)>=0) {
4578            model->setBestSolution(CBC_STRONGSOL,
4579                                   newObjectiveValue,
4580                                   solver->getColSolution()) ;
4581            model->setLastHeuristic(NULL);
4582            model->incrementUsed(solver->getColSolution());
4583            cutoff =model->getCutoff();
4584            if (newObjectiveValue >= cutoff)    //  *new* cutoff
4585              objectiveChange = 1.0e100 ;
4586          }
4587        }
4588      } else if (iStatus==1) {
4589        objectiveChange = 1.0e100 ;
4590      } else {
4591        // Can't say much as we did not finish
4592        choice.finishedUp = false ;
4593      }
4594      choice.upMovement = objectiveChange ;
4595     
4596      // restore bounds
4597      for ( j=0;j<numberColumns;j++) {
4598        if (saveLower[j] != lower[j])
4599          solver->setColLower(j,saveLower[j]);
4600        if (saveUpper[j] != upper[j])
4601          solver->setColUpper(j,saveUpper[j]);
4602      }
4603    }
4604    // If objective goes above certain amount we can set bound
4605    int jInt = back[iColumn];
4606    newLower[jInt]=upperValue;
4607    if (choice.finishedDown)
4608      objLower[jInt]=choice.downMovement+objectiveValue_;
4609    else
4610      objLower[jInt]=objectiveValue_;
4611    newUpper[jInt]=lowerValue;
4612    if (choice.finishedUp)
4613      objUpper[jInt]=choice.upMovement+objectiveValue_;
4614    else
4615      objUpper[jInt]=objectiveValue_;
4616    objMin = CoinMin(CoinMin(objLower[jInt],objUpper[jInt]),objMin);
4617    /*
4618      End of evaluation for this candidate variable. Possibilities are:
4619      * Both sides below cutoff; this variable is a candidate for branching.
4620      * Both sides infeasible or above the objective cutoff: no further action
4621      here. Break from the evaluation loop and assume the node will be purged
4622      by the caller.
4623      * One side below cutoff: Install the branch (i.e., fix the variable). Break
4624      from the evaluation loop and assume the node will be reoptimised by the
4625      caller.
4626    */
4627    if (choice.upMovement<1.0e100) {
4628      if(choice.downMovement<1.0e100) {
4629        objMax = CoinMax(CoinMax(objLower[jInt],objUpper[jInt]),objMax);
4630        // In case solution coming in was odd
4631        choice.upMovement = CoinMax(0.0,choice.upMovement);
4632        choice.downMovement = CoinMax(0.0,choice.downMovement);
4633        // feasible -
4634        model->messageHandler()->message(CBC_STRONG,*model->messagesPointer())
4635          << iObject << iColumn
4636          <<choice.downMovement<<choice.numIntInfeasDown
4637          <<choice.upMovement<<choice.numIntInfeasUp
4638          <<value
4639          <<CoinMessageEol;
4640      } else {
4641        // up feasible, down infeasible
4642        anyAction=-1;
4643        if (!satisfied)
4644          needResolve=true;
4645        choice.fix=1;
4646        numberToFix++;
4647        saveLower[iColumn]=upperValue;
4648        solver->setColLower(iColumn,upperValue);
4649      }
4650    } else {
4651      if(choice.downMovement<1.0e100) {
4652        // down feasible, up infeasible
4653        anyAction=-1;
4654        if (!satisfied)
4655          needResolve=true;
4656        choice.fix=-1;
4657        numberToFix++;
4658        saveUpper[iColumn]=lowerValue;
4659        solver->setColUpper(iColumn,lowerValue);
4660      } else {
4661        // neither side feasible
4662        anyAction=-2;
4663        printf("Both infeasible for choice %d sequence %d\n",i,
4664               model->object(choice.objectNumber)->columnNumber());
4665        delete ws;
4666        ws=NULL;
4667        //solver->writeMps("bad");
4668        numberToFix=-1;
4669        delete choice.possibleBranch;
4670        choice.possibleBranch=NULL;
4671        break;
4672      }
4673    }
4674    delete choice.possibleBranch;
4675    if (numberIterationsAllowed<=0)
4676      break;
4677    //printf("obj %d, col %d, down %g up %g value %g\n",iObject,iColumn,
4678    //     choice.downMovement,choice.upMovement,value);
4679  }
4680  printf("Best possible solution %g, can fix more if solution of %g found - looked at %d variables in %d iterations\n",
4681         objMin,objMax,iDo,model->numberAnalyzeIterations()-numberIterationsAllowed);
4682  model->setNumberAnalyzeIterations(numberIterationsAllowed);
4683  // Delete the snapshot
4684  solver->unmarkHotStart();
4685  // back to normal
4686  solver->setHintParam(OsiDoInBranchAndCut,true,OsiHintDo,NULL) ;
4687  solver->setIntParam(OsiMaxNumIterationHotStart,saveLimit);
4688  // restore basis
4689  solver->setWarmStart(ws);
4690  delete ws;
4691   
4692  delete [] sort;
4693  delete [] whichObject;
4694  delete [] saveLower;
4695  delete [] saveUpper;
4696  delete [] back;
4697  // restore solution
4698  solver->setColSolution(saveSolution);
4699# ifdef COIN_HAS_CLP
4700  if (osiclp) 
4701    osiclp->setSpecialOptions(saveClpOptions);
4702# endif
4703  model->reserveCurrentSolution(saveSolution);
4704  delete [] saveSolution;
4705  if (needResolve)
4706    solver->resolve();
4707  return numberToFix;
4708}
4709
4710
4711CbcNode::CbcNode(const CbcNode & rhs) 
4712{ 
4713#ifdef CHECK_NODE
4714  printf("CbcNode %x Constructor from rhs %x\n",this,&rhs);
4715#endif
4716  if (rhs.nodeInfo_)
4717    nodeInfo_ = rhs.nodeInfo_->clone();
4718  else
4719    nodeInfo_=NULL;
4720  objectiveValue_=rhs.objectiveValue_;
4721  guessedObjectiveValue_ = rhs.guessedObjectiveValue_;
4722  sumInfeasibilities_ = rhs.sumInfeasibilities_;
4723  if (rhs.branch_)
4724    branch_=rhs.branch_->clone();
4725  else
4726    branch_=NULL;
4727  depth_ = rhs.depth_;
4728  numberUnsatisfied_ = rhs.numberUnsatisfied_;
4729  nodeNumber_ = rhs.nodeNumber_;
4730  state_ = rhs.state_;
4731  if (nodeInfo_)
4732    assert ((state_&2)!=0);
4733  else
4734    assert ((state_&2)==0);
4735}
4736
4737CbcNode &
4738CbcNode::operator=(const CbcNode & rhs)
4739{
4740  if (this != &rhs) {
4741    delete nodeInfo_;
4742    if (rhs.nodeInfo_)
4743      nodeInfo_ = rhs.nodeInfo_->clone();
4744    else
4745      nodeInfo_ = NULL;
4746    objectiveValue_=rhs.objectiveValue_;
4747    guessedObjectiveValue_ = rhs.guessedObjectiveValue_;
4748    sumInfeasibilities_ = rhs.sumInfeasibilities_;
4749    if (rhs.branch_)
4750      branch_=rhs.branch_->clone();
4751    else
4752      branch_=NULL,
4753    depth_ = rhs.depth_;
4754    numberUnsatisfied_ = rhs.numberUnsatisfied_;
4755    nodeNumber_ = rhs.nodeNumber_;
4756    state_ = rhs.state_;
4757    if (nodeInfo_)
4758      assert ((state_&2)!=0);
4759    else
4760      assert ((state_&2)==0);
4761  }
4762  return *this;
4763}
4764CbcNode::~CbcNode ()
4765{
4766#ifdef CHECK_NODE
4767  if (nodeInfo_) {
4768    printf("CbcNode %x Destructor nodeInfo %x (%d)\n",
4769         this,nodeInfo_,nodeInfo_->numberPointingToThis());
4770    //assert(nodeInfo_->numberPointingToThis()>=0);
4771  } else {
4772    printf("CbcNode %x Destructor nodeInfo %x (?)\n",
4773         this,nodeInfo_);
4774  }
4775#endif
4776  if (nodeInfo_) {
4777    // was if (nodeInfo_&&(state_&2)!=0) {
4778    nodeInfo_->nullOwner();
4779    int numberToDelete=nodeInfo_->numberBranchesLeft();
4780    //    CbcNodeInfo * parent = nodeInfo_->parent();
4781    //assert (nodeInfo_->numberPointingToThis()>0);
4782    if (nodeInfo_->decrement(numberToDelete)==0||(state_&2)==0) {
4783      if ((state_&2)==0) 
4784        nodeInfo_->nullParent();
4785      delete nodeInfo_;
4786    } else {
4787      //printf("node %x nodeinfo %x parent %x\n",this,nodeInfo_,nodeInfo_->parent());
4788      // anyway decrement parent
4789      //if (parent)
4790      ///parent->decrement(1);
4791    }
4792  }
4793  delete branch_;
4794}
4795// Decrement  active cut counts
4796void 
4797CbcNode::decrementCuts(int change)
4798{
4799  if (nodeInfo_)
4800    assert ((state_&2)!=0);
4801  else
4802    assert ((state_&2)==0);
4803  if(nodeInfo_) {
4804    nodeInfo_->decrementCuts(change);
4805  }
4806}
4807void 
4808CbcNode::decrementParentCuts(CbcModel * model, int change)
4809{
4810  if (nodeInfo_)
4811    assert ((state_&2)!=0);
4812  else
4813    assert ((state_&2)==0);
4814  if(nodeInfo_) {
4815    nodeInfo_->decrementParentCuts(model, change);
4816  }
4817}
4818
4819/*
4820  Initialize reference counts (numberPointingToThis, numberBranchesLeft_)
4821  in the attached nodeInfo_.
4822*/
4823void
4824CbcNode::initializeInfo()
4825{
4826  assert(nodeInfo_ && branch_) ;
4827  nodeInfo_->initializeInfo(branch_->numberBranches());
4828  assert ((state_&2)!=0);
4829  assert (nodeInfo_->numberBranchesLeft()==
4830          branch_->numberBranchesLeft());
4831}
4832// Nulls out node info
4833void 
4834CbcNode::nullNodeInfo()
4835{
4836  nodeInfo_=NULL;
4837  // say not active
4838  state_ &= ~2;
4839}
4840
4841int
4842CbcNode::branch(OsiSolverInterface * solver)
4843{
4844  double changeInGuessed;
4845  assert (nodeInfo_->numberBranchesLeft()==
4846          branch_->numberBranchesLeft());
4847  if (!solver)
4848    changeInGuessed=branch_->branch();
4849  else
4850    changeInGuessed=branch_->branch(solver);
4851  guessedObjectiveValue_+= changeInGuessed;
4852  //#define PRINTIT
4853#ifdef PRINTIT
4854  int numberLeft = nodeInfo_->numberBranchesLeft();
4855  CbcNodeInfo * parent = nodeInfo_->parent();
4856  int parentNodeNumber = -1;
4857  //CbcBranchingObject * object1 = branch_->object_;
4858  //OsiObject * object = object1->
4859  //int sequence = object->columnNumber);
4860  int id=-1;
4861  double value=0.0;
4862  if (branch_) {
4863    id = branch_->variable();
4864    value = branch_->value();
4865  }
4866  printf("id %d value %g objvalue %g\n",id,value,objectiveValue_);
4867  if (parent)
4868    parentNodeNumber = parent->nodeNumber();
4869  printf("Node number %d, %s, way %d, depth %d, parent node number %d\n",
4870         nodeInfo_->nodeNumber(),(numberLeft==2) ? "leftBranch" : "rightBranch",
4871         way(),depth_,parentNodeNumber);
4872#endif
4873  return nodeInfo_->branchedOn();
4874}
4875/* Active arm of the attached OsiBranchingObject.
4876 
4877   In the simplest instance, coded -1 for the down arm of the branch, +1 for
4878   the up arm. But see OsiBranchingObject::way()
4879     Use nodeInfo--.numberBranchesLeft_ to see how active
4880*/
4881int 
4882CbcNode::way() const
4883{
4884  if (branch_) {
4885    CbcBranchingObject * obj =
4886      dynamic_cast <CbcBranchingObject *>(branch_) ;
4887    if (obj) {
4888      return obj->way();
4889    } else {
4890      OsiTwoWayBranchingObject * obj2 =
4891      dynamic_cast <OsiTwoWayBranchingObject *>(branch_) ;
4892      assert (obj2);
4893      return obj2->way();
4894    }
4895  } else {
4896    return 0;
4897  }
4898}
4899/* Create a branching object for the node
4900
4901    The routine scans the object list of the model and selects a set of
4902    unsatisfied objects as candidates for branching. The candidates are
4903    evaluated, and an appropriate branch object is installed.
4904
4905    The numberPassesLeft is decremented to stop fixing one variable each time
4906    and going on and on (e.g. for stock cutting, air crew scheduling)
4907
4908    If evaluation determines that an object is monotone or infeasible,
4909    the routine returns immediately. In the case of a monotone object,
4910    the branch object has already been called to modify the model.
4911
4912    Return value:
4913    <ul>
4914      <li>  0: A branching object has been installed
4915      <li> -1: A monotone object was discovered
4916      <li> -2: An infeasible object was discovered
4917    </ul>
4918    Branch state:
4919    <ul>
4920      <li> -1: start
4921      <li> -1: A monotone object was discovered
4922      <li> -2: An infeasible object was discovered
4923    </ul>
4924*/
4925int 
4926CbcNode::chooseOsiBranch (CbcModel * model,
4927                          CbcNode * lastNode,
4928                          OsiBranchingInformation * usefulInfo,
4929                          int branchState)
4930{
4931  int returnStatus=0;
4932  if (lastNode)
4933    depth_ = lastNode->depth_+1;
4934  else
4935    depth_ = 0;
4936  OsiSolverInterface * solver = model->solver();
4937  objectiveValue_ = solver->getObjValue()*solver->getObjSense();
4938  usefulInfo->objectiveValue_ = objectiveValue_;
4939  usefulInfo->depth_ = depth_;
4940  const double * saveInfoSol = usefulInfo->solution_;
4941  double * saveSolution = new double[solver->getNumCols()];
4942  memcpy(saveSolution,solver->getColSolution(),solver->getNumCols()*sizeof(double));
4943  usefulInfo->solution_ = saveSolution;
4944  OsiChooseVariable * choose = model->branchingMethod()->chooseMethod();
4945  int numberUnsatisfied=-1;
4946  if (branchState<0) {
4947    // initialize
4948    // initialize sum of "infeasibilities"
4949    sumInfeasibilities_ = 0.0;
4950    numberUnsatisfied = choose->setupList(usefulInfo,true);
4951    numberUnsatisfied_ = numberUnsatisfied;
4952    branchState=0;
4953    if (numberUnsatisfied_<0) {
4954      // infeasible
4955      delete [] saveSolution;
4956      return -2;
4957    }
4958  }
4959  // unset best
4960  int best=-1;
4961  choose->setBestObjectIndex(-1);
4962  if (numberUnsatisfied) {
4963    if (branchState>0||!choose->numberOnList()) {
4964      // we need to return at once - don't do strong branching or anything
4965      if (choose->numberOnList()||!choose->numberStrong()) {
4966        best = choose->candidates()[0];
4967        choose->setBestObjectIndex(best);
4968      } else {
4969        // nothing on list - need to try again - keep any solution
4970        numberUnsatisfied = choose->setupList(usefulInfo, false);
4971        numberUnsatisfied_ = numberUnsatisfied;
4972        if (numberUnsatisfied) {
4973          best = choose->candidates()[0];
4974          choose->setBestObjectIndex(best);
4975        }
4976      }
4977    } else {
4978      // carry on with strong branching or whatever
4979      int returnCode = choose->chooseVariable(solver, usefulInfo,true);
4980      // update number of strong iterations etc
4981      model->incrementStrongInfo(choose->numberStrongDone(),choose->numberStrongIterations(),
4982                                 returnCode==-1 ? 0:choose->numberStrongFixed(),returnCode==-1);
4983      if (returnCode>1) {
4984        // has fixed some
4985        returnStatus=-1;
4986      } else if (returnCode==-1) {
4987        // infeasible
4988        returnStatus=-2;
4989      } else if (returnCode==0) {
4990        // normal
4991        returnStatus=0;
4992        numberUnsatisfied=1;
4993      } else {
4994        // ones on list satisfied - double check
4995        numberUnsatisfied = choose->setupList(usefulInfo, false);
4996        numberUnsatisfied_ = numberUnsatisfied;
4997        if (numberUnsatisfied) {
4998          best = choose->candidates()[0];
4999          choose->setBestObjectIndex(best);
5000        }
5001      }
5002    }
5003  } 
5004  delete branch_;
5005  branch_ = NULL;
5006  guessedObjectiveValue_ = COIN_DBL_MAX;//objectiveValue_; // for now
5007  if (!returnStatus) {
5008    if (numberUnsatisfied) {
5009      // create branching object
5010      const OsiObject * obj = model->solver()->object(choose->bestObjectIndex());
5011      //const OsiSolverInterface * solver = usefulInfo->solver_;
5012      branch_ = obj->createBranch(model->solver(),usefulInfo,obj->whichWay());
5013    }
5014  }
5015  usefulInfo->solution_=saveInfoSol;
5016  delete [] saveSolution;
5017  // may have got solution
5018  if (choose->goodSolution()
5019      &&model->problemFeasibility()->feasible(model,-1)>=0) {
5020    // yes
5021    double objValue = choose->goodObjectiveValue();
5022    model->setBestSolution(CBC_STRONGSOL,
5023                                     objValue,
5024                                     choose->goodSolution()) ;
5025    model->setLastHeuristic(NULL);
5026    model->incrementUsed(choose->goodSolution());
5027    choose->clearGoodSolution();
5028  }
5029  return returnStatus;
5030}
5031int 
5032CbcNode::chooseClpBranch (CbcModel * model,
5033                       CbcNode * lastNode)
5034{ 
5035  assert(lastNode);
5036  depth_ = lastNode->depth_+1;
5037  delete branch_;
5038  branch_=NULL;
5039  OsiSolverInterface * solver = model->solver();
5040  //double saveObjectiveValue = solver->getObjValue();
5041  //double objectiveValue = CoinMax(solver->getObjSense()*saveObjectiveValue,objectiveValue_);
5042  const double * lower = solver->getColLower();
5043  const double * upper = solver->getColUpper();
5044  // point to useful information
5045  OsiBranchingInformation usefulInfo = model->usefulInformation();
5046  // and modify
5047  usefulInfo.depth_=depth_;
5048  int i;
5049  //bool beforeSolution = model->getSolutionCount()==0;
5050  int numberObjects = model->numberObjects();
5051  int numberColumns = model->getNumCols();
5052  double * saveUpper = new double[numberColumns];
5053  double * saveLower = new double[numberColumns];
5054
5055  // Save solution in case heuristics need good solution later
5056 
5057  double * saveSolution = new double[numberColumns];
5058  memcpy(saveSolution,solver->getColSolution(),numberColumns*sizeof(double));
5059  model->reserveCurrentSolution(saveSolution);
5060  for (i=0;i<numberColumns;i++) {
5061    saveLower[i] = lower[i];
5062    saveUpper[i] = upper[i];
5063  }
5064  // Save basis
5065  CoinWarmStart * ws = solver->getWarmStart();
5066  numberUnsatisfied_ = 0;
5067  // initialize sum of "infeasibilities"
5068  sumInfeasibilities_ = 0.0;
5069  // Note looks as if off end (hidden one)
5070  OsiObject * object = model->modifiableObject(numberObjects);
5071  CbcGeneralDepth * thisOne = dynamic_cast <CbcGeneralDepth *> (object);
5072  assert (thisOne);
5073  OsiClpSolverInterface * clpSolver
5074    = dynamic_cast<OsiClpSolverInterface *> (solver);
5075  assert (clpSolver);
5076  ClpSimplex * simplex = clpSolver->getModelPtr();
5077  int preferredWay;
5078  double infeasibility = object->infeasibility(&usefulInfo,preferredWay);
5079  if (thisOne->whichSolution()>=0) {
5080    ClpNode * nodeInfo = thisOne->nodeInfo(thisOne->whichSolution());
5081    nodeInfo->applyNode(simplex,2);
5082    simplex->dual();
5083    double cutoff = model->getCutoff();
5084    bool goodSolution=true;
5085    if (simplex->status()) {
5086      simplex->writeMps("bad7.mps",2);
5087      if (nodeInfo->objectiveValue()>cutoff-1.0e-2)
5088        goodSolution=false;
5089      else
5090        assert (!simplex->status());
5091    }
5092    if (goodSolution) {
5093      double newObjectiveValue = solver->getObjSense()*solver->getObjValue();
5094      // See if integer solution
5095      int numInf;
5096      int numInf2;
5097      bool gotSol = model->feasibleSolution(numInf,numInf2);
5098      if (!gotSol) {
5099        printf("numinf %d\n",numInf);
5100        double * sol = simplex->primalColumnSolution();
5101        for (int i=0;i<numberColumns;i++) {
5102          if (simplex->isInteger(i)) {
5103            double value = floor(sol[i]+0.5);
5104            if (fabs(value-sol[i])>1.0e-7) {
5105              printf("%d value %g\n",i,sol[i]);
5106              if (fabs(value-sol[i])<1.0e-3) {
5107                sol[i]=value;
5108              }
5109            }
5110          }
5111        }
5112        simplex->writeMps("bad8.mps",2);
5113        bool gotSol = model->feasibleSolution(numInf,numInf2);
5114        if (!gotSol) 
5115          assert (gotSol);
5116      }
5117      model->setBestSolution(CBC_STRONGSOL,
5118                             newObjectiveValue,
5119                             solver->getColSolution()) ;
5120      model->setLastHeuristic(NULL);
5121      model->incrementUsed(solver->getColSolution());
5122    }
5123  }
5124  // restore bounds
5125  { for (int j=0;j<numberColumns;j++) {
5126      if (saveLower[j] != lower[j])
5127        solver->setColLower(j,saveLower[j]);
5128      if (saveUpper[j] != upper[j])
5129        solver->setColUpper(j,saveUpper[j]);
5130    }
5131  }
5132  // restore basis
5133  solver->setWarmStart(ws);
5134  delete ws;
5135  int anyAction;
5136  //#define CHECK_PATH
5137#ifdef CHECK_PATH
5138extern int gotGoodNode_Z;
5139 if (gotGoodNode_Z>=0)
5140   printf("good node %d %g\n",gotGoodNode_Z,infeasibility);
5141#endif
5142  if (infeasibility>0.0) {
5143    if (infeasibility==COIN_DBL_MAX) {
5144      anyAction=-2; // infeasible
5145    } else {
5146      branch_=thisOne->createBranch(preferredWay);
5147      // Set to firat one (and change when re-pushing)
5148      CbcGeneralBranchingObject * branch = 
5149        dynamic_cast <CbcGeneralBranchingObject *> (branch_);
5150      branch->state(objectiveValue_,sumInfeasibilities_,
5151                    numberUnsatisfied_,0);
5152      branch->setNode(this);
5153      anyAction=0;
5154    }
5155  } else {
5156    anyAction=-1;
5157  }
5158#ifdef CHECK_PATH
5159 gotGoodNode_Z=-1;
5160#endif
5161  // Set guessed solution value
5162  guessedObjectiveValue_ = objectiveValue_+1.0e-5;
5163  delete [] saveLower;
5164  delete [] saveUpper;
5165 
5166  // restore solution
5167  solver->setColSolution(saveSolution);
5168  delete [] saveSolution;
5169  return anyAction;
5170}
5171/* Double checks in case node can change its mind!
5172   Returns objective value
5173   Can change objective etc */
5174double
5175CbcNode::checkIsCutoff(double cutoff)
5176{
5177  branch_->checkIsCutoff(cutoff);
5178  return objectiveValue_;
5179}
Note: See TracBrowser for help on using the repository browser.