source: branches/devel/Cbc/src/CbcModel.cpp @ 606

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

for lou and callable cbcmain

  • Property svn:eol-style set to native
  • Property svn:keywords set to Author Date Id Revision
File size: 303.3 KB
Line 
1// Copyright (C) 2002, International Business Machines
2// Corporation and others.  All Rights Reserved.
3#if defined(_MSC_VER)
4// Turn off compiler warning about long names
5#  pragma warning(disable:4786)
6#endif
7
8#include "CbcConfig.h"
9
10#include <string>
11//#define CBC_DEBUG 1
12//#define CHECK_CUT_COUNTS
13//#define CHECK_NODE_FULL
14//#define NODE_LOG
15//#define GLOBAL_CUTS_JUST_POINTERS
16#include <cassert>
17#include <cmath>
18#include <cfloat>
19
20#ifdef COIN_HAS_CLP
21// include Presolve from Clp
22#include "ClpPresolve.hpp"
23#include "OsiClpSolverInterface.hpp"
24#endif
25
26#include "CbcEventHandler.hpp"
27
28#include "OsiSolverInterface.hpp"
29#include "OsiAuxInfo.hpp"
30#include "OsiSolverBranch.hpp"
31#include "OsiChooseVariable.hpp"
32#include "CoinWarmStartBasis.hpp"
33#include "CoinPackedMatrix.hpp"
34#include "CoinHelperFunctions.hpp"
35#include "CbcBranchActual.hpp"
36#include "CbcBranchDynamic.hpp"
37#include "CbcHeuristic.hpp"
38#include "CbcModel.hpp"
39#include "CbcTreeLocal.hpp"
40#include "CbcStatistics.hpp"
41#include "CbcStrategy.hpp"
42#include "CbcMessage.hpp"
43#include "OsiRowCut.hpp"
44#include "OsiColCut.hpp"
45#include "OsiRowCutDebugger.hpp"
46#include "OsiCuts.hpp"
47#include "CbcCountRowCut.hpp"
48#include "CbcCutGenerator.hpp"
49#include "CbcFeasibilityBase.hpp"
50// include Probing
51#include "CglProbing.hpp"
52// include preprocessing
53#include "CglPreProcess.hpp"
54#include "CglDuplicateRow.hpp"
55#include "CglStored.hpp"
56
57#include "CoinTime.hpp"
58#include "CoinMpsIO.hpp"
59
60#include "CbcCompareActual.hpp"
61#include "CbcTree.hpp"
62/* Various functions local to CbcModel.cpp */
63
64namespace {
65
66//-------------------------------------------------------------------
67// Returns the greatest common denominator of two
68// positive integers, a and b, found using Euclid's algorithm
69//-------------------------------------------------------------------
70static int gcd(int a, int b) 
71{
72  int remainder = -1;
73  // make sure a<=b (will always remain so)
74  if(a > b) {
75    // Swap a and b
76    int temp = a;
77    a = b;
78    b = temp;
79  }
80  // if zero then gcd is nonzero (zero may occur in rhs of packed)
81  if (!a) {
82    if (b) {
83      return b;
84    } else {
85      printf("**** gcd given two zeros!!\n");
86      abort();
87    }
88  }
89  while (remainder) {
90    remainder = b % a;
91    b = a;
92    a = remainder;
93  }
94  return b;
95}
96
97
98
99#ifdef CHECK_NODE_FULL
100
101/*
102  Routine to verify that tree linkage is correct. The invariant that is tested
103  is
104
105  reference count = (number of actual references) + (number of branches left)
106
107  The routine builds a set of paired arrays, info and count, by traversing the
108  tree. Each CbcNodeInfo is recorded in info, and the number of times it is
109  referenced (via the parent field) is recorded in count. Then a final check is
110  made to see if the numberPointingToThis_ field agrees.
111*/
112
113void verifyTreeNodes (const CbcTree * branchingTree, const CbcModel &model)
114
115{if (model.getNodeCount()==661) return;  printf("*** CHECKING tree after %d nodes\n",model.getNodeCount()) ;
116 
117  int j ;
118  int nNodes = branchingTree->size() ;
119# define MAXINFO 1000
120  int *count = new int [MAXINFO] ;
121  CbcNodeInfo **info = new CbcNodeInfo*[MAXINFO] ;
122  int nInfo = 0 ;
123/*
124  Collect all CbcNodeInfo objects in info, by starting from each live node and
125  traversing back to the root. Nodes in the live set should have unexplored
126  branches remaining.
127
128  TODO: The `while (nodeInfo)' loop could be made to break on reaching a
129        common ancester (nodeInfo is found in info[k]). Alternatively, the
130        check could change to signal an error if nodeInfo is not found above a
131        common ancestor.
132*/
133  for (j = 0 ; j < nNodes ; j++)
134  { CbcNode *node = branchingTree->nodePointer(j) ;
135  if (!node)
136    continue;
137    CbcNodeInfo *nodeInfo = node->nodeInfo() ; 
138    int change = node->nodeInfo()->numberBranchesLeft() ;
139    assert(change) ;
140    while (nodeInfo)
141    { int k ;
142      for (k = 0 ; k < nInfo ; k++)
143      { if (nodeInfo == info[k]) break ; }
144      if (k == nInfo)
145      { assert(nInfo < MAXINFO) ;
146        nInfo++ ;
147        info[k] = nodeInfo ;
148        count[k] = 0 ; }
149      nodeInfo = nodeInfo->parent() ; } }
150/*
151  Walk the info array. For each nodeInfo, look up its parent in info and
152  increment the corresponding count.
153*/
154  for (j = 0 ; j < nInfo ; j++)
155  { CbcNodeInfo *nodeInfo = info[j] ;
156    nodeInfo = nodeInfo->parent() ;
157    if (nodeInfo)
158    { int k ;
159      for (k = 0 ; k < nInfo ; k++)
160      { if (nodeInfo == info[k]) break ; }
161      assert (k < nInfo) ;
162      count[k]++ ; } }
163/*
164  Walk the info array one more time and check that the invariant holds. The
165  number of references (numberPointingToThis()) should equal the sum of the
166  number of actual references (held in count[]) plus the number of potential
167  references (unexplored branches, numberBranchesLeft()).
168*/
169  for (j = 0;j < nInfo;j++) {
170    CbcNodeInfo * nodeInfo = info[j] ;
171    if (nodeInfo) {
172      int k ;
173      for (k = 0;k < nInfo;k++)
174        if (nodeInfo == info[k])
175          break ;
176      printf("Nodeinfo %x - %d left, %d count\n",
177             nodeInfo,
178             nodeInfo->numberBranchesLeft(),
179             nodeInfo->numberPointingToThis()) ;
180      assert(nodeInfo->numberPointingToThis() ==
181             count[k]+nodeInfo->numberBranchesLeft()) ; } }
182
183  delete [] count ;
184  delete [] info ;
185 
186  return ; }
187
188#endif  /* CHECK_NODE_FULL */
189
190
191
192#ifdef CHECK_CUT_COUNTS
193
194/*
195  Routine to verify that cut reference counts are correct.
196*/
197void verifyCutCounts (const CbcTree * branchingTree, CbcModel &model)
198
199{ printf("*** CHECKING cuts after %d nodes\n",model.getNodeCount()) ;
200
201  int j ;
202  int nNodes = branchingTree->size() ;
203
204/*
205  cut.tempNumber_ exists for the purpose of doing this verification. Clear it
206  in all cuts. We traverse the tree by starting from each live node and working
207  back to the root. At each CbcNodeInfo, check for cuts.
208*/
209  for (j = 0 ; j < nNodes ; j++)
210  { CbcNode *node = branchingTree->nodePointer(j) ;
211    CbcNodeInfo * nodeInfo = node->nodeInfo() ;
212    assert (node->nodeInfo()->numberBranchesLeft()) ;
213    while (nodeInfo)
214    { int k ;
215      for (k = 0 ; k < nodeInfo->numberCuts() ; k++)
216      { CbcCountRowCut *cut = nodeInfo->cuts()[k] ;
217        if (cut) cut->tempNumber_ = 0; }
218      nodeInfo = nodeInfo->parent() ; } }
219/*
220  Walk the live set again, this time collecting the list of cuts in use at each
221  node. addCuts1 will collect the cuts in model.addedCuts_. Take into account
222  that when we recreate the basis for a node, we compress out the slack cuts.
223*/
224  for (j = 0 ; j < nNodes ; j++)
225  { CoinWarmStartBasis *debugws = model.getEmptyBasis() ;
226    CbcNode *node = branchingTree->nodePointer(j) ;
227    CbcNodeInfo *nodeInfo = node->nodeInfo(); 
228    int change = node->nodeInfo()->numberBranchesLeft() ;
229    printf("Node %d %x (info %x) var %d way %d obj %g",j,node,
230           node->nodeInfo(),node->columnNumber(),node->way(),
231           node->objectiveValue()) ;
232
233    model.addCuts1(node,debugws) ;
234
235    int i ;
236    int numberRowsAtContinuous = model.numberRowsAtContinuous() ;
237    CbcCountRowCut **addedCuts = model.addedCuts() ;
238    for (i = 0 ; i < model.currentNumberCuts() ; i++)
239    { CoinWarmStartBasis::Status status = 
240        debugws->getArtifStatus(i+numberRowsAtContinuous) ;
241      if (status != CoinWarmStartBasis::basic && addedCuts[i])
242      { addedCuts[i]->tempNumber_ += change ; } }
243
244    while (nodeInfo)
245    { nodeInfo = nodeInfo->parent() ;
246      if (nodeInfo) printf(" -> %x",nodeInfo); }
247    printf("\n") ;
248    delete debugws ; }
249/*
250  The moment of truth: We've tallied up the references by direct scan of the
251  search tree. Check for agreement with the count in the cut.
252
253  TODO: Rewrite to check and print mismatch only when tempNumber_ == 0?
254*/
255  for (j = 0 ; j < nNodes ; j++)
256  { CbcNode *node = branchingTree->nodePointer(j) ;
257    CbcNodeInfo *nodeInfo = node->nodeInfo(); 
258    while (nodeInfo)
259    { int k ;
260      for (k = 0 ; k < nodeInfo->numberCuts() ; k++)
261      { CbcCountRowCut *cut = nodeInfo->cuts()[k] ;
262        if (cut && cut->tempNumber_ >= 0)
263        { if (cut->tempNumber_ != cut->numberPointingToThis()) 
264            printf("mismatch %x %d %x %d %d\n",nodeInfo,k,
265                    cut,cut->tempNumber_,cut->numberPointingToThis()) ;
266          else
267            printf("   match %x %d %x %d %d\n", nodeInfo,k,
268                   cut,cut->tempNumber_,cut->numberPointingToThis()) ;
269          cut->tempNumber_ = -1 ; } }
270      nodeInfo = nodeInfo->parent() ; } }
271
272  return ; }
273
274#endif /* CHECK_CUT_COUNTS */
275
276
277//#define CHECK_CUT_SIZE
278#ifdef CHECK_CUT_SIZE
279
280/*
281  Routine to verify that cut reference counts are correct.
282*/
283void verifyCutSize (const CbcTree * branchingTree, CbcModel &model)
284{ 
285
286  int j ;
287  int nNodes = branchingTree->size() ;
288  int totalCuts=0;
289
290/*
291  cut.tempNumber_ exists for the purpose of doing this verification. Clear it
292  in all cuts. We traverse the tree by starting from each live node and working
293  back to the root. At each CbcNodeInfo, check for cuts.
294*/
295  for (j = 0 ; j < nNodes ; j++) {
296    CbcNode *node = branchingTree->nodePointer(j) ;
297    CbcNodeInfo * nodeInfo = node->nodeInfo() ;
298    assert (node->nodeInfo()->numberBranchesLeft()) ;
299    while (nodeInfo) {
300      totalCuts += nodeInfo->numberCuts();
301      nodeInfo = nodeInfo->parent() ;
302    }
303  }
304  printf("*** CHECKING cuts (size) after %d nodes - %d cuts\n",model.getNodeCount(),totalCuts) ;
305  return ;
306}
307
308#endif /* CHECK_CUT_SIZE */
309
310}
311
312 /* End unnamed namespace for CbcModel.cpp */
313
314
315static double trueIncrement=0.0;
316void 
317CbcModel::analyzeObjective ()
318/*
319  Try to find a minimum change in the objective function. The first scan
320  checks that there are no continuous variables with non-zero coefficients,
321  and grabs the largest objective coefficient associated with an unfixed
322  integer variable. The second scan attempts to scale up the objective
323  coefficients to a point where they are sufficiently close to integer that
324  we can pretend they are integer, and calculate a gcd over the coefficients
325  of interest. This will be the minimum increment for the scaled coefficients.
326  The final action is to scale the increment back for the original coefficients
327  and install it, if it's better than the existing value.
328
329  John's note: We could do better than this.
330
331  John's second note - apologies for changing s to z
332*/
333{ const double *objective = getObjCoefficients() ;
334  const double *lower = getColLower() ;
335  const double *upper = getColUpper() ;
336/*
337  Take a first scan to see if there are unfixed continuous variables in the
338  objective.  If so, the minimum objective change could be arbitrarily small.
339  Also pick off the maximum coefficient of an unfixed integer variable.
340
341  If the objective is found to contain only integer variables, set the
342  fathoming discipline to strict.
343*/
344  double maximumCost = 0.0 ;
345  trueIncrement=0.0;
346  bool possibleMultiple = true ;
347  int iColumn ;
348  int numberColumns = getNumCols() ;
349  for (iColumn = 0 ; iColumn < numberColumns ; iColumn++)
350  { if (upper[iColumn] > lower[iColumn]+1.0e-8)
351    { if (isInteger(iColumn)) 
352        maximumCost = CoinMax(maximumCost,fabs(objective[iColumn])) ;
353      else if (objective[iColumn]) 
354        possibleMultiple = false ; } }
355  setIntParam(CbcModel::CbcFathomDiscipline,possibleMultiple) ;
356/*
357  If a nontrivial increment is possible, try and figure it out. We're looking
358  for gcd(c<j>) for all c<j> that are coefficients of unfixed integer
359  variables. Since the c<j> might not be integers, try and inflate them
360  sufficiently that they look like integers (and we'll deflate the gcd
361  later).
362
363  2520.0 is used as it is a nice multiple of 2,3,5,7
364*/
365    if (possibleMultiple&&maximumCost)
366    { int increment = 0 ;
367      double multiplier = 2520.0 ;
368      while (10.0*multiplier*maximumCost < 1.0e8)
369        multiplier *= 10.0 ;
370    int bigIntegers = 0; // Count of large costs which are integer
371    for (iColumn = 0 ; iColumn < numberColumns ; iColumn++) {
372      if (upper[iColumn] > lower[iColumn]+1.0e-8) {
373        if (isInteger(iColumn)&&objective[iColumn]) {
374          double value = fabs(objective[iColumn])*multiplier ;
375          if (value <2.1e9) {
376            int nearest = (int) floor(value+0.5) ;
377            if (fabs(value-floor(value+0.5)) > 1.0e-8)
378              { increment = 0 ;
379              break ; }
380            else if (!increment)
381              { increment = nearest ; }
382            else
383              { increment = gcd(increment,nearest) ; }
384          } else {
385            // large value - may still be multiple of 1.0
386            value = fabs(objective[iColumn]);
387            if (fabs(value-floor(value+0.5)) > 1.0e-8) {
388              increment=0;
389              break;
390            } else {
391              bigIntegers++;
392            }
393          }
394        }
395      }
396    }
397
398/*
399  If the increment beats the current value for objective change, install it.
400*/
401      if (increment)
402      { double value = increment ;
403        double cutoff = getDblParam(CbcModel::CbcCutoffIncrement) ;
404        if (bigIntegers) {
405          // allow for 1.0
406          increment = gcd(increment,(int) multiplier);
407          value = increment;
408        }
409        value /= multiplier ;
410        trueIncrement=CoinMax(cutoff,value);;
411        if (value*0.999 > cutoff)
412        { messageHandler()->message(CBC_INTEGERINCREMENT,
413                                          messages())
414            << value << CoinMessageEol ;
415          setDblParam(CbcModel::CbcCutoffIncrement,value*0.999) ; } } }
416
417  return ; 
418}
419
420
421/**
422  \todo
423  Normally, it looks like we enter here from command dispatch in the main
424  routine, after calling the solver for an initial solution
425  (CbcModel::initialSolve, which simply calls the solver's initialSolve
426  routine.) The first thing we do is call resolve. Presumably there are
427  circumstances where this is nontrivial? There's also a call from
428  CbcModel::originalModel (tied up with integer presolve), which should be
429  checked.
430
431*/
432
433/*
434  The overall flow can be divided into three stages:
435    * Prep: Check that the lp relaxation remains feasible at the root. If so,
436      do all the setup for B&C.
437    * Process the root node: Generate cuts, apply heuristics, and in general do
438      the best we can to resolve the problem without B&C.
439    * Do B&C search until we hit a limit or exhaust the search tree.
440 
441  Keep in mind that in general there is no node in the search tree that
442  corresponds to the active subproblem. The active subproblem is represented
443  by the current state of the model,  of the solver, and of the constraint
444  system held by the solver.
445*/
446
447void CbcModel::branchAndBound(int doStatistics) 
448
449{
450/*
451  Capture a time stamp before we start.
452*/
453  dblParam_[CbcStartSeconds] = CoinCpuTime();
454  strongInfo_[0]=0;
455  strongInfo_[1]=0;
456  strongInfo_[2]=0;
457  numberStrongIterations_ = 0;
458  // original solver (only set if pre-processing)
459  OsiSolverInterface * originalSolver=NULL;
460  int numberOriginalObjects=numberObjects_;
461  OsiObject ** originalObject = NULL;
462  // Set up strategies
463  if (strategy_) {
464    // May do preprocessing
465    originalSolver = solver_;
466    strategy_->setupOther(*this);
467    if (strategy_->preProcessState()) {
468      // pre-processing done
469      if (strategy_->preProcessState()<0) {
470        // infeasible
471        handler_->message(CBC_INFEAS,messages_)<< CoinMessageEol ;
472        status_ = 0 ;
473        secondaryStatus_ = 1;
474        originalContinuousObjective_ = COIN_DBL_MAX;
475        return ; 
476      } else if (numberObjects_&&object_) {
477        numberOriginalObjects=numberObjects_;
478        // redo sequence
479        numberIntegers_=0;
480        int numberColumns = getNumCols();
481        int nOrig = originalSolver->getNumCols();
482        CglPreProcess * process = strategy_->process();
483        assert (process);
484        const int * originalColumns = process->originalColumns();
485        // allow for cliques etc
486        nOrig = CoinMax(nOrig,originalColumns[numberColumns-1]+1);
487        // try and redo debugger
488        OsiRowCutDebugger * debugger = const_cast<OsiRowCutDebugger *> (solver_->getRowCutDebuggerAlways());
489        if (debugger)
490          debugger->redoSolution(numberColumns,originalColumns);
491        originalObject = object_;
492        // object number or -1
493        int * temp = new int[nOrig];
494        int iColumn;
495        for (iColumn=0;iColumn<nOrig;iColumn++) 
496          temp[iColumn]=-1;
497        int iObject;
498        int nNonInt=0;
499        for (iObject=0;iObject<numberOriginalObjects;iObject++) {
500          iColumn = originalObject[iObject]->columnNumber();
501          if (iColumn<0) {
502            nNonInt++;
503          } else {
504            temp[iColumn]=iObject;
505          }
506        }
507        int numberNewIntegers=0;
508        int numberOldIntegers=0;
509        int numberOldOther=0;
510        for (iColumn=0;iColumn<numberColumns;iColumn++) {
511          int jColumn = originalColumns[iColumn];
512          if (temp[jColumn]>=0) {
513            int iObject= temp[jColumn];
514            CbcSimpleInteger * obj =
515              dynamic_cast <CbcSimpleInteger *>(originalObject[iObject]) ;
516            if (obj) 
517              numberOldIntegers++;
518            else
519              numberOldOther++;
520          } else if (isInteger(iColumn)) {
521            numberNewIntegers++;
522          }
523        }
524        /*
525          Allocate an array to hold the indices of the integer variables.
526          Make a large enough array for all objects
527        */
528        numberObjects_= numberNewIntegers+numberOldIntegers+numberOldOther+nNonInt;
529        object_ = new OsiObject * [numberObjects_];
530        integerVariable_ = new int [numberNewIntegers+numberOldIntegers];
531        /*
532          Walk the variables again, filling in the indices and creating objects for
533          the integer variables. Initially, the objects hold the index and upper &
534          lower bounds.
535        */
536        numberIntegers_=0;
537        for (iColumn=0;iColumn<numberColumns;iColumn++) {
538          int jColumn = originalColumns[iColumn];
539          if (temp[jColumn]>=0) {
540            int iObject= temp[jColumn];
541            CbcSimpleInteger * obj =
542              dynamic_cast <CbcSimpleInteger *>(originalObject[iObject]) ;
543            if (obj) {
544              object_[numberIntegers_] = originalObject[iObject]->clone();
545              // redo ids etc
546              object_[numberIntegers_]->resetSequenceEtc(numberColumns,originalColumns);
547              integerVariable_[numberIntegers_++]=iColumn;
548            }
549          } else if (isInteger(iColumn)) {
550            object_[numberIntegers_] =
551              new CbcSimpleInteger(this,iColumn);
552            integerVariable_[numberIntegers_++]=iColumn;
553          }
554        }
555        numberObjects_=numberIntegers_;
556        // Now append other column stuff
557        for (iColumn=0;iColumn<numberColumns;iColumn++) {
558          int jColumn = originalColumns[iColumn];
559          if (temp[jColumn]>=0) {
560            int iObject= temp[jColumn];
561            CbcSimpleInteger * obj =
562              dynamic_cast <CbcSimpleInteger *>(originalObject[iObject]) ;
563            if (!obj) {
564              object_[numberObjects_] = originalObject[iObject]->clone();
565              // redo ids etc
566              CbcObject * obj =
567              dynamic_cast <CbcObject *>(object_[numberObjects_]) ;
568              assert (obj);
569              obj->redoSequenceEtc(this,numberColumns,originalColumns);
570              numberObjects_++;
571            }
572          }
573        }
574        // now append non column stuff
575        for (iObject=0;iObject<numberOriginalObjects;iObject++) {
576          iColumn = originalObject[iObject]->columnNumber();
577          if (iColumn<0) {
578            object_[numberObjects_] = originalObject[iObject]->clone();
579            // redo ids etc
580            CbcObject * obj =
581              dynamic_cast <CbcObject *>(object_[numberObjects_]) ;
582            assert (obj);
583            obj->redoSequenceEtc(this,numberColumns,originalColumns);
584            numberObjects_++;
585          }
586        }
587        delete [] temp;
588        if (!numberObjects_)
589          handler_->message(CBC_NOINT,messages_) << CoinMessageEol ;
590      } else {
591        int numberColumns = getNumCols();
592        CglPreProcess * process = strategy_->process();
593        assert (process);
594        const int * originalColumns = process->originalColumns();
595        // try and redo debugger
596        OsiRowCutDebugger * debugger = const_cast<OsiRowCutDebugger *> (solver_->getRowCutDebuggerAlways());
597        if (debugger)
598          debugger->redoSolution(numberColumns,originalColumns);
599      }
600    } else {
601      //no preprocessing
602      originalSolver=NULL;
603    }
604    strategy_->setupCutGenerators(*this);
605    strategy_->setupHeuristics(*this);
606    // Set strategy print level to models
607    strategy_->setupPrinting(*this,handler_->logLevel());
608  }
609  eventHappened_=false;
610  CbcEventHandler *eventHandler = getEventHandler() ;
611  // set up for probing
612  probingInfo_ = new CglTreeProbingInfo(solver_);
613
614  // Try for dominated columns
615  if ((specialOptions_&64)!=0) {
616    CglDuplicateRow dupcuts(solver_);
617    dupcuts.setMode(2);
618    CglStored * storedCuts = dupcuts.outDuplicates(solver_);
619    addCutGenerator(storedCuts,1,"StoredCuts from dominated");
620  }
621  if (!nodeCompare_)
622    nodeCompare_=new CbcCompareDefault();;
623  // See if hot start wanted
624  CbcCompareBase * saveCompare = NULL;
625  if (hotstartSolution_) {
626    if (strategy_&&strategy_->preProcessState()>0) {
627      CglPreProcess * process = strategy_->process();
628      assert (process);
629      int n = solver_->getNumCols();
630      const int * originalColumns = process->originalColumns();
631      // columns should be in order ... but
632      double * tempS = new double[n];
633      for (int i=0;i<n;i++) {
634        int iColumn = originalColumns[i];
635        tempS[i]=hotstartSolution_[iColumn];
636      }
637      delete [] hotstartSolution_;
638      hotstartSolution_=tempS;
639      if (hotstartPriorities_) {
640        int * tempP = new int [n];
641        for (int i=0;i<n;i++) {
642          int iColumn = originalColumns[i];
643          tempP[i]=hotstartPriorities_[iColumn];
644        }
645        delete [] hotstartPriorities_;
646        hotstartPriorities_=tempP;
647      }
648    }
649    saveCompare = nodeCompare_;
650    // depth first
651    nodeCompare_ = new CbcCompareDepth();
652  }
653  if (!problemFeasibility_)
654    problemFeasibility_=new CbcFeasibilityBase();
655# ifdef CBC_DEBUG
656  std::string problemName ;
657  solver_->getStrParam(OsiProbName,problemName) ;
658  printf("Problem name - %s\n",problemName.c_str()) ;
659  solver_->setHintParam(OsiDoReducePrint,false,OsiHintDo,0) ;
660# endif
661/*
662  Assume we're done, and see if we're proven wrong.
663*/
664  status_ = 0 ;
665  secondaryStatus_ = 0;
666  phase_=0;
667/*
668  Scan the variables, noting the integer variables. Create an
669  CbcSimpleInteger object for each integer variable.
670*/
671  findIntegers(false) ;
672  // If dynamic pseudo costs then do
673  if (numberBeforeTrust_)
674    convertToDynamic();
675  // Set up char array to say if integer
676  delete [] integerInfo_;
677  {
678    int n = solver_->getNumCols();
679    integerInfo_ = new char [n];
680    for (int i=0;i<n;i++) {
681      if (solver_->isInteger(i))
682        integerInfo_[i]=1;
683      else
684        integerInfo_[i]=0;
685    }
686  }
687  if (preferredWay_) {
688    // set all unset ones
689    for (int iObject = 0 ; iObject < numberObjects_ ; iObject++) {
690      CbcObject * obj =
691        dynamic_cast <CbcObject *>(object_[iObject]) ;
692      if (obj&&!obj->preferredWay())
693        obj->setPreferredWay(preferredWay_);
694    }
695  } 
696/*
697  Ensure that objects on the lists of OsiObjects, heuristics, and cut
698  generators attached to this model all refer to this model.
699*/
700  synchronizeModel() ;
701  if (!solverCharacteristics_) {
702    OsiBabSolver * solverCharacteristics = dynamic_cast<OsiBabSolver *> (solver_->getAuxiliaryInfo());
703    if (solverCharacteristics) {
704      solverCharacteristics_ = solverCharacteristics;
705    } else {
706      // replace in solver
707      OsiBabSolver defaultC;
708      solver_->setAuxiliaryInfo(&defaultC);
709      solverCharacteristics_ = dynamic_cast<OsiBabSolver *> (solver_->getAuxiliaryInfo());
710    }
711  }
712  solverCharacteristics_->setSolver(solver_);
713  // Set so we can tell we are in initial phase in resolve
714  continuousObjective_ = -COIN_DBL_MAX ;
715/*
716  Solve the relaxation.
717
718  Apparently there are circumstances where this will be non-trivial --- i.e.,
719  we've done something since initialSolve that's trashed the solution to the
720  continuous relaxation.
721*/
722  bool feasible;
723  // If NLP then we assume already solved outside branchAndbound
724  if (!solverCharacteristics_->solverType()||solverCharacteristics_->solverType()==4) {
725    feasible=resolve(NULL,0) != 0 ;
726  } else {
727    // pick up given status
728    feasible = (solver_->isProvenOptimal() &&
729                !solver_->isDualObjectiveLimitReached()) ;
730  }
731  if (problemFeasibility_->feasible(this,0)<0) {
732    feasible=false; // pretend infeasible
733  }
734/*
735  If the linear relaxation of the root is infeasible, bail out now. Otherwise,
736  continue with processing the root node.
737*/
738  if (!feasible) {
739    status_ = 0 ;
740    if (!solver_->isProvenDualInfeasible()) {
741      handler_->message(CBC_INFEAS,messages_)<< CoinMessageEol ;
742      secondaryStatus_ = 1;
743    } else {
744      handler_->message(CBC_UNBOUNDED,messages_)<< CoinMessageEol ;
745      secondaryStatus_ = 7;
746    }
747    originalContinuousObjective_ = COIN_DBL_MAX;
748    solverCharacteristics_ = NULL;
749    return ;
750  }
751  // Convert to Osi if wanted
752  bool useOsiBranching=false;
753  //OsiBranchingInformation * persistentInfo = NULL;
754  if (branchingMethod_&&branchingMethod_->chooseMethod()) {
755    useOsiBranching=true;
756    //persistentInfo = new OsiBranchingInformation(solver_);
757    if (numberOriginalObjects) {
758      for (int iObject = 0 ; iObject < numberObjects_ ; iObject++) {
759        CbcObject * obj =
760          dynamic_cast <CbcObject *>(object_[iObject]) ;
761        if (obj) {
762          CbcSimpleInteger * obj2 =
763            dynamic_cast <CbcSimpleInteger *>(obj) ;
764          if (obj2) {
765            // back to Osi land
766            object_[iObject]=obj2->osiObject();
767            delete obj;
768          } else {
769            OsiSimpleInteger * obj3 =
770              dynamic_cast <OsiSimpleInteger *>(obj) ;
771            if (!obj3) {
772              OsiSOS * obj4 =
773                dynamic_cast <OsiSOS *>(obj) ;
774              if (!obj4) {
775                CbcSOS * obj5 =
776                  dynamic_cast <CbcSOS *>(obj) ;
777                if (obj5) {
778                  // back to Osi land
779                  object_[iObject]=obj5->osiObject(solver_);
780                } else {
781                  printf("Code up CbcObject type in Osi land\n");
782                  abort();
783                }
784              }
785            }
786          }
787        }
788      }
789      // and add to solver if none
790      if (!solver_->numberObjects()) {
791        solver_->addObjects(numberObjects_,object_);
792      } else {
793        if (solver_->numberObjects()!=numberOriginalObjects) {
794          printf("should have trapped that solver has objects before\n");
795          abort();
796        }
797      }
798    } else {
799      // do from solver
800      deleteObjects(false);
801      solver_->findIntegersAndSOS(false);
802      numberObjects_=solver_->numberObjects();
803      object_ = solver_->objects();
804      ownObjects_ = false;
805    }
806    branchingMethod_->chooseMethod()->setSolver(solver_);
807  }
808  // take off heuristics if have to
809  if (numberHeuristics_) {
810    int numberOdd=0;
811    for (int i=0;i<numberObjects_;i++) {
812      if (!object_[i]->canDoHeuristics()) 
813        numberOdd++;
814    }
815    if (numberOdd) {
816      int k=0;
817      for (int i=0;i<numberHeuristics_;i++) {
818        if (!heuristic_[i]->canDealWithOdd())
819          delete heuristic_[i];
820        else
821          heuristic_[k++]=heuristic_[i];
822      }
823      if (!k) {
824        delete [] heuristic_;
825        heuristic_=NULL;
826      }
827      numberHeuristics_=k;
828      handler_->message(CBC_HEURISTICS_OFF,messages_)<< numberOdd<<CoinMessageEol ;
829    }
830  }
831  // Save objective (just so user can access it)
832  originalContinuousObjective_ = solver_->getObjValue();
833  bestPossibleObjective_=originalContinuousObjective_;
834  sumChangeObjective1_=0.0;
835  sumChangeObjective2_=0.0;
836/*
837  OsiRowCutDebugger knows an optimal answer for a subset of MIP problems.
838  Assuming it recognises the problem, when called upon it will check a cut to
839  see if it cuts off the optimal answer.
840*/
841  // If debugger exists set specialOptions_ bit
842  if (solver_->getRowCutDebuggerAlways())
843    specialOptions_ |= 1;
844
845# ifdef CBC_DEBUG
846  if ((specialOptions_&1)==0)
847    solver_->activateRowCutDebugger(problemName.c_str()) ;
848  if (solver_->getRowCutDebuggerAlways())
849    specialOptions_ |= 1;
850# endif
851
852/*
853  Begin setup to process a feasible root node.
854*/
855  bestObjective_ = CoinMin(bestObjective_,1.0e50) ;
856  numberSolutions_ = 0 ;
857  stateOfSearch_=0;
858  numberHeuristicSolutions_ = 0 ;
859  // Everything is minimization
860  { 
861    // needed to sync cutoffs
862    double value ;
863    solver_->getDblParam(OsiDualObjectiveLimit,value) ;
864    dblParam_[CbcCurrentCutoff]= value * solver_->getObjSense();
865  }
866  double cutoff=getCutoff() ;
867  double direction = solver_->getObjSense() ;
868  dblParam_[CbcOptimizationDirection]=direction;
869  if (cutoff < 1.0e20&&direction<0.0)
870    messageHandler()->message(CBC_CUTOFF_WARNING1,
871                                    messages())
872                                      << cutoff << -cutoff << CoinMessageEol ;
873  if (cutoff > bestObjective_)
874    cutoff = bestObjective_ ;
875  setCutoff(cutoff) ;
876/*
877  We probably already have a current solution, but just in case ...
878*/
879  int numberColumns = getNumCols() ;
880  if (!currentSolution_)
881    currentSolution_ = new double[numberColumns] ;
882  testSolution_ = currentSolution_;
883  /* Tell solver we are in Branch and Cut
884     Could use last parameter for subtle differences */
885  solver_->setHintParam(OsiDoInBranchAndCut,true,OsiHintDo,NULL) ;
886#ifdef COIN_HAS_CLP
887 {
888   OsiClpSolverInterface * clpSolver
889     = dynamic_cast<OsiClpSolverInterface *> (solver_);
890   if (clpSolver) {
891#define CLP_QUICK_OPTIONS
892#ifdef CLP_QUICK_OPTIONS
893     // Try and re-use regions   
894     ClpSimplex * simplex = clpSolver->getModelPtr();
895     simplex->setPersistenceFlag(1);
896     clpSolver->deleteScaleFactors();
897     int value=131072;
898     clpSolver->setSpecialOptions(clpSolver->specialOptions()|value);
899     simplex->setSpecialOptions(simplex->specialOptions()|value);
900     //if (simplex->numberRows()<50)
901     //simplex->setAlphaAccuracy(1.0);
902     //clpSolver->setSpecialOptions((clpSolver->specialOptions()&~128)|65536);
903#endif
904     if ((specialOptions_&32)==0) {
905       ClpSimplex * clpSimplex = clpSolver->getModelPtr();
906       // take off names
907       clpSimplex->dropNames();
908     }
909     // no crunch if mostly continuous
910     int numberColumns = solver_->getNumCols();
911     if (numberColumns>1000&&numberIntegers_*4<numberColumns)
912       clpSolver->setSpecialOptions(clpSolver->specialOptions()&(~1));
913   }
914 }
915#endif
916/*
917  Create a copy of the solver, thus capturing the original (root node)
918  constraint system (aka the continuous system).
919*/
920  continuousSolver_ = solver_->clone() ;
921
922  numberRowsAtContinuous_ = getNumRows() ;
923/*
924  Check the objective to see if we can deduce a nontrivial increment. If
925  it's better than the current value for CbcCutoffIncrement, it'll be
926  installed.
927*/
928  if(solverCharacteristics_->reducedCostsAccurate())
929    analyzeObjective() ;
930/*
931  Set up for cut generation. addedCuts_ holds the cuts which are relevant for
932  the active subproblem. whichGenerator will be used to record the generator
933  that produced a given cut.
934*/
935  maximumWhich_ = 1000 ;
936  delete [] whichGenerator_;
937  whichGenerator_ = new int[maximumWhich_] ;
938  memset(whichGenerator_,0,maximumWhich_*sizeof(int));
939  int currentNumberCuts = 0 ;
940  maximumNumberCuts_ = 0 ;
941  currentNumberCuts_ = 0 ;
942  delete [] addedCuts_ ;
943  addedCuts_ = NULL ;
944/*
945  Set up an empty heap and associated data structures to hold the live set
946  (problems which require further exploration).
947*/
948  tree_->setComparison(*nodeCompare_) ;
949/*
950  Used to record the path from a node to the root of the search tree, so that
951  we can then traverse from the root to the node when restoring a subproblem.
952*/
953  maximumDepth_ = 10 ;
954  delete [] walkback_ ;
955  walkback_ = new CbcNodeInfo * [maximumDepth_] ;
956/*
957  Used to generate bound edits for CbcPartialNodeInfo.
958*/
959  double * lowerBefore = new double [numberColumns] ;
960  double * upperBefore = new double [numberColumns] ;
961/*
962 
963  Generate cuts at the root node and reoptimise. solveWithCuts does the heavy
964  lifting. It will iterate a generate/reoptimise loop (including reduced cost
965  fixing) until no cuts are generated, the change in objective falls off,  or
966  the limit on the number of rounds of cut generation is exceeded.
967
968  At the end of all this, any cuts will be recorded in cuts and also
969  installed in the solver's constraint system. We'll have reoptimised, and
970  removed any slack cuts (numberOldActiveCuts_ and numberNewCuts_ have been
971  adjusted accordingly).
972
973  Tell cut generators they can be a bit more aggressive at root node
974
975  TODO: Why don't we make a copy of the solution after solveWithCuts?
976  TODO: If numberUnsatisfied == 0, don't we have a solution?
977*/
978  phase_=1;
979  int iCutGenerator;
980  for (iCutGenerator = 0;iCutGenerator<numberCutGenerators_;iCutGenerator++) {
981    CglCutGenerator * generator = generator_[iCutGenerator]->generator();
982    generator->setAggressiveness(generator->getAggressiveness()+100);
983  }
984  OsiCuts cuts ;
985  int anyAction = -1 ;
986  numberOldActiveCuts_ = 0 ;
987  numberNewCuts_ = 0 ;
988  // Array to mark solution
989  delete [] usedInSolution_;
990  usedInSolution_ = new int[numberColumns];
991  CoinZeroN(usedInSolution_,numberColumns);
992/*
993  For printing totals and for CbcNode (numberNodes_)
994*/
995  numberIterations_ = 0 ;
996  numberNodes_ = 0 ;
997  numberNodes2_ = 0 ;
998  maximumStatistics_=0;
999  maximumDepthActual_=0;
1000  numberDJFixed_=0.0;
1001  statistics_ = NULL;
1002  // Do on switch
1003  if (doStatistics) {
1004    maximumStatistics_=10000;
1005    statistics_ = new CbcStatistics * [maximumStatistics_];
1006    memset(statistics_,0,maximumStatistics_*sizeof(CbcStatistics *));
1007  }
1008
1009  { int iObject ;
1010    int preferredWay ;
1011    int numberUnsatisfied = 0 ;
1012    memcpy(currentSolution_,solver_->getColSolution(),
1013           numberColumns*sizeof(double)) ;
1014    // point to useful information
1015    OsiBranchingInformation usefulInfo=usefulInformation();
1016
1017    for (iObject = 0 ; iObject < numberObjects_ ; iObject++)
1018    { double infeasibility =
1019          object_[iObject]->infeasibility(&usefulInfo,preferredWay) ;
1020      if (infeasibility ) numberUnsatisfied++ ; }
1021    // replace solverType
1022    if(solverCharacteristics_->tryCuts())  {
1023
1024      if (numberUnsatisfied)   {
1025        feasible = solveWithCuts(cuts,maximumCutPassesAtRoot_,
1026                                 NULL);
1027      } else if (solverCharacteristics_->solutionAddsCuts()||
1028                 solverCharacteristics_->alwaysTryCutsAtRootNode()) {
1029        // may generate cuts and turn the solution
1030        //to an infeasible one
1031        feasible = solveWithCuts(cuts, 1,
1032                                 NULL);
1033#if 0
1034        currentNumberCuts_ = cuts.sizeRowCuts();
1035        if (currentNumberCuts_ >= maximumNumberCuts_) {
1036          maximumNumberCuts_ = currentNumberCuts;
1037          delete [] addedCuts_;
1038          addedCuts_ = new CbcCountRowCut * [maximumNumberCuts_];
1039        }
1040#endif
1041      }
1042    }
1043    // check extra info on feasibility
1044    if (!solverCharacteristics_->mipFeasible())
1045      feasible = false;
1046  }
1047  // make cut generators less aggressive
1048  for (iCutGenerator = 0;iCutGenerator<numberCutGenerators_;iCutGenerator++) {
1049    CglCutGenerator * generator = generator_[iCutGenerator]->generator();
1050    generator->setAggressiveness(generator->getAggressiveness()-100);
1051  }
1052  currentNumberCuts_ = numberNewCuts_ ;
1053/*
1054  We've taken the continuous relaxation as far as we can. Time to branch.
1055  The first order of business is to actually create a node. chooseBranch
1056  currently uses strong branching to evaluate branch object candidates,
1057  unless forced back to simple branching. If chooseBranch concludes that a
1058  branching candidate is monotone (anyAction == -1) or infeasible (anyAction
1059  == -2) when forced to integer values, it returns here immediately.
1060
1061  Monotone variables trigger a call to resolve(). If the problem remains
1062  feasible, try again to choose a branching variable. At the end of the loop,
1063  resolved == true indicates that some variables were fixed.
1064
1065  Loss of feasibility will result in the deletion of newNode.
1066*/
1067
1068  bool resolved = false ;
1069  CbcNode *newNode = NULL ;
1070  numberFixedAtRoot_=0;
1071  numberFixedNow_=0;
1072  int numberIterationsAtContinuous = numberIterations_;
1073  //solverCharacteristics_->setSolver(solver_);
1074  if (feasible) {
1075    newNode = new CbcNode ;
1076    // Set objective value (not so obvious if NLP etc)
1077    setObjectiveValue(newNode,NULL);
1078    anyAction = -1 ;
1079    // To make depth available we may need a fake node
1080    CbcNode fakeNode;
1081    if (!currentNode_) {
1082      // Not true if sub trees assert (!numberNodes_);
1083      currentNode_=&fakeNode;
1084    }
1085    phase_=3;
1086    // only allow 1000 passes
1087    int numberPassesLeft=1000;
1088    // This is first crude step
1089    if (numberAnalyzeIterations_) {
1090      delete [] analyzeResults_;
1091      analyzeResults_ = new double [4*numberIntegers_];
1092      numberFixedAtRoot_=newNode->analyze(this,analyzeResults_);
1093      if (numberFixedAtRoot_>0) {
1094        printf("%d fixed by analysis\n",numberFixedAtRoot_);
1095        setPointers(solver_);
1096        numberFixedNow_ = numberFixedAtRoot_;
1097      } else if (numberFixedAtRoot_<0) {
1098        printf("analysis found to be infeasible\n");
1099        anyAction=-2;
1100        delete newNode ;
1101        newNode = NULL ;
1102        feasible = false ;
1103      }
1104    }
1105#if 0
1106    while (anyAction == -1)
1107    {
1108      // Set objective value (not so obvious if NLP etc)
1109      setObjectiveValue(newNode,NULL);
1110      if (numberBeforeTrust_==0 ) {
1111        anyAction = newNode->chooseBranch(this,NULL,numberPassesLeft) ;
1112      } else {
1113        OsiSolverBranch * branches=NULL;
1114        anyAction = newNode->chooseDynamicBranch(this,NULL,branches,numberPassesLeft) ;
1115        if (anyAction==-3)
1116          anyAction = newNode->chooseBranch(this,NULL,numberPassesLeft) ; // dynamic did nothing
1117      }
1118      numberPassesLeft--;
1119      if (anyAction == -1)
1120      {
1121        //if (solverCharacteristics_->solverType()!=4)
1122          feasible = resolve(NULL,10) != 0 ;
1123          //else
1124          //feasible = solveWithCuts(cuts,maximumCutPassesAtRoot_,NULL);
1125      if (problemFeasibility_->feasible(this,0)<0) {
1126        feasible=false; // pretend infeasible
1127      }
1128      reducedCostFix() ;
1129      resolved = true ;
1130#       ifdef CBC_DEBUG
1131      printf("Resolve (root) as something fixed, Obj value %g %d rows\n",
1132             solver_->getObjValue(),
1133             solver_->getNumRows()) ;
1134#       endif
1135      if (!feasible) anyAction = -2 ; }
1136      if (anyAction == -2||newNode->objectiveValue() >= cutoff)
1137      { delete newNode ;
1138        newNode = NULL ;
1139        feasible = false ; } }
1140#else
1141  OsiSolverBranch * branches = NULL;
1142  anyAction = chooseBranch(newNode, numberPassesLeft, NULL, cuts,resolved,
1143                           NULL,NULL,NULL,branches);
1144  if (anyAction == -2||newNode->objectiveValue() >= cutoff) {
1145    if (anyAction != -2) {
1146      // zap parent nodeInfo
1147#ifdef COIN_DEVELOP
1148      printf("zapping CbcNodeInfo %x\n",newNode->nodeInfo()->parent());
1149#endif
1150      if (newNode->nodeInfo())
1151        newNode->nodeInfo()->nullParent();
1152    }
1153    delete newNode ;
1154    newNode = NULL ;
1155    feasible = false ;
1156  }
1157#endif
1158  }
1159/*
1160  At this point, the root subproblem is infeasible or fathomed by bound
1161  (newNode == NULL), or we're live with an objective value that satisfies the
1162  current objective cutoff.
1163*/
1164  assert (!newNode || newNode->objectiveValue() <= cutoff) ;
1165  // Save address of root node as we don't want to delete it
1166  CbcNode * rootNode = newNode;
1167  // initialize for print out
1168  int lastDepth=0;
1169  int lastUnsatisfied=0;
1170  if (newNode)
1171    lastUnsatisfied=newNode->numberUnsatisfied();
1172/*
1173  The common case is that the lp relaxation is feasible but doesn't satisfy
1174  integrality (i.e., newNode->branchingObject(), indicating we've been able to
1175  select a branching variable). Remove any cuts that have gone slack due to
1176  forcing monotone variables. Then tack on an CbcFullNodeInfo object and full
1177  basis (via createInfo()) and stash the new cuts in the nodeInfo (via
1178  addCuts()). If, by some miracle, we have an integral solution at the root
1179  (newNode->branchingObject() is NULL), takeOffCuts() will ensure that the solver holds
1180  a valid solution for use by setBestSolution().
1181*/
1182  CoinWarmStartBasis *lastws = 0 ;
1183  if (feasible && newNode->branchingObject())
1184  { if (resolved)
1185    { takeOffCuts(cuts,false,NULL) ;
1186#     ifdef CHECK_CUT_COUNTS
1187      { printf("Number of rows after chooseBranch fix (root)"
1188               "(active only) %d\n",
1189                numberRowsAtContinuous_+numberNewCuts_+numberOldActiveCuts_) ;
1190        const CoinWarmStartBasis* debugws =
1191          dynamic_cast <const CoinWarmStartBasis*>(solver_->getWarmStart()) ;
1192        debugws->print() ;
1193        delete debugws ; }
1194#     endif
1195    }
1196  //newNode->createInfo(this,NULL,NULL,NULL,NULL,0,0) ;
1197    newNode->nodeInfo()->addCuts(cuts,
1198                                 newNode->numberBranches(),whichGenerator_) ;
1199    if (lastws) delete lastws ;
1200    lastws = dynamic_cast<CoinWarmStartBasis*>(solver_->getWarmStart()) ;
1201  }
1202/*
1203  Continuous data to be used later
1204*/
1205  continuousObjective_ = solver_->getObjValue()*solver_->getObjSense();
1206  continuousInfeasibilities_ = 0 ;
1207  if (newNode)
1208  { continuousObjective_ = newNode->objectiveValue() ;
1209    delete [] continuousSolution_;
1210    continuousSolution_ = CoinCopyOfArray(solver_->getColSolution(),
1211                                             numberColumns);
1212    continuousInfeasibilities_ = newNode->numberUnsatisfied() ; }
1213/*
1214  Bound may have changed so reset in objects
1215*/
1216  { int i ;
1217    for (i = 0;i < numberObjects_;i++)
1218      object_[i]->resetBounds(solver_) ; }
1219  stoppedOnGap_ = false ;
1220/*
1221  Feasible? Then we should have either a live node prepped for future
1222  expansion (indicated by variable() >= 0), or (miracle of miracles) an
1223  integral solution at the root node.
1224
1225  initializeInfo sets the reference counts in the nodeInfo object.  Since
1226  this node is still live, push it onto the heap that holds the live set.
1227*/
1228  double bestValue = 0.0 ;
1229  if (newNode) {
1230    bestValue = newNode->objectiveValue();
1231    if (newNode->branchingObject()) {
1232      newNode->initializeInfo() ;
1233      tree_->push(newNode) ;
1234      if (statistics_) {
1235        if (numberNodes2_==maximumStatistics_) {
1236          maximumStatistics_ = 2*maximumStatistics_;
1237          CbcStatistics ** temp = new CbcStatistics * [maximumStatistics_];
1238          memset(temp,0,maximumStatistics_*sizeof(CbcStatistics *));
1239          memcpy(temp,statistics_,numberNodes2_*sizeof(CbcStatistics *));
1240          delete [] statistics_;
1241          statistics_=temp;
1242        }
1243        assert (!statistics_[numberNodes2_]);
1244        statistics_[numberNodes2_]=new CbcStatistics(newNode);
1245      }
1246      numberNodes2_++;
1247#     ifdef CHECK_NODE
1248      printf("Node %x on tree\n",newNode) ;
1249#     endif
1250    } else {
1251      // continuous is integer
1252      double objectiveValue = newNode->objectiveValue();
1253      setBestSolution(CBC_SOLUTION,objectiveValue,
1254                      solver_->getColSolution()) ;
1255      delete newNode ;
1256      newNode = NULL ;
1257    }
1258  }
1259
1260  if (printFrequency_ <= 0) {
1261    printFrequency_ = 1000 ;
1262    if (getNumCols() > 2000)
1263      printFrequency_ = 100 ;
1264  }
1265  /*
1266    It is possible that strong branching fixes one variable and then the code goes round
1267    again and again.  This can take too long.  So we need to warn user - just once.
1268  */
1269  numberLongStrong_=0;
1270/*
1271  At last, the actual branch-and-cut search loop, which will iterate until
1272  the live set is empty or we hit some limit (integrality gap, time, node
1273  count, etc.). The overall flow is to rebuild a subproblem, reoptimise using
1274  solveWithCuts(), choose a branching pattern with chooseBranch(), and finally
1275  add the node to the live set.
1276
1277  The first action is to winnow the live set to remove nodes which are worse
1278  than the current objective cutoff.
1279*/
1280  while (!tree_->empty()) {
1281#ifdef BONMIN
1282    assert(!solverCharacteristics_->solutionAddsCuts() || solverCharacteristics_->mipFeasible());
1283#endif
1284    if (cutoff > getCutoff()) {
1285      double newCutoff = getCutoff();
1286      if (analyzeResults_) {
1287        // see if we could fix any (more)
1288        int n=0;
1289        double * newLower = analyzeResults_;
1290        double * objLower = newLower+numberIntegers_;
1291        double * newUpper = objLower+numberIntegers_;
1292        double * objUpper = newUpper+numberIntegers_;
1293        for (int i=0;i<numberIntegers_;i++) {
1294          if (objLower[i]>newCutoff) {
1295            n++;
1296            if (objUpper[i]>newCutoff) {
1297              newCutoff = -COIN_DBL_MAX;
1298              break;
1299            }
1300          } else if (objUpper[i]>newCutoff) {
1301            n++;
1302          }
1303        }
1304        if (newCutoff==-COIN_DBL_MAX) {
1305          printf("Root analysis says finished\n");
1306        } else if (n>numberFixedNow_) {
1307          printf("%d more fixed by analysis - now %d\n",n-numberFixedNow_,n);
1308          numberFixedNow_=n;
1309        }
1310      }
1311      if (eventHandler) {
1312        if (!eventHandler->event(CbcEventHandler::solution)) {
1313          eventHappened_=true; // exit
1314        }
1315      }
1316      // Do from deepest
1317      tree_->cleanTree(this, newCutoff,bestPossibleObjective_) ;
1318      nodeCompare_->newSolution(this) ;
1319      nodeCompare_->newSolution(this,continuousObjective_,
1320                                continuousInfeasibilities_) ;
1321      tree_->setComparison(*nodeCompare_) ;
1322      if (tree_->empty())
1323        break; // finished
1324    }
1325    cutoff = getCutoff() ;
1326/*
1327  Periodic activities: Opportunities to
1328    + tweak the nodeCompare criteria,
1329    + check if we've closed the integrality gap enough to quit,
1330    + print a summary line to let the user know we're working
1331*/
1332    if ((numberNodes_%1000) == 0) {
1333      bool redoTree=nodeCompare_->every1000Nodes(this, numberNodes_) ;
1334#ifdef CHECK_CUT_SIZE
1335      verifyCutSize (tree_, *this);
1336#endif
1337      // redo tree if wanted
1338      if (redoTree)
1339        tree_->setComparison(*nodeCompare_) ;
1340    }
1341    if (saveCompare&&!hotstartSolution_) {
1342      // hotstart switched off
1343      delete nodeCompare_; // off depth first
1344      nodeCompare_=saveCompare;
1345      saveCompare=NULL;
1346      // redo tree
1347      tree_->setComparison(*nodeCompare_) ;
1348    }
1349    if ((numberNodes_%printFrequency_) == 0) {
1350      int j ;
1351      int nNodes = tree_->size() ;
1352      bestPossibleObjective_ = 1.0e100 ;
1353      for (j = 0;j < nNodes;j++) {
1354        CbcNode * node = tree_->nodePointer(j) ;
1355        if (node&&node->objectiveValue() < bestPossibleObjective_)
1356          bestPossibleObjective_ = node->objectiveValue() ;
1357      }
1358      if (!intParam_[CbcPrinting]) {
1359        messageHandler()->message(CBC_STATUS,messages())
1360          << numberNodes_<< nNodes<< bestObjective_<< bestPossibleObjective_
1361          <<getCurrentSeconds()
1362          << CoinMessageEol ;
1363      } else {
1364        messageHandler()->message(CBC_STATUS2,messages())
1365          << numberNodes_<< nNodes<< bestObjective_<< bestPossibleObjective_
1366          <<lastDepth<<lastUnsatisfied<<numberIterations_
1367          <<getCurrentSeconds()
1368          << CoinMessageEol ;
1369      }
1370      if (!eventHandler->event(CbcEventHandler::treeStatus)) {
1371        eventHappened_=true; // exit
1372      }
1373    }
1374    // If no solution but many nodes - signal change in strategy
1375    if (numberNodes_>2*numberObjects_+1000&&stateOfSearch_!=2)
1376      stateOfSearch_=3;
1377    // See if can stop on gap
1378    double testGap = CoinMax(dblParam_[CbcAllowableGap],
1379                             CoinMax(fabs(bestObjective_),fabs(bestPossibleObjective_))
1380                             *dblParam_[CbcAllowableFractionGap]);
1381    if (bestObjective_-bestPossibleObjective_ < testGap && getCutoffIncrement()>=0.0) {
1382      stoppedOnGap_ = true ;
1383    }
1384
1385#   ifdef CHECK_NODE_FULL
1386    verifyTreeNodes(tree_,*this) ;
1387#   endif
1388#   ifdef CHECK_CUT_COUNTS
1389    verifyCutCounts(tree_,*this) ;
1390#   endif
1391
1392/*
1393  Now we come to the meat of the loop. To create the active subproblem, we'll
1394  pop the most promising node in the live set, rebuild the subproblem it
1395  represents, and then execute the current arm of the branch to create the
1396  active subproblem.
1397*/
1398    CbcNode *node = tree_->bestNode(cutoff) ;
1399    // Possible one on tree worse than cutoff
1400    if (!node)
1401      continue;
1402    currentNode_=node; // so can be accessed elsewhere
1403#ifdef CBC_DEBUG
1404    printf("%d unsat, way %d, obj %g est %g\n",
1405           node->numberUnsatisfied(),node->way(),node->objectiveValue(),
1406           node->guessedObjectiveValue());
1407#endif
1408    // Save clone in branching decision
1409    if(branchingMethod_)
1410      branchingMethod_->saveBranchingObject(node->modifiableBranchingObject());
1411    bool nodeOnTree=false; // Node has been popped
1412    // Say not on optimal path
1413    bool onOptimalPath=false;
1414#   ifdef CHECK_NODE
1415    printf("Node %x popped from tree - %d left, %d count\n",node,
1416           node->nodeInfo()->numberBranchesLeft(),
1417           node->nodeInfo()->numberPointingToThis()) ;
1418    printf("\tdepth = %d, z =  %g, unsat = %d, var = %d.\n",
1419           node->depth(),node->objectiveValue(),
1420           node->numberUnsatisfied(),
1421           node->columnNumber()) ;
1422#   endif
1423    lastDepth=node->depth();
1424    lastUnsatisfied=node->numberUnsatisfied();
1425
1426/*
1427  Rebuild the subproblem for this node:  Call addCuts() to adjust the model
1428  to recreate the subproblem for this node (set proper variable bounds, add
1429  cuts, create a basis).  This may result in the problem being fathomed by
1430  bound or infeasibility. Returns 1 if node is fathomed.
1431  Execute the current arm of the branch: If the problem survives, save the
1432  resulting variable bounds and call branch() to modify variable bounds
1433  according to the current arm of the branching object. If we're processing
1434  the final arm of the branching object, flag the node for removal from the
1435  live set.
1436*/
1437    CbcNodeInfo * nodeInfo = node->nodeInfo() ;
1438    newNode = NULL ;
1439    if (!addCuts(node,lastws,numberFixedNow_>numberFixedAtRoot_))
1440    { int i ;
1441      const double * lower = getColLower() ;
1442      const double * upper = getColUpper() ;
1443      for (i = 0 ; i < numberColumns ; i++)
1444      { lowerBefore[i]= lower[i] ;
1445        upperBefore[i]= upper[i] ; }
1446      if ((solverCharacteristics_->extraCharacteristics()&2)!=0) {
1447        solverCharacteristics_->setBeforeLower(lowerBefore);
1448        solverCharacteristics_->setBeforeUpper(upperBefore);
1449      }
1450      bool deleteNode ;
1451      if (messageHandler()->logLevel()>2)
1452        node->modifiableBranchingObject()->print();
1453      int retCode;
1454      if (!useOsiBranching) 
1455        retCode = node->branch(NULL); // old way
1456      else
1457        retCode = node->branch(solver_); // new way
1458      if (retCode)
1459      { 
1460        // set nodenumber correctly
1461        node->nodeInfo()->setNodeNumber(numberNodes2_);
1462        tree_->push(node) ;
1463        if (statistics_) {
1464          if (numberNodes2_==maximumStatistics_) {
1465            maximumStatistics_ = 2*maximumStatistics_;
1466            CbcStatistics ** temp = new CbcStatistics * [maximumStatistics_];
1467            memset(temp,0,maximumStatistics_*sizeof(CbcStatistics *));
1468            memcpy(temp,statistics_,numberNodes2_*sizeof(CbcStatistics *));
1469            delete [] statistics_;
1470            statistics_=temp;
1471          }
1472          assert (!statistics_[numberNodes2_]);
1473          statistics_[numberNodes2_]=new CbcStatistics(node);
1474        }
1475        numberNodes2_++;
1476        nodeOnTree=true; // back on tree
1477        deleteNode = false ;
1478#       ifdef CHECK_NODE
1479        printf("Node %x pushed back on tree - %d left, %d count\n",node,
1480               nodeInfo->numberBranchesLeft(),
1481               nodeInfo->numberPointingToThis()) ;
1482#       endif
1483      }
1484      else
1485      { deleteNode = true ;
1486      if (!nodeInfo->numberBranchesLeft())
1487        nodeInfo->allBranchesGone(); // can clean up
1488      }
1489
1490      if ((specialOptions_&1)!=0) {
1491        /*
1492          This doesn't work as intended --- getRowCutDebugger will return null
1493          unless the current feasible solution region includes the optimal solution
1494          that RowCutDebugger knows. There's no way to tell inactive from off the
1495          optimal path.
1496        */
1497        const OsiRowCutDebugger *debugger = solver_->getRowCutDebugger() ;
1498        if (debugger) {
1499          onOptimalPath=true;
1500          printf("On optimal path\n") ;
1501        }
1502      }
1503/*
1504  Reoptimize, possibly generating cuts and/or using heuristics to find
1505  solutions.  Cut reference counts are unaffected unless we lose feasibility,
1506  in which case solveWithCuts() will make the adjustment.
1507*/
1508      phase_=2;
1509      cuts = OsiCuts() ;
1510      currentNumberCuts = solver_->getNumRows()-numberRowsAtContinuous_ ;
1511      int saveNumber = numberIterations_;
1512      if(solverCharacteristics_->solutionAddsCuts()) {
1513        int returnCode=resolve(node ? node->nodeInfo() : NULL,1);
1514        feasible = returnCode != 0;
1515        if (feasible) {
1516          int iObject ;
1517          int preferredWay ;
1518          int numberUnsatisfied = 0 ;
1519          memcpy(currentSolution_,solver_->getColSolution(),
1520                 numberColumns*sizeof(double)) ;
1521          // point to useful information
1522          OsiBranchingInformation usefulInfo=usefulInformation();
1523         
1524          for (iObject = 0 ; iObject < numberObjects_ ; iObject++) {
1525            double infeasibility =
1526              object_[iObject]->infeasibility(&usefulInfo,preferredWay) ;
1527            if (infeasibility ) numberUnsatisfied++ ;
1528          }
1529          if (returnCode>0) {
1530            if (numberUnsatisfied)   {
1531              feasible = solveWithCuts(cuts,maximumCutPasses_,node);
1532            } else {
1533              // may generate cuts and turn the solution
1534              //to an infeasible one
1535              feasible = solveWithCuts(cuts, 1,
1536                                       node);
1537#if 0
1538              currentNumberCuts_ = cuts.sizeRowCuts();
1539              if (currentNumberCuts_ >= maximumNumberCuts_) {
1540                maximumNumberCuts_ = currentNumberCuts;
1541                delete [] addedCuts_;
1542                addedCuts_ = new CbcCountRowCut * [maximumNumberCuts_];
1543              }
1544#endif
1545            }
1546          }
1547          // check extra info on feasibility
1548          if (!solverCharacteristics_->mipFeasible()) {
1549            feasible = false;
1550            solverCharacteristics_->setMipBound(-COIN_DBL_MAX);
1551          }
1552        }
1553      } else {
1554        // normal
1555        //int zzzzzz=0;
1556        //if (zzzzzz)
1557        //solver_->writeMps("before");
1558        feasible = solveWithCuts(cuts,maximumCutPasses_,node);
1559      }
1560      if ((specialOptions_&1)!=0&&onOptimalPath) {
1561        if (!solver_->getRowCutDebugger()) {
1562          // dj fix did something???
1563          solver_->writeMps("infeas2");
1564          assert (solver_->getRowCutDebugger()) ;
1565        }
1566      }
1567      if (statistics_) {
1568        assert (numberNodes2_);
1569        assert (statistics_[numberNodes2_-1]);
1570        assert (statistics_[numberNodes2_-1]->node()==numberNodes2_-1);
1571        statistics_[numberNodes2_-1]->endOfBranch(numberIterations_-saveNumber,
1572                                               feasible ? solver_->getObjValue()
1573                                               : COIN_DBL_MAX);
1574      }
1575/*
1576  Check for abort on limits: node count, solution count, time, integrality gap.
1577*/
1578      double totalTime = CoinCpuTime()-dblParam_[CbcStartSeconds] ;
1579      if (numberNodes_ < intParam_[CbcMaxNumNode] &&
1580          numberSolutions_ < intParam_[CbcMaxNumSol] &&
1581          totalTime < dblParam_[CbcMaximumSeconds] &&
1582          !stoppedOnGap_&&!eventHappened_) 
1583      {
1584/*
1585  Are we still feasible? If so, create a node and do the work to attach a
1586  branching object, reoptimising as needed if chooseBranch() identifies
1587  monotone objects.
1588
1589  Finally, attach a partial nodeInfo object and store away any cuts that we
1590  created back in solveWithCuts. addCuts() will initialise the reference
1591  counts for these new cuts.
1592
1593  This next test can be problematic if we've discovered an
1594  alternate equivalent answer and subsequently fathom the solution
1595  known to the row cut debugger due to bounds.
1596*/
1597        if (onOptimalPath) {
1598          bool objLim = solver_->isDualObjectiveLimitReached() ;
1599          if (!feasible && !objLim) {
1600            printf("infeas2\n");
1601            solver_->writeMps("infeas");
1602            CoinWarmStartBasis *slack =
1603              dynamic_cast<CoinWarmStartBasis *>(solver_->getEmptyWarmStart()) ;
1604            solver_->setWarmStart(slack);
1605            delete slack ;
1606            solver_->setHintParam(OsiDoReducePrint,false,OsiHintDo,0) ;
1607            solver_->initialSolve();
1608            assert (!solver_->isProvenOptimal());
1609          }
1610          assert (feasible || objLim);
1611        }
1612        bool checkingNode=false;
1613        if (feasible) {
1614          newNode = new CbcNode ;
1615          // Set objective value (not so obvious if NLP etc)
1616          setObjectiveValue(newNode,node);
1617          anyAction =-1 ;
1618          resolved = false ;
1619          if (newNode->objectiveValue() >= getCutoff()) 
1620            anyAction=-2;
1621          // only allow at most a few passes
1622          int numberPassesLeft=5;
1623          checkingNode=true;
1624#if 0
1625          while (anyAction == -1)
1626          {
1627            // Set objective value (not so obvious if NLP etc)
1628            setObjectiveValue(newNode,node);
1629            if (numberBeforeTrust_==0 ) {
1630              anyAction = newNode->chooseBranch(this,node,numberPassesLeft) ;
1631            } else {
1632              OsiSolverBranch * branches=NULL;
1633              anyAction = newNode->chooseDynamicBranch(this,node,branches,numberPassesLeft) ;
1634              if (anyAction==-3)
1635                anyAction = newNode->chooseBranch(this,node,numberPassesLeft) ; // dynamic did nothing
1636            }
1637            if (solverCharacteristics_ &&
1638                solverCharacteristics_->solutionAddsCuts() && // we are in some OA based bab
1639                feasible && (newNode->numberUnsatisfied()==0) //solution has become integer feasible during strong branching
1640                )
1641              {
1642                //in the present case we need to check here integer infeasibility if the node is not fathomed we will have to do the loop
1643                // again
1644                std::cout<<solver_<<std::endl;
1645                resolve(solver_);
1646                double objval = solver_->getObjValue();
1647                setBestSolution(CBC_SOLUTION, objval,
1648                                solver_->getColSolution()) ;
1649                lastHeuristic_ = NULL;
1650                int easy=2;
1651                if (!solverCharacteristics_->mipFeasible())//did we prove that the node could be pruned?
1652                  feasible = false;
1653                // Reset the bound now
1654                solverCharacteristics_->setMipBound(-COIN_DBL_MAX);
1655               
1656               
1657                //numberPassesLeft++;
1658                solver_->setHintParam(OsiDoInBranchAndCut,true,OsiHintDo,&easy) ;
1659                feasible &= resolve(node ? node->nodeInfo() : NULL,11) != 0 ;
1660                solver_->setHintParam(OsiDoInBranchAndCut,true,OsiHintDo,NULL) ;
1661                resolved = true ;
1662                if (problemFeasibility_->feasible(this,0)<0) {
1663                  feasible=false; // pretend infeasible
1664                }
1665                if(feasible)
1666                  anyAction = -1;
1667                else
1668                  anyAction = -2;
1669              }
1670/*
1671  Yep, false positives for sure. And no easy way to distinguish honest
1672  infeasibility from `found a solution and tightened objective target.'
1673
1674            if (onOptimalPath)
1675              assert (anyAction!=-2); // can be useful but gives false positives on strong
1676*/
1677            numberPassesLeft--;
1678            if (numberPassesLeft<=-1) {
1679              if (!numberLongStrong_)
1680                messageHandler()->message(CBC_WARNING_STRONG,
1681                                          messages()) << CoinMessageEol ;
1682              numberLongStrong_++;
1683            }
1684            if (anyAction == -1)
1685            {
1686              // can do quick optimality check
1687              int easy=2;
1688              solver_->setHintParam(OsiDoInBranchAndCut,true,OsiHintDo,&easy) ;
1689              //if (solverCharacteristics_->solverType()!=4)
1690                feasible = resolve(node ? node->nodeInfo() : NULL,11) != 0 ;
1691                //else
1692                //feasible = solveWithCuts(cuts,maximumCutPasses_,node);
1693              solver_->setHintParam(OsiDoInBranchAndCut,true,OsiHintDo,NULL) ;
1694              resolved = true ;
1695              if (problemFeasibility_->feasible(this,0)<0) {
1696                feasible=false; // pretend infeasible
1697              }
1698              if (feasible)
1699              {
1700                // Set objective value (not so obvious if NLP etc)
1701                setObjectiveValue(newNode,node);
1702                reducedCostFix() ;
1703                if (newNode->objectiveValue() >= getCutoff())
1704                  anyAction=-2;
1705              }
1706              else
1707              { anyAction = -2 ; } } }
1708          if (anyAction >= 0)
1709          { if (resolved)
1710            { bool needValidSolution = (newNode->branchingObject() == NULL) ;
1711              takeOffCuts(cuts,needValidSolution,NULL) ;
1712#             ifdef CHECK_CUT_COUNTS
1713              { printf("Number of rows after chooseBranch fix (node)"
1714                       "(active only) %d\n",
1715                        numberRowsAtContinuous_+numberNewCuts_+
1716                        numberOldActiveCuts_) ;
1717                const CoinWarmStartBasis* debugws =
1718                  dynamic_cast<const CoinWarmStartBasis*>
1719                    (solver_->getWarmStart()) ;
1720                debugws->print() ;
1721                delete debugws ; }
1722#             endif
1723            }
1724            newNode->createInfo(this,node,lastws,lowerBefore,upperBefore,
1725                                numberOldActiveCuts_,numberNewCuts_) ;
1726            if (newNode->numberUnsatisfied()) {
1727              maximumDepthActual_ = CoinMax(maximumDepthActual_,newNode->depth());
1728              newNode->initializeInfo() ;
1729              newNode->nodeInfo()->addCuts(cuts,newNode->numberBranches(),
1730                                           whichGenerator_) ; } } }
1731        else {
1732          anyAction = -2 ;
1733          // Reset bound anyway (no harm if not odd)
1734          solverCharacteristics_->setMipBound(-COIN_DBL_MAX);
1735          node->nodeInfo()->decrement();
1736        }
1737        // May have slipped through i.e. anyAction == 0 and objective above cutoff
1738        // I think this will screw up cut reference counts if executed.
1739        // We executed addCuts just above. (lh)
1740        if ( anyAction >=0 ) {
1741          assert (newNode);
1742          if (newNode->objectiveValue() >= getCutoff()) {
1743            anyAction = -2; // say bad after all
1744#ifdef COIN_DEVELOP
1745            printf("zapping2 CbcNodeInfo %x\n",newNode->nodeInfo()->parent());
1746#endif
1747            // zap parent nodeInfo
1748            if (newNode->nodeInfo())
1749              newNode->nodeInfo()->nullParent();
1750          }
1751        }
1752/*
1753  If we end up infeasible, we can delete the new node immediately. Since this
1754  node won't be needing the cuts we collected, decrement the reference counts.
1755  If we are feasible, then we'll be placing this node into the live set, so
1756  increment the reference count in the current (parent) nodeInfo.
1757*/
1758        if (anyAction == -2)
1759        { delete newNode ;
1760          newNode = NULL ;
1761          // say strong doing well
1762          if (checkingNode)
1763            setSpecialOptions(specialOptions_|8);
1764          for (i = 0 ; i < currentNumberCuts_ ; i++)
1765          { if (addedCuts_[i])
1766            { if (!addedCuts_[i]->decrement(1))
1767                delete addedCuts_[i] ; } } }
1768        else
1769        { nodeInfo->increment() ;
1770        if ((numberNodes_%20)==0) {
1771          // say strong not doing as well
1772          setSpecialOptions(specialOptions_&~8);
1773        }
1774#else
1775        OsiSolverBranch * branches=NULL;
1776        // point to useful information
1777        anyAction = chooseBranch(newNode, numberPassesLeft,node, cuts,resolved,
1778                                 lastws, lowerBefore, upperBefore, branches);
1779/*
1780  If we end up infeasible, we can delete the new node immediately. Since this
1781  node won't be needing the cuts we collected, decrement the reference counts.
1782  If we are feasible, then we'll be placing this node into the live set, so
1783  increment the reference count in the current (parent) nodeInfo.
1784*/
1785        if (anyAction == -2)
1786          { delete newNode ;
1787          newNode = NULL ;
1788          // say strong doing well
1789          if (checkingNode)
1790            setSpecialOptions(specialOptions_|8);
1791          for (i = 0 ; i < currentNumberCuts_ ; i++)
1792            { if (addedCuts_[i])
1793              { if (!addedCuts_[i]->decrement(1))
1794                delete addedCuts_[i] ; } } }
1795        else
1796          { nodeInfo->increment() ;
1797          if ((numberNodes_%20)==0) {
1798            // say strong not doing as well
1799            setSpecialOptions(specialOptions_&~8);
1800          }
1801        }
1802#endif
1803        }
1804/*
1805  At this point, there are three possibilities:
1806    * newNode is live and will require further branching to resolve
1807      (variable() >= 0). Increment the cut reference counts by
1808      numberBranches() to allow for use by children of this node, and
1809      decrement by 1 because we've executed one arm of the branch of our
1810      parent (consuming one reference). Before we push newNode onto the
1811      search tree, try for a heuristic solution.
1812    * We have a solution, in which case newNode is non-null but we have no
1813      branching variable. Decrement the cut counts and save the solution.
1814    * The node was found to be infeasible, in which case it's already been
1815      deleted, and newNode is null.
1816*/
1817        if (!eventHandler->event(CbcEventHandler::node)) {
1818          eventHappened_=true; // exit
1819        }
1820        assert (!newNode || newNode->objectiveValue() <= getCutoff()) ;
1821        if (statistics_) {
1822          assert (numberNodes2_);
1823          assert (statistics_[numberNodes2_-1]);
1824          assert (statistics_[numberNodes2_-1]->node()==numberNodes2_-1);
1825          if (newNode)
1826            statistics_[numberNodes2_-1]->updateInfeasibility(newNode->numberUnsatisfied());
1827          else
1828            statistics_[numberNodes2_-1]->sayInfeasible();
1829        }
1830        if (newNode) {
1831          if (newNode->branchingObject() == NULL&&solverCharacteristics_->solverType()==4) {
1832            // need to check if any cuts would do anything
1833            OsiCuts theseCuts;
1834            // reset probing info
1835            if (probingInfo_)
1836              probingInfo_->initializeFixing();
1837            for (int i = 0;i<numberCutGenerators_;i++) {
1838              bool generate = generator_[i]->normal();
1839              // skip if not optimal and should be (maybe a cut generator has fixed variables)
1840              if (generator_[i]->needsOptimalBasis()&&!solver_->basisIsAvailable())
1841                generate=false;
1842              if (!generator_[i]->mustCallAgain())
1843                generate=false; // only special cuts
1844              if (generate) {
1845                generator_[i]->generateCuts(theseCuts,true,NULL) ;
1846                int numberRowCutsAfter = theseCuts.sizeRowCuts() ;
1847                if (numberRowCutsAfter) {
1848                  // need dummy branch
1849                  newNode->setBranchingObject(new CbcDummyBranchingObject(this));
1850                  newNode->nodeInfo()->initializeInfo(1);
1851                  break;
1852                }
1853              }
1854            }
1855          }
1856          if (newNode->branchingObject())
1857          { handler_->message(CBC_BRANCH,messages_)
1858               << numberNodes_<< newNode->objectiveValue()
1859               << newNode->numberUnsatisfied()<< newNode->depth()
1860               << CoinMessageEol ;
1861            // Increment cut counts (taking off current)
1862            int numberLeft = newNode->numberBranches() ;
1863            for (i = 0;i < currentNumberCuts_;i++)
1864            { if (addedCuts_[i])
1865              {
1866#               ifdef CHECK_CUT_COUNTS
1867                printf("Count on cut %x increased by %d\n",addedCuts_[i],
1868                        numberLeft-1) ;
1869#               endif
1870                addedCuts_[i]->increment(numberLeft-1) ; } }
1871
1872            double estValue = newNode->guessedObjectiveValue() ;
1873            int found = -1 ;
1874            // no - overhead on small problems solver_->resolve() ;     // double check current optimal
1875            // assert (!solver_->getIterationCount());
1876            double * newSolution = new double [numberColumns] ;
1877            double heurValue = getCutoff() ;
1878            int iHeur ;
1879            for (iHeur = 0 ; iHeur < numberHeuristics_ ; iHeur++)
1880            { double saveValue = heurValue ;
1881              int ifSol = heuristic_[iHeur]->solution(heurValue,newSolution) ;
1882              if (ifSol > 0) {
1883                // new solution found
1884                found = iHeur ;
1885                incrementUsed(newSolution);
1886              }
1887              else
1888              if (ifSol < 0)    // just returning an estimate
1889              { estValue = CoinMin(heurValue,estValue) ;
1890                heurValue = saveValue ; } }
1891            if (found >= 0) {
1892              setBestSolution(CBC_ROUNDING,heurValue,newSolution) ;
1893              lastHeuristic_ = heuristic_[found];
1894            }
1895            delete [] newSolution ;
1896            newNode->setGuessedObjectiveValue(estValue) ;
1897            tree_->push(newNode) ;
1898            if (statistics_) {
1899              if (numberNodes2_==maximumStatistics_) {
1900                maximumStatistics_ = 2*maximumStatistics_;
1901                CbcStatistics ** temp = new CbcStatistics * [maximumStatistics_];
1902                memset(temp,0,maximumStatistics_*sizeof(CbcStatistics *));
1903                memcpy(temp,statistics_,numberNodes2_*sizeof(CbcStatistics *));
1904                delete [] statistics_;
1905                statistics_=temp;
1906              }
1907              assert (!statistics_[numberNodes2_]);
1908              statistics_[numberNodes2_]=new CbcStatistics(newNode);
1909            }
1910            numberNodes2_++;
1911#           ifdef CHECK_NODE
1912            printf("Node %x pushed on tree c\n",newNode) ;
1913#           endif
1914          }
1915          else
1916          { 
1917            if(solverCharacteristics_ && //we may be in a non standard bab
1918               solverCharacteristics_->solutionAddsCuts()// we are in some kind of OA based bab.
1919               )
1920              {
1921                std::cerr<<"You should never get here"<<std::endl;
1922                throw CoinError("Nodes should not be fathomed on integer infeasibility in this setting",
1923                                "branchAndBound","CbcModel") ;
1924              }
1925            for (i = 0 ; i < currentNumberCuts_ ; i++)
1926            { if (addedCuts_[i])
1927              { if (!addedCuts_[i]->decrement(1))
1928                  delete addedCuts_[i] ; } }
1929          double objectiveValue = newNode->objectiveValue();
1930            setBestSolution(CBC_SOLUTION,objectiveValue,
1931                            solver_->getColSolution()) ;
1932            lastHeuristic_ = NULL;
1933            incrementUsed(solver_->getColSolution());
1934            //assert(nodeInfo->numberPointingToThis() <= 2) ;
1935            // avoid accidental pruning, if newNode was final branch arm
1936            nodeInfo->increment();
1937            delete newNode ;
1938            nodeInfo->decrement() ; } }
1939/*
1940  This node has been completely expanded and can be removed from the live
1941  set.
1942*/
1943        if (deleteNode) 
1944          delete node ; }
1945/*
1946  End of the non-abort actions. The next block of code is executed if we've
1947  aborted because we hit one of the limits. Clean up by deleting the live set
1948  and break out of the node processing loop. Note that on an abort, node may
1949  have been pushed back onto the tree for further processing, in which case
1950  it'll be deleted in cleanTree. We need to check.
1951*/
1952      else
1953      { 
1954        tree_->cleanTree(this,-COIN_DBL_MAX,bestPossibleObjective_) ;
1955        delete nextRowCut_;
1956        // We need to get rid of node if is has already been popped from tree
1957        if (!nodeOnTree&&!stoppedOnGap_&&node!=rootNode)
1958          delete node;
1959        if (stoppedOnGap_)
1960        { messageHandler()->message(CBC_GAP,messages())
1961            << bestObjective_-bestPossibleObjective_
1962            << dblParam_[CbcAllowableGap]
1963            << dblParam_[CbcAllowableFractionGap]*100.0
1964            << CoinMessageEol ;
1965        secondaryStatus_ = 2;
1966          status_ = 0 ; }
1967        else
1968        if (isNodeLimitReached())
1969        { handler_->message(CBC_MAXNODES,messages_) << CoinMessageEol ;
1970        secondaryStatus_ = 3;
1971          status_ = 1 ; }
1972        else
1973        if (totalTime >= dblParam_[CbcMaximumSeconds])
1974        { handler_->message(CBC_MAXTIME,messages_) << CoinMessageEol ; 
1975        secondaryStatus_ = 4;
1976          status_ = 1 ; }
1977        else
1978        if (eventHappened_)
1979        { handler_->message(CBC_EVENT,messages_) << CoinMessageEol ; 
1980        secondaryStatus_ = 5;
1981          status_ = 5 ; }
1982        else
1983        { handler_->message(CBC_MAXSOLS,messages_) << CoinMessageEol ;
1984        secondaryStatus_ = 6;
1985          status_ = 1 ; }
1986        break ; }
1987/*
1988  Delete cuts to get back to the original system.
1989
1990  I'm thinking this is redundant --- the call to addCuts that conditions entry
1991  to this code block also performs this action.
1992*/
1993      int numberToDelete = getNumRows()-numberRowsAtContinuous_ ;
1994      if (numberToDelete)
1995      { int * delRows = new int[numberToDelete] ;
1996        int i ;
1997        for (i = 0 ; i < numberToDelete ; i++)
1998        { delRows[i] = i+numberRowsAtContinuous_ ; }
1999        solver_->deleteRows(numberToDelete,delRows) ;
2000        delete [] delRows ; } }
2001/*
2002  This node fathomed when addCuts atttempted to revive it. Toss it.
2003*/
2004    else
2005      { delete node ; } }
2006/*
2007  That's it, we've exhausted the search tree, or broken out of the loop because
2008  we hit some limit on evaluation.
2009
2010  We may have got an intelligent tree so give it one more chance
2011*/
2012  // Tell solver we are not in Branch and Cut
2013  solver_->setHintParam(OsiDoInBranchAndCut,false,OsiHintDo,NULL) ;
2014  tree_->endSearch();
2015  //  If we did any sub trees - did we give up on any?
2016  if ( numberStoppedSubTrees_)
2017    status_=1;
2018  if (!status_) {
2019    // Set best possible unless stopped on gap
2020    if(secondaryStatus_ != 2)
2021      bestPossibleObjective_=bestObjective_;
2022    handler_->message(CBC_END_GOOD,messages_)
2023      << bestObjective_ << numberIterations_ << numberNodes_<<getCurrentSeconds()
2024      << CoinMessageEol ;
2025  } else {
2026    handler_->message(CBC_END,messages_)
2027      << bestObjective_ <<bestPossibleObjective_
2028      << numberIterations_ << numberNodes_<<getCurrentSeconds()
2029      << CoinMessageEol ;
2030  }
2031  if (numberStrongIterations_)
2032    handler_->message(CBC_STRONG_STATS,messages_)
2033      << strongInfo_[0] << numberStrongIterations_ << strongInfo_[2]
2034      << strongInfo_[1] << CoinMessageEol ;
2035  handler_->message(CBC_OTHER_STATS,messages_)
2036      << maximumDepthActual_
2037      << numberDJFixed_ << CoinMessageEol ;
2038  if (statistics_) {
2039    // report in some way
2040    int * lookup = new int[numberObjects_];
2041    int i;
2042    for (i=0;i<numberObjects_;i++) 
2043      lookup[i]=-1;
2044    bool goodIds=true;
2045    for (i=0;i<numberObjects_;i++) {
2046      int iColumn = object_[i]->columnNumber();
2047      if(iColumn>=0&&iColumn<numberColumns) {
2048        if (lookup[i]==-1) {
2049          lookup[i]=iColumn;
2050        } else {
2051          goodIds=false;
2052          break;
2053        }
2054      } else {
2055        goodIds=false;
2056        break;
2057      }
2058    }
2059    if (!goodIds) {
2060      delete [] lookup;
2061      lookup=NULL;
2062    }
2063    if (doStatistics==3) {
2064      printf("  node parent depth column   value                    obj      inf\n");
2065      for ( i=0;i<numberNodes2_;i++) {
2066        statistics_[i]->print(lookup);
2067      }
2068    }
2069    if (doStatistics>1) {
2070      // Find last solution
2071      int k;
2072      for (k=numberNodes2_-1;k>=0;k--) {
2073        if (statistics_[k]->endingObjective()!=COIN_DBL_MAX&&
2074            !statistics_[k]->endingInfeasibility())
2075          break;
2076      }
2077      if (k>=0) {
2078        int depth=statistics_[k]->depth();
2079        int * which = new int[depth+1];
2080        for (i=depth;i>=0;i--) {
2081          which[i]=k;
2082          k=statistics_[k]->parentNode();
2083        }
2084        printf("  node parent depth column   value                    obj      inf\n");
2085        for (i=0;i<=depth;i++) {
2086          statistics_[which[i]]->print(lookup);
2087        }
2088        delete [] which;
2089      }
2090    }
2091    // now summary
2092    int maxDepth=0;
2093    double averageSolutionDepth=0.0;
2094    int numberSolutions=0;
2095    double averageCutoffDepth=0.0;
2096    double averageSolvedDepth=0.0;
2097    int numberCutoff=0;
2098    int numberDown=0;
2099    int numberFirstDown=0;
2100    double averageInfDown=0.0;
2101    double averageObjDown=0.0;
2102    int numberCutoffDown=0;
2103    int numberUp=0;
2104    int numberFirstUp=0;
2105    double averageInfUp=0.0;
2106    double averageObjUp=0.0;
2107    int numberCutoffUp=0;
2108    double averageNumberIterations1=0.0;
2109    double averageValue=0.0;
2110    for ( i=0;i<numberNodes2_;i++) {
2111      int depth =  statistics_[i]->depth(); 
2112      int way =  statistics_[i]->way(); 
2113      double value = statistics_[i]->value(); 
2114      double startingObjective =  statistics_[i]->startingObjective(); 
2115      int startingInfeasibility = statistics_[i]->startingInfeasibility(); 
2116      double endingObjective = statistics_[i]->endingObjective(); 
2117      int endingInfeasibility = statistics_[i]->endingInfeasibility(); 
2118      maxDepth = CoinMax(depth,maxDepth);
2119      // Only for completed
2120      averageNumberIterations1 += statistics_[i]->numberIterations();
2121      averageValue += value;
2122      if (endingObjective!=COIN_DBL_MAX&&!endingInfeasibility) {
2123        numberSolutions++;
2124        averageSolutionDepth += depth;
2125      }
2126      if (endingObjective==COIN_DBL_MAX) {
2127        numberCutoff++;
2128        averageCutoffDepth += depth;
2129        if (way<0) {
2130          numberDown++;
2131          numberCutoffDown++;
2132          if (way==-1)
2133            numberFirstDown++;
2134        } else {
2135          numberUp++;
2136          numberCutoffUp++;
2137          if (way==1)
2138            numberFirstUp++;
2139        }
2140      } else {
2141        averageSolvedDepth += depth;
2142        if (way<0) {
2143          numberDown++;
2144          averageInfDown += startingInfeasibility-endingInfeasibility;
2145          averageObjDown += endingObjective-startingObjective;
2146          if (way==-1)
2147            numberFirstDown++;
2148        } else {
2149          numberUp++;
2150          averageInfUp += startingInfeasibility-endingInfeasibility;
2151          averageObjUp += endingObjective-startingObjective;
2152          if (way==1)
2153            numberFirstUp++;
2154        }
2155      }
2156    }
2157    // Now print
2158    if (numberSolutions)
2159      averageSolutionDepth /= (double) numberSolutions;
2160    int numberSolved = numberNodes2_-numberCutoff;
2161    double averageNumberIterations2=numberIterations_-averageNumberIterations1
2162      -numberIterationsAtContinuous;
2163    if(numberCutoff) {
2164      averageCutoffDepth /= (double) numberCutoff;
2165      averageNumberIterations2 /= (double) numberCutoff;
2166    }
2167    if (numberNodes2_) 
2168      averageValue /= (double) numberNodes2_;
2169    if (numberSolved) {
2170      averageNumberIterations1 /= (double) numberSolved;
2171      averageSolvedDepth /= (double) numberSolved;
2172    }
2173    printf("%d solution(s) were found (by branching) at an average depth of %g\n",
2174           numberSolutions,averageSolutionDepth);
2175    printf("average value of variable being branched on was %g\n",
2176           averageValue);
2177    printf("%d nodes were cutoff at an average depth of %g with iteration count of %g\n",
2178           numberCutoff,averageCutoffDepth,averageNumberIterations2);
2179    printf("%d nodes were solved at an average depth of %g with iteration count of %g\n",
2180           numberSolved,averageSolvedDepth,averageNumberIterations1);
2181    if (numberDown) {
2182      averageInfDown /= (double) numberDown;
2183      averageObjDown /= (double) numberDown;
2184    }
2185    printf("Down %d nodes (%d first, %d second) - %d cutoff, rest decrease numinf %g increase obj %g\n",
2186           numberDown,numberFirstDown,numberDown-numberFirstDown,numberCutoffDown,
2187           averageInfDown,averageObjDown);
2188    if (numberUp) {
2189      averageInfUp /= (double) numberUp;
2190      averageObjUp /= (double) numberUp;
2191    }
2192    printf("Up %d nodes (%d first, %d second) - %d cutoff, rest decrease numinf %g increase obj %g\n",
2193           numberUp,numberFirstUp,numberUp-numberFirstUp,numberCutoffUp,
2194           averageInfUp,averageObjUp);
2195    for ( i=0;i<numberNodes2_;i++) 
2196      delete statistics_[i];
2197    delete [] statistics_;
2198    statistics_=NULL;
2199    maximumStatistics_=0;
2200    delete [] lookup;
2201  }
2202/*
2203  If we think we have a solution, restore and confirm it with a call to
2204  setBestSolution().  We need to reset the cutoff value so as not to fathom
2205  the solution on bounds.  Note that calling setBestSolution( ..., true)
2206  leaves the continuousSolver_ bounds vectors fixed at the solution value.
2207
2208  Running resolve() here is a failsafe --- setBestSolution has already
2209  reoptimised using the continuousSolver_. If for some reason we fail to
2210  prove optimality, run the problem again after instructing the solver to
2211  tell us more.
2212
2213  If all looks good, replace solver_ with continuousSolver_, so that the
2214  outside world will be able to obtain information about the solution using
2215  public methods.
2216*/
2217  if (bestSolution_&&(solverCharacteristics_->solverType()<2||solverCharacteristics_->solverType()==4)) 
2218  { setCutoff(1.0e50) ; // As best solution should be worse than cutoff
2219    phase_=5;
2220    double increment = getDblParam(CbcModel::CbcCutoffIncrement) ;
2221    if ((specialOptions_&4)!=0)
2222      bestObjective_ += 100.0*increment+1.0e-3; // only set if we are going to solve
2223    setBestSolution(CBC_END_SOLUTION,bestObjective_,bestSolution_,true) ;
2224    continuousSolver_->resolve() ;
2225    if (!continuousSolver_->isProvenOptimal())
2226    { continuousSolver_->messageHandler()->setLogLevel(2) ;
2227      continuousSolver_->initialSolve() ; }
2228    delete solver_ ;
2229    // above deletes solverCharacteristics_
2230    solverCharacteristics_ = NULL;
2231    solver_ = continuousSolver_ ;
2232    setPointers(solver_);
2233    continuousSolver_ = NULL ; }
2234/*
2235  Clean up dangling objects. continuousSolver_ may already be toast.
2236*/
2237  delete lastws ;
2238  delete [] whichGenerator_ ;
2239  whichGenerator_=NULL;
2240  delete [] lowerBefore ;
2241  delete [] upperBefore ;
2242  delete [] walkback_ ;
2243  walkback_ = NULL ;
2244  delete [] addedCuts_ ;
2245  addedCuts_ = NULL ;
2246  //delete persistentInfo;
2247  // Get rid of characteristics
2248  solverCharacteristics_=NULL;
2249  if (continuousSolver_)
2250  { delete continuousSolver_ ;
2251    continuousSolver_ = NULL ; }
2252/*
2253  Destroy global cuts by replacing with an empty OsiCuts object.
2254*/
2255  globalCuts_= OsiCuts() ;
2256  if (!bestSolution_) {
2257    // make sure lp solver is infeasible
2258    int numberColumns = solver_->getNumCols();
2259    const double * columnLower = solver_->getColLower();
2260    int iColumn;
2261    for (iColumn=0;iColumn<numberColumns;iColumn++) {
2262      if (solver_->isInteger(iColumn)) 
2263        solver_->setColUpper(iColumn,columnLower[iColumn]);
2264    }
2265    solver_->initialSolve();
2266  }
2267  if (strategy_&&strategy_->preProcessState()>0) {
2268    // undo preprocessing
2269    CglPreProcess * process = strategy_->process();
2270    assert (process);
2271    int n = originalSolver->getNumCols();
2272    if (bestSolution_) {
2273      delete [] bestSolution_;
2274      bestSolution_ = new double [n];
2275      process->postProcess(*solver_);
2276    }
2277    strategy_->deletePreProcess();
2278    // Solution now back in originalSolver
2279    delete solver_;
2280    solver_=originalSolver;
2281    if (bestSolution_)
2282      memcpy(bestSolution_,solver_->getColSolution(),n*sizeof(double));
2283    // put back original objects if there were any
2284    if (originalObject) {
2285      int iColumn;
2286      assert (ownObjects_);
2287      for (iColumn=0;iColumn<numberObjects_;iColumn++) 
2288        delete object_[iColumn];
2289      delete [] object_;
2290      numberObjects_ = numberOriginalObjects;
2291      object_=originalObject;
2292      delete [] integerVariable_;
2293      numberIntegers_=0;
2294      for (iColumn=0;iColumn<n;iColumn++) {
2295        if (solver_->isInteger(iColumn))
2296          numberIntegers_++;
2297      }
2298      integerVariable_ = new int[numberIntegers_];
2299      numberIntegers_=0;
2300      for (iColumn=0;iColumn<n;iColumn++) {
2301        if (solver_->isInteger(iColumn))
2302            integerVariable_[numberIntegers_++]=iColumn;
2303      }
2304    }
2305  }
2306#ifdef CLP_QUICK_OPTIONS
2307  {
2308    OsiClpSolverInterface * clpSolver
2309      = dynamic_cast<OsiClpSolverInterface *> (solver_);
2310    if (clpSolver) {
2311      // Try and re-use regions   
2312      ClpSimplex * simplex = clpSolver->getModelPtr();
2313      simplex->setPersistenceFlag(0);
2314      clpSolver->deleteScaleFactors();
2315      clpSolver->setSpecialOptions(clpSolver->specialOptions()&(~131072));
2316      simplex->setSpecialOptions(simplex->specialOptions()&(~131072));
2317      simplex->setAlphaAccuracy(-1.0);
2318      //clpSolver->setSpecialOptions((clpSolver->specialOptions()&~128)|65536);
2319    }
2320  }
2321#endif
2322  return ;
2323 }
2324
2325
2326// Solve the initial LP relaxation
2327void 
2328CbcModel::initialSolve() 
2329{
2330  assert (solver_);
2331  assert (!solverCharacteristics_);
2332  OsiBabSolver * solverCharacteristics = dynamic_cast<OsiBabSolver *> (solver_->getAuxiliaryInfo());
2333  if (solverCharacteristics) {
2334    solverCharacteristics_ = solverCharacteristics;
2335  } else {
2336    // replace in solver
2337    OsiBabSolver defaultC;
2338    solver_->setAuxiliaryInfo(&defaultC);
2339    solverCharacteristics_ = dynamic_cast<OsiBabSolver *> (solver_->getAuxiliaryInfo());
2340  }
2341  solverCharacteristics_->setSolver(solver_);
2342  solver_->setHintParam(OsiDoInBranchAndCut,true,OsiHintDo,NULL) ;
2343  solver_->initialSolve();
2344  solver_->setHintParam(OsiDoInBranchAndCut,false,OsiHintDo,NULL) ;
2345  // But set up so Jon Lee will be happy
2346  status_=-1;
2347  secondaryStatus_ = -1;
2348  originalContinuousObjective_ = solver_->getObjValue()*solver_->getObjSense();
2349  delete [] continuousSolution_;
2350  continuousSolution_ = CoinCopyOfArray(solver_->getColSolution(),
2351                                             solver_->getNumCols());
2352  setPointers(solver_);
2353  solverCharacteristics_ = NULL;
2354}
2355
2356/*! \brief Get an empty basis object
2357
2358  Return an empty CoinWarmStartBasis object with the requested capacity,
2359  appropriate for the current solver. The object is cloned from the object
2360  cached as emptyWarmStart_. If there is no cached object, the routine
2361  queries the solver for a warm start object, empties it, and caches the
2362  result.
2363*/
2364
2365CoinWarmStartBasis *CbcModel::getEmptyBasis (int ns, int na) const
2366
2367{ CoinWarmStartBasis *emptyBasis ;
2368/*
2369  Acquire an empty basis object, if we don't yet have one.
2370*/
2371  if (emptyWarmStart_ == 0)
2372  { if (solver_ == 0)
2373    { throw CoinError("Cannot construct basis without solver!",
2374                      "getEmptyBasis","CbcModel") ; }
2375    emptyBasis =
2376        dynamic_cast<CoinWarmStartBasis *>(solver_->getEmptyWarmStart()) ;
2377    if (emptyBasis == 0)
2378    { throw CoinError(
2379        "Solver does not appear to use a basis-oriented warm start.",
2380        "getEmptyBasis","CbcModel") ; }
2381    emptyBasis->setSize(0,0) ;
2382    emptyWarmStart_ = dynamic_cast<CoinWarmStart *>(emptyBasis) ; }
2383/*
2384  Clone the empty basis object, resize it as requested, and return.
2385*/
2386  emptyBasis = dynamic_cast<CoinWarmStartBasis *>(emptyWarmStart_->clone()) ;
2387  assert(emptyBasis) ;
2388  if (ns != 0 || na != 0) emptyBasis->setSize(ns,na) ;
2389
2390  return (emptyBasis) ; }
2391   
2392
2393/** Default Constructor
2394
2395  Creates an empty model without an associated solver.
2396*/
2397CbcModel::CbcModel() 
2398
2399:
2400  solver_(NULL),
2401  ourSolver_(true),
2402  continuousSolver_(NULL),
2403  referenceSolver_(NULL),
2404  defaultHandler_(true),
2405  emptyWarmStart_(NULL),
2406  bestObjective_(COIN_DBL_MAX),
2407  bestPossibleObjective_(COIN_DBL_MAX),
2408  sumChangeObjective1_(0.0),
2409  sumChangeObjective2_(0.0),
2410  bestSolution_(NULL),
2411  currentSolution_(NULL),
2412  testSolution_(NULL),
2413  minimumDrop_(1.0e-4),
2414  numberSolutions_(0),
2415  stateOfSearch_(0),
2416  hotstartSolution_(NULL),
2417  hotstartPriorities_(NULL),
2418  numberHeuristicSolutions_(0),
2419  numberNodes_(0),
2420  numberNodes2_(0),
2421  numberIterations_(0),
2422  status_(-1),
2423  secondaryStatus_(-1),
2424  numberIntegers_(0),
2425  numberRowsAtContinuous_(0),
2426  maximumNumberCuts_(0),
2427  phase_(0),
2428  currentNumberCuts_(0),
2429  maximumDepth_(0),
2430  walkback_(NULL),
2431  addedCuts_(NULL),
2432  nextRowCut_(NULL),
2433  currentNode_(NULL),
2434  integerVariable_(NULL),
2435  integerInfo_(NULL),
2436  continuousSolution_(NULL),
2437  usedInSolution_(NULL),
2438  specialOptions_(0),
2439  subTreeModel_(NULL),
2440  numberStoppedSubTrees_(0),
2441  presolve_(0),
2442  numberStrong_(5),
2443  numberBeforeTrust_(0),
2444  numberPenalties_(20),
2445  penaltyScaleFactor_(3.0),
2446  numberAnalyzeIterations_(0),
2447  analyzeResults_(NULL),
2448  numberInfeasibleNodes_(0),
2449  problemType_(0),
2450  printFrequency_(0),
2451  numberCutGenerators_(0),
2452  generator_(NULL),
2453  virginGenerator_(NULL),
2454  numberHeuristics_(0),
2455  heuristic_(NULL),
2456  lastHeuristic_(NULL),
2457  eventHandler_(0),
2458  numberObjects_(0),
2459  object_(NULL),
2460  ownObjects_(true),
2461  originalColumns_(NULL),
2462  howOftenGlobalScan_(1),
2463  numberGlobalViolations_(0),
2464  continuousObjective_(COIN_DBL_MAX),
2465  originalContinuousObjective_(COIN_DBL_MAX),
2466  continuousInfeasibilities_(INT_MAX),
2467  maximumCutPassesAtRoot_(20),
2468  maximumCutPasses_(10),
2469  preferredWay_(0),
2470  currentPassNumber_(0),
2471  maximumWhich_(1000),
2472  whichGenerator_(NULL),
2473  maximumStatistics_(0),
2474  statistics_(NULL),
2475  maximumDepthActual_(0),
2476  numberDJFixed_(0.0),
2477  probingInfo_(NULL),
2478  numberFixedAtRoot_(0),
2479  numberFixedNow_(0),
2480  stoppedOnGap_(false),
2481  eventHappened_(false),
2482  numberLongStrong_(0),
2483  numberOldActiveCuts_(0),
2484  numberNewCuts_(0),
2485  sizeMiniTree_(0),
2486  searchStrategy_(-1),
2487  numberStrongIterations_(0),
2488  resolveAfterTakeOffCuts_(true)
2489{
2490  memset(intParam_,0,sizeof(intParam_));
2491  intParam_[CbcMaxNumNode] = 2147483647;
2492  intParam_[CbcMaxNumSol] = 9999999;
2493  intParam_[CbcFathomDiscipline] = 0;
2494
2495  dblParam_[CbcIntegerTolerance] = 1e-6;
2496  dblParam_[CbcInfeasibilityWeight] = 0.0;
2497  dblParam_[CbcCutoffIncrement] = 1e-5;
2498  dblParam_[CbcAllowableGap] = 1.0e-10;
2499  dblParam_[CbcAllowableFractionGap] = 0.0;
2500  dblParam_[CbcMaximumSeconds] = 1.0e100;
2501  dblParam_[CbcCurrentCutoff] = 1.0e100;
2502  dblParam_[CbcOptimizationDirection] = 1.0;
2503  dblParam_[CbcCurrentObjectiveValue] = 1.0e100;
2504  dblParam_[CbcCurrentMinimizationObjectiveValue] = 1.0e100;
2505  dblParam_[CbcStartSeconds] = 0.0;
2506  strongInfo_[0]=0;
2507  strongInfo_[1]=0;
2508  strongInfo_[2]=0;
2509  solverCharacteristics_ = NULL;
2510  nodeCompare_=new CbcCompareDefault();;
2511  problemFeasibility_=new CbcFeasibilityBase();
2512  tree_= new CbcTree();
2513  branchingMethod_=NULL;
2514  cutModifier_=NULL;
2515  strategy_=NULL;
2516  parentModel_=NULL;
2517  cbcColLower_ = NULL;
2518  cbcColUpper_ = NULL;
2519  cbcRowLower_ = NULL;
2520  cbcRowUpper_ = NULL;
2521  cbcColSolution_ = NULL;
2522  cbcRowPrice_ = NULL;
2523  cbcReducedCost_ = NULL;
2524  cbcRowActivity_ = NULL;
2525  appData_=NULL;
2526  handler_ = new CoinMessageHandler();
2527  handler_->setLogLevel(2);
2528  messages_ = CbcMessage();
2529  eventHandler_ = new CbcEventHandler() ;
2530}
2531
2532/** Constructor from solver.
2533
2534  Creates a model complete with a clone of the solver passed as a parameter.
2535*/
2536
2537CbcModel::CbcModel(const OsiSolverInterface &rhs)
2538:
2539  continuousSolver_(NULL),
2540  referenceSolver_(NULL),
2541  defaultHandler_(true),
2542  emptyWarmStart_(NULL),
2543  bestObjective_(COIN_DBL_MAX),
2544  bestPossibleObjective_(COIN_DBL_MAX),
2545  sumChangeObjective1_(0.0),
2546  sumChangeObjective2_(0.0),
2547  minimumDrop_(1.0e-4),
2548  numberSolutions_(0),
2549  stateOfSearch_(0),
2550  hotstartSolution_(NULL),
2551  hotstartPriorities_(NULL),
2552  numberHeuristicSolutions_(0),
2553  numberNodes_(0),
2554  numberNodes2_(0),
2555  numberIterations_(0),
2556  status_(-1),
2557  secondaryStatus_(-1),
2558  numberRowsAtContinuous_(0),
2559  maximumNumberCuts_(0),
2560  phase_(0),
2561  currentNumberCuts_(0),
2562  maximumDepth_(0),
2563  walkback_(NULL),
2564  addedCuts_(NULL),
2565  nextRowCut_(NULL),
2566  currentNode_(NULL),
2567  integerInfo_(NULL),
2568  specialOptions_(0),
2569  subTreeModel_(NULL),
2570  numberStoppedSubTrees_(0),
2571  presolve_(0),
2572  numberStrong_(5),
2573  numberBeforeTrust_(0),
2574  numberPenalties_(20),
2575  penaltyScaleFactor_(3.0),
2576  numberAnalyzeIterations_(0),
2577  analyzeResults_(NULL),
2578  numberInfeasibleNodes_(0),
2579  problemType_(0),
2580  printFrequency_(0),
2581  numberCutGenerators_(0),
2582  generator_(NULL),
2583  virginGenerator_(NULL),
2584  numberHeuristics_(0),
2585  heuristic_(NULL),
2586  lastHeuristic_(NULL),
2587  eventHandler_(0),
2588  numberObjects_(0),
2589  object_(NULL),
2590  ownObjects_(true),
2591  originalColumns_(NULL),
2592  howOftenGlobalScan_(1),
2593  numberGlobalViolations_(0),
2594  continuousObjective_(COIN_DBL_MAX),
2595  originalContinuousObjective_(COIN_DBL_MAX),
2596  continuousInfeasibilities_(INT_MAX),
2597  maximumCutPassesAtRoot_(20),
2598  maximumCutPasses_(10),
2599  preferredWay_(0),
2600  currentPassNumber_(0),
2601  maximumWhich_(1000),
2602  whichGenerator_(NULL),
2603  maximumStatistics_(0),
2604  statistics_(NULL),
2605  maximumDepthActual_(0),
2606  numberDJFixed_(0.0),
2607  probingInfo_(NULL),
2608  numberFixedAtRoot_(0),
2609  numberFixedNow_(0),
2610  stoppedOnGap_(false),
2611  eventHappened_(false),
2612  numberLongStrong_(0),
2613  numberOldActiveCuts_(0),
2614  numberNewCuts_(0),
2615  sizeMiniTree_(0),
2616  searchStrategy_(-1),
2617  numberStrongIterations_(0),
2618  resolveAfterTakeOffCuts_(true)
2619{
2620  memset(intParam_,0,sizeof(intParam_));
2621  intParam_[CbcMaxNumNode] = 2147483647;
2622  intParam_[CbcMaxNumSol] = 9999999;
2623  intParam_[CbcFathomDiscipline] = 0;
2624
2625  dblParam_[CbcIntegerTolerance] = 1e-6;
2626  dblParam_[CbcInfeasibilityWeight] = 0.0;
2627  dblParam_[CbcCutoffIncrement] = 1e-5;
2628  dblParam_[CbcAllowableGap] = 1.0e-10;
2629  dblParam_[CbcAllowableFractionGap] = 0.0;
2630  dblParam_[CbcMaximumSeconds] = 1.0e100;
2631  dblParam_[CbcCurrentCutoff] = 1.0e100;
2632  dblParam_[CbcOptimizationDirection] = 1.0;
2633  dblParam_[CbcCurrentObjectiveValue] = 1.0e100;
2634  dblParam_[CbcCurrentMinimizationObjectiveValue] = 1.0e100;
2635  dblParam_[CbcStartSeconds] = 0.0;
2636  strongInfo_[0]=0;
2637  strongInfo_[1]=0;
2638  strongInfo_[2]=0;
2639  solverCharacteristics_ = NULL;
2640
2641  nodeCompare_=new CbcCompareDefault();;
2642  problemFeasibility_=new CbcFeasibilityBase();
2643  tree_= new CbcTree();
2644  branchingMethod_=NULL;
2645  cutModifier_=NULL;
2646  strategy_=NULL;
2647  parentModel_=NULL;
2648  appData_=NULL;
2649  handler_ = new CoinMessageHandler();
2650  handler_->setLogLevel(2);
2651  messages_ = CbcMessage();
2652  eventHandler_ = new CbcEventHandler() ;
2653  solver_ = rhs.clone();
2654  referenceSolver_ = solver_->clone();
2655  ourSolver_ = true ;
2656  cbcColLower_ = NULL;
2657  cbcColUpper_ = NULL;
2658  cbcRowLower_ = NULL;
2659  cbcRowUpper_ = NULL;
2660  cbcColSolution_ = NULL;
2661  cbcRowPrice_ = NULL;
2662  cbcReducedCost_ = NULL;
2663  cbcRowActivity_ = NULL;
2664
2665  // Initialize solution and integer variable vectors
2666  bestSolution_ = NULL; // to say no solution found
2667  numberIntegers_=0;
2668  int numberColumns = solver_->getNumCols();
2669  int iColumn;
2670  if (numberColumns) {
2671    // Space for current solution
2672    currentSolution_ = new double[numberColumns];
2673    continuousSolution_ = new double[numberColumns];
2674    usedInSolution_ = new int[numberColumns];
2675    for (iColumn=0;iColumn<numberColumns;iColumn++) {
2676      if( solver_->isInteger(iColumn)) 
2677        numberIntegers_++;
2678    }
2679  } else {
2680    // empty model
2681    currentSolution_=NULL;
2682    continuousSolution_=NULL;
2683    usedInSolution_=NULL;
2684  }
2685  testSolution_=currentSolution_;
2686  if (numberIntegers_) {
2687    integerVariable_ = new int [numberIntegers_];
2688    numberIntegers_=0;
2689    for (iColumn=0;iColumn<numberColumns;iColumn++) {
2690      if( solver_->isInteger(iColumn)) 
2691        integerVariable_[numberIntegers_++]=iColumn;
2692    }
2693  } else {
2694    integerVariable_ = NULL;
2695  }
2696}
2697
2698/*
2699  Assign a solver to the model (model assumes ownership)
2700
2701  The integer variable vector is initialized if it's not already present.
2702  If deleteSolver then current solver deleted (if model owned)
2703
2704  Assuming ownership matches usage in OsiSolverInterface
2705  (cf. assignProblem, loadProblem).
2706
2707  TODO: What to do about solver parameters? A simple copy likely won't do it,
2708        because the SI must push the settings into the underlying solver. In
2709        the context of switching solvers in cbc, this means that command line
2710        settings will get lost. Stash the command line somewhere and reread it
2711        here, maybe?
2712 
2713  TODO: More generally, how much state should be transferred from the old
2714        solver to the new solver? Best perhaps to see how usage develops.
2715        What's done here mimics the CbcModel(OsiSolverInterface) constructor.
2716*/
2717void
2718CbcModel::assignSolver(OsiSolverInterface *&solver, bool deleteSolver)
2719
2720{
2721  // resize best solution if exists
2722  if (bestSolution_) {
2723    int nOld = solver_->getNumCols();
2724    int nNew = solver->getNumCols();
2725    if (nNew>nOld) {
2726      double * temp = new double[nNew];
2727      memcpy(temp,bestSolution_,nOld*sizeof(double));
2728      memset(temp+nOld,0,(nNew-nOld)*sizeof(double));
2729      delete [] bestSolution_;
2730      bestSolution_=temp;
2731    }
2732  }
2733  // Keep the current message level for solver (if solver exists)
2734  if (solver_)
2735    solver->messageHandler()->setLogLevel(solver_->messageHandler()->logLevel()) ;
2736
2737  if (ourSolver_&&deleteSolver) delete solver_ ;
2738  solver_ = solver;
2739  solver = NULL ;
2740  ourSolver_ = true ;
2741/*
2742  Basis information is solver-specific.
2743*/
2744  if (emptyWarmStart_)
2745  { delete emptyWarmStart_  ;
2746    emptyWarmStart_ = 0 ; }
2747/*
2748  Initialize integer variable vector.
2749*/
2750  numberIntegers_=0;
2751  int numberColumns = solver_->getNumCols();
2752  int iColumn;
2753  for (iColumn=0;iColumn<numberColumns;iColumn++) {
2754    if( solver_->isInteger(iColumn)) 
2755      numberIntegers_++;
2756  }
2757  if (numberIntegers_) {
2758    delete [] integerVariable_;
2759    integerVariable_ = new int [numberIntegers_];
2760    numberIntegers_=0;
2761    for (iColumn=0;iColumn<numberColumns;iColumn++) {
2762      if( solver_->isInteger(iColumn)) 
2763        integerVariable_[numberIntegers_++]=iColumn;
2764    }
2765  } else {
2766    integerVariable_ = NULL;
2767  }
2768
2769  return ;
2770}
2771
2772// Copy constructor.
2773
2774CbcModel::CbcModel(const CbcModel & rhs, bool noTree)
2775:
2776  continuousSolver_(NULL),
2777  referenceSolver_(NULL),
2778  defaultHandler_(rhs.defaultHandler_),
2779  emptyWarmStart_(NULL),
2780  bestObjective_(rhs.bestObjective_),
2781  bestPossibleObjective_(rhs.bestPossibleObjective_),
2782  sumChangeObjective1_(rhs.sumChangeObjective1_),
2783  sumChangeObjective2_(rhs.sumChangeObjective2_),
2784  minimumDrop_(rhs.minimumDrop_),
2785  numberSolutions_(rhs.numberSolutions_),
2786  stateOfSearch_(rhs.stateOfSearch_),
2787  numberHeuristicSolutions_(rhs.numberHeuristicSolutions_),
2788  numberNodes_(rhs.numberNodes_),
2789  numberNodes2_(rhs.numberNodes2_),
2790  numberIterations_(rhs.numberIterations_),
2791  status_(rhs.status_),
2792  secondaryStatus_(rhs.secondaryStatus_),
2793  specialOptions_(rhs.specialOptions_),
2794  subTreeModel_(rhs.subTreeModel_),
2795  numberStoppedSubTrees_(rhs.numberStoppedSubTrees_),
2796  presolve_(rhs.presolve_),
2797  numberStrong_(rhs.numberStrong_),
2798  numberBeforeTrust_(rhs.numberBeforeTrust_),
2799  numberPenalties_(rhs.numberPenalties_),
2800  penaltyScaleFactor_(rhs.penaltyScaleFactor_),
2801  numberAnalyzeIterations_(rhs.numberAnalyzeIterations_),
2802  analyzeResults_(NULL),
2803  numberInfeasibleNodes_(rhs.numberInfeasibleNodes_),
2804  problemType_(rhs.problemType_),
2805  printFrequency_(rhs.printFrequency_),
2806  howOftenGlobalScan_(rhs.howOftenGlobalScan_),
2807  numberGlobalViolations_(rhs.numberGlobalViolations_),
2808  continuousObjective_(rhs.continuousObjective_),
2809  originalContinuousObjective_(rhs.originalContinuousObjective_),
2810  continuousInfeasibilities_(rhs.continuousInfeasibilities_),
2811  maximumCutPassesAtRoot_(rhs.maximumCutPassesAtRoot_),
2812  maximumCutPasses_( rhs.maximumCutPasses_),
2813  preferredWay_(rhs.preferredWay_),
2814  currentPassNumber_(rhs.currentPassNumber_),
2815  maximumWhich_(rhs.maximumWhich_),
2816  whichGenerator_(NULL),
2817  maximumStatistics_(0),
2818  statistics_(NULL),
2819  maximumDepthActual_(0),
2820  numberDJFixed_(0.0),
2821  probingInfo_(NULL),
2822  numberFixedAtRoot_(rhs.numberFixedAtRoot_),
2823  numberFixedNow_(rhs.numberFixedNow_),
2824  stoppedOnGap_(rhs.stoppedOnGap_),
2825  eventHappened_(rhs.eventHappened_),
2826  numberLongStrong_(rhs.numberLongStrong_),
2827  numberOldActiveCuts_(rhs.numberOldActiveCuts_),
2828  numberNewCuts_(rhs.numberNewCuts_),
2829  sizeMiniTree_(rhs.sizeMiniTree_),
2830  searchStrategy_(rhs.searchStrategy_),
2831  numberStrongIterations_(rhs.numberStrongIterations_),
2832  resolveAfterTakeOffCuts_(rhs.resolveAfterTakeOffCuts_)
2833{
2834  memcpy(intParam_,rhs.intParam_,sizeof(intParam_));
2835  memcpy(dblParam_,rhs.dblParam_,sizeof(dblParam_));
2836  strongInfo_[0]=rhs.strongInfo_[0];
2837  strongInfo_[1]=rhs.strongInfo_[1];
2838  strongInfo_[2]=rhs.strongInfo_[2];
2839  solverCharacteristics_ = NULL;
2840  if (rhs.emptyWarmStart_) emptyWarmStart_ = rhs.emptyWarmStart_->clone() ;
2841  if (defaultHandler_) {
2842    handler_ = new CoinMessageHandler();
2843    handler_->setLogLevel(2);
2844  } else {
2845    handler_ = rhs.handler_;
2846  }
2847  messageHandler()->setLogLevel(rhs.messageHandler()->logLevel());
2848  numberCutGenerators_ = rhs.numberCutGenerators_;
2849  if (numberCutGenerators_) {
2850    generator_ = new CbcCutGenerator * [numberCutGenerators_];
2851    virginGenerator_ = new CbcCutGenerator * [numberCutGenerators_];
2852    int i;
2853    for (i=0;i<numberCutGenerators_;i++) {
2854      generator_[i]=new CbcCutGenerator(*rhs.generator_[i]);
2855      virginGenerator_[i]=new CbcCutGenerator(*rhs.virginGenerator_[i]);
2856    }
2857  } else {
2858    generator_=NULL;
2859    virginGenerator_=NULL;
2860  }
2861  if (!noTree)
2862    globalCuts_ = rhs.globalCuts_;
2863  numberHeuristics_ = rhs.numberHeuristics_;
2864  if (numberHeuristics_) {
2865    heuristic_ = new CbcHeuristic * [numberHeuristics_];
2866    int i;
2867    for (i=0;i<numberHeuristics_;i++) {
2868      heuristic_[i]=rhs.heuristic_[i]->clone();
2869    }
2870  } else {
2871    heuristic_=NULL;
2872  }
2873  lastHeuristic_ = NULL;
2874  if (rhs.eventHandler_)
2875  { eventHandler_ = new CbcEventHandler(*rhs.eventHandler_) ; }
2876  else
2877  { eventHandler_ = NULL ; }
2878  ownObjects_ = rhs.ownObjects_;
2879  if (ownObjects_) {
2880    numberObjects_=rhs.numberObjects_;
2881    if (numberObjects_) {
2882      object_ = new OsiObject * [numberObjects_];
2883      int i;
2884      for (i=0;i<numberObjects_;i++) 
2885        object_[i]=(rhs.object_[i])->clone();
2886    } else {
2887      object_=NULL;
2888    }
2889  } else {
2890    // assume will be redone
2891    numberObjects_=0;
2892    object_=NULL;
2893  }
2894  if (rhs.referenceSolver_)
2895    referenceSolver_ = rhs.referenceSolver_->clone();
2896  else
2897    referenceSolver_=NULL;
2898  if (!noTree||!rhs.continuousSolver_)
2899    solver_ = rhs.solver_->clone();
2900  else
2901    solver_ = rhs.continuousSolver_->clone();
2902  if (rhs.originalColumns_) {
2903    int numberColumns = solver_->getNumCols();
2904    originalColumns_= new int [numberColumns];
2905    memcpy(originalColumns_,rhs.originalColumns_,numberColumns*sizeof(int));
2906  } else {
2907    originalColumns_=NULL;
2908  }
2909  if (maximumWhich_&&rhs.whichGenerator_)
2910    whichGenerator_ = CoinCopyOfArray(rhs.whichGenerator_,maximumWhich_);
2911  nodeCompare_=rhs.nodeCompare_->clone();
2912  problemFeasibility_=rhs.problemFeasibility_->clone();
2913  tree_= rhs.tree_->clone();
2914  if (rhs.branchingMethod_)
2915    branchingMethod_=rhs.branchingMethod_->clone();
2916  else
2917    branchingMethod_=NULL;
2918  if (rhs.cutModifier_)
2919    cutModifier_=rhs.cutModifier_->clone();
2920  else
2921    cutModifier_=NULL;
2922  cbcColLower_ = NULL;
2923  cbcColUpper_ = NULL;
2924  cbcRowLower_ = NULL;
2925  cbcRowUpper_ = NULL;
2926  cbcColSolution_ = NULL;
2927  cbcRowPrice_ = NULL;
2928  cbcReducedCost_ = NULL;
2929  cbcRowActivity_ = NULL;
2930  if (rhs.strategy_)
2931    strategy_=rhs.strategy_->clone();
2932  else
2933    strategy_=NULL;
2934  parentModel_=rhs.parentModel_;
2935  appData_=rhs.appData_;
2936  messages_ = rhs.messages_;
2937  ourSolver_ = true ;
2938  messageHandler()->setLogLevel(rhs.messageHandler()->logLevel());
2939  numberIntegers_=rhs.numberIntegers_;
2940  if (numberIntegers_) {
2941    integerVariable_ = new int [numberIntegers_];
2942    memcpy(integerVariable_,rhs.integerVariable_,numberIntegers_*sizeof(int));
2943    integerInfo_ = CoinCopyOfArray(rhs.integerInfo_,solver_->getNumCols());
2944  } else {
2945    integerVariable_ = NULL;
2946    integerInfo_=NULL;
2947  }
2948  if (rhs.hotstartSolution_) {
2949    int numberColumns = solver_->getNumCols();
2950    hotstartSolution_ = CoinCopyOfArray(rhs.hotstartSolution_,numberColumns);
2951    hotstartPriorities_ = CoinCopyOfArray(rhs.hotstartPriorities_,numberColumns);
2952  } else {
2953    hotstartSolution_ = NULL;
2954    hotstartPriorities_ =NULL;
2955  }
2956  if (rhs.bestSolution_&&!noTree) {
2957    int numberColumns = solver_->getNumCols();
2958    bestSolution_ = new double[numberColumns];
2959    memcpy(bestSolution_,rhs.bestSolution_,numberColumns*sizeof(double));
2960  } else {
2961    bestSolution_=NULL;
2962  }
2963  if (!noTree) {
2964    int numberColumns = solver_->getNumCols();
2965    currentSolution_ = CoinCopyOfArray(rhs.currentSolution_,numberColumns);
2966    continuousSolution_ = CoinCopyOfArray(rhs.continuousSolution_,numberColumns);
2967    usedInSolution_ = CoinCopyOfArray(rhs.usedInSolution_,numberColumns);
2968  } else {
2969    currentSolution_=NULL;
2970    continuousSolution_=NULL;
2971    usedInSolution_=NULL;
2972  }
2973  testSolution_=currentSolution_;
2974  numberRowsAtContinuous_ = rhs.numberRowsAtContinuous_;
2975  maximumNumberCuts_=rhs.maximumNumberCuts_;
2976  phase_ = rhs.phase_;
2977  currentNumberCuts_=rhs.currentNumberCuts_;
2978  maximumDepth_= rhs.maximumDepth_;
2979  if (noTree) {
2980    bestObjective_ = COIN_DBL_MAX;
2981    numberSolutions_ =0;
2982    stateOfSearch_= 0;
2983    numberHeuristicSolutions_=0;
2984    numberNodes_=0;
2985    numberNodes2_=0;
2986    numberIterations_=0;
2987    status_=0;
2988    subTreeModel_=NULL;
2989    numberStoppedSubTrees_=0;
2990    continuousObjective_=COIN_DBL_MAX;
2991    originalContinuousObjective_=COIN_DBL_MAX;
2992    continuousInfeasibilities_=INT_MAX;
2993    maximumNumberCuts_=0;
2994    tree_->cleanTree(this,-COIN_DBL_MAX,bestPossibleObjective_);
2995    bestPossibleObjective_ = COIN_DBL_MAX;
2996  }
2997  // These are only used as temporary arrays so need not be filled
2998  if (maximumNumberCuts_) {
2999    addedCuts_ = new CbcCountRowCut * [maximumNumberCuts_];
3000  } else {
3001    addedCuts_ = NULL;
3002  }
3003  nextRowCut_ = NULL;
3004  currentNode_ = NULL;
3005  if (maximumDepth_)
3006    walkback_ = new CbcNodeInfo * [maximumDepth_];
3007  else
3008    walkback_ = NULL;
3009  synchronizeModel();
3010}
3011 
3012// Assignment operator
3013CbcModel & 
3014CbcModel::operator=(const CbcModel& rhs)
3015{
3016  if (this!=&rhs) {
3017    gutsOfDestructor();
3018    if (defaultHandler_)
3019    { delete handler_;
3020      handler_ = NULL; }
3021    defaultHandler_ = rhs.defaultHandler_;
3022    if (defaultHandler_)
3023    { handler_ = new CoinMessageHandler();
3024      handler_->setLogLevel(2); }
3025    else
3026    { handler_ = rhs.handler_; }
3027    messages_ = rhs.messages_;
3028    messageHandler()->setLogLevel(rhs.messageHandler()->logLevel());
3029
3030    delete solver_;
3031    if (rhs.solver_)
3032    { solver_ = rhs.solver_->clone() ; }
3033    else
3034    { solver_ = 0 ; }
3035    ourSolver_ = true ;
3036    delete continuousSolver_ ;
3037    if (rhs.continuousSolver_)
3038    { solver_ = rhs.continuousSolver_->clone() ; }
3039    else
3040    { continuousSolver_ = 0 ; }
3041
3042    if (rhs.referenceSolver_)
3043    { solver_ = rhs.referenceSolver_->clone() ; }
3044    else
3045    { referenceSolver_ = NULL ; }
3046
3047    delete emptyWarmStart_ ;
3048    if (rhs.emptyWarmStart_)
3049    { emptyWarmStart_ = rhs.emptyWarmStart_->clone() ; }
3050    else
3051    { emptyWarmStart_ = 0 ; }
3052
3053    bestObjective_ = rhs.bestObjective_;
3054    bestPossibleObjective_=rhs.bestPossibleObjective_;
3055    sumChangeObjective1_=rhs.sumChangeObjective1_;
3056    sumChangeObjective2_=rhs.sumChangeObjective2_;
3057    delete [] bestSolution_;
3058    if (rhs.bestSolution_) {
3059      int numberColumns = rhs.getNumCols();
3060      bestSolution_ = new double[numberColumns];
3061      memcpy(bestSolution_,rhs.bestSolution_,numberColumns*sizeof(double));
3062    } else {
3063      bestSolution_=NULL;
3064    }
3065    int numberColumns = solver_->getNumCols();
3066    currentSolution_ = CoinCopyOfArray(rhs.currentSolution_,numberColumns);
3067    continuousSolution_ = CoinCopyOfArray(rhs.continuousSolution_,numberColumns);
3068    usedInSolution_ = CoinCopyOfArray(rhs.usedInSolution_,numberColumns);
3069    testSolution_=currentSolution_;
3070    minimumDrop_ = rhs.minimumDrop_;
3071    numberSolutions_=rhs.numberSolutions_;
3072    stateOfSearch_= rhs.stateOfSearch_;
3073    numberHeuristicSolutions_=rhs.numberHeuristicSolutions_;
3074    numberNodes_ = rhs.numberNodes_;
3075    numberNodes2_ = rhs.numberNodes2_;
3076    numberIterations_ = rhs.numberIterations_;
3077    status_ = rhs.status_;
3078    secondaryStatus_ = rhs.secondaryStatus_;
3079    specialOptions_ = rhs.specialOptions_;
3080    subTreeModel_ = rhs.subTreeModel_;
3081    numberStoppedSubTrees_ = rhs.numberStoppedSubTrees_;
3082    presolve_ = rhs.presolve_;
3083    numberStrong_ = rhs.numberStrong_;
3084    numberBeforeTrust_ = rhs.numberBeforeTrust_;
3085    numberPenalties_ = rhs.numberPenalties_;
3086    penaltyScaleFactor_ = rhs.penaltyScaleFactor_;
3087    numberAnalyzeIterations_ = rhs.numberAnalyzeIterations_;
3088    delete [] analyzeResults_;
3089    analyzeResults_ = NULL;
3090    numberInfeasibleNodes_ = rhs.numberInfeasibleNodes_;
3091    problemType_ = rhs.problemType_;
3092    printFrequency_ = rhs.printFrequency_;
3093    howOftenGlobalScan_=rhs.howOftenGlobalScan_;
3094    numberGlobalViolations_=rhs.numberGlobalViolations_;
3095    continuousObjective_=rhs.continuousObjective_;
3096    originalContinuousObjective_ = rhs.originalContinuousObjective_;
3097    continuousInfeasibilities_ = rhs.continuousInfeasibilities_;
3098    maximumCutPassesAtRoot_ = rhs.maximumCutPassesAtRoot_;
3099    maximumCutPasses_ = rhs.maximumCutPasses_;
3100    preferredWay_ = rhs.preferredWay_;
3101    currentPassNumber_ = rhs.currentPassNumber_;
3102    memcpy(intParam_,rhs.intParam_,sizeof(intParam_));
3103    memcpy(dblParam_,rhs.dblParam_,sizeof(dblParam_));
3104    globalCuts_ = rhs.globalCuts_;
3105    int i;
3106    for (i=0;i<numberCutGenerators_;i++) {
3107      delete generator_[i];
3108      delete virginGenerator_[i];
3109    }
3110    delete [] generator_;
3111    delete [] virginGenerator_;
3112    delete [] heuristic_;
3113    maximumWhich_ = rhs.maximumWhich_;
3114    delete [] whichGenerator_;
3115    whichGenerator_ = NULL;
3116    if (maximumWhich_&&rhs.whichGenerator_)
3117      whichGenerator_ = CoinCopyOfArray(rhs.whichGenerator_,maximumWhich_);
3118    for (i=0;i<maximumStatistics_;i++)
3119      delete statistics_[i];
3120    delete [] statistics_;
3121    maximumStatistics_ = 0;
3122    statistics_ = NULL;
3123    delete probingInfo_;
3124    probingInfo_=NULL;
3125    numberFixedAtRoot_ = rhs.numberFixedAtRoot_;
3126    numberFixedNow_ = rhs.numberFixedNow_;
3127    stoppedOnGap_ = rhs.stoppedOnGap_;
3128    eventHappened_ = rhs.eventHappened_;
3129    numberLongStrong_ = rhs.numberLongStrong_;
3130    numberOldActiveCuts_ = rhs.numberOldActiveCuts_;
3131    numberNewCuts_ = rhs.numberNewCuts_;
3132    resolveAfterTakeOffCuts_=rhs.resolveAfterTakeOffCuts_;
3133    sizeMiniTree_ = rhs.sizeMiniTree_;
3134    searchStrategy_ = rhs.searchStrategy_;
3135    numberStrongIterations_ = rhs.numberStrongIterations_;
3136    strongInfo_[0]=rhs.strongInfo_[0];
3137    strongInfo_[1]=rhs.strongInfo_[1];
3138    strongInfo_[2]=rhs.strongInfo_[2];
3139    solverCharacteristics_ = NULL;
3140    lastHeuristic_ = NULL;
3141    numberCutGenerators_ = rhs.numberCutGenerators_;
3142    if (numberCutGenerators_) {
3143      generator_ = new CbcCutGenerator * [numberCutGenerators_];
3144      virginGenerator_ = new CbcCutGenerator * [numberCutGenerators_];
3145      int i;
3146      for (i=0;i<numberCutGenerators_;i++) {
3147        generator_[i]=new CbcCutGenerator(*rhs.generator_[i]);
3148        virginGenerator_[i]=new CbcCutGenerator(*rhs.virginGenerator_[i]);
3149      }
3150    } else {
3151      generator_=NULL;
3152      virginGenerator_=NULL;
3153    }
3154    numberHeuristics_ = rhs.numberHeuristics_;
3155    if (numberHeuristics_) {
3156      heuristic_ = new CbcHeuristic * [numberHeuristics_];
3157      memcpy(heuristic_,rhs.heuristic_,
3158             numberHeuristics_*sizeof(CbcHeuristic *));
3159    } else {
3160      heuristic_=NULL;
3161    }
3162    lastHeuristic_ = NULL;
3163    if (eventHandler_)
3164      delete eventHandler_ ;
3165    if (rhs.eventHandler_)
3166    { eventHandler_ = new CbcEventHandler(*rhs.eventHandler_) ; }
3167    else
3168    { eventHandler_ = NULL ; }
3169    if (ownObjects_) {
3170      for (i=0;i<numberObjects_;i++)
3171        delete object_[i];
3172      delete [] object_;
3173      numberObjects_=rhs.numberObjects_;
3174      if (numberObjects_) {
3175        object_ = new OsiObject * [numberObjects_];
3176        int i;
3177        for (i=0;i<numberObjects_;i++) 
3178          object_[i]=(rhs.object_[i])->clone();
3179      } else {
3180        object_=NULL;
3181    }
3182    } else {
3183      // assume will be redone
3184      numberObjects_=0;
3185      object_=NULL;
3186    }
3187    delete [] originalColumns_;
3188    if (rhs.originalColumns_) {
3189      int numberColumns = rhs.getNumCols();
3190      originalColumns_= new int [numberColumns];
3191      memcpy(originalColumns_,rhs.originalColumns_,numberColumns*sizeof(int));
3192    } else {
3193      originalColumns_=NULL;
3194    }
3195    nodeCompare_=rhs.nodeCompare_->clone();
3196    problemFeasibility_=rhs.problemFeasibility_->clone();
3197    delete tree_;
3198    tree_= rhs.tree_->clone();
3199    if (rhs.branchingMethod_)
3200      branchingMethod_=rhs.branchingMethod_->clone();
3201    else
3202      branchingMethod_=NULL;
3203    if (rhs.cutModifier_)
3204      cutModifier_=rhs.cutModifier_->clone();
3205    else
3206      cutModifier_=NULL;
3207    delete strategy_;
3208    if (rhs.strategy_)
3209      strategy_=rhs.strategy_->clone();
3210    else
3211      strategy_=NULL;
3212    parentModel_=rhs.parentModel_;
3213    appData_=rhs.appData_;
3214
3215    delete [] integerVariable_;
3216    numberIntegers_=rhs.numberIntegers_;
3217    if (numberIntegers_) {
3218      integerVariable_ = new int [numberIntegers_];
3219      memcpy(integerVariable_,rhs.integerVariable_,
3220             numberIntegers_*sizeof(int));
3221      integerInfo_ = CoinCopyOfArray(rhs.integerInfo_,solver_->getNumCols());
3222    } else {
3223      integerVariable_ = NULL;
3224      integerInfo_=NULL;
3225    }
3226    if (rhs.hotstartSolution_) {
3227      int numberColumns = solver_->getNumCols();
3228      hotstartSolution_ = CoinCopyOfArray(rhs.hotstartSolution_,numberColumns);
3229      hotstartPriorities_ = CoinCopyOfArray(rhs.hotstartPriorities_,numberColumns);
3230    } else {
3231      hotstartSolution_ = NULL;
3232      hotstartPriorities_ =NULL;
3233    }
3234    numberRowsAtContinuous_ = rhs.numberRowsAtContinuous_;
3235    maximumNumberCuts_=rhs.maximumNumberCuts_;
3236    phase_ = rhs.phase_;
3237    currentNumberCuts_=rhs.currentNumberCuts_;
3238    maximumDepth_= rhs.maximumDepth_;
3239    delete [] addedCuts_;
3240    delete [] walkback_;
3241    // These are only used as temporary arrays so need not be filled
3242    if (maximumNumberCuts_) {
3243      addedCuts_ = new CbcCountRowCut * [maximumNumberCuts_];
3244    } else {
3245      addedCuts_ = NULL;
3246    }
3247    nextRowCut_ = NULL;
3248    currentNode_ = NULL;
3249    if (maximumDepth_)
3250      walkback_ = new CbcNodeInfo * [maximumDepth_];
3251    else
3252      walkback_ = NULL;
3253    synchronizeModel();
3254    cbcColLower_ = NULL;
3255    cbcColUpper_ = NULL;
3256    cbcRowLower_ = NULL;
3257    cbcRowUpper_ = NULL;
3258    cbcColSolution_ = NULL;
3259    cbcRowPrice_ = NULL;
3260    cbcReducedCost_ = NULL;
3261    cbcRowActivity_ = NULL;
3262  }
3263  return *this;
3264}
3265 
3266// Destructor
3267CbcModel::~CbcModel ()
3268{
3269  if (defaultHandler_) {
3270    delete handler_;
3271    handler_ = NULL;
3272  }
3273  delete tree_;
3274  tree_=NULL;
3275  if (ourSolver_) delete solver_;
3276  gutsOfDestructor();
3277  delete eventHandler_ ;
3278  eventHandler_ = NULL ;
3279}
3280// Clears out as much as possible (except solver)
3281void 
3282CbcModel::gutsOfDestructor()
3283{
3284  delete referenceSolver_;
3285  referenceSolver_=NULL;
3286  int i;
3287  for (i=0;i<numberCutGenerators_;i++) {
3288    delete generator_[i];
3289    delete virginGenerator_[i];
3290  }
3291  delete [] generator_;
3292  delete [] virginGenerator_;
3293  generator_=NULL;
3294  virginGenerator_=NULL;
3295  for (i=0;i<numberHeuristics_;i++)
3296    delete heuristic_[i];
3297  delete [] heuristic_;
3298  heuristic_=NULL;
3299  delete nodeCompare_;
3300  nodeCompare_=NULL;
3301  delete problemFeasibility_;
3302  problemFeasibility_=NULL;
3303  delete [] originalColumns_;
3304  originalColumns_=NULL;
3305  delete strategy_;
3306  gutsOfDestructor2();
3307}
3308// Clears out enough to reset CbcModel
3309void 
3310CbcModel::gutsOfDestructor2()
3311{
3312  delete [] integerInfo_;
3313  integerInfo_=NULL;
3314  delete [] integerVariable_;
3315  integerVariable_=NULL;
3316  int i;
3317  if (ownObjects_) {
3318    for (i=0;i<numberObjects_;i++)
3319      delete object_[i];
3320    delete [] object_;
3321  }
3322  ownObjects_=true;
3323  object_=NULL;
3324  numberIntegers_=0;
3325  numberObjects_=0;
3326  delete emptyWarmStart_ ;
3327  emptyWarmStart_ =NULL;
3328  delete continuousSolver_;
3329  continuousSolver_=NULL;
3330  delete [] bestSolution_;
3331  bestSolution_=NULL;
3332  delete [] currentSolution_;
3333  currentSolution_=NULL;
3334  delete [] continuousSolution_;
3335  continuousSolution_=NULL;
3336  solverCharacteristics_=NULL;
3337  delete [] usedInSolution_;
3338  usedInSolution_ = NULL;
3339  testSolution_=NULL;
3340  lastHeuristic_ = NULL;
3341  delete [] addedCuts_;
3342  addedCuts_=NULL;
3343  nextRowCut_ = NULL;
3344  currentNode_ = NULL;
3345  delete [] walkback_;
3346  walkback_=NULL;
3347  delete [] whichGenerator_;
3348  whichGenerator_ = NULL;
3349  for (i=0;i<maximumStatistics_;i++)
3350    delete statistics_[i];
3351  delete [] statistics_;
3352  statistics_=NULL;
3353  maximumDepthActual_ = 0;
3354  numberDJFixed_ =0.0;
3355  delete probingInfo_;
3356  probingInfo_ = NULL;
3357  maximumStatistics_=0;
3358  delete [] analyzeResults_;
3359  analyzeResults_=NULL;
3360  // Below here is whatever consensus is
3361  ourSolver_=true;
3362  bestObjective_=COIN_DBL_MAX;
3363  bestPossibleObjective_=COIN_DBL_MAX;
3364  sumChangeObjective1_=0.0;
3365  sumChangeObjective2_=0.0;
3366  numberSolutions_=0;
3367  stateOfSearch_=0;
3368  delete [] hotstartSolution_;
3369  hotstartSolution_=NULL;
3370  delete [] hotstartPriorities_;
3371  hotstartPriorities_=NULL;
3372  numberHeuristicSolutions_=0;
3373  numberNodes_=0;
3374  numberNodes2_=0;
3375  numberIterations_=0;
3376  status_=-1;
3377  secondaryStatus_=-1;
3378  maximumNumberCuts_=0;
3379  phase_=0;
3380  currentNumberCuts_=0;
3381  maximumDepth_=0;
3382  nextRowCut_=NULL;
3383  currentNode_=NULL;
3384  delete branchingMethod_;
3385  branchingMethod_=NULL;
3386  delete cutModifier_;
3387  cutModifier_=NULL;
3388  // clear out tree
3389  if (tree_&&tree_->size())
3390    tree_->cleanTree(this, -1.0e100,bestPossibleObjective_) ;
3391  subTreeModel_=NULL;
3392  numberStoppedSubTrees_=0;
3393  numberInfeasibleNodes_=0;
3394  numberGlobalViolations_=0;
3395  continuousObjective_=0.0;
3396  originalContinuousObjective_=0.0;
3397  continuousInfeasibilities_=0;
3398  numberFixedAtRoot_=0;
3399  numberFixedNow_=0;
3400  stoppedOnGap_=false;
3401  eventHappened_=false;
3402  numberLongStrong_=0;
3403  numberOldActiveCuts_=0;
3404  numberNewCuts_=0;
3405  searchStrategy_=-1;
3406  numberStrongIterations_=0;
3407  // Parameters which need to be reset
3408  dblParam_[CbcCutoffIncrement] = 1e-5;
3409  dblParam_[CbcCurrentCutoff] = 1.0e100;
3410  dblParam_[CbcCurrentObjectiveValue] = 1.0e100;
3411  dblParam_[CbcCurrentMinimizationObjectiveValue] = 1.0e100;
3412}
3413// Save a copy of the current solver so can be reset to
3414void 
3415CbcModel::saveReferenceSolver()
3416{
3417  delete referenceSolver_;
3418  referenceSolver_= solver_->clone();
3419}
3420
3421// Uses a copy of reference solver to be current solver
3422void 
3423CbcModel::resetToReferenceSolver()
3424{
3425  delete solver_;
3426  solver_ = referenceSolver_->clone();
3427  // clear many things
3428  gutsOfDestructor2();
3429  // Reset cutoff
3430  // Solvers know about direction
3431  double direction = solver_->getObjSense();
3432  double value;
3433  solver_->getDblParam(OsiDualObjectiveLimit,value); 
3434  setCutoff(value*direction);
3435}
3436
3437// Are there a numerical difficulties?
3438bool 
3439CbcModel::isAbandoned() const
3440{
3441  return status_ == 2;
3442}
3443// Is optimality proven?
3444bool 
3445CbcModel::isProvenOptimal() const
3446{
3447  if (!status_ && bestObjective_<1.0e30)
3448    return true;
3449  else
3450    return false;
3451}
3452// Is  infeasiblity proven (or none better than cutoff)?
3453bool 
3454CbcModel::isProvenInfeasible() const
3455{
3456  if (!status_ && bestObjective_>=1.0e30)
3457    return true;
3458  else
3459    return false;
3460}
3461// Was continuous solution unbounded
3462bool 
3463CbcModel::isContinuousUnbounded() const
3464{
3465  if (!status_ && secondaryStatus_==7)
3466    return true;
3467  else
3468    return false;
3469}
3470// Was continuous solution unbounded
3471bool 
3472CbcModel::isProvenDualInfeasible() const
3473{
3474  if (!status_ && secondaryStatus_==7)
3475    return true;
3476  else
3477    return false;
3478}
3479// Node limit reached?
3480bool 
3481CbcModel::isNodeLimitReached() const
3482{
3483  return numberNodes_ >= intParam_[CbcMaxNumNode];
3484}
3485// Time limit reached?
3486bool 
3487CbcModel::isSecondsLimitReached() const
3488{
3489  if (status_==1&&secondaryStatus_==4)
3490    return true;
3491  else
3492    return false;
3493}
3494// Solution limit reached?
3495bool 
3496CbcModel::isSolutionLimitReached() const
3497{
3498  return numberSolutions_ >= intParam_[CbcMaxNumSol];
3499}
3500// Set language
3501void 
3502CbcModel::newLanguage(CoinMessages::Language language)
3503{
3504  messages_ = CbcMessage(language);
3505}
3506void 
3507CbcModel::setNumberStrong(int number)
3508{
3509  if (number<0)
3510    numberStrong_=0;
3511   else
3512    numberStrong_=number;
3513}
3514void 
3515CbcModel::setNumberBeforeTrust(int number)
3516{
3517  if (number<-1) {
3518    numberBeforeTrust_=0;
3519  } else {
3520    numberBeforeTrust_=number;
3521    //numberStrong_ = CoinMax(numberStrong_,1);
3522  }
3523}
3524void 
3525CbcModel::setNumberPenalties(int number)
3526{
3527  if (number<=0) {
3528    numberPenalties_=0;
3529  } else {
3530    numberPenalties_=number;
3531  }
3532}
3533void 
3534CbcModel::setPenaltyScaleFactor(double value)
3535{
3536  if (value<=0) {
3537    penaltyScaleFactor_=3.0;
3538  } else {
3539    penaltyScaleFactor_=value;
3540  }
3541}
3542void 
3543CbcModel::setHowOftenGlobalScan(int number)
3544{
3545  if (number<-1)
3546    howOftenGlobalScan_=0;
3547   else
3548    howOftenGlobalScan_=number;
3549}
3550
3551// Add one generator
3552void 
3553CbcModel::addCutGenerator(CglCutGenerator * generator,
3554                          int howOften, const char * name,
3555                          bool normal, bool atSolution,
3556                          bool whenInfeasible,int howOftenInSub,
3557                          int whatDepth, int whatDepthInSub)
3558{
3559  CbcCutGenerator ** temp = generator_;
3560  generator_ = new CbcCutGenerator * [numberCutGenerators_+1];
3561  memcpy(generator_,temp,numberCutGenerators_*sizeof(CbcCutGenerator *));
3562  delete[] temp ;
3563  generator_[numberCutGenerators_]= 
3564    new CbcCutGenerator(this,generator, howOften, name,
3565                        normal,atSolution,whenInfeasible,howOftenInSub,
3566                        whatDepth, whatDepthInSub);
3567  // and before any cahnges
3568  temp = virginGenerator_;
3569  virginGenerator_ = new CbcCutGenerator * [numberCutGenerators_+1];
3570  memcpy(virginGenerator_,temp,numberCutGenerators_*sizeof(CbcCutGenerator *));
3571  delete[] temp ;
3572  virginGenerator_[numberCutGenerators_++]= 
3573    new CbcCutGenerator(this,generator, howOften, name,
3574                        normal,atSolution,whenInfeasible,howOftenInSub,
3575                        whatDepth, whatDepthInSub);
3576                                                         
3577}
3578// Add one heuristic
3579void 
3580CbcModel::addHeuristic(CbcHeuristic * generator)
3581{
3582  CbcHeuristic ** temp = heuristic_;
3583  heuristic_ = new CbcHeuristic * [numberHeuristics_+1];
3584  memcpy(heuristic_,temp,numberHeuristics_*sizeof(CbcHeuristic *));
3585  delete [] temp;
3586  heuristic_[numberHeuristics_++]=generator->clone();
3587}
3588
3589/*
3590  The last subproblem handled by the solver is not necessarily related to the
3591  one being recreated, so the first action is to remove all cuts from the
3592  constraint system.  Next, traverse the tree from node to the root to
3593  determine the basis size required for this subproblem and create an empty
3594  basis with the right capacity.  Finally, traverse the tree from root to
3595  node, adjusting bounds in the constraint system, adjusting the basis, and
3596  collecting the cuts that must be added to the constraint system.
3597  applyToModel does the heavy lifting.
3598
3599  addCuts1 is used in contexts where all that's desired is the list of cuts:
3600  the node is already fathomed, and we're collecting cuts so that we can
3601  adjust reference counts as we prune nodes. Arguably the two functions
3602  should be separated. The culprit is applyToModel, which performs cut
3603  collection and model adjustment.
3604
3605  Certainly in the contexts where all we need is a list of cuts, there's no
3606  point in passing in a valid basis --- an empty basis will do just fine.
3607*/
3608void CbcModel::addCuts1 (CbcNode * node, CoinWarmStartBasis *&lastws)
3609{ int i;
3610  int nNode=0;
3611  int numberColumns = getNumCols();
3612  CbcNodeInfo * nodeInfo = node->nodeInfo();
3613
3614/*
3615  Remove all cuts from the constraint system.
3616  (original comment includes ``see note below for later efficiency'', but
3617  the reference isn't clear to me).
3618*/
3619  int currentNumberCuts = solver_->getNumRows()-numberRowsAtContinuous_;
3620  int *which = new int[currentNumberCuts];
3621  for (i = 0 ; i < currentNumberCuts ; i++)
3622    which[i] = i+numberRowsAtContinuous_;
3623  solver_->deleteRows(currentNumberCuts,which);
3624  delete [] which;
3625/*
3626  Accumulate the path from node to the root in walkback_, and accumulate a
3627  cut count in currentNumberCuts.
3628
3629  original comment: when working then just unwind until where new node joins
3630  old node (for cuts?)
3631*/
3632  currentNumberCuts = 0;
3633  while (nodeInfo) {
3634    //printf("nNode = %d, nodeInfo = %x\n",nNode,nodeInfo);
3635    walkback_[nNode++]=nodeInfo;
3636    currentNumberCuts += nodeInfo->numberCuts() ;
3637    nodeInfo = nodeInfo->parent() ;
3638    if (nNode==maximumDepth_) {
3639      maximumDepth_ *= 2;
3640      CbcNodeInfo ** temp = new CbcNodeInfo * [maximumDepth_];
3641      for (i=0;i<nNode;i++) 
3642        temp[i] = walkback_[i];
3643      delete [] walkback_;
3644      walkback_ = temp;
3645    }
3646  }
3647/*
3648  Create an empty basis with sufficient capacity for the constraint system
3649  we'll construct: original system plus cuts. Make sure we have capacity to
3650  record those cuts in addedCuts_.
3651
3652  The method of adjusting the basis at a FullNodeInfo object (the root, for
3653  example) is to use a copy constructor to duplicate the basis held in the
3654  nodeInfo, then resize it and return the new basis object. Guaranteed,
3655  lastws will point to a different basis when it returns. We pass in a basis
3656  because we need the parameter to return the allocated basis, and it's an
3657  easy way to pass in the size. But we take a hit for memory allocation.
3658*/
3659  currentNumberCuts_=currentNumberCuts;
3660  if (currentNumberCuts >= maximumNumberCuts_) {
3661    maximumNumberCuts_ = currentNumberCuts;
3662    delete [] addedCuts_;
3663    addedCuts_ = new CbcCountRowCut * [maximumNumberCuts_];
3664  }
3665  lastws->setSize(numberColumns,numberRowsAtContinuous_+currentNumberCuts);
3666/*
3667  This last bit of code traverses the path collected in walkback_ from the
3668  root back to node. At the end of the loop,
3669   * lastws will be an appropriate basis for node;
3670   * variable bounds in the constraint system will be set to be correct for
3671     node; and
3672   * addedCuts_ will be set to a list of cuts that need to be added to the
3673     constraint system at node.
3674  applyToModel does all the heavy lifting.
3675*/
3676  currentNumberCuts=0;
3677  while (nNode) {
3678    --nNode;
3679    walkback_[nNode]->applyToModel(this,lastws,addedCuts_,currentNumberCuts);
3680  }
3681}
3682
3683/*
3684  adjustCuts might be a better name: If the node is feasible, we sift through
3685  the cuts collected by addCuts1, add the ones that are tight and omit the
3686  ones that are loose. If the node is infeasible, we just adjust the
3687  reference counts to reflect that we're about to prune this node and its
3688  descendants.
3689*/
3690int CbcModel::addCuts (CbcNode *node, CoinWarmStartBasis *&lastws,bool canFix)
3691{
3692/*
3693  addCuts1 performs step 1 of restoring the subproblem at this node; see the
3694  comments there.
3695*/
3696  addCuts1(node,lastws);
3697  int i;
3698  int numberColumns = getNumCols();
3699  CbcNodeInfo * nodeInfo = node->nodeInfo();
3700  double cutoff = getCutoff() ;
3701  int currentNumberCuts=currentNumberCuts_;
3702  if (canFix) {
3703    bool feasible=true;
3704    const double *lower = solver_->getColLower() ;
3705    const double *upper = solver_->getColUpper() ;
3706    double * newLower = analyzeResults_;
3707    double * objLower = newLower+numberIntegers_;
3708    double * newUpper = objLower+numberIntegers_;
3709    double * objUpper = newUpper+numberIntegers_;
3710    int n=0;
3711    for (i=0;i<numberIntegers_;i++) {
3712      int iColumn = integerVariable_[i];
3713      bool changed=false;
3714      double lo = 0.0;
3715      double up = 0.0;
3716      if (objLower[i]>cutoff) {
3717        lo = lower[iColumn];
3718        up = upper[iColumn];
3719        if (lo<newLower[i]) {
3720          lo = newLower[i];
3721          solver_->setColLower(iColumn,lo);
3722          changed=true;
3723          n++;
3724        }
3725        if (objUpper[i]>cutoff) {
3726          if (up>newUpper[i]) {
3727            up = newUpper[i];
3728            solver_->setColUpper(iColumn,up);
3729            changed=true;
3730            n++;
3731          }
3732        }
3733      } else if (objUpper[i]>cutoff) {
3734        lo = lower[iColumn];
3735        up = upper[iColumn];
3736        if (up>newUpper[i]) {
3737          up = newUpper[i];
3738          solver_->setColUpper(iColumn,up);
3739          changed=true;
3740          n++;
3741        }
3742      }
3743      if (changed&&lo>up) {
3744        feasible=false;
3745        break;
3746      }
3747    }
3748    if (!feasible) {
3749      printf("analysis says node infeas\n");
3750      cutoff=-COIN_DBL_MAX;
3751    }
3752  }
3753/*
3754  If the node can't be fathomed by bound, reinstall tight cuts in the
3755  constraint system. Even if there are no cuts, we'll want to set the
3756  reconstructed basis in the solver.
3757*/
3758  if (node->objectiveValue() < cutoff)
3759  { 
3760#   ifdef CBC_CHECK_BASIS
3761    printf("addCuts: expanded basis; rows %d+%d\n",
3762           numberRowsAtContinuous_,currentNumberCuts);
3763    lastws->print();
3764#   endif
3765/*
3766  Adjust the basis and constraint system so that we retain only active cuts.
3767  There are three steps:
3768    1) Scan the basis. Sort the cuts into effective cuts to be kept and
3769       loose cuts to be dropped.
3770    2) Drop the loose cuts and resize the basis to fit.
3771    3) Install the tight cuts in the constraint system (applyRowCuts) and
3772       and install the basis (setWarmStart).
3773  Use of compressRows conveys we're compressing the basis and not just
3774  tweaking the artificialStatus_ array.
3775*/
3776    if (currentNumberCuts > 0) {
3777      int numberToAdd = 0;
3778      const OsiRowCut **addCuts;
3779      int numberToDrop = 0 ;
3780      int *cutsToDrop ;
3781      addCuts = new const OsiRowCut* [currentNumberCuts];
3782      cutsToDrop = new int[currentNumberCuts] ;
3783      assert (currentNumberCuts+numberRowsAtContinuous_<=lastws->getNumArtificial());
3784      for (i=0;i<currentNumberCuts;i++) {
3785        CoinWarmStartBasis::Status status = 
3786          lastws->getArtifStatus(i+numberRowsAtContinuous_);
3787        if (addedCuts_[i] &&
3788            (status != CoinWarmStartBasis::basic ||
3789             addedCuts_[i]->effectiveness()==COIN_DBL_MAX)) {
3790#         ifdef CHECK_CUT_COUNTS
3791          printf("Using cut %d %x as row %d\n",i,addedCuts_[i],
3792                 numberRowsAtContinuous_+numberToAdd);
3793#         endif
3794          addCuts[numberToAdd++] = new OsiRowCut(*addedCuts_[i]);
3795        } else {
3796#         ifdef CHECK_CUT_COUNTS
3797          printf("Dropping cut %d %x\n",i,addedCuts_[i]);
3798#         endif
3799          addedCuts_[i]=NULL;
3800          cutsToDrop[numberToDrop++] = numberRowsAtContinuous_+i ;
3801        }
3802      }
3803      int numberRowsNow=numberRowsAtContinuous_+numberToAdd;
3804      lastws->compressRows(numberToDrop,cutsToDrop) ;
3805      lastws->resize(numberRowsNow,numberColumns);
3806      solver_->applyRowCuts(numberToAdd,addCuts);
3807#     ifdef CBC_CHECK_BASIS
3808      printf("addCuts: stripped basis; rows %d + %d\n",
3809             numberRowsAtContinuous_,numberToAdd);
3810      lastws->print();
3811#     endif
3812      for (i=0;i<numberToAdd;i++)
3813        delete addCuts[i];
3814      delete [] addCuts;
3815      delete [] cutsToDrop ;
3816    }
3817/*
3818  Set the basis in the solver.
3819*/
3820    solver_->setWarmStart(lastws);
3821#if 0
3822    if ((numberNodes_%printFrequency_)==0) {
3823      printf("Objective %g, depth %d, unsatisfied %d\n",
3824             node->objectiveValue(),
3825             node->depth(),node->numberUnsatisfied());
3826    }
3827#endif
3828/*
3829  Clean up and we're out of here.
3830*/
3831    numberNodes_++;
3832    return 0;
3833  } 
3834/*
3835  This node has been fathomed by bound as we try to revive it out of the live
3836  set. Adjust the cut reference counts to reflect that we no longer need to
3837  explore the remaining branch arms, hence they will no longer reference any
3838  cuts. Cuts whose reference count falls to zero are deleted.
3839*/
3840  else
3841  { int i;
3842    int numberLeft = nodeInfo->numberBranchesLeft();
3843    for (i = 0 ; i < currentNumberCuts ; i++)
3844    { if (addedCuts_[i])
3845      { if (!addedCuts_[i]->decrement(numberLeft))
3846        { delete addedCuts_[i];
3847          addedCuts_[i] = NULL; } } }
3848    return 1 ; }
3849}
3850
3851
3852/*
3853  Perform reduced cost fixing on integer variables.
3854
3855  The variables in question are already nonbasic at bound. We're just nailing
3856  down the current situation.
3857*/
3858
3859int CbcModel::reducedCostFix ()
3860
3861{
3862  if(!solverCharacteristics_->reducedCostsAccurate())
3863    return 0; //NLP
3864  double cutoff = getCutoff() ;
3865  double direction = solver_->getObjSense() ;
3866  double gap = cutoff - solver_->getObjValue()*direction ;
3867  double tolerance;
3868  solver_->getDblParam(OsiDualTolerance,tolerance) ;
3869  if (gap<=0.0)
3870    return 0;
3871  gap += 100.0*tolerance;
3872  double integerTolerance = getDblParam(CbcIntegerTolerance) ;
3873
3874  const double *lower = solver_->getColLower() ;
3875  const double *upper = solver_->getColUpper() ;
3876  const double *solution = solver_->getColSolution() ;
3877  const double *reducedCost = solver_->getReducedCost() ;
3878
3879  int numberFixed = 0 ;
3880
3881# ifdef COIN_HAS_CLP
3882  OsiClpSolverInterface * clpSolver
3883    = dynamic_cast<OsiClpSolverInterface *> (solver_);
3884  ClpSimplex * clpSimplex=NULL;
3885  if (clpSolver) 
3886    clpSimplex = clpSolver->getModelPtr();
3887# endif
3888  for (int i = 0 ; i < numberIntegers_ ; i++)
3889  { int iColumn = integerVariable_[i] ;
3890    double djValue = direction*reducedCost[iColumn] ;
3891    if (upper[iColumn]-lower[iColumn] > integerTolerance)
3892    { if (solution[iColumn] < lower[iColumn]+integerTolerance && djValue > gap)
3893      { solver_->setColUpper(iColumn,lower[iColumn]) ;
3894#ifdef COIN_HAS_CLP
3895      // may just have been fixed before
3896      if (clpSimplex)
3897        assert(clpSimplex->getColumnStatus(iColumn)==ClpSimplex::atLowerBound||
3898               clpSimplex->getColumnStatus(iColumn)==ClpSimplex::isFixed);
3899#endif
3900      numberFixed++ ; }
3901      else
3902      if (solution[iColumn] > upper[iColumn]-integerTolerance && -djValue > gap)
3903      { solver_->setColLower(iColumn,upper[iColumn]) ;
3904#ifdef COIN_HAS_CLP
3905      // may just have been fixed before
3906      if (clpSimplex)
3907        assert(clpSimplex->getColumnStatus(iColumn)==ClpSimplex::atUpperBound||
3908               clpSimplex->getColumnStatus(iColumn)==ClpSimplex::isFixed);
3909#endif
3910      numberFixed++ ; } } }
3911  numberDJFixed_ += numberFixed;
3912 
3913  return numberFixed; }
3914
3915// Collect coding to replace whichGenerator
3916void
3917CbcModel::resizeWhichGenerator(int numberNow, int numberAfter)
3918{
3919  if (numberAfter > maximumWhich_) {
3920    maximumWhich_ = CoinMax(maximumWhich_*2+100,numberAfter) ;
3921    int * temp = new int[2*maximumWhich_] ;
3922    memcpy(temp,whichGenerator_,numberNow*sizeof(int)) ;
3923    delete [] whichGenerator_ ;
3924    whichGenerator_ = temp ;
3925    memset(whichGenerator_+numberNow,0,(maximumWhich_-numberNow)*sizeof(int));
3926  }
3927}
3928
3929/** Solve the model using cuts
3930
3931  This version takes off redundant cuts from node.
3932  Returns true if feasible.
3933
3934  \todo
3935  Why do I need to resolve the problem? What has been done between the last
3936  relaxation and calling solveWithCuts?
3937
3938  If numberTries == 0 then user did not want any cuts.
3939*/
3940
3941bool 
3942CbcModel::solveWithCuts (OsiCuts &cuts, int numberTries, CbcNode *node)
3943/*
3944  Parameters:
3945    numberTries: (i) the maximum number of iterations for this round of cut
3946                     generation; if negative then we don't mind if drop is tiny.
3947   
3948    cuts:       (o) all cuts generated in this round of cut generation
3949
3950    node: (i)     So we can update dynamic pseudo costs
3951*/
3952                       
3953
3954{
3955# ifdef COIN_HAS_CLP
3956  OsiClpSolverInterface * clpSolver
3957    = dynamic_cast<OsiClpSolverInterface *> (solver_);
3958  int saveClpOptions=0;
3959  if (clpSolver) 
3960    saveClpOptions = clpSolver->specialOptions();
3961# endif
3962  bool feasible = true ;
3963  int lastNumberCuts = 0 ;
3964  double lastObjective = -1.0e100 ;
3965  int violated = 0 ;
3966  int numberRowsAtStart = solver_->getNumRows() ;
3967  //printf("solver had %d rows\n",numberRowsAtStart);
3968  int numberColumns = solver_->getNumCols() ;
3969  CoinBigIndex numberElementsAtStart = solver_->getNumElements();
3970
3971  numberOldActiveCuts_ = numberRowsAtStart-numberRowsAtContinuous_ ;
3972  numberNewCuts_ = 0 ;
3973
3974  bool onOptimalPath = false ;
3975  const OsiRowCutDebugger *debugger = NULL;
3976  if ((specialOptions_&1)!=0) {
3977    /*
3978      See OsiRowCutDebugger for details. In a nutshell, make sure that current
3979      variable values do not conflict with a known optimal solution. (Obviously
3980      this can be fooled when there are multiple solutions.)
3981    */
3982    debugger = solver_->getRowCutDebugger() ;
3983    if (debugger) 
3984      onOptimalPath = (debugger->onOptimalPath(*solver_)) ;
3985  }
3986  OsiCuts slackCuts;
3987/*
3988  Resolve the problem. If we've lost feasibility, might as well bail out right
3989  after the debug stuff. The resolve will also refresh cached copies of the
3990  solver solution (cbcColLower_, ...) held by CbcModel.
3991*/
3992  double objectiveValue = solver_->getObjValue()*solver_->getObjSense();
3993  if (node)
3994    objectiveValue= node->objectiveValue();
3995  int returnCode = resolve(node ? node->nodeInfo() : NULL,1);
3996#ifdef COIN_DEVELOP
3997  //if (!solver_->getIterationCount()&&solver_->isProvenOptimal())
3998  //printf("zero iterations on first solve of branch\n");
3999#endif
4000  if (node&&!node->nodeInfo()->numberBranchesLeft())
4001    node->nodeInfo()->allBranchesGone(); // can clean up
4002  feasible = returnCode  != 0 ;
4003  if (returnCode<0)
4004    numberTries=0;
4005  if (problemFeasibility_->feasible(this,0)<0) {
4006    feasible=false; // pretend infeasible
4007  }
4008 
4009  // Update branching information if wanted
4010  if(node &&branchingMethod_)
4011    branchingMethod_->updateInformation(solver_,node);
4012
4013#ifdef CBC_DEBUG
4014  if (feasible)
4015  { printf("Obj value %g (%s) %d rows\n",solver_->getObjValue(),
4016           (solver_->isProvenOptimal())?"proven":"unproven",
4017           solver_->getNumRows()) ; }
4018 
4019  else
4020  { printf("Infeasible %d rows\n",solver_->getNumRows()) ; }
4021#endif
4022  if ((specialOptions_&1)!=0) {
4023/*
4024  If the RowCutDebugger said we were compatible with the optimal solution,
4025  and now we're suddenly infeasible, we might be confused. Then again, we
4026  may have fathomed by bound, heading for a rediscovery of an optimal solution.
4027*/
4028    if (onOptimalPath && !solver_->isDualObjectiveLimitReached()) {
4029      if (!feasible) {
4030        solver_->writeMps("infeas");
4031        CoinWarmStartBasis *slack =
4032          dynamic_cast<CoinWarmStartBasis *>(solver_->getEmptyWarmStart()) ;
4033        solver_->setWarmStart(slack);
4034        delete slack ;
4035        solver_->setHintParam(OsiDoReducePrint,false,OsiHintDo,0) ;
4036        solver_->initialSolve();
4037      }
4038      assert(feasible) ;
4039    }
4040  }
4041
4042  if (!feasible) {
4043    numberInfeasibleNodes_++;
4044# ifdef COIN_HAS_CLP
4045    if (clpSolver) 
4046    clpSolver->setSpecialOptions(saveClpOptions);
4047# endif
4048    return (false) ;
4049  }
4050  sumChangeObjective1_ += solver_->getObjValue()*solver_->getObjSense()
4051    - objectiveValue ;
4052  if ( CoinCpuTime()-dblParam_[CbcStartSeconds] > dblParam_[CbcMaximumSeconds] )
4053    numberTries=0; // exit
4054  //if ((numberNodes_%100)==0)
4055  //printf("XXa sum obj changed by %g\n",sumChangeObjective1_);
4056  objectiveValue = solver_->getObjValue()*solver_->getObjSense();
4057  // Return at once if numberTries zero
4058  if (!numberTries) {
4059    cuts=OsiCuts();
4060    numberNewCuts_=0;
4061# ifdef COIN_HAS_CLP
4062    if (clpSolver) 
4063    clpSolver->setSpecialOptions(saveClpOptions);
4064# endif
4065    return true;
4066  }
4067/*
4068  Do reduced cost fixing.
4069*/
4070  reducedCostFix() ;
4071/*
4072  Set up for at most numberTries rounds of cut generation. If numberTries is
4073  negative, we'll ignore the minimumDrop_ cutoff and keep generating cuts for
4074  the specified number of rounds.
4075*/
4076  double minimumDrop = minimumDrop_ ;
4077  if (numberTries<0)
4078  { numberTries = -numberTries ;
4079    minimumDrop = -1.0 ; }
4080/*
4081  Is it time to scan the cuts in order to remove redundant cuts? If so, set
4082  up to do it.
4083*/
4084# define SCANCUTS 100 
4085  int *countColumnCuts = NULL ;
4086  // Always accumulate row cut counts
4087  int * countRowCuts =new int[numberCutGenerators_+numberHeuristics_] ;
4088  memset(countRowCuts,0,
4089         (numberCutGenerators_+numberHeuristics_)*sizeof(int)) ;
4090  bool fullScan = false ;
4091  if ((numberNodes_%SCANCUTS) == 0)
4092  { fullScan = true ;
4093    countColumnCuts = new int[numberCutGenerators_+numberHeuristics_] ;
4094    memset(countColumnCuts,0,
4095           (numberCutGenerators_+numberHeuristics_)*sizeof(int)) ; }
4096
4097  double direction = solver_->getObjSense() ;
4098  double startObjective = solver_->getObjValue()*direction ;
4099
4100  currentPassNumber_ = 0 ;
4101  double primalTolerance = 1.0e-7 ;
4102  // We may need to keep going on
4103  bool keepGoing=false;
4104/*
4105  Begin cut generation loop. Cuts generated during each iteration are
4106  collected in theseCuts. The loop can be divided into four phases:
4107   1) Prep: Fix variables using reduced cost. In the first iteration only,
4108      consider scanning globalCuts_ and activating any applicable cuts.
4109   2) Cut Generation: Call each generator and heuristic registered in the
4110      generator_ and heuristic_ arrays. Newly generated global cuts are
4111      copied to globalCuts_ at this time.
4112   3) Cut Installation and Reoptimisation: Install column and row cuts in
4113      the solver. Copy row cuts to cuts (parameter). Reoptimise.
4114   4) Cut Purging: takeOffCuts() removes inactive cuts from the solver, and
4115      does the necessary bookkeeping in the model.
4116*/
4117  do
4118  { currentPassNumber_++ ;
4119    numberTries-- ;
4120    if (numberTries<0&&keepGoing) {
4121      // switch off all normal ones
4122      for (int i = 0;i<numberCutGenerators_;i++) {
4123        if (!generator_[i]->mustCallAgain())
4124          generator_[i]->setSwitchedOff(true);
4125      }
4126    }
4127    keepGoing=false;
4128    OsiCuts theseCuts ;
4129/*
4130  Scan previously generated global column and row cuts to see if any are
4131  useful.
4132*/
4133    int numberViolated=0;
4134    if (currentPassNumber_ == 1 && howOftenGlobalScan_ > 0 &&
4135        (numberNodes_%howOftenGlobalScan_) == 0)
4136    { int numberCuts = globalCuts_.sizeColCuts() ;
4137      int i;
4138      // possibly extend whichGenerator
4139      resizeWhichGenerator(numberViolated, numberViolated+numberCuts);
4140      for ( i = 0 ; i < numberCuts ; i++)
4141      { OsiColCut *thisCut = globalCuts_.colCutPtr(i) ;
4142        if (thisCut->violated(cbcColSolution_)>primalTolerance) {
4143          printf("Global cut added - violation %g\n",
4144                 thisCut->violated(cbcColSolution_)) ;
4145          whichGenerator_[numberViolated++]=-1;
4146#ifndef GLOBAL_CUTS_JUST_POINTERS
4147          theseCuts.insert(*thisCut) ;
4148#else
4149          theseCuts.insert(thisCut) ;
4150#endif
4151        }
4152      }
4153      numberCuts = globalCuts_.sizeRowCuts() ;
4154      // possibly extend whichGenerator
4155      resizeWhichGenerator(numberViolated, numberViolated+numberCuts);
4156      for ( i = 0;i<numberCuts;i++) {
4157        OsiRowCut * thisCut = globalCuts_.rowCutPtr(i) ;
4158        if (thisCut->violated(cbcColSolution_)>primalTolerance) {
4159          //printf("Global cut added - violation %g\n",
4160          // thisCut->violated(cbcColSolution_)) ;
4161          whichGenerator_[numberViolated++]=-1;
4162#ifndef GLOBAL_CUTS_JUST_POINTERS
4163          theseCuts.insert(*thisCut) ;
4164#else
4165          theseCuts.insert(thisCut) ;
4166#endif
4167        }
4168      }
4169      numberGlobalViolations_+=numberViolated;
4170    }
4171/*
4172  Generate new cuts (global and/or local) and/or apply heuristics.  If
4173  CglProbing is used, then it should be first as it can fix continuous
4174  variables.
4175
4176  At present, CglProbing is the only case where generateCuts will return
4177  true. generateCuts actually modifies variable bounds in the solver when
4178  CglProbing indicates that it can fix a variable. Reoptimisation is required
4179  to take full advantage.
4180
4181  The need to resolve here should only happen after a heuristic solution.
4182  (Note default OSI implementation of optimalBasisIsAvailable always returns
4183  false.)
4184*/
4185    if (solverCharacteristics_->warmStart()&&
4186        !solver_->optimalBasisIsAvailable()) {
4187      //printf("XXXXYY no opt basis\n");
4188      resolve(node ? node->nodeInfo() : NULL,3);
4189    }
4190    if (nextRowCut_) {
4191      // branch was a cut - add it
4192      theseCuts.insert(*nextRowCut_);
4193      if (handler_->logLevel()>1)
4194        nextRowCut_->print();
4195      const OsiRowCut * cut=nextRowCut_;
4196      double lb = cut->lb();
4197      double ub = cut->ub();
4198      int n=cut->row().getNumElements();
4199      const int * column = cut->row().getIndices();
4200      const double * element = cut->row().getElements();
4201      double sum=0.0;
4202      for (int i=0;i<n;i++) {
4203        int iColumn = column[i];
4204        double value = element[i];
4205        //if (cbcColSolution_[iColumn]>1.0e-7)
4206        //printf("value of %d is %g\n",iColumn,cbcColSolution_[iColumn]);
4207        sum += value * cbcColSolution_[iColumn];
4208      }
4209      delete nextRowCut_;
4210      nextRowCut_=NULL;
4211      if (handler_->logLevel()>1)
4212        printf("applying branch cut, sum is %g, bounds %g %g\n",sum,lb,ub);
4213      // possibly extend whichGenerator
4214      resizeWhichGenerator(numberViolated, numberViolated+1);
4215      // set whichgenerator (also serves as marker to say don't delete0
4216      whichGenerator_[numberViolated++]=-2;
4217    }
4218    double * newSolution = new double [numberColumns] ;
4219    double heuristicValue = getCutoff() ;
4220    int found = -1; // no solution found
4221
4222    // reset probing info
4223    if (probingInfo_)
4224      probingInfo_->initializeFixing();
4225    for (int i = 0;i<numberCutGenerators_+numberHeuristics_;i++) {
4226      int numberRowCutsBefore = theseCuts.sizeRowCuts() ;
4227      int numberColumnCutsBefore = theseCuts.sizeColCuts() ;
4228      if (i<numberCutGenerators_) {
4229        bool generate = generator_[i]->normal();
4230        // skip if not optimal and should be (maybe a cut generator has fixed variables)
4231        if (generator_[i]->needsOptimalBasis()&&!solver_->basisIsAvailable())
4232          generate=false;
4233        if (generator_[i]->switchedOff())
4234          generate=false;;
4235        if (generate) {
4236          bool mustResolve = 
4237            generator_[i]->generateCuts(theseCuts,fullScan,node) ;
4238          int numberRowCutsAfter = theseCuts.sizeRowCuts() ;
4239          if(numberRowCutsBefore < numberRowCutsAfter &&
4240             generator_[i]->mustCallAgain())
4241            keepGoing=true; // say must go round
4242          // Check last cut to see if infeasible
4243          if(numberRowCutsBefore < numberRowCutsAfter) {
4244            const OsiRowCut * thisCut = theseCuts.rowCutPtr(numberRowCutsAfter-1) ;
4245            if (thisCut->lb()>thisCut->ub()) {
4246              feasible = false; // sub-problem is infeasible
4247              break;
4248            }
4249          }
4250#ifdef CBC_DEBUG
4251          {
4252            int numberRowCutsAfter = theseCuts.sizeRowCuts() ;
4253            int k ;
4254            for (k = numberRowCutsBefore;k<numberRowCutsAfter;k++) {
4255              OsiRowCut thisCut = theseCuts.rowCut(k) ;
4256              /* check size of elements.
4257                 We can allow smaller but this helps debug generators as it
4258                 is unsafe to have small elements */
4259              int n=thisCut.row().getNumElements();
4260              const int * column = thisCut.row().getIndices();
4261              const double * element = thisCut.row().getElements();
4262              //assert (n);
4263              for (int i=0;i<n;i++) {
4264                double value = element[i];
4265                assert(fabs(value)>1.0e-12&&fabs(value)<1.0e20);
4266              }
4267            }
4268          }
4269#endif
4270          if (mustResolve) {
4271            int returncode = resolve(node ? node->nodeInfo() : NULL,2);
4272            feasible = returnCode  != 0 ;
4273            if (returncode<0)
4274              numberTries=0;
4275            if ((specialOptions_&1)!=0) {
4276              debugger = solver_->getRowCutDebugger() ;
4277              if (debugger) 
4278                onOptimalPath = (debugger->onOptimalPath(*solver_)) ;
4279              else
4280                onOptimalPath=false;
4281              if (onOptimalPath && !solver_->isDualObjectiveLimitReached())
4282                assert(feasible) ;
4283            }
4284            if (!feasible)
4285              break ;
4286          }
4287        }
4288      } else {
4289        // see if heuristic will do anything
4290        double saveValue = heuristicValue ;
4291        int ifSol = 
4292          heuristic_[i-numberCutGenerators_]->solution(heuristicValue,
4293                                                       newSolution,
4294                                                       theseCuts) ;
4295        if (ifSol>0) {
4296          // better solution found
4297          found = i-numberCutGenerators_ ;
4298          incrementUsed(newSolution);
4299        } else if (ifSol<0) {
4300          heuristicValue = saveValue ;
4301        }
4302      }
4303      int numberRowCutsAfter = theseCuts.sizeRowCuts() ;
4304      int numberColumnCutsAfter = theseCuts.sizeColCuts() ;
4305
4306      if ((specialOptions_&1)!=0) {
4307        if (onOptimalPath) {
4308          int k ;
4309          for (k = numberRowCutsBefore;k<numberRowCutsAfter;k++) {
4310            OsiRowCut thisCut = theseCuts.rowCut(k) ;
4311            if(debugger->invalidCut(thisCut)) {
4312              solver_->writeMps("badCut");
4313            }
4314            assert(!debugger->invalidCut(thisCut)) ;
4315          }
4316        }
4317      }
4318/*
4319  The cut generator/heuristic has done its thing, and maybe it generated some
4320  cuts and/or a new solution.  Do a bit of bookkeeping: load
4321  whichGenerator[i] with the index of the generator responsible for a cut,
4322  and place cuts flagged as global in the global cut pool for the model.
4323
4324  lastNumberCuts is the sum of cuts added in previous iterations; it's the
4325  offset to the proper starting position in whichGenerator.
4326*/
4327      int numberBefore =
4328            numberRowCutsBefore+numberColumnCutsBefore+lastNumberCuts ;
4329      int numberAfter =
4330            numberRowCutsAfter+numberColumnCutsAfter+lastNumberCuts ;
4331      // possibly extend whichGenerator
4332      resizeWhichGenerator(numberBefore, numberAfter);
4333      int j ;
4334      if (fullScan) {
4335        // counts
4336        countColumnCuts[i] += numberColumnCutsAfter-numberColumnCutsBefore ;
4337      }
4338      countRowCuts[i] += numberRowCutsAfter-numberRowCutsBefore ;
4339     
4340      for (j = numberRowCutsBefore;j<numberRowCutsAfter;j++) {
4341        whichGenerator_[numberBefore++] = i ;
4342        const OsiRowCut * thisCut = theseCuts.rowCutPtr(j) ;
4343        if (thisCut->lb()>thisCut->ub())
4344          violated=-2; // sub-problem is infeasible
4345        if (thisCut->globallyValid()) {
4346          // add to global list
4347          OsiRowCut newCut(*thisCut);
4348          newCut.setGloballyValid(2);
4349          newCut.mutableRow().setTestForDuplicateIndex(false);
4350          globalCuts_.insert(newCut) ;
4351        }
4352      }
4353      for (j = numberColumnCutsBefore;j<numberColumnCutsAfter;j++) {
4354        whichGenerator_[numberBefore++] = i ;
4355        const OsiColCut * thisCut = theseCuts.colCutPtr(j) ;
4356        if (thisCut->globallyValid()) {
4357          // add to global list
4358          OsiColCut newCut(*thisCut);
4359          newCut.setGloballyValid(2);
4360          globalCuts_.insert(newCut) ;
4361        }
4362      }
4363    }
4364    // If at root node and first pass do heuristics without cuts
4365    if (!numberNodes_&&currentPassNumber_==1) {
4366      // Save number solutions
4367      int saveNumberSolutions = numberSolutions_;
4368      int saveNumberHeuristicSolutions = numberHeuristicSolutions_;
4369      // make sure bounds tight
4370      { 
4371        for (int i = 0;i < numberObjects_;i++)
4372          object_[i]->resetBounds(solver_) ;
4373      }
4374      for (int i = 0;i<numberHeuristics_;i++) {
4375        // see if heuristic will do anything
4376        double saveValue = heuristicValue ;
4377        int ifSol = heuristic_[i]->solution(heuristicValue,
4378                                            newSolution);
4379        if (ifSol>0) {
4380          // better solution found
4381          found = i ;
4382          incrementUsed(newSolution);
4383          // increment number of solutions so other heuristics can test
4384          numberSolutions_++;
4385          numberHeuristicSolutions_++;
4386        } else {
4387          heuristicValue = saveValue ;
4388        }
4389      }
4390      // Restore number solutions
4391      numberSolutions_ = saveNumberSolutions;
4392      numberHeuristicSolutions_ = saveNumberHeuristicSolutions;
4393    }
4394    // Add in any violated saved cuts
4395    if (!theseCuts.sizeRowCuts()&&!theseCuts.sizeColCuts()) {
4396      int numberOld = theseCuts.sizeRowCuts();
4397      int numberCuts = slackCuts.sizeRowCuts() ;
4398      int i;
4399      // possibly extend whichGenerator
4400      resizeWhichGenerator(numberOld, numberOld+numberCuts);
4401      for ( i = 0;i<numberCuts;i++) {
4402        const OsiRowCut * thisCut = slackCuts.rowCutPtr(i) ;
4403        if (thisCut->violated(cbcColSolution_)>100.0*primalTolerance) {
4404          if (messageHandler()->logLevel()>2)
4405            printf("Old cut added - violation %g\n",
4406                   thisCut->violated(cbcColSolution_)) ;
4407          whichGenerator_[numberOld++]=-1;
4408          theseCuts.insert(*thisCut) ;
4409        }
4410      }
4411    }
4412/*
4413  End of the loop to exercise each generator/heuristic.
4414
4415  Did any of the heuristics turn up a new solution? Record it before we free
4416  the vector.
4417*/
4418    if (found >= 0) { 
4419      phase_=4;
4420      incrementUsed(newSolution);
4421      setBestSolution(CBC_ROUNDING,heuristicValue,newSolution) ;
4422      lastHeuristic_ = heuristic_[found];
4423      CbcTreeLocal * tree
4424          = dynamic_cast<CbcTreeLocal *> (tree_);
4425      if (tree)
4426        tree->passInSolution(bestSolution_,heuristicValue);
4427    }
4428    delete [] newSolution ;
4429
4430#if 0
4431    // switch on to get all cuts printed
4432    theseCuts.printCuts() ;
4433#endif
4434    int numberColumnCuts = theseCuts.sizeColCuts() ;
4435    int numberRowCuts = theseCuts.sizeRowCuts() ;
4436    if (violated>=0)
4437      violated = numberRowCuts + numberColumnCuts ;
4438/*
4439  Apply column cuts (aka bound tightening). This may be partially redundant
4440  for column cuts returned by CglProbing, as generateCuts installs bounds
4441  from CglProbing when it determines it can fix a variable.
4442
4443  TODO: Looks like the use of violated has evolved. The value set above is
4444        completely ignored. All that's left is violated == -1 indicates some
4445        cut is violated, violated == -2 indicates infeasibility. Only
4446        infeasibility warrants exceptional action.
4447
4448  TODO: Strikes me that this code will fail to detect infeasibility, because
4449        the breaks escape the inner loops but the outer loop keeps going.
4450        Infeasibility in an early cut will be overwritten if a later cut is
4451        merely violated.
4452*/
4453    if (numberColumnCuts) {
4454
4455#ifdef CBC_DEBUG
4456      double * oldLower = new double [numberColumns] ;
4457      double * oldUpper = new double [numberColumns] ;
4458      memcpy(oldLower,cbcColLower_,numberColumns*sizeof(double)) ;
4459      memcpy(oldUpper,cbcColUpper_,numberColumns*sizeof(double)) ;
4460#endif
4461
4462      double integerTolerance = getDblParam(CbcIntegerTolerance) ;
4463      for (int i = 0;i<numberColumnCuts;i++) {
4464        const OsiColCut * thisCut = theseCuts.colCutPtr(i) ;
4465        const CoinPackedVector & lbs = thisCut->lbs() ;
4466        const CoinPackedVector & ubs = thisCut->ubs() ;
4467        int j ;
4468        int n ;
4469        const int * which ;
4470        const double * values ;
4471        n = lbs.getNumElements() ;
4472        which = lbs.getIndices() ;
4473        values = lbs.getElements() ;
4474        for (j = 0;j<n;j++) {
4475          int iColumn = which[j] ;
4476          double value = cbcColSolution_[iColumn] ;
4477#if CBC_DEBUG>1
4478          printf("%d %g %g %g %g\n",iColumn,oldLower[iColumn],
4479                 cbcColSolution_[iColumn],oldUpper[iColumn],values[j]) ;
4480#endif
4481          solver_->setColLower(iColumn,values[j]) ;
4482          if (value<values[j]-integerTolerance)
4483            violated = -1 ;
4484          if (values[j]>cbcColUpper_[iColumn]+integerTolerance) {
4485            // infeasible
4486            violated = -2 ;
4487            break ;
4488          }
4489        }
4490        n = ubs.getNumElements() ;
4491        which = ubs.getIndices() ;
4492        values = ubs.getElements() ;
4493        for (j = 0;j<n;j++) {
4494          int iColumn = which[j] ;
4495          double value = cbcColSolution_[iColumn] ;
4496#if CBC_DEBUG>1
4497          printf("%d %g %g %g %g\n",iColumn,oldLower[iColumn],
4498                 cbcColSolution_[iColumn],oldUpper[iColumn],values[j]) ;
4499#endif
4500          solver_->setColUpper(iColumn,values[j]) ;
4501          if (value>values[j]+integerTolerance)
4502            violated = -1 ;
4503          if (values[j]<cbcColLower_[iColumn]-integerTolerance) {
4504            // infeasible
4505            violated = -2 ;
4506            break ;
4507          }
4508        }
4509      }
4510#ifdef CBC_DEBUG
4511      delete [] oldLower ;
4512      delete [] oldUpper ;
4513#endif
4514    }
4515/*
4516  End installation of column cuts. The break here escapes the numberTries
4517  loop.
4518*/
4519    if (violated == -2||!feasible) {
4520      // infeasible
4521      feasible = false ;
4522      violated = -2;
4523      if (!numberNodes_) 
4524        messageHandler()->message(CBC_INFEAS,
4525                                  messages())
4526                                    << CoinMessageEol ;
4527      break ;
4528    }
4529/*
4530  Now apply the row (constraint) cuts. This is a bit more work because we need
4531  to obtain and augment the current basis.
4532
4533  TODO: Why do this work, if there are no row cuts? The current basis will do
4534        just fine.
4535*/
4536    int numberRowsNow = solver_->getNumRows() ;
4537    assert(numberRowsNow == numberRowsAtStart+lastNumberCuts) ;
4538    int numberToAdd = theseCuts.sizeRowCuts() ;
4539    numberNewCuts_ = lastNumberCuts+numberToAdd ;
4540/*
4541  Now actually add the row cuts and reoptimise.
4542
4543  Install the cuts in the solver using applyRowCuts and
4544  augment the basis with the corresponding slack. We also add each row cut to
4545  the set of row cuts (cuts.insert()) supplied as a parameter. The new basis
4546  must be set with setWarmStart().
4547
4548  TODO: Seems to me the original code could allocate addCuts with size 0, if
4549        numberRowCuts was 0 and numberColumnCuts was nonzero. That might
4550        explain the memory fault noted in the comment by AJK.  Unfortunately,
4551        just commenting out the delete[] results in massive memory leaks. Try
4552        a revision to separate the row cut case. Why do we need addCuts at
4553        all? A typing issue, apparently: OsiCut vs. OsiRowCut.
4554
4555  TODO: It looks to me as if numberToAdd and numberRowCuts are identical at
4556        this point. Confirm & get rid of one of them.
4557
4558  TODO: Any reason why the three loops can't be consolidated?
4559*/
4560    if (numberRowCuts > 0 || numberColumnCuts > 0)
4561    { if (numberToAdd > 0)
4562      { int i ;
4563        // Faster to add all at once
4564        const OsiRowCut ** addCuts = new const OsiRowCut * [numberToAdd] ;
4565        for (i = 0 ; i < numberToAdd ; i++)
4566        { addCuts[i] = &theseCuts.rowCut(i) ; }
4567        solver_->applyRowCuts(numberToAdd,addCuts) ;
4568        // AJK this caused a memory fault on Win32
4569        // may be okay now with ** form
4570        delete [] addCuts ;
4571        for (i = 0 ; i < numberToAdd ; i++)
4572        { cuts.insert(theseCuts.rowCut(i)) ; }
4573        CoinWarmStartBasis * basis = dynamic_cast<CoinWarmStartBasis*>(solver_->getWarmStart()) ;
4574        assert(basis != NULL); // make sure not volume
4575/* dylp bug
4576
4577  Consistent size used by OsiDylp as sanity check. Implicit resize seen
4578  as an error. Hence this call to resize is necessary.
4579*/
4580        basis->resize(numberRowsAtStart+numberNewCuts_,numberColumns) ;
4581        for (i = 0 ; i < numberToAdd ; i++) 
4582        { basis->setArtifStatus(numberRowsNow+i,
4583                                 CoinWarmStartBasis::basic) ; }
4584        if (solver_->setWarmStart(basis) == false)
4585        { throw CoinError("Fail setWarmStart() after cut installation.",
4586                          "solveWithCuts","CbcModel") ; }
4587        delete basis;
4588      }
4589      feasible = resolve(node ? node->nodeInfo() : NULL,2) ;
4590      if ( CoinCpuTime()-dblParam_[CbcStartSeconds] > dblParam_[CbcMaximumSeconds] )
4591        numberTries=0; // exit
4592#     ifdef CBC_DEBUG
4593      printf("Obj value after cuts %g %d rows\n",solver_->getObjValue(),
4594              solver_->getNumRows()) ;
4595      if (onOptimalPath && !solver_->isDualObjectiveLimitReached())
4596        assert(feasible) ;
4597#     endif
4598    }
4599/*
4600  No cuts. Cut short the cut generation (numberTries) loop.
4601*/
4602    else
4603    { numberTries = 0 ; }
4604/*
4605  If the problem is still feasible, first, call takeOffCuts() to remove cuts
4606  that are now slack. takeOffCuts() will call the solver to reoptimise if
4607  that's needed to restore a valid solution.
4608
4609  Next, see if we should quit due to diminishing returns:
4610    * we've tried three rounds of cut generation and we're getting
4611      insufficient improvement in the objective; or
4612    * we generated no cuts; or
4613    * the solver declared optimality with 0 iterations after we added the
4614      cuts generated in this round.
4615  If we decide to keep going, prep for the next iteration.
4616
4617  It just seems more safe to tell takeOffCuts() to call resolve(), even if
4618  we're not continuing cut generation. Otherwise code executed between here
4619  and final disposition of the node will need to be careful not to access the
4620  lp solution. It can happen that we lose feasibility in takeOffCuts ---
4621  numerical jitters when the cutoff bound is epsilon less than the current
4622  best, and we're evaluating an alternative optimum.
4623
4624  TODO: After successive rounds of code motion, there seems no need to
4625        distinguish between the two checks for aborting the cut generation
4626        loop. Confirm and clean up.
4627*/
4628    if (feasible)
4629    { int cutIterations = solver_->getIterationCount() ;
4630      if (numberOldActiveCuts_+numberNewCuts_) {
4631        OsiCuts * saveCuts = node ? NULL : &slackCuts;
4632        takeOffCuts(cuts,resolveAfterTakeOffCuts_,saveCuts) ;
4633        if (solver_->isDualObjectiveLimitReached()&&resolveAfterTakeOffCuts_)
4634          { feasible = false ;
4635#       ifdef CBC_DEBUG
4636          double z = solver_->getObjValue() ;
4637          double cut = getCutoff() ;
4638          printf("Lost feasibility by %g in takeOffCuts; z = %g, cutoff = %g\n",
4639                 z-cut,z,cut) ;
4640#       endif
4641          }
4642      }
4643      if (feasible) {
4644        numberRowsAtStart = numberOldActiveCuts_+numberRowsAtContinuous_ ;
4645        lastNumberCuts = numberNewCuts_ ;
4646        if (direction*solver_->getObjValue() < lastObjective+minimumDrop &&
4647            currentPassNumber_ >= 3 && !keepGoing)
4648          { numberTries = 0 ; }
4649        if (numberRowCuts+numberColumnCuts == 0 || cutIterations == 0)
4650          { break ; }
4651        if (numberTries > 0) {
4652          reducedCostFix() ;
4653          lastObjective = direction*solver_->getObjValue() ;
4654        }
4655      }
4656    }
4657/*
4658  We've lost feasibility --- this node won't be referencing the cuts we've
4659  been collecting, so decrement the reference counts.
4660*/
4661    if (!feasible)
4662    { int i ;
4663      for (i = 0;i<currentNumberCuts_;i++) {
4664        // take off node
4665        if (addedCuts_[i]) {
4666          if (!addedCuts_[i]->decrement())
4667            delete addedCuts_[i] ;
4668          addedCuts_[i] = NULL ;
4669        }
4670      }
4671      numberTries = 0 ;
4672    }
4673  } while (numberTries>0||keepGoing) ;
4674  {
4675    // switch on
4676    for (int i = 0;i<numberCutGenerators_;i++) 
4677      generator_[i]->setSwitchedOff(false);
4678  }
4679  //check feasibility.
4680  //If solution seems to be integer feasible calling setBestSolution
4681  //will eventually add extra global cuts which we need to install at
4682  //the nodes
4683
4684
4685  if(feasible&&solverCharacteristics_->solutionAddsCuts())  //check integer feasibility
4686  {
4687    bool integerFeasible = true;
4688    const double * save = testSolution_;
4689    testSolution_ = solver_->getColSolution();
4690    // point to useful information
4691    OsiBranchingInformation usefulInfo=usefulInformation();
4692    for (int i=0;i<numberObjects_ && integerFeasible;i++)
4693    {
4694      int preferredWay;
4695      double infeasibility = object_[i]->infeasibility(&usefulInfo,preferredWay);
4696      if(infeasibility)
4697        integerFeasible = false;
4698    }
4699    testSolution_ = save;
4700    if(integerFeasible)//update
4701    {
4702      double objValue = solver_->getObjValue();
4703      int numberGlobalBefore = globalCuts_.sizeRowCuts();
4704      // SOLUTION2 so won't up cutoff or print message
4705      setBestSolution(CBC_SOLUTION2, objValue, 
4706                      solver_->getColSolution(),0);
4707      int numberGlobalAfter = globalCuts_.sizeRowCuts();
4708      int numberToAdd = numberGlobalAfter - numberGlobalBefore;
4709      if(numberToAdd > 0)
4710        //We have added some cuts say they are tight at that node
4711        //Basis and lp should already have been updated
4712      {
4713        feasible = (solver_->isProvenOptimal() &&
4714                    !solver_->isDualObjectiveLimitReached()) ;
4715        if(feasible)
4716        {
4717          int numberCuts = numberNewCuts_ = cuts.sizeRowCuts();
4718      // possibly extend whichGenerator
4719          resizeWhichGenerator(numberCuts, numberToAdd+numberCuts);
4720         
4721          for (int i = numberGlobalBefore ; i < numberGlobalAfter ; i++)
4722            {       
4723              whichGenerator_[numberNewCuts_++]=-1;
4724#ifndef GLOBAL_CUTS_JUST_POINTERS
4725              cuts.insert(globalCuts_.rowCut(i)) ; 
4726#else
4727              OsiRowCut * rowCutPointer = globalCuts_.rowCutPtr(i);
4728              cuts.insert(rowCutPointer) ; 
4729#endif
4730            }
4731          numberNewCuts_ = lastNumberCuts+numberToAdd;
4732          //now take off the cuts which are not tight anymore
4733          takeOffCuts(cuts,resolveAfterTakeOffCuts_, NULL) ;
4734          if (solver_->isDualObjectiveLimitReached()&&resolveAfterTakeOffCuts_)
4735            { feasible = false ;}
4736        }
4737        if(!feasible) //node will be fathomed
4738          { 
4739          for (int i = 0;i<currentNumberCuts_;i++) 
4740            {
4741            // take off node
4742            if (addedCuts_[i])
4743              {
4744                if (!addedCuts_[i]->decrement())
4745                  delete addedCuts_[i] ;
4746                addedCuts_[i] = NULL ;
4747              }
4748            }
4749        }
4750      }
4751    }
4752  }
4753/*
4754  Reduced cost fix at end. Must also check feasible, in case we've popped out
4755  because a generator indicated we're infeasible.
4756*/
4757  if (feasible && solver_->isProvenOptimal())
4758    reducedCostFix();
4759  // If at root node do heuristics
4760  if (!numberNodes_) {
4761    // First see if any cuts are slack
4762    int numberRows = solver_->getNumRows();
4763    int numberAdded = numberRows-numberRowsAtContinuous_;
4764    if (numberAdded) {
4765      CoinWarmStartBasis * basis = dynamic_cast<CoinWarmStartBasis*>(solver_->getWarmStart()) ;
4766      assert(basis != NULL);
4767      int * added = new int[numberAdded];
4768      int nDelete = 0;
4769      for (int j=numberRowsAtContinuous_;j<numberRows;j++) {
4770        if (basis->getArtifStatus(j)==CoinWarmStartBasis::basic) {
4771          //printf("%d slack!\n",j);
4772          added[nDelete++]=j;
4773        }
4774      }
4775      if (nDelete) {
4776        solver_->deleteRows(nDelete,added);
4777      }
4778      delete [] added;
4779      delete basis ;
4780    }
4781    // mark so heuristics can tell
4782    int savePass=currentPassNumber_;
4783    currentPassNumber_=999999;
4784    double * newSolution = new double [numberColumns] ;
4785    double heuristicValue = getCutoff() ;
4786    int found = -1; // no solution found
4787    for (int i = 0;i<numberHeuristics_;i++) {
4788      // see if heuristic will do anything
4789      double saveValue = heuristicValue ;
4790      int ifSol = heuristic_[i]->solution(heuristicValue,
4791                                          newSolution);
4792      if (ifSol>0) {
4793        // better solution found
4794        found = i ;
4795        incrementUsed(newSolution);
4796      } else {
4797        heuristicValue = saveValue ;
4798      }
4799    }
4800    currentPassNumber_=savePass;
4801    if (found >= 0) { 
4802      phase_=4;
4803      incrementUsed(newSolution);
4804      setBestSolution(CBC_ROUNDING,heuristicValue,newSolution) ;
4805      lastHeuristic_ = heuristic_[found];
4806    }
4807    delete [] newSolution ;
4808  }
4809  // Up change due to cuts
4810  if (feasible)
4811    sumChangeObjective2_ += solver_->getObjValue()*solver_->getObjSense()
4812      - objectiveValue ;
4813  //if ((numberNodes_%100)==0)
4814  //printf("XXb sum obj changed by %g\n",sumChangeObjective2_);
4815/*
4816  End of cut generation loop.
4817
4818  Now, consider if we want to disable or adjust the frequency of use for any
4819  of the cut generators. If the client specified a positive number for
4820  howOften, it will never change. If the original value was negative, it'll
4821  be converted to 1000000+|howOften|, and this value will be adjusted each
4822  time fullScan is true. Actual cut generation is performed every
4823  howOften%1000000 nodes; the 1000000 offset is just a convenient way to
4824  specify that the frequency is adjustable.
4825
4826  During cut generation, we recorded the number of cuts produced by each
4827  generator for this node. For all cuts, whichGenerator records the generator
4828  that produced a cut.
4829
4830  TODO: All this should probably be hidden in a method of the CbcCutGenerator
4831  class.
4832
4833  TODO: Can the loop that scans over whichGenerator to accumulate per generator
4834        counts be replaced by values in countRowCuts and countColumnCuts?
4835*/
4836#ifdef NODE_LOG
4837  int fatherNum = (node == NULL) ? -1 : node->nodeNumber();
4838  double value =  (node == NULL) ? -1 : node->branchingObject()->value();
4839  string bigOne = (solver_->getIterationCount() > 30) ? "*******" : "";
4840  string way = (node == NULL) ? "" : (node->branchingObject()->way())==1 ? "Down" : "Up";
4841  std::cout<<"Node "<<numberNodes_<<", father "<<fatherNum<<", #iterations "<<solver_->getIterationCount()<<", sol value : "<<solver_->getObjValue()<<std::endl;
4842#endif
4843  if (fullScan&&numberCutGenerators_) {
4844    /* If cuts just at root node then it will probably be faster to
4845       update matrix and leave all in */
4846    int willBeCutsInTree=0;
4847    double thisObjective = solver_->getObjValue()*direction ;
4848    // get sizes
4849    int numberRowsAdded = solver_->getNumRows()-numberRowsAtStart;
4850    CoinBigIndex numberElementsAdded =  solver_->getNumElements()-numberElementsAtStart ;
4851    double densityOld = ((double) numberElementsAtStart)/((double) numberRowsAtStart);
4852    double densityNew = numberRowsAdded ? ((double) (numberElementsAdded))/((double) numberRowsAdded)
4853      : 0.0;
4854    if (!numberNodes_) {
4855      if (numberRowsAdded)
4856        handler_->message(CBC_CUTS_STATS,messages_)
4857          <<numberRowsAdded
4858          <<densityNew
4859          <<CoinMessageEol ;
4860      if (thisObjective-startObjective<1.0e-5&&numberElementsAdded>0.2*numberElementsAtStart)
4861        willBeCutsInTree=-1;
4862    }
4863    if ((numberRowsAdded>100+0.5*numberRowsAtStart
4864         ||numberElementsAdded>0.5*numberElementsAtStart)
4865        &&(densityNew>200.0&&numberRowsAdded>100&&densityNew>2.0*densityOld)) {
4866      // much bigger
4867      //if (thisObjective-startObjective<0.1*fabs(startObjective)+1.0e-5)
4868      //willBeCutsInTree=-1;
4869      //printf("Cuts will be taken off , %d rows added with density %g\n",
4870      //     numberRowsAdded,densityNew);
4871    }
4872    if (densityNew>100.0&&numberRowsAdded>2&&densityNew>2.0*densityOld) {
4873      //if (thisObjective-startObjective<0.1*fabs(startObjective)+1.0e-5)
4874      //willBeCutsInTree=-2;
4875      //printf("Density says no cuts ? , %d rows added with density %g\n",
4876      //     numberRowsAdded,densityNew);
4877    }
4878    // Root node or every so often - see what to turn off
4879    int i ;
4880    for (i = 0;i<numberCutGenerators_;i++) {
4881      int howOften = generator_[i]->howOften() ;
4882      if (howOften>-90) 
4883        willBeCutsInTree=0;
4884    }
4885    if (!numberNodes_)
4886      handler_->message(CBC_ROOT,messages_)
4887        <<numberNewCuts_
4888        <<startObjective<<thisObjective
4889        <<currentPassNumber_
4890        <<CoinMessageEol ;
4891    int * count = new int[numberCutGenerators_] ;
4892    memset(count,0,numberCutGenerators_*sizeof(int)) ;
4893    int numberActiveGenerators=0;
4894    for (i = 0;i<numberNewCuts_;i++) {
4895      int iGenerator = whichGenerator_[i];
4896      if (iGenerator>=0&&iGenerator<numberCutGenerators_)
4897        count[iGenerator]++ ;
4898    }
4899    double totalCuts = 0.0 ;
4900    //#define JUST_ACTIVE
4901    for (i = 0;i<numberCutGenerators_;i++) {
4902      if (countRowCuts[i]||countColumnCuts[i])
4903        numberActiveGenerators++;
4904#ifdef JUST_ACTIVE
4905      totalCuts += count[i] + 5.0*countColumnCuts[i] ;
4906#else
4907      totalCuts += countRowCuts[i] + 5.0*countColumnCuts[i] ;
4908#endif
4909    }
4910    int iProbing=-1;
4911    double small = (0.5* totalCuts) /
4912      ((double) numberActiveGenerators) ;
4913    for (i = 0;i<numberCutGenerators_;i++) {
4914      int howOften = generator_[i]->howOften() ;
4915      /*  Probing can be set to just do column cuts in treee.
4916          But if doing good then leave as on */
4917      bool probingWasOnBut=false;
4918      CglProbing * probing = dynamic_cast<CglProbing*>(generator_[i]->generator());
4919      if (probing&&!numberNodes_) {
4920        iProbing = i;
4921        if (probing->rowCuts()==-3) {
4922          probingWasOnBut=true;
4923          howOften = -98;
4924          probing->setRowCuts(3);
4925        }
4926      }
4927      if (willBeCutsInTree<0&&howOften==-98)
4928        howOften =-99;
4929      if (howOften==-98&&generator_[i]->switchOffIfLessThan()<0) {
4930        if (thisObjective-startObjective<0.005*fabs(startObjective)+1.0e-5)
4931          howOften=-99; // switch off
4932        if (thisObjective-startObjective<0.1*fabs(startObjective)+1.0e-5
4933            &&5*solver_->getNumRows()<solver_->getNumCols())
4934          howOften=-99; // switch off
4935      }
4936      if (howOften<-99)
4937        continue ;
4938      if (howOften<0||howOften >= 1000000) {
4939        if( !numberNodes_) {
4940          // If small number switch mostly off
4941#ifdef JUST_ACTIVE
4942          double thisCuts = count[i] + 5.0*countColumnCuts[i] ;
4943#else
4944          double thisCuts = countRowCuts[i] + 5.0*countColumnCuts[i] ;
4945#endif
4946          if (!thisCuts||howOften == -99) {
4947            if (howOften == -99||howOften == -98) 
4948              howOften = -100 ;
4949            else
4950              howOften = 1000000+SCANCUTS; // wait until next time
4951          } else if (thisCuts<small) {
4952            if (howOften!=1&&!probingWasOnBut) {
4953              if (generator_[i]->whatDepth()<0||howOften!=-1) {
4954                int k = (int) sqrt(small/thisCuts) ;
4955                if (howOften!=-98)
4956                  howOften = k+1000000 ;
4957                else
4958                  howOften=-100;
4959              } else {
4960                howOften=1;
4961              }
4962            } else {
4963              howOften=1;
4964              // allow cuts
4965              probingWasOnBut=false;
4966            }
4967          } else {
4968            if (thisObjective-startObjective<0.1*fabs(startObjective)+1.0e-5&&generator_[i]->whatDepth()<0)
4969              generator_[i]->setWhatDepth(5);
4970            howOften = 1+1000000 ;
4971          }
4972        }
4973        // If cuts useless switch off
4974        if (numberNodes_>=100000&&sumChangeObjective1_>2.0e2*(sumChangeObjective2_+1.0e-12)) {
4975          howOften = 1000000+SCANCUTS; // wait until next time
4976          //printf("switch off cut %d due to lack of use\n",i);
4977        }
4978      }
4979      if (!numberNodes_) {
4980        if (probingWasOnBut&&howOften==-100) {
4981          probing->setRowCuts(-3);
4982          howOften=1;
4983        }
4984        if (howOften==1)
4985          generator_[i]->setWhatDepth(1);
4986       
4987        if (howOften>=0&&generator_[i]->generator()->mayGenerateRowCutsInTree())
4988          willBeCutsInTree=1;
4989      }
4990      generator_[i]->setHowOften(howOften) ;
4991      if (howOften>=1000000&&howOften<2000000&&0) {
4992        // Go to depth
4993        int bias=1;
4994        if (howOften==1+1000000)
4995          generator_[i]->setWhatDepth(bias+1);
4996        else if (howOften<=10+1000000)
4997          generator_[i]->setWhatDepth(bias+2);
4998        else
4999          generator_[i]->setWhatDepth(bias+1000);
5000      }
5001      int newFrequency = generator_[i]->howOften()%1000000 ;
5002      // increment cut counts
5003      generator_[i]->incrementNumberCutsInTotal(countRowCuts[i]);
5004      generator_[i]->incrementNumberCutsActive(count[i]);
5005      if (handler_->logLevel()>1||!numberNodes_) {
5006        handler_->message(CBC_GENERATOR,messages_)
5007          <<i
5008          <<generator_[i]->cutGeneratorName()
5009          <<countRowCuts[i]<<count[i]
5010          <<countColumnCuts[i];
5011        handler_->printing(!numberNodes_&&generator_[i]->timing())
5012          <<generator_[i]->timeInCutGenerator();
5013        handler_->message()
5014          <<newFrequency
5015          <<CoinMessageEol ;
5016      }
5017    } 
5018    delete [] count ;
5019    if( !numberNodes_) {
5020      // save statistics
5021      for (i = 0;i<numberCutGenerators_;i++) {
5022        generator_[i]->setNumberCutsAtRoot(generator_[i]->numberCutsInTotal());
5023        generator_[i]->setNumberActiveCutsAtRoot(generator_[i]->numberCutsActive());
5024      }
5025      // decide on pseudo cost strategy
5026      int howOften = iProbing>=0 ? generator_[iProbing]->howOften() : 0;
5027      if ((howOften %1000000)!=1) 
5028        howOften = 0;
5029      //if (howOften) {
5030      //CglProbing * probing = dynamic_cast<CglProbing*>(generator_[iProbing]->generator());
5031      //}
5032      howOften = 0;
5033      if (howOften) {
5034        printf("** method 1\n");
5035        //CglProbing * probing = dynamic_cast<CglProbing*>(generator_[iProbing]->generator());
5036        generator_[iProbing]->setWhatDepth(1);
5037        // could set no row cuts
5038        //if (thisObjective-startObjective<0.001*fabs(startObjective)+1.0e-5)
5039        // probing->setRowCuts(0);
5040        for (int i=0;i<numberObjects_;i++) {
5041          CbcSimpleIntegerDynamicPseudoCost * obj =
5042            dynamic_cast <CbcSimpleIntegerDynamicPseudoCost *>(object_[i]) ;
5043          if (obj) 
5044            obj->setMethod(1);
5045        }
5046      }
5047      if (willBeCutsInTree==-2)
5048        willBeCutsInTree=0;
5049      if( willBeCutsInTree<=0) {
5050        // Take off cuts
5051        cuts = OsiCuts();
5052        numberNewCuts_=0;
5053        if (!willBeCutsInTree) {
5054          // update size of problem
5055          numberRowsAtContinuous_ = solver_->getNumRows() ;
5056        } else {
5057          // take off cuts
5058          int numberRows = solver_->getNumRows();
5059          int numberAdded = numberRows-numberRowsAtContinuous_;
5060          if (numberAdded) {
5061            int * added = new int[numberAdded];
5062            for (int i=0;i<numberAdded;i++)
5063              added[i]=i+numberRowsAtContinuous_;
5064            solver_->deleteRows(numberAdded,added);
5065            delete [] added;
5066            // resolve so optimal
5067            resolve(solver_);
5068          }
5069        }
5070#ifdef COIN_HAS_CLP
5071        OsiClpSolverInterface * clpSolver
5072          = dynamic_cast<OsiClpSolverInterface *> (solver_);
5073        if (clpSolver) {
5074          // Maybe solver might like to know only column bounds will change
5075          //int options = clpSolver->specialOptions();
5076          //clpSolver->setSpecialOptions(options|128);
5077          clpSolver->synchronizeModel();
5078        }
5079#endif
5080      } else {
5081#ifdef COIN_HAS_CLP
5082        OsiClpSolverInterface * clpSolver
5083          = dynamic_cast<OsiClpSolverInterface *> (solver_);
5084        if (clpSolver) {
5085        // make sure factorization can't carry over
5086          int options = clpSolver->specialOptions();
5087          clpSolver->setSpecialOptions(options&(~8));
5088        }
5089#endif
5090      }
5091    }
5092  } else {
5093#ifdef COIN_HAS_CLP
5094    OsiClpSolverInterface * clpSolver
5095      = dynamic_cast<OsiClpSolverInterface *> (solver_);
5096    if (clpSolver) {
5097      // Maybe solver might like to know only column bounds will change
5098      //int options = clpSolver->specialOptions();
5099      //clpSolver->setSpecialOptions(options|128);
5100      clpSolver->synchronizeModel();
5101    }
5102#endif
5103    if (numberCutGenerators_) {
5104      int i;
5105      // add to counts anyway
5106      for (i = 0;i<numberCutGenerators_;i++) 
5107        generator_[i]->incrementNumberCutsInTotal(countRowCuts[i]);
5108      // What if not feasible as cuts may have helped
5109      if (feasible) {
5110        for (i = 0;i<numberNewCuts_;i++) {
5111          int iGenerator = whichGenerator_[i];
5112          if (iGenerator>=0)
5113            generator_[iGenerator]->incrementNumberCutsActive();
5114        }
5115      }
5116    }
5117  }
5118
5119  delete [] countRowCuts ;
5120  delete [] countColumnCuts ;
5121
5122
5123#ifdef CHECK_CUT_COUNTS
5124  if (feasible)
5125  { 
5126    CoinWarmStartBasis * basis = dynamic_cast<CoinWarmStartBasis*>(solver_->getWarmStart()) ;
5127    printf("solveWithCuts: Number of rows at end (only active cuts) %d\n",
5128           numberRowsAtContinuous_+numberNewCuts_+numberOldActiveCuts_) ;
5129    basis->print() ; delete basis;}
5130#endif
5131#ifdef CBC_DEBUG
5132  if (onOptimalPath && !solver_->isDualObjectiveLimitReached())
5133    assert(feasible) ;
5134#endif
5135# ifdef COIN_HAS_CLP
5136  if (clpSolver) 
5137    clpSolver->setSpecialOptions(saveClpOptions);
5138# endif
5139
5140  return feasible ; }
5141
5142
5143/*
5144  Remove slack cuts. We obtain a basis and scan it. Cuts with basic slacks
5145  are purged. If any cuts are purged, resolve() is called to restore the
5146  solution held in the solver.  If resolve() pivots, there's the possibility
5147  that a slack may be pivoted in (trust me :-), so the process iterates.
5148  Setting allowResolve to false will suppress reoptimisation (but see note
5149  below).
5150
5151  At the level of the solver's constraint system, loose cuts are really
5152  deleted.  There's an implicit assumption that deleteRows will also update
5153  the active basis in the solver.
5154
5155  At the level of nodes and models, it's more complicated.
5156
5157  New cuts exist only in the collection of cuts passed as a parameter. They
5158  are deleted from the collection and that's the end of them.
5159
5160  Older cuts have made it into addedCuts_. Two separate actions are needed.
5161  The reference count for the CbcCountRowCut object is decremented. If this
5162  count falls to 0, the node which owns the cut is located, the reference to
5163  the cut is removed, and then the cut object is destroyed (courtesy of the
5164  CbcCountRowCut destructor). We also need to set the addedCuts_ entry to
5165  NULL. This is important so that when it comes time to generate basis edits
5166  we can tell this cut was dropped from the basis during processing of the
5167  node.
5168
5169  NOTE: In general, it's necessary to call resolve() after purging slack
5170        cuts.  Deleting constraints constitutes a change in the problem, and
5171        an OSI is not required to maintain a valid solution when the problem
5172        is changed. But ... it's really useful to maintain the active basis,
5173        and the OSI is supposed to do that. (Yes, it's splitting hairs.) In
5174        some places, it's possible to know that the solution will never be
5175        consulted after this call, only the basis.  (E.g., this routine is
5176        called as a last act before generating info to place the node in the
5177        live set.) For such use, set allowResolve to false.
5178 
5179  TODO: No real harm would be done if we just ignored the rare occasion when
5180        the call to resolve() pivoted a slack back into the basis. It's a
5181        minor inefficiency, at worst. But it does break assertions which
5182        check that there are no loose cuts in the basis. It might be better
5183        to remove the assertions.
5184*/
5185
5186void
5187CbcModel::takeOffCuts (OsiCuts &newCuts,
5188                       bool allowResolve, OsiCuts * saveCuts)
5189
5190{ // int resolveIterations = 0 ;
5191  int firstOldCut = numberRowsAtContinuous_ ;
5192  int totalNumberCuts = numberNewCuts_+numberOldActiveCuts_ ;
5193  int *solverCutIndices = new int[totalNumberCuts] ;
5194  int *newCutIndices = new int[numberNewCuts_] ;
5195  const CoinWarmStartBasis* ws ;
5196  CoinWarmStartBasis::Status status ;
5197  bool needPurge = true ;
5198/*
5199  The outer loop allows repetition of purge in the event that reoptimisation
5200  changes the basis. To start an iteration, clear the deletion counts and grab
5201  the current basis.
5202*/
5203  while (needPurge)
5204  { int numberNewToDelete = 0 ;
5205    int numberOldToDelete = 0 ;
5206    int i ;
5207    ws = dynamic_cast<const CoinWarmStartBasis*>(solver_->getWarmStart()) ;
5208/*
5209  Scan the basis entries of the old cuts generated prior to this round of cut
5210  generation.  Loose cuts are `removed' by decrementing their reference count
5211  and setting the addedCuts_ entry to NULL. (If the reference count falls to
5212  0, they're really deleted.  See CbcModel and CbcCountRowCut doc'n for
5213  principles of cut handling.)
5214*/
5215    int oldCutIndex = 0 ;
5216    for (i = 0 ; i < numberOldActiveCuts_ ; i++)
5217    { status = ws->getArtifStatus(i+firstOldCut) ;
5218      while (!addedCuts_[oldCutIndex]) oldCutIndex++ ;
5219      assert(oldCutIndex < currentNumberCuts_) ;
5220      // always leave if from nextRowCut_
5221      if (status == CoinWarmStartBasis::basic&&
5222          addedCuts_[oldCutIndex]->effectiveness()!=COIN_DBL_MAX)
5223      { solverCutIndices[numberOldToDelete++] = i+firstOldCut ;
5224        if (saveCuts) {
5225          // send to cut pool
5226          OsiRowCut * slackCut = addedCuts_[oldCutIndex];
5227          if (slackCut->effectiveness()!=-1.234) {
5228            slackCut->setEffectiveness(-1.234);
5229            saveCuts->insert(*slackCut);
5230          }
5231        }
5232        if (addedCuts_[oldCutIndex]->decrement() == 0)
5233          delete addedCuts_[oldCutIndex] ;
5234        addedCuts_[oldCutIndex] = NULL ;
5235        oldCutIndex++ ; }
5236      else
5237      { oldCutIndex++ ; } }
5238/*
5239  Scan the basis entries of the new cuts generated with this round of cut
5240  generation.  At this point, newCuts is the only record of the new cuts, so
5241  when we delete loose cuts from newCuts, they're really gone. newCuts is a
5242  vector, so it's most efficient to compress it (eraseRowCut) from back to
5243  front.
5244*/
5245    int firstNewCut = firstOldCut+numberOldActiveCuts_ ;
5246    int k = 0 ;
5247    for (i = 0 ; i < numberNewCuts_ ; i++)
5248    { status = ws->getArtifStatus(i+firstNewCut) ;
5249      if (status == CoinWarmStartBasis::basic&&whichGenerator_[i]!=-2)
5250      { solverCutIndices[numberNewToDelete+numberOldToDelete] = i+firstNewCut ;
5251        newCutIndices[numberNewToDelete++] = i ; }
5252      else
5253      { // save which generator did it
5254        whichGenerator_[k++] = whichGenerator_[i] ; } }
5255    delete ws ;
5256    for (i = numberNewToDelete-1 ; i >= 0 ; i--)
5257    { int iCut = newCutIndices[i] ;
5258      if (saveCuts) {
5259        // send to cut pool
5260        OsiRowCut * slackCut = newCuts.rowCutPtr(iCut);
5261        if (slackCut->effectiveness()!=-1.234) {
5262          slackCut->setEffectiveness(-1.234);
5263          saveCuts->insert(*slackCut);
5264        }
5265      }
5266      newCuts.eraseRowCut(iCut) ; }
5267/*
5268  Did we delete anything? If so, delete the cuts from the constraint system
5269  held in the solver and reoptimise unless we're forbidden to do so. If the
5270  call to resolve() results in pivots, there's the possibility we again have
5271  basic slacks. Repeat the purging loop.
5272*/
5273    if (numberNewToDelete+numberOldToDelete > 0)
5274    { solver_->deleteRows(numberNewToDelete+numberOldToDelete,
5275                          solverCutIndices) ;
5276      numberNewCuts_ -= numberNewToDelete ;
5277      numberOldActiveCuts_ -= numberOldToDelete ;
5278#     ifdef CBC_DEBUG
5279      printf("takeOffCuts: purged %d+%d cuts\n", numberOldToDelete,
5280             numberNewToDelete );
5281#     endif
5282      if (allowResolve)
5283      { 
5284        phase_=3;
5285        // can do quick optimality check
5286        int easy=2;
5287        solver_->setHintParam(OsiDoInBranchAndCut,true,OsiHintDo,&easy) ;
5288        resolve(solver_) ;
5289        setPointers(solver_);
5290        solver_->setHintParam(OsiDoInBranchAndCut,true,OsiHintDo,NULL) ;
5291        if (solver_->getIterationCount() == 0)
5292        { needPurge = false ; }
5293#       ifdef CBC_DEBUG
5294        else
5295          { printf( "Repeating purging loop. %d iters.\n",
5296                    solver_->getIterationCount()); }
5297#       endif
5298      }
5299      else
5300      { needPurge = false ; } }
5301    else
5302    { needPurge = false ; } }
5303/*
5304  Clean up and return.
5305*/
5306  delete [] solverCutIndices ;
5307  delete [] newCutIndices ;
5308}
5309
5310
5311
5312int
5313CbcModel::resolve(CbcNodeInfo * parent, int whereFrom)
5314{
5315  // We may have deliberately added in violated cuts - check to avoid message
5316  int iRow;
5317  int numberRows = solver_->getNumRows();
5318  const double * rowLower = solver_->getRowLower();
5319  const double * rowUpper = solver_->getRowUpper();
5320  bool feasible=true;
5321  for (iRow= numberRowsAtContinuous_;iRow<numberRows;iRow++) {
5322    if (rowLower[iRow]>rowUpper[iRow]+1.0e-8)
5323      feasible=false;
5324  }
5325  // Can't happen if strong branching as would have been found before
5326  if (!numberStrong_&&numberObjects_>numberIntegers_) {
5327    int iColumn;
5328    int numberColumns = solver_->getNumCols();
5329    const double * columnLower = solver_->getColLower();
5330    const double * columnUpper = solver_->getColUpper();
5331    for (iColumn= 0;iColumn<numberColumns;iColumn++) {
5332      if (columnLower[iColumn]>columnUpper[iColumn]+1.0e-5)
5333        feasible=false;
5334    }
5335  }
5336#if 0
5337  {
5338    int iColumn;
5339    int numberColumns = solver_->getNumCols();
5340    const double * columnLower = solver_->getColLower();
5341    const double * columnUpper = solver_->getColUpper();
5342    const double * sol = solver_->getColSolution();
5343    bool feasible=true;
5344    double xx[37];
5345    memset(xx,0,sizeof(xx));
5346    xx[6]=xx[7]=xx[9]=xx[18]=xx[26]=xx[36]=1.0;
5347    for (iColumn= 0;iColumn<numberColumns;iColumn++) {
5348      if (solver_->isInteger(iColumn)) {
5349        //printf("yy %d at %g (bounds %g, %g)\n",iColumn,
5350        //     sol[iColumn],columnLower[iColumn],columnUpper[iColumn]);
5351        if (columnLower[iColumn]>xx[iColumn]||
5352            columnUpper[iColumn]<xx[iColumn])
5353          feasible=false;
5354        //solver_->setColLower(iColumn,CoinMin(xx[iColumn],columnUpper[iColumn]));
5355        //solver_->setColUpper(iColumn,CoinMax(xx[iColumn],columnLower[iColumn]));
5356      }
5357    }
5358    if (feasible)
5359      printf("On Feasible path\n");
5360  }
5361#endif
5362/*
5363  Reoptimize. Consider the possibility that we should fathom on bounds. But be
5364  careful --- where the objective takes on integral values, we may want to keep
5365  a solution where the objective is right on the cutoff.
5366*/
5367  if (feasible)
5368    {
5369      resolve(solver_) ;
5370      numberIterations_ += solver_->getIterationCount() ;
5371      feasible = (solver_->isProvenOptimal() &&
5372                  !solver_->isDualObjectiveLimitReached()) ;
5373    }
5374  if (0&&feasible) {
5375    const double * lb = solver_->getColLower();
5376    const double * ub = solver_->getColUpper();
5377    const double * x = solver_->getColSolution();
5378    const double * dj = solver_->getReducedCost();
5379    int numberColumns = solver_->getNumCols();
5380    for (int i=0;i<numberColumns;i++) {
5381      if (dj[i]>1.0e-4&&ub[i]-lb[i]>1.0e-4&&x[i]>lb[i]+1.0e-4)
5382        printf("error %d %g %g %g %g\n",i,dj[i],lb[i],x[i],ub[i]);
5383      if (dj[i]<-1.0e-4&&ub[i]-lb[i]>1.0e-4&&x[i]<ub[i]-1.0e-4)
5384        printf("error %d %g %g %g %g\n",i,dj[i],lb[i],x[i],ub[i]);
5385    }
5386  } 
5387  if (false&&!feasible&& continuousObjective_ <-1.0e30) {
5388    // at root node - double double check
5389    bool saveTakeHint;
5390    OsiHintStrength saveStrength;
5391    solver_->getHintParam(OsiDoDualInResolve,saveTakeHint,saveStrength);
5392    if (saveTakeHint||saveStrength==OsiHintIgnore) {
5393      solver_->setHintParam(OsiDoDualInResolve,false,OsiHintDo) ;
5394      resolve(solver_);
5395      solver_->setHintParam(OsiDoDualInResolve,saveTakeHint,saveStrength);
5396      numberIterations_ += solver_->getIterationCount() ;
5397      feasible = solver_->isProvenOptimal();
5398      //      solver_->writeMps("infeas");
5399    }
5400  }
5401  if (cutModifier_&&feasible&&!solverCharacteristics_->solutionAddsCuts()) {
5402    //double increment = getDblParam(CbcModel::CbcCutoffIncrement) ;
5403    double cutoff ;
5404    solver_->getDblParam(OsiDualObjectiveLimit,cutoff) ;
5405    double distance = fabs(cutoff-solver_->getObjValue());
5406    if (distance<10.0*trueIncrement) {
5407      double offset;
5408      solver_->getDblParam(OsiObjOffset,offset);
5409      double objFixedValue = -offset;
5410      double objValue=0.0;
5411      double direction = solver_->getObjSense();
5412      const double * solution = solver_->getColSolution();
5413      const double * objective = solver_->getObjCoefficients();
5414      const double * columnLower = solver_->getColLower();
5415      const double * columnUpper = solver_->getColUpper();
5416      int numberColumns = solver_->getNumCols();
5417      int increment = 0 ;
5418      double multiplier = 1.0/trueIncrement;
5419      int bigIntegers = 0; // Count of large costs which are integer
5420      for (int iColumn=0;iColumn<numberColumns;iColumn++) {
5421        double value = solution[iColumn];
5422        // make sure clean
5423        value = CoinMin(value,columnUpper[iColumn]);
5424        value = CoinMax(value,columnLower[iColumn]);
5425        double cost = direction * objective[iColumn];
5426        if (cost) {
5427          if (columnLower[iColumn]<columnUpper[iColumn]) {
5428            objValue += value*cost;
5429            value = fabs(cost)*multiplier ;
5430            if (value <2.1e9) {
5431              int nearest = (int) floor(value+0.5) ;
5432              assert (fabs(value-floor(value+0.5)) < 1.0e-8);
5433              if (!increment) 
5434                increment = nearest ; 
5435              else
5436                increment = gcd(increment,nearest) ;
5437            } else {
5438              // large value - may still be multiple of 1.0
5439              value = fabs(objective[iColumn]);
5440              assert(fabs(value-floor(value+0.5)) < 1.0e-8);
5441              bigIntegers++;
5442            }
5443          } else {
5444            // fixed
5445            objFixedValue += value*cost;
5446          }
5447        }
5448      }
5449      if (increment) {
5450        double value = increment ;
5451        value /= multiplier ;
5452        if (value>trueIncrement) {
5453          double x = objValue/value;
5454          x = ceil(x-1.0e-5);
5455          x *= value;
5456          //printf("fixed %g, variable %g -> %g, sum %g - cutoff %g\n",
5457          // objFixedValue,objValue,x,x+objFixedValue,cutoff);
5458          x += objFixedValue;
5459          if (x>cutoff + 1.0e-5*fabs(cutoff)+1.0e-5) {
5460            //printf("Node cutoff\n");
5461            feasible=false;
5462          }
5463        } else {
5464          value = trueIncrement;
5465          double x = objValue/value;
5466          x = ceil(x-1.0e-5);
5467          x *= value;
5468          x += objFixedValue;
5469          if (x>cutoff + 1.0e-5*fabs(cutoff)+1.0e-5) {
5470            //printf("Node cutoff\n");
5471            feasible=false;
5472          }
5473        }
5474      }
5475    }
5476  }
5477
5478  setPointers(solver_);
5479  int returnStatus = feasible ? 1 : 0;
5480  if (strategy_) {
5481    // user can play clever tricks here
5482    int status = strategy_->status(this,parent,whereFrom);
5483    if (status>=0) {
5484      if (status==0)
5485        returnStatus = 1;
5486      else if(status==1)
5487        returnStatus=-1;
5488      else
5489        returnStatus=0;
5490    }
5491  }
5492  return returnStatus ;
5493}
5494
5495
5496/* Set up objects.  Only do ones whose length is in range.
5497   If makeEquality true then a new model may be returned if
5498   modifications had to be made, otherwise "this" is returned.
5499
5500   Could use Probing at continuous to extend objects
5501*/
5502CbcModel * 
5503CbcModel::findCliques(bool makeEquality,
5504                      int atLeastThisMany, int lessThanThis, int defaultValue)
5505{
5506  // No objects are allowed to exist
5507  assert(numberObjects_==numberIntegers_||!numberObjects_);
5508  CoinPackedMatrix matrixByRow(*solver_->getMatrixByRow());
5509  int numberRows = solver_->getNumRows();
5510  int numberColumns = solver_->getNumCols();
5511
5512  // We may want to add columns
5513  int numberSlacks=0;
5514  int * rows = new int[numberRows];
5515  double * element =new double[numberRows];
5516
5517  int iRow;
5518
5519  findIntegers(true);
5520  numberObjects_=numberIntegers_;
5521
5522  int numberCliques=0;
5523  OsiObject ** object = new OsiObject * [numberRows];
5524  int * which = new int[numberIntegers_];
5525  char * type = new char[numberIntegers_];
5526  int * lookup = new int[numberColumns];
5527  int i;
5528  for (i=0;i<numberColumns;i++) 
5529    lookup[i]=-1;
5530  for (i=0;i<numberIntegers_;i++) 
5531    lookup[integerVariable_[i]]=i;
5532
5533  // Statistics
5534  int totalP1=0,totalM1=0;
5535  int numberBig=0,totalBig=0;
5536  int numberFixed=0;
5537
5538  // Row copy
5539  const double * elementByRow = matrixByRow.getElements();
5540  const int * column = matrixByRow.getIndices();
5541  const CoinBigIndex * rowStart = matrixByRow.getVectorStarts();
5542  const int * rowLength = matrixByRow.getVectorLengths();
5543
5544  // Column lengths for slacks
5545  const int * columnLength = solver_->getMatrixByCol()->getVectorLengths();
5546
5547  const double * lower = getColLower();
5548  const double * upper = getColUpper();
5549  const double * rowLower = getRowLower();
5550  const double * rowUpper = getRowUpper();
5551
5552  for (iRow=0;iRow<numberRows;iRow++) {
5553    int numberP1=0, numberM1=0;
5554    int j;
5555    double upperValue=rowUpper[iRow];
5556    double lowerValue=rowLower[iRow];
5557    bool good=true;
5558    int slack = -1;
5559    for (j=rowStart[iRow];j<rowStart[iRow]+rowLength[iRow];j++) {
5560      int iColumn = column[j];
5561      int iInteger=lookup[iColumn];
5562      if (upper[iColumn]-lower[iColumn]<1.0e-8) {
5563        // fixed
5564        upperValue -= lower[iColumn]*elementByRow[j];
5565        lowerValue -= lower[iColumn]*elementByRow[j];
5566        continue;
5567      } else if (upper[iColumn]!=1.0||lower[iColumn]!=0.0) {
5568        good = false;
5569        break;
5570      } else if (iInteger<0) {
5571        good = false;
5572        break;
5573      } else {
5574        if (columnLength[iColumn]==1)
5575          slack = iInteger;
5576      }
5577      if (fabs(elementByRow[j])!=1.0) {
5578        good=false;
5579        break;
5580      } else if (elementByRow[j]>0.0) {
5581        which[numberP1++]=iInteger;
5582      } else {
5583        numberM1++;
5584        which[numberIntegers_-numberM1]=iInteger;
5585      }
5586    }
5587    int iUpper = (int) floor(upperValue+1.0e-5);
5588    int iLower = (int) ceil(lowerValue-1.0e-5);
5589    int state=0;
5590    if (upperValue<1.0e6) {
5591      if (iUpper==1-numberM1)
5592        state=1;
5593      else if (iUpper==-numberM1)
5594        state=2;
5595      else if (iUpper<-numberM1)
5596        state=3;
5597    }
5598    if (!state&&lowerValue>-1.0e6) {
5599      if (-iLower==1-numberP1)
5600        state=-1;
5601      else if (-iLower==-numberP1)
5602        state=-2;
5603      else if (-iLower<-numberP1)
5604        state=-3;
5605    }
5606    if (good&&state) {
5607      if (abs(state)==3) {
5608        // infeasible
5609        numberObjects_=-1;
5610        break;
5611      } else if (abs(state)==2) {
5612        // we can fix all
5613        numberFixed += numberP1+numberM1;
5614        if (state>0) {
5615          // fix all +1 at 0, -1 at 1
5616          for (i=0;i<numberP1;i++)
5617            solver_->setColUpper(integerVariable_[which[i]],0.0);
5618          for (i=0;i<numberM1;i++)
5619            solver_->setColLower(integerVariable_[which[numberIntegers_-i-1]],
5620                                 1.0);
5621        } else {
5622          // fix all +1 at 1, -1 at 0
5623          for (i=0;i<numberP1;i++)
5624            solver_->setColLower(integerVariable_[which[i]],1.0);
5625          for (i=0;i<numberM1;i++)
5626            solver_->setColUpper(integerVariable_[which[numberIntegers_-i-1]],
5627                                 0.0);
5628        }
5629      } else {
5630        int length = numberP1+numberM1;
5631        if (length >= atLeastThisMany&&length<lessThanThis) {
5632          // create object
5633          bool addOne=false;
5634          int objectType;
5635          if (iLower==iUpper) {
5636            objectType=1;
5637          } else {
5638            if (makeEquality) {
5639              objectType=1;
5640              element[numberSlacks]=state;
5641              rows[numberSlacks++]=iRow;
5642              addOne=true;
5643            } else {
5644              objectType=0;
5645            }
5646          }
5647          if (state>0) {
5648            totalP1 += numberP1;
5649            totalM1 += numberM1;
5650            for (i=0;i<numberP1;i++)
5651              type[i]=1;
5652            for (i=0;i<numberM1;i++) {
5653              which[numberP1]=which[numberIntegers_-i-1];
5654              type[numberP1++]=0;
5655            }
5656          } else {
5657            totalP1 += numberM1;
5658            totalM1 += numberP1;
5659            for (i=0;i<numberP1;i++)
5660              type[i]=0;
5661            for (i=0;i<numberM1;i++) {
5662              which[numberP1]=which[numberIntegers_-i-1];
5663              type[numberP1++]=1;
5664            }
5665          }
5666          if (addOne) {
5667            // add in slack
5668            which[numberP1]=numberIntegers_+numberSlacks-1;
5669            slack = numberP1;
5670            type[numberP1++]=1;
5671          } else if (slack >= 0) {
5672            for (i=0;i<numberP1;i++) {
5673              if (which[i]==slack) {
5674                slack=i;
5675              }
5676            }
5677          }
5678          object[numberCliques] = new CbcClique(this,objectType,numberP1,
5679                                              which,type,
5680                                               1000000+numberCliques,slack);
5681          numberCliques++;
5682        } else if (numberP1+numberM1 >= lessThanThis) {
5683          // too big
5684          numberBig++;
5685          totalBig += numberP1+numberM1;
5686        }
5687      }
5688    }
5689  }
5690  delete [] which;
5691  delete [] type;
5692  delete [] lookup;
5693  if (numberCliques<0) {
5694    printf("*** Problem infeasible\n");
5695  } else {
5696    if (numberCliques)
5697      printf("%d cliques of average size %g found, %d P1, %d M1\n",
5698             numberCliques,
5699             ((double)(totalP1+totalM1))/((double) numberCliques),
5700             totalP1,totalM1);
5701    else
5702      printf("No cliques found\n");
5703    if (numberBig)
5704      printf("%d large cliques ( >= %d) found, total %d\n",
5705             numberBig,lessThanThis,totalBig);
5706    if (numberFixed)
5707      printf("%d variables fixed\n",numberFixed);
5708  }
5709  if (numberCliques>0&&numberSlacks&&makeEquality) {
5710    printf("adding %d integer slacks\n",numberSlacks);
5711    // add variables to make equality rows
5712    int * temp = new int[numberIntegers_+numberSlacks];
5713    memcpy(temp,integerVariable_,numberIntegers_*sizeof(int));
5714    // Get new model
5715    CbcModel * newModel = new CbcModel(*this);
5716    OsiSolverInterface * newSolver = newModel->solver();
5717    for (i=0;i<numberSlacks;i++) {
5718      temp[i+numberIntegers_]=i+numberColumns;
5719      int iRow = rows[i];
5720      double value = element[i];
5721      double lowerValue = 0.0;
5722      double upperValue = 1.0;
5723      double objValue  = 0.0;
5724      CoinPackedVector column(1,&iRow,&value);
5725      newSolver->addCol(column,lowerValue,upperValue,objValue);
5726      // set integer
5727      newSolver->setInteger(numberColumns+i);
5728      if (value >0)
5729        newSolver->setRowLower(iRow,rowUpper[iRow]);
5730      else
5731        newSolver->setRowUpper(iRow,rowLower[iRow]);
5732    }
5733    // replace list of integers
5734    for (i=0;i<newModel->numberObjects_;i++)
5735      delete newModel->object_[i];
5736    newModel->numberObjects_ = 0;
5737    delete [] newModel->object_;
5738    newModel->object_=NULL;
5739    newModel->findIntegers(true); //Set up all integer objects
5740    for (i=0;i<numberIntegers_;i++) {
5741      newModel->modifiableObject(i)->setPriority(object_[i]->priority());
5742    }
5743    if (originalColumns_) {
5744      // old model had originalColumns
5745      delete [] newModel->originalColumns_;
5746      newModel->originalColumns_ = new int[numberColumns+numberSlacks];
5747      memcpy(newModel->originalColumns_,originalColumns_,numberColumns*sizeof(int));
5748      // mark as not in previous model
5749      for (i=numberColumns;i<numberColumns+numberSlacks;i++)
5750        newModel->originalColumns_[i]=-1;
5751    }
5752    delete [] rows;
5753    delete [] element;
5754    newModel->addObjects(numberCliques,object);
5755    assert (ownObjects_);
5756    for (;i<numberCliques;i++) 
5757      delete object[i];
5758    delete [] object;
5759    newModel->synchronizeModel();
5760    return newModel;
5761  } else {
5762    assert (ownObjects_);
5763    if (numberCliques>0) {
5764      addObjects(numberCliques,object);
5765      for (;i<numberCliques;i++) 
5766        delete object[i];
5767      synchronizeModel();
5768    }
5769    delete [] object;
5770    delete [] rows;
5771    delete [] element;
5772    return this;
5773  }
5774}
5775// Fill in useful estimates
5776void 
5777CbcModel::pseudoShadow(double * down, double * up)
5778{
5779  // Column copy of matrix
5780  const double * element = solver_->getMatrixByCol()->getElements();
5781  const int * row = solver_->getMatrixByCol()->getIndices();
5782  const CoinBigIndex * columnStart = solver_->getMatrixByCol()->getVectorStarts();
5783  const int * columnLength = solver_->getMatrixByCol()->getVectorLengths();
5784  const double *objective = solver_->getObjCoefficients() ;
5785  int numberColumns = solver_->getNumCols() ;
5786  double direction = solver_->getObjSense();
5787  int iColumn;
5788  const double * dual = cbcRowPrice_;
5789  down = new double[numberColumns];
5790  up = new double[numberColumns];
5791  double upSum=1.0e-20;
5792  double downSum = 1.0e-20;
5793  int numberIntegers=0;
5794  for (iColumn=0;iColumn<numberColumns;iColumn++) {
5795    CoinBigIndex start = columnStart[iColumn];
5796    CoinBigIndex end = start + columnLength[iColumn];
5797    double upValue = 0.0;
5798    double downValue = 0.0;
5799    double value = direction*objective[iColumn];
5800    if (value) {
5801      if (value>0.0)
5802        upValue += value;
5803      else
5804        downValue -= value;
5805    }
5806    for (CoinBigIndex j=start;j<end;j++) {
5807      int iRow = row[j];
5808      value = -dual[iRow];
5809      if (value) {
5810        value *= element[j];
5811        if (value>0.0)
5812          upValue += value;
5813        else
5814          downValue -= value;
5815      }
5816    }
5817    // use dj if bigger
5818    double dj = cbcReducedCost_[iColumn];
5819    upValue = CoinMax(upValue,dj);
5820    downValue = CoinMax(downValue,-dj);
5821    up[iColumn]=upValue;
5822    down[iColumn]=downValue;
5823    if (solver_->isInteger(iColumn)) {
5824      if (!numberNodes_)
5825        printf("%d - dj %g up %g down %g cost %g\n",
5826               iColumn,dj,upValue,downValue,objective[iColumn]);
5827      upSum += upValue;
5828      downSum += downValue;
5829      numberIntegers++;
5830    }
5831  }
5832  if (numberIntegers) {
5833    double smallDown = 0.01*(downSum/((double) numberIntegers));
5834    double smallUp = 0.01*(upSum/((double) numberIntegers));
5835    for (int i=0;i<numberObjects_;i++) {
5836      CbcSimpleIntegerDynamicPseudoCost * obj1 =
5837        dynamic_cast <CbcSimpleIntegerDynamicPseudoCost *>(object_[i]) ;
5838      if (obj1) {
5839        iColumn = obj1->columnNumber();
5840        double upPseudoCost = obj1->upDynamicPseudoCost();
5841        double saveUp = upPseudoCost;
5842        upPseudoCost = CoinMax(upPseudoCost,smallUp);
5843        upPseudoCost = CoinMax(upPseudoCost,up[iColumn]);
5844        upPseudoCost = CoinMax(upPseudoCost,0.1*down[iColumn]);
5845        obj1->setUpDynamicPseudoCost(upPseudoCost);
5846        if (upPseudoCost>saveUp&&!numberNodes_)
5847          printf("For %d up went from %g to %g\n",
5848                 iColumn,saveUp,upPseudoCost);
5849        double downPseudoCost = obj1->downDynamicPseudoCost();
5850        double saveDown = downPseudoCost;
5851        downPseudoCost = CoinMax(downPseudoCost,smallDown);
5852        downPseudoCost = CoinMax(downPseudoCost,down[iColumn]);
5853        downPseudoCost = CoinMax(downPseudoCost,0.1*down[iColumn]);
5854        obj1->setDownDynamicPseudoCost(downPseudoCost);
5855        if (downPseudoCost>saveDown&&!numberNodes_)
5856          printf("For %d down went from %g to %g\n",
5857                 iColumn,saveDown,downPseudoCost);
5858      }
5859    }
5860  }
5861  delete [] down;
5862  delete [] up;
5863}
5864
5865/*
5866  Set branching priorities.
5867
5868  Setting integer priorities looks pretty robust; the call to findIntegers
5869  makes sure that SimpleInteger objects are in place. Setting priorities for
5870  other objects is entirely dependent on their existence, and the routine may
5871  quietly fail in several directions.
5872*/
5873
5874void 
5875CbcModel::passInPriorities (const int * priorities,
5876                            bool ifObject)
5877{
5878  findIntegers(false);
5879  int i;
5880  if (priorities) {
5881    int i0=0;
5882    int i1=numberObjects_-1;
5883    if (ifObject) {
5884      for (i=numberIntegers_;i<numberObjects_;i++) {
5885        object_[i]->setPriority(priorities[i-numberIntegers_]);
5886      }
5887      i0=numberIntegers_;
5888    } else {
5889      for (i=0;i<numberIntegers_;i++) {
5890        object_[i]->setPriority(priorities[i]);
5891      }
5892      i1=numberIntegers_-1;
5893    }
5894    messageHandler()->message(CBC_PRIORITY,
5895                              messages())
5896                                << i0<<i1<<numberObjects_ << CoinMessageEol ;
5897  }
5898}
5899
5900// Delete all object information
5901void 
5902CbcModel::deleteObjects(bool getIntegers)
5903{
5904  if (ownObjects_) {
5905    int i;
5906    for (i=0;i<numberObjects_;i++)
5907      delete object_[i];
5908    delete [] object_;
5909  }
5910  object_ =