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

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

possible undefined reference on windows?

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