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

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

changes to try and improve performance

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