source: branches/reengAnn/Cbc/src/CbcNode.cpp @ 1265

Last change on this file since 1265 was 1265, checked in by lou, 10 years ago

BUILDABLE. Fold in additional comments (in varying amount) for about 20 cbc
source files.

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