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

Last change on this file since 607 was 607, checked in by lou, 13 years ago

Add names for heuristic, analogous to names for cut generators.

  • Property svn:eol-style set to native
  • Property svn:keywords set to Author Date Id Revision
File size: 303.4 KB
Line 
1// Copyright (C) 2002, International Business Machines
2// Corporation and others.  All Rights Reserved.
3#if defined(_MSC_VER)
4// Turn off compiler warning about long names
5#  pragma warning(disable:4786)
6#endif
7
8#include "CbcConfig.h"
9
10#include <string>
11//#define CBC_DEBUG 1
12//#define CHECK_CUT_COUNTS
13//#define CHECK_NODE_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, const char *name)
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  if (name)
3588  { heuristic_[numberHeuristics_]->setHeuristicName(name) ; }
3589  numberHeuristics_++ ;
3590}
3591
3592/*
3593  The last subproblem handled by the solver is not necessarily related to the
3594  one being recreated, so the first action is to remove all cuts from the
3595  constraint system.  Next, traverse the tree from node to the root to
3596  determine the basis size required for this subproblem and create an empty
3597  basis with the right capacity.  Finally, traverse the tree from root to
3598  node, adjusting bounds in the constraint system, adjusting the basis, and
3599  collecting the cuts that must be added to the constraint system.
3600  applyToModel does the heavy lifting.
3601
3602  addCuts1 is used in contexts where all that's desired is the list of cuts:
3603  the node is already fathomed, and we're collecting cuts so that we can
3604  adjust reference counts as we prune nodes. Arguably the two functions
3605  should be separated. The culprit is applyToModel, which performs cut
3606  collection and model adjustment.
3607
3608  Certainly in the contexts where all we need is a list of cuts, there's no
3609  point in passing in a valid basis --- an empty basis will do just fine.
3610*/
3611void CbcModel::addCuts1 (CbcNode * node, CoinWarmStartBasis *&lastws)
3612{ int i;
3613  int nNode=0;
3614  int numberColumns = getNumCols();
3615  CbcNodeInfo * nodeInfo = node->nodeInfo();
3616
3617/*
3618  Remove all cuts from the constraint system.
3619  (original comment includes ``see note below for later efficiency'', but
3620  the reference isn't clear to me).
3621*/
3622  int currentNumberCuts = solver_->getNumRows()-numberRowsAtContinuous_;
3623  int *which = new int[currentNumberCuts];
3624  for (i = 0 ; i < currentNumberCuts ; i++)
3625    which[i] = i+numberRowsAtContinuous_;
3626  solver_->deleteRows(currentNumberCuts,which);
3627  delete [] which;
3628/*
3629  Accumulate the path from node to the root in walkback_, and accumulate a
3630  cut count in currentNumberCuts.
3631
3632  original comment: when working then just unwind until where new node joins
3633  old node (for cuts?)
3634*/
3635  currentNumberCuts = 0;
3636  while (nodeInfo) {
3637    //printf("nNode = %d, nodeInfo = %x\n",nNode,nodeInfo);
3638    walkback_[nNode++]=nodeInfo;
3639    currentNumberCuts += nodeInfo->numberCuts() ;
3640    nodeInfo = nodeInfo->parent() ;
3641    if (nNode==maximumDepth_) {
3642      maximumDepth_ *= 2;
3643      CbcNodeInfo ** temp = new CbcNodeInfo * [maximumDepth_];
3644      for (i=0;i<nNode;i++) 
3645        temp[i] = walkback_[i];
3646      delete [] walkback_;
3647      walkback_ = temp;
3648    }
3649  }
3650/*
3651  Create an empty basis with sufficient capacity for the constraint system
3652  we'll construct: original system plus cuts. Make sure we have capacity to
3653  record those cuts in addedCuts_.
3654
3655  The method of adjusting the basis at a FullNodeInfo object (the root, for
3656  example) is to use a copy constructor to duplicate the basis held in the
3657  nodeInfo, then resize it and return the new basis object. Guaranteed,
3658  lastws will point to a different basis when it returns. We pass in a basis
3659  because we need the parameter to return the allocated basis, and it's an
3660  easy way to pass in the size. But we take a hit for memory allocation.
3661*/
3662  currentNumberCuts_=currentNumberCuts;
3663  if (currentNumberCuts >= maximumNumberCuts_) {
3664    maximumNumberCuts_ = currentNumberCuts;
3665    delete [] addedCuts_;
3666    addedCuts_ = new CbcCountRowCut * [maximumNumberCuts_];
3667  }
3668  lastws->setSize(numberColumns,numberRowsAtContinuous_+currentNumberCuts);
3669/*
3670  This last bit of code traverses the path collected in walkback_ from the
3671  root back to node. At the end of the loop,
3672   * lastws will be an appropriate basis for node;
3673   * variable bounds in the constraint system will be set to be correct for
3674     node; and
3675   * addedCuts_ will be set to a list of cuts that need to be added to the
3676     constraint system at node.
3677  applyToModel does all the heavy lifting.
3678*/
3679  currentNumberCuts=0;
3680  while (nNode) {
3681    --nNode;
3682    walkback_[nNode]->applyToModel(this,lastws,addedCuts_,currentNumberCuts);
3683  }
3684}
3685
3686/*
3687  adjustCuts might be a better name: If the node is feasible, we sift through
3688  the cuts collected by addCuts1, add the ones that are tight and omit the
3689  ones that are loose. If the node is infeasible, we just adjust the
3690  reference counts to reflect that we're about to prune this node and its
3691  descendants.
3692*/
3693int CbcModel::addCuts (CbcNode *node, CoinWarmStartBasis *&lastws,bool canFix)
3694{
3695/*
3696  addCuts1 performs step 1 of restoring the subproblem at this node; see the
3697  comments there.
3698*/
3699  addCuts1(node,lastws);
3700  int i;
3701  int numberColumns = getNumCols();
3702  CbcNodeInfo * nodeInfo = node->nodeInfo();
3703  double cutoff = getCutoff() ;
3704  int currentNumberCuts=currentNumberCuts_;
3705  if (canFix) {
3706    bool feasible=true;
3707    const double *lower = solver_->getColLower() ;
3708    const double *upper = solver_->getColUpper() ;
3709    double * newLower = analyzeResults_;
3710    double * objLower = newLower+numberIntegers_;
3711    double * newUpper = objLower+numberIntegers_;
3712    double * objUpper = newUpper+numberIntegers_;
3713    int n=0;
3714    for (i=0;i<numberIntegers_;i++) {
3715      int iColumn = integerVariable_[i];
3716      bool changed=false;
3717      double lo = 0.0;
3718      double up = 0.0;
3719      if (objLower[i]>cutoff) {
3720        lo = lower[iColumn];
3721        up = upper[iColumn];
3722        if (lo<newLower[i]) {
3723          lo = newLower[i];
3724          solver_->setColLower(iColumn,lo);
3725          changed=true;
3726          n++;
3727        }
3728        if (objUpper[i]>cutoff) {
3729          if (up>newUpper[i]) {
3730            up = newUpper[i];
3731            solver_->setColUpper(iColumn,up);
3732            changed=true;
3733            n++;
3734          }
3735        }
3736      } else if (objUpper[i]>cutoff) {
3737        lo = lower[iColumn];
3738        up = upper[iColumn];
3739        if (up>newUpper[i]) {
3740          up = newUpper[i];
3741          solver_->setColUpper(iColumn,up);
3742          changed=true;
3743          n++;
3744        }
3745      }
3746      if (changed&&lo>up) {
3747        feasible=false;
3748        break;
3749      }
3750    }
3751    if (!feasible) {
3752      printf("analysis says node infeas\n");
3753      cutoff=-COIN_DBL_MAX;
3754    }
3755  }
3756/*
3757  If the node can't be fathomed by bound, reinstall tight cuts in the
3758  constraint system. Even if there are no cuts, we'll want to set the
3759  reconstructed basis in the solver.
3760*/
3761  if (node->objectiveValue() < cutoff)
3762  { 
3763#   ifdef CBC_CHECK_BASIS
3764    printf("addCuts: expanded basis; rows %d+%d\n",
3765           numberRowsAtContinuous_,currentNumberCuts);
3766    lastws->print();
3767#   endif
3768/*
3769  Adjust the basis and constraint system so that we retain only active cuts.
3770  There are three steps:
3771    1) Scan the basis. Sort the cuts into effective cuts to be kept and
3772       loose cuts to be dropped.
3773    2) Drop the loose cuts and resize the basis to fit.
3774    3) Install the tight cuts in the constraint system (applyRowCuts) and
3775       and install the basis (setWarmStart).
3776  Use of compressRows conveys we're compressing the basis and not just
3777  tweaking the artificialStatus_ array.
3778*/
3779    if (currentNumberCuts > 0) {
3780      int numberToAdd = 0;
3781      const OsiRowCut **addCuts;
3782      int numberToDrop = 0 ;
3783      int *cutsToDrop ;
3784      addCuts = new const OsiRowCut* [currentNumberCuts];
3785      cutsToDrop = new int[currentNumberCuts] ;
3786      assert (currentNumberCuts+numberRowsAtContinuous_<=lastws->getNumArtificial());
3787      for (i=0;i<currentNumberCuts;i++) {
3788        CoinWarmStartBasis::Status status = 
3789          lastws->getArtifStatus(i+numberRowsAtContinuous_);
3790        if (addedCuts_[i] &&
3791            (status != CoinWarmStartBasis::basic ||
3792             addedCuts_[i]->effectiveness()==COIN_DBL_MAX)) {
3793#         ifdef CHECK_CUT_COUNTS
3794          printf("Using cut %d %x as row %d\n",i,addedCuts_[i],
3795                 numberRowsAtContinuous_+numberToAdd);
3796#         endif
3797          addCuts[numberToAdd++] = new OsiRowCut(*addedCuts_[i]);
3798        } else {
3799#         ifdef CHECK_CUT_COUNTS
3800          printf("Dropping cut %d %x\n",i,addedCuts_[i]);
3801#         endif
3802          addedCuts_[i]=NULL;
3803          cutsToDrop[numberToDrop++] = numberRowsAtContinuous_+i ;
3804        }
3805      }
3806      int numberRowsNow=numberRowsAtContinuous_+numberToAdd;
3807      lastws->compressRows(numberToDrop,cutsToDrop) ;
3808      lastws->resize(numberRowsNow,numberColumns);
3809      solver_->applyRowCuts(numberToAdd,addCuts);
3810#     ifdef CBC_CHECK_BASIS
3811      printf("addCuts: stripped basis; rows %d + %d\n",
3812             numberRowsAtContinuous_,numberToAdd);
3813      lastws->print();
3814#     endif
3815      for (i=0;i<numberToAdd;i++)
3816        delete addCuts[i];
3817      delete [] addCuts;
3818      delete [] cutsToDrop ;
3819    }
3820/*
3821  Set the basis in the solver.
3822*/
3823    solver_->setWarmStart(lastws);
3824#if 0
3825    if ((numberNodes_%printFrequency_)==0) {
3826      printf("Objective %g, depth %d, unsatisfied %d\n",
3827             node->objectiveValue(),
3828             node->depth(),node->numberUnsatisfied());
3829    }
3830#endif
3831/*
3832  Clean up and we're out of here.
3833*/
3834    numberNodes_++;
3835    return 0;
3836  } 
3837/*
3838  This node has been fathomed by bound as we try to revive it out of the live
3839  set. Adjust the cut reference counts to reflect that we no longer need to
3840  explore the remaining branch arms, hence they will no longer reference any
3841  cuts. Cuts whose reference count falls to zero are deleted.
3842*/
3843  else
3844  { int i;
3845    int numberLeft = nodeInfo->numberBranchesLeft();
3846    for (i = 0 ; i < currentNumberCuts ; i++)
3847    { if (addedCuts_[i])
3848      { if (!addedCuts_[i]->decrement(numberLeft))
3849        { delete addedCuts_[i];
3850          addedCuts_[i] = NULL; } } }
3851    return 1 ; }
3852}
3853
3854
3855/*
3856  Perform reduced cost fixing on integer variables.
3857
3858  The variables in question are already nonbasic at bound. We're just nailing
3859  down the current situation.
3860*/
3861
3862int CbcModel::reducedCostFix ()
3863
3864{
3865  if(!solverCharacteristics_->reducedCostsAccurate())
3866    return 0; //NLP
3867  double cutoff = getCutoff() ;
3868  double direction = solver_->getObjSense() ;
3869  double gap = cutoff - solver_->getObjValue()*direction ;
3870  double tolerance;
3871  solver_->getDblParam(OsiDualTolerance,tolerance) ;
3872  if (gap<=0.0)
3873    return 0;
3874  gap += 100.0*tolerance;
3875  double integerTolerance = getDblParam(CbcIntegerTolerance) ;
3876
3877  const double *lower = solver_->getColLower() ;
3878  const double *upper = solver_->getColUpper() ;
3879  const double *solution = solver_->getColSolution() ;
3880  const double *reducedCost = solver_->getReducedCost() ;
3881
3882  int numberFixed = 0 ;
3883
3884# ifdef COIN_HAS_CLP
3885  OsiClpSolverInterface * clpSolver
3886    = dynamic_cast<OsiClpSolverInterface *> (solver_);
3887  ClpSimplex * clpSimplex=NULL;
3888  if (clpSolver) 
3889    clpSimplex = clpSolver->getModelPtr();
3890# endif
3891  for (int i = 0 ; i < numberIntegers_ ; i++)
3892  { int iColumn = integerVariable_[i] ;
3893    double djValue = direction*reducedCost[iColumn] ;
3894    if (upper[iColumn]-lower[iColumn] > integerTolerance)
3895    { if (solution[iColumn] < lower[iColumn]+integerTolerance && djValue > gap)
3896      { solver_->setColUpper(iColumn,lower[iColumn]) ;
3897#ifdef COIN_HAS_CLP
3898      // may just have been fixed before
3899      if (clpSimplex)
3900        assert(clpSimplex->getColumnStatus(iColumn)==ClpSimplex::atLowerBound||
3901               clpSimplex->getColumnStatus(iColumn)==ClpSimplex::isFixed);
3902#endif
3903      numberFixed++ ; }
3904      else
3905      if (solution[iColumn] > upper[iColumn]-integerTolerance && -djValue > gap)
3906      { solver_->setColLower(iColumn,upper[iColumn]) ;
3907#ifdef COIN_HAS_CLP
3908      // may just have been fixed before
3909      if (clpSimplex)
3910        assert(clpSimplex->getColumnStatus(iColumn)==ClpSimplex::atUpperBound||
3911               clpSimplex->getColumnStatus(iColumn)==ClpSimplex::isFixed);
3912#endif
3913      numberFixed++ ; } } }
3914  numberDJFixed_ += numberFixed;
3915 
3916  return numberFixed; }
3917
3918// Collect coding to replace whichGenerator
3919void
3920CbcModel::resizeWhichGenerator(int numberNow, int numberAfter)
3921{
3922  if (numberAfter > maximumWhich_) {
3923    maximumWhich_ = CoinMax(maximumWhich_*2+100,numberAfter) ;
3924    int * temp = new int[2*maximumWhich_] ;
3925    memcpy(temp,whichGenerator_,numberNow*sizeof(int)) ;
3926    delete [] whichGenerator_ ;
3927    whichGenerator_ = temp ;
3928    memset(whichGenerator_+numberNow,0,(maximumWhich_-numberNow)*sizeof(int));
3929  }
3930}
3931
3932/** Solve the model using cuts
3933
3934  This version takes off redundant cuts from node.
3935  Returns true if feasible.
3936
3937  \todo
3938  Why do I need to resolve the problem? What has been done between the last
3939  relaxation and calling solveWithCuts?
3940
3941  If numberTries == 0 then user did not want any cuts.
3942*/
3943
3944bool 
3945CbcModel::solveWithCuts (OsiCuts &cuts, int numberTries, CbcNode *node)
3946/*
3947  Parameters:
3948    numberTries: (i) the maximum number of iterations for this round of cut
3949                     generation; if negative then we don't mind if drop is tiny.
3950   
3951    cuts:       (o) all cuts generated in this round of cut generation
3952
3953    node: (i)     So we can update dynamic pseudo costs
3954*/
3955                       
3956
3957{
3958# ifdef COIN_HAS_CLP
3959  OsiClpSolverInterface * clpSolver
3960    = dynamic_cast<OsiClpSolverInterface *> (solver_);
3961  int saveClpOptions=0;
3962  if (clpSolver) 
3963    saveClpOptions = clpSolver->specialOptions();
3964# endif
3965  bool feasible = true ;
3966  int lastNumberCuts = 0 ;
3967  double lastObjective = -1.0e100 ;
3968  int violated = 0 ;
3969  int numberRowsAtStart = solver_->getNumRows() ;
3970  //printf("solver had %d rows\n",numberRowsAtStart);
3971  int numberColumns = solver_->getNumCols() ;
3972  CoinBigIndex numberElementsAtStart = solver_->getNumElements();
3973
3974  numberOldActiveCuts_ = numberRowsAtStart-numberRowsAtContinuous_ ;
3975  numberNewCuts_ = 0 ;
3976
3977  bool onOptimalPath = false ;
3978  const OsiRowCutDebugger *debugger = NULL;
3979  if ((specialOptions_&1)!=0) {
3980    /*
3981      See OsiRowCutDebugger for details. In a nutshell, make sure that current
3982      variable values do not conflict with a known optimal solution. (Obviously
3983      this can be fooled when there are multiple solutions.)
3984    */
3985    debugger = solver_->getRowCutDebugger() ;
3986    if (debugger) 
3987      onOptimalPath = (debugger->onOptimalPath(*solver_)) ;
3988  }
3989  OsiCuts slackCuts;
3990/*
3991  Resolve the problem. If we've lost feasibility, might as well bail out right
3992  after the debug stuff. The resolve will also refresh cached copies of the
3993  solver solution (cbcColLower_, ...) held by CbcModel.
3994*/
3995  double objectiveValue = solver_->getObjValue()*solver_->getObjSense();
3996  if (node)
3997    objectiveValue= node->objectiveValue();
3998  int returnCode = resolve(node ? node->nodeInfo() : NULL,1);
3999#ifdef COIN_DEVELOP
4000  //if (!solver_->getIterationCount()&&solver_->isProvenOptimal())
4001  //printf("zero iterations on first solve of branch\n");
4002#endif
4003  if (node&&!node->nodeInfo()->numberBranchesLeft())
4004    node->nodeInfo()->allBranchesGone(); // can clean up
4005  feasible = returnCode  != 0 ;
4006  if (returnCode<0)
4007    numberTries=0;
4008  if (problemFeasibility_->feasible(this,0)<0) {
4009    feasible=false; // pretend infeasible
4010  }
4011 
4012  // Update branching information if wanted
4013  if(node &&branchingMethod_)
4014    branchingMethod_->updateInformation(solver_,node);
4015
4016#ifdef CBC_DEBUG
4017  if (feasible)
4018  { printf("Obj value %g (%s) %d rows\n",solver_->getObjValue(),
4019           (solver_->isProvenOptimal())?"proven":"unproven",
4020           solver_->getNumRows()) ; }
4021 
4022  else
4023  { printf("Infeasible %d rows\n",solver_->getNumRows()) ; }
4024#endif
4025  if ((specialOptions_&1)!=0) {
4026/*
4027  If the RowCutDebugger said we were compatible with the optimal solution,
4028  and now we're suddenly infeasible, we might be confused. Then again, we
4029  may have fathomed by bound, heading for a rediscovery of an optimal solution.
4030*/
4031    if (onOptimalPath && !solver_->isDualObjectiveLimitReached()) {
4032      if (!feasible) {
4033        solver_->writeMps("infeas");
4034        CoinWarmStartBasis *slack =
4035          dynamic_cast<CoinWarmStartBasis *>(solver_->getEmptyWarmStart()) ;
4036        solver_->setWarmStart(slack);
4037        delete slack ;
4038        solver_->setHintParam(OsiDoReducePrint,false,OsiHintDo,0) ;
4039        solver_->initialSolve();
4040      }
4041      assert(feasible) ;
4042    }
4043  }
4044
4045  if (!feasible) {
4046    numberInfeasibleNodes_++;
4047# ifdef COIN_HAS_CLP
4048    if (clpSolver) 
4049    clpSolver->setSpecialOptions(saveClpOptions);
4050# endif
4051    return (false) ;
4052  }
4053  sumChangeObjective1_ += solver_->getObjValue()*solver_->getObjSense()
4054    - objectiveValue ;
4055  if ( CoinCpuTime()-dblParam_[CbcStartSeconds] > dblParam_[CbcMaximumSeconds] )
4056    numberTries=0; // exit
4057  //if ((numberNodes_%100)==0)
4058  //printf("XXa sum obj changed by %g\n",sumChangeObjective1_);
4059  objectiveValue = solver_->getObjValue()*solver_->getObjSense();
4060  // Return at once if numberTries zero
4061  if (!numberTries) {
4062    cuts=OsiCuts();
4063    numberNewCuts_=0;
4064# ifdef COIN_HAS_CLP
4065    if (clpSolver) 
4066    clpSolver->setSpecialOptions(saveClpOptions);
4067# endif
4068    return true;
4069  }
4070/*
4071  Do reduced cost fixing.
4072*/
4073  reducedCostFix() ;
4074/*
4075  Set up for at most numberTries rounds of cut generation. If numberTries is
4076  negative, we'll ignore the minimumDrop_ cutoff and keep generating cuts for
4077  the specified number of rounds.
4078*/
4079  double minimumDrop = minimumDrop_ ;
4080  if (numberTries<0)
4081  { numberTries = -numberTries ;
4082    minimumDrop = -1.0 ; }
4083/*
4084  Is it time to scan the cuts in order to remove redundant cuts? If so, set
4085  up to do it.
4086*/
4087# define SCANCUTS 100 
4088  int *countColumnCuts = NULL ;
4089  // Always accumulate row cut counts
4090  int * countRowCuts =new int[numberCutGenerators_+numberHeuristics_] ;
4091  memset(countRowCuts,0,
4092         (numberCutGenerators_+numberHeuristics_)*sizeof(int)) ;
4093  bool fullScan = false ;
4094  if ((numberNodes_%SCANCUTS) == 0)
4095  { fullScan = true ;
4096    countColumnCuts = new int[numberCutGenerators_+numberHeuristics_] ;
4097    memset(countColumnCuts,0,
4098           (numberCutGenerators_+numberHeuristics_)*sizeof(int)) ; }
4099
4100  double direction = solver_->getObjSense() ;
4101  double startObjective = solver_->getObjValue()*direction ;
4102
4103  currentPassNumber_ = 0 ;
4104  double primalTolerance = 1.0e-7 ;
4105  // We may need to keep going on
4106  bool keepGoing=false;
4107/*
4108  Begin cut generation loop. Cuts generated during each iteration are
4109  collected in theseCuts. The loop can be divided into four phases:
4110   1) Prep: Fix variables using reduced cost. In the first iteration only,
4111      consider scanning globalCuts_ and activating any applicable cuts.
4112   2) Cut Generation: Call each generator and heuristic registered in the
4113      generator_ and heuristic_ arrays. Newly generated global cuts are
4114      copied to globalCuts_ at this time.
4115   3) Cut Installation and Reoptimisation: Install column and row cuts in
4116      the solver. Copy row cuts to cuts (parameter). Reoptimise.
4117   4) Cut Purging: takeOffCuts() removes inactive cuts from the solver, and
4118      does the necessary bookkeeping in the model.
4119*/
4120  do
4121  { currentPassNumber_++ ;
4122    numberTries-- ;
4123    if (numberTries<0&&keepGoing) {
4124      // switch off all normal ones
4125      for (int i = 0;i<numberCutGenerators_;i++) {
4126        if (!generator_[i]->mustCallAgain())
4127          generator_[i]->setSwitchedOff(true);
4128      }
4129    }
4130    keepGoing=false;
4131    OsiCuts theseCuts ;
4132/*
4133  Scan previously generated global column and row cuts to see if any are
4134  useful.
4135*/
4136    int numberViolated=0;
4137    if (currentPassNumber_ == 1 && howOftenGlobalScan_ > 0 &&
4138        (numberNodes_%howOftenGlobalScan_) == 0)
4139    { int numberCuts = globalCuts_.sizeColCuts() ;
4140      int i;
4141      // possibly extend whichGenerator
4142      resizeWhichGenerator(numberViolated, numberViolated+numberCuts);
4143      for ( i = 0 ; i < numberCuts ; i++)
4144      { OsiColCut *thisCut = globalCuts_.colCutPtr(i) ;
4145        if (thisCut->violated(cbcColSolution_)>primalTolerance) {
4146          printf("Global cut added - violation %g\n",
4147                 thisCut->violated(cbcColSolution_)) ;
4148          whichGenerator_[numberViolated++]=-1;
4149#ifndef GLOBAL_CUTS_JUST_POINTERS
4150          theseCuts.insert(*thisCut) ;
4151#else
4152          theseCuts.insert(thisCut) ;
4153#endif
4154        }
4155      }
4156      numberCuts = globalCuts_.sizeRowCuts() ;
4157      // possibly extend whichGenerator
4158      resizeWhichGenerator(numberViolated, numberViolated+numberCuts);
4159      for ( i = 0;i<numberCuts;i++) {
4160        OsiRowCut * thisCut = globalCuts_.rowCutPtr(i) ;
4161        if (thisCut->violated(cbcColSolution_)>primalTolerance) {
4162          //printf("Global cut added - violation %g\n",
4163          // thisCut->violated(cbcColSolution_)) ;
4164          whichGenerator_[numberViolated++]=-1;
4165#ifndef GLOBAL_CUTS_JUST_POINTERS
4166          theseCuts.insert(*thisCut) ;
4167#else
4168          theseCuts.insert(thisCut) ;
4169#endif
4170        }
4171      }
4172      numberGlobalViolations_+=numberViolated;
4173    }
4174/*
4175  Generate new cuts (global and/or local) and/or apply heuristics.  If
4176  CglProbing is used, then it should be first as it can fix continuous
4177  variables.
4178
4179  At present, CglProbing is the only case where generateCuts will return
4180  true. generateCuts actually modifies variable bounds in the solver when
4181  CglProbing indicates that it can fix a variable. Reoptimisation is required
4182  to take full advantage.
4183
4184  The need to resolve here should only happen after a heuristic solution.
4185  (Note default OSI implementation of optimalBasisIsAvailable always returns
4186  false.)
4187*/
4188    if (solverCharacteristics_->warmStart()&&
4189        !solver_->optimalBasisIsAvailable()) {
4190      //printf("XXXXYY no opt basis\n");
4191      resolve(node ? node->nodeInfo() : NULL,3);
4192    }
4193    if (nextRowCut_) {
4194      // branch was a cut - add it
4195      theseCuts.insert(*nextRowCut_);
4196      if (handler_->logLevel()>1)
4197        nextRowCut_->print();
4198      const OsiRowCut * cut=nextRowCut_;
4199      double lb = cut->lb();
4200      double ub = cut->ub();
4201      int n=cut->row().getNumElements();
4202      const int * column = cut->row().getIndices();
4203      const double * element = cut->row().getElements();
4204      double sum=0.0;
4205      for (int i=0;i<n;i++) {
4206        int iColumn = column[i];
4207        double value = element[i];
4208        //if (cbcColSolution_[iColumn]>1.0e-7)
4209        //printf("value of %d is %g\n",iColumn,cbcColSolution_[iColumn]);
4210        sum += value * cbcColSolution_[iColumn];
4211      }
4212      delete nextRowCut_;
4213      nextRowCut_=NULL;
4214      if (handler_->logLevel()>1)
4215        printf("applying branch cut, sum is %g, bounds %g %g\n",sum,lb,ub);
4216      // possibly extend whichGenerator
4217      resizeWhichGenerator(numberViolated, numberViolated+1);
4218      // set whichgenerator (also serves as marker to say don't delete0
4219      whichGenerator_[numberViolated++]=-2;
4220    }
4221    double * newSolution = new double [numberColumns] ;
4222    double heuristicValue = getCutoff() ;
4223    int found = -1; // no solution found
4224
4225    // reset probing info
4226    if (probingInfo_)
4227      probingInfo_->initializeFixing();
4228    for (int i = 0;i<numberCutGenerators_+numberHeuristics_;i++) {
4229      int numberRowCutsBefore = theseCuts.sizeRowCuts() ;
4230      int numberColumnCutsBefore = theseCuts.sizeColCuts() ;
4231      if (i<numberCutGenerators_) {
4232        bool generate = generator_[i]->normal();
4233        // skip if not optimal and should be (maybe a cut generator has fixed variables)
4234        if (generator_[i]->needsOptimalBasis()&&!solver_->basisIsAvailable())
4235          generate=false;
4236        if (generator_[i]->switchedOff())
4237          generate=false;;
4238        if (generate) {
4239          bool mustResolve = 
4240            generator_[i]->generateCuts(theseCuts,fullScan,node) ;
4241          int numberRowCutsAfter = theseCuts.sizeRowCuts() ;
4242          if(numberRowCutsBefore < numberRowCutsAfter &&
4243             generator_[i]->mustCallAgain())
4244            keepGoing=true; // say must go round
4245          // Check last cut to see if infeasible
4246          if(numberRowCutsBefore < numberRowCutsAfter) {
4247            const OsiRowCut * thisCut = theseCuts.rowCutPtr(numberRowCutsAfter-1) ;
4248            if (thisCut->lb()>thisCut->ub()) {
4249              feasible = false; // sub-problem is infeasible
4250              break;
4251            }
4252          }
4253#ifdef CBC_DEBUG
4254          {
4255            int numberRowCutsAfter = theseCuts.sizeRowCuts() ;
4256            int k ;
4257            for (k = numberRowCutsBefore;k<numberRowCutsAfter;k++) {
4258              OsiRowCut thisCut = theseCuts.rowCut(k) ;
4259              /* check size of elements.
4260                 We can allow smaller but this helps debug generators as it
4261                 is unsafe to have small elements */
4262              int n=thisCut.row().getNumElements();
4263              const int * column = thisCut.row().getIndices();
4264              const double * element = thisCut.row().getElements();
4265              //assert (n);
4266              for (int i=0;i<n;i++) {
4267                double value = element[i];
4268                assert(fabs(value)>1.0e-12&&fabs(value)<1.0e20);
4269              }
4270            }
4271          }
4272#endif
4273          if (mustResolve) {
4274            int returncode = resolve(node ? node->nodeInfo() : NULL,2);
4275            feasible = returnCode  != 0 ;
4276            if (returncode<0)
4277              numberTries=0;
4278            if ((specialOptions_&1)!=0) {
4279              debugger = solver_->getRowCutDebugger() ;
4280              if (debugger) 
4281                onOptimalPath = (debugger->onOptimalPath(*solver_)) ;
4282              else
4283                onOptimalPath=false;
4284              if (onOptimalPath && !solver_->isDualObjectiveLimitReached())
4285                assert(feasible) ;
4286            }
4287            if (!feasible)
4288              break ;
4289          }
4290        }
4291      } else {
4292        // see if heuristic will do anything
4293        double saveValue = heuristicValue ;
4294        int ifSol = 
4295          heuristic_[i-numberCutGenerators_]->solution(heuristicValue,
4296                                                       newSolution,
4297                                                       theseCuts) ;
4298        if (ifSol>0) {
4299          // better solution found
4300          found = i-numberCutGenerators_ ;
4301          incrementUsed(newSolution);
4302        } else if (ifSol<0) {
4303          heuristicValue = saveValue ;
4304        }
4305      }
4306      int numberRowCutsAfter = theseCuts.sizeRowCuts() ;
4307      int numberColumnCutsAfter = theseCuts.sizeColCuts() ;
4308
4309      if ((specialOptions_&1)!=0) {
4310        if (onOptimalPath) {
4311          int k ;
4312          for (k = numberRowCutsBefore;k<numberRowCutsAfter;k++) {
4313            OsiRowCut thisCut = theseCuts.rowCut(k) ;
4314            if(debugger->invalidCut(thisCut)) {
4315              solver_->writeMps("badCut");
4316            }
4317            assert(!debugger->invalidCut(thisCut)) ;
4318          }
4319        }
4320      }
4321/*
4322  The cut generator/heuristic has done its thing, and maybe it generated some
4323  cuts and/or a new solution.  Do a bit of bookkeeping: load
4324  whichGenerator[i] with the index of the generator responsible for a cut,
4325  and place cuts flagged as global in the global cut pool for the model.
4326
4327  lastNumberCuts is the sum of cuts added in previous iterations; it's the
4328  offset to the proper starting position in whichGenerator.
4329*/
4330      int numberBefore =
4331            numberRowCutsBefore+numberColumnCutsBefore+lastNumberCuts ;
4332      int numberAfter =
4333            numberRowCutsAfter+numberColumnCutsAfter+lastNumberCuts ;
4334      // possibly extend whichGenerator
4335      resizeWhichGenerator(numberBefore, numberAfter);
4336      int j ;
4337      if (fullScan) {
4338        // counts
4339        countColumnCuts[i] += numberColumnCutsAfter-numberColumnCutsBefore ;
4340      }
4341      countRowCuts[i] += numberRowCutsAfter-numberRowCutsBefore ;
4342     
4343      for (j = numberRowCutsBefore;j<numberRowCutsAfter;j++) {
4344        whichGenerator_[numberBefore++] = i ;
4345        const OsiRowCut * thisCut = theseCuts.rowCutPtr(j) ;
4346        if (thisCut->lb()>thisCut->ub())
4347          violated=-2; // sub-problem is infeasible
4348        if (thisCut->globallyValid()) {
4349          // add to global list
4350          OsiRowCut newCut(*thisCut);
4351          newCut.setGloballyValid(2);
4352          newCut.mutableRow().setTestForDuplicateIndex(false);
4353          globalCuts_.insert(newCut) ;
4354        }
4355      }
4356      for (j = numberColumnCutsBefore;j<numberColumnCutsAfter;j++) {
4357        whichGenerator_[numberBefore++] = i ;
4358        const OsiColCut * thisCut = theseCuts.colCutPtr(j) ;
4359        if (thisCut->globallyValid()) {
4360          // add to global list
4361          OsiColCut newCut(*thisCut);
4362          newCut.setGloballyValid(2);
4363          globalCuts_.insert(newCut) ;
4364        }
4365      }
4366    }
4367    // If at root node and first pass do heuristics without cuts
4368    if (!numberNodes_&&currentPassNumber_==1) {
4369      // Save number solutions
4370      int saveNumberSolutions = numberSolutions_;
4371      int saveNumberHeuristicSolutions = numberHeuristicSolutions_;
4372      // make sure bounds tight
4373      { 
4374        for (int i = 0;i < numberObjects_;i++)
4375          object_[i]->resetBounds(solver_) ;
4376      }
4377      for (int i = 0;i<numberHeuristics_;i++) {
4378        // see if heuristic will do anything
4379        double saveValue = heuristicValue ;
4380        int ifSol = heuristic_[i]->solution(heuristicValue,
4381                                            newSolution);
4382        if (ifSol>0) {
4383          // better solution found
4384          found = i ;
4385          incrementUsed(newSolution);
4386          // increment number of solutions so other heuristics can test
4387          numberSolutions_++;
4388          numberHeuristicSolutions_++;
4389        } else {
4390          heuristicValue = saveValue ;
4391        }
4392      }
4393      // Restore number solutions
4394      numberSolutions_ = saveNumberSolutions;
4395      numberHeuristicSolutions_ = saveNumberHeuristicSolutions;
4396    }
4397    // Add in any violated saved cuts
4398    if (!theseCuts.sizeRowCuts()&&!theseCuts.sizeColCuts()) {
4399      int numberOld = theseCuts.sizeRowCuts();
4400      int numberCuts = slackCuts.sizeRowCuts() ;
4401      int i;
4402      // possibly extend whichGenerator
4403      resizeWhichGenerator(numberOld, numberOld+numberCuts);
4404      for ( i = 0;i<numberCuts;i++) {
4405        const OsiRowCut * thisCut = slackCuts.rowCutPtr(i) ;
4406        if (thisCut->violated(cbcColSolution_)>100.0*primalTolerance) {
4407          if (messageHandler()->logLevel()>2)
4408            printf("Old cut added - violation %g\n",
4409                   thisCut->violated(cbcColSolution_)) ;
4410          whichGenerator_[numberOld++]=-1;
4411          theseCuts.insert(*thisCut) ;
4412        }
4413      }
4414    }
4415/*
4416  End of the loop to exercise each generator/heuristic.
4417
4418  Did any of the heuristics turn up a new solution? Record it before we free
4419  the vector.
4420*/
4421    if (found >= 0) { 
4422      phase_=4;
4423      incrementUsed(newSolution);
4424      setBestSolution(CBC_ROUNDING,heuristicValue,newSolution) ;
4425      lastHeuristic_ = heuristic_[found];
4426      CbcTreeLocal * tree
4427          = dynamic_cast<CbcTreeLocal *> (tree_);
4428      if (tree)
4429        tree->passInSolution(bestSolution_,heuristicValue);
4430    }
4431    delete [] newSolution ;
4432
4433#if 0
4434    // switch on to get all cuts printed
4435    theseCuts.printCuts() ;
4436#endif
4437    int numberColumnCuts = theseCuts.sizeColCuts() ;
4438    int numberRowCuts = theseCuts.sizeRowCuts() ;
4439    if (violated>=0)
4440      violated = numberRowCuts + numberColumnCuts ;
4441/*
4442  Apply column cuts (aka bound tightening). This may be partially redundant
4443  for column cuts returned by CglProbing, as generateCuts installs bounds
4444  from CglProbing when it determines it can fix a variable.
4445
4446  TODO: Looks like the use of violated has evolved. The value set above is
4447        completely ignored. All that's left is violated == -1 indicates some
4448        cut is violated, violated == -2 indicates infeasibility. Only
4449        infeasibility warrants exceptional action.
4450
4451  TODO: Strikes me that this code will fail to detect infeasibility, because
4452        the breaks escape the inner loops but the outer loop keeps going.
4453        Infeasibility in an early cut will be overwritten if a later cut is
4454        merely violated.
4455*/
4456    if (numberColumnCuts) {
4457
4458#ifdef CBC_DEBUG
4459      double * oldLower = new double [numberColumns] ;
4460      double * oldUpper = new double [numberColumns] ;
4461      memcpy(oldLower,cbcColLower_,numberColumns*sizeof(double)) ;
4462      memcpy(oldUpper,cbcColUpper_,numberColumns*sizeof(double)) ;
4463#endif
4464
4465      double integerTolerance = getDblParam(CbcIntegerTolerance) ;
4466      for (int i = 0;i<numberColumnCuts;i++) {
4467        const OsiColCut * thisCut = theseCuts.colCutPtr(i) ;
4468        const CoinPackedVector & lbs = thisCut->lbs() ;
4469        const CoinPackedVector & ubs = thisCut->ubs() ;
4470        int j ;
4471        int n ;
4472        const int * which ;
4473        const double * values ;
4474        n = lbs.getNumElements() ;
4475        which = lbs.getIndices() ;
4476        values = lbs.getElements() ;
4477        for (j = 0;j<n;j++) {
4478          int iColumn = which[j] ;
4479          double value = cbcColSolution_[iColumn] ;
4480#if CBC_DEBUG>1
4481          printf("%d %g %g %g %g\n",iColumn,oldLower[iColumn],
4482                 cbcColSolution_[iColumn],oldUpper[iColumn],values[j]) ;
4483#endif
4484          solver_->setColLower(iColumn,values[j]) ;
4485          if (value<values[j]-integerTolerance)
4486            violated = -1 ;
4487          if (values[j]>cbcColUpper_[iColumn]+integerTolerance) {
4488            // infeasible
4489            violated = -2 ;
4490            break ;
4491          }
4492        }
4493        n = ubs.getNumElements() ;
4494        which = ubs.getIndices() ;
4495        values = ubs.getElements() ;
4496        for (j = 0;j<n;j++) {
4497          int iColumn = which[j] ;
4498          double value = cbcColSolution_[iColumn] ;
4499#if CBC_DEBUG>1
4500          printf("%d %g %g %g %g\n",iColumn,oldLower[iColumn],
4501                 cbcColSolution_[iColumn],oldUpper[iColumn],values[j]) ;
4502#endif
4503          solver_->setColUpper(iColumn,values[j]) ;
4504          if (value>values[j]+integerTolerance)
4505            violated = -1 ;
4506          if (values[j]<cbcColLower_[iColumn]-integerTolerance) {
4507            // infeasible
4508            violated = -2 ;
4509            break ;
4510          }
4511        }
4512      }
4513#ifdef CBC_DEBUG
4514      delete [] oldLower ;
4515      delete [] oldUpper ;
4516#endif
4517    }
4518/*
4519  End installation of column cuts. The break here escapes the numberTries
4520  loop.
4521*/
4522    if (violated == -2||!feasible) {
4523      // infeasible
4524      feasible = false ;
4525      violated = -2;
4526      if (!numberNodes_) 
4527        messageHandler()->message(CBC_INFEAS,
4528                                  messages())
4529                                    << CoinMessageEol ;
4530      break ;
4531    }
4532/*
4533  Now apply the row (constraint) cuts. This is a bit more work because we need
4534  to obtain and augment the current basis.
4535
4536  TODO: Why do this work, if there are no row cuts? The current basis will do
4537        just fine.
4538*/
4539    int numberRowsNow = solver_->getNumRows() ;
4540    assert(numberRowsNow == numberRowsAtStart+lastNumberCuts) ;
4541    int numberToAdd = theseCuts.sizeRowCuts() ;
4542    numberNewCuts_ = lastNumberCuts+numberToAdd ;
4543/*
4544  Now actually add the row cuts and reoptimise.
4545
4546  Install the cuts in the solver using applyRowCuts and
4547  augment the basis with the corresponding slack. We also add each row cut to
4548  the set of row cuts (cuts.insert()) supplied as a parameter. The new basis
4549  must be set with setWarmStart().
4550
4551  TODO: Seems to me the original code could allocate addCuts with size 0, if
4552        numberRowCuts was 0 and numberColumnCuts was nonzero. That might
4553        explain the memory fault noted in the comment by AJK.  Unfortunately,
4554        just commenting out the delete[] results in massive memory leaks. Try
4555        a revision to separate the row cut case. Why do we need addCuts at
4556        all? A typing issue, apparently: OsiCut vs. OsiRowCut.
4557
4558  TODO: It looks to me as if numberToAdd and numberRowCuts are identical at
4559        this point. Confirm & get rid of one of them.
4560
4561  TODO: Any reason why the three loops can't be consolidated?
4562*/
4563    if (numberRowCuts > 0 || numberColumnCuts > 0)
4564    { if (numberToAdd > 0)
4565      { int i ;
4566        // Faster to add all at once
4567        const OsiRowCut ** addCuts = new const OsiRowCut * [numberToAdd] ;
4568        for (i = 0 ; i < numberToAdd ; i++)
4569        { addCuts[i] = &theseCuts.rowCut(i) ; }
4570        solver_->applyRowCuts(numberToAdd,addCuts) ;
4571        // AJK this caused a memory fault on Win32
4572        // may be okay now with ** form
4573        delete [] addCuts ;
4574        for (i = 0 ; i < numberToAdd ; i++)
4575        { cuts.insert(theseCuts.rowCut(i)) ; }
4576        CoinWarmStartBasis * basis = dynamic_cast<CoinWarmStartBasis*>(solver_->getWarmStart()) ;
4577        assert(basis != NULL); // make sure not volume
4578/* dylp bug
4579
4580  Consistent size used by OsiDylp as sanity check. Implicit resize seen
4581  as an error. Hence this call to resize is necessary.
4582*/
4583        basis->resize(numberRowsAtStart+numberNewCuts_,numberColumns) ;
4584        for (i = 0 ; i < numberToAdd ; i++) 
4585        { basis->setArtifStatus(numberRowsNow+i,
4586                                 CoinWarmStartBasis::basic) ; }
4587        if (solver_->setWarmStart(basis) == false)
4588        { throw CoinError("Fail setWarmStart() after cut installation.",
4589                          "solveWithCuts","CbcModel") ; }
4590        delete basis;
4591      }
4592      feasible = resolve(node ? node->nodeInfo() : NULL,2) ;
4593      if ( CoinCpuTime()-dblParam_[CbcStartSeconds] > dblParam_[CbcMaximumSeconds] )
4594        numberTries=0; // exit
4595#     ifdef CBC_DEBUG
4596      printf("Obj value after cuts %g %d rows\n",solver_->getObjValue(),
4597              solver_->getNumRows()) ;
4598      if (onOptimalPath && !solver_->isDualObjectiveLimitReached())
4599        assert(feasible) ;
4600#     endif
4601    }
4602/*
4603  No cuts. Cut short the cut generation (numberTries) loop.
4604*/
4605    else
4606    { numberTries = 0 ; }
4607/*
4608  If the problem is still feasible, first, call takeOffCuts() to remove cuts
4609  that are now slack. takeOffCuts() will call the solver to reoptimise if
4610  that's needed to restore a valid solution.
4611
4612  Next, see if we should quit due to diminishing returns:
4613    * we've tried three rounds of cut generation and we're getting
4614      insufficient improvement in the objective; or
4615    * we generated no cuts; or
4616    * the solver declared optimality with 0 iterations after we added the
4617      cuts generated in this round.
4618  If we decide to keep going, prep for the next iteration.
4619
4620  It just seems more safe to tell takeOffCuts() to call resolve(), even if
4621  we're not continuing cut generation. Otherwise code executed between here
4622  and final disposition of the node will need to be careful not to access the
4623  lp solution. It can happen that we lose feasibility in takeOffCuts ---
4624  numerical jitters when the cutoff bound is epsilon less than the current
4625  best, and we're evaluating an alternative optimum.
4626
4627  TODO: After successive rounds of code motion, there seems no need to
4628        distinguish between the two checks for aborting the cut generation
4629        loop. Confirm and clean up.
4630*/
4631    if (feasible)
4632    { int cutIterations = solver_->getIterationCount() ;
4633      if (numberOldActiveCuts_+numberNewCuts_) {
4634        OsiCuts * saveCuts = node ? NULL : &slackCuts;
4635        takeOffCuts(cuts,resolveAfterTakeOffCuts_,saveCuts) ;
4636        if (solver_->isDualObjectiveLimitReached()&&resolveAfterTakeOffCuts_)
4637          { feasible = false ;
4638#       ifdef CBC_DEBUG
4639          double z = solver_->getObjValue() ;
4640          double cut = getCutoff() ;
4641          printf("Lost feasibility by %g in takeOffCuts; z = %g, cutoff = %g\n",
4642                 z-cut,z,cut) ;
4643#       endif
4644          }
4645      }
4646      if (feasible) {
4647        numberRowsAtStart = numberOldActiveCuts_+numberRowsAtContinuous_ ;
4648        lastNumberCuts = numberNewCuts_ ;
4649        if (direction*solver_->getObjValue() < lastObjective+minimumDrop &&
4650            currentPassNumber_ >= 3 && !keepGoing)
4651          { numberTries = 0 ; }
4652        if (numberRowCuts+numberColumnCuts == 0 || cutIterations == 0)
4653          { break ; }
4654        if (numberTries > 0) {
4655          reducedCostFix() ;
4656          lastObjective = direction*solver_->getObjValue() ;
4657        }
4658      }
4659    }
4660/*
4661  We've lost feasibility --- this node won't be referencing the cuts we've
4662  been collecting, so decrement the reference counts.
4663*/
4664    if (!feasible)
4665    { int i ;
4666      for (i = 0;i<currentNumberCuts_;i++) {
4667        // take off node
4668        if (addedCuts_[i]) {
4669          if (!addedCuts_[i]->decrement())
4670            delete addedCuts_[i] ;
4671          addedCuts_[i] = NULL ;
4672        }
4673      }
4674      numberTries = 0 ;
4675    }
4676  } while (numberTries>0||keepGoing) ;
4677  {
4678    // switch on
4679    for (int i = 0;i<numberCutGenerators_;i++) 
4680      generator_[i]->setSwitchedOff(false);
4681  }
4682  //check feasibility.
4683  //If solution seems to be integer feasible calling setBestSolution
4684  //will eventually add extra global cuts which we need to install at
4685  //the nodes
4686
4687
4688  if(feasible&&solverCharacteristics_->solutionAddsCuts())  //check integer feasibility
4689  {
4690    bool integerFeasible = true;
4691    const double * save = testSolution_;
4692    testSolution_ = solver_->getColSolution();
4693    // point to useful information
4694    OsiBranchingInformation usefulInfo=usefulInformation();
4695    for (int i=0;i<numberObjects_ && integerFeasible;i++)
4696    {
4697      int preferredWay;
4698      double infeasibility = object_[i]->infeasibility(&usefulInfo,preferredWay);
4699      if(infeasibility)
4700        integerFeasible = false;
4701    }
4702    testSolution_ = save;
4703    if(integerFeasible)//update
4704    {
4705      double objValue = solver_->getObjValue();
4706      int numberGlobalBefore = globalCuts_.sizeRowCuts();
4707      // SOLUTION2 so won't up cutoff or print message
4708      setBestSolution(CBC_SOLUTION2, objValue, 
4709                      solver_->getColSolution(),0);
4710      int numberGlobalAfter = globalCuts_.sizeRowCuts();
4711      int numberToAdd = numberGlobalAfter - numberGlobalBefore;
4712      if(numberToAdd > 0)
4713        //We have added some cuts say they are tight at that node
4714        //Basis and lp should already have been updated
4715      {
4716        feasible = (solver_->isProvenOptimal() &&
4717                    !solver_->isDualObjectiveLimitReached()) ;
4718        if(feasible)
4719        {
4720          int numberCuts = numberNewCuts_ = cuts.sizeRowCuts();
4721      // possibly extend whichGenerator
4722          resizeWhichGenerator(numberCuts, numberToAdd+numberCuts);
4723         
4724          for (int i = numberGlobalBefore ; i < numberGlobalAfter ; i++)
4725            {       
4726              whichGenerator_[numberNewCuts_++]=-1;
4727#ifndef GLOBAL_CUTS_JUST_POINTERS
4728              cuts.insert(globalCuts_.rowCut(i)) ; 
4729#else
4730              OsiRowCut * rowCutPointer = globalCuts_.rowCutPtr(i);
4731              cuts.insert(rowCutPointer) ; 
4732#endif
4733            }
4734          numberNewCuts_ = lastNumberCuts+numberToAdd;
4735          //now take off the cuts which are not tight anymore
4736          takeOffCuts(cuts,resolveAfterTakeOffCuts_, NULL) ;
4737          if (solver_->isDualObjectiveLimitReached()&&resolveAfterTakeOffCuts_)
4738            { feasible = false ;}
4739        }
4740        if(!feasible) //node will be fathomed
4741          { 
4742          for (int i = 0;i<currentNumberCuts_;i++) 
4743            {
4744            // take off node
4745            if (addedCuts_[i])
4746              {
4747                if (!addedCuts_[i]->decrement())
4748                  delete addedCuts_[i] ;
4749                addedCuts_[i] = NULL ;
4750              }
4751            }
4752        }
4753      }
4754    }
4755  }
4756/*
4757  Reduced cost fix at end. Must also check feasible, in case we've popped out
4758  because a generator indicated we're infeasible.
4759*/
4760  if (feasible && solver_->isProvenOptimal())
4761    reducedCostFix();
4762  // If at root node do heuristics
4763  if (!numberNodes_) {
4764    // First see if any cuts are slack
4765    int numberRows = solver_->getNumRows();
4766    int numberAdded = numberRows-numberRowsAtContinuous_;
4767    if (numberAdded) {
4768      CoinWarmStartBasis * basis = dynamic_cast<CoinWarmStartBasis*>(solver_->getWarmStart()) ;
4769      assert(basis != NULL);
4770      int * added = new int[numberAdded];
4771      int nDelete = 0;
4772      for (int j=numberRowsAtContinuous_;j<numberRows;j++) {
4773        if (basis->getArtifStatus(j)==CoinWarmStartBasis::basic) {
4774          //printf("%d slack!\n",j);
4775          added[nDelete++]=j;
4776        }
4777      }
4778      if (nDelete) {
4779        solver_->deleteRows(nDelete,added);
4780      }
4781      delete [] added;
4782      delete basis ;
4783    }
4784    // mark so heuristics can tell
4785    int savePass=currentPassNumber_;
4786    currentPassNumber_=999999;
4787    double * newSolution = new double [numberColumns] ;
4788    double heuristicValue = getCutoff() ;
4789    int found = -1; // no solution found
4790    for (int i = 0;i<numberHeuristics_;i++) {
4791      // see if heuristic will do anything
4792      double saveValue = heuristicValue ;
4793      int ifSol = heuristic_[i]->solution(heuristicValue,
4794                                          newSolution);
4795      if (ifSol>0) {
4796        // better solution found
4797        found = i ;
4798        incrementUsed(newSolution);
4799      } else {
4800        heuristicValue = saveValue ;
4801      }
4802    }
4803    currentPassNumber_=savePass;
4804    if (found >= 0) { 
4805      phase_=4;
4806      incrementUsed(newSolution);
4807      setBestSolution(CBC_ROUNDING,heuristicValue,newSolution) ;
4808      lastHeuristic_ = heuristic_[found];
4809    }
4810    delete [] newSolution ;
4811  }
4812  // Up change due to cuts
4813  if (feasible)
4814    sumChangeObjective2_ += solver_->getObjValue()*solver_->getObjSense()
4815      - objectiveValue ;
4816  //if ((numberNodes_%100)==0)
4817  //printf("XXb sum obj changed by %g\n",sumChangeObjective2_);
4818/*
4819  End of cut generation loop.
4820
4821  Now, consider if we want to disable or adjust the frequency of use for any
4822  of the cut generators. If the client specified a positive number for
4823  howOften, it will never change. If the original value was negative, it'll
4824  be converted to 1000000+|howOften|, and this value will be adjusted each
4825  time fullScan is true. Actual cut generation is performed every
4826  howOften%1000000 nodes; the 1000000 offset is just a convenient way to
4827  specify that the frequency is adjustable.
4828
4829  During cut generation, we recorded the number of cuts produced by each
4830  generator for this node. For all cuts, whichGenerator records the generator
4831  that produced a cut.
4832
4833  TODO: All this should probably be hidden in a method of the CbcCutGenerator
4834  class.
4835
4836  TODO: Can the loop that scans over whichGenerator to accumulate per generator
4837        counts be replaced by values in countRowCuts and countColumnCuts?
4838*/
4839#ifdef NODE_LOG
4840  int fatherNum = (node == NULL) ? -1 : node->nodeNumber();
4841  double value =  (node == NULL) ? -1 : node->branchingObject()->value();
4842  string bigOne = (solver_->getIterationCount() > 30) ? "*******" : "";
4843  string way = (node == NULL) ? "" : (node->branchingObject()->way())==1 ? "Down" : "Up";
4844  std::cout<<"Node "<<numberNodes_<<", father "<<fatherNum<<", #iterations "<<solver_->getIterationCount()<<", sol value : "<<solver_->getObjValue()<<std::endl;
4845#endif
4846  if (fullScan&&numberCutGenerators_) {
4847    /* If cuts just at root node then it will probably be faster to
4848       update matrix and leave all in */
4849    int willBeCutsInTree=0;
4850    double thisObjective = solver_->getObjValue()*direction ;
4851    // get sizes
4852    int numberRowsAdded = solver_->getNumRows()-numberRowsAtStart;
4853    CoinBigIndex numberElementsAdded =  solver_->getNumElements()-numberElementsAtStart ;
4854    double densityOld = ((double) numberElementsAtStart)/((double) numberRowsAtStart);
4855    double densityNew = numberRowsAdded ? ((double) (numberElementsAdded))/((double) numberRowsAdded)
4856      : 0.0;
4857    if (!numberNodes_) {
4858      if (numberRowsAdded)
4859        handler_->message(CBC_CUTS_STATS,messages_)
4860          <<numberRowsAdded
4861          <<densityNew
4862          <<CoinMessageEol ;
4863      if (thisObjective-startObjective<1.0e-5&&numberElementsAdded>0.2*numberElementsAtStart)
4864        willBeCutsInTree=-1;
4865    }
4866    if ((numberRowsAdded>100+0.5*numberRowsAtStart
4867         ||numberElementsAdded>0.5*numberElementsAtStart)
4868        &&(densityNew>200.0&&numberRowsAdded>100&&densityNew>2.0*densityOld)) {
4869      // much bigger
4870      //if (thisObjective-startObjective<0.1*fabs(startObjective)+1.0e-5)
4871      //willBeCutsInTree=-1;
4872      //printf("Cuts will be taken off , %d rows added with density %g\n",
4873      //     numberRowsAdded,densityNew);
4874    }
4875    if (densityNew>100.0&&numberRowsAdded>2&&densityNew>2.0*densityOld) {
4876      //if (thisObjective-startObjective<0.1*fabs(startObjective)+1.0e-5)
4877      //willBeCutsInTree=-2;
4878      //printf("Density says no cuts ? , %d rows added with density %g\n",
4879      //     numberRowsAdded,densityNew);
4880    }
4881    // Root node or every so often - see what to turn off
4882    int i ;
4883    for (i = 0;i<numberCutGenerators_;i++) {
4884      int howOften = generator_[i]->howOften() ;
4885      if (howOften>-90) 
4886        willBeCutsInTree=0;
4887    }
4888    if (!numberNodes_)
4889      handler_->message(CBC_ROOT,messages_)
4890        <<numberNewCuts_
4891        <<startObjective<<thisObjective
4892        <<currentPassNumber_
4893        <<CoinMessageEol ;
4894    int * count = new int[numberCutGenerators_] ;
4895    memset(count,0,numberCutGenerators_*sizeof(int)) ;
4896    int numberActiveGenerators=0;
4897    for (i = 0;i<numberNewCuts_;i++) {
4898      int iGenerator = whichGenerator_[i];
4899      if (iGenerator>=0&&iGenerator<numberCutGenerators_)
4900        count[iGenerator]++ ;
4901    }
4902    double totalCuts = 0.0 ;
4903    //#define JUST_ACTIVE
4904    for (i = 0;i<numberCutGenerators_;i++) {
4905      if (countRowCuts[i]||countColumnCuts[i])
4906        numberActiveGenerators++;
4907#ifdef JUST_ACTIVE
4908      totalCuts += count[i] + 5.0*countColumnCuts[i] ;
4909#else
4910      totalCuts += countRowCuts[i] + 5.0*countColumnCuts[i] ;
4911#endif
4912    }
4913    int iProbing=-1;
4914    double small = (0.5* totalCuts) /
4915      ((double) numberActiveGenerators) ;
4916    for (i = 0;i<numberCutGenerators_;i++) {
4917      int howOften = generator_[i]->howOften() ;
4918      /*  Probing can be set to just do column cuts in treee.
4919          But if doing good then leave as on */
4920      bool probingWasOnBut=false;
4921      CglProbing * probing = dynamic_cast<CglProbing*>(generator_[i]->generator());
4922      if (probing&&!numberNodes_) {
4923        iProbing = i;
4924        if (probing->rowCuts()==-3) {
4925          probingWasOnBut=true;
4926          howOften = -98;
4927          probing->setRowCuts(3);
4928        }
4929      }
4930      if (willBeCutsInTree<0&&howOften==-98)
4931        howOften =-99;
4932      if (howOften==-98&&generator_[i]->switchOffIfLessThan()<0) {
4933        if (thisObjective-startObjective<0.005*fabs(startObjective)+1.0e-5)
4934          howOften=-99; // switch off
4935        if (thisObjective-startObjective<0.1*fabs(startObjective)+1.0e-5
4936            &&5*solver_->getNumRows()<solver_->getNumCols())
4937          howOften=-99; // switch off
4938      }
4939      if (howOften<-99)
4940        continue ;
4941      if (howOften<0||howOften >= 1000000) {
4942        if( !numberNodes_) {
4943          // If small number switch mostly off
4944#ifdef JUST_ACTIVE
4945          double thisCuts = count[i] + 5.0*countColumnCuts[i] ;
4946#else
4947          double thisCuts = countRowCuts[i] + 5.0*countColumnCuts[i] ;
4948#endif
4949          if (!thisCuts||howOften == -99) {
4950            if (howOften == -99||howOften == -98) 
4951              howOften = -100 ;
4952            else
4953              howOften = 1000000+SCANCUTS; // wait until next time
4954          } else if (thisCuts<small) {
4955            if (howOften!=1&&!probingWasOnBut) {
4956              if (generator_[i]->whatDepth()<0||howOften!=-1) {
4957                int k = (int) sqrt(small/thisCuts) ;
4958                if (howOften!=-98)
4959                  howOften = k+1000000 ;
4960                else
4961                  howOften=-100;
4962              } else {
4963                howOften=1;
4964              }
4965            } else {
4966              howOften=1;
4967              // allow cuts
4968              probingWasOnBut=false;
4969            }
4970          } else {
4971            if (thisObjective-startObjective<0.1*fabs(startObjective)+1.0e-5&&generator_[i]->whatDepth()<0)
4972              generator_[i]->setWhatDepth(5);
4973            howOften = 1+1000000 ;
4974          }
4975        }
4976        // If cuts useless switch off
4977        if (numberNodes_>=100000&&sumChangeObjective1_>2.0e2*(sumChangeObjective2_+1.0e-12)) {
4978          howOften = 1000000+SCANCUTS; // wait until next time
4979          //printf("switch off cut %d due to lack of use\n",i);
4980        }
4981      }
4982      if (!numberNodes_) {
4983        if (probingWasOnBut&&howOften==-100) {
4984          probing->setRowCuts(-3);
4985          howOften=1;
4986        }
4987        if (howOften==1)
4988          generator_[i]->setWhatDepth(1);
4989       
4990        if (howOften>=0&&generator_[i]->generator()->mayGenerateRowCutsInTree())
4991          willBeCutsInTree=1;
4992      }
4993      generator_[i]->setHowOften(howOften) ;
4994      if (howOften>=1000000&&howOften<2000000&&0) {
4995        // Go to depth
4996        int bias=1;
4997        if (howOften==1+1000000)
4998          generator_[i]->setWhatDepth(bias+1);
4999        else if (howOften<=10+1000000)
5000          generator_[i]->setWhatDepth(bias+2);
5001        else
5002          generator_[i]->setWhatDepth(bias+1000);
5003      }
5004      int newFrequency = generator_[i]->howOften()%1000000 ;
5005      // increment cut counts
5006      generator_[i]->incrementNumberCutsInTotal(countRowCuts[i]);
5007      generator_[i]->incrementNumberCutsActive(count[i]);
5008      if (handler_->logLevel()>1||!numberNodes_) {
5009        handler_->message(CBC_GENERATOR,messages_)
5010          <<i
5011          <<generator_[i]->cutGeneratorName()
5012          <<countRowCuts[i]<<count[i]
5013          <<countColumnCuts[i];
5014        handler_->printing(!numberNodes_&&generator_[i]->timing())
5015          <<generator_[i]->timeInCutGenerator();
5016        handler_->message()
5017          <<newFrequency
5018          <<CoinMessageEol ;
5019      }
5020    } 
5021    delete [] count ;
5022    if( !numberNodes_) {
5023      // save statistics
5024      for (i = 0;i<numberCutGenerators_;i++) {
5025        generator_[i]->setNumberCutsAtRoot(generator_[i]->numberCutsInTotal());
5026        generator_[i]->setNumberActiveCutsAtRoot(generator_[i]->numberCutsActive());
5027      }
5028      // decide on pseudo cost strategy
5029      int howOften = iProbing>=0 ? generator_[iProbing]->howOften() : 0;
5030      if ((howOften %1000000)!=1) 
5031        howOften = 0;
5032      //if (howOften) {
5033      //CglProbing * probing = dynamic_cast<CglProbing*>(generator_[iProbing]->generator());
5034      //}
5035      howOften = 0;
5036      if (howOften) {
5037        printf("** method 1\n");
5038        //CglProbing * probing = dynamic_cast<CglProbing*>(generator_[iProbing]->generator());
5039        generator_[iProbing]->setWhatDepth(1);
5040        // could set no row cuts
5041        //if (thisObjective-startObjective<0.001*fabs(startObjective)+1.0e-5)
5042        // probing->setRowCuts(0);
5043        for (int i=0;i<numberObjects_;i++) {
5044          CbcSimpleIntegerDynamicPseudoCost * obj =
5045            dynamic_cast <CbcSimpleIntegerDynamicPseudoCost *>(object_[i]) ;
5046          if (obj) 
5047            obj->setMethod(1);
5048        }
5049      }
5050      if (willBeCutsInTree==-2)
5051        willBeCutsInTree=0;
5052      if( willBeCutsInTree<=0) {
5053        // Take off cuts
5054        cuts = OsiCuts();
5055        numberNewCuts_=0;
5056        if (!willBeCutsInTree) {
5057          // update size of problem
5058          numberRowsAtContinuous_ = solver_->getNumRows() ;
5059        } else {
5060          // take off cuts
5061          int numberRows = solver_->getNumRows();
5062          int numberAdded = numberRows-numberRowsAtContinuous_;
5063          if (numberAdded) {
5064            int * added = new int[numberAdded];
5065            for (int i=0;i<numberAdded;i++)
5066              added[i]=i+numberRowsAtContinuous_;
5067            solver_->deleteRows(numberAdded,added);
5068            delete [] added;
5069            // resolve so optimal
5070            resolve(solver_);
5071          }
5072        }
5073#ifdef COIN_HAS_CLP
5074        OsiClpSolverInterface * clpSolver
5075          = dynamic_cast<OsiClpSolverInterface *> (solver_);
5076        if (clpSolver) {
5077          // Maybe solver might like to know only column bounds will change
5078          //int options = clpSolver->specialOptions();
5079          //clpSolver->setSpecialOptions(options|128);
5080          clpSolver->synchronizeModel();
5081        }
5082#endif
5083      } else {
5084#ifdef COIN_HAS_CLP
5085        OsiClpSolverInterface * clpSolver
5086          = dynamic_cast<OsiClpSolverInterface *> (solver_);
5087        if (clpSolver) {
5088        // make sure factorization can't carry over
5089          int options = clpSolver->specialOptions();
5090          clpSolver->setSpecialOptions(options&(~8));
5091        }
5092#endif
5093      }
5094    }
5095  } else {
5096#ifdef COIN_HAS_CLP
5097    OsiClpSolverInterface * clpSolver
5098      = dynamic_cast<OsiClpSolverInterface *> (solver_);
5099    if (clpSolver) {
5100      // Maybe solver might like to know only column bounds will change
5101      //int options = clpSolver->specialOptions();
5102      //clpSolver->setSpecialOptions(options|128);
5103      clpSolver->synchronizeModel();
5104    }
5105#endif
5106    if (numberCutGenerators_) {
5107      int i;
5108      // add to counts anyway
5109      for (i = 0;i<numberCutGenerators_;i++) 
5110        generator_[i]->incrementNumberCutsInTotal(countRowCuts[i]);
5111      // What if not feasible as cuts may have helped
5112      if (feasible) {
5113        for (i = 0;i<numberNewCuts_;i++) {
5114          int iGenerator = whichGenerator_[i];
5115          if (iGenerator>=0)
5116            generator_[iGenerator]->incrementNumberCutsActive();
5117        }
5118      }
5119    }
5120  }
5121
5122  delete [] countRowCuts ;
5123  delete [] countColumnCuts ;
5124
5125
5126#ifdef CHECK_CUT_COUNTS
5127  if (feasible)
5128  { 
5129    CoinWarmStartBasis * basis = dynamic_cast<CoinWarmStartBasis*>(solver_->getWarmStart()) ;
5130    printf("solveWithCuts: Number of rows at end (only active cuts) %d\n",
5131           numberRowsAtContinuous_+numberNewCuts_+numberOldActiveCuts_) ;
5132    basis->print() ; delete basis;}
5133#endif
5134#ifdef CBC_DEBUG
5135  if (onOptimalPath && !solver_->isDualObjectiveLimitReached())
5136    assert(feasible) ;
5137#endif
5138# ifdef COIN_HAS_CLP
5139  if (clpSolver) 
5140    clpSolver->setSpecialOptions(saveClpOptions);
5141# endif
5142
5143  return feasible ; }
5144
5145
5146/*
5147  Remove slack cuts. We obtain a basis and scan it. Cuts with basic slacks
5148  are purged. If any cuts are purged, resolve() is called to restore the
5149  solution held in the solver.  If resolve() pivots, there's the possibility
5150  that a slack may be pivoted in (trust me :-), so the process iterates.
5151  Setting allowResolve to false will suppress reoptimisation (but see note
5152  below).
5153
5154  At the level of the solver's constraint system, loose cuts are really
5155  deleted.  There's an implicit assumption that deleteRows will also update
5156  the active basis in the solver.
5157
5158  At the level of nodes and models, it's more complicated.
5159
5160  New cuts exist only in the collection of cuts passed as a parameter. They
5161  are deleted from the collection and that's the end of them.
5162
5163  Older cuts have made it into addedCuts_. Two separate actions are needed.
5164  The reference count for the CbcCountRowCut object is decremented. If this
5165  count falls to 0, the node which owns the cut is located, the reference to
5166  the cut is removed, and then the cut object is destroyed (courtesy of the
5167  CbcCountRowCut destructor). We also need to set the addedCuts_ entry to
5168  NULL. This is important so that when it comes time to generate basis edits
5169  we can tell this cut was dropped from the basis during processing of the
5170  node.
5171
5172  NOTE: In general, it's necessary to call resolve() after purging slack
5173        cuts.  Deleting constraints constitutes a change in the problem, and
5174        an OSI is not required to maintain a valid solution when the problem
5175        is changed. But ... it's really useful to maintain the active basis,
5176        and the OSI is supposed to do that. (Yes, it's splitting hairs.) In
5177        some places, it's possible to know that the solution will never be
5178        consulted after this call, only the basis.  (E.g., this routine is
5179        called as a last act before generating info to place the node in the
5180        live set.) For such use, set allowResolve to false.
5181 
5182  TODO: No real harm would be done if we just ignored the rare occasion when
5183        the call to resolve() pivoted a slack back into the basis. It's a
5184        minor inefficiency, at worst. But it does break assertions which
5185        check that there are no loose cuts in the basis. It might be better
5186        to remove the assertions.
5187*/
5188
5189void
5190CbcModel::takeOffCuts (OsiCuts &newCuts,
5191                       bool allowResolve, OsiCuts * saveCuts)
5192
5193{ // int resolveIterations = 0 ;
5194  int firstOldCut = numberRowsAtContinuous_ ;
5195  int totalNumberCuts = numberNewCuts_+numberOldActiveCuts_ ;
5196  int *solverCutIndices = new int[totalNumberCuts] ;
5197  int *newCutIndices = new int[numberNewCuts_] ;
5198  const CoinWarmStartBasis* ws ;
5199  CoinWarmStartBasis::Status status ;
5200  bool needPurge = true ;
5201/*
5202  The outer loop allows repetition of purge in the event that reoptimisation
5203  changes the basis. To start an iteration, clear the deletion counts and grab
5204  the current basis.
5205*/
5206  while (needPurge)
5207  { int numberNewToDelete = 0 ;
5208    int numberOldToDelete = 0 ;
5209    int i ;
5210    ws = dynamic_cast<const CoinWarmStartBasis*>(solver_->getWarmStart()) ;
5211/*
5212  Scan the basis entries of the old cuts generated prior to this round of cut
5213  generation.  Loose cuts are `removed' by decrementing their reference count
5214  and setting the addedCuts_ entry to NULL. (If the reference count falls to
5215  0, they're really deleted.  See CbcModel and CbcCountRowCut doc'n for
5216  principles of cut handling.)
5217*/
5218    int oldCutIndex = 0 ;
5219    for (i = 0 ; i < numberOldActiveCuts_ ; i++)
5220    { status = ws->getArtifStatus(i+firstOldCut) ;
5221      while (!addedCuts_[oldCutIndex]) oldCutIndex++ ;
5222      assert(oldCutIndex < currentNumberCuts_) ;
5223      // always leave if from nextRowCut_
5224      if (status == CoinWarmStartBasis::basic&&
5225          addedCuts_[oldCutIndex]->effectiveness()!=COIN_DBL_MAX)
5226      { solverCutIndices[numberOldToDelete++] = i+firstOldCut ;
5227        if (saveCuts) {
5228          // send to cut pool
5229          OsiRowCut * slackCut = addedCuts_[oldCutIndex];
5230          if (slackCut->effectiveness()!=-1.234) {
5231            slackCut->setEffectiveness(-1.234);
5232            saveCuts->insert(*slackCut);
5233          }
5234        }
5235        if (addedCuts_[oldCutIndex]->decrement() == 0)
5236          delete addedCuts_[oldCutIndex] ;
5237        addedCuts_[oldCutIndex] = NULL ;
5238        oldCutIndex++ ; }
5239      else
5240      { oldCutIndex++ ; } }
5241/*
5242  Scan the basis entries of the new cuts generated with this round of cut
5243  generation.  At this point, newCuts is the only record of the new cuts, so
5244  when we delete loose cuts from newCuts, they're really gone. newCuts is a
5245  vector, so it's most efficient to compress it (eraseRowCut) from back to
5246  front.
5247*/
5248    int firstNewCut = firstOldCut+numberOldActiveCuts_ ;
5249    int k = 0 ;
5250    for (i = 0 ; i < numberNewCuts_ ; i++)
5251    { status = ws->getArtifStatus(i+firstNewCut) ;
5252      if (status == CoinWarmStartBasis::basic&&whichGenerator_[i]!=-2)
5253      { solverCutIndices[numberNewToDelete+numberOldToDelete] = i+firstNewCut ;
5254        newCutIndices[numberNewToDelete++] = i ; }
5255      else
5256      { // save which generator did it
5257        whichGenerator_[k++] = whichGenerator_[i] ; } }
5258    delete ws ;
5259    for (i = numberNewToDelete-1 ; i >= 0 ; i--)
5260    { int iCut = newCutIndices[i] ;
5261      if (saveCuts) {
5262        // send to cut pool
5263        OsiRowCut * slackCut = newCuts.rowCutPtr(iCut);
5264        if (slackCut->effectiveness()!=-1.234) {
5265          slackCut->setEffectiveness(-1.234);
5266          saveCuts->insert(*slackCut);
5267        }
5268      }
5269      newCuts.eraseRowCut(iCut) ; }
5270/*
5271  Did we delete anything? If so, delete the cuts from the constraint system
5272  held in the solver and reoptimise unless we're forbidden to do so. If the
5273  call to resolve() results in pivots, there's the possibility we again have
5274  basic slacks. Repeat the purging loop.
5275*/
5276    if (numberNewToDelete+numberOldToDelete > 0)
5277    { solver_->deleteRows(numberNewToDelete+numberOldToDelete,
5278                          solverCutIndices) ;
5279      numberNewCuts_ -= numberNewToDelete ;
5280      numberOldActiveCuts_ -= numberOldToDelete ;
5281#     ifdef CBC_DEBUG
5282      printf("takeOffCuts: purged %d+%d cuts\n", numberOldToDelete,
5283             numberNewToDelete );
5284#     endif
5285      if (allowResolve)
5286      { 
5287        phase_=3;
5288        // can do quick optimality check
5289        int easy=2;
5290        solver_->setHintParam(OsiDoInBranchAndCut,true,OsiHintDo,&easy) ;
5291        resolve(solver_) ;
5292        setPointers(solver_);
5293        solver_->setHintParam(OsiDoInBranchAndCut,true,OsiHintDo,NULL) ;
5294        if (solver_->getIterationCount() == 0)
5295        { needPurge = false ; }
5296#       ifdef CBC_DEBUG
5297        else
5298          { printf( "Repeating purging loop. %d iters.\n",
5299                    solver_->getIterationCount()); }
5300#       endif
5301      }
5302      else
5303      { needPurge = false ; } }
5304    else
5305    { needPurge = false ; } }
5306/*
5307  Clean up and return.
5308*/
5309  delete [] solverCutIndices ;
5310  delete [] newCutIndices ;
5311}
5312
5313
5314
5315int
5316CbcModel::resolve(CbcNodeInfo * parent, int whereFrom)
5317{
5318  // We may have deliberately added in violated cuts - check to avoid message
5319  int iRow;
5320  int numberRows = solver_->getNumRows();
5321  const double * rowLower = solver_->getRowLower();
5322  const double * rowUpper = solver_->getRowUpper();
5323  bool feasible=true;
5324  for (iRow= numberRowsAtContinuous_;iRow<numberRows;iRow++) {
5325    if (rowLower[iRow]>rowUpper[iRow]+1.0e-8)
5326      feasible=false;
5327  }
5328  // Can't happen if strong branching as would have been found before
5329  if (!numberStrong_&&numberObjects_>numberIntegers_) {
5330    int iColumn;
5331    int numberColumns = solver_->getNumCols();
5332    const double * columnLower = solver_->getColLower();
5333    const double * columnUpper = solver_->getColUpper();
5334    for (iColumn= 0;iColumn<numberColumns;iColumn++) {
5335      if (columnLower[iColumn]>columnUpper[iColumn]+1.0e-5)
5336        feasible=false;
5337    }
5338  }
5339#if 0
5340  {
5341    int iColumn;
5342    int numberColumns = solver_->getNumCols();
5343    const double * columnLower = solver_->getColLower();
5344    const double * columnUpper = solver_->getColUpper();
5345    const double * sol = solver_->getColSolution();
5346    bool feasible=true;
5347    double xx[37];
5348    memset(xx,0,sizeof(xx));
5349    xx[6]=xx[7]=xx[9]=xx[18]=xx[26]=xx[36]=1.0;
5350    for (iColumn= 0;iColumn<numberColumns;iColumn++) {
5351      if (solver_->isInteger(iColumn)) {
5352        //printf("yy %d at %g (bounds %g, %g)\n",iColumn,
5353        //     sol[iColumn],columnLower[iColumn],columnUpper[iColumn]);
5354        if (columnLower[iColumn]>xx[iColumn]||
5355            columnUpper[iColumn]<xx[iColumn])
5356          feasible=false;
5357        //solver_->setColLower(iColumn,CoinMin(xx[iColumn],columnUpper[iColumn]));
5358        //solver_->setColUpper(iColumn,CoinMax(xx[iColumn],columnLower[iColumn]));
5359      }
5360    }
5361    if (feasible)
5362      printf("On Feasible path\n");
5363  }
5364#endif
5365/*
5366  Reoptimize. Consider the possibility that we should fathom on bounds. But be
5367  careful --- where the objective takes on integral values, we may want to keep
5368  a solution where the objective is right on the cutoff.
5369*/
5370  if (feasible)
5371    {
5372      resolve(solver_) ;
5373      numberIterations_ += solver_->getIterationCount() ;
5374      feasible = (solver_->isProvenOptimal() &&
5375                  !solver_->isDualObjectiveLimitReached()) ;
5376    }
5377  if (0&&feasible) {
5378    const double * lb = solver_->getColLower();
5379    const double * ub = solver_->getColUpper();
5380    const double * x = solver_->getColSolution();
5381    const double * dj = solver_->getReducedCost();
5382    int numberColumns = solver_->getNumCols();
5383    for (int i=0;i<numberColumns;i++) {
5384      if (dj[i]>1.0e-4&&ub[i]-lb[i]>1.0e-4&&x[i]>lb[i]+1.0e-4)
5385        printf("error %d %g %g %g %g\n",i,dj[i],lb[i],x[i],ub[i]);
5386      if (dj[i]<-1.0e-4&&ub[i]-lb[i]>1.0e-4&&x[i]<ub[i]-1.0e-4)
5387        printf("error %d %g %g %g %g\n",i,dj[i],lb[i],x[i],ub[i]);
5388    }
5389  } 
5390  if (false&&!feasible&& continuousObjective_ <-1.0e30) {
5391    // at root node - double double check
5392    bool saveTakeHint;
5393    OsiHintStrength saveStrength;
5394    solver_->getHintParam(OsiDoDualInResolve,saveTakeHint,saveStrength);
5395    if (saveTakeHint||saveStrength==OsiHintIgnore) {
5396      solver_->setHintParam(OsiDoDualInResolve,false,OsiHintDo) ;
5397      resolve(solver_);
5398      solver_->setHintParam(OsiDoDualInResolve,saveTakeHint,saveStrength);
5399      numberIterations_ += solver_->getIterationCount() ;
5400      feasible = solver_->isProvenOptimal();
5401      //      solver_->writeMps("infeas");
5402    }
5403  }
5404  if (cutModifier_&&feasible&&!solverCharacteristics_->solutionAddsCuts()) {
5405    //double increment = getDblParam(CbcModel::CbcCutoffIncrement) ;
5406    double cutoff ;
5407    solver_->getDblParam(OsiDualObjectiveLimit,cutoff) ;
5408    double distance = fabs(cutoff-solver_->getObjValue());
5409    if (distance<10.0*trueIncrement) {
5410      double offset;
5411      solver_->getDblParam(OsiObjOffset,offset);
5412      double objFixedValue = -offset;
5413      double objValue=0.0;
5414      double direction = solver_->getObjSense();
5415      const double * solution = solver_->getColSolution();
5416      const double * objective = solver_->getObjCoefficients();
5417      const double * columnLower = solver_->getColLower();
5418      const double * columnUpper = solver_->getColUpper();
5419      int numberColumns = solver_->getNumCols();
5420      int increment = 0 ;
5421      double multiplier = 1.0/trueIncrement;
5422      int bigIntegers = 0; // Count of large costs which are integer
5423      for (int iColumn=0;iColumn<numberColumns;iColumn++) {
5424        double value = solution[iColumn];
5425        // make sure clean
5426        value = CoinMin(value,columnUpper[iColumn]);
5427        value = CoinMax(value,columnLower[iColumn]);
5428        double cost = direction * objective[iColumn];
5429        if (cost) {
5430          if (columnLower[iColumn]<columnUpper[iColumn]) {
5431            objValue += value*cost;
5432            value = fabs(cost)*multiplier ;
5433            if (value <2.1e9) {
5434              int nearest = (int) floor(value+0.5) ;
5435              assert (fabs(value-floor(value+0.5)) < 1.0e-8);
5436              if (!increment) 
5437                increment = nearest ; 
5438              else
5439                increment = gcd(increment,nearest) ;
5440            } else {
5441              // large value - may still be multiple of 1.0
5442              value = fabs(objective[iColumn]);
5443              assert(fabs(value-floor(value+0.5)) < 1.0e-8);
5444              bigIntegers++;
5445            }
5446          } else {
5447            // fixed
5448            objFixedValue += value*cost;
5449          }
5450        }
5451      }
5452      if (increment) {
5453        double value = increment ;
5454        value /= multiplier ;
5455        if (value>trueIncrement) {
5456          double x = objValue/value;
5457          x = ceil(x-1.0e-5);
5458          x *= value;
5459          //printf("fixed %g, variable %g -> %g, sum %g - cutoff %g\n",
5460          // objFixedValue,objValue,x,x+objFixedValue,cutoff);
5461          x += objFixedValue;
5462          if (x>cutoff + 1.0e-5*fabs(cutoff)+1.0e-5) {
5463            //printf("Node cutoff\n");
5464            feasible=false;
5465          }
5466        } else {
5467          value = trueIncrement;
5468          double x = objValue/value;
5469          x = ceil(x-1.0e-5);
5470          x *= value;
5471          x += objFixedValue;
5472          if (x>cutoff + 1.0e-5*fabs(cutoff)+1.0e-5) {
5473            //printf("Node cutoff\n");
5474            feasible=false;
5475          }
5476        }
5477      }
5478    }
5479  }
5480
5481  setPointers(solver_);
5482  int returnStatus = feasible ? 1 : 0;
5483  if (strategy_) {
5484    // user can play clever tricks here
5485    int status = strategy_->status(this,parent,whereFrom);
5486    if (status>=0) {
5487      if (status==0)
5488        returnStatus = 1;
5489      else if(status==1)
5490        returnStatus=-1;
5491      else
5492        returnStatus=0;
5493    }
5494  }
5495  return returnStatus ;
5496}
5497
5498
5499/* Set up objects.  Only do ones whose length is in range.
5500   If makeEquality true then a new model may be returned if
5501   modifications had to be made, otherwise "this" is returned.
5502
5503   Could use Probing at continuous to extend objects
5504*/
5505CbcModel * 
5506CbcModel::findCliques(bool makeEquality,
5507                      int atLeastThisMany, int lessThanThis, int defaultValue)
5508{
5509  // No objects are allowed to exist
5510  assert(numberObjects_==numberIntegers_||!numberObjects_);
5511  CoinPackedMatrix matrixByRow(*solver_->getMatrixByRow());
5512  int numberRows = solver_->getNumRows();
5513  int numberColumns = solver_->getNumCols();
5514
5515  // We may want to add columns
5516  int numberSlacks=0;
5517  int * rows = new int[numberRows];
5518  double * element =new double[numberRows];
5519
5520  int iRow;
5521
5522  findIntegers(true);
5523  numberObjects_=numberIntegers_;
5524
5525  int numberCliques=0;
5526  OsiObject ** object = new OsiObject * [numberRows];
5527  int * which = new int[numberIntegers_];
5528  char * type = new char[numberIntegers_];
5529  int * lookup = new int[numberColumns];
5530  int i;
5531  for (i=0;i<numberColumns;i++) 
5532    lookup[i]=-1;
5533  for (i=0;i<numberIntegers_;i++) 
5534    lookup[integerVariable_[i]]=i;
5535
5536  // Statistics
5537  int totalP1=0,totalM1=0;
5538  int numberBig=0,totalBig=0;
5539  int numberFixed=0;
5540
5541  // Row copy
5542  const double * elementByRow = matrixByRow.getElements();
5543  const int * column = matrixByRow.getIndices();
5544  const CoinBigIndex * rowStart = matrixByRow.getVectorStarts();
5545  const int * rowLength = matrixByRow.getVectorLengths();
5546
5547  // Column lengths for slacks
5548  const int * columnLength = solver_->getMatrixByCol()->getVectorLengths();
5549
5550  const double * lower = getColLower();
5551  const double * upper = getColUpper();
5552  const double * rowLower = getRowLower();
5553  const double * rowUpper = getRowUpper();
5554
5555  for (iRow=0;iRow<numberRows;iRow++) {
5556    int numberP1=0, numberM1=0;
5557    int j;
5558    double upperValue=rowUpper[iRow];
5559    double lowerValue=rowLower[iRow];
5560    bool good=true;
5561    int slack = -1;
5562    for (j=rowStart[iRow];j<rowStart[iRow]+rowLength[iRow];j++) {
5563      int iColumn = column[j];
5564      int iInteger=lookup[iColumn];
5565      if (upper[iColumn]-lower[iColumn]<1.0e-8) {
5566        // fixed
5567        upperValue -= lower[iColumn]*elementByRow[j];
5568        lowerValue -= lower[iColumn]*elementByRow[j];
5569        continue;
5570      } else if (upper[iColumn]!=1.0||lower[iColumn]!=0.0) {
5571        good = false;
5572        break;
5573      } else if (iInteger<0) {
5574        good = false;
5575        break;
5576      } else {
5577        if (columnLength[iColumn]==1)
5578          slack = iInteger;
5579      }
5580      if (fabs(elementByRow[j])!=1.0) {
5581        good=false;
5582        break;
5583      } else if (elementByRow[j]>0.0) {
5584        which[numberP1++]=iInteger;
5585      } else {
5586        numberM1++;
5587        which[numberIntegers_-numberM1]=iInteger;
5588      }
5589    }
5590    int iUpper = (int) floor(upperValue+1.0e-5);
5591    int iLower = (int) ceil(lowerValue-1.0e-5);
5592    int state=0;
5593    if (upperValue<1.0e6) {
5594      if (iUpper==1-numberM1)
5595        state=1;
5596      else if (iUpper==-numberM1)
5597        state=2;
5598      else if (iUpper<-numberM1)
5599        state=3;
5600    }
5601    if (!state&&lowerValue>-1.0e6) {
5602      if (-iLower==1-numberP1)
5603        state=-1;
5604      else if (-iLower==-numberP1)
5605        state=-2;
5606      else if (-iLower<-numberP1)
5607        state=-3;
5608    }
5609    if (good&&state) {
5610      if (abs(state)==3) {
5611        // infeasible
5612        numberObjects_=-1;
5613        break;
5614      } else if (abs(state)==2) {
5615        // we can fix all
5616        numberFixed += numberP1+numberM1;
5617        if (state>0) {
5618          // fix all +1 at 0, -1 at 1
5619          for (i=0;i<numberP1;i++)
5620            solver_->setColUpper(integerVariable_[which[i]],0.0);
5621          for (i=0;i<numberM1;i++)
5622            solver_->setColLower(integerVariable_[which[numberIntegers_-i-1]],
5623                                 1.0);
5624        } else {
5625          // fix all +1 at 1, -1 at 0
5626          for (i=0;i<numberP1;i++)
5627            solver_->setColLower(integerVariable_[which[i]],1.0);
5628          for (i=0;i<numberM1;i++)
5629            solver_->setColUpper(integerVariable_[which[numberIntegers_-i-1]],
5630                                 0.0);
5631        }
5632      } else {
5633        int length = numberP1+numberM1;
5634        if (length >= atLeastThisMany&&length<lessThanThis) {
5635          // create object
5636          bool addOne=false;
5637          int objectType;
5638          if (iLower==iUpper) {
5639            objectType=1;
5640          } else {
5641            if (makeEquality) {
5642              objectType=1;
5643              element[numberSlacks]=state;
5644              rows[numberSlacks++]=iRow;
5645              addOne=true;
5646            } else {
5647              objectType=0;
5648            }
5649          }
5650          if (state>0) {
5651            totalP1 += numberP1;
5652            totalM1 += numberM1;
5653            for (i=0;i<numberP1;i++)
5654              type[i]=1;
5655            for (i=0;i<numberM1;i++) {
5656              which[numberP1]=which[numberIntegers_-i-1];
5657              type[numberP1++]=0;
5658            }
5659          } else {
5660            totalP1 += numberM1;
5661            totalM1 += numberP1;
5662            for (i=0;i<numberP1;i++)
5663              type[i]=0;
5664            for (i=0;i<numberM1;i++) {
5665              which[numberP1]=which[numberIntegers_-i-1];
5666              type[numberP1++]=1;
5667            }
5668          }
5669          if (addOne) {
5670            // add in slack
5671            which[numberP1]=numberIntegers_+numberSlacks-1;
5672            slack = numberP1;
5673            type[numberP1++]=1;
5674          } else if (slack >= 0) {
5675            for (i=0;i<numberP1;i++) {
5676              if (which[i]==slack) {
5677                slack=i;
5678              }
5679            }
5680          }
5681          object[numberCliques] = new CbcClique(this,objectType,numberP1,
5682                                              which,type,
5683                                               1000000+numberCliques,slack);
5684          numberCliques++;
5685        } else if (numberP1+numberM1 >= lessThanThis) {
5686          // too big
5687          numberBig++;
5688          totalBig += numberP1+numberM1;
5689        }
5690      }
5691    }
5692  }
5693  delete [] which;
5694  delete [] type;
5695  delete [] lookup;
5696  if (numberCliques<0) {
5697    printf("*** Problem infeasible\n");
5698  } else {
5699    if (numberCliques)
5700      printf("%d cliques of average size %g found, %d P1, %d M1\n",
5701             numberCliques,
5702             ((double)(totalP1+totalM1))/((double) numberCliques),
5703             totalP1,totalM1);
5704    else
5705      printf("No cliques found\n");
5706    if (numberBig)
5707      printf("%d large cliques ( >= %d) found, total %d\n",
5708             numberBig,lessThanThis,totalBig);
5709    if (numberFixed)
5710      printf("%d variables fixed\n",numberFixed);
5711  }
5712  if (numberCliques>0&&numberSlacks&&makeEquality) {
5713    printf("adding %d integer slacks\n",numberSlacks);
5714    // add variables to make equality rows
5715    int * temp = new int[numberIntegers_+numberSlacks];
5716    memcpy(temp,integerVariable_,numberIntegers_*sizeof(int));
5717    // Get new model
5718    CbcModel * newModel = new CbcModel(*this);
5719    OsiSolverInterface * newSolver = newModel->solver();
5720    for (i=0;i<numberSlacks;i++) {
5721      temp[i+numberIntegers_]=i+numberColumns;
5722      int iRow = rows[i];
5723      double value = element[i];
5724      double lowerValue = 0.0;
5725      double upperValue = 1.0;
5726      double objValue  = 0.0;
5727      CoinPackedVector column(1,&iRow,&value);
5728      newSolver->addCol(column,lowerValue,upperValue,objValue);
5729      // set integer
5730      newSolver->setInteger(numberColumns+i);
5731      if (value >0)
5732        newSolver->setRowLower(iRow,rowUpper[iRow]);
5733      else
5734        newSolver->setRowUpper(iRow,rowLower[iRow]);
5735    }
5736    // replace list of integers
5737    for (i=0;i<newModel->numberObjects_;i++)
5738      delete newModel->object_[i];
5739    newModel->numberObjects_ = 0;
5740    delete [] newModel->object_;
5741    newModel->object_=NULL;
5742    newModel->findIntegers(true); //Set up all integer objects
5743    for (i=0;i<numberIntegers_;i++) {
5744      newModel->modifiableObject(i)->setPriority(object_[i]->priority());
5745    }
5746    if (originalColumns_) {
5747      // old model had originalColumns
5748      delete [] newModel->originalColumns_;
5749      newModel->originalColumns_ = new int[numberColumns+numberSlacks];
5750      memcpy(newModel->originalColumns_,originalColumns_,numberColumns*sizeof(int));
5751      // mark as not in previous model
5752      for (i=numberColumns;i<numberColumns+numberSlacks;i++)
5753        newModel->originalColumns_[i]=-1;
5754    }
5755    delete [] rows;
5756    delete [] element;
5757    newModel->addObjects(numberCliques,object);
5758    assert (ownObjects_);
5759    for (;i<numberCliques;i++) 
5760      delete object[i];
5761    delete [] object;
5762    newModel->synchronizeModel();
5763    return newModel;
5764  } else {
5765    assert (ownObjects_);
5766    if (numberCliques>0) {
5767      addObjects(numberCliques,object);
5768      for (;i<numberCliques;i++) 
5769        delete object[i];
5770      synchronizeModel();
5771    }
5772    delete [] object;
5773    delete [] rows;
5774    delete [] element;
5775    return this;
5776  }
5777}
5778// Fill in useful estimates
5779void 
5780CbcModel::pseudoShadow(double * down, double * up)
5781{
5782  // Column copy of matrix
5783  const double * element = solver_->getMatrixByCol()->getElements();
5784  const int * row = solver_->getMatrixByCol()->getIndices();
5785  const CoinBigIndex * columnStart = solver_->getMatrixByCol()->getVectorStarts();
5786  const int * columnLength = solver_->getMatrixByCol()->getVectorLengths();
5787  const double *objective = solver_->getObjCoefficients() ;
5788  int numberColumns = solver_->getNumCols() ;
5789  double direction = solver_->getObjSense();
5790  int iColumn;
5791  const double * dual = cbcRowPrice_;
5792  down = new double[numberColumns];
5793  up = new double[numberColumns];
5794  double upSum=1.0e-20;
5795  double downSum = 1.0e-20;
5796  int numberIntegers=0;
5797  for (iColumn=0;iColumn<numberColumns;iColumn++) {
5798    CoinBigIndex start = columnStart[iColumn];
5799    CoinBigIndex end = start + columnLength[iColumn];
5800    double upValue = 0.0;
5801    double downValue = 0.0;
5802    double value = direction*objective[iColumn];
5803    if (value) {
5804      if (value>0.0)
5805        upValue += value;
5806      else
5807        downValue -= value;
5808    }
5809    for (CoinBigIndex j=start;j<end;j++) {
5810      int iRow = row[j];
5811      value = -dual[iRow];
5812      if (value) {
5813        value *= element[j];
5814        if (value>0.0)
5815          upValue += value;
5816        else
5817          downValue -= value;
5818      }
5819    }
5820    // use dj if bigger
5821    double dj = cbcReducedCost_[iColumn];
5822    upValue = CoinMax(upValue,dj);
5823    downValue = CoinMax(downValue,-dj);
5824    up[iColumn]=upValue;
5825    down[iColumn]=downValue;
5826    if (solver_->isInteger(iColumn)) {
5827      if (!numberNodes_)
5828        printf("%d - dj %g up %g down %g cost %g\n",
5829               iColumn,dj,upValue,downValue,objective[iColumn]);
5830      upSum += upValue;
5831      downSum += downValue;
5832      numberIntegers++;
5833    }
5834  }
5835  if (numberIntegers) {
5836    double smallDown = 0.01*(downSum/((double) numberIntegers));
5837    double smallUp = 0.01*(upSum/((double) numberIntegers));
5838    for (int i=0;i<numberObjects_;i++) {
5839      CbcSimpleIntegerDynamicPseudoCost * obj1 =
5840        dynamic_cast <CbcSimpleIntegerDynamicPseudoCost *>(object_[i]) ;
5841      if (obj1) {
5842        iColumn = obj1->columnNumber();
5843        double upPseudoCost = obj1->upDynamicPseudoCost();
5844        double saveUp = upPseudoCost;
5845        upPseudoCost = CoinMax(upPseudoCost,smallUp);
5846        upPseudoCost = CoinMax(upPseudoCost,up[iColumn]);
5847        upPseudoCost = CoinMax(upPseudoCost,0.1*down[iColumn]);
5848        obj1->setUpDynamicPseudoCost(upPseudoCost);
5849        if (upPseudoCost>saveUp&&!numberNodes_)
5850          printf("For %d up went from %g to %g\n",
5851                 iColumn,saveUp,upPseudoCost);
5852        double downPseudoCost = obj1->downDynamicPseudoCost();
5853        double saveDown = downPseudoCost;
5854        downPseudoCost = CoinMax(downPseudoCost,smallDown);
5855        downPseudoCost = CoinMax(downPseudoCost,down[iColumn]);
5856        downPseudoCost = CoinMax(downPseudoCost,0.1*down[iColumn]);
5857        obj1->setDownDynamicPseudoCost(downPseudoCost);
5858        if (downPseudoCost>saveDown&&!numberNodes_)
5859          printf("For %d down went from %g to %g\n",
5860                 iColumn,saveDown,downPseudoCost);
5861      }
5862    }
5863  }
5864  delete [] down;
5865  delete [] up;
5866}
5867
5868/*
5869  Set branching priorities.
5870
5871  Setting integer priorities looks pretty robust; the call to findIntegers
5872  makes sure that SimpleInteger objects are in place. Setting priorities for
5873  other objects is entirely dependent on their existence, and the routine may
5874  quietly fail in several directions.
5875*/
5876
5877void 
5878CbcModel::passInPriorities (const int * priorities,
5879                            bool ifObject)
5880{
5881  findIntegers(false);
5882  int i;
5883  if (priorities) {
5884    int i0=0;
5885    int i1=numberObjects_-1;
5886    if (ifObject) {
5887      for (i=numberIntegers_;i<numberObjects_;i++) {
5888        object_[i]->setPriority(priorities[i-numberIntegers_]);
5889      }
5890      i0=numberIntegers_;
5891    } else {
5892      for (i=0;i<numberIntegers_;i++) {
5893        object_[i]->setPriority(priorities[i]);
5894      }
5895      i1=numberIntegers_-1;
5896    }
5897    messageHandler()->message(CBC_PRIORITY,
5898                              messages())
5899                                << i0<<i1<<numberObjects_ << CoinMessageEol ;
5900  }
5901}
5902
5903// Delete all object information
5904void 
5905CbcModel::deleteObjects(bool getIntegers)
5906{
5907  if (ownObjects_) {
5908    int i;
5909    for