source: stable/2.0/Cbc/src/CbcModel.cpp @ 844

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

make useSolution work better

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