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

Last change on this file since 940 was 940, checked in by forrest, 14 years ago

for my experiments

  • Property svn:eol-style set to native
  • Property svn:keywords set to Author Date Id Revision
File size: 449.0 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/* Most of copy constructor
4789      mode - 0 copy but don't delete before
4790             1 copy and delete before
4791             2 copy and delete before (but use virgin generators)
4792*/
4793void 
4794CbcModel::gutsOfCopy(const CbcModel & rhs,int mode)
4795{
4796  minimumDrop_ = rhs.minimumDrop_;
4797  specialOptions_ = rhs.specialOptions_;
4798  numberStrong_ = rhs.numberStrong_;
4799  numberBeforeTrust_ = rhs.numberBeforeTrust_;
4800  numberPenalties_ = rhs.numberPenalties_;
4801  printFrequency_ = rhs.printFrequency_;
4802# ifdef COIN_HAS_CLP
4803  fastNodeDepth_ = rhs.fastNodeDepth_;
4804#endif
4805  howOftenGlobalScan_ = rhs.howOftenGlobalScan_;
4806  maximumCutPassesAtRoot_ = rhs.maximumCutPassesAtRoot_;
4807  maximumCutPasses_ =  rhs.maximumCutPasses_;
4808  preferredWay_ = rhs.preferredWay_;
4809  resolveAfterTakeOffCuts_ = rhs.resolveAfterTakeOffCuts_;
4810  numberThreads_ = rhs.numberThreads_;
4811  threadMode_ = rhs.threadMode_;
4812  memcpy(intParam_,rhs.intParam_,sizeof(intParam_));
4813  memcpy(dblParam_,rhs.dblParam_,sizeof(dblParam_));
4814  int i;
4815  if (mode) {
4816    for (i=0;i<numberCutGenerators_;i++) {
4817      delete generator_[i];
4818      delete virginGenerator_[i];
4819    }
4820    delete [] generator_;
4821    delete [] virginGenerator_;
4822    for (i=0;i<numberHeuristics_;i++) {
4823      delete heuristic_[i];
4824    }
4825    delete [] heuristic_;
4826    delete eventHandler_;
4827    delete branchingMethod_;
4828  }
4829  numberCutGenerators_ = rhs.numberCutGenerators_;
4830  if (numberCutGenerators_) {
4831    generator_ = new CbcCutGenerator * [numberCutGenerators_];
4832    virginGenerator_ = new CbcCutGenerator * [numberCutGenerators_];
4833    int i;
4834    for (i=0;i<numberCutGenerators_;i++) {
4835      if (mode<2)
4836        generator_[i]=new CbcCutGenerator(*rhs.generator_[i]);
4837      else
4838        generator_[i]=new CbcCutGenerator(*rhs.virginGenerator_[i]);
4839      virginGenerator_[i]=new CbcCutGenerator(*rhs.virginGenerator_[i]);
4840    }
4841  } else {
4842    generator_=NULL;
4843    virginGenerator_=NULL;
4844  }
4845  numberHeuristics_ = rhs.numberHeuristics_;
4846  if (numberHeuristics_) {
4847    heuristic_ = new CbcHeuristic * [numberHeuristics_];
4848    int i;
4849    for (i=0;i<numberHeuristics_;i++) {
4850      heuristic_[i]=rhs.heuristic_[i]->clone();
4851    }
4852  } else {
4853    heuristic_=NULL;
4854  }
4855  if (rhs.eventHandler_)
4856    eventHandler_ = rhs.eventHandler_->clone() ; 
4857  else
4858    eventHandler_ = NULL ; 
4859  if (rhs.branchingMethod_)
4860    branchingMethod_=rhs.branchingMethod_->clone();
4861  else
4862    branchingMethod_=NULL;
4863  messageHandler()->setLogLevel(rhs.messageHandler()->logLevel());
4864  synchronizeModel();
4865}
4866// Move status, nodes etc etc across
4867void 
4868CbcModel::moveInfo(const CbcModel & rhs)
4869{
4870  bestObjective_ = rhs.bestObjective_;
4871  bestPossibleObjective_=rhs.bestPossibleObjective_;
4872  numberSolutions_=rhs.numberSolutions_;
4873  numberHeuristicSolutions_=rhs.numberHeuristicSolutions_;
4874  numberNodes_ = rhs.numberNodes_;
4875  numberNodes2_ = rhs.numberNodes2_;
4876  numberIterations_ = rhs.numberIterations_;
4877  status_ = rhs.status_;
4878  secondaryStatus_ = rhs.secondaryStatus_;
4879  numberStoppedSubTrees_ = rhs.numberStoppedSubTrees_;
4880  numberInfeasibleNodes_ = rhs.numberInfeasibleNodes_;
4881  continuousObjective_=rhs.continuousObjective_;
4882  originalContinuousObjective_ = rhs.originalContinuousObjective_;
4883  continuousInfeasibilities_ = rhs.continuousInfeasibilities_;
4884  numberFixedAtRoot_ = rhs.numberFixedAtRoot_;
4885  numberFixedNow_ = rhs.numberFixedNow_;
4886  stoppedOnGap_ = rhs.stoppedOnGap_;
4887  eventHappened_ = rhs.eventHappened_;
4888  numberLongStrong_ = rhs.numberLongStrong_;
4889  numberStrongIterations_ = rhs.numberStrongIterations_;
4890  strongInfo_[0]=rhs.strongInfo_[0];
4891  strongInfo_[1]=rhs.strongInfo_[1];
4892  strongInfo_[2]=rhs.strongInfo_[2];
4893  numberRowsAtContinuous_ = rhs.numberRowsAtContinuous_;
4894  maximumDepth_= rhs.maximumDepth_;
4895}
4896// Save a copy of the current solver so can be reset to
4897void 
4898CbcModel::saveReferenceSolver()
4899{
4900  delete referenceSolver_;
4901  referenceSolver_= solver_->clone();
4902}
4903
4904// Uses a copy of reference solver to be current solver
4905void 
4906CbcModel::resetToReferenceSolver()
4907{
4908  delete solver_;
4909  solver_ = referenceSolver_->clone();
4910  // clear many things
4911  gutsOfDestructor2();
4912  // Reset cutoff
4913  // Solvers know about direction
4914  double direction = solver_->getObjSense();
4915  double value;
4916  solver_->getDblParam(OsiDualObjectiveLimit,value); 
4917  setCutoff(value*direction);
4918}
4919
4920// Are there a numerical difficulties?
4921bool 
4922CbcModel::isAbandoned() const
4923{
4924  return status_ == 2;
4925}
4926// Is optimality proven?
4927bool 
4928CbcModel::isProvenOptimal() const
4929{
4930  if (!status_ && bestObjective_<1.0e30)
4931    return true;
4932  else
4933    return false;
4934}
4935// Is  infeasiblity proven (or none better than cutoff)?
4936bool 
4937CbcModel::isProvenInfeasible() const
4938{
4939  if (!status_ && bestObjective_>=1.0e30)
4940    return true;
4941  else
4942    return false;
4943}
4944// Was continuous solution unbounded
4945bool 
4946CbcModel::isContinuousUnbounded() const
4947{
4948  if (!status_ && secondaryStatus_==7)
4949    return true;
4950  else
4951    return false;
4952}
4953// Was continuous solution unbounded
4954bool 
4955CbcModel::isProvenDualInfeasible() const
4956{
4957  if (!status_ && secondaryStatus_==7)
4958    return true;
4959  else
4960    return false;
4961}
4962// Node limit reached?
4963bool 
4964CbcModel::isNodeLimitReached() const
4965{
4966  return numberNodes_ >= intParam_[CbcMaxNumNode];
4967}
4968// Time limit reached?
4969bool 
4970CbcModel::isSecondsLimitReached() const
4971{
4972  if (status_==1&&secondaryStatus_==4)
4973    return true;
4974  else
4975    return false;
4976}
4977// Solution limit reached?
4978bool 
4979CbcModel::isSolutionLimitReached() const
4980{
4981  return numberSolutions_ >= intParam_[CbcMaxNumSol];
4982}
4983// Set language
4984void 
4985CbcModel::newLanguage(CoinMessages::Language language)
4986{
4987  messages_ = CbcMessage(language);
4988}
4989void 
4990CbcModel::setNumberStrong(int number)
4991{
4992  if (number<0)
4993    numberStrong_=0;
4994   else
4995    numberStrong_=number;
4996}
4997void 
4998CbcModel::setNumberBeforeTrust(int number)
4999{
5000  if (number<-3) {
5001    numberBeforeTrust_=0;
5002  } else {
5003    numberBeforeTrust_=number;
5004    //numberStrong_ = CoinMax(numberStrong_,1);
5005  }
5006}
5007void 
5008CbcModel::setNumberPenalties(int number)
5009{
5010  if (number<=0) {
5011    numberPenalties_=0;
5012  } else {
5013    numberPenalties_=number;
5014  }
5015}
5016void 
5017CbcModel::setPenaltyScaleFactor(double value)
5018{
5019  if (value<=0) {
5020    penaltyScaleFactor_=3.0;
5021  } else {
5022    penaltyScaleFactor_=value;
5023  }
5024}
5025void 
5026CbcModel::setHowOftenGlobalScan(int number)
5027{
5028  if (number<-1)
5029    howOftenGlobalScan_=0;
5030   else
5031    howOftenGlobalScan_=number;
5032}
5033
5034// Add one generator
5035void 
5036CbcModel::addCutGenerator(CglCutGenerator * generator,
5037                          int howOften, const char * name,
5038                          bool normal, bool atSolution,
5039                          bool whenInfeasible,int howOftenInSub,
5040                          int whatDepth, int whatDepthInSub)
5041{
5042  CbcCutGenerator ** temp = generator_;
5043  generator_ = new CbcCutGenerator * [numberCutGenerators_+1];
5044  memcpy(generator_,temp,numberCutGenerators_*sizeof(CbcCutGenerator *));
5045  delete[] temp ;
5046  generator_[numberCutGenerators_]= 
5047    new CbcCutGenerator(this,generator, howOften, name,
5048                        normal,atSolution,whenInfeasible,howOftenInSub,
5049                        whatDepth, whatDepthInSub);
5050  // and before any changes
5051  temp = virginGenerator_;
5052  virginGenerator_ = new CbcCutGenerator * [numberCutGenerators_+1];
5053  memcpy(virginGenerator_,temp,numberCutGenerators_*sizeof(CbcCutGenerator *));
5054  delete[] temp ;
5055  virginGenerator_[numberCutGenerators_++]= 
5056    new CbcCutGenerator(this,generator, howOften, name,
5057                        normal,atSolution,whenInfeasible,howOftenInSub,
5058                        whatDepth, whatDepthInSub);
5059                                                         
5060}
5061// Add one heuristic
5062void 
5063CbcModel::addHeuristic(CbcHeuristic * generator, const char *name,
5064                       int before)
5065{
5066  CbcHeuristic ** temp = heuristic_;
5067  heuristic_ = new CbcHeuristic * [numberHeuristics_+1];
5068  memcpy(heuristic_,temp,numberHeuristics_*sizeof(CbcHeuristic *));
5069  delete [] temp;
5070  int where;
5071  if (before<0||before>=numberHeuristics_) {
5072    where=numberHeuristics_;
5073  } else {
5074    // move up
5075    for (int i=numberHeuristics_;i>before;i--) 
5076      heuristic_[i]=heuristic_[i-1];
5077    where=before;
5078  }
5079  heuristic_[where]=generator->clone();
5080  if (name)
5081    heuristic_[where]->setHeuristicName(name) ; 
5082  heuristic_[where]->setSeed(987654321+where);
5083  numberHeuristics_++ ;
5084}
5085
5086/*
5087  The last subproblem handled by the solver is not necessarily related to the
5088  one being recreated, so the first action is to remove all cuts from the
5089  constraint system.  Next, traverse the tree from node to the root to
5090  determine the basis size required for this subproblem and create an empty
5091  basis with the right capacity.  Finally, traverse the tree from root to
5092  node, adjusting bounds in the constraint system, adjusting the basis, and
5093  collecting the cuts that must be added to the constraint system.
5094  applyToModel does the heavy lifting.
5095
5096  addCuts1 is used in contexts where all that's desired is the list of cuts:
5097  the node is already fathomed, and we're collecting cuts so that we can
5098  adjust reference counts as we prune nodes. Arguably the two functions
5099  should be separated. The culprit is applyToModel, which performs cut
5100  collection and model adjustment.
5101
5102  Certainly in the contexts where all we need is a list of cuts, there's no
5103  point in passing in a valid basis --- an empty basis will do just fine.
5104*/
5105void CbcModel::addCuts1 (CbcNode * node, CoinWarmStartBasis *&lastws)
5106{ int i;
5107  int nNode=0;
5108  int numberColumns = getNumCols();
5109  CbcNodeInfo * nodeInfo = node->nodeInfo();
5110
5111/*
5112  Remove all cuts from the constraint system.
5113  (original comment includes ``see note below for later efficiency'', but
5114  the reference isn't clear to me).
5115*/
5116  int currentNumberCuts ;
5117  if (numberThreads_<=0) {
5118    solver_->restoreBaseModel(numberRowsAtContinuous_);
5119  } else {
5120    // *** Fix later
5121    currentNumberCuts = solver_->getNumRows()-numberRowsAtContinuous_;
5122    int *which = new int[currentNumberCuts];
5123    for (i = 0 ; i < currentNumberCuts ; i++)
5124      which[i] = i+numberRowsAtContinuous_;
5125    solver_->deleteRows(currentNumberCuts,which);
5126    delete [] which;
5127  }
5128/*
5129  Accumulate the path from node to the root in walkback_, and accumulate a
5130  cut count in currentNumberCuts.
5131
5132  original comment: when working then just unwind until where new node joins
5133  old node (for cuts?)
5134*/
5135  currentNumberCuts = 0;
5136  while (nodeInfo) {
5137    //printf("nNode = %d, nodeInfo = %x\n",nNode,nodeInfo);
5138    walkback_[nNode++]=nodeInfo;
5139    currentNumberCuts += nodeInfo->numberCuts() ;
5140    nodeInfo = nodeInfo->parent() ;
5141    if (nNode==maximumDepth_) {
5142      maximumDepth_ *= 2;
5143      CbcNodeInfo ** temp = new CbcNodeInfo * [maximumDepth_];
5144      for (i=0;i<nNode;i++) 
5145        temp[i] = walkback_[i];
5146      delete [] walkback_;
5147      walkback_ = temp;
5148    }
5149  }
5150/*
5151  Create an empty basis with sufficient capacity for the constraint system
5152  we'll construct: original system plus cuts. Make sure we have capacity to
5153  record those cuts in addedCuts_.
5154
5155  The method of adjusting the basis at a FullNodeInfo object (the root, for
5156  example) is to use a copy constructor to duplicate the basis held in the
5157  nodeInfo, then resize it and return the new basis object. Guaranteed,
5158  lastws will point to a different basis when it returns. We pass in a basis
5159  because we need the parameter to return the allocated basis, and it's an
5160  easy way to pass in the size. But we take a hit for memory allocation.
5161*/
5162  currentNumberCuts_=currentNumberCuts;
5163  if (currentNumberCuts > maximumNumberCuts_) {
5164    maximumNumberCuts_ = currentNumberCuts;
5165    delete [] addedCuts_;
5166    addedCuts_ = new CbcCountRowCut * [maximumNumberCuts_];
5167  }
5168  lastws->setSize(numberColumns,numberRowsAtContinuous_+currentNumberCuts);
5169/*
5170  This last bit of code traverses the path collected in walkback_ from the
5171  root back to node. At the end of the loop,
5172   * lastws will be an appropriate basis for node;
5173   * variable bounds in the constraint system will be set to be correct for
5174     node; and
5175   * addedCuts_ will be set to a list of cuts that need to be added to the
5176     constraint system at node.
5177  applyToModel does all the heavy lifting.
5178*/
5179  currentNumberCuts=0;
5180  //#define CBC_PRINT2
5181#ifdef CBC_PRINT2
5182  printf("Starting bounds at node %d\n",numberNodes_);
5183#endif
5184  while (nNode) {
5185    --nNode;
5186    walkback_[nNode]->applyToModel(this,lastws,addedCuts_,currentNumberCuts);
5187  }
5188  if (!lastws->fullBasis()) {
5189    printf("******* bad basis\n");
5190    int numberRows = lastws->getNumArtificial();
5191    int i;
5192    for (i=0;i<numberRows;i++)
5193      lastws->setArtifStatus(i,CoinWarmStartBasis::basic);
5194    int numberColumns = lastws->getNumStructural();
5195    for (i=0;i<numberColumns;i++) {
5196      if (lastws->getStructStatus(i)==CoinWarmStartBasis::basic)
5197        lastws->setStructStatus(i,CoinWarmStartBasis::atLowerBound);
5198    }
5199#if 0
5200  } else {
5201    // OPTION - take off slack cuts
5202    // First see if any cuts are slack
5203    int numberAdded = currentNumberCuts;
5204    if (saveNode<2&&false) {
5205      printf("nNode %d cuts %d\n",saveNode,currentNumberCuts);
5206      for (int i=0;i<currentNumberCuts;i++)
5207        addedCuts_[i]->print();
5208    }
5209    if (numberAdded&&saveNode<5&&!parentModel_) {
5210#if 0
5211      currentNumberCuts=0;
5212      for (int j=numberRowsAtContinuous_;
5213           j<numberAdded+numberRowsAtContinuous_;j++) {
5214        CoinWarmStartBasis::Status status = lastws->getArtifStatus(j);
5215        if (status!=CoinWarmStartBasis::basic) {
5216          lastws->setArtifStatus(currentNumberCuts+numberRowsAtContinuous_,
5217                                 status);
5218          addedCuts_[currentNumberCuts++]=addedCuts_[j-numberRowsAtContinuous_];
5219        }
5220      }
5221      if (currentNumberCuts<numberAdded) {
5222        printf("deleting %d rows\n",numberAdded-currentNumberCuts);
5223        lastws->resize(currentNumberCuts+numberRowsAtContinuous_,
5224                       lastws->getNumStructural());
5225        currentNumberCuts_=currentNumberCuts;
5226      }
5227#else
5228      int nDelete=0;
5229      for (int j=numberRowsAtContinuous_;
5230           j<numberAdded+numberRowsAtContinuous_;j++) {
5231        CoinWarmStartBasis::Status status = lastws->getArtifStatus(j);
5232        if (status==CoinWarmStartBasis::basic)
5233          nDelete++;
5234      }
5235      if (nDelete)
5236        printf("depth %d can delete %d\n",saveNode-1,nDelete);
5237#endif
5238    }
5239#endif
5240  }
5241  if (0) {
5242    int numberDebugValues=18;
5243    double * debugValues = new double[numberDebugValues];
5244    CoinZeroN(debugValues,numberDebugValues);
5245    debugValues[1]=6.0;
5246    debugValues[3]=60.0;
5247    debugValues[4]=6.0;
5248    debugValues[6]=60.0;
5249    debugValues[7]=16.0;
5250    debugValues[9]=70.0;
5251    debugValues[10]=7.0;
5252    debugValues[12]=70.0;
5253    debugValues[13]=12.0;
5254    debugValues[15]=75.0;
5255    int nBad=0;
5256    for (int j=0;j<numberColumns;j++) {
5257      if (integerInfo_[j]) {
5258        if(solver_->getColLower()[j]>debugValues[j]||
5259           solver_->getColUpper()[j]<debugValues[j]) {
5260          printf("** (%g) ** ",debugValues[j]);
5261          nBad++;
5262        }
5263        printf("%d bounds %g %g\n",j,solver_->getColLower()[j],solver_->getColUpper()[j]);
5264      }
5265    }
5266    if (nBad)
5267      printf("%d BAD\n",nBad);
5268    else
5269      printf("OKAY\n");
5270    delete [] debugValues;
5271  }
5272}
5273
5274/*
5275  adjustCuts might be a better name: If the node is feasible, we sift through
5276  the cuts collected by addCuts1, add the ones that are tight and omit the
5277  ones that are loose. If the node is infeasible, we just adjust the
5278  reference counts to reflect that we're about to prune this node and its
5279  descendants.
5280*/
5281int CbcModel::addCuts (CbcNode *node, CoinWarmStartBasis *&lastws,bool canFix)
5282{
5283/*
5284  addCuts1 performs step 1 of restoring the subproblem at this node; see the
5285  comments there.
5286*/
5287  addCuts1(node,lastws);
5288  int i;
5289  int numberColumns = getNumCols();
5290  if (solver_->getNumRows()>maximumRows_) {
5291    maximumRows_ = solver_->getNumRows();
5292    workingBasis_.resize(maximumRows_,numberColumns);
5293  }
5294  CbcNodeInfo * nodeInfo = node->nodeInfo();
5295  double cutoff = getCutoff() ;
5296  int currentNumberCuts=currentNumberCuts_;
5297  if (canFix) {
5298    bool feasible=true;
5299    const double *lower = solver_->getColLower() ;
5300    const double *upper = solver_->getColUpper() ;
5301    double * newLower = analyzeResults_;
5302    double * objLower = newLower+numberIntegers_;
5303    double * newUpper = objLower+numberIntegers_;
5304    double * objUpper = newUpper+numberIntegers_;
5305    int n=0;
5306    for (i=0;i<numberIntegers_;i++) {
5307      int iColumn = integerVariable_[i];
5308      bool changed=false;
5309      double lo = 0.0;
5310      double up = 0.0;
5311      if (objLower[i]>cutoff) {
5312        lo = lower[iColumn];
5313        up = upper[iColumn];
5314        if (lo<newLower[i]) {
5315          lo = newLower[i];
5316          solver_->setColLower(iColumn,lo);
5317          changed=true;
5318          n++;
5319        }
5320        if (objUpper[i]>cutoff) {
5321          if (up>newUpper[i]) {
5322            up = newUpper[i];
5323            solver_->setColUpper(iColumn,up);
5324            changed=true;
5325            n++;
5326          }
5327        }
5328      } else if (objUpper[i]>cutoff) {
5329        lo = lower[iColumn];
5330        up = upper[iColumn];
5331        if (up>newUpper[i]) {
5332          up = newUpper[i];
5333          solver_->setColUpper(iColumn,up);
5334          changed=true;
5335          n++;
5336        }
5337      }
5338      if (changed&&lo>up) {
5339        feasible=false;
5340        break;
5341      }
5342    }
5343    if (!feasible) {
5344      printf("analysis says node infeas\n");
5345      cutoff=-COIN_DBL_MAX;
5346    }
5347  }
5348/*
5349  If the node can't be fathomed by bound, reinstall tight cuts in the
5350  constraint system. Even if there are no cuts, we'll want to set the
5351  reconstructed basis in the solver.
5352*/
5353  if (node->objectiveValue() < cutoff||numberThreads_)
5354  { 
5355    //#   define CBC_CHECK_BASIS
5356#   ifdef CBC_CHECK_BASIS
5357    printf("addCuts: expanded basis; rows %d+%d\n",
5358           numberRowsAtContinuous_,currentNumberCuts);
5359    lastws->print();
5360#   endif
5361/*
5362  Adjust the basis and constraint system so that we retain only active cuts.
5363  There are three steps:
5364    1) Scan the basis. Sort the cuts into effective cuts to be kept and
5365       loose cuts to be dropped.
5366    2) Drop the loose cuts and resize the basis to fit.
5367    3) Install the tight cuts in the constraint system (applyRowCuts) and
5368       and install the basis (setWarmStart).
5369  Use of compressRows conveys we're compressing the basis and not just
5370  tweaking the artificialStatus_ array.
5371*/
5372    if (currentNumberCuts > 0) {
5373      int numberToAdd = 0;
5374      const OsiRowCut **addCuts;
5375      int numberToDrop = 0 ;
5376      int *cutsToDrop ;
5377      addCuts = new const OsiRowCut* [currentNumberCuts];
5378      cutsToDrop = new int[currentNumberCuts] ;
5379      assert (currentNumberCuts+numberRowsAtContinuous_<=lastws->getNumArtificial());
5380      for (i=0;i<currentNumberCuts;i++) {
5381        CoinWarmStartBasis::Status status = 
5382          lastws->getArtifStatus(i+numberRowsAtContinuous_);
5383        if (addedCuts_[i] &&
5384            (status != CoinWarmStartBasis::basic ||
5385             addedCuts_[i]->effectiveness()==COIN_DBL_MAX)) {
5386#         ifdef CHECK_CUT_COUNTS
5387          printf("Using cut %d %x as row %d\n",i,addedCuts_[i],
5388                 numberRowsAtContinuous_+numberToAdd);
5389#         endif
5390          addCuts[numberToAdd++] = addedCuts_[i];
5391        } else {
5392#         ifdef CHECK_CUT_COUNTS
5393          printf("Dropping cut %d %x\n",i,addedCuts_[i]);
5394#         endif
5395          addedCuts_[i]=NULL;
5396          cutsToDrop[numberToDrop++] = numberRowsAtContinuous_+i ;
5397        }
5398      }
5399      int numberRowsNow=numberRowsAtContinuous_+numberToAdd;
5400      lastws->compressRows(numberToDrop,cutsToDrop) ;
5401      lastws->resize(numberRowsNow,numberColumns);
5402      solver_->applyRowCuts(numberToAdd,addCuts);
5403#     ifdef CBC_CHECK_BASIS
5404      printf("addCuts: stripped basis; rows %d + %d\n",
5405             numberRowsAtContinuous_,numberToAdd);
5406      lastws->print();
5407#     endif
5408      //for (i=0;i<numberToAdd;i++)
5409      //delete addCuts[i];
5410      delete [] addCuts;
5411      delete [] cutsToDrop ;
5412    }
5413/*
5414  Set the basis in the solver.
5415*/
5416    solver_->setWarmStart(lastws);
5417#if 0
5418    if ((numberNodes_%printFrequency_)==0) {
5419      printf("Objective %g, depth %d, unsatisfied %d\n",
5420             node->objectiveValue(),
5421             node->depth(),node->numberUnsatisfied());
5422    }
5423#endif
5424/*
5425  Clean up and we're out of here.
5426*/
5427    numberNodes_++;
5428    return 0;
5429  } 
5430/*
5431  This node has been fathomed by bound as we try to revive it out of the live
5432  set. Adjust the cut reference counts to reflect that we no longer need to
5433  explore the remaining branch arms, hence they will no longer reference any
5434  cuts. Cuts whose reference count falls to zero are deleted. 
5435*/
5436  else
5437  { int i;
5438    if (currentNumberCuts) {
5439#ifndef CBC_DETERMINISTIC_THREAD
5440      lockThread();
5441#endif
5442      int numberLeft = nodeInfo->numberBranchesLeft();
5443      for (i = 0 ; i < currentNumberCuts ; i++)
5444        { if (addedCuts_[i])
5445          { if (!addedCuts_[i]->decrement(numberLeft))
5446            { delete addedCuts_[i];
5447            addedCuts_[i] = NULL; } } }
5448#ifndef CBC_DETERMINISTIC_THREAD
5449      unlockThread();
5450#endif
5451    }
5452    return 1 ; }
5453}
5454
5455
5456/*
5457  Perform reduced cost fixing on integer variables.
5458
5459  The variables in question are already nonbasic at bound. We're just nailing
5460  down the current situation.
5461*/
5462
5463int CbcModel::reducedCostFix ()
5464
5465{
5466  if(!solverCharacteristics_->reducedCostsAccurate())
5467    return 0; //NLP
5468  double cutoff = getCutoff() ;
5469  double direction = solver_->getObjSense() ;
5470  double gap = cutoff - solver_->getObjValue()*direction ;
5471  double tolerance;
5472  solver_->getDblParam(OsiDualTolerance,tolerance) ;
5473  if (gap<=0.0)
5474    gap = tolerance; //return 0;
5475  gap += 100.0*tolerance;
5476  double integerTolerance = getDblParam(CbcIntegerTolerance) ;
5477
5478  const double *lower = solver_->getColLower() ;
5479  const double *upper = solver_->getColUpper() ;
5480  const double *solution = solver_->getColSolution() ;
5481  const double *reducedCost = solver_->getReducedCost() ;
5482
5483  int numberFixed = 0 ;
5484
5485# ifdef COIN_HAS_CLP
5486  OsiClpSolverInterface * clpSolver
5487    = dynamic_cast<OsiClpSolverInterface *> (solver_);
5488  ClpSimplex * clpSimplex=NULL;
5489  if (clpSolver) 
5490    clpSimplex = clpSolver->getModelPtr();
5491# endif
5492  for (int i = 0 ; i < numberIntegers_ ; i++) {
5493    int iColumn = integerVariable_[i] ;
5494    double djValue = direction*reducedCost[iColumn] ;
5495    if (upper[iColumn]-lower[iColumn] > integerTolerance) {
5496      if (solution[iColumn] < lower[iColumn]+integerTolerance && djValue > gap) {
5497#ifdef COIN_HAS_CLP
5498        // may just have been fixed before
5499        if (clpSimplex) {
5500          if (clpSimplex->getColumnStatus(iColumn)==ClpSimplex::basic) {
5501#ifdef COIN_DEVELOP
5502            printf("DJfix %d has status of %d, dj of %g gap %g, bounds %g %g\n",
5503                   iColumn,clpSimplex->getColumnStatus(iColumn),
5504                   djValue,gap,lower[iColumn],upper[iColumn]);
5505#endif
5506          } else {         
5507            assert(clpSimplex->getColumnStatus(iColumn)==ClpSimplex::atLowerBound||
5508                   clpSimplex->getColumnStatus(iColumn)==ClpSimplex::isFixed);
5509          }
5510        }
5511#endif
5512        solver_->setColUpper(iColumn,lower[iColumn]) ;
5513        numberFixed++ ;
5514      } else if (solution[iColumn] > upper[iColumn]-integerTolerance && -djValue > gap) {
5515#ifdef COIN_HAS_CLP
5516        // may just have been fixed before
5517        if (clpSimplex) {
5518          if (clpSimplex->getColumnStatus(iColumn)==ClpSimplex::basic) {
5519#ifdef COIN_DEVELOP
5520            printf("DJfix %d has status of %d, dj of %g gap %g, bounds %g %g\n",
5521                   iColumn,clpSimplex->getColumnStatus(iColumn),
5522                   djValue,gap,lower[iColumn],upper[iColumn]);
5523#endif
5524          } else {         
5525            assert(clpSimplex->getColumnStatus(iColumn)==ClpSimplex::atUpperBound||
5526                   clpSimplex->getColumnStatus(iColumn)==ClpSimplex::isFixed);
5527          }
5528        }
5529#endif
5530        solver_->setColLower(iColumn,upper[iColumn]) ;
5531        numberFixed++ ;
5532      }
5533    }
5534  }
5535  numberDJFixed_ += numberFixed;
5536  return numberFixed; }
5537
5538// Collect coding to replace whichGenerator
5539void
5540CbcModel::resizeWhichGenerator(int numberNow, int numberAfter)
5541{
5542  if (numberAfter > maximumWhich_) {
5543    maximumWhich_ = CoinMax(maximumWhich_*2+100,numberAfter) ;
5544    int * temp = new int[2*maximumWhich_] ;
5545    memcpy(temp,whichGenerator_,numberNow*sizeof(int)) ;
5546    delete [] whichGenerator_ ;
5547    whichGenerator_ = temp ;
5548    memset(whichGenerator_+numberNow,0,(maximumWhich_-numberNow)*sizeof(int));
5549  }
5550}
5551
5552/** Solve the model using cuts
5553
5554  This version takes off redundant cuts from node.
5555  Returns true if feasible.
5556
5557  \todo
5558  Why do I need to resolve the problem? What has been done between the last
5559  relaxation and calling solveWithCuts?
5560
5561  If numberTries == 0 then user did not want any cuts.
5562*/
5563
5564bool 
5565CbcModel::solveWithCuts (OsiCuts &cuts, int numberTries, CbcNode *node)
5566/*
5567  Parameters:
5568    numberTries: (i) the maximum number of iterations for this round of cut
5569                     generation; if negative then we don't mind if drop is tiny.
5570   
5571    cuts:       (o) all cuts generated in this round of cut generation
5572
5573    node: (i)     So we can update dynamic pseudo costs
5574*/
5575                       
5576
5577{
5578# ifdef COIN_HAS_CLP
5579  OsiClpSolverInterface * clpSolver
5580    = dynamic_cast<OsiClpSolverInterface *> (solver_);
5581  int saveClpOptions=0;
5582  if (clpSolver) 
5583    saveClpOptions = clpSolver->specialOptions();
5584# endif
5585  //solver_->writeMps("saved");
5586#ifdef CBC_THREAD
5587  CbcModel ** threadModel = NULL;
5588  Coin_pthread_t * threadId = NULL;
5589  pthread_cond_t condition_main;
5590  pthread_mutex_t condition_mutex;
5591  pthread_mutex_t * mutex2 = NULL;
5592  pthread_cond_t * condition2 = NULL;
5593  threadStruct * threadInfo = NULL;
5594  void * saveMutex = NULL;
5595  if (numberThreads_&&(threadMode_&2)!=0&&!numberNodes_) {
5596    threadId = new Coin_pthread_t [numberThreads_];
5597    pthread_cond_init(&condition_main,NULL);
5598    pthread_mutex_init(&condition_mutex,NULL);
5599    threadModel = new CbcModel * [numberThreads_];
5600    threadInfo = new threadStruct [numberThreads_+1];
5601    mutex2 = new pthread_mutex_t [numberThreads_];
5602    condition2 = new pthread_cond_t [numberThreads_];
5603    saveMutex = mutex_;
5604    for (int i=0;i<numberThreads_;i++) {
5605      pthread_mutex_init(mutex2+i,NULL);
5606      pthread_cond_init(condition2+i,NULL);
5607      threadId[i].status=0;
5608      threadModel[i]=new CbcModel;
5609      threadModel[i]->generator_ = new CbcCutGenerator * [1];
5610      delete threadModel[i]->solver_;
5611      threadModel[i]->solver_=NULL;
5612      threadModel[i]->numberThreads_=numberThreads_;
5613      mutex_ = (void *) (threadInfo+i);
5614      threadInfo[i].thisModel=(CbcModel *) threadModel[i];
5615      threadInfo[i].baseModel=this;
5616      threadInfo[i].threadIdOfBase.thr=pthread_self();
5617      threadInfo[i].mutex2=mutex2+i;
5618      threadInfo[i].condition2=condition2+i;
5619      threadInfo[i].returnCode=-1;
5620      pthread_create(&(threadId[i].thr),NULL,doCutsThread,threadInfo+i);
5621      threadId[i].status = 1;
5622    }
5623    // Do a partial one for base model
5624    threadInfo[numberThreads_].baseModel=this;
5625    mutex_ = (void *) (threadInfo+numberThreads_);
5626    threadInfo[numberThreads_].condition2=&condition_main;
5627    threadInfo[numberThreads_].mutex2=&condition_mutex;
5628  }
5629#endif
5630  bool feasible = true ;
5631  int lastNumberCuts = 0 ;
5632  int violated = 0 ;
5633  int numberRowsAtStart = solver_->getNumRows() ;
5634  //printf("solver had %d rows\n",numberRowsAtStart);
5635  int numberColumns = solver_->getNumCols() ;
5636  CoinBigIndex numberElementsAtStart = solver_->getNumElements();
5637
5638  numberOldActiveCuts_ = numberRowsAtStart-numberRowsAtContinuous_ ;
5639  numberNewCuts_ = 0 ;
5640
5641  bool onOptimalPath = false ;
5642  const OsiRowCutDebugger *debugger = NULL;
5643  if ((specialOptions_&1)!=0) {
5644    /*
5645      See OsiRowCutDebugger for details. In a nutshell, make sure that current
5646      variable values do not conflict with a known optimal solution. (Obviously
5647      this can be fooled when there are multiple solutions.)
5648    */
5649    debugger = solver_->getRowCutDebugger() ;
5650    if (debugger) 
5651      onOptimalPath = (debugger->onOptimalPath(*solver_)) ;
5652  }
5653  OsiCuts slackCuts;
5654/*
5655  Resolve the problem. If we've lost feasibility, might as well bail out right
5656  after the debug stuff. The resolve will also refresh cached copies of the
5657  solver solution (cbcColLower_, ...) held by CbcModel.
5658*/
5659  double objectiveValue = solver_->getObjValue()*solver_->getObjSense();
5660  if (node)
5661    objectiveValue= node->objectiveValue();
5662  int returnCode = resolve(node ? node->nodeInfo() : NULL,1);
5663#if COIN_DEVELOP>1
5664  //if (!solver_->getIterationCount()&&solver_->isProvenOptimal())
5665  //printf("zero iterations on first solve of branch\n");
5666#endif
5667  double lastObjective = solver_->getObjValue()*solver_->getObjSense();
5668  //double firstObjective = lastObjective+1.0e-8+1.0e-12*fabs(lastObjective);
5669  if (node&&node->nodeInfo()&&!node->nodeInfo()->numberBranchesLeft())
5670    node->nodeInfo()->allBranchesGone(); // can clean up
5671  feasible = returnCode  != 0 ;
5672  if (returnCode<0)
5673    numberTries=0;
5674  if (problemFeasibility_->feasible(this,0)<0) {
5675    feasible=false; // pretend infeasible
5676  }
5677 
5678#if NEW_UPDATE_OBJECT==0
5679  // Update branching information if wanted
5680  if(node &&branchingMethod_)
5681    branchingMethod_->updateInformation(solver_,node);
5682#elif NEW_UPDATE_OBJECT<2
5683  // Update branching information if wanted
5684  if(node &&branchingMethod_) {
5685    OsiBranchingObject * bobj = node->modifiableBranchingObject();
5686    CbcBranchingObject * cbcobj = dynamic_cast<CbcBranchingObject *> (bobj);
5687    if (cbcobj) {
5688      CbcObject * object = cbcobj->object();
5689      CbcObjectUpdateData update = object->createUpdateInformation(solver_,node,cbcobj);
5690      object->updateInformation(update);
5691    } else {
5692      branchingMethod_->updateInformation(solver_,node);
5693    }
5694  }
5695#else
5696  // Update branching information if wanted
5697  if(node &&branchingMethod_) {
5698    OsiBranchingObject * bobj = node->modifiableBranchingObject();
5699    CbcBranchingObject * cbcobj = dynamic_cast<CbcBranchingObject *> (bobj);
5700    if (cbcobj&&cbcobj->object()) {
5701      CbcObject * object = cbcobj->object();
5702      CbcObjectUpdateData update = object->createUpdateInformation(solver_,node,cbcobj);
5703      // have to compute object number as not saved
5704      CbcSimpleInteger * simpleObject =
5705          dynamic_cast <CbcSimpleInteger *>(object) ;
5706      int iObject;
5707      int iColumn = simpleObject->columnNumber();
5708      for (iObject = 0 ; iObject < numberObjects_ ; iObject++) {
5709        simpleObject =
5710          dynamic_cast <CbcSimpleInteger *>(object_[iObject]) ;
5711        if (simpleObject->columnNumber()==iColumn)
5712          break;
5713      }
5714      assert (iObject<numberObjects_);
5715      update.objectNumber_ = iObject;
5716      addUpdateInformation(update);
5717    } else {
5718      OsiIntegerBranchingObject * obj = dynamic_cast<OsiIntegerBranchingObject *> (bobj);
5719      if (obj) {
5720        const OsiObject * object = obj->originalObject();
5721        // have to compute object number as not saved
5722        int iObject;
5723        int iColumn = object->columnNumber();
5724        for (iObject = 0 ; iObject < numberObjects_ ; iObject++) {
5725          if (object_[iObject]->columnNumber()==iColumn)
5726            break;
5727        }
5728        assert (iObject<numberObjects_);
5729        int branch = obj->firstBranch();
5730        if (obj->branchIndex()==2)
5731          branch = 1-branch;
5732        assert (branch==0||branch==1);
5733        double originalValue=node->objectiveValue();
5734        double objectiveValue = solver_->getObjValue()*solver_->getObjSense();
5735        double changeInObjective = CoinMax(0.0,objectiveValue-originalValue);
5736        int iStatus = (feasible) ? 0 : 0;
5737        double value = obj->value();
5738        double movement;
5739        if (branch)
5740          movement = ceil(value)-value;
5741        else
5742          movement = value -floor(value);
5743#if 0
5744        // OUT as much too complicated - we are not at a natural hotstart place
5745        OsiBranchingInformation usefulInfo=usefulInformation();
5746        // hotInfo is meant for BEFORE a branch so we need to fool
5747        // was much simpler with alternate method
5748        double save[3];
5749        save[0]=usefulInfo.lower_[iColumn];
5750        save[1]=usefulInfo.solution_[iColumn];
5751        save[2]=usefulInfo.upper_[iColumn];
5752        usefulInfo.lower_[iColumn]=floor(value);
5753        usefulInfo.solution_[iColumn]=value;
5754        usefulInfo.upper_[iColumn]=ceil(value);
5755        OsiHotInfo hotInfo(solver_,&usefulInfo,&object,0);
5756        usefulInfo.lower_[iColumn]=save[0];
5757        usefulInfo.solution_[iColumn]=save[1];
5758        usefulInfo.upper_[iColumn]=save[2];
5759        if (branch) {
5760          hotInfo.setUpStatus(iStatus);
5761          hotInfo.setUpChange(changeInObjective);
5762          //object->setUpEstimate(movement);
5763        } else {
5764          hotInfo.setDownStatus(iStatus);
5765          hotInfo.setDownChange(changeInObjective);
5766          //object->setDownEstimate(movement);
5767        }
5768        branchingMethod_->chooseMethod()->updateInformation(&usefulInfo,branch,&hotInfo);
5769#else
5770        branchingMethod_->chooseMethod()->updateInformation(iObject,branch,changeInObjective,
5771                                                            movement,iStatus);
5772#endif
5773      }
5774    }
5775  }
5776#endif
5777
5778#ifdef CBC_DEBUG
5779  if (feasible)
5780  { printf("Obj value %g (%s) %d rows\n",solver_->getObjValue(),
5781           (solver_->isProvenOptimal())?"proven":"unproven",
5782           solver_->getNumRows()) ; }
5783 
5784  else
5785  { printf("Infeasible %d rows\n",solver_->getNumRows()) ; }
5786#endif
5787  if ((specialOptions_&1)!=0) {
5788/*
5789  If the RowCutDebugger said we were compatible with the optimal solution,
5790  and now we're suddenly infeasible, we might be confused. Then again, we
5791  may have fathomed by bound, heading for a rediscovery of an optimal solution.
5792*/
5793    if (onOptimalPath && !solver_->isDualObjectiveLimitReached()) {
5794      if (!feasible) {
5795        solver_->writeMpsNative("infeas.mps",NULL,NULL,2);
5796        solver_->getRowCutDebuggerAlways()->printOptimalSolution(*solver_);
5797        CoinWarmStartBasis *slack =
5798          dynamic_cast<CoinWarmStartBasis *>(solver_->getEmptyWarmStart()) ;
5799        solver_->setWarmStart(slack);
5800        delete slack ;
5801        solver_->setHintParam(OsiDoReducePrint,false,OsiHintDo,0) ;
5802        solver_->initialSolve();
5803      }
5804      assert(feasible) ;
5805    }
5806  }
5807
5808  if (!feasible) {
5809    numberInfeasibleNodes_++;
5810# ifdef COIN_HAS_CLP
5811    if (clpSolver) 
5812    clpSolver->setSpecialOptions(saveClpOptions);
5813# endif
5814    return (false) ;
5815  }
5816  sumChangeObjective1_ += solver_->getObjValue()*solver_->getObjSense()
5817    - objectiveValue ;
5818  if ( getCurrentSeconds() > dblParam_[CbcMaximumSeconds] )
5819    numberTries=0; // exit
5820  //if ((numberNodes_%100)==0)
5821  //printf("XXa sum obj changed by %g\n",sumChangeObjective1_);
5822  objectiveValue = solver_->getObjValue()*solver_->getObjSense();
5823  // Return at once if numberTries zero
5824  if (!numberTries) {
5825    cuts=OsiCuts();
5826    numberNewCuts_=0;
5827# ifdef COIN_HAS_CLP
5828    if (clpSolver) 
5829    clpSolver->setSpecialOptions(saveClpOptions);
5830# endif
5831    return true;
5832  }
5833/*
5834  Do reduced cost fixing.
5835*/
5836  int xxxxxx=0;
5837  if(xxxxxx)
5838  solver_->resolve();
5839  reducedCostFix() ;
5840/*
5841  Set up for at most numberTries rounds of cut generation. If numberTries is
5842  negative, we'll ignore the minimumDrop_ cutoff and keep generating cuts for
5843  the specified number of rounds.
5844*/
5845  double minimumDrop = minimumDrop_ ;
5846  if (numberTries<0)
5847  { numberTries = -numberTries ;
5848    minimumDrop = -1.0 ; }
5849/*
5850  Is it time to scan the cuts in order to remove redundant cuts? If so, set
5851  up to do it.
5852*/
5853  int fullScan = 0 ;
5854  if ((numberNodes_%SCANCUTS) == 0||(specialOptions_&256)!=0)
5855  { fullScan = 1 ;
5856    if (!numberNodes_||(specialOptions_&256)!=0)
5857      fullScan=2;
5858    specialOptions_ &= ~256; // mark as full scan done
5859  }
5860
5861  double direction = solver_->getObjSense() ;
5862  double startObjective = solver_->getObjValue()*direction ;
5863
5864  currentPassNumber_ = 0 ;
5865  double primalTolerance = 1.0e-7 ;
5866  // We may need to keep going on
5867  bool keepGoing=false;
5868/*
5869  Begin cut generation loop. Cuts generated during each iteration are
5870  collected in theseCuts. The loop can be divided into four phases:
5871   1) Prep: Fix variables using reduced cost. In the first iteration only,
5872      consider scanning globalCuts_ and activating any applicable cuts.
5873   2) Cut Generation: Call each generator and heuristic registered in the
5874      generator_ and heuristic_ arrays. Newly generated global cuts are
5875      copied to globalCuts_ at this time.
5876   3) Cut Installation and Reoptimisation: Install column and row cuts in
5877      the solver. Copy row cuts to cuts (parameter). Reoptimise.
5878   4) Cut Purging: takeOffCuts() removes inactive cuts from the solver, and
5879      does the necessary bookkeeping in the model.
5880*/
5881  do
5882  { currentPassNumber_++ ;
5883    numberTries--