source: stable/2.2/Cbc/src/CbcNode.cpp @ 1094

Last change on this file since 1094 was 1094, checked in by forrest, 10 years ago

for nonlinear

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