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

Last change on this file since 1061 was 1061, checked in by forrest, 13 years ago

out printf

  • Property svn:eol-style set to native
  • Property svn:keywords set to Author Date Id Revision
File size: 179.4 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 = 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 = NULL; //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      //#define CBC_NODE7 3
2711#if CBC_NODE7>1
2712      double * activeWeight = NULL;
2713      if (usefulInfo.columnLength_) {
2714        double tolerance=1.0e-6;
2715        int numberRows = solver->getNumRows();
2716        activeWeight = new double [numberRows];
2717        for (int iRow = 0;iRow<numberRows;iRow++) {
2718          // could use pi to see if active or activity
2719#if 1
2720          if (usefulInfo.rowActivity_[iRow]>usefulInfo.rowUpper_[iRow]-tolerance
2721              ||usefulInfo.rowActivity_[iRow]<usefulInfo.rowLower_[iRow]+tolerance) {
2722            activeWeight[iRow]=0.0;
2723          } else {
2724            activeWeight[iRow]=-1.0;
2725          }
2726#else
2727          if (fabs(usefulInfo.pi_[iRow])>1.0e-6) {
2728            activeWeight[iRow]=0.0;
2729          } else {
2730            activeWeight[iRow]=-1.0;
2731          }
2732#endif
2733        }
2734        for (int iColumn=0;iColumn<numberColumns;iColumn++) {
2735          if (solver->isInteger(iColumn)) {
2736            double value = usefulInfo.solution_[iColumn];
2737            if (fabs(value-floor(value+0.5))>1.0e-6) {
2738              CoinBigIndex start = usefulInfo.columnStart_[iColumn];
2739              CoinBigIndex end = start + usefulInfo.columnLength_[iColumn];
2740              for (CoinBigIndex j=start;j<end;j++) {
2741                int iRow = usefulInfo.row_[j];
2742                if (activeWeight[iRow]>=0.0)
2743                  activeWeight[iRow] += 1.0;
2744              }
2745            }
2746          }
2747        }
2748        for (int iRow = 0;iRow<numberRows;iRow++) {
2749          if (activeWeight[iRow]>0.0) {
2750            // could use pi
2751            activeWeight[iRow] = 1.0/activeWeight[iRow];
2752#if 0
2753            activeWeight[iRow] *= fabs(usefulInfo.pi_[iRow]);
2754#endif
2755          }
2756        }
2757      }
2758#endif
2759      /* Problem type as set by user or found by analysis.  This will be extended
2760         0 - not known
2761         1 - Set partitioning <=
2762         2 - Set partitioning ==
2763         3 - Set covering
2764         4 - all +- 1 or all +1 and odd
2765      */
2766      int problemType = model->problemType();
2767#define PRINT_STUFF -1
2768      for (i=0;i<numberObjects;i++) {
2769        OsiObject * object = model->modifiableObject(i);
2770        CbcSimpleIntegerDynamicPseudoCost * dynamicObject =
2771          dynamic_cast <CbcSimpleIntegerDynamicPseudoCost *>(object) ;
2772        int preferredWay;
2773        double infeasibility = object->infeasibility(&usefulInfo,preferredWay);
2774        int priorityLevel = object->priority();
2775#define ZERO_ONE 0
2776#define ZERO_FAKE 1.0e20;
2777#if ZERO_ONE==1
2778        // branch on 0-1 first (temp)
2779        if (fabs(saveSolution[dynamicObject->columnNumber()])<1.0)
2780          priorityLevel--;
2781#endif
2782#if ZERO_ONE==2
2783        if (fabs(saveSolution[dynamicObject->columnNumber()])<1.0)
2784          infeasibility *= ZERO_FAKE;
2785#endif
2786        if (infeasibility) {
2787          int iColumn = numberColumns+i;
2788          bool gotDown=false;
2789          int numberThisDown = 0;
2790          bool gotUp=false;
2791          int numberThisUp = 0;
2792          // check branching method
2793          if (dynamicObject) {
2794            if (branchingMethod!=dynamicObject->method()) {
2795              if (branchingMethod==-1)
2796                branchingMethod = dynamicObject->method();
2797              else
2798                branchingMethod = 100;
2799            }
2800            iColumn = dynamicObject->columnNumber();
2801            gotDown=false;
2802            numberThisDown = dynamicObject->numberTimesDown();
2803            if (numberThisDown>=numberBeforeTrust)
2804              gotDown=true;
2805            gotUp=false;
2806            numberThisUp = dynamicObject->numberTimesUp();
2807            if (numberThisUp>=numberBeforeTrust)
2808              gotUp=true;
2809#ifdef CBC_NODE7
2810            if (usefulInfo.columnLength_&&iColumn>=0&&iColumn<numberColumns) {
2811              CoinBigIndex start = usefulInfo.columnStart_[iColumn];
2812              CoinBigIndex end = start + usefulInfo.columnLength_[iColumn];
2813              double upValue = 0.0;
2814              double downValue = 0.0;
2815#if CBC_NODE7==1||CBC_NODE7==3
2816              double value = usefulInfo.direction_*usefulInfo.objective_[iColumn];
2817              // Bias cost
2818              if (value>0.0)
2819                upValue += 1.5*value;
2820              else
2821                downValue -= 1.5*value;
2822              for (CoinBigIndex j=start;j<end;j++) {
2823                int iRow = usefulInfo.row_[j];
2824                value = -usefulInfo.pi_[iRow];
2825                if (value) {
2826                  value *= usefulInfo.elementByColumn_[j];
2827#if CBC_NODE7==3
2828                  value *= activeWeight[iRow];
2829#endif
2830                  if (value>0.0)
2831                    upValue += value;
2832                  else
2833                    downValue -= value;
2834                }
2835              }
2836              value = usefulInfo.solution_[iColumn];
2837              value -= floor(value);
2838              downValue *= value;
2839              upValue *= (1.0-value);
2840#if 0
2841              printf("%d inf %g ord %g %g shadow %g %g\n",
2842                     iColumn,infeasibility,
2843                     object->downEstimate(),object->upEstimate(),
2844                     downValue,upValue);
2845#endif
2846              //if (!numberThisDown||!numberThisUp)
2847              //if (!gotDown||!gotUp)
2848#if 0
2849              infeasibility = CoinMax(CoinMax(downValue,upValue),infeasibility);
2850#else
2851              infeasibility = CoinMax(downValue,upValue);
2852#endif
2853#elif CBC_NODE7==2
2854              double value=0.0;
2855              for (CoinBigIndex j=start;j<end;j++) {
2856                int iRow = usefulInfo.row_[j];
2857                if (activeWeight[iRow])
2858                  value += fabs(usefulInfo.elementByColumn_[j])*activeWeight[iRow];
2859              }
2860              double solValue = usefulInfo.solution_[iColumn];
2861              solValue -= floor(solValue);
2862              downValue = value*solValue;
2863              upValue = value*(1.0-solValue);
2864#if 0
2865              printf("%d inf %g ord %g %g shadow %g %g\n",
2866                     iColumn,infeasibility,
2867                     object->downEstimate(),object->upEstimate(),
2868                     downValue,upValue);
2869#endif
2870              //if (!numberThisDown||!numberThisUp)
2871              //if (!gotDown||!gotUp)
2872              infeasibility = CoinMax(CoinMax(downValue,upValue),infeasibility);
2873              infeasibility = CoinMax(downValue,upValue);
2874              infeasibility = value;
2875#endif
2876            }
2877#endif
2878            //double gap = saveUpper[iColumn]-saveLower[iColumn];
2879            // Give precedence to ones with gap of 1.0
2880            //assert(gap>0.0);
2881            //infeasibility /= CoinMin(gap,100.0);
2882            if (!depth_&&false) {
2883              // try closest to 0.5
2884              double part =saveSolution[iColumn]-floor(saveSolution[iColumn]);
2885              infeasibility = fabs(0.5-part);
2886            }
2887            if (problemType>0&&problemType<4&&false) {
2888              // try closest to 0.5
2889              double part =saveSolution[iColumn]-floor(saveSolution[iColumn]);
2890              infeasibility = 0.5-fabs(0.5-part);
2891            }
2892            if (probingInfo) {
2893              int iSeq = backward[iColumn];
2894              assert (iSeq>=0);
2895              infeasibility = 1.0 + (toZero[iSeq+1]-toZero[iSeq])+ 
2896                5.0*CoinMin(toOne[iSeq]-toZero[iSeq],toZero[iSeq+1]-toOne[iSeq]);
2897              if (toZero[iSeq+1]>toZero[iSeq]) {
2898                numberUnsatisProbed++;
2899              } else {
2900                numberUnsatisNotProbed++;
2901              }
2902            }
2903            if ((numberNodes%PRINT_STUFF)==0&&PRINT_STUFF>0)
2904              printf("%d down %d %g up %d %g - infeas %g\n",
2905                     i,numberThisDown,object->downEstimate(),numberThisUp,object->upEstimate(),
2906                     infeasibility);
2907          } else {
2908            // see if SOS
2909            CbcSOS * sosObject =
2910              dynamic_cast <CbcSOS *>(object) ;
2911            assert (sosObject);
2912            gotDown=false;
2913            numberThisDown = sosObject->numberTimesDown();
2914            if (numberThisDown>=numberBeforeTrust)
2915              gotDown=true;
2916            gotUp=false;
2917            numberThisUp = sosObject->numberTimesUp();
2918            if (numberThisUp>=numberBeforeTrust)
2919              gotUp=true;
2920          }
2921          // Increase estimated degradation to solution
2922          estimatedDegradation += CoinMin(object->upEstimate(),object->downEstimate());
2923          downEstimate[i]=object->downEstimate();
2924          upEstimate[i]=object->upEstimate();
2925          numberUnsatisfied_++;
2926          sumInfeasibilities_ += infeasibility;
2927          // Better priority? Flush choices.
2928          if (priorityLevel<bestPriority) {
2929            numberToDo=0;
2930            bestPriority = priorityLevel;
2931            iBestGot=-1;
2932            best=0.0;
2933            numberNotTrusted=0;
2934          } else if (priorityLevel>bestPriority) {
2935            continue;
2936          }
2937          if (!gotUp||!gotDown)
2938            numberNotTrusted++;
2939          // Check for suitability based on infeasibility.
2940          if ((gotDown&&gotUp)&&numberStrong>0) {
2941            sort[numberToDo]=-infeasibility;
2942            if (infeasibility>best) {
2943              best=infeasibility;
2944              iBestGot=numberToDo;
2945            }
2946          } else {
2947            objectMark[neededPenalties]=numberToDo;
2948            which[neededPenalties++]=iColumn;
2949            sort[numberToDo]=-10.0*infeasibility;
2950            if (!(numberThisUp+numberThisDown))
2951              sort[numberToDo] *= 100.0; // make even more likely
2952            if (iColumn<numberColumns) {
2953              double part =saveSolution[iColumn]-floor(saveSolution[iColumn]);
2954              if (1.0-fabs(part-0.5)>bestNot) {
2955                iBestNot=numberToDo;
2956                bestNot = 1.0-fabs(part-0.5);
2957              }
2958            } else {
2959              // SOS
2960              if (-sort[numberToDo]>bestNot) {
2961                iBestNot=numberToDo;
2962                bestNot = -sort[numberToDo];
2963              }
2964            }
2965          }
2966          if (model->messageHandler()->logLevel()>3) { 
2967            printf("%d (%d) down %d %g up %d %g - infeas %g - sort %g solution %g\n",
2968                   i,iColumn,numberThisDown,object->downEstimate(),numberThisUp,object->upEstimate(),
2969                   infeasibility,sort[numberToDo],saveSolution[iColumn]);
2970          }
2971          whichObject[numberToDo++]=i;
2972        } else {
2973          // for debug
2974          downEstimate[i]=-1.0;
2975          upEstimate[i]=-1.0;
2976        }
2977      }
2978#if CBC_NODE7>1
2979      delete [] activeWeight;
2980#endif
2981      if (numberUnsatisfied_) {
2982        if (probingInfo&&false)
2983          printf("nunsat %d, %d probed, %d other 0-1\n",numberUnsatisfied_,
2984                 numberUnsatisProbed,numberUnsatisNotProbed);
2985        // some infeasibilities - go to next steps
2986        break;
2987      } else if (!iPass) {
2988        // may just need resolve
2989        solver->resolve();
2990        memcpy(saveSolution,solver->getColSolution(),numberColumns*sizeof(double));
2991        model->reserveCurrentSolution(saveSolution);
2992        if (!solver->isProvenOptimal()) {
2993          // infeasible
2994          anyAction=-2;
2995          break;
2996        }
2997      } else if (iPass==1) {
2998        // looks like a solution - get paranoid
2999        bool roundAgain=false;
3000        // get basis
3001        CoinWarmStartBasis * ws = dynamic_cast<CoinWarmStartBasis*>(solver->getWarmStart());
3002        if (!ws)
3003          break;
3004        double tolerance;
3005        solver->getDblParam(OsiPrimalTolerance,tolerance);
3006        for (i=0;i<numberColumns;i++) {
3007          double value = saveSolution[i];
3008          if (value<lower[i]-tolerance) {
3009            saveSolution[i]=lower[i];
3010            roundAgain=true;
3011            ws->setStructStatus(i,CoinWarmStartBasis::atLowerBound);
3012          } else if (value>upper[i]+tolerance) {
3013            saveSolution[i]=upper[i];
3014            roundAgain=true;
3015            ws->setStructStatus(i,CoinWarmStartBasis::atUpperBound);
3016          } 
3017        }
3018        if (roundAgain) {
3019          // restore basis
3020          solver->setWarmStart(ws);
3021          solver->setColSolution(saveSolution);
3022          delete ws;
3023          bool takeHint;
3024          OsiHintStrength strength;
3025          solver->getHintParam(OsiDoDualInResolve,takeHint,strength);
3026          solver->setHintParam(OsiDoDualInResolve,false,OsiHintDo) ;
3027          solver->resolve();
3028          solver->setHintParam(OsiDoDualInResolve,takeHint,strength) ;
3029          memcpy(saveSolution,solver->getColSolution(),numberColumns*sizeof(double));
3030          model->reserveCurrentSolution(saveSolution);
3031          if (!solver->isProvenOptimal()) {
3032            // infeasible
3033            anyAction=-2;
3034            break;
3035          }
3036        } else {
3037          delete ws;
3038          break;
3039        }
3040      }
3041    }
3042    if (anyAction==-2) {
3043      break;
3044    }
3045    bool solveAll=false; // set true to say look at all even if some fixed (experiment)
3046    solveAll=true;
3047    // skip if solution
3048    if (!numberUnsatisfied_)
3049      break;
3050    //bool skipAll = (numberBeforeTrust>20&&numberNodes>20000&&numberNotTrusted==0);
3051    bool skipAll = numberNotTrusted==0||numberToDo==1;
3052    bool doneHotStart=false;
3053    int searchStrategy = saveSearchStrategy>=0 ? (saveSearchStrategy%10) : -1;
3054    if (0) {
3055      int numberIterations = model->getIterationCount();
3056      int numberStrongIterations = model->numberStrongIterations();
3057      printf("node %d saveSearch %d search %d - its %d strong %d\n",
3058             numberNodes,saveSearchStrategy,searchStrategy,
3059             numberIterations,numberStrongIterations);
3060    }
3061#ifndef CBC_WEAK_STRONG
3062    if (((numberNodes%20)==0&&searchStrategy!=2)||(model->specialOptions()&8)!=0)
3063      skipAll=false;
3064#endif
3065    if (!newWay) {
3066    // 10 up always use %10, 20 up as 10 and allow penalties
3067    // But adjust depending on ratio of iterations
3068    if (searchStrategy>0&&saveSearchStrategy<10) {
3069      if (numberBeforeTrust>=5&&numberBeforeTrust<=10) {
3070        if (searchStrategy!=2) {
3071          if (depth_>5) {
3072            int numberIterations = model->getIterationCount();
3073            int numberStrongIterations = model->numberStrongIterations();
3074            if (numberStrongIterations>numberIterations+10000) {
3075              searchStrategy=2;
3076              //skipAll=true;
3077            } else if (numberStrongIterations*4+1000<numberIterations||depth_<5) {
3078              searchStrategy=3;
3079              skipAll=false;
3080            }
3081          } else {
3082            searchStrategy=3;
3083            skipAll=false;
3084          }
3085        } else {
3086          //skipAll=true;
3087        }
3088      }
3089    }
3090    } else {
3091    // But adjust depending on ratio of iterations
3092    if (saveSearchStrategy<0) {
3093      // unset
3094      if ((numberNodes%20)==0||(model->specialOptions()&8)!=0) {
3095        // Do numberStrong
3096        searchStrategy=3;
3097      } else if (depth_<5) {
3098        // Do numberStrong
3099        searchStrategy=2;
3100      } else {
3101        int numberIterations = model->getIterationCount();
3102        int numberStrongIterations = model->numberStrongIterations();
3103        int numberRows = solver->getNumRows();
3104        if (numberStrongIterations>numberIterations+CoinMin(10000,10*numberRows)) {
3105          // off
3106          searchStrategy=0;
3107        } else if (numberStrongIterations*4+1000<numberIterations) {
3108          // Do numberStrong if not trusted
3109          searchStrategy=2;
3110        } else {
3111          searchStrategy=1;
3112        }
3113      }
3114    }
3115    if (searchStrategy<3&&(!numberNotTrusted||!searchStrategy))
3116      skipAll=true;
3117    else
3118      skipAll=false;
3119    }
3120    // worth trying if too many times
3121    // Save basis
3122    CoinWarmStart * ws = NULL;
3123    // save limit
3124    int saveLimit=0;
3125    solver->getIntParam(OsiMaxNumIterationHotStart,saveLimit);
3126    if (!skipAll) {
3127      ws = solver->getWarmStart();
3128      int limit=100;
3129#if 0
3130      int averageBranchIterations = model->getIterationCount()/(model->getNodeCount()+1);
3131      if (numberNodes)
3132        limit = CoinMin(CoinMax(limit,2*averageBranchIterations),500);
3133      else
3134        limit = 500;
3135#endif
3136      if ((!saveStateOfSearch||searchStrategy>3)&&saveLimit<limit&&saveLimit==100)
3137        solver->setIntParam(OsiMaxNumIterationHotStart,limit); 
3138    }
3139    // Say which one will be best
3140    int whichChoice=0;
3141    int bestChoice;
3142    if (iBestGot>=0)
3143      bestChoice=iBestGot;
3144    else
3145      bestChoice=iBestNot;
3146    assert (bestChoice>=0);
3147    // If we have hit max time don't do strong branching
3148    bool hitMaxTime = ( CoinCpuTime()-model->getDblParam(CbcModel::CbcStartSeconds) > 
3149                        model->getDblParam(CbcModel::CbcMaximumSeconds));
3150    // also give up if we are looping round too much
3151    if (hitMaxTime||numberPassesLeft<=0||(!numberNotTrusted&&false)||branchingMethod==11) {
3152      int iObject = whichObject[bestChoice];
3153      OsiObject * object = model->modifiableObject(iObject);
3154      int preferredWay;
3155      object->infeasibility(&usefulInfo,preferredWay);
3156      CbcSimpleInteger * obj =
3157        dynamic_cast <CbcSimpleInteger *>(object) ;
3158      if (obj) {
3159        branch_=obj->createBranch(solver,&usefulInfo,preferredWay);
3160      } else {
3161        CbcObject * obj =
3162          dynamic_cast <CbcObject *>(object) ;
3163        assert (obj);
3164        branch_=obj->createBranch(preferredWay);
3165      }
3166      {
3167        CbcBranchingObject * branchObj =
3168          dynamic_cast <CbcBranchingObject *>(branch_) ;
3169        assert (branchObj);
3170        branchObj->way(preferredWay);
3171      }
3172      delete ws;
3173      ws=NULL;
3174      break;
3175    } else {
3176      // say do fast
3177      int easy=1;
3178      if (!skipAll)
3179        solver->setHintParam(OsiDoInBranchAndCut,true,OsiHintDo,&easy) ;
3180      int iDo;
3181#ifdef RANGING
3182      if (useRanging) {
3183      if ((skipAll&&numberBeforeTrust&&saveSearchStrategy<20)||saveSearchStrategy<10)
3184        numberPenalties=0;
3185      {
3186        // off penalties if too much
3187        double needed = neededPenalties;
3188        needed *= numberRows;
3189        if (needed>1.0e6&&numberNodes&&saveSearchStrategy<20) {
3190          numberPenalties=0;
3191          neededPenalties=0;
3192        }
3193      }
3194#     ifdef COIN_HAS_CLP
3195      if (osiclp&&numberPenalties&&neededPenalties) {
3196        xPen += neededPenalties;
3197        which--;
3198        which[0]=neededPenalties;
3199        osiclp->passInRanges(which);
3200        // Mark hot start and get ranges
3201        if (kPass) {
3202          // until can work out why solution can go funny
3203          int save = osiclp->specialOptions();
3204          osiclp->setSpecialOptions(save|256);
3205          solver->markHotStart();
3206          osiclp->setSpecialOptions(save);
3207        } else {
3208          solver->markHotStart();
3209        }
3210        doneHotStart=true;
3211        xMark++;
3212        kPass++;
3213        osiclp->passInRanges(NULL);
3214        const double * downCost=osiclp->upRange();
3215        const double * upCost=osiclp->downRange();
3216        //printf("numberTodo %d needed %d numberpenalties %d\n",numberToDo,neededPenalties,numberPenalties);
3217        double invTrust = 1.0/((double) numberBeforeTrust);
3218        for (int i=0;i<neededPenalties;i++) {
3219          int j = objectMark[i];
3220          int iObject = whichObject[j];
3221          OsiObject * object = model->modifiableObject(iObject);
3222          CbcSimpleIntegerDynamicPseudoCost * dynamicObject =
3223            dynamic_cast <CbcSimpleIntegerDynamicPseudoCost *>(object) ;
3224          int iSequence=dynamicObject->columnNumber();
3225          double value = saveSolution[iSequence];
3226          value -= floor(value);
3227          double upPenalty = CoinMin(upCost[i],1.0e110)*(1.0-value);
3228          double downPenalty = CoinMin(downCost[i],1.0e110)*value;
3229          if (!numberBeforeTrust) {
3230            // override
3231            downEstimate[iObject]=downPenalty;
3232            upEstimate[iObject]=upPenalty;
3233          } else {
3234            int numberThisDown = dynamicObject->numberTimesDown();
3235            if (numberThisDown<numberBeforeTrust) {
3236              double fraction = ((double) numberThisDown)*invTrust;
3237              downEstimate[iObject] = fraction*downEstimate[iObject]+(1.0-fraction)*downPenalty;
3238            }
3239            int numberThisUp = dynamicObject->numberTimesUp();
3240            if (numberThisUp<numberBeforeTrust) {
3241              double fraction = ((double) numberThisUp)*invTrust;
3242              upEstimate[iObject] = fraction*upEstimate[iObject]+(1.0-fraction)*upPenalty;
3243            }
3244          }
3245          sort[j] = - CoinMin(downEstimate[iObject],upEstimate[iObject]);
3246#ifdef CBC_WEAK_STRONG
3247          sort[j] -= 1.0e10; // make more likely to be chosen
3248#endif
3249          //if ((numberNodes%PRINT_STUFF)==0&&PRINT_STUFF>0)
3250          if (!numberNodes)
3251            printf("%d pen down ps %g -> %g up ps %g -> %g\n",
3252                   iObject,downCost[i],downPenalty,upCost[i],upPenalty);
3253        }
3254      } else
3255#     endif     /* COIN_HAS_CLP */
3256      {
3257        if (!skipAll) {
3258          // Mark hot start
3259          solver->markHotStart();
3260          doneHotStart=true;
3261          xMark++;
3262          //if (solver->isProvenPrimalInfeasible())
3263          //printf("**** Hot start says node infeasible\n");
3264        }
3265        // make sure best will be first
3266        if (iBestGot>=0)
3267          sort[iBestGot]=-1.0e120;
3268      }
3269      } else {
3270#endif          /* RANGING */
3271      if (!skipAll) {
3272        // Mark hot start
3273        doneHotStart=true;
3274        solver->markHotStart();
3275        xMark++;
3276      }
3277      // make sure best will be first
3278      if (iBestGot>=0)
3279        sort[iBestGot]=-COIN_DBL_MAX;
3280#ifdef RANGING
3281      }
3282#endif          /* RANGING */
3283      // Actions 0 - exit for repeat, 1 resolve and try old choice,2 exit for continue
3284#define ACTION 0
3285#if ACTION<2
3286      if (anyAction)
3287        numberToDo=0; // skip as we will be trying again
3288#endif
3289      // Sort
3290      CoinSort_2(sort,sort+numberToDo,whichObject);
3291      // Change in objective opposite infeasible
3292      double worstFeasible=0.0;
3293      // Just first if strong off
3294      if (!numberStrong)
3295        numberToDo=CoinMin(numberToDo,1);
3296#if NODE_NEW
3297      if (searchStrategy==2)
3298        numberToDo=CoinMin(numberToDo,20);
3299#endif
3300      iDo=0;
3301      int saveLimit2;
3302      solver->getIntParam(OsiMaxNumIterationHotStart,saveLimit2);
3303      bool doQuickly = false; // numberToDo>2*numberStrong;
3304      if (searchStrategy==2)
3305        doQuickly=true;
3306      //printf("todo %d, strong %d\n",numberToDo,numberStrong);
3307      int numberTest=numberNotTrusted>0 ? numberStrong : (numberStrong+1)/2;
3308      int numberTest2 = 2*numberStrong;
3309      //double distanceToCutoff2 = model->getCutoff()-objectiveValue_;
3310      if (!newWay) {
3311      if (searchStrategy==3) {
3312        // Previously decided we need strong
3313        doQuickly=false;
3314        numberTest = numberStrong;
3315        //numberTest2 = 1000000;
3316      }
3317      //if (searchStrategy<0||searchStrategy==1)
3318        //numberTest2 = 1000000;
3319#if 0
3320      if (numberBeforeTrust>20&&(numberNodes>20000||(numberNodes>200&&numberNotTrusted==0))) {
3321        if ((numberNodes%20)!=0) {
3322          numberTest=0;
3323          doQuickly=true;
3324        }
3325      }
3326#else
3327      // Try nearly always off
3328      if (searchStrategy<2) {
3329        if ((numberNodes%20)!=0) {
3330          if ((model->specialOptions()&8)==0) {
3331            numberTest=0;
3332            doQuickly=true;
3333          }
3334        } else {
3335          doQuickly=false;
3336          numberTest=2*numberStrong;
3337          skipAll=false;
3338        }
3339      } else if (searchStrategy!=3) {
3340        doQuickly=true;
3341        numberTest=numberStrong;
3342      }
3343#endif
3344      if (depth_<8&&numberStrong) {
3345        if (searchStrategy!=2) {
3346          doQuickly=false;
3347          int numberRows = solver->getNumRows();
3348          // whether to do this or not is important - think
3349          if (numberRows<300||numberRows+numberColumns<2500) {
3350            if (depth_<7)
3351              numberStrong = CoinMin(3*numberStrong,numberToDo);
3352            if (!depth_) 
3353              numberStrong=CoinMin(6*numberStrong,numberToDo);
3354          }
3355          numberTest=numberStrong;
3356          skipAll=false;
3357        }
3358        //model->setStateOfSearch(2); // use min min
3359      }
3360      // could adjust using average iterations per branch
3361      // double average = ((double)model->getIterationCount())/
3362      //((double) model->getNodeCount()+1.0);
3363      // if too many and big then just do 10 its
3364      if (!skipAll&&saveStateOfSearch) {
3365        //if (numberNotTrusted>3*numberStrong&&numberRows>250&&numberColumns>1000&&saveLimit==100)
3366          // off solver->setIntParam(OsiMaxNumIterationHotStart,10);
3367      }
3368      // make negative for test
3369      distanceToCutoff = - distanceToCutoff;
3370      if (numberObjects>-100) {
3371        // larger
3372        distanceToCutoff *= 100.0;
3373      }
3374        distanceToCutoff = -COIN_DBL_MAX;
3375      // Do at least 5 strong
3376      if (numberColumns<1000&&(depth_<15||numberNodes<1000000))
3377        numberTest = CoinMax(numberTest,5);
3378      if ((model->specialOptions()&8)==0) {
3379        if (skipAll) {
3380          numberTest=0;
3381          doQuickly=true;
3382        }
3383      } else {
3384        // do 5 as strong is fixing
3385        numberTest = CoinMax(numberTest,5);
3386      }
3387      } else {
3388        if (!depth_&&numberToDo<200&&solver->getNumRows()<2000&&
3389            numberColumns<2000&&false) {
3390          numberStrong = numberToDo;
3391          doQuickly = false;
3392        }
3393      int numberTest=numberNotTrusted>0 ? numberStrong : (numberStrong+1)/2;
3394      int numberTest2 = 2*numberStrong;
3395      if (searchStrategy>=3) {
3396        // Previously decided we need strong
3397        doQuickly=false;
3398        if (depth_<7)
3399          numberStrong *=3;
3400        if (!depth_) 
3401          numberStrong=CoinMin(6*numberStrong,numberToDo);
3402        numberTest = numberStrong;
3403        numberTest2 *= 2;
3404      } else if (searchStrategy==2||(searchStrategy==1&&depth_<6)) {
3405        numberStrong *=2;
3406        if (!depth_) 
3407          numberStrong=CoinMin(2*numberStrong,numberToDo);
3408        numberTest = numberStrong;
3409      } else if (searchStrategy==1&&numberNotTrusted) {
3410        numberTest = numberStrong;
3411      } else {
3412        numberTest=0;
3413        skipAll=true;
3414      }
3415      distanceToCutoff=model->getCutoff()-objectiveValue_;
3416      // make negative for test
3417      distanceToCutoff = - distanceToCutoff;
3418      if (numberObjects>-100) {
3419        // larger
3420        distanceToCutoff *= 100.0;
3421      }
3422      distanceToCutoff = -COIN_DBL_MAX;
3423      if (skipAll) {
3424        numberTest=0;
3425        doQuickly=true;
3426      }
3427      }
3428      // temp - always switch off
3429      if (0) {
3430        int numberIterations = model->getIterationCount();
3431        int numberStrongIterations = model->numberStrongIterations();
3432        int numberRows = solver->getNumRows();
3433        if (numberStrongIterations>numberIterations+CoinMin(100,10*numberRows)&&depth_>=5) {
3434          skipAll=true;
3435          newWay=false;
3436          numberTest=0;
3437          doQuickly=true;
3438        }
3439      }
3440#if 0
3441      // temp - always switch on
3442      if (0) {
3443        int numberIterations = model->getIterationCount();
3444        int numberStrongIterations = model->numberStrongIterations();
3445        if (2*numberStrongIterations<numberIterations||depth_<=5) {
3446          skipAll=false;
3447          newWay=false;
3448          numberTest=CoinMax(numberTest,numberStrong);
3449          doQuickly=false;
3450        }
3451      }
3452#endif
3453      px[0]=numberTest;
3454      px[1]=numberTest2;
3455      px[2]= doQuickly ? 1 : -1;
3456      px[3]=numberStrong;
3457      if (!newWay) {
3458        if (numberColumns>8*solver->getNumRows()&&false) {
3459          printf("skipAll %c doQuickly %c numberTest %d numberTest2 %d numberNot %d\n",
3460                 skipAll ? 'Y' : 'N',doQuickly ? 'Y' : 'N',numberTest,numberTest2,numberNotTrusted);
3461          numberTest = CoinMin(numberTest,model->numberStrong());
3462          numberTest2 = CoinMin(numberTest2,model->numberStrong());
3463          printf("new test,test2 %d %d\n",numberTest,numberTest2);
3464        }
3465      }
3466      //printf("skipAll %c doQuickly %c numberTest %d numberTest2 %d numberNot %d\n",
3467      //   skipAll ? 'Y' : 'N',doQuickly ? 'Y' : 'N',numberTest,numberTest2,numberNotTrusted);
3468      // See if we want mini tree
3469      bool wantMiniTree=false;
3470      if (model->sizeMiniTree()&&depth_>7&&saveStateOfSearch>0)
3471        wantMiniTree=true;
3472      numberMini=0;
3473      //if (skipAll&&numberTest==0&&doQuickly)
3474      //numberToDo = 1; // trust previous stuff
3475      bool couldChooseFirst = false ; //(skipAll&&numberTest==0&&doQuickly);
3476      //skipAll=false;
3477      int realMaxHotIterations=999999;
3478#if 0
3479      if (saveSearchStrategy==-1&&!model->parentModel()&&
3480          !depth_&&saveLimit==100) {
3481        realMaxHotIterations=saveLimit;
3482        saveLimit2=200;
3483        solver->setIntParam(OsiMaxNumIterationHotStart,saveLimit2);
3484      }
3485#endif
3486      for ( iDo=0;iDo<numberToDo;iDo++) {
3487        int iObject = whichObject[iDo];
3488        OsiObject * object = model->modifiableObject(iObject);
3489        CbcSimpleIntegerDynamicPseudoCost * dynamicObject =
3490          dynamic_cast <CbcSimpleIntegerDynamicPseudoCost *>(object) ;
3491        int iColumn = dynamicObject ? dynamicObject->columnNumber() : numberColumns+iObject;
3492        int preferredWay;
3493        double infeasibility = object->infeasibility(&usefulInfo,preferredWay);
3494        // may have become feasible
3495        if (!infeasibility)
3496          continue;
3497        CbcSimpleInteger * obj =
3498          dynamic_cast <CbcSimpleInteger *>(object) ;
3499        if (obj) {
3500          if (choiceObject) {
3501            obj->fillCreateBranch(choiceObject,&usefulInfo,preferredWay);
3502            choiceObject->setObject(dynamicObject);
3503          } else {
3504            choice.possibleBranch=obj->createBranch(solver,&usefulInfo,preferredWay);
3505          }
3506        } else {
3507          CbcObject * obj =
3508            dynamic_cast <CbcObject *>(object) ;
3509          assert (obj);
3510          choice.possibleBranch=obj->createBranch(preferredWay);
3511        }
3512        // Save which object it was
3513        choice.objectNumber=iObject;
3514        choice.numIntInfeasUp = numberUnsatisfied_;
3515        choice.numIntInfeasDown = numberUnsatisfied_;
3516        choice.upMovement = upEstimate[iObject];
3517        choice.downMovement = downEstimate[iObject];
3518        assert (choice.upMovement>=0.0);
3519        assert (choice.downMovement>=0.0);
3520        choice.fix=0; // say not fixed
3521        double maxChange = 0.5*(choice.upMovement+choice.downMovement);
3522        maxChange = CoinMin(choice.upMovement,choice.downMovement);
3523        maxChange = CoinMax(choice.upMovement,choice.downMovement);
3524        if (searchStrategy==2)
3525          maxChange = COIN_DBL_MAX;
3526        //maxChange *= 5.0;
3527        if (dynamicObject&&dynamicObject->method()==1)
3528          maxChange *= 0.1; // probing
3529        // see if can skip strong branching
3530        int canSkip = choice.possibleBranch->fillStrongInfo(choice);
3531#if 0
3532        if (!newWay) {
3533          if ((maxChange>distanceToCutoff2)&&(!doQuickly||(numberTest>0&&searchStrategy!=2)))
3534          canSkip=0;
3535        } else {
3536        if (skipAll)
3537          canSkip=1;
3538        else if (numberTest>0&&searchStrategy>=3)
3539          canSkip=0;
3540        }
3541        if (!numberBeforeTrust) {
3542          canSkip=1;
3543        }
3544        if (sort[iDo]<distanceToCutoff)
3545          canSkip=0;
3546        if (((numberTest2<=0&&numberTest<=0)||skipAll)&&sort[iDo]>distanceToCutoff) {
3547          canSkip=1; // always skip
3548          if (iDo>20) {
3549            if (!choiceObject) {
3550              delete choice.possibleBranch;
3551              choice.possibleBranch=NULL;
3552            }
3553            break; // give up anyway
3554          }
3555        }
3556#else
3557        if (((numberTest2<=0&&numberTest<=0)||skipAll)&&sort[iDo]>distanceToCutoff) {
3558          //canSkip=1; // always skip
3559          if (iDo>20) {
3560            if (!choiceObject) {
3561              delete choice.possibleBranch;
3562              choice.possibleBranch=NULL;
3563            }
3564            break; // give up anyway
3565          }
3566        }
3567#endif
3568        if (model->messageHandler()->logLevel()>3&&numberBeforeTrust&&dynamicObject) 
3569          dynamicObject->print(1,choice.possibleBranch->value());
3570        // was if (!canSkip)
3571        if (newWay)
3572        numberTest2--;
3573        if (!canSkip) {
3574          //#ifndef RANGING
3575          if (!doneHotStart) {
3576            // Mark hot start
3577            doneHotStart=true;
3578            solver->markHotStart();
3579            xMark++;
3580          }
3581          //#endif
3582          assert (!couldChooseFirst);
3583          numberTest--;
3584          if (!newWay)
3585          numberTest2--;
3586          // just do a few
3587#if NODE_NEW == 5  || NODE_NEW == 2
3588          //if (canSkip)
3589          if (searchStrategy==2)
3590            solver->setIntParam(OsiMaxNumIterationHotStart,10); 
3591#endif
3592          double objectiveChange ;
3593          double newObjectiveValue=1.0e100;
3594          int j;
3595          // status is 0 finished, 1 infeasible and other
3596          int iStatus;
3597          if (0) {
3598            CbcDynamicPseudoCostBranchingObject * cbcobj = dynamic_cast<CbcDynamicPseudoCostBranchingObject *> (choice.possibleBranch);
3599            if (cbcobj) {
3600              CbcSimpleIntegerDynamicPseudoCost * object = cbcobj->object();
3601              printf("strong %d ",iDo);
3602              object->print(1,0.5);
3603            }
3604          }
3605          /*
3606            Try the down direction first. (Specify the initial branching alternative as
3607            down with a call to way(-1). Each subsequent call to branch() performs the
3608            specified branch and advances the branch object state to the next branch
3609            alternative.)
3610          */
3611          choice.possibleBranch->way(-1) ;
3612#if NEW_UPDATE_OBJECT==0
3613          decision->saveBranchingObject( choice.possibleBranch);
3614#endif
3615          choice.possibleBranch->branch() ;
3616          solver->solveFromHotStart() ;
3617          bool needHotStartUpdate=false;
3618          numberStrongDone++;
3619          numberStrongIterations += solver->getIterationCount();
3620          /*
3621            We now have an estimate of objective degradation that we can use for strong
3622            branching. If we're over the cutoff, the variable is monotone up.
3623            If we actually made it to optimality, check for a solution, and if we have
3624            a good one, call setBestSolution to process it. Note that this may reduce the
3625            cutoff, so we check again to see if we can declare this variable monotone.
3626          */
3627          if (solver->isProvenOptimal())
3628            iStatus=0; // optimal
3629          else if (solver->isIterationLimitReached()
3630                   &&!solver->isDualObjectiveLimitReached())
3631            iStatus=2; // unknown
3632          else
3633            iStatus=1; // infeasible
3634          if (iStatus!=2&&solver->getIterationCount()>
3635              realMaxHotIterations)
3636            numberUnfinished++;
3637          newObjectiveValue = solver->getObjSense()*solver->getObjValue();
3638          choice.numItersDown = solver->getIterationCount();
3639          objectiveChange = CoinMax(newObjectiveValue  - objectiveValue_,0.0);
3640          // Update branching information if wanted
3641#if NEW_UPDATE_OBJECT==0
3642          decision->updateInformation( solver,this);
3643#elif NEW_UPDATE_OBJECT<2
3644          CbcBranchingObject * cbcobj = dynamic_cast<CbcBranchingObject *> (choice.possibleBranch);
3645          if (cbcobj) {
3646            CbcObject * object = cbcobj->object();
3647            assert (object);
3648            CbcObjectUpdateData update = object->createUpdateInformation(solver,this,cbcobj);
3649            object->updateInformation(update);
3650          } else {
3651            decision->updateInformation( solver,this);
3652          }
3653#else
3654          CbcBranchingObject * cbcobj = dynamic_cast<CbcBranchingObject *> (choice.possibleBranch);
3655          if (cbcobj) {
3656            CbcObject * object = cbcobj->object();
3657            assert (object) ;
3658            CbcObjectUpdateData update = object->createUpdateInformation(solver,this,cbcobj);
3659            update.objectNumber_ = choice.objectNumber;
3660            model->addUpdateInformation(update);
3661          } else {
3662            decision->updateInformation( solver,this);
3663          }
3664#endif
3665          if (!iStatus) {
3666            choice.finishedDown = true ;
3667            if (newObjectiveValue>=cutoff) {
3668              objectiveChange = 1.0e100; // say infeasible
3669              numberStrongInfeasible++;
3670            } else {
3671              // See if integer solution
3672              if (model->feasibleSolution(choice.numIntInfeasDown,
3673                                          choice.numObjInfeasDown)
3674                  &&model->problemFeasibility()->feasible(model,-1)>=0) {
3675                if (auxiliaryInfo->solutionAddsCuts()) {
3676                  needHotStartUpdate=true;
3677                  solver->unmarkHotStart();
3678                }
3679                model->setBestSolution(CBC_STRONGSOL,
3680                                       newObjectiveValue,
3681                                       solver->getColSolution()) ;
3682                if (needHotStartUpdate) {
3683                  solver->resolve();
3684                  newObjectiveValue = solver->getObjSense()*solver->getObjValue();
3685                  objectiveChange = CoinMax(newObjectiveValue  - objectiveValue_,0.0);
3686                  model->feasibleSolution(choice.numIntInfeasDown,
3687                                          choice.numObjInfeasDown);
3688                }
3689                model->setLastHeuristic(NULL);
3690                model->incrementUsed(solver->getColSolution());
3691                cutoff =model->getCutoff();
3692                if (newObjectiveValue >= cutoff) {      //  *new* cutoff
3693                  objectiveChange = 1.0e100 ;
3694                  numberStrongInfeasible++;
3695                }
3696              }
3697            }
3698          } else if (iStatus==1) {
3699            objectiveChange = 1.0e100 ;
3700            numberStrongInfeasible++;
3701          } else {
3702            // Can't say much as we did not finish
3703            choice.finishedDown = false ;
3704            numberUnfinished++;
3705          }
3706          choice.downMovement = objectiveChange ;
3707         
3708          // restore bounds
3709          for ( j=0;j<numberColumns;j++) {
3710            if (saveLower[j] != lower[j])
3711              solver->setColLower(j,saveLower[j]);
3712            if (saveUpper[j] != upper[j])
3713              solver->setColUpper(j,saveUpper[j]);
3714          }
3715          if(needHotStartUpdate) {
3716            needHotStartUpdate = false;
3717            solver->resolve();
3718            //we may again have an integer feasible solution
3719            int numberIntegerInfeasibilities;
3720            int numberObjectInfeasibilities;
3721            if (model->feasibleSolution(
3722                                        numberIntegerInfeasibilities,
3723                                        numberObjectInfeasibilities)) {
3724#ifdef BONMIN
3725              //In this case node has become integer feasible, let us exit the loop
3726              std::cout<<"Node has become integer feasible"<<std::endl;
3727              numberUnsatisfied_ = 0;
3728              break;
3729#endif
3730              double objValue = solver->getObjValue();
3731              model->setBestSolution(CBC_STRONGSOL,
3732                                     objValue,
3733                                     solver->getColSolution()) ;
3734              solver->resolve();
3735              cutoff =model->getCutoff();
3736            }
3737            solver->markHotStart();
3738          }
3739          //printf("Down on %d, status is %d, obj %g its %d cost %g finished %d inf %d infobj %d\n",
3740          //printf("Down on %d, status is %d, obj %g its %d cost %g finished %d inf %d infobj %d\n",
3741          //     choice.objectNumber,iStatus,newObjectiveValue,choice.numItersDown,
3742          //     choice.downMovement,choice.finishedDown,choice.numIntInfeasDown,
3743          //     choice.numObjInfeasDown);
3744         
3745          // repeat the whole exercise, forcing the variable up
3746#if NEW_UPDATE_OBJECT==0
3747          decision->saveBranchingObject( choice.possibleBranch);
3748#endif
3749          choice.possibleBranch->branch();
3750          solver->solveFromHotStart() ;
3751          numberStrongDone++;
3752          numberStrongIterations += solver->getIterationCount();
3753          /*
3754            We now have an estimate of objective degradation that we can use for strong
3755            branching. If we're over the cutoff, the variable is monotone up.
3756            If we actually made it to optimality, check for a solution, and if we have
3757            a good one, call setBestSolution to process it. Note that this may reduce the
3758            cutoff, so we check again to see if we can declare this variable monotone.
3759          */
3760          if (solver->isProvenOptimal())
3761            iStatus=0; // optimal
3762          else if (solver->isIterationLimitReached()
3763                   &&!solver->isDualObjectiveLimitReached())
3764            iStatus=2; // unknown
3765          else
3766            iStatus=1; // infeasible
3767          if (iStatus!=2&&solver->getIterationCount()>
3768              realMaxHotIterations)
3769            numberUnfinished++;
3770          newObjectiveValue = solver->getObjSense()*solver->getObjValue();
3771          choice.numItersUp = solver->getIterationCount();
3772          objectiveChange = CoinMax(newObjectiveValue  - objectiveValue_,0.0);
3773          // Update branching information if wanted
3774#if NEW_UPDATE_OBJECT==0
3775          decision->updateInformation( solver,this);
3776#elif NEW_UPDATE_OBJECT<2
3777          cbcobj = dynamic_cast<CbcBranchingObject *> (choice.possibleBranch);
3778          if (cbcobj) {
3779            CbcObject * object = cbcobj->object();
3780            CbcObjectUpdateData update = object->createUpdateInformation(solver,this,cbcobj);
3781            object->updateInformation(update);
3782          } else {
3783            decision->updateInformation( solver,this);
3784          }
3785#else
3786          cbcobj = dynamic_cast<CbcBranchingObject *> (choice.possibleBranch);
3787          if (cbcobj) {
3788            CbcObject * object = cbcobj->object();
3789            assert (object) ;
3790            CbcObjectUpdateData update = object->createUpdateInformation(solver,this,cbcobj);
3791            update.objectNumber_ = choice.objectNumber;
3792            model->addUpdateInformation(update);
3793          } else {
3794            decision->updateInformation( solver,this);
3795          }
3796#endif
3797          if (!iStatus) {
3798            choice.finishedUp = true ;
3799            if (newObjectiveValue>=cutoff) {
3800              objectiveChange = 1.0e100; // say infeasible
3801              numberStrongInfeasible++;
3802            } else {
3803              // See if integer solution
3804              if (model->feasibleSolution(choice.numIntInfeasUp,
3805                                          choice.numObjInfeasUp)
3806                  &&model->problemFeasibility()->feasible(model,-1)>=0) {
3807#ifdef BONMIN
3808                std::cout<<"Node has become integer feasible"<<std::endl;
3809                numberUnsatisfied_ = 0;
3810                break;
3811#endif
3812                if (auxiliaryInfo->solutionAddsCuts()) {
3813                  needHotStartUpdate=true;
3814                  solver->unmarkHotStart();
3815                }
3816                model->setBestSolution(CBC_STRONGSOL,
3817                                       newObjectiveValue,
3818                                       solver->getColSolution()) ;
3819                if (needHotStartUpdate) {
3820                  solver->resolve();
3821                  newObjectiveValue = solver->getObjSense()*solver->getObjValue();
3822                  objectiveChange = CoinMax(newObjectiveValue  - objectiveValue_,0.0);
3823                  model->feasibleSolution(choice.numIntInfeasDown,
3824                                          choice.numObjInfeasDown);
3825                }
3826                model->setLastHeuristic(NULL);
3827                model->incrementUsed(solver->getColSolution());
3828                cutoff =model->getCutoff();
3829                if (newObjectiveValue >= cutoff) {      //  *new* cutoff
3830                  objectiveChange = 1.0e100 ;
3831                  numberStrongInfeasible++;
3832                }
3833              }
3834            }
3835          } else if (iStatus==1) {
3836            objectiveChange = 1.0e100 ;
3837            numberStrongInfeasible++;
3838          } else {
3839            // Can't say much as we did not finish
3840            choice.finishedUp = false ;
3841            numberUnfinished++;
3842          }
3843          choice.upMovement = objectiveChange ;
3844         
3845          // restore bounds
3846          for ( j=0;j<numberColumns;j++) {
3847            if (saveLower[j] != lower[j])
3848              solver->setColLower(j,saveLower[j]);
3849            if (saveUpper[j] != upper[j])
3850              solver->setColUpper(j,saveUpper[j]);
3851          }
3852          if(needHotStartUpdate) {
3853            needHotStartUpdate = false;
3854            solver->resolve();
3855            //we may again have an integer feasible solution
3856            int numberIntegerInfeasibilities;
3857            int numberObjectInfeasibilities;
3858            if (model->feasibleSolution(
3859                                        numberIntegerInfeasibilities,
3860                                        numberObjectInfeasibilities)) {
3861              double objValue = solver->getObjValue();
3862              model->setBestSolution(CBC_STRONGSOL,
3863                                     objValue,
3864                                     solver->getColSolution()) ;
3865              solver->resolve();
3866              cutoff =model->getCutoff();
3867            }
3868            solver->markHotStart();
3869          }
3870         
3871          //printf("Up on %d, status is %d, obj %g its %d cost %g finished %d inf %d infobj %d\n",
3872          //     choice.objectNumber,iStatus,newObjectiveValue,choice.numItersUp,
3873          //     choice.upMovement,choice.finishedUp,choice.numIntInfeasUp,
3874          //     choice.numObjInfeasUp);
3875        }
3876   
3877        solver->setIntParam(OsiMaxNumIterationHotStart,saveLimit2); 
3878        /*
3879          End of evaluation for this candidate variable. Possibilities are:
3880          * Both sides below cutoff; this variable is a candidate for branching.
3881          * Both sides infeasible or above the objective cutoff: no further action
3882          here. Break from the evaluation loop and assume the node will be purged
3883          by the caller.
3884          * One side below cutoff: Install the branch (i.e., fix the variable). Break
3885          from the evaluation loop and assume the node will be reoptimised by the
3886          caller.
3887        */
3888        // reset
3889        choice.possibleBranch->resetNumberBranchesLeft();
3890        if (choice.upMovement<1.0e100) {
3891          if(choice.downMovement<1.0e100) {
3892            // In case solution coming in was odd
3893            choice.upMovement = CoinMax(0.0,choice.upMovement);
3894            choice.downMovement = CoinMax(0.0,choice.downMovement);
3895            if (couldChooseFirst)
3896              printf("candidate %d up %g down %g sort %g\n",iDo,choice.upMovement,choice.downMovement,sort[iDo]);
3897#if ZERO_ONE==2
3898            // branch on 0-1 first (temp)
3899            if (fabs(choice.possibleBranch->value())<1.0) {
3900              choice.upMovement *= ZERO_FAKE;
3901              choice.downMovement *= ZERO_FAKE;
3902            }
3903#endif
3904            // feasible - see which best
3905            if (!canSkip) {
3906              if (iColumn==-46) {
3907                printf("sort %g downest %g upest %g ",sort[iDo],downEstimate[iObject],
3908                     upEstimate[iObject]);
3909                printf("downMove %g upMove %g value %g current pseudo %g %g\n",
3910                       choice.downMovement,choice.upMovement,choice.possibleBranch->value(),
3911                       dynamicObject->downDynamicPseudoCost(),dynamicObject->upDynamicPseudoCost());
3912              }
3913              if (model->messageHandler()->logLevel()>3) 
3914                printf("sort %g downest %g upest %g ",sort[iDo],downEstimate[iObject],
3915                     upEstimate[iObject]);
3916              model->messageHandler()->message(CBC_STRONG,*model->messagesPointer())
3917                << iObject << iColumn
3918                <<choice.downMovement<<choice.numIntInfeasDown
3919                <<choice.upMovement<<choice.numIntInfeasUp
3920                <<choice.possibleBranch->value()
3921                <<CoinMessageEol;
3922            }
3923            //if (!stateOfSearch)
3924            //choice.numIntInfeasDown=99999; // temp fudge
3925            if (wantMiniTree)
3926              decision->setBestCriterion(-1.0);
3927            double bestCriterion = -1.0;
3928            //double gap = saveUpper[iColumn]-saveLower[iColumn];
3929            // Give precedence to ones with gap of 1.0
3930            //assert(gap>0.0);
3931            double factor = 1.0; //changeFactor/CoinMin(gap,100.0);
3932            int betterWay;
3933            {
3934              CbcBranchingObject * branchObj =
3935                dynamic_cast <CbcBranchingObject *>(branch_) ;
3936              if (branch_)
3937                assert (branchObj);
3938              betterWay = decision->betterBranch(choice.possibleBranch,
3939                                                     branchObj,
3940                                                     choice.upMovement*factor,
3941                                                     choice.numIntInfeasUp ,
3942                                                     choice.downMovement*factor,
3943                                                     choice.numIntInfeasDown );
3944            }
3945            if (wantMiniTree) {
3946              double criterion = decision->getBestCriterion();
3947              sort[numberMini]=-criterion;
3948              whichObject[numberMini++]=whichObject[iDo];
3949              assert (betterWay);
3950              if (criterion>bestCriterion) 
3951                bestCriterion=criterion;
3952              else
3953                betterWay=0;
3954            }
3955            if (iDo>=changeStrategy) {
3956              // make less likely
3957              changeStrategy+=numberStrong;
3958              changeFactor *= 0.9;
3959            }
3960            if (betterWay) {
3961              // C) create branching object
3962              if (choiceObject) {
3963                delete branch_;
3964                branch_ = choice.possibleBranch->clone();
3965              } else {
3966                delete branch_;
3967                branch_ = choice.possibleBranch;
3968                choice.possibleBranch=NULL;
3969              }
3970              {
3971                CbcBranchingObject * branchObj =
3972                  dynamic_cast <CbcBranchingObject *>(branch_) ;
3973                assert (branchObj);
3974                //branchObj->way(preferredWay);
3975                branchObj->way(betterWay);
3976              }
3977              if (couldChooseFirst)
3978                printf("choosing %d way %d\n",iDo,betterWay);
3979              bestChoice = choice.objectNumber;
3980              whichChoice = iDo;
3981              if (numberStrong<=1) {
3982                delete ws;
3983                ws=NULL;
3984                break;
3985              }
3986            } else {
3987              if (!choiceObject) {
3988                delete choice.possibleBranch;
3989                choice.possibleBranch=NULL;
3990              }
3991              if (iDo>=2*numberStrong) {
3992                delete ws;
3993                ws=NULL;
3994                break;
3995              }
3996              if (!dynamicObject||dynamicObject->numberTimesUp()>1) {
3997                if (iDo-whichChoice>=numberStrong) {
3998                  if (!choiceObject) {
3999                    delete choice.possibleBranch;
4000                    choice.possibleBranch=NULL;
4001                  }
4002                  break; // give up
4003                }
4004              } else {
4005                if (iDo-whichChoice>=2*numberStrong) {
4006                  delete ws;
4007                  ws=NULL;
4008                  if (!choiceObject) {
4009                    delete choice.possibleBranch;
4010                    choice.possibleBranch=NULL;
4011                  }
4012                  break; // give up
4013                }
4014              }
4015            }
4016          } else {
4017            // up feasible, down infeasible
4018            anyAction=-1;
4019            worstFeasible = CoinMax(worstFeasible,choice.upMovement);
4020            model->messageHandler()->message(CBC_STRONG,*model->messagesPointer())
4021              << iObject << iColumn
4022              <<choice.downMovement<<choice.numIntInfeasDown
4023              <<choice.upMovement<<choice.numIntInfeasUp
4024              <<choice.possibleBranch->value()
4025              <<CoinMessageEol;
4026            //printf("Down infeasible for choice %d sequence %d\n",i,
4027            // model->object(choice.objectNumber)->columnNumber());
4028            if (!solveAll) {
4029              choice.possibleBranch->way(1);
4030              choice.possibleBranch->branch();
4031              if (!choiceObject) {
4032                delete choice.possibleBranch;
4033                choice.possibleBranch=NULL;
4034              }
4035              delete ws;
4036              ws=NULL;
4037              break;
4038            } else {
4039              choice.fix=1;
4040              fixObject[numberToFix++]=choice;
4041#define FIXNOW
4042#ifdef FIXNOW
4043#if 0
4044              double value = ceil(saveSolution[iColumn]);
4045              saveLower[iColumn]=value;
4046              solver->setColLower(iColumn,value);
4047#else
4048              choice.possibleBranch->fix(solver,saveLower,saveUpper,1);
4049#endif
4050#endif
4051              if (!choiceObject) {
4052                choice.possibleBranch=NULL;
4053              } else {
4054                choiceObject = new CbcDynamicPseudoCostBranchingObject(*choiceObject);
4055                choice.possibleBranch=choiceObject;
4056              }
4057#ifdef FIXNOW
4058              assert(doneHotStart);
4059              solver->unmarkHotStart();
4060              solver->resolve();
4061              solver->markHotStart();
4062              // may be infeasible (if other way stopped on iterations)
4063              if (!solver->isProvenOptimal()) {
4064                // neither side feasible
4065                anyAction=-2;
4066                if (!choiceObject) {
4067                  delete choice.possibleBranch;
4068                  choice.possibleBranch=NULL;
4069                }
4070                //printf("Both infeasible for choice %d sequence %d\n",i,
4071                // model->object(choice.objectNumber)->columnNumber());
4072                delete ws;
4073                ws=NULL;
4074                break;
4075              }
4076#endif
4077            }
4078          }
4079        } else {
4080          if(choice.downMovement<1.0e100) {
4081            // down feasible, up infeasible
4082            anyAction=-1;
4083            worstFeasible = CoinMax(worstFeasible,choice.downMovement);
4084            model->messageHandler()->message(CBC_STRONG,*model->messagesPointer())
4085              << iObject << iColumn
4086              <<choice.downMovement<<choice.numIntInfeasDown
4087              <<choice.upMovement<<choice.numIntInfeasUp
4088              <<choice.possibleBranch->value()
4089              <<CoinMessageEol;
4090            //printf("Up infeasible for choice %d sequence %d\n",i,
4091            // model->object(choice.objectNumber)->columnNumber());
4092            if (!solveAll) {
4093              choice.possibleBranch->way(-1);
4094              choice.possibleBranch->branch();
4095              if (!choiceObject) {
4096                delete choice.possibleBranch;
4097                choice.possibleBranch=NULL;
4098              }
4099              delete ws;
4100              ws=NULL;
4101              break;
4102            } else {
4103              choice.fix=-1;
4104              fixObject[numberToFix++]=choice;
4105#ifdef FIXNOW
4106#if 0
4107              double value = floor(saveSolution[iColumn]);
4108              saveUpper[iColumn]=value;
4109              solver->setColUpper(iColumn,value);
4110#else
4111              choice.possibleBranch->fix(solver,saveLower,saveUpper,-1);
4112#endif
4113#endif
4114              if (!choiceObject) {
4115                choice.possibleBranch=NULL;
4116              } else {
4117                choiceObject = new CbcDynamicPseudoCostBranchingObject(*choiceObject);
4118                choice.possibleBranch=choiceObject;
4119              }
4120#ifdef FIXNOW
4121              assert(doneHotStart);
4122              solver->unmarkHotStart();
4123              solver->resolve();
4124              solver->markHotStart();
4125              // may be infeasible (if other way stopped on iterations)
4126              if (!solver->isProvenOptimal()) {
4127                // neither side feasible
4128                anyAction=-2;
4129                if (!choiceObject) {
4130                  delete choice.possibleBranch;
4131                  choice.possibleBranch=NULL;
4132                }
4133                //printf("Both infeasible for choice %d sequence %d\n",i,
4134                // model->object(choice.objectNumber)->columnNumber());
4135                delete ws;
4136                ws=NULL;
4137                break;
4138              }
4139#endif
4140            }
4141          } else {
4142            // neither side feasible
4143            anyAction=-2;
4144            if (!choiceObject) {
4145              delete choice.possibleBranch;
4146              choice.possibleBranch=NULL;
4147            }
4148            //printf("Both infeasible for choice %d sequence %d\n",i,
4149            // model->object(choice.objectNumber)->columnNumber());
4150            delete ws;
4151            ws=NULL;
4152            break;
4153          }
4154        }
4155        // Check max time
4156        hitMaxTime = ( CoinCpuTime()-model->getDblParam(CbcModel::CbcStartSeconds) > 
4157                       model->getDblParam(CbcModel::CbcMaximumSeconds));
4158        if (hitMaxTime) {
4159          // make sure rest are fast
4160          doQuickly=true;
4161          for ( int jDo=iDo+1;jDo<numberToDo;jDo++) {
4162            int iObject = whichObject[iDo];
4163            OsiObject * object = model->modifiableObject(iObject);
4164            CbcSimpleIntegerDynamicPseudoCost * dynamicObject =
4165              dynamic_cast <CbcSimpleIntegerDynamicPseudoCost *>(object) ;
4166            if (dynamicObject) 
4167              dynamicObject->setNumberBeforeTrust(0);
4168          }
4169          numberTest=0;
4170          distanceToCutoff=-COIN_DBL_MAX;
4171        }
4172        if (!choiceObject) {
4173          delete choice.possibleBranch;
4174        }
4175      }
4176      double averageChange = model->sumChangeObjective()/((double) model->getNodeCount());
4177      if (depth_<10||worstFeasible>0.2*averageChange) 
4178        solveAll=false;
4179      if (model->messageHandler()->logLevel()>3||false) { 
4180        if (anyAction==-2) {
4181          printf("infeasible\n");
4182        } else if(anyAction==-1) {
4183          if (!solveAll)
4184            printf("%d fixed\n",numberToFix);
4185          else
4186            printf("%d fixed AND choosing %d iDo %d iChosenWhen %d numberToDo %d\n",numberToFix,bestChoice,
4187                   iDo,whichChoice,numberToDo);
4188        } else {
4189          int iObject = whichObject[whichChoice];
4190          OsiObject * object = model->modifiableObject(iObject);
4191          CbcSimpleIntegerDynamicPseudoCost * dynamicObject =
4192            dynamic_cast <CbcSimpleIntegerDynamicPseudoCost *>(object) ;
4193          if (dynamicObject) {
4194            int iColumn = dynamicObject->columnNumber();
4195            printf("choosing %d (column %d) iChosenWhen %d numberToDo %d\n",bestChoice,
4196                   iColumn,whichChoice,numberToDo);
4197          }
4198        }
4199      }
4200      if (doneHotStart) {
4201        // Delete the snapshot
4202        solver->unmarkHotStart();
4203        // back to normal
4204        solver->setHintParam(OsiDoInBranchAndCut,true,OsiHintDo,NULL) ;
4205        // restore basis
4206        solver->setWarmStart(ws);
4207      }
4208      solver->setIntParam(OsiMaxNumIterationHotStart,saveLimit);
4209      // Unless infeasible we will carry on
4210      // But we could fix anyway
4211      if (numberToFix&&!hitMaxTime) {
4212        if (anyAction==-2) {
4213          // take off
4214          for (i = 0 ; i < numberToFix ; i++) {
4215            delete fixObject[i].possibleBranch;
4216          }
4217        } else {
4218          // apply and take off
4219          for (i = 0 ; i < numberToFix ; i++) {
4220#ifndef FIXNOW
4221            fixObject[i].possibleBranch->way(fixObject[i].fix) ;
4222            fixObject[i].possibleBranch->branch() ;
4223#endif
4224            delete fixObject[i].possibleBranch;
4225          }
4226          bool feasible=true;
4227#if ACTION <2
4228          if (solveAll) {
4229            // can do quick optimality check
4230            int easy=2;
4231            solver->setHintParam(OsiDoInBranchAndCut,true,OsiHintDo,&easy) ;
4232            solver->resolve() ;
4233            solver->setHintParam(OsiDoInBranchAndCut,true,OsiHintDo,NULL) ;
4234            feasible = solver->isProvenOptimal();
4235            if (feasible) {
4236              anyAction=0;
4237              numberMini=0;
4238              memcpy(saveSolution,solver->getColSolution(),numberColumns*sizeof(double));
4239              model->reserveCurrentSolution(saveSolution);
4240              memcpy(saveLower,solver->getColLower(),numberColumns*sizeof(double));
4241              memcpy(saveUpper,solver->getColUpper(),numberColumns*sizeof(double));
4242              model->setPointers(solver);
4243              // See if candidate still possible
4244              if (branch_) {
4245                const OsiObject * object = model->object(bestChoice);
4246                int preferredWay;
4247                double infeasibility = object->infeasibility(&usefulInfo,preferredWay);
4248                if (!infeasibility) {
4249                  // take out
4250                  delete branch_;
4251                  branch_=NULL;
4252                } else {
4253                  CbcBranchingObject * branchObj =
4254                    dynamic_cast <CbcBranchingObject *>(branch_) ;
4255                  assert (branchObj);
4256                  branchObj->way(preferredWay);
4257                }
4258              }
4259            } else {
4260              anyAction=-2;
4261              finished=true;
4262            }
4263          }
4264#endif
4265          // If  fixed then round again
4266          if (!branch_&&anyAction!=-2) {
4267            finished=false;
4268          }
4269          // If these in then different action
4270#if ACTION == 1
4271          if (!anyAction)
4272            anyAction=-1;
4273          finished=true;
4274#endif
4275        }
4276      }
4277      delete ws;
4278    }
4279  }
4280  if (model->messageHandler()->logLevel()>2) 
4281    printf("%d strong, %d iters, %d pen, %d mark, %d fixed, action %d nnott %d nt %d, %d dq %s ns %d\n",
4282         numberStrongDone,numberStrongIterations,xPen,xMark,
4283           numberToFix,anyAction,numberNotTrusted,px[0],px[1],px[2]>0 ? "y" : "n",px[3]);
4284  // update number of strong iterations etc
4285  model->incrementStrongInfo(numberStrongDone,numberStrongIterations,
4286                             anyAction==-2 ? 0:numberToFix,anyAction==-2);
4287  if (!newWay) {
4288  if (((model->searchStrategy()+1)%1000)==0) {
4289#ifndef COIN_DEVELOP
4290    if (solver->messageHandler()->logLevel()>1)
4291#endif
4292      printf("%d strong, %d iters, %d inf, %d not finished, %d not trusted\n",
4293             numberStrongDone,numberStrongIterations,numberStrongInfeasible,numberUnfinished,
4294             numberNotTrusted);
4295    // decide what to do
4296    int strategy=1;
4297    if (((numberUnfinished*4>numberStrongDone&&
4298         numberStrongInfeasible*40<numberStrongDone)||
4299         numberStrongInfeasible<0)&&model->numberStrong()<10&&model->numberBeforeTrust()<20&&model->numberObjects()>1000) {
4300      strategy=2;
4301#ifdef COIN_DEVELOP
4302      //if (model->logLevel()>1)
4303        printf("going to strategy 2\n");
4304#endif
4305#if NODE_NEW==1
4306      // Basically off
4307      model->setNumberStrong(0);
4308      model->setNumberBeforeTrust(-2);
4309#elif NODE_NEW==2
4310      // Weaken
4311      model->setNumberStrong(2);
4312      model->setNumberBeforeTrust(1);
4313#elif NODE_NEW==3
4314      // Off and ranging
4315      model->setNumberStrong(0);
4316      model->setNumberBeforeTrust(-1);
4317#elif NODE_NEW==4
4318      // Off and pseudo shadow prices
4319      model->setNumberStrong(0);
4320      model->setNumberBeforeTrust(-2);
4321#elif NODE_NEW==5
4322      // Just fewer iterations
4323#endif
4324      model->synchronizeNumberBeforeTrust();
4325    }
4326    if (numberNodes)
4327      strategy=1;  // should only happen after hot start
4328    if (model->searchStrategy()<999)
4329      model->setSearchStrategy(strategy);
4330  }
4331  if (model->searchStrategy()==1&&numberNodes>500&&numberNodes<-510) {
4332#ifndef COIN_DEVELOP
4333    if (solver->messageHandler()->logLevel()>1)
4334#endif
4335      printf("after %d nodes - %d strong, %d iters, %d inf, %d not finished, %d not trusted\n",
4336             numberNodes,numberStrongDone,numberStrongIterations,numberStrongInfeasible,numberUnfinished,
4337             numberNotTrusted);
4338    // decide what to do
4339    if (numberUnfinished*10>numberStrongDone+1||
4340        !numberStrongInfeasible) {
4341      //if (model->logLevel()>1)
4342        printf("going to strategy 2\n");
4343#if NODE_NEW==1
4344      // Basically off
4345      model->setNumberStrong(0);
4346      model->setNumberBeforeTrust(-2);
4347#elif NODE_NEW==2
4348      // Weaken
4349      model->setNumberStrong(2);
4350      model->setNumberBeforeTrust(1);
4351#elif NODE_NEW==3
4352      // Off and ranging
4353      model->setNumberStrong(0);
4354      model->setNumberBeforeTrust(-1);
4355#elif NODE_NEW==4
4356      // Off and pseudo shadow prices
4357      model->setNumberStrong(0);
4358      model->setNumberBeforeTrust(-2);
4359#elif NODE_NEW==5
4360      // Just fewer iterations
4361#endif
4362      model->synchronizeNumberBeforeTrust();
4363      model->setSearchStrategy(2);
4364    }
4365  }
4366  }
4367  //if (numberToFix&&depth_<5)
4368  //printf("%d fixed by strong at depth %d\n",numberToFix,depth_);
4369  // Set guessed solution value
4370  guessedObjectiveValue_ = objectiveValue_+estimatedDegradation;
4371 
4372  // Get collection of branches if mini tree wanted
4373  if (anyAction==0&&numberMini&&numberMini>1) {
4374    // Sort
4375    CoinSort_2(sort,sort+numberMini,whichObject);
4376    delete branch_;
4377    branch_=NULL;
4378    numberMini = CoinMin(numberMini,model->sizeMiniTree());
4379    anyAction=numberMini;
4380    branches = new OsiSolverBranch[numberMini];
4381    for (int iDo=0;iDo<numberMini;iDo++) {
4382      int iObject = whichObject[iDo];
4383      OsiObject * object = model->modifiableObject(iObject);
4384      CbcSimpleInteger * obj =
4385        dynamic_cast <CbcSimpleInteger *>(object) ;
4386      OsiSolverBranch * oneBranch;
4387      if (obj) {
4388        oneBranch = obj->solverBranch(solver,&usefulInfo);
4389      } else {
4390        CbcObject * obj =
4391          dynamic_cast <CbcObject *>(object) ;
4392        assert (obj);
4393        oneBranch = obj->solverBranch();
4394      }
4395      branches[iDo]=*oneBranch;
4396      delete oneBranch;
4397    }
4398  }
4399/*
4400  Cleanup, then we're finished
4401*/
4402  if (!model->branchingMethod())
4403    delete decision;
4404
4405  delete choiceObject;
4406  delete [] fixObject;
4407  delete [] sort;
4408  delete [] whichObject;
4409  delete [] objectMark;
4410  delete [] saveLower;
4411  delete [] saveUpper;
4412  delete [] upEstimate;
4413  delete [] downEstimate;
4414# ifdef COIN_HAS_CLP
4415  if (osiclp) 
4416    osiclp->setSpecialOptions(saveClpOptions);
4417# endif
4418  // restore solution
4419  solver->setColSolution(saveSolution);
4420  model->reserveCurrentSolution(saveSolution);
4421  delete [] saveSolution;
4422  model->setStateOfSearch(saveStateOfSearch);
4423  model->setLogLevel(saveLogLevel);
4424  // delete extra regions
4425  if (usefulInfo.usefulRegion_) {
4426    delete [] usefulInfo.usefulRegion_;
4427    delete [] usefulInfo.indexRegion_;
4428    delete [] usefulInfo.pi_;
4429    usefulInfo.usefulRegion_ = NULL;
4430    usefulInfo.indexRegion_ = NULL;
4431    usefulInfo.pi_ = NULL;
4432  }
4433  return anyAction;
4434}
4435int CbcNode::analyze (CbcModel *model, double * results)
4436{
4437  int i;
4438  int numberIterationsAllowed = model->numberAnalyzeIterations();
4439  OsiSolverInterface * solver = model->solver();
4440  objectiveValue_ = solver->getObjSense()*solver->getObjValue();
4441  double cutoff =model->getCutoff();
4442  const double * lower = solver->getColLower();
4443  const double * upper = solver->getColUpper();
4444  const double * dj = solver->getReducedCost();
4445  int numberObjects = model->numberObjects();
4446  int numberColumns = model->getNumCols();
4447  // Initialize arrays
4448  int numberIntegers = model->numberIntegers();
4449  int * back = new int[numberColumns];
4450  const int * integerVariable = model->integerVariable();
4451  for (i=0;i<numberColumns;i++) 
4452    back[i]=-1;
4453  // What results is
4454  double * newLower = results;
4455  double * objLower = newLower+numberIntegers;
4456  double * newUpper = objLower+numberIntegers;
4457  double * objUpper = newUpper+numberIntegers;
4458  for (i=0;i<numberIntegers;i++) {
4459    int iColumn = integerVariable[i];
4460    back[iColumn]=i;
4461    newLower[i]=0.0;
4462    objLower[i]=-COIN_DBL_MAX;
4463    newUpper[i]=0.0;
4464    objUpper[i]=-COIN_DBL_MAX;
4465  }
4466  double * saveUpper = new double[numberColumns];
4467  double * saveLower = new double[numberColumns];
4468  int anyAction=0;
4469  // Save solution in case heuristics need good solution later
4470 
4471  double * saveSolution = new double[numberColumns];
4472  memcpy(saveSolution,solver->getColSolution(),numberColumns*sizeof(double));
4473  model->reserveCurrentSolution(saveSolution);
4474  for (i=0;i<numberColumns;i++) {
4475    saveLower[i] = lower[i];
4476    saveUpper[i] = upper[i];
4477  }
4478  // Get arrays to sort
4479  double * sort = new double[numberObjects];
4480  int * whichObject = new int[numberObjects];
4481  int numberToFix=0;
4482  int numberToDo=0;
4483  double integerTolerance = 
4484    model->getDblParam(CbcModel::CbcIntegerTolerance);
4485  // point to useful information
4486  OsiBranchingInformation usefulInfo = model->usefulInformation();
4487  // and modify
4488  usefulInfo.depth_=depth_;
4489     
4490  // compute current state
4491  int numberObjectInfeasibilities; // just odd ones
4492  int numberIntegerInfeasibilities;
4493  model->feasibleSolution(
4494                          numberIntegerInfeasibilities,
4495                          numberObjectInfeasibilities);
4496# ifdef COIN_HAS_CLP
4497  OsiClpSolverInterface * osiclp = dynamic_cast< OsiClpSolverInterface*> (solver);
4498  int saveClpOptions=0;
4499  bool fastIterations = (model->specialOptions()&8)!=0;
4500  if (osiclp&&fastIterations) {
4501    // for faster hot start
4502    saveClpOptions = osiclp->specialOptions();
4503    osiclp->setSpecialOptions(saveClpOptions|8192);
4504  }
4505# else
4506  bool fastIterations = false ;
4507# endif
4508  /*
4509    Scan for branching objects that indicate infeasibility. Choose candidates
4510    using priority as the first criteria, then integer infeasibility.
4511   
4512    The algorithm is to fill the array with a set of good candidates (by
4513    infeasibility) with priority bestPriority.  Finding a candidate with
4514    priority better (less) than bestPriority flushes the choice array. (This
4515    serves as initialization when the first candidate is found.)
4516   
4517  */
4518  numberToDo=0;
4519  for (i=0;i<numberObjects;i++) {
4520    OsiObject * object = model->modifiableObject(i);
4521    CbcSimpleIntegerDynamicPseudoCost * dynamicObject =
4522      dynamic_cast <CbcSimpleIntegerDynamicPseudoCost *>(object) ;
4523    if(!dynamicObject)
4524      continue;
4525    int preferredWay;
4526    double infeasibility = object->infeasibility(&usefulInfo,preferredWay);
4527    int iColumn = dynamicObject->columnNumber();
4528    if (saveUpper[iColumn]==saveLower[iColumn])
4529      continue;
4530    if (infeasibility)
4531      sort[numberToDo]=-1.0e10-infeasibility;
4532    else
4533      sort[numberToDo]=-fabs(dj[iColumn]);
4534    whichObject[numberToDo++]=i;
4535  }
4536  // Save basis
4537  CoinWarmStart * ws = solver->getWarmStart();
4538  int saveLimit;
4539  solver->getIntParam(OsiMaxNumIterationHotStart,saveLimit);
4540  int targetIterations = CoinMax(500,numberIterationsAllowed/numberObjects);
4541  if (saveLimit<targetIterations)
4542    solver->setIntParam(OsiMaxNumIterationHotStart,targetIterations); 
4543  // Mark hot start
4544  solver->markHotStart();
4545  // Sort
4546  CoinSort_2(sort,sort+numberToDo,whichObject);
4547  //double distanceToCutoff=model->getCutoff()-objectiveValue_;
4548  double * currentSolution = model->currentSolution();
4549  double objMin = 1.0e50;
4550  double objMax = -1.0e50;
4551  bool needResolve=false;
4552  int iDo;
4553  for (iDo=0;iDo<numberToDo;iDo++) {
4554    CbcStrongInfo choice;
4555    int iObject = whichObject[iDo];
4556    OsiObject * object = model->modifiableObject(iObject);
4557    CbcSimpleIntegerDynamicPseudoCost * dynamicObject =
4558      dynamic_cast <CbcSimpleIntegerDynamicPseudoCost *>(object) ;
4559    int iColumn = dynamicObject->columnNumber();
4560    int preferredWay;
4561    object->infeasibility(&usefulInfo,preferredWay);
4562    double value = currentSolution[iColumn];
4563    double nearest = floor(value+0.5);
4564    double lowerValue = floor(value);
4565    bool satisfied=false;
4566    if (fabs(value-nearest)<=integerTolerance||value<saveLower[iColumn]||value>saveUpper[iColumn]) {
4567      satisfied=true;
4568      double newValue;
4569      if (nearest<saveUpper[iColumn]) {
4570        newValue = nearest + 1.0001*integerTolerance;
4571        lowerValue = nearest;
4572      } else {
4573        newValue = nearest - 1.0001*integerTolerance;
4574        lowerValue = nearest-1;
4575      }
4576      currentSolution[iColumn]=newValue;
4577    }
4578    double upperValue = lowerValue+1.0;
4579    CbcSimpleInteger * obj =
4580      dynamic_cast <CbcSimpleInteger *>(object) ;
4581    if (obj) {
4582      choice.possibleBranch=obj->createBranch(solver,&usefulInfo,preferredWay);
4583    } else {
4584      CbcObject * obj =
4585        dynamic_cast <CbcObject *>(object) ;
4586      assert (obj);
4587      choice.possibleBranch=obj->createBranch(preferredWay);
4588    }
4589    currentSolution[iColumn]=value;
4590    // Save which object it was
4591    choice.objectNumber=iObject;
4592    choice.numIntInfeasUp = numberUnsatisfied_;
4593    choice.numIntInfeasDown = numberUnsatisfied_;
4594    choice.downMovement = 0.0;
4595    choice.upMovement = 0.0;
4596    choice.numItersDown = 0;
4597    choice.numItersUp = 0;
4598    choice.fix=0; // say not fixed
4599    double objectiveChange ;
4600    double newObjectiveValue=1.0e100;
4601    int j;
4602    // status is 0 finished, 1 infeasible and other
4603    int iStatus;
4604    /*
4605      Try the down direction first. (Specify the initial branching alternative as
4606      down with a call to way(-1). Each subsequent call to branch() performs the
4607      specified branch and advances the branch object state to the next branch
4608      alternative.)
4609    */
4610    choice.possibleBranch->way(-1) ;
4611    choice.possibleBranch->branch() ;
4612    if (fabs(value-lowerValue)>integerTolerance) {
4613      solver->solveFromHotStart() ;
4614      /*
4615        We now have an estimate of objective degradation that we can use for strong
4616        branching. If we're over the cutoff, the variable is monotone up.
4617        If we actually made it to optimality, check for a solution, and if we have
4618        a good one, call setBestSolution to process it. Note that this may reduce the
4619        cutoff, so we check again to see if we can declare this variable monotone.
4620      */
4621      if (solver->isProvenOptimal())
4622        iStatus=0; // optimal
4623      else if (solver->isIterationLimitReached()
4624               &&!solver->isDualObjectiveLimitReached())
4625        iStatus=2; // unknown
4626      else
4627        iStatus=1; // infeasible
4628      newObjectiveValue = solver->getObjSense()*solver->getObjValue();
4629      choice.numItersDown = solver->getIterationCount();
4630      numberIterationsAllowed -= choice.numItersDown;
4631      objectiveChange = newObjectiveValue  - objectiveValue_;
4632      if (!iStatus) {
4633        choice.finishedDown = true ;
4634        if (newObjectiveValue>=cutoff) {
4635          objectiveChange = 1.0e100; // say infeasible
4636        } else {
4637          // See if integer solution
4638          if (model->feasibleSolution(choice.numIntInfeasDown,
4639                                      choice.numObjInfeasDown)
4640              &&model->problemFeasibility()->feasible(model,-1)>=0) {
4641            model->setBestSolution(CBC_STRONGSOL,
4642                                   newObjectiveValue,
4643                                   solver->getColSolution()) ;
4644            model->setLastHeuristic(NULL);
4645            model->incrementUsed(solver->getColSolution());
4646            cutoff =model->getCutoff();
4647            if (newObjectiveValue >= cutoff)    //  *new* cutoff
4648              objectiveChange = 1.0e100 ;
4649          }
4650        }
4651      } else if (iStatus==1) {
4652        objectiveChange = 1.0e100 ;
4653      } else {
4654        // Can't say much as we did not finish
4655        choice.finishedDown = false ;
4656      }
4657      choice.downMovement = objectiveChange ;
4658    }
4659    // restore bounds
4660    for ( j=0;j<numberColumns;j++) {
4661      if (saveLower[j] != lower[j])
4662        solver->setColLower(j,saveLower[j]);
4663      if (saveUpper[j] != upper[j])
4664        solver->setColUpper(j,saveUpper[j]);
4665    }
4666    // repeat the whole exercise, forcing the variable up
4667    choice.possibleBranch->branch();
4668    if (fabs(value-upperValue)>integerTolerance) {
4669      solver->solveFromHotStart() ;
4670      /*
4671        We now have an estimate of objective degradation that we can use for strong
4672        branching. If we're over the cutoff, the variable is monotone up.
4673        If we actually made it to optimality, check for a solution, and if we have
4674        a good one, call setBestSolution to process it. Note that this may reduce the
4675        cutoff, so we check again to see if we can declare this variable monotone.
4676      */
4677      if (solver->isProvenOptimal())
4678        iStatus=0; // optimal
4679      else if (solver->isIterationLimitReached()
4680               &&!solver->isDualObjectiveLimitReached())
4681        iStatus=2; // unknown
4682      else
4683        iStatus=1; // infeasible
4684      newObjectiveValue = solver->getObjSense()*solver->getObjValue();
4685      choice.numItersUp = solver->getIterationCount();
4686      numberIterationsAllowed -= choice.numItersUp;
4687      objectiveChange = newObjectiveValue  - objectiveValue_;
4688      if (!iStatus) {
4689        choice.finishedUp = true ;
4690        if (newObjectiveValue>=cutoff) {
4691          objectiveChange = 1.0e100; // say infeasible
4692        } else {
4693          // See if integer solution
4694          if (model->feasibleSolution(choice.numIntInfeasUp,
4695                                      choice.numObjInfeasUp)
4696              &&model->problemFeasibility()->feasible(model,-1)>=0) {
4697            model->setBestSolution(CBC_STRONGSOL,
4698                                   newObjectiveValue,
4699                                   solver->getColSolution()) ;
4700            model->setLastHeuristic(NULL);
4701            model->incrementUsed(solver->getColSolution());
4702            cutoff =model->getCutoff();
4703            if (newObjectiveValue >= cutoff)    //  *new* cutoff
4704              objectiveChange = 1.0e100 ;
4705          }
4706        }
4707      } else if (iStatus==1) {
4708        objectiveChange = 1.0e100 ;
4709      } else {
4710        // Can't say much as we did not finish
4711        choice.finishedUp = false ;
4712      }
4713      choice.upMovement = objectiveChange ;
4714     
4715      // restore bounds
4716      for ( j=0;j<numberColumns;j++) {
4717        if (saveLower[j] != lower[j])
4718          solver->setColLower(j,saveLower[j]);
4719        if (saveUpper[j] != upper[j])
4720          solver->setColUpper(j,saveUpper[j]);
4721      }
4722    }
4723    // If objective goes above certain amount we can set bound
4724    int jInt = back[iColumn];
4725    newLower[jInt]=upperValue;
4726    if (choice.finishedDown)
4727      objLower[jInt]=choice.downMovement+objectiveValue_;
4728    else
4729      objLower[jInt]=objectiveValue_;
4730    newUpper[jInt]=lowerValue;
4731    if (choice.finishedUp)
4732      objUpper[jInt]=choice.upMovement+objectiveValue_;
4733    else
4734      objUpper[jInt]=objectiveValue_;
4735    objMin = CoinMin(CoinMin(objLower[jInt],objUpper[jInt]),objMin);
4736    /*
4737      End of evaluation for this candidate variable. Possibilities are:
4738      * Both sides below cutoff; this variable is a candidate for branching.
4739      * Both sides infeasible or above the objective cutoff: no further action
4740      here. Break from the evaluation loop and assume the node will be purged
4741      by the caller.
4742      * One side below cutoff: Install the branch (i.e., fix the variable). Break
4743      from the evaluation loop and assume the node will be reoptimised by the
4744      caller.
4745    */
4746    if (choice.upMovement<1.0e100) {
4747      if(choice.downMovement<1.0e100) {
4748        objMax = CoinMax(CoinMax(objLower[jInt],objUpper[jInt]),objMax);
4749        // In case solution coming in was odd
4750        choice.upMovement = CoinMax(0.0,choice.upMovement);
4751        choice.downMovement = CoinMax(0.0,choice.downMovement);
4752        // feasible -
4753        model->messageHandler()->message(CBC_STRONG,*model->messagesPointer())
4754          << iObject << iColumn
4755          <<choice.downMovement<<choice.numIntInfeasDown
4756          <<choice.upMovement<<choice.numIntInfeasUp
4757          <<value
4758          <<CoinMessageEol;
4759      } else {
4760        // up feasible, down infeasible
4761        anyAction=-1;
4762        if (!satisfied)
4763          needResolve=true;
4764        choice.fix=1;
4765        numberToFix++;
4766        saveLower[iColumn]=upperValue;
4767        solver->setColLower(iColumn,upperValue);
4768      }
4769    } else {
4770      if(choice.downMovement<1.0e100) {
4771        // down feasible, up infeasible
4772        anyAction=-1;
4773        if (!satisfied)
4774          needResolve=true;
4775        choice.fix=-1;
4776        numberToFix++;
4777        saveUpper[iColumn]=lowerValue;
4778        solver->setColUpper(iColumn,lowerValue);
4779      } else {
4780        // neither side feasible
4781        anyAction=-2;
4782        printf("Both infeasible for choice %d sequence %d\n",i,
4783               model->object(choice.objectNumber)->columnNumber());
4784        delete ws;
4785        ws=NULL;
4786        //solver->writeMps("bad");
4787        numberToFix=-1;
4788        delete choice.possibleBranch;
4789        choice.possibleBranch=NULL;
4790        break;
4791      }
4792    }
4793    delete choice.possibleBranch;
4794    if (numberIterationsAllowed<=0)
4795      break;
4796    //printf("obj %d, col %d, down %g up %g value %g\n",iObject,iColumn,
4797    //     choice.downMovement,choice.upMovement,value);
4798  }
4799  printf("Best possible solution %g, can fix more if solution of %g found - looked at %d variables in %d iterations\n",
4800         objMin,objMax,iDo,model->numberAnalyzeIterations()-numberIterationsAllowed);
4801  model->setNumberAnalyzeIterations(numberIterationsAllowed);
4802  // Delete the snapshot
4803  solver->unmarkHotStart();
4804  // back to normal
4805  solver->setHintParam(OsiDoInBranchAndCut,true,OsiHintDo,NULL) ;
4806  solver->setIntParam(OsiMaxNumIterationHotStart,saveLimit);
4807  // restore basis
4808  solver->setWarmStart(ws);
4809  delete ws;
4810   
4811  delete [] sort;
4812  delete [] whichObject;
4813  delete [] saveLower;
4814  delete [] saveUpper;
4815  delete [] back;
4816  // restore solution
4817  solver->setColSolution(saveSolution);
4818# ifdef COIN_HAS_CLP
4819  if (osiclp) 
4820    osiclp->setSpecialOptions(saveClpOptions);
4821# endif
4822  model->reserveCurrentSolution(saveSolution);
4823  delete [] saveSolution;
4824  if (needResolve)
4825    solver->resolve();
4826  return numberToFix;
4827}
4828
4829
4830CbcNode::CbcNode(const CbcNode & rhs) 
4831{ 
4832#ifdef CHECK_NODE
4833  printf("CbcNode %x Constructor from rhs %x\n",this,&rhs);
4834#endif
4835  if (rhs.nodeInfo_)
4836    nodeInfo_ = rhs.nodeInfo_->clone();
4837  else
4838    nodeInfo_=NULL;
4839  objectiveValue_=rhs.objectiveValue_;
4840  guessedObjectiveValue_ = rhs.guessedObjectiveValue_;
4841  sumInfeasibilities_ = rhs.sumInfeasibilities_;
4842  if (rhs.branch_)
4843    branch_=rhs.branch_->clone();
4844  else
4845    branch_=NULL;
4846  depth_ = rhs.depth_;
4847  numberUnsatisfied_ = rhs.numberUnsatisfied_;
4848  nodeNumber_ = rhs.nodeNumber_;
4849  state_ = rhs.state_;
4850  if (nodeInfo_)
4851    assert ((state_&2)!=0);
4852  else
4853    assert ((state_&2)==0);
4854}
4855
4856CbcNode &
4857CbcNode::operator=(const CbcNode & rhs)
4858{
4859  if (this != &rhs) {
4860    delete nodeInfo_;
4861    if (rhs.nodeInfo_)
4862      nodeInfo_ = rhs.nodeInfo_->clone();
4863    else
4864      nodeInfo_ = NULL;
4865    objectiveValue_=rhs.objectiveValue_;
4866    guessedObjectiveValue_ = rhs.guessedObjectiveValue_;
4867    sumInfeasibilities_ = rhs.sumInfeasibilities_;
4868    if (rhs.branch_)
4869      branch_=rhs.branch_->clone();
4870    else
4871      branch_=NULL,
4872    depth_ = rhs.depth_;
4873    numberUnsatisfied_ = rhs.numberUnsatisfied_;
4874    nodeNumber_ = rhs.nodeNumber_;
4875    state_ = rhs.state_;
4876    if (nodeInfo_)
4877      assert ((state_&2)!=0);
4878    else
4879      assert ((state_&2)==0);
4880  }
4881  return *this;
4882}
4883CbcNode::~CbcNode ()
4884{
4885#ifdef CHECK_NODE
4886  if (nodeInfo_) {
4887    printf("CbcNode %x Destructor nodeInfo %x (%d)\n",
4888         this,nodeInfo_,nodeInfo_->numberPointingToThis());
4889    //assert(nodeInfo_->numberPointingToThis()>=0);
4890  } else {
4891    printf("CbcNode %x Destructor nodeInfo %x (?)\n",
4892         this,nodeInfo_);
4893  }
4894#endif
4895  if (nodeInfo_) {
4896    // was if (nodeInfo_&&(state_&2)!=0) {
4897    nodeInfo_->nullOwner();
4898    int numberToDelete=nodeInfo_->numberBranchesLeft();
4899    //    CbcNodeInfo * parent = nodeInfo_->parent();
4900    //assert (nodeInfo_->numberPointingToThis()>0);
4901    if (nodeInfo_->decrement(numberToDelete)==0||(state_&2)==0) {
4902      if ((state_&2)==0) 
4903        nodeInfo_->nullParent();
4904      delete nodeInfo_;
4905    } else {
4906      //printf("node %x nodeinfo %x parent %x\n",this,nodeInfo_,nodeInfo_->parent());
4907      // anyway decrement parent
4908      //if (parent)
4909      ///parent->decrement(1);
4910    }
4911  }
4912  delete branch_;
4913}
4914// Decrement  active cut counts
4915void 
4916CbcNode::decrementCuts(int change)
4917{
4918  if (nodeInfo_)
4919    assert ((state_&2)!=0);
4920  else
4921    assert ((state_&2)==0);
4922  if(nodeInfo_) {
4923    nodeInfo_->decrementCuts(change);
4924  }
4925}
4926void 
4927CbcNode::decrementParentCuts(CbcModel * model, int change)
4928{
4929  if (nodeInfo_)
4930    assert ((state_&2)!=0);
4931  else
4932    assert ((state_&2)==0);
4933  if(nodeInfo_) {
4934    nodeInfo_->decrementParentCuts(model, change);
4935  }
4936}
4937
4938/*
4939  Initialize reference counts (numberPointingToThis, numberBranchesLeft_)
4940  in the attached nodeInfo_.
4941*/
4942void
4943CbcNode::initializeInfo()
4944{
4945  assert(nodeInfo_ && branch_) ;
4946  nodeInfo_->initializeInfo(branch_->numberBranches());
4947  assert ((state_&2)!=0);
4948  assert (nodeInfo_->numberBranchesLeft()==
4949          branch_->numberBranchesLeft());
4950}
4951// Nulls out node info
4952void 
4953CbcNode::nullNodeInfo()
4954{
4955  nodeInfo_=NULL;
4956  // say not active
4957  state_ &= ~2;
4958}
4959
4960int
4961CbcNode::branch(OsiSolverInterface * solver)
4962{
4963  double changeInGuessed;
4964  assert (nodeInfo_->numberBranchesLeft()==
4965          branch_->numberBranchesLeft());
4966  if (!solver)
4967    changeInGuessed=branch_->branch();
4968  else
4969    changeInGuessed=branch_->branch(solver);
4970  guessedObjectiveValue_+= changeInGuessed;
4971  //#define PRINTIT
4972#ifdef PRINTIT
4973  int numberLeft = nodeInfo_->numberBranchesLeft();
4974  CbcNodeInfo * parent = nodeInfo_->parent();
4975  int parentNodeNumber = -1;
4976  //CbcBranchingObject * object1 = branch_->object_;
4977  //OsiObject * object = object1->
4978  //int sequence = object->columnNumber);
4979  int id=-1;
4980  double value=0.0;
4981  if (branch_) {
4982    id = branch_->variable();
4983    value = branch_->value();
4984  }
4985  printf("id %d value %g objvalue %g\n",id,value,objectiveValue_);
4986  if (parent)
4987    parentNodeNumber = parent->nodeNumber();
4988  printf("Node number %d, %s, way %d, depth %d, parent node number %d\n",
4989         nodeInfo_->nodeNumber(),(numberLeft==2) ? "leftBranch" : "rightBranch",
4990         way(),depth_,parentNodeNumber);
4991#endif
4992  return nodeInfo_->branchedOn();
4993}
4994/* Active arm of the attached OsiBranchingObject.
4995 
4996   In the simplest instance, coded -1 for the down arm of the branch, +1 for
4997   the up arm. But see OsiBranchingObject::way()
4998     Use nodeInfo--.numberBranchesLeft_ to see how active
4999*/
5000int 
5001CbcNode::way() const
5002{
5003  if (branch_) {
5004    CbcBranchingObject * obj =
5005      dynamic_cast <CbcBranchingObject *>(branch_) ;
5006    if (obj) {
5007      return obj->way();
5008    } else {
5009      OsiTwoWayBranchingObject * obj2 =
5010      dynamic_cast <OsiTwoWayBranchingObject *>(branch_) ;
5011      assert (obj2);
5012      return obj2->way();
5013    }
5014  } else {
5015    return 0;
5016  }
5017}
5018/* Create a branching object for the node
5019
5020    The routine scans the object list of the model and selects a set of
5021    unsatisfied objects as candidates for branching. The candidates are
5022    evaluated, and an appropriate branch object is installed.
5023
5024    The numberPassesLeft is decremented to stop fixing one variable each time
5025    and going on and on (e.g. for stock cutting, air crew scheduling)
5026
5027    If evaluation determines that an object is monotone or infeasible,
5028    the routine returns immediately. In the case of a monotone object,
5029    the branch object has already been called to modify the model.
5030
5031    Return value:
5032    <ul>
5033      <li>  0: A branching object has been installed
5034      <li> -1: A monotone object was discovered
5035      <li> -2: An infeasible object was discovered
5036    </ul>
5037    Branch state:
5038    <ul>
5039      <li> -1: start
5040      <li> -1: A monotone object was discovered
5041      <li> -2: An infeasible object was discovered
5042    </ul>
5043*/
5044int 
5045CbcNode::chooseOsiBranch (CbcModel * model,
5046                          CbcNode * lastNode,
5047                          OsiBranchingInformation * usefulInfo,
5048                          int branchState)
5049{
5050  int returnStatus=0;
5051  if (lastNode)
5052    depth_ = lastNode->depth_+1;
5053  else
5054    depth_ = 0;
5055  OsiSolverInterface * solver = model->solver();
5056  objectiveValue_ = solver->getObjValue()*solver->getObjSense();
5057  usefulInfo->objectiveValue_ = objectiveValue_;
5058  usefulInfo->depth_ = depth_;
5059  const double * saveInfoSol = usefulInfo->solution_;
5060  double * saveSolution = new double[solver->getNumCols()];
5061  memcpy(saveSolution,solver->getColSolution(),solver->getNumCols()*sizeof(double));
5062  usefulInfo->solution_ = saveSolution;
5063  OsiChooseVariable * choose = model->branchingMethod()->chooseMethod();
5064  int numberUnsatisfied=-1;
5065  if (branchState<0) {
5066    // initialize
5067    // initialize sum of "infeasibilities"
5068    sumInfeasibilities_ = 0.0;
5069    numberUnsatisfied = choose->setupList(usefulInfo,true);
5070    numberUnsatisfied_ = numberUnsatisfied;
5071    branchState=0;
5072    if (numberUnsatisfied_<0) {
5073      // infeasible
5074      delete [] saveSolution;
5075      return -2;
5076    }
5077  }
5078  // unset best
5079  int best=-1;
5080  choose->setBestObjectIndex(-1);
5081  if (numberUnsatisfied) {
5082    if (branchState>0||!choose->numberOnList()) {
5083      // we need to return at once - don't do strong branching or anything
5084      if (choose->numberOnList()||!choose->numberStrong()) {
5085        best = choose->candidates()[0];
5086        choose->setBestObjectIndex(best);
5087      } else {
5088        // nothing on list - need to try again - keep any solution
5089        numberUnsatisfied = choose->setupList(usefulInfo, false);
5090        numberUnsatisfied_ = numberUnsatisfied;
5091        if (numberUnsatisfied) {
5092          best = choose->candidates()[0];
5093          choose->setBestObjectIndex(best);
5094        }
5095      }
5096    } else {
5097      // carry on with strong branching or whatever
5098      int returnCode = choose->chooseVariable(solver, usefulInfo,true);
5099      // update number of strong iterations etc
5100      model->incrementStrongInfo(choose->numberStrongDone(),choose->numberStrongIterations(),
5101                                 returnCode==-1 ? 0:choose->numberStrongFixed(),returnCode==-1);
5102      if (returnCode>1) {
5103        // has fixed some
5104        returnStatus=-1;
5105      } else if (returnCode==-1) {
5106        // infeasible
5107        returnStatus=-2;
5108      } else if (returnCode==0) {
5109        // normal
5110        returnStatus=0;
5111        numberUnsatisfied=1;
5112      } else {
5113        // ones on list satisfied - double check
5114        numberUnsatisfied = choose->setupList(usefulInfo, false);
5115        numberUnsatisfied_ = numberUnsatisfied;
5116        if (numberUnsatisfied) {
5117          best = choose->candidates()[0];
5118          choose->setBestObjectIndex(best);
5119        }
5120      }
5121    }
5122  } 
5123  delete branch_;
5124  branch_ = NULL;
5125  guessedObjectiveValue_ = COIN_DBL_MAX;//objectiveValue_; // for now
5126  if (!returnStatus) {
5127    if (numberUnsatisfied) {
5128      // create branching object
5129      const OsiObject * obj = model->solver()->object(choose->bestObjectIndex());
5130      //const OsiSolverInterface * solver = usefulInfo->solver_;
5131      branch_ = obj->createBranch(model->solver(),usefulInfo,obj->whichWay());
5132    }
5133  }
5134  usefulInfo->solution_=saveInfoSol;
5135  delete [] saveSolution;
5136  // may have got solution
5137  if (choose->goodSolution()
5138      &&model->problemFeasibility()->feasible(model,-1)>=0) {
5139    // yes
5140    double objValue = choose->goodObjectiveValue();
5141    model->setBestSolution(CBC_STRONGSOL,
5142                                     objValue,
5143                                     choose->goodSolution()) ;
5144    model->setLastHeuristic(NULL);
5145    model->incrementUsed(choose->goodSolution());
5146    choose->clearGoodSolution();
5147  }
5148  return returnStatus;
5149}
5150int 
5151CbcNode::chooseClpBranch (CbcModel * model,
5152                       CbcNode * lastNode)
5153{ 
5154  assert(lastNode);
5155  depth_ = lastNode->depth_+1;
5156  delete branch_;
5157  branch_=NULL;
5158  OsiSolverInterface * solver = model->solver();
5159  //double saveObjectiveValue = solver->getObjValue();
5160  //double objectiveValue = CoinMax(solver->getObjSense()*saveObjectiveValue,objectiveValue_);
5161  const double * lower = solver->getColLower();
5162  const double * upper = solver->getColUpper();
5163  // point to useful information
5164  OsiBranchingInformation usefulInfo = model->usefulInformation();
5165  // and modify
5166  usefulInfo.depth_=depth_;
5167  int i;
5168  //bool beforeSolution = model->getSolutionCount()==0;
5169  int numberObjects = model->numberObjects();
5170  int numberColumns = model->getNumCols();
5171  double * saveUpper = new double[numberColumns];
5172  double * saveLower = new double[numberColumns];
5173
5174  // Save solution in case heuristics need good solution later
5175 
5176  double * saveSolution = new double[numberColumns];
5177  memcpy(saveSolution,solver->getColSolution(),numberColumns*sizeof(double));
5178  model->reserveCurrentSolution(saveSolution);
5179  for (i=0;i<numberColumns;i++) {
5180    saveLower[i] = lower[i];
5181    saveUpper[i] = upper[i];
5182  }
5183  // Save basis
5184  CoinWarmStart * ws = solver->getWarmStart();
5185  numberUnsatisfied_ = 0;
5186  // initialize sum of "infeasibilities"
5187  sumInfeasibilities_ = 0.0;
5188  // Note looks as if off end (hidden one)
5189  OsiObject * object = model->modifiableObject(numberObjects);
5190  CbcGeneralDepth * thisOne = dynamic_cast <CbcGeneralDepth *> (object);
5191  assert (thisOne);
5192  OsiClpSolverInterface * clpSolver
5193    = dynamic_cast<OsiClpSolverInterface *> (solver);
5194  assert (clpSolver);
5195  ClpSimplex * simplex = clpSolver->getModelPtr();
5196  int preferredWay;
5197  double infeasibility = object->infeasibility(&usefulInfo,preferredWay);
5198  if (thisOne->whichSolution()>=0) {
5199    ClpNode * nodeInfo = thisOne->nodeInfo(thisOne->whichSolution());
5200    nodeInfo->applyNode(simplex,2);
5201    simplex->dual();
5202    double cutoff = model->getCutoff();
5203    bool goodSolution=true;
5204    if (simplex->status()) {
5205      simplex->writeMps("bad7.mps",2);
5206      if (nodeInfo->objectiveValue()>cutoff-1.0e-2)
5207        goodSolution=false;
5208      else
5209        assert (!simplex->status());
5210    }
5211    if (goodSolution) {
5212      double newObjectiveValue = solver->getObjSense()*solver->getObjValue();
5213      // See if integer solution
5214      int numInf;
5215      int numInf2;
5216      bool gotSol = model->feasibleSolution(numInf,numInf2);
5217      if (!gotSol) {
5218        printf("numinf %d\n",numInf);
5219        double * sol = simplex->primalColumnSolution();
5220        for (int i=0;i<numberColumns;i++) {
5221          if (simplex->isInteger(i)) {
5222            double value = floor(sol[i]+0.5);
5223            if (fabs(value-sol[i])>1.0e-7) {
5224              printf("%d value %g\n",i,sol[i]);
5225              if (fabs(value-sol[i])<1.0e-3) {
5226                sol[i]=value;
5227              }
5228            }
5229          }
5230        }
5231        simplex->writeMps("bad8.mps",2);
5232        bool gotSol = model->feasibleSolution(numInf,numInf2);
5233        if (!gotSol) 
5234          assert (gotSol);
5235      }
5236      model->setBestSolution(CBC_STRONGSOL,
5237                             newObjectiveValue,
5238                             solver->getColSolution()) ;
5239      model->setLastHeuristic(NULL);
5240      model->incrementUsed(solver->getColSolution());
5241    }
5242  }
5243  // restore bounds
5244  { for (int j=0;j<numberColumns;j++) {
5245      if (saveLower[j] != lower[j])
5246        solver->setColLower(j,saveLower[j]);
5247      if (saveUpper[j] != upper[j])
5248        solver->setColUpper(j,saveUpper[j]);
5249    }
5250  }
5251  // restore basis
5252  solver->setWarmStart(ws);
5253  delete ws;
5254  int anyAction;
5255  //#define CHECK_PATH
5256#ifdef CHECK_PATH
5257extern int gotGoodNode_Z;
5258 if (gotGoodNode_Z>=0)
5259   printf("good node %d %g\n",gotGoodNode_Z,infeasibility);
5260#endif
5261  if (infeasibility>0.0) {
5262    if (infeasibility==COIN_DBL_MAX) {
5263      anyAction=-2; // infeasible
5264    } else {
5265      branch_=thisOne->createBranch(preferredWay);
5266      // Set to firat one (and change when re-pushing)
5267      CbcGeneralBranchingObject * branch = 
5268        dynamic_cast <CbcGeneralBranchingObject *> (branch_);
5269      branch->state(objectiveValue_,sumInfeasibilities_,
5270                    numberUnsatisfied_,0);
5271      branch->setNode(this);
5272      anyAction=0;
5273    }
5274  } else {
5275    anyAction=-1;
5276  }
5277#ifdef CHECK_PATH
5278 gotGoodNode_Z=-1;
5279#endif
5280  // Set guessed solution value
5281  guessedObjectiveValue_ = objectiveValue_+1.0e-5;
5282  delete [] saveLower;
5283  delete [] saveUpper;
5284 
5285  // restore solution
5286  solver->setColSolution(saveSolution);
5287  delete [] saveSolution;
5288  return anyAction;
5289}
5290/* Double checks in case node can change its mind!
5291   Returns objective value
5292   Can change objective etc */
5293double
5294CbcNode::checkIsCutoff(double cutoff)
5295{
5296  branch_->checkIsCutoff(cutoff);
5297  return objectiveValue_;
5298}
Note: See TracBrowser for help on using the repository browser.