source: trunk/Cbc/src/CbcModel.cpp @ 687

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

add threaded to trunk

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