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

Last change on this file since 747 was 747, checked in by forrest, 12 years ago

max times and fix bug in fpump

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