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

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

random number generator to CbcModel? and allow cut pruning on size

  • Property svn:eol-style set to native
  • Property svn:keywords set to Author Date Id Revision
File size: 461.1 KB
Line 
1// Copyright (C) 2002, International Business Machines
2// Corporation and others.  All Rights Reserved.
3#if defined(_MSC_VER)
4// Turn off compiler warning about long names
5#  pragma warning(disable:4786)
6#endif
7//#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  currentDepth_=0;
1511  workingBasis_.resize(maximumRows_,numberColumns);
1512/*
1513  Set up an empty heap and associated data structures to hold the live set
1514  (problems which require further exploration).
1515*/
1516  tree_->setComparison(*nodeCompare_) ;
1517/*
1518  Used to record the path from a node to the root of the search tree, so that
1519  we can then traverse from the root to the node when restoring a subproblem.
1520*/
1521  maximumDepth_ = 10 ;
1522  delete [] walkback_ ;
1523  walkback_ = new CbcNodeInfo * [maximumDepth_] ;
1524/*
1525  Used to generate bound edits for CbcPartialNodeInfo.
1526*/
1527  double * lowerBefore = new double [numberColumns] ;
1528  double * upperBefore = new double [numberColumns] ;
1529/*
1530 
1531  Generate cuts at the root node and reoptimise. solveWithCuts does the heavy
1532  lifting. It will iterate a generate/reoptimise loop (including reduced cost
1533  fixing) until no cuts are generated, the change in objective falls off,  or
1534  the limit on the number of rounds of cut generation is exceeded.
1535
1536  At the end of all this, any cuts will be recorded in cuts and also
1537  installed in the solver's constraint system. We'll have reoptimised, and
1538  removed any slack cuts (numberOldActiveCuts_ and numberNewCuts_ have been
1539  adjusted accordingly).
1540
1541  Tell cut generators they can be a bit more aggressive at root node
1542
1543  TODO: Why don't we make a copy of the solution after solveWithCuts?
1544  TODO: If numberUnsatisfied == 0, don't we have a solution?
1545*/
1546  phase_=1;
1547  int iCutGenerator;
1548  for (iCutGenerator = 0;iCutGenerator<numberCutGenerators_;iCutGenerator++) {
1549    CglCutGenerator * generator = generator_[iCutGenerator]->generator();
1550    generator->setAggressiveness(generator->getAggressiveness()+100);
1551  }
1552  OsiCuts cuts ;
1553  int anyAction = -1 ;
1554  numberOldActiveCuts_ = 0 ;
1555  numberNewCuts_ = 0 ;
1556  // Array to mark solution
1557  delete [] usedInSolution_;
1558  usedInSolution_ = new int[numberColumns];
1559  CoinZeroN(usedInSolution_,numberColumns);
1560/*
1561  For printing totals and for CbcNode (numberNodes_)
1562*/
1563  numberIterations_ = 0 ;
1564  numberNodes_ = 0 ;
1565  numberNodes2_ = 0 ;
1566  maximumStatistics_=0;
1567  maximumDepthActual_=0;
1568  numberDJFixed_=0.0;
1569  // Do heuristics
1570  doHeuristicsAtRoot();
1571  if ( intParam_[CbcMaxNumNode] < 0)
1572    eventHappened_=true; // stop as fast as possible
1573  statistics_ = NULL;
1574  // Do on switch
1575  if (doStatistics>0&&doStatistics<100) {
1576    maximumStatistics_=10000;
1577    statistics_ = new CbcStatistics * [maximumStatistics_];
1578    memset(statistics_,0,maximumStatistics_*sizeof(CbcStatistics *));
1579  }
1580
1581  { int iObject ;
1582    int preferredWay ;
1583    int numberUnsatisfied = 0 ;
1584    delete [] currentSolution_;
1585    currentSolution_ = new double [numberColumns];
1586    testSolution_ = currentSolution_;
1587    memcpy(currentSolution_,solver_->getColSolution(),
1588           numberColumns*sizeof(double)) ;
1589    // point to useful information
1590    OsiBranchingInformation usefulInfo=usefulInformation();
1591
1592    for (iObject = 0 ; iObject < numberObjects_ ; iObject++)
1593    { double infeasibility =
1594          object_[iObject]->infeasibility(&usefulInfo,preferredWay) ;
1595      if (infeasibility ) numberUnsatisfied++ ; }
1596    // replace solverType
1597    if(solverCharacteristics_->tryCuts())  {
1598
1599      if (numberUnsatisfied)   {
1600        // User event
1601        if (!eventHappened_)
1602          feasible = solveWithCuts(cuts,maximumCutPassesAtRoot_,
1603                                   NULL);
1604        else
1605          feasible=false;
1606      } else if (solverCharacteristics_->solutionAddsCuts()||
1607                 solverCharacteristics_->alwaysTryCutsAtRootNode()) {
1608        // may generate cuts and turn the solution
1609        //to an infeasible one
1610        feasible = solveWithCuts(cuts, 1,
1611                                 NULL);
1612      }
1613    }
1614    // check extra info on feasibility
1615    if (!solverCharacteristics_->mipFeasible())
1616      feasible = false;
1617  }
1618  // make cut generators less aggressive
1619  for (iCutGenerator = 0;iCutGenerator<numberCutGenerators_;iCutGenerator++) {
1620    CglCutGenerator * generator = generator_[iCutGenerator]->generator();
1621    generator->setAggressiveness(generator->getAggressiveness()-100);
1622  }
1623  currentNumberCuts_ = numberNewCuts_ ;
1624  // See if can stop on gap
1625  stoppedOnGap_ = false ;
1626  bestPossibleObjective_ = solver_->getObjValue()*solver_->getObjSense();
1627  double testGap = CoinMax(dblParam_[CbcAllowableGap],
1628                           CoinMax(fabs(bestObjective_),fabs(bestPossibleObjective_))
1629                           *dblParam_[CbcAllowableFractionGap]);
1630  if (bestObjective_-bestPossibleObjective_ < testGap && getCutoffIncrement()>=0.0) {
1631    if (bestPossibleObjective_<getCutoff())
1632      stoppedOnGap_ = true ;
1633    feasible = false;
1634  }
1635  // User event
1636  if (eventHappened_)
1637    feasible=false;
1638  if(fastNodeDepth_>=0&&!parentModel_) {
1639    // add in a general depth object
1640    CbcObject * obj = 
1641      new CbcGeneralDepth(this,fastNodeDepth_);
1642    addObjects(1,&obj);
1643    // mark as done
1644    fastNodeDepth_ += 1000000;
1645    delete obj;
1646    // fake number of objects
1647    numberObjects_--;
1648  }
1649#ifdef COIN_HAS_CLP
1650#ifdef NO_CRUNCH
1651 {
1652   OsiClpSolverInterface * clpSolver
1653     = dynamic_cast<OsiClpSolverInterface *> (solver_);
1654   if (clpSolver&&!parentModel_) {
1655     ClpSimplex * clpSimplex = clpSolver->getModelPtr();
1656     clpSimplex->setSpecialOptions(clpSimplex->specialOptions()|131072);
1657     //clpSimplex->startPermanentArrays();
1658     clpSimplex->setPersistenceFlag(2);
1659   }
1660 }
1661#endif
1662#endif
1663 // Save copy of solver
1664 OsiSolverInterface * saveSolver = NULL;
1665 if (!parentModel_&&(specialOptions_&(512+2048))!=0)
1666   saveSolver = solver_->clone();
1667 if (!parentModel_&&(specialOptions_&2048)!=0) {
1668   // See if worth trying reduction
1669   bool tryNewSearch=solverCharacteristics_->reducedCostsAccurate();
1670   int numberColumns = getNumCols();
1671   if (tryNewSearch) {
1672     double cutoff = getCutoff() ;
1673     saveSolver->resolve();
1674     double direction = saveSolver->getObjSense() ;
1675     double gap = cutoff - saveSolver->getObjValue()*direction ;
1676     double tolerance;
1677     saveSolver->getDblParam(OsiDualTolerance,tolerance) ;
1678     if (gap<=0.0)
1679       gap = tolerance; 
1680     gap += 100.0*tolerance;
1681     double integerTolerance = getDblParam(CbcIntegerTolerance) ;
1682     
1683     const double *lower = saveSolver->getColLower() ;
1684     const double *upper = saveSolver->getColUpper() ;
1685     const double *solution = saveSolver->getColSolution() ;
1686     const double *reducedCost = saveSolver->getReducedCost() ;
1687     
1688     int numberFixed = 0 ;
1689     int numberFixed2=0;
1690     for (int i = 0 ; i < numberIntegers_ ; i++) {
1691       int iColumn = integerVariable_[i] ;
1692       double djValue = direction*reducedCost[iColumn] ;
1693       if (upper[iColumn]-lower[iColumn] > integerTolerance) {
1694         if (solution[iColumn] < lower[iColumn]+integerTolerance && djValue > gap) {
1695           saveSolver->setColUpper(iColumn,lower[iColumn]) ;
1696           numberFixed++ ;
1697         } else if (solution[iColumn] > upper[iColumn]-integerTolerance && -djValue > gap) {
1698           saveSolver->setColLower(iColumn,upper[iColumn]) ;
1699           numberFixed++ ;
1700         }
1701       } else {
1702         numberFixed2++;
1703       }
1704     }
1705#ifdef COIN_DEVELOP
1706     if ((specialOptions_&1)!=0) {
1707       const OsiRowCutDebugger *debugger = saveSolver->getRowCutDebugger() ;
1708       if (debugger) { 
1709         printf("Contains optimal\n") ;
1710         saveSolver->writeMps("reduced");
1711       } else {
1712         abort();
1713       }
1714     }
1715     printf("Restart could fix %d integers (%d already fixed)\n",
1716            numberFixed+numberFixed2,numberFixed2);
1717#endif
1718     numberFixed += numberFixed2;
1719     if (numberFixed*20<numberColumns)
1720       tryNewSearch=false;
1721   }
1722   if (tryNewSearch) {
1723     // back to solver without cuts?
1724#if 0
1725     OsiSolverInterface * solver2 = continuousSolver_->clone();
1726#else
1727     OsiSolverInterface * solver2 = saveSolver->clone();
1728#endif
1729     const double *lower = saveSolver->getColLower() ;
1730     const double *upper = saveSolver->getColUpper() ;
1731     for (int i = 0 ; i < numberIntegers_ ; i++) {
1732       int iColumn = integerVariable_[i] ;
1733       solver2->setColLower(iColumn,lower[iColumn]);
1734       solver2->setColUpper(iColumn,upper[iColumn]);
1735     }
1736     // swap
1737     delete saveSolver;
1738     saveSolver=solver2;
1739     double * newSolution = new double[numberColumns];
1740     double objectiveValue=cutoff;
1741     CbcSerendipity heuristic(*this);
1742     if (bestSolution_)
1743       heuristic.setInputSolution(bestSolution_,bestObjective_);
1744     heuristic.setFractionSmall(0.9);
1745     heuristic.setFeasibilityPumpOptions(1008013);
1746     // Use numberNodes to say how many are original rows
1747     heuristic.setNumberNodes(continuousSolver_->getNumRows());
1748     if (continuousSolver_->getNumRows()<
1749         solver_->getNumRows())
1750       printf("%d rows added ZZZZZ\n",
1751              solver_->getNumRows()-continuousSolver_->getNumRows());
1752     int returnCode= heuristic.smallBranchAndBound(saveSolver,
1753                                                   -1,newSolution,
1754                                                   objectiveValue,
1755                                                   cutoff,"Reduce");
1756     if (returnCode==-1) {
1757#ifdef COIN_DEVELOP
1758       printf("Restart - not small enough to do search after fixing\n");
1759#endif
1760       delete [] newSolution;
1761     } else {
1762       assert (returnCode>=0);
1763       if ((returnCode&1)!=0) {
1764         // increment number of solutions so other heuristics can test
1765         numberSolutions_++;
1766         numberHeuristicSolutions_++;
1767         lastHeuristic_ = NULL;
1768         setBestSolution(CBC_ROUNDING,objectiveValue,newSolution) ;
1769       }
1770       delete [] newSolution;
1771       feasible=false; // stop search
1772     }
1773   } 
1774 }
1775/*
1776  We've taken the continuous relaxation as far as we can. Time to branch.
1777  The first order of business is to actually create a node. chooseBranch
1778  currently uses strong branching to evaluate branch object candidates,
1779  unless forced back to simple branching. If chooseBranch concludes that a
1780  branching candidate is monotone (anyAction == -1) or infeasible (anyAction
1781  == -2) when forced to integer values, it returns here immediately.
1782
1783  Monotone variables trigger a call to resolve(). If the problem remains
1784  feasible, try again to choose a branching variable. At the end of the loop,
1785  resolved == true indicates that some variables were fixed.
1786
1787  Loss of feasibility will result in the deletion of newNode.
1788*/
1789
1790  bool resolved = false ;
1791  CbcNode *newNode = NULL ;
1792  numberFixedAtRoot_=0;
1793  numberFixedNow_=0;
1794  int numberIterationsAtContinuous = numberIterations_;
1795  //solverCharacteristics_->setSolver(solver_);
1796  if (feasible) {
1797    if (probingInfo_) {
1798      int number01 = probingInfo_->numberIntegers();
1799      //const fixEntry * entry = probingInfo_->fixEntries();
1800      const int * toZero = probingInfo_->toZero();
1801      //const int * toOne = probingInfo_->toOne();
1802      //const int * integerVariable = probingInfo_->integerVariable();
1803      if (toZero[number01]) {
1804        for (int i = 0;i<numberCutGenerators_;i++) {
1805          CglFakeClique * clique = dynamic_cast<CglFakeClique*>(generator_[i]->generator());
1806          if (clique) {
1807            OsiSolverInterface * fakeSolver = probingInfo_->analyze(*solver_,1);
1808            if (fakeSolver) {
1809              printf("Probing fake solver has %d rows\n",fakeSolver->getNumRows());
1810              //if (fakeSolver)
1811              //fakeSolver->writeMps("bad");
1812              if (generator_[i]->numberCutsInTotal())
1813                generator_[i]->setHowOften(1);
1814            }
1815            clique->assignSolver(fakeSolver);
1816            //stored->setProbingInfo(probingInfo_);
1817            break;
1818          }
1819        }
1820      }
1821      delete probingInfo_;
1822      probingInfo_=NULL;
1823    }
1824    newNode = new CbcNode ;
1825    // Set objective value (not so obvious if NLP etc)
1826    setObjectiveValue(newNode,NULL);
1827    anyAction = -1 ;
1828    // To make depth available we may need a fake node
1829    CbcNode fakeNode;
1830    if (!currentNode_) {
1831      // Not true if sub trees assert (!numberNodes_);
1832      currentNode_=&fakeNode;
1833    }
1834    phase_=3;
1835    // only allow 1000 passes
1836    int numberPassesLeft=1000;
1837    // This is first crude step
1838    if (numberAnalyzeIterations_) {
1839      delete [] analyzeResults_;
1840      analyzeResults_ = new double [4*numberIntegers_];
1841      numberFixedAtRoot_=newNode->analyze(this,analyzeResults_);
1842      if (numberFixedAtRoot_>0) {
1843        printf("%d fixed by analysis\n",numberFixedAtRoot_);
1844        setPointers(solver_);
1845        numberFixedNow_ = numberFixedAtRoot_;
1846      } else if (numberFixedAtRoot_<0) {
1847        printf("analysis found to be infeasible\n");
1848        anyAction=-2;
1849        delete newNode ;
1850        newNode = NULL ;
1851        feasible = false ;
1852      }
1853    }
1854    OsiSolverBranch * branches = NULL;
1855    anyAction = chooseBranch(newNode, numberPassesLeft, NULL, cuts,resolved,
1856                             NULL,NULL,NULL,branches);
1857    if (anyAction == -2||newNode->objectiveValue() >= cutoff) {
1858      if (anyAction != -2) {
1859        // zap parent nodeInfo
1860#ifdef COIN_DEVELOP
1861        printf("zapping CbcNodeInfo %x\n",newNode->nodeInfo()->parent());
1862#endif
1863        if (newNode->nodeInfo())
1864          newNode->nodeInfo()->nullParent();
1865      }
1866      delete newNode ;
1867      newNode = NULL ;
1868      feasible = false ;
1869    }
1870  }
1871/*
1872  At this point, the root subproblem is infeasible or fathomed by bound
1873  (newNode == NULL), or we're live with an objective value that satisfies the
1874  current objective cutoff.
1875*/
1876  assert (!newNode || newNode->objectiveValue() <= cutoff) ;
1877  // Save address of root node as we don't want to delete it
1878  // initialize for print out
1879  int lastDepth=0;
1880  int lastUnsatisfied=0;
1881  if (newNode)
1882    lastUnsatisfied=newNode->numberUnsatisfied();
1883/*
1884  The common case is that the lp relaxation is feasible but doesn't satisfy
1885  integrality (i.e., newNode->branchingObject(), indicating we've been able to
1886  select a branching variable). Remove any cuts that have gone slack due to
1887  forcing monotone variables. Then tack on an CbcFullNodeInfo object and full
1888  basis (via createInfo()) and stash the new cuts in the nodeInfo (via
1889  addCuts()). If, by some miracle, we have an integral solution at the root
1890  (newNode->branchingObject() is NULL), takeOffCuts() will ensure that the solver holds
1891  a valid solution for use by setBestSolution().
1892*/
1893  CoinWarmStartBasis *lastws = NULL ;
1894  if (feasible && newNode->branchingObject())
1895  { if (resolved)
1896    { takeOffCuts(cuts,false,NULL) ;
1897#     ifdef CHECK_CUT_COUNTS
1898      { printf("Number of rows after chooseBranch fix (root)"
1899               "(active only) %d\n",
1900                numberRowsAtContinuous_+numberNewCuts_+numberOldActiveCuts_) ;
1901        const CoinWarmStartBasis* debugws =
1902          dynamic_cast <const CoinWarmStartBasis*>(solver_->getWarmStart()) ;
1903        debugws->print() ;
1904        delete debugws ; }
1905#     endif
1906    }
1907  //newNode->createInfo(this,NULL,NULL,NULL,NULL,0,0) ;
1908    //newNode->nodeInfo()->addCuts(cuts,
1909    //                   newNode->numberBranches(),whichGenerator_) ;
1910    if (lastws) delete lastws ;
1911    lastws = dynamic_cast<CoinWarmStartBasis*>(solver_->getWarmStart()) ;
1912  }
1913/*
1914  Continuous data to be used later
1915*/
1916  continuousObjective_ = solver_->getObjValue()*solver_->getObjSense();
1917  continuousInfeasibilities_ = 0 ;
1918  if (newNode)
1919  { continuousObjective_ = newNode->objectiveValue() ;
1920    delete [] continuousSolution_;
1921    continuousSolution_ = CoinCopyOfArray(solver_->getColSolution(),
1922                                             numberColumns);
1923    continuousInfeasibilities_ = newNode->numberUnsatisfied() ; }
1924/*
1925  Bound may have changed so reset in objects
1926*/
1927  { int i ;
1928    for (i = 0;i < numberObjects_;i++)
1929      object_[i]->resetBounds(solver_) ; }
1930/*
1931  Feasible? Then we should have either a live node prepped for future
1932  expansion (indicated by variable() >= 0), or (miracle of miracles) an
1933  integral solution at the root node.
1934
1935  initializeInfo sets the reference counts in the nodeInfo object.  Since
1936  this node is still live, push it onto the heap that holds the live set.
1937*/
1938  double bestValue = 0.0 ;
1939  if (newNode) {
1940    bestValue = newNode->objectiveValue();
1941    if (newNode->branchingObject()) {
1942      newNode->initializeInfo() ;
1943      tree_->push(newNode) ;
1944      if (statistics_) {
1945        if (numberNodes2_==maximumStatistics_) {
1946          maximumStatistics_ = 2*maximumStatistics_;
1947          CbcStatistics ** temp = new CbcStatistics * [maximumStatistics_];
1948          memset(temp,0,maximumStatistics_*sizeof(CbcStatistics *));
1949          memcpy(temp,statistics_,numberNodes2_*sizeof(CbcStatistics *));
1950          delete [] statistics_;
1951          statistics_=temp;
1952        }
1953        assert (!statistics_[numberNodes2_]);
1954        statistics_[numberNodes2_]=new CbcStatistics(newNode,this);
1955      }
1956      numberNodes2_++;
1957#     ifdef CHECK_NODE
1958      printf("Node %x on tree\n",newNode) ;
1959#     endif
1960    } else {
1961      // continuous is integer
1962      double objectiveValue = newNode->objectiveValue();
1963      setBestSolution(CBC_SOLUTION,objectiveValue,
1964                      solver_->getColSolution()) ;
1965      delete newNode ;
1966      newNode = NULL ;
1967    }
1968  }
1969
1970  if (printFrequency_ <= 0) {
1971    printFrequency_ = 1000 ;
1972    if (getNumCols() > 2000)
1973      printFrequency_ = 100 ;
1974  }
1975  /*
1976    It is possible that strong branching fixes one variable and then the code goes round
1977    again and again.  This can take too long.  So we need to warn user - just once.
1978  */
1979  numberLongStrong_=0;
1980  double totalTime = 0.0;
1981  CbcNode * createdNode=NULL;
1982#ifdef CBC_THREAD
1983  CbcModel ** threadModel = NULL;
1984  Coin_pthread_t * threadId = NULL;
1985  int * threadCount = NULL;
1986  pthread_mutex_t mutex;
1987  pthread_cond_t condition_main;
1988  pthread_mutex_t condition_mutex;
1989  pthread_mutex_t * mutex2 = NULL;
1990  pthread_cond_t * condition2 = NULL;
1991  threadStruct * threadInfo = NULL;
1992#ifdef CBC_NORMAL_THREAD
1993  bool locked=false;
1994#endif
1995  int threadStats[6];
1996#ifdef CBC_DETERMINISTIC_THREAD
1997  int defaultParallelIterations=500;
1998  int defaultParallelNodes=10;
1999#endif
2000  memset(threadStats,0,sizeof(threadStats));
2001  double timeWaiting=0.0;
2002  // For now just one model
2003  if (numberThreads_) {
2004    nodeCompare_->sayThreaded(); // need to use addresses
2005    threadId = new Coin_pthread_t [numberThreads_];
2006    threadCount = new int [numberThreads_];
2007    CoinZeroN(threadCount,numberThreads_);
2008    pthread_mutex_init(&mutex,NULL);
2009    pthread_cond_init(&condition_main,NULL);
2010    pthread_mutex_init(&condition_mutex,NULL);
2011    threadModel = new CbcModel * [numberThreads_+1];
2012    threadInfo = new threadStruct [numberThreads_+1];
2013    mutex2 = new pthread_mutex_t [numberThreads_];
2014    condition2 = new pthread_cond_t [numberThreads_];
2015#ifdef CBC_DETERMINISTIC_THREAD
2016    // May need for deterministic
2017    saveObjects=new OsiObject * [numberObjects_];
2018    for (int i=0;i<numberObjects_;i++) {
2019      saveObjects[i] = object_[i]->clone();
2020    }
2021#endif
2022    // we don't want a strategy object
2023    CbcStrategy * saveStrategy = strategy_;
2024    strategy_ = NULL;
2025    for (int i=0;i<numberThreads_;i++) {
2026      pthread_mutex_init(mutex2+i,NULL);
2027      pthread_cond_init(condition2+i,NULL);
2028      threadId[i].status=0;
2029      threadInfo[i].baseModel=this;
2030      threadModel[i]=new CbcModel(*this);
2031#ifdef COIN_HAS_CLP
2032      // Solver may need to know about model
2033      CbcModel * thisModel = threadModel[i];
2034      CbcOsiSolver * solver =
2035              dynamic_cast<CbcOsiSolver *>(thisModel->solver()) ;
2036      if (solver)
2037        solver->setCbcModel(thisModel);
2038#endif
2039      mutex_ = (void *) (threadInfo+i);
2040      threadModel[i]->moveToModel(this,-1);
2041      threadInfo[i].thisModel=threadModel[i];
2042      threadInfo[i].node=NULL;
2043      threadInfo[i].createdNode=NULL;
2044      threadInfo[i].threadIdOfBase.thr=pthread_self();
2045      threadInfo[i].mutex=&mutex;
2046      threadInfo[i].mutex2=mutex2+i;
2047      threadInfo[i].condition2=condition2+i;
2048      threadInfo[i].returnCode=-1;
2049      threadInfo[i].timeLocked=0.0;
2050      threadInfo[i].timeWaitingToLock=0.0;
2051      threadInfo[i].timeWaitingToStart=0.0;
2052      threadInfo[i].timeInThread=0.0;
2053      threadInfo[i].numberTimesLocked=0;
2054      threadInfo[i].numberTimesUnlocked=0;
2055      threadInfo[i].numberTimesWaitingToStart=0;
2056      threadInfo[i].dantzigState=0; // 0 unset, -1 waiting to be set, 1 set
2057      threadInfo[i].locked=false;
2058#if CBC_THREAD_DEBUG
2059      threadInfo[i].threadNumber=i+2;
2060#endif
2061#ifdef CBC_DETERMINISTIC_THREAD
2062      threadInfo[i].delNode = NULL;
2063      threadInfo[i].maxDeleteNode=0;
2064      threadInfo[i].nDeleteNode=0;
2065      threadInfo[i].nodesThisTime=0;
2066      threadInfo[i].iterationsThisTime=0;
2067#endif
2068      pthread_create(&(threadId[i].thr),NULL,doNodesThread,threadInfo+i);
2069      threadId[i].status = 1;
2070    }
2071    strategy_ = saveStrategy;
2072    // Do a partial one for base model
2073    threadInfo[numberThreads_].baseModel=this;
2074    threadModel[numberThreads_]=this;
2075    mutex_ = (void *) (threadInfo+numberThreads_);
2076    threadInfo[numberThreads_].node=NULL;
2077    threadInfo[numberThreads_].mutex=&mutex;
2078    threadInfo[numberThreads_].condition2=&condition_main;
2079    threadInfo[numberThreads_].mutex2=&condition_mutex;
2080    threadInfo[numberThreads_].timeLocked=0.0;
2081    threadInfo[numberThreads_].timeWaitingToLock=0.0;
2082    threadInfo[numberThreads_].numberTimesLocked=0;
2083    threadInfo[numberThreads_].numberTimesUnlocked=0;
2084    threadInfo[numberThreads_].locked=false;
2085#if CBC_THREAD_DEBUG
2086    threadInfo[numberThreads_].threadNumber=1;
2087#endif
2088  }
2089#endif
2090#ifdef COIN_HAS_CLP
2091  {
2092    OsiClpSolverInterface * clpSolver
2093      = dynamic_cast<OsiClpSolverInterface *> (solver_);
2094    if (clpSolver&&!parentModel_) {
2095      clpSolver->computeLargestAway();
2096    }
2097  }
2098#endif
2099/*
2100  At last, the actual branch-and-cut search loop, which will iterate until
2101  the live set is empty or we hit some limit (integrality gap, time, node
2102  count, etc.). The overall flow is to rebuild a subproblem, reoptimise using
2103  solveWithCuts(), choose a branching pattern with chooseBranch(), and finally
2104  add the node to the live set.
2105
2106  The first action is to winnow the live set to remove nodes which are worse
2107  than the current objective cutoff.
2108*/
2109  if (solver_->getRowCutDebuggerAlways()) {
2110    OsiRowCutDebugger * debuggerX = const_cast<OsiRowCutDebugger *> (solver_->getRowCutDebuggerAlways());
2111    const OsiRowCutDebugger *debugger = solver_->getRowCutDebugger() ;
2112    if (!debugger) {
2113      // infeasible!!
2114      printf("before search\n");
2115      const double * lower = solver_->getColLower();
2116      const double * upper = solver_->getColUpper();
2117      const double * solution = debuggerX->optimalSolution();
2118      int numberColumns = solver_->getNumCols();
2119      for (int i=0;i<numberColumns;i++) {
2120        if (solver_->isInteger(i)) {
2121          if (solution[i]<lower[i]-1.0e-6||solution[i]>upper[i]+1.0e-6)
2122            printf("**** ");
2123          printf("%d %g <= %g <= %g\n",
2124                 i,lower[i],solution[i],upper[i]);
2125        }
2126      }
2127      //abort();
2128    }
2129  }
2130#ifdef CBC_DETERMINISTIC_THREAD
2131#define MAX_DEL_NODE 1
2132  CbcNode * delNode[MAX_DEL_NODE+1];
2133  int nDeleteNode=0;
2134  bool goneParallel=false;
2135#endif
2136  // For Printing etc when parallel
2137  int lastEvery1000=0;
2138  int lastPrintEvery=0;
2139  while (true) {
2140#ifdef CBC_NORMAL_THREAD
2141    if (!locked) {
2142      lockThread();
2143      locked=true;
2144    }
2145#endif
2146#ifdef COIN_HAS_CLP
2147    // Possible change of pivot method
2148    if(!savePivotMethod&&!parentModel_) {
2149      OsiClpSolverInterface * clpSolver
2150        = dynamic_cast<OsiClpSolverInterface *> (solver_);
2151      if (clpSolver&&numberNodes_>=100&&numberNodes_<200) {
2152        if (numberIterations_<numberNodes_*20) {
2153          ClpSimplex * simplex = clpSolver->getModelPtr();
2154          ClpDualRowPivot * pivotMethod=simplex->dualRowPivot();
2155          ClpDualRowDantzig * pivot =
2156            dynamic_cast< ClpDualRowDantzig*>(pivotMethod);
2157          if (!pivot) {
2158            savePivotMethod = pivotMethod->clone(true);
2159            ClpDualRowDantzig dantzig;
2160            simplex->setDualRowPivotAlgorithm(dantzig);
2161#ifdef COIN_DEVELOP
2162            printf("%d node, %d iterations ->Dantzig\n",numberNodes_,
2163                   numberIterations_);
2164#endif
2165#ifdef CBC_THREAD
2166            for (int i=0;i<numberThreads_;i++) {
2167              threadInfo[i].dantzigState=-1;
2168            }
2169#endif
2170          }
2171        }
2172      }
2173    }
2174#endif
2175    if (tree_->empty()) {
2176#ifdef CBC_NORMAL_THREAD
2177      if (numberThreads_) {
2178#ifdef COIN_DEVELOP
2179        printf("empty\n");
2180#endif
2181        // may still be outstanding nodes
2182        int iThread;
2183        for (iThread=0;iThread<numberThreads_;iThread++) {
2184          if (threadId[iThread].status) {
2185            if (threadInfo[iThread].returnCode==0) 
2186              break;
2187          }
2188        }
2189        if (iThread<numberThreads_) {
2190#ifdef COIN_DEVELOP
2191          printf("waiting for thread %d code 0\n",iThread);
2192#endif
2193#ifndef CBC_DETERMINISTIC_THREAD
2194          unlockThread();
2195#endif
2196          locked = false;
2197          pthread_cond_signal(threadInfo[iThread].condition2); // unlock in case
2198          while (true) {
2199            pthread_mutex_lock(&condition_mutex);
2200            struct timespec absTime;
2201            clock_gettime(CLOCK_REALTIME,&absTime);
2202            double time = absTime.tv_sec+1.0e-9*absTime.tv_nsec;
2203            absTime.tv_nsec += 1000000; // millisecond
2204            if (absTime.tv_nsec>=1000000000) {
2205              absTime.tv_nsec -= 1000000000;
2206              absTime.tv_sec++;
2207            }
2208            pthread_cond_timedwait(&condition_main,&condition_mutex,&absTime);
2209            clock_gettime(CLOCK_REALTIME,&absTime);
2210            double time2 = absTime.tv_sec+1.0e-9*absTime.tv_nsec;
2211            timeWaiting += time2-time;
2212            pthread_mutex_unlock(&condition_mutex);
2213            if (threadInfo[iThread].returnCode!=0) 
2214              break;
2215            pthread_cond_signal(threadInfo[iThread].condition2); // unlock
2216          }
2217          threadModel[iThread]->moveToModel(this,1);
2218          assert (threadInfo[iThread].returnCode==1);
2219          if (threadInfo[iThread].dantzigState==-1) {
2220            // 0 unset, -1 waiting to be set, 1 set
2221            threadInfo[iThread].dantzigState=1;
2222            CbcModel * model = threadInfo[iThread].thisModel;
2223            OsiClpSolverInterface * clpSolver2
2224              = dynamic_cast<OsiClpSolverInterface *> (model->solver());
2225            assert (clpSolver2);
2226            ClpSimplex * simplex2 = clpSolver2->getModelPtr();
2227            ClpDualRowDantzig dantzig;
2228            simplex2->setDualRowPivotAlgorithm(dantzig);
2229          }
2230          // say available
2231          threadInfo[iThread].returnCode=-1;
2232          threadStats[4]++;
2233#ifdef COIN_DEVELOP
2234          printf("thread %d code now -1\n",iThread);
2235#endif
2236          continue;
2237        } else {
2238#ifdef COIN_DEVELOP
2239          printf("no threads at code 0 \n");
2240#endif
2241          // now check if any have just finished
2242          for (iThread=0;iThread<numberThreads_;iThread++) {
2243            if (threadId[iThread].status) {
2244              if (threadInfo[iThread].returnCode==1) 
2245                break;
2246            }
2247          }
2248          if (iThread<numberThreads_) {
2249#ifndef CBC_DETERMINISTIC_THREAD
2250            unlockThread();
2251#endif
2252            locked = false;
2253            threadModel[iThread]->moveToModel(this,1);
2254            assert (threadInfo[iThread].returnCode==1);
2255            // say available
2256            threadInfo[iThread].returnCode=-1;
2257            threadStats[4]++;
2258#ifdef COIN_DEVELOP
2259            printf("thread %d code now -1\n",iThread);
2260#endif
2261            continue;
2262          }
2263        }
2264        if (!tree_->empty()) {
2265#ifdef COIN_DEVELOP
2266          printf("tree not empty!!!!!!\n");
2267#endif
2268          continue;
2269        }
2270        for (iThread=0;iThread<numberThreads_;iThread++) {
2271          if (threadId[iThread].status) {
2272            if (threadInfo[iThread].returnCode!=-1) { 
2273              printf("bad end of tree\n");
2274              abort();
2275            }
2276          }
2277        }
2278#ifdef COIN_DEVELOP
2279        printf("finished ************\n");
2280#endif
2281      }
2282#ifndef CBC_DETERMINISTIC_THREAD
2283      unlockThread();
2284#endif
2285      locked=false; // not needed as break
2286#endif
2287      break;
2288    }
2289#ifdef CBC_NORMAL_THREAD
2290    unlockThread();
2291    locked = false;
2292#endif
2293    // If done 100 nodes see if worth trying reduction
2294    if (numberNodes_==100&&saveSolver) {
2295      bool tryNewSearch=solverCharacteristics_->reducedCostsAccurate();
2296      int numberColumns = getNumCols();
2297      if (tryNewSearch) {
2298        double cutoff = getCutoff() ;
2299        saveSolver->resolve();
2300        double direction = saveSolver->getObjSense() ;
2301        double gap = cutoff - saveSolver->getObjValue()*direction ;
2302        double tolerance;
2303        saveSolver->getDblParam(OsiDualTolerance,tolerance) ;
2304        if (gap<=0.0)
2305          gap = tolerance; 
2306        gap += 100.0*tolerance;
2307        double integerTolerance = getDblParam(CbcIntegerTolerance) ;
2308       
2309        const double *lower = saveSolver->getColLower() ;
2310        const double *upper = saveSolver->getColUpper() ;
2311        const double *solution = saveSolver->getColSolution() ;
2312        const double *reducedCost = saveSolver->getReducedCost() ;
2313       
2314        int numberFixed = 0 ;
2315        int numberFixed2=0;
2316        printf("gap %g\n",gap);
2317        for (int i = 0 ; i < numberIntegers_ ; i++) {
2318          int iColumn = integerVariable_[i] ;
2319          double djValue = direction*reducedCost[iColumn] ;
2320          if (upper[iColumn]-lower[iColumn] > integerTolerance) {
2321            if (solution[iColumn] < lower[iColumn]+integerTolerance && djValue > gap) {
2322              //printf("%d to lb on dj of %g - bounds %g %g\n",
2323              //     iColumn,djValue,lower[iColumn],upper[iColumn]);
2324              saveSolver->setColUpper(iColumn,lower[iColumn]) ;
2325              numberFixed++ ;
2326            } else if (solution[iColumn] > upper[iColumn]-integerTolerance && -djValue > gap) {
2327              //printf("%d to ub on dj of %g - bounds %g %g\n",
2328              //     iColumn,djValue,lower[iColumn],upper[iColumn]);
2329              saveSolver->setColLower(iColumn,upper[iColumn]) ;
2330              numberFixed++ ;
2331            }
2332          } else {
2333            //printf("%d has dj of %g - already fixed to %g\n",
2334            //     iColumn,djValue,lower[iColumn]);
2335            numberFixed2++;
2336          }
2337        }
2338#ifdef COIN_DEVELOP
2339        if ((specialOptions_&1)!=0) {
2340          const OsiRowCutDebugger *debugger = saveSolver->getRowCutDebugger() ;
2341          if (debugger) { 
2342            printf("Contains optimal\n") ;
2343            saveSolver->writeMps("reduced");
2344          } else {
2345            abort();
2346          }
2347        }
2348        printf("Restart could fix %d integers (%d already fixed)\n",
2349               numberFixed+numberFixed2,numberFixed2);
2350#endif
2351        numberFixed += numberFixed2;
2352        if (numberFixed*10<numberColumns)
2353          tryNewSearch=false;
2354      }
2355      if (tryNewSearch) {
2356        // back to solver without cuts?
2357#if 1
2358        OsiSolverInterface * solver2 = continuousSolver_->clone();
2359#else
2360        OsiSolverInterface * solver2 = saveSolver->clone();
2361#endif
2362        const double *lower = saveSolver->getColLower() ;
2363        const double *upper = saveSolver->getColUpper() ;
2364        for (int i = 0 ; i < numberIntegers_ ; i++) {
2365          int iColumn = integerVariable_[i] ;
2366          solver2->setColLower(iColumn,lower[iColumn]);
2367          solver2->setColUpper(iColumn,upper[iColumn]);
2368        }
2369        // swap
2370        delete saveSolver;
2371        saveSolver=solver2;
2372        double * newSolution = new double[numberColumns];
2373        double objectiveValue=cutoff;
2374        CbcSerendipity heuristic(*this);
2375        if (bestSolution_)
2376          heuristic.setInputSolution(bestSolution_,bestObjective_);
2377        heuristic.setFractionSmall(0.6);
2378        heuristic.setFeasibilityPumpOptions(1008013);
2379        int returnCode= heuristic.smallBranchAndBound(saveSolver,
2380                                                      -1,newSolution,
2381                                                      objectiveValue,
2382                                                      cutoff,"Reduce");
2383        if (returnCode==-1) {
2384#ifdef COIN_DEVELOP
2385          printf("Restart - not small enough to do search after fixing\n");
2386#endif
2387          delete [] newSolution;
2388        } else {
2389          assert (returnCode>=0);
2390          if ((returnCode&1)!=0) {
2391            // increment number of solutions so other heuristics can test
2392            numberSolutions_++;
2393            numberHeuristicSolutions_++;
2394            lastHeuristic_ = NULL;
2395            setBestSolution(CBC_ROUNDING,objectiveValue,newSolution) ;
2396          }
2397          delete [] newSolution;
2398          break;
2399        }
2400      } 
2401      delete saveSolver;
2402      saveSolver=NULL;
2403    }
2404/*
2405  Check for abort on limits: node count, solution count, time, integrality gap.
2406*/
2407    totalTime = getCurrentSeconds() ;
2408    double maxSeconds = getMaximumSeconds();
2409    if (parentModel_)
2410      maxSeconds=CoinMin(maxSeconds,parentModel_->getMaximumSeconds());
2411    if (!(numberNodes_ < intParam_[CbcMaxNumNode] &&
2412          numberSolutions_ < intParam_[CbcMaxNumSol] &&
2413          totalTime < maxSeconds &&
2414          !stoppedOnGap_&&!eventHappened_&&(maximumNumberIterations_<0||
2415                                            numberIterations_<maximumNumberIterations_))) {
2416      // out of loop
2417      break;
2418    }
2419#ifdef BONMIN
2420    assert(!solverCharacteristics_->solutionAddsCuts() || solverCharacteristics_->mipFeasible());
2421#endif
2422    if (cutoff > getCutoff()) {
2423      double newCutoff = getCutoff();
2424      if (analyzeResults_) {
2425        // see if we could fix any (more)
2426        int n=0;
2427        double * newLower = analyzeResults_;
2428        double * objLower = newLower+numberIntegers_;
2429        double * newUpper = objLower+numberIntegers_;
2430        double * objUpper = newUpper+numberIntegers_;
2431        for (int i=0;i<numberIntegers_;i++) {
2432          if (objLower[i]>newCutoff) {
2433            n++;
2434            if (objUpper[i]>newCutoff) {
2435              newCutoff = -COIN_DBL_MAX;
2436              break;
2437            }
2438          } else if (objUpper[i]>newCutoff) {
2439            n++;
2440          }
2441        }
2442        if (newCutoff==-COIN_DBL_MAX) {
2443          printf("Root analysis says finished\n");
2444        } else if (n>numberFixedNow_) {
2445          printf("%d more fixed by analysis - now %d\n",n-numberFixedNow_,n);
2446          numberFixedNow_=n;
2447        }
2448      }
2449      if (eventHandler) {
2450        if (!eventHandler->event(CbcEventHandler::solution)) {
2451          eventHappened_=true; // exit
2452        }
2453      }
2454#ifndef CBC_DETERMINISTIC_THREAD
2455      lockThread();
2456#endif
2457      // Do from deepest
2458      tree_->cleanTree(this, newCutoff,bestPossibleObjective_) ;
2459      nodeCompare_->newSolution(this) ;
2460      nodeCompare_->newSolution(this,continuousObjective_,
2461                                continuousInfeasibilities_) ;
2462      tree_->setComparison(*nodeCompare_) ;
2463      if (tree_->empty()) {
2464#ifndef CBC_DETERMINISTIC_THREAD
2465        unlockThread();
2466#endif
2467        // For threads we need to check further
2468        //break; // finished
2469        continue;
2470      }
2471#ifndef CBC_DETERMINISTIC_THREAD
2472      unlockThread();
2473#endif
2474    }
2475    cutoff = getCutoff() ;
2476/*
2477    Periodic activities: Opportunities to
2478    + tweak the nodeCompare criteria,
2479    + check if we've closed the integrality gap enough to quit,
2480    + print a summary line to let the user know we're working
2481*/
2482    if (numberNodes_>=lastEvery1000) {
2483#ifndef CBC_DETERMINISTIC_THREAD
2484      lockThread();
2485#endif
2486#ifdef COIN_HAS_CLP
2487      // Possible change of pivot method
2488      if(!savePivotMethod&&!parentModel_) {
2489        OsiClpSolverInterface * clpSolver
2490          = dynamic_cast<OsiClpSolverInterface *> (solver_);
2491        if (clpSolver&&numberNodes_>=1000&&numberNodes_<2000) {
2492          if (numberIterations_<numberNodes_*20) {
2493            ClpSimplex * simplex = clpSolver->getModelPtr();
2494            ClpDualRowPivot * pivotMethod=simplex->dualRowPivot();
2495            ClpDualRowDantzig * pivot =
2496              dynamic_cast< ClpDualRowDantzig*>(pivotMethod);
2497            if (!pivot) {
2498              savePivotMethod = pivotMethod->clone(true);
2499              ClpDualRowDantzig dantzig;
2500              simplex->setDualRowPivotAlgorithm(dantzig);
2501#ifdef COIN_DEVELOP
2502              printf("%d node, %d iterations ->Dantzig\n",numberNodes_,
2503                     numberIterations_);
2504#endif
2505#ifdef CBC_THREAD
2506              for (int i=0;i<numberThreads_;i++) {
2507                threadInfo[i].dantzigState=-1;
2508              }
2509#endif
2510            }
2511          }
2512        }
2513      }
2514#endif
2515      lastEvery1000 = numberNodes_ + 1000;
2516      bool redoTree=nodeCompare_->every1000Nodes(this, numberNodes_) ;
2517#ifdef CHECK_CUT_SIZE
2518      verifyCutSize (tree_, *this);
2519#endif
2520      // redo tree if wanted
2521      if (redoTree)
2522        tree_->setComparison(*nodeCompare_) ;
2523#ifndef CBC_DETERMINISTIC_THREAD
2524      unlockThread();
2525#endif
2526    }
2527    if (saveCompare&&!hotstartSolution_) {
2528      // hotstart switched off
2529      delete nodeCompare_; // off depth first
2530      nodeCompare_=saveCompare;
2531      saveCompare=NULL;
2532      // redo tree
2533#ifndef CBC_DETERMINISTIC_THREAD
2534      lockThread();
2535#endif
2536      tree_->setComparison(*nodeCompare_) ;
2537#ifndef CBC_DETERMINISTIC_THREAD
2538      unlockThread();
2539#endif
2540    }
2541    if (numberNodes_>=lastPrintEvery) {
2542      lastPrintEvery = numberNodes_ + printFrequency_;
2543#ifdef CBC_INSTRUMENT
2544      if (0) {
2545        printf("==Start instrument\n");
2546        for (int iObject=0;iObject<numberObjects_;iObject++) {
2547          CbcSimpleIntegerDynamicPseudoCost * obj =
2548            dynamic_cast <CbcSimpleIntegerDynamicPseudoCost *>(object_[iObject]) ;
2549          if (obj)
2550            obj->print();
2551        }
2552        printf("==End instrument\n");
2553      }
2554#endif
2555#ifndef CBC_DETERMINISTIC_THREAD
2556      lockThread();
2557#endif
2558      int nNodes = tree_->size() ;
2559
2560      //MODIF PIERRE
2561      bestPossibleObjective_ = tree_->getBestPossibleObjective();
2562#ifndef CBC_DETERMINISTIC_THREAD
2563      unlockThread();
2564#endif
2565      if (!intParam_[CbcPrinting]) {
2566        messageHandler()->message(CBC_STATUS,messages())
2567          << numberNodes_<< nNodes<< bestObjective_<< bestPossibleObjective_
2568          <<getCurrentSeconds()
2569          << CoinMessageEol ;
2570      } else if (intParam_[CbcPrinting]==1) {
2571        messageHandler()->message(CBC_STATUS2,messages())
2572          << numberNodes_<< nNodes<< bestObjective_<< bestPossibleObjective_
2573          <<lastDepth<<lastUnsatisfied<<numberIterations_
2574          <<getCurrentSeconds()
2575          << CoinMessageEol ;
2576      } else if (!numberExtraIterations_) {
2577        messageHandler()->message(CBC_STATUS2,messages())
2578          << numberNodes_<< nNodes<< bestObjective_<< bestPossibleObjective_
2579          <<lastDepth<<lastUnsatisfied<<numberIterations_
2580          <<getCurrentSeconds()
2581          << CoinMessageEol ;
2582      } else {
2583        messageHandler()->message(CBC_STATUS3,messages())
2584          << numberNodes_<<numberExtraNodes_<< nNodes<< bestObjective_<< bestPossibleObjective_
2585          <<lastDepth<<lastUnsatisfied<<numberIterations_<<numberExtraIterations_
2586          <<getCurrentSeconds()
2587          << CoinMessageEol ;
2588      }
2589      if (!eventHandler->event(CbcEventHandler::treeStatus)) {
2590        eventHappened_=true; // exit
2591      }
2592    }
2593    // See if can stop on gap
2594    double testGap = CoinMax(dblParam_[CbcAllowableGap],
2595                             CoinMax(fabs(bestObjective_),fabs(bestPossibleObjective_))
2596                             *dblParam_[CbcAllowableFractionGap]);
2597    if (bestObjective_-bestPossibleObjective_ < testGap && getCutoffIncrement()>=0.0) {
2598      stoppedOnGap_ = true ;
2599    }
2600   
2601#ifdef CHECK_NODE_FULL
2602    verifyTreeNodes(tree_,*this) ;
2603#   endif
2604#   ifdef CHECK_CUT_COUNTS
2605    verifyCutCounts(tree_,*this) ;
2606#   endif
2607/*
2608  Now we come to the meat of the loop. To create the active subproblem, we'll
2609  pop the most promising node in the live set, rebuild the subproblem it
2610  represents, and then execute the current arm of the branch to create the
2611  active subproblem.
2612*/
2613#ifdef BACK_TO_OLD_WAY //old way without threads #ifndef CBC_THREAD
2614    CbcNode *node = tree_->bestNode(cutoff) ;
2615    // Possible one on tree worse than cutoff
2616    if (!node||node->objectiveValue()>cutoff)
2617      continue;
2618    int currentNumberCuts = 0 ;
2619    currentNode_=node; // so can be accessed elsewhere
2620#ifdef CBC_DEBUG
2621    printf("%d unsat, way %d, obj %g est %g\n",
2622           node->numberUnsatisfied(),node->way(),node->objectiveValue(),
2623           node->guessedObjectiveValue());
2624#endif
2625#if NEW_UPDATE_OBJECT==0
2626    // Save clone in branching decision
2627    if(branchingMethod_)
2628      branchingMethod_->saveBranchingObject(node->modifiableBranchingObject());
2629#endif
2630    // Say not on optimal path
2631    bool onOptimalPath=false;
2632#   ifdef CHECK_NODE
2633    printf("Node %x popped from tree - %d left, %d count\n",node,
2634           node->nodeInfo()->numberBranchesLeft(),
2635           node->nodeInfo()->numberPointingToThis()) ;
2636    printf("\tdepth = %d, z =  %g, unsat = %d, var = %d.\n",
2637           node->depth(),node->objectiveValue(),
2638           node->numberUnsatisfied(),
2639           node->columnNumber()) ;
2640#   endif
2641    lastDepth=node->depth();
2642    lastUnsatisfied=node->numberUnsatisfied();
2643
2644/*
2645  Rebuild the subproblem for this node:  Call addCuts() to adjust the model
2646  to recreate the subproblem for this node (set proper variable bounds, add
2647  cuts, create a basis).  This may result in the problem being fathomed by
2648  bound or infeasibility. Returns 1 if node is fathomed.
2649  Execute the current arm of the branch: If the problem survives, save the
2650  resulting variable bounds and call branch() to modify variable bounds
2651  according to the current arm of the branching object. If we're processing
2652  the final arm of the branching object, flag the node for removal from the
2653  live set.
2654*/
2655    CbcNodeInfo * nodeInfo = node->nodeInfo() ;
2656    newNode = NULL ;
2657    int branchesLeft=0;
2658    if (!addCuts(node,lastws,numberFixedNow_>numberFixedAtRoot_))
2659    { int i ;
2660      const double * lower = getColLower() ;
2661      const double * upper = getColUpper() ;
2662      for (i = 0 ; i < numberColumns ; i++)
2663      { lowerBefore[i]= lower[i] ;
2664        upperBefore[i]= upper[i] ; }
2665      if ((solverCharacteristics_->extraCharacteristics()&2)!=0) {
2666        solverCharacteristics_->setBeforeLower(lowerBefore);
2667        solverCharacteristics_->setBeforeUpper(upperBefore);
2668      }
2669      if (messageHandler()->logLevel()>2)
2670        node->modifiableBranchingObject()->print();
2671      if (!useOsiBranching) 
2672        branchesLeft = node->branch(NULL); // old way
2673      else
2674        branchesLeft = node->branch(solver_); // new way
2675      if (branchesLeft) {
2676        // set nodenumber correctly
2677        node->nodeInfo()->setNodeNumber(numberNodes2_);
2678        tree_->push(node) ;
2679        if (statistics_) {
2680          if (numberNodes2_==maximumStatistics_) {
2681            maximumStatistics_ = 2*maximumStatistics_;
2682            CbcStatistics ** temp = new CbcStatistics * [maximumStatistics_];
2683            memset(temp,0,maximumStatistics_*sizeof(CbcStatistics *));
2684            memcpy(temp,statistics_,numberNodes2_*sizeof(CbcStatistics *));
2685            delete [] statistics_;
2686            statistics_=temp;
2687          }
2688          assert (!statistics_[numberNodes2_]);
2689          statistics_[numberNodes2_]=new CbcStatistics(node,this);
2690        }
2691        numberNodes2_++;
2692        //nodeOnTree=true; // back on tree
2693        //deleteNode = false ;
2694#       ifdef CHECK_NODE
2695        printf("Node %x pushed back on tree - %d left, %d count\n",node,
2696               nodeInfo->numberBranchesLeft(),
2697               nodeInfo->numberPointingToThis()) ;
2698#       endif
2699      } else {
2700        //deleteNode = true ;
2701        if (!nodeInfo->numberBranchesLeft())
2702          nodeInfo->allBranchesGone(); // can clean up
2703      }
2704      if ((specialOptions_&1)!=0) {
2705        /*
2706          This doesn't work as intended --- getRowCutDebugger will return null
2707          unless the current feasible solution region includes the optimal solution
2708          that RowCutDebugger knows. There's no way to tell inactive from off the
2709          optimal path.
2710        */
2711        const OsiRowCutDebugger *debugger = solver_->getRowCutDebugger() ;
2712        if (debugger) {
2713          onOptimalPath=true;
2714          printf("On optimal path\n") ;
2715        }
2716      }
2717     
2718/*
2719  Reoptimize, possibly generating cuts and/or using heuristics to find
2720  solutions.  Cut reference counts are unaffected unless we lose feasibility,
2721  in which case solveWithCuts() will make the adjustment.
2722*/
2723      phase_=2;
2724      cuts = OsiCuts() ;
2725      currentNumberCuts = solver_->getNumRows()-numberRowsAtContinuous_ ;
2726      int saveNumber = numberIterations_;
2727      if(solverCharacteristics_->solutionAddsCuts()) {
2728        int returnCode=resolve(node ? node->nodeInfo() : NULL,1);
2729        feasible = returnCode != 0;
2730        if (feasible) {
2731          int iObject ;
2732          int preferredWay ;
2733          int numberUnsatisfied = 0 ;
2734          memcpy(currentSolution_,solver_->getColSolution(),
2735                 numberColumns*sizeof(double)) ;
2736          // point to useful information
2737          OsiBranchingInformation usefulInfo=usefulInformation();
2738         
2739          for (iObject = 0 ; iObject < numberObjects_ ; iObject++) {
2740            double infeasibility =
2741              object_[iObject]->infeasibility(&usefulInfo,preferredWay) ;
2742            if (infeasibility ) numberUnsatisfied++ ;
2743          }
2744          if (returnCode>0) {
2745            if (numberUnsatisfied)   {
2746              feasible = solveWithCuts(cuts,maximumCutPasses_,node);
2747            } else {
2748              // may generate cuts and turn the solution
2749              //to an infeasible one
2750              feasible = solveWithCuts(cuts, 1,
2751                                       node);
2752#if 0
2753              currentNumberCuts_ = cuts.sizeRowCuts();
2754              if (currentNumberCuts_ >= maximumNumberCuts_) {
2755                maximumNumberCuts_ = currentNumberCuts;
2756                delete [] addedCuts_;
2757                addedCuts_ = new CbcCountRowCut * [maximumNumberCuts_];
2758              }
2759#endif
2760            }
2761          }
2762          // check extra info on feasibility
2763          if (!solverCharacteristics_->mipFeasible()) {
2764            feasible = false;
2765            solverCharacteristics_->setMipBound(-COIN_DBL_MAX);
2766          }
2767        }
2768      } else {
2769        // normal
2770        //int zzzzzz=0;
2771        //if (zzzzzz)
2772        //solver_->writeMps("before");
2773        feasible = solveWithCuts(cuts,maximumCutPasses_,node);
2774      }
2775      if ((specialOptions_&1)!=0&&onOptimalPath) {
2776        if (!solver_->getRowCutDebugger()) {
2777          if (solver_->getRowCutDebuggerAlways()->optimalValue()<
2778              getCutoff()-1.0e-5) {
2779            // dj fix did something???
2780            solver_->writeMpsNative("infeas2.mps",NULL,NULL,2);
2781            solver_->getRowCutDebuggerAlways()->printOptimalSolution(*solver_);
2782            assert (solver_->getRowCutDebugger()) ;
2783          }
2784        }
2785      }
2786      if (statistics_) {
2787        assert (numberNodes2_);
2788        assert (statistics_[numberNodes2_-1]);
2789        assert (statistics_[numberNodes2_-1]->node()==numberNodes2_-1);
2790        statistics_[numberNodes2_-1]->endOfBranch(numberIterations_-saveNumber,
2791                                               feasible ? solver_->getObjValue()
2792                                               : COIN_DBL_MAX);
2793      }
2794/*
2795  Are we still feasible? If so, create a node and do the work to attach a
2796  branching object, reoptimising as needed if chooseBranch() identifies
2797  monotone objects.
2798
2799  Finally, attach a partial nodeInfo object and store away any cuts that we
2800  created back in solveWithCuts. addCuts() will initialise the reference
2801  counts for these new cuts.
2802
2803  This next test can be problematic if we've discovered an
2804  alternate equivalent answer and subsequently fathom the solution
2805  known to the row cut debugger due to bounds.
2806*/
2807        if (onOptimalPath) {
2808          bool objLim = solver_->isDualObjectiveLimitReached() ;
2809          if (!feasible && !objLim) {
2810            printf("infeas2\n");
2811            solver_->writeMpsNative("infeas.mps",NULL,NULL,2);
2812            solver_->getRowCutDebuggerAlways()->printOptimalSolution(*solver_);
2813            CoinWarmStartBasis *slack =
2814              dynamic_cast<CoinWarmStartBasis *>(solver_->getEmptyWarmStart()) ;
2815            solver_->setWarmStart(slack);
2816            delete slack ;
2817            solver_->setHintParam(OsiDoReducePrint,false,OsiHintDo,0) ;
2818            solver_->initialSolve();
2819            assert (!solver_->isProvenOptimal());
2820          }
2821          assert (feasible || objLim);
2822        }
2823        bool checkingNode=false;
2824        if (feasible) {
2825          newNode = new CbcNode ;//Regular node of the tree
2826          // Set objective value (not so obvious if NLP etc)
2827          setObjectiveValue(newNode,node);
2828          anyAction =-1 ;
2829          resolved = false ;
2830          if (newNode->objectiveValue() >= getCutoff()) 
2831            anyAction=-2;
2832          // only allow at most a few passes
2833          int numberPassesLeft=5;
2834          checkingNode=true;
2835        OsiSolverBranch * branches=NULL;
2836        // point to useful information
2837        anyAction = chooseBranch(newNode, numberPassesLeft,node, cuts,resolved,
2838                                 lastws, lowerBefore, upperBefore, branches);
2839/*
2840  If we end up infeasible, we can delete the new node immediately. Since this
2841  node won't be needing the cuts we collected, decrement the reference counts.
2842  If we are feasible, then we'll be placing this node into the live set, so
2843  increment the reference count in the current (parent) nodeInfo.
2844*/
2845        if (anyAction == -2)
2846          { delete newNode ;
2847          newNode = NULL ;
2848          // say strong doing well
2849          if (checkingNode)
2850            setSpecialOptions(specialOptions_|8);
2851          for (i = 0 ; i < currentNumberCuts_ ; i++)
2852            { if (addedCuts_[i])
2853              { if (!addedCuts_[i]->decrement(1))
2854                delete addedCuts_[i] ; } } }
2855        else
2856          { nodeInfo->increment() ;
2857          if ((numberNodes_%20)==0) {
2858            // say strong not doing as well
2859            setSpecialOptions(specialOptions_&~8);
2860          }
2861        }
2862        }
2863/*
2864  At this point, there are three possibilities:
2865    * newNode is live and will require further branching to resolve
2866      (variable() >= 0). Increment the cut reference counts by
2867      numberBranches() to allow for use by children of this node, and
2868      decrement by 1 because we've executed one arm of the branch of our
2869      parent (consuming one reference). Before we push newNode onto the
2870      search tree, try for a heuristic solution.
2871    * We have a solution, in which case newNode is non-null but we have no
2872      branching variable. Decrement the cut counts and save the solution.
2873    * The node was found to be infeasible, in which case it's already been
2874      deleted, and newNode is null.
2875*/
2876        if (!eventHandler->event(CbcEventHandler::node)) {
2877          eventHappened_=true; // exit
2878        }
2879        assert (!newNode || newNode->objectiveValue() <= getCutoff()) ;
2880        if (statistics_) {
2881          assert (numberNodes2_);
2882          assert (statistics_[numberNodes2_-1]);
2883          assert (statistics_[numberNodes2_-1]->node()==numberNodes2_-1);
2884          if (newNode)
2885            statistics_[numberNodes2_-1]->updateInfeasibility(newNode->numberUnsatisfied());
2886          else
2887            statistics_[numberNodes2_-1]->sayInfeasible();
2888        }
2889        if (newNode) {
2890          if (newNode->branchingObject() == NULL&&solverCharacteristics_->solverType()==4) {
2891            // need to check if any cuts would do anything
2892            OsiCuts theseCuts;
2893            // reset probing info
2894            //if (probingInfo_)
2895            //probingInfo_->initializeFixing();
2896            for (int i = 0;i<numberCutGenerators_;i++) {
2897              bool generate = generator_[i]->normal();
2898              // skip if not optimal and should be (maybe a cut generator has fixed variables)
2899              if (generator_[i]->needsOptimalBasis()&&!solver_->basisIsAvailable())
2900                generate=false;
2901              if (!generator_[i]->mustCallAgain())
2902                generate=false; // only special cuts
2903              if (generate) {
2904                generator_[i]->generateCuts(theseCuts,1,solver_,NULL) ;
2905                int numberRowCutsAfter = theseCuts.sizeRowCuts() ;
2906                if (numberRowCutsAfter) {
2907                  // need dummy branch
2908                  newNode->setBranchingObject(new CbcDummyBranchingObject(this));
2909                  newNode->nodeInfo()->initializeInfo(1);
2910                  break;
2911                }
2912              }
2913            }
2914          }
2915          if (newNode->branchingObject())
2916          { handler_->message(CBC_BRANCH,messages_)
2917               << numberNodes_<< newNode->objectiveValue()
2918               << newNode->numberUnsatisfied()<< newNode->depth()
2919               << CoinMessageEol ;
2920            // Increment cut counts (taking off current)
2921            int numberLeft = newNode->numberBranches() ;
2922            for (i = 0;i < currentNumberCuts_;i++)
2923            { if (addedCuts_[i])
2924              {
2925#               ifdef CHECK_CUT_COUNTS
2926                printf("Count on cut %x increased by %d\n",addedCuts_[i],
2927                        numberLeft-1) ;
2928#               endif
2929                addedCuts_[i]->increment(numberLeft-1) ; } }
2930
2931            double estValue = newNode->guessedObjectiveValue() ;
2932            int found = -1 ;
2933            // no - overhead on small problems solver_->resolve() ;     // double check current optimal
2934            // assert (!solver_->getIterationCount());
2935            double * newSolution = new double [numberColumns] ;
2936            double heurValue = getCutoff() ;
2937            int iHeur ;
2938            for (iHeur = 0 ; iHeur < numberHeuristics_ ; iHeur++) {
2939              double saveValue = heurValue ;
2940              int ifSol = heuristic_[iHeur]->solution(heurValue,newSolution) ;
2941              if (ifSol > 0) {
2942                // new solution found
2943                heuristic_[iHeur]->incrementNumberSolutionsFound();
2944                found = iHeur ;
2945                incrementUsed(newSolution);
2946                lastHeuristic_ = heuristic_[found];
2947                setBestSolution(CBC_ROUNDING,heurValue,newSolution) ;
2948              } else if (ifSol < 0) {
2949                // just returning an estimate
2950                estValue = CoinMin(heurValue,estValue) ;
2951                heurValue = saveValue ;
2952              }
2953            }
2954            delete [] newSolution ;
2955            newNode->setGuessedObjectiveValue(estValue) ;
2956            tree_->push(newNode) ;
2957            if (statistics_) {
2958              if (numberNodes2_==maximumStatistics_) {
2959                maximumStatistics_ = 2*maximumStatistics_;
2960                CbcStatistics ** temp = new CbcStatistics * [maximumStatistics_];
2961                memset(temp,0,maximumStatistics_*sizeof(CbcStatistics *));
2962                memcpy(temp,statistics_,numberNodes2_*sizeof(CbcStatistics *));
2963                delete [] statistics_;
2964                statistics_=temp;
2965              }
2966              assert (!statistics_[numberNodes2_]);
2967              statistics_[numberNodes2_]=new CbcStatistics(newNode,this);
2968            }
2969            numberNodes2_++;
2970#           ifdef CHECK_NODE
2971            printf("Node %x pushed on tree c\n",newNode) ;
2972#           endif
2973          }
2974          else
2975          { 
2976            if(solverCharacteristics_ && //we may be in a non standard bab
2977               solverCharacteristics_->solutionAddsCuts()// we are in some kind of OA based bab.
2978               )
2979              {
2980                std::cerr<<"You should never get here"<<std::endl;
2981                throw CoinError("Nodes should not be fathomed on integer infeasibility in this setting",
2982                                "branchAndBound","CbcModel") ;
2983              }
2984            for (i = 0 ; i < currentNumberCuts_ ; i++)
2985            { if (addedCuts_[i])
2986              { if (!addedCuts_[i]->decrement(1))
2987                  delete addedCuts_[i] ; } }
2988          double objectiveValue = newNode->objectiveValue();
2989            setBestSolution(CBC_SOLUTION,objectiveValue,
2990                            solver_->getColSolution()) ;
2991            lastHeuristic_ = NULL;
2992            incrementUsed(solver_->getColSolution());
2993            //assert(nodeInfo->numberPointingToThis() <= 2) ;
2994            // avoid accidental pruning, if newNode was final branch arm
2995            nodeInfo->increment();
2996            delete newNode ;
2997            nodeInfo->decrement() ; } }
2998/*
2999  This node has been completely expanded and can be removed from the live
3000  set.
3001*/
3002      if (branchesLeft)
3003      { 
3004      }
3005      else
3006      { 
3007        if (!nodeInfo->numberBranchesLeft())
3008          nodeInfo->allBranchesGone(); // can clean up
3009        delete node ; }
3010    } else {
3011      // add cuts found to be infeasible (on bound)!
3012      abort();
3013      delete node;
3014    }
3015/*
3016  Delete cuts to get back to the original system.
3017
3018  I'm thinking this is redundant --- the call to addCuts that conditions entry
3019  to this code block also performs this action.
3020*/
3021      int numberToDelete = getNumRows()-numberRowsAtContinuous_ ;
3022      if (numberToDelete)
3023      { int * delRows = new int[numberToDelete] ;
3024        int i ;
3025        for (i = 0 ; i < numberToDelete ; i++)
3026        { delRows[i] = i+numberRowsAtContinuous_ ; }
3027        solver_->deleteRows(numberToDelete,delRows) ;
3028        delete [] delRows ; }
3029#else // end of not CBC_THREAD
3030#ifndef CBC_DETERMINISTIC_THREAD
3031      CbcNode *node = tree_->bestNode(cutoff) ;
3032      // Possible one on tree worse than cutoff
3033      if (!node||node->objectiveValue()>cutoff) 
3034        continue;
3035    if (!numberThreads_) {
3036#else
3037      if (!numberThreads_||(tree_->size()<5*numberThreads_&&!goneParallel)) {
3038      CbcNode *node = tree_->bestNode(cutoff) ;
3039      // Possible one on tree worse than cutoff
3040      if (!node||node->objectiveValue()>cutoff)
3041        continue;
3042#endif
3043      doOneNode(this,node,createdNode);
3044#ifdef CBC_DETERMINISTIC_THREAD
3045        assert (createdNode);
3046        if (!createdNode->active()) {
3047          //if (createdNode->nodeInfo()) {
3048          //createdNode->nodeInfo()->throwAway();
3049          //}
3050          delete createdNode;
3051          createdNode=NULL;
3052        } else {
3053          // Say one more pointing to this
3054          node->nodeInfo()->increment() ;
3055          tree_->push(createdNode) ;
3056        }
3057        //if (node) {
3058        //assert (node->active());
3059        if (node->active()) {
3060          assert (node->nodeInfo());
3061          if (node->nodeInfo()->numberBranchesLeft()) {
3062            tree_->push(node) ;
3063          } else {
3064            node->setActive(false);
3065          }
3066        } else {
3067          if (node->nodeInfo()) {
3068            if (!node->nodeInfo()->numberBranchesLeft())
3069              node->nodeInfo()->allBranchesGone(); // can clean up
3070            // So will delete underlying stuff
3071            node->setActive(true);
3072          }
3073          delNode[nDeleteNode++]=node;
3074          node=NULL;
3075        } 
3076        if (nDeleteNode>=MAX_DEL_NODE) {
3077          for (int i=0;i<nDeleteNode;i++) {
3078            //printf("trying to del %d %x\n",i,delNode[i]);
3079            delete delNode[i];
3080            //printf("done to del %d %x\n",i,delNode[i]);
3081          }
3082          nDeleteNode=0;
3083        }
3084#endif
3085      } else {
3086#ifdef CBC_NORMAL_THREAD
3087        threadStats[0]++;
3088        //need to think
3089        int iThread;
3090        // Start one off if any available
3091        for (iThread=0;iThread<numberThreads_;iThread++) {
3092          if (threadInfo[iThread].returnCode==-1) 
3093            break;
3094        }
3095        if (iThread<numberThreads_) {
3096          threadInfo[iThread].node=node;
3097          assert (threadInfo[iThread].returnCode==-1);
3098          // say in use
3099          threadInfo[iThread].returnCode=0;
3100          threadModel[iThread]->moveToModel(this,0);
3101          pthread_cond_signal(threadInfo[iThread].condition2); // unlock
3102          threadCount[iThread]++;
3103        }
3104#ifndef CBC_DETERMINISTIC_THREAD
3105        lockThread();
3106#endif
3107        locked=true;
3108        // see if any finished
3109        for (iThread=0;iThread<numberThreads_;iThread++) {
3110          if (threadInfo[iThread].returnCode>0) 
3111            break;
3112        }
3113#ifndef CBC_DETERMINISTIC_THREAD
3114        unlockThread();
3115#endif
3116        locked=false;
3117        if (iThread<numberThreads_) {
3118          threadModel[iThread]->moveToModel(this,1);
3119          assert (threadInfo[iThread].returnCode==1);
3120          // say available
3121          threadInfo[iThread].returnCode=-1;
3122          // carry on
3123          threadStats[3]++;
3124        } else {
3125          // Start one off if any available
3126          for (iThread=0;iThread<numberThreads_;iThread++) {
3127            if (threadInfo[iThread].returnCode==-1) 
3128              break;
3129          }
3130          if (iThread<numberThreads_) {
3131#ifndef CBC_DETERMINISTIC_THREAD
3132            lockThread();
3133#endif
3134            locked=true;
3135            // If any on tree get
3136            if (!tree_->empty()) {
3137              //node = tree_->bestNode(cutoff) ;
3138              //assert (node);
3139              threadStats[1]++;
3140              continue; // ** get another node
3141            }
3142#ifndef CBC_DETERMINISTIC_THREAD
3143            unlockThread();
3144#endif
3145            locked=false;
3146          }
3147          // wait (for debug could sleep and use test)
3148          bool finished=false;
3149          while (!finished) {
3150            pthread_mutex_lock(&condition_mutex);
3151            struct timespec absTime;
3152            clock_gettime(CLOCK_REALTIME,&absTime);
3153            double time = absTime.tv_sec+1.0e-9*absTime.tv_nsec;
3154            absTime.tv_nsec += 1000000; // millisecond
3155            if (absTime.tv_nsec>=1000000000) {
3156              absTime.tv_nsec -= 1000000000;
3157              absTime.tv_sec++;
3158            }
3159            pthread_cond_timedwait(&condition_main,&condition_mutex,&absTime);
3160            clock_gettime(CLOCK_REALTIME,&absTime);
3161            double time2 = absTime.tv_sec+1.0e-9*absTime.tv_nsec;
3162            timeWaiting += time2-time;
3163            pthread_mutex_unlock(&condition_mutex);
3164            for (iThread=0;iThread<numberThreads_;iThread++) {
3165              if (threadInfo[iThread].returnCode>0) {
3166                finished=true;
3167                break;
3168              } else if (threadInfo[iThread].returnCode==0) {
3169                pthread_cond_signal(threadInfo[iThread].condition2); // unlock
3170              }
3171            }
3172          }
3173          assert (iThread<numberThreads_);
3174          threadModel[iThread]->moveToModel(this,1);
3175          node = threadInfo[iThread].node;
3176          threadInfo[iThread].node=NULL;
3177          assert (threadInfo[iThread].returnCode==1);
3178          // say available
3179          threadInfo[iThread].returnCode=-1;
3180          // carry on
3181          threadStats[2]++;
3182        }
3183#else
3184        // Deterministic parallel
3185#ifndef CBC_DETERMINISTIC_THREAD
3186        abort();
3187#endif
3188#ifdef CBC_THREAD
3189        int saveTreeSize = tree_->size();
3190        goneParallel=true;
3191        int nAffected=splitModel(numberThreads_,threadModel,defaultParallelNodes);
3192        int saveTreeSize2 = tree_->size();
3193        int iThread;
3194        // do all until finished
3195        for (iThread=0;iThread<numberThreads_;iThread++) {
3196          // obviously tune
3197          threadInfo[iThread].nDeleteNode=defaultParallelIterations;
3198        }
3199        // Save current state
3200        int iObject;
3201        for (iObject=0;iObject<numberObjects_;iObject++) {
3202          saveObjects[iObject]->updateBefore(object_[iObject]);
3203        }
3204        for (iThread=0;iThread<numberThreads_;iThread++) {
3205          threadInfo[iThread].returnCode=0;
3206          pthread_cond_signal(threadInfo[iThread].condition2); // unlock
3207#if 0
3208          //wait!!
3209          bool finished=false;
3210          while (!finished) {
3211            pthread_mutex_lock(&condition_mutex);
3212            struct timespec absTime;
3213            clock_gettime(CLOCK_REALTIME,&absTime);
3214            double time = absTime.tv_sec+1.0e-9*absTime.tv_nsec;
3215            absTime.tv_nsec += 1000000; // millisecond
3216            if (absTime.tv_nsec>=1000000000) {
3217              absTime.tv_nsec -= 1000000000;
3218              absTime.tv_sec++;
3219            }
3220            pthread_cond_timedwait(&condition_main,&condition_mutex,&absTime);
3221            clock_gettime(CLOCK_REALTIME,&absTime);
3222            double time2 = absTime.tv_sec+1.0e-9*absTime.tv_nsec;
3223            timeWaiting += time2-time;
3224            pthread_mutex_unlock(&condition_mutex);
3225            finished=true;
3226            if (threadInfo[iThread].returnCode<=0) {
3227              finished=false;
3228            }
3229          }
3230#endif
3231        }
3232        // wait
3233        bool finished=false;
3234        while (!finished) {
3235          pthread_mutex_lock(&condition_mutex);
3236          struct timespec absTime;
3237          clock_gettime(CLOCK_REALTIME,&absTime);
3238          double time = absTime.tv_sec+1.0e-9*absTime.tv_nsec;
3239          absTime.tv_nsec += 1000000; // millisecond
3240          if (absTime.tv_nsec>=1000000000) {
3241            absTime.tv_nsec -= 1000000000;
3242            absTime.tv_sec++;
3243          }
3244          pthread_cond_timedwait(&condition_main,&condition_mutex,&absTime);
3245          clock_gettime(CLOCK_REALTIME,&absTime);
3246          double time2 = absTime.tv_sec+1.0e-9*absTime.tv_nsec;
3247          timeWaiting += time2-time;
3248          pthread_mutex_unlock(&condition_mutex);
3249          finished=true;
3250          for (iThread=0;iThread<numberThreads_;iThread++) {
3251            if (threadInfo[iThread].returnCode<=0) {
3252              finished=false;
3253            }
3254          }
3255        }
3256        // Unmark marked
3257        for (int i=0;i<nAffected;i++) {
3258          walkback_[i]->unmark();
3259        }
3260        assert (saveTreeSize2 == tree_->size());
3261        if (0) { 
3262          // put back cut counts
3263          for (int i=0;i<nAffected;i++) {
3264            walkback_[i]->decrementCuts(1000000);
3265          }
3266        }
3267#ifndef NDEBUG
3268        for (iObject=0;iObject<numberObjects_;iObject++) {
3269          CbcSimpleIntegerDynamicPseudoCost * obj =
3270            dynamic_cast <CbcSimpleIntegerDynamicPseudoCost *>(object_[iObject]) ;
3271          CbcSimpleIntegerDynamicPseudoCost * obj2 =
3272            dynamic_cast <CbcSimpleIntegerDynamicPseudoCost *>(saveObjects[iObject]) ;
3273          assert (obj->same(obj2));
3274        }
3275#endif
3276        int iModel;
3277        double scaleFactor=1.0;
3278        for (iModel=0;iModel<numberThreads_;iModel++) {
3279          //printf("model %d tree size %d\n",iModel,threadModel[iModel]->tree_->size());
3280          if (saveTreeSize>4*numberThreads_*defaultParallelNodes) {
3281            if (!threadModel[iModel]->tree_->size()) {
3282              scaleFactor *= 1.05;
3283            }
3284          }
3285          threadModel[iModel]->moveToModel(this,11);
3286          // Update base model
3287          OsiObject ** threadObject = threadModel[iModel]->object_;
3288          for (iObject=0;iObject<numberObjects_;iObject++) {
3289            object_[iObject]->updateAfter(threadObject[iObject],saveObjects[iObject]);
3290          }
3291        }
3292        if (scaleFactor!=1.0) {
3293          int newNumber = (int) (defaultParallelNodes * scaleFactor+0.5001);
3294          if (newNumber*2<defaultParallelIterations) {
3295            printf("Changing tree size from %d to %d\n",
3296                   defaultParallelNodes,newNumber);
3297            defaultParallelNodes = newNumber;
3298          }
3299        }
3300        printf("Tree sizes %d %d %d - affected %d\n",saveTreeSize,saveTreeSize2,tree_->size(),nAffected);
3301        // later remember random may not be thread neutral
3302#endif
3303#endif
3304      }
3305      //lastDepth=node->depth();
3306      //lastUnsatisfied=node->numberUnsatisfied();
3307#endif // end of CBC_THREAD
3308  }
3309#ifdef CBC_DETERMINISTIC_THREAD
3310  if (nDeleteNode) {
3311    for (int i=0;i<nDeleteNode;i++) {
3312      delete delNode[i];
3313    }
3314    nDeleteNode=0;
3315  }
3316#endif
3317#ifdef CBC_THREAD
3318  if (numberThreads_) {
3319    //printf("stats ");
3320    //for (unsigned int j=0;j<sizeof(threadStats)/sizeof(int);j++)
3321    //printf("%d ",threadStats[j]);
3322    //printf("\n");
3323    int i;
3324    // Seems to be bug in CoinCpu on Linux - does threads as well despite documentation
3325    double time=0.0;
3326    for (i=0;i<numberThreads_;i++) 
3327      time += threadInfo[i].timeInThread;
3328    bool goodTimer = time<(getCurrentSeconds());
3329    //bool stopped = (!(numberNodes_ < intParam_[CbcMaxNumNode] &&
3330    //        numberSolutions_ < intParam_[CbcMaxNumSol] &&
3331    //        totalTime < dblParam_[CbcMaximumSeconds] &&
3332    //        !stoppedOnGap_&&!eventHappened_));
3333    for (i=0;i<numberThreads_;i++) {
3334      while (threadInfo[i].returnCode==0) {
3335        pthread_cond_signal(threadInfo[i].condition2); // unlock
3336        pthread_mutex_lock(&condition_mutex);
3337        struct timespec absTime;
3338        clock_gettime(CLOCK_REALTIME,&absTime);
3339        absTime.tv_nsec += 1000000; // millisecond
3340        if (absTime.tv_nsec>=1000000000) {
3341          absTime.tv_nsec -= 1000000000;
3342          absTime.tv_sec++;
3343        }
3344        pthread_cond_timedwait(&condition_main,&condition_mutex,&absTime);
3345        clock_gettime(CLOCK_REALTIME,&absTime);
3346        pthread_mutex_unlock(&condition_mutex);
3347      }
3348      pthread_cond_signal(threadInfo[i].condition2); // unlock
3349      pthread_mutex_lock(&condition_mutex); // not sure necessary but have had one hang on interrupt
3350      threadModel[i]->numberThreads_=0; // say exit
3351#ifdef CBC_DETERMINISTIC_THREAD
3352      delete [] threadInfo[i].delNode;
3353#endif
3354      threadInfo[i].returnCode=0;
3355      pthread_mutex_unlock(&condition_mutex);
3356      pthread_cond_signal(threadInfo[i].condition2); // unlock
3357      //if (!stopped)
3358      //pthread_join(threadId[i],NULL);
3359      int returnCode;
3360      returnCode=pthread_join(threadId[i].thr,NULL);
3361      threadId[i].status = 0;
3362      assert (!returnCode);
3363        //else
3364        //pthread_kill(threadId[i]); // kill rather than try and synchronize
3365      threadModel[i]->moveToModel(this,2);
3366      pthread_mutex_destroy (threadInfo[i].mutex2);
3367      pthread_cond_destroy (threadInfo[i].condition2);
3368      assert (threadInfo[i].numberTimesLocked==threadInfo[i].numberTimesUnlocked);
3369      handler_->message(CBC_THREAD_STATS,messages_)
3370        <<"Thread";
3371      handler_->printing(true)
3372        <<i<<threadCount[i]<<threadInfo[i].timeWaitingToStart;
3373      handler_->printing(goodTimer)<<threadInfo[i].timeInThread;
3374      handler_->printing(false)<<0.0;
3375      handler_->printing(true)<<threadInfo[i].numberTimesLocked
3376        <<threadInfo[i].timeLocked<<threadInfo[i].timeWaitingToLock
3377        <<CoinMessageEol;
3378    }
3379    assert (threadInfo[numberThreads_].numberTimesLocked==threadInfo[numberThreads_].numberTimesUnlocked);
3380    handler_->message(CBC_THREAD_STATS,messages_)
3381      <<"Main thread";
3382    handler_->printing(false)<<0<<0<<0.0;
3383    handler_->printing(false)<<0.0;
3384    handler_->printing(true)<<timeWaiting;
3385    handler_->printing(true)<<threadInfo[numberThreads_].numberTimesLocked
3386      <<threadInfo[numberThreads_].timeLocked<<threadInfo[numberThreads_].timeWaitingToLock
3387      <<CoinMessageEol;
3388    pthread_mutex_destroy (&mutex);
3389    pthread_cond_destroy (&condition_main);
3390    pthread_mutex_destroy (&condition_mutex);
3391    // delete models (here in case some point to others)
3392    for (i=0;i<numberThreads_;i++) {
3393      delete threadModel[i];
3394    }
3395    delete [] mutex2;
3396    delete [] condition2;
3397    delete [] threadId;
3398    delete [] threadInfo;
3399    delete [] threadModel;
3400    delete [] threadCount;
3401    mutex_=NULL;
3402    // adjust time to allow for children on some systems
3403    dblParam_[CbcStartSeconds] -= CoinCpuTimeJustChildren();
3404  }
3405#endif
3406/*
3407  End of the non-abort actions. The next block of code is executed if we've
3408  aborted because we hit one of the limits. Clean up by deleting the live set
3409  and break out of the node processing loop. Note that on an abort, node may
3410  have been pushed back onto the tree for further processing, in which case
3411  it'll be deleted in cleanTree. We need to check.
3412*/
3413    if (!(numberNodes_ < intParam_[CbcMaxNumNode] &&
3414        numberSolutions_ < intParam_[CbcMaxNumSol] &&
3415        totalTime < dblParam_[CbcMaximumSeconds] &&
3416        !stoppedOnGap_&&!eventHappened_)) {
3417      if (tree_->size()) {
3418        double dummyBest;
3419        tree_->cleanTree(this,-COIN_DBL_MAX,dummyBest) ;
3420      }
3421      delete nextRowCut_;
3422      if (stoppedOnGap_)
3423        { messageHandler()->message(CBC_GAP,messages())
3424          << bestObjective_-bestPossibleObjective_
3425          << dblParam_[CbcAllowableGap]
3426          << dblParam_[CbcAllowableFractionGap]*100.0
3427          << CoinMessageEol ;
3428        secondaryStatus_ = 2;
3429        status_ = 0 ; }
3430        else
3431          if (isNodeLimitReached())
3432            { handler_->message(CBC_MAXNODES,messages_) << CoinMessageEol ;
3433            secondaryStatus_ = 3;
3434            status_ = 1 ; }
3435          else
3436        if (totalTime >= dblParam_[CbcMaximumSeconds])
3437          { handler_->message(CBC_MAXTIME,messages_) << CoinMessageEol ; 
3438          secondaryStatus_ = 4;
3439          status_ = 1 ; }
3440        else
3441          if (eventHappened_)
3442            { handler_->message(CBC_EVENT,messages_) << CoinMessageEol ; 
3443            secondaryStatus_ = 5;
3444            status_ = 5 ; }
3445          else
3446            { handler_->message(CBC_MAXSOLS,messages_) << CoinMessageEol ;
3447            secondaryStatus_ = 6;
3448            status_ = 1 ; }
3449    }
3450/*
3451  That's it, we've exhausted the search tree, or broken out of the loop because
3452  we hit some limit on evaluation.
3453
3454  We may have got an intelligent tree so give it one more chance
3455*/
3456  // Tell solver we are not in Branch and Cut
3457  solver_->setHintParam(OsiDoInBranchAndCut,false,OsiHintDo,NULL) ;
3458  tree_->endSearch();
3459  //  If we did any sub trees - did we give up on any?
3460  if ( numberStoppedSubTrees_)
3461    status_=1;
3462  if (!status_) {
3463    // Set best possible unless stopped on gap
3464    if(secondaryStatus_ != 2)
3465      bestPossibleObjective_=bestObjective_;
3466    handler_->message(CBC_END_GOOD,messages_)
3467      << bestObjective_ << numberIterations_ << numberNodes_<<getCurrentSeconds()
3468      << CoinMessageEol ;
3469  } else {
3470    handler_->message(CBC_END,messages_)
3471      << bestObjective_ <<bestPossibleObjective_
3472      << numberIterations_ << numberNodes_<<getCurrentSeconds()
3473      << CoinMessageEol ;
3474  }
3475  if (numberStrongIterations_)
3476    handler_->message(CBC_STRONG_STATS,messages_)
3477      << strongInfo_[0] << numberStrongIterations_ << strongInfo_[2]
3478      << strongInfo_[1] << CoinMessageEol ;
3479  if (!numberExtraNodes_) 
3480    handler_->message(CBC_OTHER_STATS,messages_)
3481      << maximumDepthActual_
3482      << numberDJFixed_ << CoinMessageEol ;
3483  else
3484    handler_->message(CBC_OTHER_STATS2,messages_)
3485      << maximumDepthActual_
3486      << numberDJFixed_ << numberExtraNodes_<<numberExtraIterations_
3487      <<CoinMessageEol ;
3488  if (doStatistics==100) {
3489    for (int i=0;i<numberObjects_;i++) {
3490      CbcSimpleIntegerDynamicPseudoCost * obj =
3491        dynamic_cast <CbcSimpleIntegerDynamicPseudoCost *>(object_[i]) ;
3492      if (obj)
3493        obj->print();
3494    }
3495  }
3496  if (statistics_) {
3497    // report in some way
3498    int * lookup = new int[numberObjects_];
3499    int i;
3500    for (i=0;i<numberObjects_;i++) 
3501      lookup[i]=-1;
3502    bool goodIds=true;
3503    for (i=0;i<numberObjects_;i++) {
3504      int iColumn = object_[i]->columnNumber();
3505      if(iColumn>=0&&iColumn<numberColumns) {
3506        if (lookup[i]==-1) {
3507          lookup[i]=iColumn;
3508        } else {
3509          goodIds=false;
3510          break;
3511        }
3512      } else {
3513        goodIds=false;
3514        break;
3515      }
3516    }
3517    if (!goodIds) {
3518      delete [] lookup;
3519      lookup=NULL;
3520    }
3521    if (doStatistics==3) {
3522      printf("  node parent depth column   value                    obj      inf\n");
3523      for ( i=0;i<numberNodes2_;i++) {
3524        statistics_[i]->print(lookup);
3525      }
3526    }
3527    if (doStatistics>1) {
3528      // Find last solution
3529      int k;
3530      for (k=numberNodes2_-1;k>=0;k--) {
3531        if (statistics_[k]->endingObjective()!=COIN_DBL_MAX&&
3532            !statistics_[k]->endingInfeasibility())
3533          break;
3534      }
3535      if (k>=0) {
3536        int depth=statistics_[k]->depth();
3537        int * which = new int[depth+1];
3538        for (i=depth;i>=0;i--) {
3539          which[i]=k;
3540          k=statistics_[k]->parentNode();
3541        }
3542        printf("  node parent depth column   value                    obj      inf\n");
3543        for (i=0;i<=depth;i++) {
3544          statistics_[which[i]]->print(lookup);
3545        }
3546        delete [] which;
3547      }
3548    }
3549    // now summary
3550    int maxDepth=0;
3551    double averageSolutionDepth=0.0;
3552    int numberSolutions=0;
3553    double averageCutoffDepth=0.0;
3554    double averageSolvedDepth=0.0;
3555    int numberCutoff=0;
3556    int numberDown=0;
3557    int numberFirstDown=0;
3558    double averageInfDown=0.0;
3559    double averageObjDown=0.0;
3560    int numberCutoffDown=0;
3561    int numberUp=0;
3562    int numberFirstUp=0;
3563    double averageInfUp=0.0;
3564    double averageObjUp=0.0;
3565    int numberCutoffUp=0;
3566    double averageNumberIterations1=0.0;
3567    double averageValue=0.0;
3568    for ( i=0;i<numberNodes2_;i++) {
3569      int depth =  statistics_[i]->depth(); 
3570      int way =  statistics_[i]->way(); 
3571      double value = statistics_[i]->value(); 
3572      double startingObjective =  statistics_[i]->startingObjective(); 
3573      int startingInfeasibility = statistics_[i]->startingInfeasibility(); 
3574      double endingObjective = statistics_[i]->endingObjective(); 
3575      int endingInfeasibility = statistics_[i]->endingInfeasibility(); 
3576      maxDepth = CoinMax(depth,maxDepth);
3577      // Only for completed
3578      averageNumberIterations1 += statistics_[i]->numberIterations();
3579      averageValue += value;
3580      if (endingObjective!=COIN_DBL_MAX&&!endingInfeasibility) {
3581        numberSolutions++;
3582        averageSolutionDepth += depth;
3583      }
3584      if (endingObjective==COIN_DBL_MAX) {
3585        numberCutoff++;
3586        averageCutoffDepth += depth;
3587        if (way<0) {
3588          numberDown++;
3589          numberCutoffDown++;
3590          if (way==-1)
3591            numberFirstDown++;
3592        } else {
3593          numberUp++;
3594          numberCutoffUp++;
3595          if (way==1)
3596            numberFirstUp++;
3597        }
3598      } else {
3599        averageSolvedDepth += depth;
3600        if (way<0) {
3601          numberDown++;
3602          averageInfDown += startingInfeasibility-endingInfeasibility;
3603          averageObjDown += endingObjective-startingObjective;
3604          if (way==-1)
3605            numberFirstDown++;
3606        } else {
3607          numberUp++;
3608          averageInfUp += startingInfeasibility-endingInfeasibility;
3609          averageObjUp += endingObjective-startingObjective;
3610          if (way==1)
3611            numberFirstUp++;
3612        }
3613      }
3614    }
3615    // Now print
3616    if (numberSolutions)
3617      averageSolutionDepth /= (double) numberSolutions;
3618    int numberSolved = numberNodes2_-numberCutoff;
3619    double averageNumberIterations2=numberIterations_-averageNumberIterations1
3620      -numberIterationsAtContinuous;
3621    if(numberCutoff) {
3622      averageCutoffDepth /= (double) numberCutoff;
3623      averageNumberIterations2 /= (double) numberCutoff;
3624    }
3625    if (numberNodes2_) 
3626      averageValue /= (double) numberNodes2_;
3627    if (numberSolved) {
3628      averageNumberIterations1 /= (double) numberSolved;
3629      averageSolvedDepth /= (double) numberSolved;
3630    }
3631    printf("%d solution(s) were found (by branching) at an average depth of %g\n",
3632           numberSolutions,averageSolutionDepth);
3633    printf("average value of variable being branched on was %g\n",
3634           averageValue);
3635    printf("%d nodes were cutoff at an average depth of %g with iteration count of %g\n",
3636           numberCutoff,averageCutoffDepth,averageNumberIterations2);
3637    printf("%d nodes were solved at an average depth of %g with iteration count of %g\n",
3638           numberSolved,averageSolvedDepth,averageNumberIterations1);
3639    if (numberDown) {
3640      averageInfDown /= (double) numberDown;
3641      averageObjDown /= (double) numberDown;
3642    }
3643    printf("Down %d nodes (%d first, %d second) - %d cutoff, rest decrease numinf %g increase obj %g\n",
3644           numberDown,numberFirstDown,numberDown-numberFirstDown,numberCutoffDown,
3645           averageInfDown,averageObjDown);
3646    if (numberUp) {
3647      averageInfUp /= (double) numberUp;
3648      averageObjUp /= (double) numberUp;
3649    }
3650    printf("Up %d nodes (%d first, %d second) - %d cutoff, rest decrease numinf %g increase obj %g\n",
3651           numberUp,numberFirstUp,numberUp-numberFirstUp,numberCutoffUp,
3652           averageInfUp,averageObjUp);
3653    for ( i=0;i<numberNodes2_;i++) 
3654      delete statistics_[i];
3655    delete [] statistics_;
3656    statistics_=NULL;
3657    maximumStatistics_=0;
3658    delete [] lookup;
3659  }
3660/*
3661  If we think we have a solution, restore and confirm it with a call to
3662  setBestSolution().  We need to reset the cutoff value so as not to fathom
3663  the solution on bounds.  Note that calling setBestSolution( ..., true)
3664  leaves the continuousSolver_ bounds vectors fixed at the solution value.
3665
3666  Running resolve() here is a failsafe --- setBestSolution has already
3667  reoptimised using the continuousSolver_. If for some reason we fail to
3668  prove optimality, run the problem again after instructing the solver to
3669  tell us more.
3670
3671  If all looks good, replace solver_ with continuousSolver_, so that the
3672  outside world will be able to obtain information about the solution using
3673  public methods.
3674*/
3675  if (bestSolution_&&(solverCharacteristics_->solverType()<2||solverCharacteristics_->solverType()==4)) 
3676  { setCutoff(1.0e50) ; // As best solution should be worse than cutoff
3677    phase_=5;
3678    double increment = getDblParam(CbcModel::CbcCutoffIncrement) ;
3679    if ((specialOptions_&4)==0)
3680      bestObjective_ += 100.0*increment+1.0e-3; // only set if we are going to solve
3681    setBestSolution(CBC_END_SOLUTION,bestObjective_,bestSolution_,1) ;
3682    continuousSolver_->resolve() ;
3683    if (!continuousSolver_->isProvenOptimal())
3684    { continuousSolver_->messageHandler()->setLogLevel(2) ;
3685      continuousSolver_->initialSolve() ; }
3686    delete solver_ ;
3687    // above deletes solverCharacteristics_
3688    solverCharacteristics_ = NULL;
3689    solver_ = continuousSolver_ ;
3690    setPointers(solver_);
3691    continuousSolver_ = NULL ; }
3692/*
3693  Clean up dangling objects. continuousSolver_ may already be toast.
3694*/
3695  delete lastws ;
3696  if (saveObjects) {
3697    for (int i=0;i<numberObjects_;i++)
3698      delete saveObjects[i];
3699    delete [] saveObjects;
3700  }
3701  numberStrong_ = saveNumberStrong;
3702  numberBeforeTrust_ = saveNumberBeforeTrust;
3703  delete [] whichGenerator_ ;
3704  whichGenerator_=NULL;
3705  delete [] lowerBefore ;
3706  delete [] upperBefore ;
3707  delete [] walkback_ ;
3708  walkback_ = NULL ;
3709  delete [] addedCuts_ ;
3710  addedCuts_ = NULL ;
3711  //delete persistentInfo;
3712  // Get rid of characteristics
3713  solverCharacteristics_=NULL;
3714  if (continuousSolver_)
3715  { delete continuousSolver_ ;
3716    continuousSolver_ = NULL ; }
3717/*
3718  Destroy global cuts by replacing with an empty OsiCuts object.
3719*/
3720  globalCuts_= OsiCuts() ;
3721  if (!bestSolution_) {
3722    // make sure lp solver is infeasible
3723    int numberColumns = solver_->getNumCols();
3724    const double * columnLower = solver_->getColLower();
3725    int iColumn;
3726    for (iColumn=0;iColumn<numberColumns;iColumn++) {
3727      if (solver_->isInteger(iColumn)) 
3728        solver_->setColUpper(iColumn,columnLower[iColumn]);
3729    }
3730    solver_->initialSolve();
3731  }
3732#ifdef COIN_HAS_CLP
3733  {
3734    OsiClpSolverInterface * clpSolver
3735      = dynamic_cast<OsiClpSolverInterface *> (solver_);
3736    if (clpSolver) {
3737      // Possible restore of pivot method
3738      if(savePivotMethod) {
3739        // model may have changed
3740        savePivotMethod->setModel(NULL);
3741        clpSolver->getModelPtr()->setDualRowPivotAlgorithm(*savePivotMethod);
3742        delete savePivotMethod;
3743      }
3744      clpSolver->setLargestAway(-1.0);
3745    }
3746  }
3747#endif
3748  if(fastNodeDepth_>=1000000&&!parentModel_) {
3749    // delete object off end
3750    delete object_[numberObjects_];
3751    fastNodeDepth_ -= 1000000;
3752  }
3753  delete saveSolver;
3754#if COIN_DEVELOP>1
3755  void print_fac_stats();
3756  //if (!parentModel_)
3757  //  print_fac_stats();
3758#endif
3759  if (strategy_&&strategy_->preProcessState()>0) {
3760    // undo preprocessing
3761    CglPreProcess * process = strategy_->process();
3762    assert (process);
3763    int n = originalSolver->getNumCols();
3764    if (bestSolution_) {
3765      delete [] bestSolution_;
3766      bestSolution_ = new double [n];
3767      process->postProcess(*solver_);
3768    }
3769    strategy_->deletePreProcess();
3770    // Solution now back in originalSolver
3771    delete solver_;
3772    solver_=originalSolver;
3773    if (bestSolution_)
3774      memcpy(bestSolution_,solver_->getColSolution(),n*sizeof(double));
3775    // put back original objects if there were any
3776    if (originalObject) {
3777      int iColumn;
3778      assert (ownObjects_);
3779      for (iColumn=0;iColumn<numberObjects_;iColumn++) 
3780        delete object_[iColumn];
3781      delete [] object_;
3782      numberObjects_ = numberOriginalObjects;
3783      object_=originalObject;
3784      delete [] integerVariable_;
3785      numberIntegers_=0;
3786      for (iColumn=0;iColumn<n;iColumn++) {
3787        if (solver_->isInteger(iColumn))
3788          numberIntegers_++;
3789      }
3790      integerVariable_ = new int[numberIntegers_];
3791      numberIntegers_=0;
3792      for (iColumn=0;iColumn<n;iColumn++) {
3793        if (solver_->isInteger(iColumn))
3794            integerVariable_[numberIntegers_++]=iColumn;
3795      }
3796    }
3797  }
3798#ifdef CLP_QUICK_OPTIONS
3799  {
3800    OsiClpSolverInterface * clpSolver
3801      = dynamic_cast<OsiClpSolverInterface *> (solver_);
3802    if (clpSolver) {
3803      // Try and re-use regions   
3804      ClpSimplex * simplex = clpSolver->getModelPtr();
3805      simplex->setPersistenceFlag(0);
3806      clpSolver->deleteScaleFactors();
3807      clpSolver->setSpecialOptions(clpSolver->specialOptions()&(~131072));
3808      simplex->setSpecialOptions(simplex->specialOptions()&(~131072));
3809      simplex->setAlphaAccuracy(-1.0);
3810      //clpSolver->setSpecialOptions((clpSolver->specialOptions()&~128)|65536);
3811    }
3812  }
3813#endif
3814  return ;
3815 }
3816
3817
3818// Solve the initial LP relaxation
3819void 
3820CbcModel::initialSolve() 
3821{
3822  assert (solver_);
3823  assert (!solverCharacteristics_);
3824  OsiBabSolver * solverCharacteristics = dynamic_cast<OsiBabSolver *> (solver_->getAuxiliaryInfo());
3825  if (solverCharacteristics) {
3826    solverCharacteristics_ = solverCharacteristics;
3827  } else {
3828    // replace in solver
3829    OsiBabSolver defaultC;
3830    solver_->setAuxiliaryInfo(&defaultC);
3831    solverCharacteristics_ = dynamic_cast<OsiBabSolver *> (solver_->getAuxiliaryInfo());
3832  }
3833  solverCharacteristics_->setSolver(solver_);
3834  solver_->setHintParam(OsiDoInBranchAndCut,true,OsiHintDo,NULL) ;
3835  solver_->initialSolve();
3836  solver_->setHintParam(OsiDoInBranchAndCut,false,OsiHintDo,NULL) ;
3837  // But set up so Jon Lee will be happy
3838  status_=-1;
3839  secondaryStatus_ = -1;
3840  originalContinuousObjective_ = solver_->getObjValue()*solver_->getObjSense();
3841  delete [] continuousSolution_;
3842  continuousSolution_ = CoinCopyOfArray(solver_->getColSolution(),
3843                                             solver_->getNumCols());
3844  setPointers(solver_);
3845  solverCharacteristics_ = NULL;
3846}
3847
3848/*! \brief Get an empty basis object
3849
3850  Return an empty CoinWarmStartBasis object with the requested capacity,
3851  appropriate for the current solver. The object is cloned from the object
3852  cached as emptyWarmStart_. If there is no cached object, the routine
3853  queries the solver for a warm start object, empties it, and caches the
3854  result.
3855*/
3856
3857CoinWarmStartBasis *CbcModel::getEmptyBasis (int ns, int na) const
3858
3859{ CoinWarmStartBasis *emptyBasis ;
3860/*
3861  Acquire an empty basis object, if we don't yet have one.
3862*/
3863  if (emptyWarmStart_ == 0)
3864  { if (solver_ == 0)
3865    { throw CoinError("Cannot construct basis without solver!",
3866                      "getEmptyBasis","CbcModel") ; }
3867    emptyBasis =
3868        dynamic_cast<CoinWarmStartBasis *>(solver_->getEmptyWarmStart()) ;
3869    if (emptyBasis == 0)
3870    { throw CoinError(
3871        "Solver does not appear to use a basis-oriented warm start.",
3872        "getEmptyBasis","CbcModel") ; }
3873    emptyBasis->setSize(0,0) ;
3874    emptyWarmStart_ = dynamic_cast<CoinWarmStart *>(emptyBasis) ; }
3875/*
3876  Clone the empty basis object, resize it as requested, and return.
3877*/
3878  emptyBasis = dynamic_cast<CoinWarmStartBasis *>(emptyWarmStart_->clone()) ;
3879  assert(emptyBasis) ;
3880  if (ns != 0 || na != 0) emptyBasis->setSize(ns,na) ;
3881
3882  return (emptyBasis) ; }
3883   
3884
3885/** Default Constructor
3886
3887  Creates an empty model without an associated solver.
3888*/
3889CbcModel::CbcModel() 
3890
3891:
3892  solver_(NULL),
3893  ownership_(0x80000000),
3894  continuousSolver_(NULL),
3895  referenceSolver_(NULL),
3896  defaultHandler_(true),
3897  emptyWarmStart_(NULL),
3898  bestObjective_(COIN_DBL_MAX),
3899  bestPossibleObjective_(COIN_DBL_MAX),
3900  sumChangeObjective1_(0.0),
3901  sumChangeObjective2_(0.0),
3902  bestSolution_(NULL),
3903  currentSolution_(NULL),
3904  testSolution_(NULL),
3905  minimumDrop_(1.0e-4),
3906  numberSolutions_(0),
3907  stateOfSearch_(0),
3908  hotstartSolution_(NULL),
3909  hotstartPriorities_(NULL),
3910  numberHeuristicSolutions_(0),
3911  numberNodes_(0),
3912  numberNodes2_(0),
3913  numberIterations_(0),
3914  status_(-1),
3915  secondaryStatus_(-1),
3916  numberIntegers_(0),
3917  numberRowsAtContinuous_(0),
3918  maximumNumberCuts_(0),
3919  phase_(0),
3920  currentNumberCuts_(0),
3921  maximumDepth_(0),
3922  walkback_(NULL),
3923  addedCuts_(NULL),
3924  nextRowCut_(NULL),
3925  currentNode_(NULL),
3926  integerVariable_(NULL),
3927  integerInfo_(NULL),
3928  continuousSolution_(NULL),
3929  usedInSolution_(NULL),
3930  specialOptions_(0),
3931  subTreeModel_(NULL),
3932  numberStoppedSubTrees_(0),
3933  mutex_(NULL),
3934  presolve_(0),
3935  numberStrong_(5),
3936  numberBeforeTrust_(10),
3937  numberPenalties_(20),
3938  stopNumberIterations_(-1),
3939  penaltyScaleFactor_(3.0),
3940  numberAnalyzeIterations_(0),
3941  analyzeResults_(NULL),
3942  numberInfeasibleNodes_(0),
3943  problemType_(0),
3944  printFrequency_(0),
3945  numberCutGenerators_(0),
3946  generator_(NULL),
3947  virginGenerator_(NULL),
3948  numberHeuristics_(0),
3949  heuristic_(NULL),
3950  lastHeuristic_(NULL),
3951# ifdef COIN_HAS_CLP
3952  fastNodeDepth_(-1),
3953#endif
3954  eventHandler_(NULL),
3955  numberObjects_(0),
3956  object_(NULL),
3957  ownObjects_(true),
3958  originalColumns_(NULL),
3959  howOftenGlobalScan_(1),
3960  numberGlobalViolations_(0),
3961  numberExtraIterations_(0),
3962  numberExtraNodes_(0),
3963  continuousObjective_(COIN_DBL_MAX),
3964  originalContinuousObjective_(COIN_DBL_MAX),
3965  continuousInfeasibilities_(COIN_INT_MAX),
3966  maximumCutPassesAtRoot_(20),
3967  maximumCutPasses_(10),
3968  preferredWay_(0),
3969  currentPassNumber_(0),
3970  maximumWhich_(1000),
3971  maximumRows_(0),
3972  currentDepth_(0),
3973  whichGenerator_(NULL),
3974  maximumStatistics_(0),
3975  statistics_(NULL),
3976  maximumDepthActual_(0),
3977  numberDJFixed_(0.0),
3978  probingInfo_(NULL),
3979  numberFixedAtRoot_(0),
3980  numberFixedNow_(0),
3981  stoppedOnGap_(false),
3982  eventHappened_(false),
3983  numberLongStrong_(0),
3984  numberOldActiveCuts_(0),
3985  numberNewCuts_(0),
3986  sizeMiniTree_(0),
3987  searchStrategy_(-1),
3988  numberStrongIterations_(0),
3989  resolveAfterTakeOffCuts_(true),
3990  maximumNumberIterations_(-1),
3991#if NEW_UPDATE_OBJECT>1
3992  numberUpdateItems_(0),
3993  maximumNumberUpdateItems_(0),
3994  updateItems_(NULL),
3995#endif
3996  numberThreads_(0),
3997  threadMode_(0)
3998{
3999  memset(intParam_,0,sizeof(intParam_));
4000  intParam_[CbcMaxNumNode] = 2147483647;
4001  intParam_[CbcMaxNumSol] = 9999999;
4002  intParam_[CbcFathomDiscipline] = 0;
4003
4004  dblParam_[CbcIntegerTolerance] = 1e-6;
4005  dblParam_[CbcInfeasibilityWeight] = 0.0;
4006  dblParam_[CbcCutoffIncrement] = 1e-5;
4007  dblParam_[CbcAllowableGap] = 1.0e-10;
4008  dblParam_[CbcAllowableFractionGap] = 0.0;
4009  dblParam_[CbcMaximumSeconds] = 1.0e100;
4010  dblParam_[CbcCurrentCutoff] = 1.0e100;
4011  dblParam_[CbcOptimizationDirection] = 1.0;
4012  dblParam_[CbcCurrentObjectiveValue] = 1.0e100;
4013  dblParam_[CbcCurrentMinimizationObjectiveValue] = 1.0e100;
4014  dblParam_[CbcStartSeconds] = 0.0;
4015  strongInfo_[0]=0;
4016  strongInfo_[1]=0;
4017  strongInfo_[2]=0;
4018  solverCharacteristics_ = NULL;
4019  nodeCompare_=new CbcCompareDefault();;
4020  problemFeasibility_=new CbcFeasibilityBase();
4021  tree_= new CbcTree();
4022  branchingMethod_=NULL;
4023  cutModifier_=NULL;
4024  strategy_=NULL;
4025  parentModel_=NULL;
4026  cbcColLower_ = NULL;
4027  cbcColUpper_ = NULL;
4028  cbcRowLower_ = NULL;
4029  cbcRowUpper_ = NULL;
4030  cbcColSolution_ = NULL;
4031  cbcRowPrice_ = NULL;
4032  cbcReducedCost_ = NULL;
4033  cbcRowActivity_ = NULL;
4034  appData_=NULL;
4035  handler_ = new CoinMessageHandler();
4036  handler_->setLogLevel(2);
4037  messages_ = CbcMessage();
4038  eventHandler_ = new CbcEventHandler() ;
4039}
4040
4041/** Constructor from solver.
4042
4043  Creates a model complete with a clone of the solver passed as a parameter.
4044*/
4045
4046CbcModel::CbcModel(const OsiSolverInterface &rhs)
4047:
4048  continuousSolver_(NULL),
4049  referenceSolver_(NULL),
4050  defaultHandler_(true),
4051  emptyWarmStart_(NULL),
4052  bestObjective_(COIN_DBL_MAX),
4053  bestPossibleObjective_(COIN_DBL_MAX),
4054  sumChangeObjective1_(0.0),
4055  sumChangeObjective2_(0.0),
4056  minimumDrop_(1.0e-4),
4057  numberSolutions_(0),
4058  stateOfSearch_(0),
4059  hotstartSolution_(NULL),
4060  hotstartPriorities_(NULL),
4061  numberHeuristicSolutions_(0),
4062  numberNodes_(0),
4063  numberNodes2_(0),
4064  numberIterations_(0),
4065  status_(-1),
4066  secondaryStatus_(-1),
4067  numberRowsAtContinuous_(0),
4068  maximumNumberCuts_(0),
4069  phase_(0),
4070  currentNumberCuts_(0),
4071  maximumDepth_(0),
4072  walkback_(NULL),
4073  addedCuts_(NULL),
4074  nextRowCut_(NULL),
4075  currentNode_(NULL),
4076  integerInfo_(NULL),
4077  specialOptions_(0),
4078  subTreeModel_(NULL),
4079  numberStoppedSubTrees_(0),
4080  mutex_(NULL),
4081  presolve_(0),
4082  numberStrong_(5),
4083  numberBeforeTrust_(10),
4084  numberPenalties_(20),
4085  stopNumberIterations_(-1),
4086  penaltyScaleFactor_(3.0),
4087  numberAnalyzeIterations_(0),
4088  analyzeResults_(NULL),
4089  numberInfeasibleNodes_(0),
4090  problemType_(0),
4091  printFrequency_(0),
4092  numberCutGenerators_(0),
4093  generator_(NULL),
4094  virginGenerator_(NULL),
4095  numberHeuristics_(0),
4096  heuristic_(NULL),
4097  lastHeuristic_(NULL),
4098# ifdef COIN_HAS_CLP
4099  fastNodeDepth_(-1),
4100#endif
4101  eventHandler_(NULL),
4102  numberObjects_(0),
4103  object_(NULL),
4104  ownObjects_(true),
4105  originalColumns_(NULL),
4106  howOftenGlobalScan_(1),
4107  numberGlobalViolations_(0),
4108  numberExtraIterations_(0),
4109  numberExtraNodes_(0),
4110  continuousObjective_(COIN_DBL_MAX),
4111  originalContinuousObjective_(COIN_DBL_MAX),
4112  continuousInfeasibilities_(COIN_INT_MAX),
4113  maximumCutPassesAtRoot_(20),
4114  maximumCutPasses_(10),
4115  preferredWay_(0),
4116  currentPassNumber_(0),
4117  maximumWhich_(1000),
4118  maximumRows_(0),
4119  currentDepth_(0),
4120  whichGenerator_(NULL),
4121  maximumStatistics_(0),
4122  statistics_(NULL),
4123  maximumDepthActual_(0),
4124  numberDJFixed_(0.0),
4125  probingInfo_(NULL),
4126  numberFixedAtRoot_(0),
4127  numberFixedNow_(0),
4128  stoppedOnGap_(false),
4129  eventHappened_(false),
4130  numberLongStrong_(0),
4131  numberOldActiveCuts_(0),
4132  numberNewCuts_(0),
4133  sizeMiniTree_(0),
4134  searchStrategy_(-1),
4135  numberStrongIterations_(0),
4136  resolveAfterTakeOffCuts_(true),
4137  maximumNumberIterations_(-1),
4138#if NEW_UPDATE_OBJECT>1
4139  numberUpdateItems_(0),
4140  maximumNumberUpdateItems_(0),
4141  updateItems_(NULL),
4142#endif
4143  numberThreads_(0),
4144  threadMode_(0)
4145{
4146  memset(intParam_,0,sizeof(intParam_));
4147  intParam_[CbcMaxNumNode] = 2147483647;
4148  intParam_[CbcMaxNumSol] = 9999999;
4149  intParam_[CbcFathomDiscipline] = 0;
4150
4151  dblParam_[CbcIntegerTolerance] = 1e-6;
4152  dblParam_[CbcInfeasibilityWeight] = 0.0;
4153  dblParam_[CbcCutoffIncrement] = 1e-5;
4154  dblParam_[CbcAllowableGap] = 1.0e-10;
4155  dblParam_[CbcAllowableFractionGap] = 0.0;
4156  dblParam_[CbcMaximumSeconds] = 1.0e100;
4157  dblParam_[CbcCurrentCutoff] = 1.0e100;
4158  dblParam_[CbcOptimizationDirection] = 1.0;
4159  dblParam_[CbcCurrentObjectiveValue] = 1.0e100;
4160  dblParam_[CbcCurrentMinimizationObjectiveValue] = 1.0e100;
4161  dblParam_[CbcStartSeconds] = 0.0;
4162  strongInfo_[0]=0;
4163  strongInfo_[1]=0;
4164  strongInfo_[2]=0;
4165  solverCharacteristics_ = NULL;
4166
4167  nodeCompare_=new CbcCompareDefault();;
4168  problemFeasibility_=new CbcFeasibilityBase();
4169  tree_= new CbcTree();
4170  branchingMethod_=NULL;
4171  cutModifier_=NULL;
4172  strategy_=NULL;
4173  parentModel_=NULL;
4174  appData_=NULL;
4175  handler_ = new CoinMessageHandler();
4176  handler_->setLogLevel(2);
4177  messages_ = CbcMessage();
4178  eventHandler_ = new CbcEventHandler() ;
4179  solver_ = rhs.clone();
4180  referenceSolver_ = solver_->clone();
4181  ownership_ = 0x80000000;
4182  cbcColLower_ = NULL;
4183  cbcColUpper_ = NULL;
4184  cbcRowLower_ = NULL;
4185  cbcRowUpper_ = NULL;
4186  cbcColSolution_ = NULL;
4187  cbcRowPrice_ = NULL;
4188  cbcReducedCost_ = NULL;
4189  cbcRowActivity_ = NULL;
4190
4191  // Initialize solution and integer variable vectors
4192  bestSolution_ = NULL; // to say no solution found
4193  numberIntegers_=0;
4194  int numberColumns = solver_->getNumCols();
4195  int iColumn;
4196  if (numberColumns) {
4197    // Space for current solution
4198    currentSolution_ = new double[numberColumns];
4199    continuousSolution_ = new double[numberColumns];
4200    usedInSolution_ = new int[numberColumns];
4201    CoinZeroN(usedInSolution_,numberColumns);
4202    for (iColumn=0;iColumn<numberColumns;iColumn++) {
4203      if( solver_->isInteger(iColumn)) 
4204        numberIntegers_++;
4205    }
4206  } else {
4207    // empty model
4208    currentSolution_=NULL;
4209    continuousSolution_=NULL;
4210    usedInSolution_=NULL;
4211  }
4212  testSolution_=currentSolution_;
4213  if (numberIntegers_) {
4214    integerVariable_ = new int [numberIntegers_];
4215    numberIntegers_=0;
4216    for (iColumn=0;iColumn<numberColumns;iColumn++) {
4217      if( solver_->isInteger(iColumn)) 
4218        integerVariable_[numberIntegers_++]=iColumn;
4219    }
4220  } else {
4221    integerVariable_ = NULL;
4222  }
4223}
4224
4225/*
4226  Assign a solver to the model (model assumes ownership)
4227
4228  The integer variable vector is initialized if it's not already present.
4229  If deleteSolver then current solver deleted (if model owned)
4230
4231  Assuming ownership matches usage in OsiSolverInterface
4232  (cf. assignProblem, loadProblem).
4233
4234  TODO: What to do about solver parameters? A simple copy likely won't do it,
4235        because the SI must push the settings into the underlying solver. In
4236        the context of switching solvers in cbc, this means that command line
4237        settings will get lost. Stash the command line somewhere and reread it
4238        here, maybe?
4239 
4240  TODO: More generally, how much state should be transferred from the old
4241        solver to the new solver? Best perhaps to see how usage develops.
4242        What's done here mimics the CbcModel(OsiSolverInterface) constructor.
4243*/
4244void
4245CbcModel::assignSolver(OsiSolverInterface *&solver, bool deleteSolver)
4246
4247{
4248  // resize best solution if exists
4249  if (bestSolution_&&solver&&solver_) {
4250    int nOld = solver_->getNumCols();
4251    int nNew = solver->getNumCols();
4252    if (nNew>nOld) {
4253      double * temp = new double[nNew];
4254      memcpy(temp,bestSolution_,nOld*sizeof(double));
4255      memset(temp+nOld,0,(nNew-nOld)*sizeof(double));
4256      delete [] bestSolution_;
4257      bestSolution_=temp;
4258    }
4259  }
4260  // Keep the current message level for solver (if solver exists)
4261  if (solver_)
4262    solver->messageHandler()->setLogLevel(solver_->messageHandler()->logLevel()) ;
4263
4264  if (modelOwnsSolver()&&deleteSolver) delete solver_ ;
4265  solver_ = solver;
4266  solver = NULL ;
4267  setModelOwnsSolver(true) ;
4268/*
4269  Basis information is solver-specific.
4270*/
4271  if (emptyWarmStart_)
4272  { delete emptyWarmStart_  ;
4273    emptyWarmStart_ = 0 ; }
4274  bestSolutionBasis_ = CoinWarmStartBasis();
4275/*
4276  Initialize integer variable vector.
4277*/
4278  numberIntegers_=0;
4279  int numberColumns = solver_->getNumCols();
4280  int iColumn;
4281  for (iColumn=0;iColumn<numberColumns;iColumn++) {
4282    if( solver_->isInteger(iColumn)) 
4283      numberIntegers_++;
4284  }
4285  delete [] integerVariable_;
4286  if (numberIntegers_) {
4287    integerVariable_ = new int [numberIntegers_];
4288    numberIntegers_=0;
4289    for (iColumn=0;iColumn<numberColumns;iColumn++) {
4290      if( solver_->isInteger(iColumn)) 
4291        integerVariable_[numberIntegers_++]=iColumn;
4292    }
4293  } else {
4294    integerVariable_ = NULL;
4295  }
4296
4297  return ;
4298}
4299
4300// Copy constructor.
4301
4302CbcModel::CbcModel(const CbcModel & rhs, bool noTree)
4303:
4304  continuousSolver_(NULL),
4305  referenceSolver_(NULL),
4306  defaultHandler_(rhs.defaultHandler_),
4307  emptyWarmStart_(NULL),
4308  bestObjective_(rhs.bestObjective_),
4309  bestPossibleObjective_(rhs.bestPossibleObjective_),
4310  sumChangeObjective1_(rhs.sumChangeObjective1_),
4311  sumChangeObjective2_(rhs.sumChangeObjective2_),
4312  minimumDrop_(rhs.minimumDrop_),
4313  numberSolutions_(rhs.numberSolutions_),
4314  stateOfSearch_(rhs.stateOfSearch_),
4315  numberHeuristicSolutions_(rhs.numberHeuristicSolutions_),
4316  numberNodes_(rhs.numberNodes_),
4317  numberNodes2_(rhs.numberNodes2_),
4318  numberIterations_(rhs.numberIterations_),
4319  status_(rhs.status_),
4320  secondaryStatus_(rhs.secondaryStatus_),
4321  specialOptions_(rhs.specialOptions_),
4322  subTreeModel_(rhs.subTreeModel_),
4323  numberStoppedSubTrees_(rhs.numberStoppedSubTrees_),
4324  mutex_(NULL),
4325  presolve_(rhs.presolve_),
4326  numberStrong_(rhs.numberStrong_),
4327  numberBeforeTrust_(rhs.numberBeforeTrust_),
4328  numberPenalties_(rhs.numberPenalties_),
4329  stopNumberIterations_(rhs.stopNumberIterations_),
4330  penaltyScaleFactor_(rhs.penaltyScaleFactor_),
4331  numberAnalyzeIterations_(rhs.numberAnalyzeIterations_),
4332  analyzeResults_(NULL),
4333  numberInfeasibleNodes_(rhs.numberInfeasibleNodes_),
4334  problemType_(rhs.problemType_),
4335  printFrequency_(rhs.printFrequency_),
4336# ifdef COIN_HAS_CLP
4337  fastNodeDepth_(rhs.fastNodeDepth_),
4338#endif
4339  howOftenGlobalScan_(rhs.howOftenGlobalScan_),
4340  numberGlobalViolations_(rhs.numberGlobalViolations_),
4341  numberExtraIterations_(rhs.numberExtraIterations_),
4342  numberExtraNodes_(rhs.numberExtraNodes_),
4343  continuousObjective_(rhs.continuousObjective_),
4344  originalContinuousObjective_(rhs.originalContinuousObjective_),
4345  continuousInfeasibilities_(rhs.continuousInfeasibilities_),
4346  maximumCutPassesAtRoot_(rhs.maximumCutPassesAtRoot_),
4347  maximumCutPasses_( rhs.maximumCutPasses_),
4348  preferredWay_(rhs.preferredWay_),
4349  currentPassNumber_(rhs.currentPassNumber_),
4350  maximumWhich_(rhs.maximumWhich_),
4351  maximumRows_(0),
4352  currentDepth_(0),
4353  whichGenerator_(NULL),
4354  maximumStatistics_(0),
4355  statistics_(NULL),
4356  maximumDepthActual_(0),
4357  numberDJFixed_(0.0),
4358  probingInfo_(NULL),
4359  numberFixedAtRoot_(rhs.numberFixedAtRoot_),
4360  numberFixedNow_(rhs.numberFixedNow_),
4361  stoppedOnGap_(rhs.stoppedOnGap_),
4362  eventHappened_(rhs.eventHappened_),
4363  numberLongStrong_(rhs.numberLongStrong_),
4364  numberOldActiveCuts_(rhs.numberOldActiveCuts_),
4365  numberNewCuts_(rhs.numberNewCuts_),
4366  sizeMiniTree_(rhs.sizeMiniTree_),
4367  searchStrategy_(rhs.searchStrategy_),
4368  numberStrongIterations_(rhs.numberStrongIterations_),
4369  resolveAfterTakeOffCuts_(rhs.resolveAfterTakeOffCuts_),
4370  maximumNumberIterations_(rhs.maximumNumberIterations_),
4371#if NEW_UPDATE_OBJECT>1
4372  numberUpdateItems_(rhs.numberUpdateItems_),
4373  maximumNumberUpdateItems_(rhs.maximumNumberUpdateItems_),
4374  updateItems_(NULL),
4375#endif
4376  numberThreads_(rhs.numberThreads_),
4377  threadMode_(rhs.threadMode_)
4378{
4379  memcpy(intParam_,rhs.intParam_,sizeof(intParam_));
4380  memcpy(dblParam_,rhs.dblParam_,sizeof(dblParam_));
4381  strongInfo_[0]=rhs.strongInfo_[0];
4382  strongInfo_[1]=rhs.strongInfo_[1];
4383  strongInfo_[2]=rhs.strongInfo_[2];
4384  solverCharacteristics_ = NULL;
4385  if (rhs.emptyWarmStart_) emptyWarmStart_ = rhs.emptyWarmStart_->clone() ;
4386  if (defaultHandler_) {
4387    handler_ = new CoinMessageHandler();
4388    handler_->setLogLevel(2);
4389  } else {
4390    handler_ = rhs.handler_;
4391  }
4392  messageHandler()->setLogLevel(rhs.messageHandler()->logLevel());
4393  numberCutGenerators_ = rhs.numberCutGenerators_;
4394  if (numberCutGenerators_) {
4395    generator_ = new CbcCutGenerator * [numberCutGenerators_];
4396    virginGenerator_ = new CbcCutGenerator * [numberCutGenerators_];
4397    int i;
4398    for (i=0;i<numberCutGenerators_;i++) {
4399      generator_[i]=new CbcCutGenerator(*rhs.generator_[i]);
4400      virginGenerator_[i]=new CbcCutGenerator(*rhs.virginGenerator_[i]);
4401    }
4402  } else {
4403    generator_=NULL;
4404    virginGenerator_=NULL;
4405  }
4406  if (!noTree)
4407    globalCuts_ = rhs.globalCuts_;
4408  numberHeuristics_ = rhs.numberHeuristics_;
4409  if (numberHeuristics_) {
4410    heuristic_ = new CbcHeuristic * [numberHeuristics_];
4411    int i;
4412    for (i=0;i<numberHeuristics_;i++) {
4413      heuristic_[i]=rhs.heuristic_[i]->clone();
4414    }
4415  } else {
4416    heuristic_=NULL;
4417  }
4418  lastHeuristic_ = NULL;
4419  if (rhs.eventHandler_)
4420    { eventHandler_ = rhs.eventHandler_->clone() ; }
4421  else
4422  { eventHandler_ = NULL ; }
4423  ownObjects_ = rhs.ownObjects_;
4424  if (ownObjects_) {
4425    numberObjects_=rhs.numberObjects_;
4426    if (numberObjects_) {
4427      object_ = new OsiObject * [numberObjects_];
4428      int i;
4429      for (i=0;i<numberObjects_;i++) {
4430        object_[i]=(rhs.object_[i])->clone();
4431        CbcObject * obj = dynamic_cast <CbcObject *>(object_[i]) ;
4432        // Could be OsiObjects
4433        if (obj)
4434          obj->setModel(this);
4435      }
4436    } else {
4437      object_=NULL;
4438    }
4439  } else {
4440    // assume will be redone
4441    numberObjects_=0;
4442    object_=NULL;
4443  }
4444  if (rhs.referenceSolver_)
4445    referenceSolver_ = rhs.referenceSolver_->clone();
4446  else
4447    referenceSolver_=NULL;
4448  if (!noTree||!rhs.continuousSolver_)
4449    solver_ = rhs.solver_->clone();
4450  else
4451    solver_ = rhs.continuousSolver_->clone();
4452  if (rhs.originalColumns_) {
4453    int numberColumns = solver_->getNumCols();
4454    originalColumns_= new int [numberColumns];
4455    memcpy(originalColumns_,rhs.originalColumns_,numberColumns*sizeof(int));
4456  } else {
4457    originalColumns_=NULL;
4458  }
4459#if NEW_UPDATE_OBJECT>1
4460  if (maximumNumberUpdateItems_) {
4461    updateItems_ = new CbcObjectUpdateData [maximumNumberUpdateItems_];
4462    for (int i=0;i<maximumNumberUpdateItems_;i++)
4463      updateItems_[i] = rhs.updateItems_[i];
4464  }
4465#endif
4466  if (maximumWhich_&&rhs.whichGenerator_)
4467    whichGenerator_ = CoinCopyOfArray(rhs.whichGenerator_,maximumWhich_);
4468  nodeCompare_=rhs.nodeCompare_->clone();
4469  problemFeasibility_=rhs.problemFeasibility_->clone();
4470  tree_= rhs.tree_->clone();
4471  if (rhs.branchingMethod_)
4472    branchingMethod_=rhs.branchingMethod_->clone();
4473  else
4474    branchingMethod_=NULL;
4475  if (rhs.cutModifier_)
4476    cutModifier_=rhs.cutModifier_->clone();
4477  else
4478    cutModifier_=NULL;
4479  cbcColLower_ = NULL;
4480  cbcColUpper_ = NULL;
4481  cbcRowLower_ = NULL;
4482  cbcRowUpper_ = NULL;
4483  cbcColSolution_ = NULL;
4484  cbcRowPrice_ = NULL;
4485  cbcReducedCost_ = NULL;
4486  cbcRowActivity_ = NULL;
4487  if (rhs.strategy_)
4488    strategy_=rhs.strategy_->clone();
4489  else
4490    strategy_=NULL;
4491  parentModel_=rhs.parentModel_;
4492  appData_=rhs.appData_;
4493  messages_ = rhs.messages_;
4494  ownership_ = 0x80000000;
4495  messageHandler()->setLogLevel(rhs.messageHandler()->logLevel());
4496  numberIntegers_=rhs.numberIntegers_;
4497  randomNumberGenerator_ = rhs.randomNumberGenerator_;
4498  if (numberIntegers_) {
4499    integerVariable_ = new int [numberIntegers_];
4500    memcpy(integerVariable_,rhs.integerVariable_,numberIntegers_*sizeof(int));
4501    integerInfo_ = CoinCopyOfArray(rhs.integerInfo_,solver_->getNumCols());
4502  } else {
4503    integerVariable_ = NULL;
4504    integerInfo_=NULL;
4505  }
4506  if (rhs.hotstartSolution_) {
4507    int numberColumns = solver_->getNumCols();
4508    hotstartSolution_ = CoinCopyOfArray(rhs.hotstartSolution_,numberColumns);
4509    hotstartPriorities_ = CoinCopyOfArray(rhs.hotstartPriorities_,numberColumns);
4510  } else {
4511    hotstartSolution_ = NULL;
4512    hotstartPriorities_ =NULL;
4513  }
4514  if (rhs.bestSolution_&&!noTree) {
4515    int numberColumns = solver_->getNumCols();
4516    bestSolution_ = new double[numberColumns];
4517    memcpy(bestSolution_,rhs.bestSolution_,numberColumns*sizeof(double));
4518  } else {
4519    bestSolution_=NULL;
4520  }
4521  if (!noTree) {
4522    int numberColumns = solver_->getNumCols();
4523    // Space for current solution
4524    currentSolution_ = new double[numberColumns];
4525    continuousSolution_ = new double[numberColumns];
4526    usedInSolution_ = new int[numberColumns];
4527    CoinZeroN(usedInSolution_,numberColumns);
4528  } else {
4529    currentSolution_=NULL;
4530    continuousSolution_=NULL;
4531    usedInSolution_=NULL;
4532  }
4533  testSolution_=currentSolution_;
4534  numberRowsAtContinuous_ = rhs.numberRowsAtContinuous_;
4535  maximumNumberCuts_=rhs.maximumNumberCuts_;
4536  phase_ = rhs.phase_;
4537  currentNumberCuts_=rhs.currentNumberCuts_;
4538  maximumDepth_= rhs.maximumDepth_;
4539  if (noTree) {
4540    bestObjective_ = COIN_DBL_MAX;
4541    numberSolutions_ =0;
4542    stateOfSearch_= 0;
4543    numberHeuristicSolutions_=0;
4544    numberNodes_=0;
4545    numberNodes2_=0;
4546    numberIterations_=0;
4547    status_=0;
4548    subTreeModel_=NULL;
4549    numberStoppedSubTrees_=0;
4550    continuousObjective_=COIN_DBL_MAX;
4551    originalContinuousObjective_=COIN_DBL_MAX;
4552    continuousInfeasibilities_=COIN_INT_MAX;
4553    maximumNumberCuts_=0;
4554    tree_->cleanTree(this,-COIN_DBL_MAX,bestPossibleObjective_);
4555    bestPossibleObjective_ = COIN_DBL_MAX;
4556  }
4557  // These are only used as temporary arrays so need not be filled
4558  if (maximumNumberCuts_) {
4559    addedCuts_ = new CbcCountRowCut * [maximumNumberCuts_];
4560  } else {
4561    addedCuts_ = NULL;
4562  }
4563  bestSolutionBasis_ = rhs.bestSolutionBasis_;
4564  nextRowCut_ = NULL;
4565  currentNode_ = NULL;
4566  if (maximumDepth_)
4567    walkback_ = new CbcNodeInfo * [maximumDepth_];
4568  else
4569    walkback_ = NULL;
4570  synchronizeModel();
4571}
4572 
4573// Assignment operator
4574CbcModel & 
4575CbcModel::operator=(const CbcModel& rhs)
4576{
4577  if (this!=&rhs) {
4578    if (modelOwnsSolver()) {
4579      delete solver_;
4580      solver_=NULL;
4581    }
4582    gutsOfDestructor();
4583    if (defaultHandler_)
4584    { delete handler_;
4585      handler_ = NULL; }
4586    defaultHandler_ = rhs.defaultHandler_;
4587    if (defaultHandler_)
4588    { handler_ = new CoinMessageHandler();
4589      handler_->setLogLevel(2); }
4590    else
4591    { handler_ = rhs.handler_; }
4592    messages_ = rhs.messages_;
4593    messageHandler()->setLogLevel(rhs.messageHandler()->logLevel());
4594    if (rhs.solver_)
4595    { solver_ = rhs.solver_->clone() ; }
4596    else
4597    { solver_ = 0 ; }
4598    ownership_ = 0x80000000;
4599    delete continuousSolver_ ;
4600    if (rhs.continuousSolver_)
4601    { continuousSolver_ = rhs.continuousSolver_->clone() ; }
4602    else
4603    { continuousSolver_ = 0 ; }
4604    delete referenceSolver_;
4605    if (rhs.referenceSolver_)
4606    { referenceSolver_ = rhs.referenceSolver_->clone() ; }
4607    else
4608    { referenceSolver_ = NULL ; }
4609
4610    delete emptyWarmStart_ ;
4611    if (rhs.emptyWarmStart_)
4612    { emptyWarmStart_ = rhs.emptyWarmStart_->clone() ; }
4613    else
4614    { emptyWarmStart_ = 0 ; }
4615
4616    bestObjective_ = rhs.bestObjective_;
4617    bestPossibleObjective_=rhs.bestPossibleObjective_;
4618    sumChangeObjective1_=rhs.sumChangeObjective1_;
4619    sumChangeObjective2_=rhs.sumChangeObjective2_;
4620    delete [] bestSolution_;
4621    if (rhs.bestSolution_) {
4622      int numberColumns = rhs.getNumCols();
4623      bestSolution_ = new double[numberColumns];
4624      memcpy(bestSolution_,rhs.bestSolution_,numberColumns*sizeof(double));
4625    } else {
4626      bestSolution_=NULL;
4627    }
4628    int numberColumns = rhs.getNumCols();
4629    if (numberColumns) {
4630      // Space for current solution
4631      currentSolution_ = new double[numberColumns];
4632      continuousSolution_ = new double[numberColumns];
4633      usedInSolution_ = new int[numberColumns];
4634      CoinZeroN(usedInSolution_,numberColumns);
4635    } else {
4636      currentSolution_=NULL;
4637      continuousSolution_=NULL;
4638      usedInSolution_=NULL;
4639    }
4640    testSolution_=currentSolution_;
4641    minimumDrop_ = rhs.minimumDrop_;
4642    numberSolutions_=rhs.numberSolutions_;
4643    stateOfSearch_= rhs.stateOfSearch_;
4644    numberHeuristicSolutions_=rhs.numberHeuristicSolutions_;
4645    numberNodes_ = rhs.numberNodes_;
4646    numberNodes2_ = rhs.numberNodes2_;
4647    numberIterations_ = rhs.numberIterations_;
4648    status_ = rhs.status_;
4649    secondaryStatus_ = rhs.secondaryStatus_;
4650    specialOptions_ = rhs.specialOptions_;
4651    subTreeModel_ = rhs.subTreeModel_;
4652    numberStoppedSubTrees_ = rhs.numberStoppedSubTrees_;
4653    mutex_ = NULL;
4654    presolve_ = rhs.presolve_;
4655    numberStrong_ = rhs.numberStrong_;
4656    numberBeforeTrust_ = rhs.numberBeforeTrust_;
4657    numberPenalties_ = rhs.numberPenalties_;
4658    stopNumberIterations_ = rhs.stopNumberIterations_;
4659    penaltyScaleFactor_ = rhs.penaltyScaleFactor_;
4660    numberAnalyzeIterations_ = rhs.numberAnalyzeIterations_;
4661    delete [] analyzeResults_;
4662    analyzeResults_ = NULL;
4663    numberInfeasibleNodes_ = rhs.numberInfeasibleNodes_;
4664    problemType_ = rhs.problemType_;
4665    printFrequency_ = rhs.printFrequency_;
4666    howOftenGlobalScan_=rhs.howOftenGlobalScan_;
4667    numberGlobalViolations_=rhs.numberGlobalViolations_;
4668    numberExtraIterations_ = rhs.numberExtraIterations_;
4669    numberExtraNodes_ = rhs.numberExtraNodes_;
4670    continuousObjective_=rhs.continuousObjective_;
4671    originalContinuousObjective_ = rhs.originalContinuousObjective_;
4672    continuousInfeasibilities_ = rhs.continuousInfeasibilities_;
4673    maximumCutPassesAtRoot_ = rhs.maximumCutPassesAtRoot_;
4674    maximumCutPasses_ = rhs.maximumCutPasses_;
4675    preferredWay_ = rhs.preferredWay_;
4676    currentPassNumber_ = rhs.currentPassNumber_;
4677    memcpy(intParam_,rhs.intParam_,sizeof(intParam_));
4678    memcpy(dblParam_,rhs.dblParam_,sizeof(dblParam_));
4679    globalCuts_ = rhs.globalCuts_;
4680    int i;
4681    for (i=0;i<numberCutGenerators_;i++) {
4682      delete generator_[i];
4683      delete virginGenerator_[i];
4684    }
4685    delete [] generator_;
4686    delete [] virginGenerator_;
4687    delete [] heuristic_;
4688    maximumWhich_ = rhs.maximumWhich_;
4689    delete [] whichGenerator_;
4690    whichGenerator_ = NULL;
4691    if (maximumWhich_&&rhs.whichGenerator_)
4692      whichGenerator_ = CoinCopyOfArray(rhs.whichGenerator_,maximumWhich_);
4693    maximumRows_=0;
4694    currentDepth_ = 0;
4695    randomNumberGenerator_ = rhs.randomNumberGenerator_;
4696    workingBasis_ = CoinWarmStartBasis();
4697    for (i=0;i<maximumStatistics_;i++)
4698      delete statistics_[i];
4699    delete [] statistics_;
4700    maximumStatistics_ = 0;
4701    statistics_ = NULL;
4702    delete probingInfo_;
4703    probingInfo_=NULL;
4704    numberFixedAtRoot_ = rhs.numberFixedAtRoot_;
4705    numberFixedNow_ = rhs.numberFixedNow_;
4706    stoppedOnGap_ = rhs.stoppedOnGap_;
4707    eventHappened_ = rhs.eventHappened_;
4708    numberLongStrong_ = rhs.numberLongStrong_;
4709    numberOldActiveCuts_ = rhs.numberOldActiveCuts_;
4710    numberNewCuts_ = rhs.numberNewCuts_;
4711    resolveAfterTakeOffCuts_=rhs.resolveAfterTakeOffCuts_;
4712    maximumNumberIterations_ = rhs.maximumNumberIterations_;
4713#if NEW_UPDATE_OBJECT>1
4714    numberUpdateItems_ = rhs.numberUpdateItems_;
4715    maximumNumberUpdateItems_ = rhs.maximumNumberUpdateItems_;
4716    delete [] updateItems_;
4717    if (maximumNumberUpdateItems_) {
4718      updateItems_ = new CbcObjectUpdateData [maximumNumberUpdateItems_];
4719      for (i=0;i<maximumNumberUpdateItems_;i++)
4720        updateItems_[i] = rhs.updateItems_[i];
4721    } else {
4722      updateItems_ = NULL;
4723    }
4724#endif
4725    numberThreads_ = rhs.numberThreads_;
4726    threadMode_ = rhs.threadMode_;
4727    sizeMiniTree_ = rhs.sizeMiniTree_;
4728    searchStrategy_ = rhs.searchStrategy_;
4729    numberStrongIterations_ = rhs.numberStrongIterations_;
4730    strongInfo_[0]=rhs.strongInfo_[0];
4731    strongInfo_[1]=rhs.strongInfo_[1];
4732    strongInfo_[2]=rhs.strongInfo_[2];
4733    solverCharacteristics_ = NULL;
4734    lastHeuristic_ = NULL;
4735    numberCutGenerators_ = rhs.numberCutGenerators_;
4736    if (numberCutGenerators_) {
4737      generator_ = new CbcCutGenerator * [numberCutGenerators_];
4738      virginGenerator_ = new CbcCutGenerator * [numberCutGenerators_];
4739      int i;
4740      for (i=0;i<numberCutGenerators_;i++) {
4741        generator_[i]=new CbcCutGenerator(*rhs.generator_[i]);
4742        virginGenerator_[i]=new CbcCutGenerator(*rhs.virginGenerator_[i]);
4743      }
4744    } else {
4745      generator_=NULL;
4746      virginGenerator_=NULL;
4747    }
4748    numberHeuristics_ = rhs.numberHeuristics_;
4749    if (numberHeuristics_) {
4750      heuristic_ = new CbcHeuristic * [numberHeuristics_];
4751      memcpy(heuristic_,rhs.heuristic_,
4752             numberHeuristics_*sizeof(CbcHeuristic *));
4753    } else {
4754      heuristic_=NULL;
4755    }
4756    lastHeuristic_ = NULL;
4757    if (eventHandler_)
4758      delete eventHandler_ ;
4759    if (rhs.eventHandler_)
4760      { eventHandler_ = rhs.eventHandler_->clone() ; }
4761    else
4762    { eventHandler_ = NULL ; }
4763# ifdef COIN_HAS_CLP
4764    fastNodeDepth_ = rhs.fastNodeDepth_;
4765#endif
4766    if (ownObjects_) {
4767      for (i=0;i<numberObjects_;i++)
4768        delete object_[i];
4769      delete [] object_;
4770      numberObjects_=rhs.numberObjects_;
4771      if (numberObjects_) {
4772        object_ = new OsiObject * [numberObjects_];
4773        int i;
4774        for (i=0;i<numberObjects_;i++) 
4775          object_[i]=(rhs.object_[i])->clone();
4776      } else {
4777        object_=NULL;
4778    }
4779    } else {
4780      // assume will be redone
4781      numberObjects_=0;
4782      object_=NULL;
4783    }
4784    delete [] originalColumns_;
4785    if (rhs.originalColumns_) {
4786      int numberColumns = rhs.getNumCols();
4787      originalColumns_= new int [numberColumns];
4788      memcpy(originalColumns_,rhs.originalColumns_,numberColumns*sizeof(int));
4789    } else {
4790      originalColumns_=NULL;
4791    }
4792    nodeCompare_=rhs.nodeCompare_->clone();
4793    problemFeasibility_=rhs.problemFeasibility_->clone();
4794    delete tree_;
4795    tree_= rhs.tree_->clone();
4796    if (rhs.branchingMethod_)
4797      branchingMethod_=rhs.branchingMethod_->clone();
4798    else
4799      branchingMethod_=NULL;
4800    if (rhs.cutModifier_)
4801      cutModifier_=rhs.cutModifier_->clone();
4802    else
4803      cutModifier_=NULL;
4804    delete strategy_;
4805    if (rhs.strategy_)
4806      strategy_=rhs.strategy_->clone();
4807    else
4808      strategy_=NULL;
4809    parentModel_=rhs.parentModel_;
4810    appData_=rhs.appData_;
4811
4812    delete [] integerVariable_;
4813    numberIntegers_=rhs.numberIntegers_;
4814    if (numberIntegers_) {
4815      integerVariable_ = new int [numberIntegers_];
4816      memcpy(integerVariable_,rhs.integerVariable_,
4817             numberIntegers_*sizeof(int));
4818      integerInfo_ = CoinCopyOfArray(rhs.integerInfo_,rhs.getNumCols());
4819    } else {
4820      integerVariable_ = NULL;
4821      integerInfo_=NULL;
4822    }
4823    if (rhs.hotstartSolution_) {
4824      int numberColumns = solver_->getNumCols();
4825      hotstartSolution_ = CoinCopyOfArray(rhs.hotstartSolution_,numberColumns);
4826      hotstartPriorities_ = CoinCopyOfArray(rhs.hotstartPriorities_,numberColumns);
4827    } else {
4828      hotstartSolution_ = NULL;
4829      hotstartPriorities_ =NULL;
4830    }
4831    numberRowsAtContinuous_ = rhs.numberRowsAtContinuous_;
4832    maximumNumberCuts_=rhs.maximumNumberCuts_;
4833    phase_ = rhs.phase_;
4834    currentNumberCuts_=rhs.currentNumberCuts_;
4835    maximumDepth_= rhs.maximumDepth_;
4836    delete [] addedCuts_;
4837    delete [] walkback_;
4838    // These are only used as temporary arrays so need not be filled
4839    if (maximumNumberCuts_) {
4840      addedCuts_ = new CbcCountRowCut * [maximumNumberCuts_];
4841    } else {
4842      addedCuts_ = NULL;
4843    }
4844    bestSolutionBasis_ = rhs.bestSolutionBasis_;
4845    nextRowCut_ = NULL;
4846    currentNode_ = NULL;
4847    if (maximumDepth_)
4848      walkback_ = new CbcNodeInfo * [maximumDepth_];
4849    else
4850      walkback_ = NULL;
4851    synchronizeModel();
4852    cbcColLower_ = NULL;
4853    cbcColUpper_ = NULL;
4854    cbcRowLower_ = NULL;
4855    cbcRowUpper_ = NULL;
4856    cbcColSolution_ = NULL;
4857    cbcRowPrice_ = NULL;
4858    cbcReducedCost_ = NULL;
4859    cbcRowActivity_ = NULL;
4860  }
4861  return *this;
4862}
4863// Destructor
4864CbcModel::~CbcModel ()
4865{
4866  if (defaultHandler_) {
4867    delete handler_;
4868    handler_ = NULL;
4869  }
4870  delete tree_;
4871  tree_=NULL;
4872  if (modelOwnsSolver()) {
4873    delete solver_;
4874    solver_ = NULL;
4875  }
4876  gutsOfDestructor();
4877  delete eventHandler_ ;
4878  eventHandler_ = NULL ;
4879}
4880// Clears out as much as possible (except solver)
4881void 
4882CbcModel::gutsOfDestructor()
4883{
4884  delete referenceSolver_;
4885  referenceSolver_=NULL;
4886  int i;
4887  for (i=0;i<numberCutGenerators_;i++) {
4888    delete generator_[i];
4889    delete virginGenerator_[i];
4890  }
4891  delete [] generator_;
4892  delete [] virginGenerator_;
4893  generator_=NULL;
4894  virginGenerator_=NULL;
4895  for (i=0;i<numberHeuristics_;i++)
4896    delete heuristic_[i];
4897  delete [] heuristic_;
4898  heuristic_=NULL;
4899  delete nodeCompare_;
4900  nodeCompare_=NULL;
4901  delete problemFeasibility_;
4902  problemFeasibility_=NULL;
4903  delete [] originalColumns_;
4904  originalColumns_=NULL;
4905  delete strategy_;
4906#if NEW_UPDATE_OBJECT>1
4907  delete [] updateItems_;
4908  updateItems_=NULL;
4909  numberUpdateItems_=0;
4910  maximumNumberUpdateItems_=0;
4911#endif
4912  gutsOfDestructor2();
4913}
4914// Clears out enough to reset CbcModel
4915void 
4916CbcModel::gutsOfDestructor2()
4917{
4918  delete [] integerInfo_;
4919  integerInfo_=NULL;
4920  delete [] integerVariable_;
4921  integerVariable_=NULL;
4922  int i;
4923  if (ownObjects_) {
4924    for (i=0;i<numberObjects_;i++)
4925      delete object_[i];
4926    delete [] object_;
4927  }
4928  ownObjects_=true;
4929  object_=NULL;
4930  numberIntegers_=0;
4931  numberObjects_=0;
4932  // Below here is whatever consensus is
4933  ownership_ = 0x80000000;
4934  delete branchingMethod_;
4935  branchingMethod_=NULL;
4936  delete cutModifier_;
4937  cutModifier_=NULL;
4938  resetModel();
4939}
4940// Clears out enough to reset CbcModel
4941void 
4942CbcModel::resetModel()
4943{
4944  delete emptyWarmStart_ ;
4945  emptyWarmStart_ =NULL;
4946  delete continuousSolver_;
4947  continuousSolver_=NULL;
4948  delete [] bestSolution_;
4949  bestSolution_=NULL;
4950  delete [] currentSolution_;
4951  currentSolution_=NULL;
4952  delete [] continuousSolution_;
4953  continuousSolution_=NULL;
4954  solverCharacteristics_=NULL;
4955  delete [] usedInSolution_;
4956  usedInSolution_ = NULL;
4957  testSolution_=NULL;
4958  lastHeuristic_ = NULL;
4959  delete [] addedCuts_;
4960  addedCuts_=NULL;
4961  nextRowCut_ = NULL;
4962  currentNode_ = NULL;
4963  delete [] walkback_;
4964  walkback_=NULL;
4965  delete [] whichGenerator_;
4966  whichGenerator_ = NULL;
4967  for (int i=0;i<maximumStatistics_;i++)
4968    delete statistics_[i];
4969  delete [] statistics_;
4970  statistics_=NULL;
4971  maximumDepthActual_ = 0;
4972  numberDJFixed_ =0.0;
4973  delete probingInfo_;
4974  probingInfo_ = NULL;
4975  maximumStatistics_=0;
4976  delete [] analyzeResults_;
4977  analyzeResults_=NULL;
4978  bestObjective_=COIN_DBL_MAX;
4979  bestPossibleObjective_=COIN_DBL_MAX;
4980  sumChangeObjective1_=0.0;
4981  sumChangeObjective2_=0.0;
4982  numberSolutions_=0;
4983  stateOfSearch_=0;
4984  delete [] hotstartSolution_;
4985  hotstartSolution_=NULL;
4986  delete [] hotstartPriorities_;
4987  hotstartPriorities_=NULL;
4988  numberHeuristicSolutions_=0;
4989  numberNodes_=0;
4990  numberNodes2_=0;
4991  numberIterations_=0;
4992  status_=-1;
4993  secondaryStatus_=-1;
4994  maximumNumberCuts_=0;
4995  phase_=0;
4996  currentNumberCuts_=0;
4997  maximumDepth_=0;
4998  nextRowCut_=NULL;
4999  currentNode_=NULL;
5000  // clear out tree
5001  if (tree_&&tree_->size())
5002    tree_->cleanTree(this, -1.0e100,bestPossibleObjective_) ;
5003  subTreeModel_=NULL;
5004  numberStoppedSubTrees_=0;
5005  numberInfeasibleNodes_=0;
5006  numberGlobalViolations_=0;
5007  numberExtraIterations_ = 0;
5008  numberExtraNodes_ = 0;
5009  continuousObjective_=0.0;
5010  originalContinuousObjective_=0.0;
5011  continuousInfeasibilities_=0;
5012  numberFixedAtRoot_=0;
5013  numberFixedNow_=0;
5014  stoppedOnGap_=false;
5015  eventHappened_=false;
5016  numberLongStrong_=0;
5017  numberOldActiveCuts_=0;
5018  numberNewCuts_=0;
5019  searchStrategy_=-1;
5020  numberStrongIterations_=0;
5021  // Parameters which need to be reset
5022  setCutoff(COIN_DBL_MAX);
5023  dblParam_[CbcCutoffIncrement] = 1e-5;
5024  dblParam_[CbcCurrentCutoff] = 1.0e100;
5025  dblParam_[CbcCurrentObjectiveValue] = 1.0e100;
5026  dblParam_[CbcCurrentMinimizationObjectiveValue] = 1.0e100;
5027}
5028/* Most of copy constructor
5029      mode - 0 copy but don't delete before
5030             1 copy and delete before
5031             2 copy and delete before (but use virgin generators)
5032*/
5033void 
5034CbcModel::gutsOfCopy(const CbcModel & rhs,int mode)
5035{
5036  minimumDrop_ = rhs.minimumDrop_;
5037  specialOptions_ = rhs.specialOptions_;
5038  numberStrong_ = rhs.numberStrong_;
5039  numberBeforeTrust_ = rhs.numberBeforeTrust_;
5040  numberPenalties_ = rhs.numberPenalties_;
5041  printFrequency_ = rhs.printFrequency_;
5042# ifdef COIN_HAS_CLP
5043  fastNodeDepth_ = rhs.fastNodeDepth_;
5044#endif
5045  howOftenGlobalScan_ = rhs.howOftenGlobalScan_;
5046  maximumCutPassesAtRoot_ = rhs.maximumCutPassesAtRoot_;
5047  maximumCutPasses_ =  rhs.maximumCutPasses_;
5048  preferredWay_ = rhs.preferredWay_;
5049  resolveAfterTakeOffCuts_ = rhs.resolveAfterTakeOffCuts_;
5050  maximumNumberIterations_ = rhs.maximumNumberIterations_;
5051  numberThreads_ = rhs.numberThreads_;
5052  threadMode_ = rhs.threadMode_;
5053  memcpy(intParam_,rhs.intParam_,sizeof(intParam_));
5054  memcpy(dblParam_,rhs.dblParam_,sizeof(dblParam_));
5055  int i;
5056  if (mode) {
5057    for (i=0;i<numberCutGenerators_;i++) {
5058      delete generator_[i];
5059      delete virginGenerator_[i];
5060    }
5061    delete [] generator_;
5062    delete [] virginGenerator_;
5063    for (i=0;i<numberHeuristics_;i++) {
5064      delete heuristic_[i];
5065    }
5066    delete [] heuristic_;
5067    delete eventHandler_;
5068    delete branchingMethod_;
5069  }
5070  numberCutGenerators_ = rhs.numberCutGenerators_;
5071  if (numberCutGenerators_) {
5072    generator_ = new CbcCutGenerator * [numberCutGenerators_];
5073    virginGenerator_ = new CbcCutGenerator * [numberCutGenerators_];
5074    int i;
5075    for (i=0;i<numberCutGenerators_;i++) {
5076      if (mode<2)
5077        generator_[i]=new CbcCutGenerator(*rhs.generator_[i]);
5078      else
5079        generator_[i]=new CbcCutGenerator(*rhs.virginGenerator_[i]);
5080      virginGenerator_[i]=new CbcCutGenerator(*rhs.virginGenerator_[i]);
5081    }
5082  } else {
5083    generator_=NULL;
5084    virginGenerator_=NULL;
5085  }
5086  numberHeuristics_ = rhs.numberHeuristics_;
5087  if (numberHeuristics_) {
5088    heuristic_ = new CbcHeuristic * [numberHeuristics_];
5089    int i;
5090    for (i=0;i<numberHeuristics_;i++) {
5091      heuristic_[i]=rhs.heuristic_[i]->clone();
5092    }
5093  } else {
5094    heuristic_=NULL;
5095  }
5096  if (rhs.eventHandler_)
5097    eventHandler_ = rhs.eventHandler_->clone() ; 
5098  else
5099    eventHandler_ = NULL ; 
5100  if (rhs.branchingMethod_)
5101    branchingMethod_=rhs.branchingMethod_->clone();
5102  else
5103    branchingMethod_=NULL;
5104  messageHandler()->setLogLevel(rhs.messageHandler()->logLevel());
5105  synchronizeModel();
5106}
5107// Move status, nodes etc etc across
5108void 
5109CbcModel::moveInfo(const CbcModel & rhs)
5110{
5111  bestObjective_ = rhs.bestObjective_;
5112  bestPossibleObjective_=rhs.bestPossibleObjective_;
5113  numberSolutions_=rhs.numberSolutions_;
5114  numberHeuristicSolutions_=rhs.numberHeuristicSolutions_;
5115  numberNodes_ = rhs.numberNodes_;
5116  numberNodes2_ = rhs.numberNodes2_;
5117  numberIterations_ = rhs.numberIterations_;
5118  status_ = rhs.status_;
5119  secondaryStatus_ = rhs.secondaryStatus_;
5120  numberStoppedSubTrees_ = rhs.numberStoppedSubTrees_;
5121  numberInfeasibleNodes_ = rhs.numberInfeasibleNodes_;
5122  continuousObjective_=rhs.continuousObjective_;
5123  originalContinuousObjective_ = rhs.originalContinuousObjective_;
5124  continuousInfeasibilities_ = rhs.continuousInfeasibilities_;
5125  numberFixedAtRoot_ = rhs.numberFixedAtRoot_;
5126  numberFixedNow_ = rhs.numberFixedNow_;
5127  stoppedOnGap_ = rhs.stoppedOnGap_;
5128  eventHappened_ = rhs.eventHappened_;
5129  numberLongStrong_ = rhs.numberLongStrong_;
5130  numberStrongIterations_ = rhs.numberStrongIterations_;
5131  strongInfo_[0]=rhs.strongInfo_[0];
5132  strongInfo_[1]=rhs.strongInfo_[1];
5133  strongInfo_[2]=rhs.strongInfo_[2];
5134  numberRowsAtContinuous_ = rhs.numberRowsAtContinuous_;
5135  maximumDepth_= rhs.maximumDepth_;
5136}
5137// Save a copy of the current solver so can be reset to
5138void 
5139CbcModel::saveReferenceSolver()
5140{
5141  delete referenceSolver_;
5142  referenceSolver_= solver_->clone();
5143}
5144
5145// Uses a copy of reference solver to be current solver
5146void 
5147CbcModel::resetToReferenceSolver()
5148{
5149  delete solver_;
5150  solver_ = referenceSolver_->clone();
5151  // clear many things
5152  gutsOfDestructor2();
5153  // Reset cutoff
5154  // Solvers know about direction
5155  double direction = solver_->getObjSense();
5156  double value;
5157  solver_->getDblParam(OsiDualObjectiveLimit,value); 
5158  setCutoff(value*direction);
5159}
5160
5161// Are there a numerical difficulties?
5162bool 
5163CbcModel::isAbandoned() const
5164{
5165  return status_ == 2;
5166}
5167// Is optimality proven?
5168bool 
5169CbcModel::isProvenOptimal() const
5170{
5171  if (!status_ && bestObjective_<1.0e30)
5172    return true;
5173  else
5174    return false;
5175}
5176// Is  infeasiblity proven (or none better than cutoff)?
5177bool 
5178CbcModel::isProvenInfeasible() const
5179{
5180  if (!status_ && bestObjective_>=1.0e30)
5181    return true;
5182  else
5183    return false;
5184}
5185// Was continuous solution unbounded
5186bool 
5187CbcModel::isContinuousUnbounded() const
5188{
5189  if (!status_ && secondaryStatus_==7)
5190    return true;
5191  else
5192    return false;
5193}
5194// Was continuous solution unbounded
5195bool 
5196CbcModel::isProvenDualInfeasible() const
5197{
5198  if (!status_ && secondaryStatus_==7)
5199    return true;
5200  else
5201    return false;
5202}
5203// Node limit reached?
5204bool 
5205CbcModel::isNodeLimitReached() const
5206{
5207  return numberNodes_ >= intParam_[CbcMaxNumNode];
5208}
5209// Time limit reached?
5210bool 
5211CbcModel::isSecondsLimitReached() const
5212{
5213  if (status_==1&&secondaryStatus_==4)
5214    return true;
5215  else
5216    return false;
5217}
5218// Solution limit reached?
5219bool 
5220CbcModel::isSolutionLimitReached() const
5221{
5222  return numberSolutions_ >= intParam_[CbcMaxNumSol];
5223}
5224// Set language
5225void 
5226CbcModel::newLanguage(CoinMessages::Language language)
5227{
5228  messages_ = CbcMessage(language);
5229}
5230void 
5231CbcModel::setNumberStrong(int number)
5232{
5233  if (number<0)
5234    numberStrong_=0;
5235   else
5236    numberStrong_=number;
5237}
5238void 
5239CbcModel::setNumberBeforeTrust(int number)
5240{
5241  if (number<-3) {
5242    numberBeforeTrust_=0;
5243  } else {
5244    numberBeforeTrust_=number;
5245    //numberStrong_ = CoinMax(numberStrong_,1);
5246  }
5247}
5248void 
5249CbcModel::setNumberPenalties(int number)
5250{
5251  if (number<=0) {
5252    numberPenalties_=0;
5253  } else {
5254    numberPenalties_=number;
5255  }
5256}
5257void 
5258CbcModel::setPenaltyScaleFactor(double value)
5259{
5260  if (value<=0) {
5261    penaltyScaleFactor_=3.0;
5262  } else {
5263    penaltyScaleFactor_=value;
5264  }
5265}
5266void 
5267CbcModel::setHowOftenGlobalScan(int number)
5268{
5269  if (number<-1)
5270    howOftenGlobalScan_=0;
5271   else
5272    howOftenGlobalScan_=number;
5273}
5274
5275// Add one generator
5276void 
5277CbcModel::addCutGenerator(CglCutGenerator * generator,
5278                          int howOften, const char * name,
5279                          bool normal, bool atSolution,
5280                          bool whenInfeasible,int howOftenInSub,
5281                          int whatDepth, int whatDepthInSub)
5282{
5283  CbcCutGenerator ** temp = generator_;
5284  generator_ = new CbcCutGenerator * [numberCutGenerators_+1];
5285  memcpy(generator_,temp,numberCutGenerators_*sizeof(CbcCutGenerator *));
5286  delete[] temp ;
5287  generator_[numberCutGenerators_]= 
5288    new CbcCutGenerator(this,generator, howOften, name,
5289                        normal,atSolution,whenInfeasible,howOftenInSub,
5290                        whatDepth, whatDepthInSub);
5291  // and before any changes
5292  temp = virginGenerator_;
5293  virginGenerator_ = new CbcCutGenerator * [numberCutGenerators_+1];
5294  memcpy(virginGenerator_,temp,numberCutGenerators_*sizeof(CbcCutGenerator *));
5295  delete[] temp ;
5296  virginGenerator_[numberCutGenerators_++]= 
5297    new CbcCutGenerator(this,generator, howOften, name,
5298                        normal,atSolution,whenInfeasible,howOftenInSub,
5299                        whatDepth, whatDepthInSub);
5300                                                         
5301}
5302// Add one heuristic
5303void 
5304CbcModel::addHeuristic(CbcHeuristic * generator, const char *name,
5305                       int before)
5306{
5307  CbcHeuristic ** temp = heuristic_;
5308  heuristic_ = new CbcHeuristic * [numberHeuristics_+1];
5309  memcpy(heuristic_,temp,numberHeuristics_*sizeof(CbcHeuristic *));
5310  delete [] temp;
5311  int where;
5312  if (before<0||before>=numberHeuristics_) {
5313    where=numberHeuristics_;
5314  } else {
5315    // move up
5316    for (int i=numberHeuristics_;i>before;i--) 
5317      heuristic_[i]=heuristic_[i-1];
5318    where=before;
5319  }
5320  heuristic_[where]=generator->clone();
5321  if (name)
5322    heuristic_[where]->setHeuristicName(name) ; 
5323  heuristic_[where]->setSeed(987654321+where);
5324  numberHeuristics_++ ;
5325}
5326
5327/*
5328  The last subproblem handled by the solver is not necessarily related to the
5329  one being recreated, so the first action is to remove all cuts from the
5330  constraint system.  Next, traverse the tree from node to the root to
5331  determine the basis size required for this subproblem and create an empty
5332  basis with the right capacity.  Finally, traverse the tree from root to
5333  node, adjusting bounds in the constraint system, adjusting the basis, and
5334  collecting the cuts that must be added to the constraint system.
5335  applyToModel does the heavy lifting.
5336
5337  addCuts1 is used in contexts where all that's desired is the list of cuts:
5338  the node is already fathomed, and we're collecting cuts so that we can
5339  adjust reference counts as we prune nodes. Arguably the two functions
5340  should be separated. The culprit is applyToModel, which performs cut
5341  collection and model adjustment.
5342
5343  Certainly in the contexts where all we need is a list of cuts, there's no
5344  point in passing in a valid basis --- an empty basis will do just fine.
5345*/
5346void CbcModel::addCuts1 (CbcNode * node, CoinWarmStartBasis *&lastws)
5347{ int i;
5348  int nNode=0;
5349  int numberColumns = getNumCols();
5350  CbcNodeInfo * nodeInfo = node->nodeInfo();
5351
5352/*
5353  Remove all cuts from the constraint system.
5354  (original comment includes ``see note below for later efficiency'', but
5355  the reference isn't clear to me).
5356*/
5357  int currentNumberCuts ;
5358  if (numberThreads_<=0) {
5359    solver_->restoreBaseModel(numberRowsAtContinuous_);
5360  } else {
5361    // *** Fix later
5362    currentNumberCuts = solver_->getNumRows()-numberRowsAtContinuous_;
5363    int *which = new int[currentNumberCuts];
5364    for (i = 0 ; i < currentNumberCuts ; i++)
5365      which[i] = i+numberRowsAtContinuous_;
5366    solver_->deleteRows(currentNumberCuts,which);
5367    delete [] which;
5368  }
5369/*
5370  Accumulate the path from node to the root in walkback_, and accumulate a
5371  cut count in currentNumberCuts.
5372
5373  original comment: when working then just unwind until where new node joins
5374  old node (for cuts?)
5375*/
5376  currentNumberCuts = 0;
5377  while (nodeInfo) {
5378    //printf("nNode = %d, nodeInfo = %x\n",nNode,nodeInfo);
5379    walkback_[nNode++]=nodeInfo;
5380    currentNumberCuts += nodeInfo->numberCuts() ;
5381    nodeInfo = nodeInfo->parent() ;
5382    if (nNode==maximumDepth_) {
5383      maximumDepth_ *= 2;
5384      CbcNodeInfo ** temp = new CbcNodeInfo * [maximumDepth_];
5385      for (i=0;i<nNode;i++) 
5386        temp[i] = walkback_[i];
5387      delete [] walkback_;
5388      walkback_ = temp;
5389    }
5390  }
5391/*
5392  Create an empty basis with sufficient capacity for the constraint system
5393  we'll construct: original system plus cuts. Make sure we have capacity to
5394  record those cuts in addedCuts_.
5395
5396  The method of adjusting the basis at a FullNodeInfo object (the root, for
5397  example) is to use a copy constructor to duplicate the basis held in the
5398  nodeInfo, then resize it and return the new basis object. Guaranteed,
5399  lastws will point to a different basis when it returns. We pass in a basis
5400  because we need the parameter to return the allocated basis, and it's an
5401  easy way to pass in the size. But we take a hit for memory allocation.
5402*/
5403  currentNumberCuts_=currentNumberCuts;
5404  if (currentNumberCuts > maximumNumberCuts_) {
5405    maximumNumberCuts_ = currentNumberCuts;
5406    delete [] addedCuts_;
5407    addedCuts_ = new CbcCountRowCut * [maximumNumberCuts_];
5408  }
5409  lastws->setSize(numberColumns,numberRowsAtContinuous_+currentNumberCuts);
5410/*
5411  This last bit of code traverses the path collected in walkback_ from the
5412  root back to node. At the end of the loop,
5413   * lastws will be an appropriate basis for node;
5414   * variable bounds in the constraint system will be set to be correct for
5415     node; and
5416   * addedCuts_ will be set to a list of cuts that need to be added to the
5417     constraint system at node.
5418  applyToModel does all the heavy lifting.
5419*/
5420  currentNumberCuts=0;
5421  //#define CBC_PRINT2
5422#ifdef CBC_PRINT2
5423  printf("Starting bounds at node %d\n",numberNodes_);
5424#endif
5425  currentDepth_=nNode;
5426  while (nNode) {
5427    --nNode;
5428    walkback_[nNode]->applyToModel(this,lastws,addedCuts_,currentNumberCuts);
5429  }
5430  if (!lastws->fullBasis()) {
5431#ifdef COIN_DEVELOP
5432    printf("******* bad basis\n");
5433#endif
5434    int numberRows = lastws->getNumArtificial();
5435    int i;
5436    for (i=0;i<numberRows;i++)
5437      lastws->setArtifStatus(i,CoinWarmStartBasis::basic);
5438    int numberColumns = lastws->getNumStructural();
5439    for (i=0;i<numberColumns;i++) {
5440      if (lastws->getStructStatus(i)==CoinWarmStartBasis::basic)
5441        lastws->setStructStatus(i,CoinWarmStartBasis::atLowerBound);
5442    }
5443#if 0
5444  } else {
5445    // OPTION - take off slack cuts
5446    // First see if any cuts are slack
5447    int numberAdded = currentNumberCuts;
5448    if (saveNode<2&&false) {
5449      printf("nNode %d cuts %d\n",saveNode,currentNumberCuts);
5450      for (int i=0;i<currentNumberCuts;i++)
5451        addedCuts_[i]->print();
5452    }
5453    if (numberAdded&&saveNode<5&&!parentModel_) {
5454#if 0
5455      currentNumberCuts=0;
5456      for (int j=numberRowsAtContinuous_;
5457           j<numberAdded+numberRowsAtContinuous_;j++) {
5458        CoinWarmStartBasis::Status status = lastws->getArtifStatus(j);
5459        if (status!=CoinWarmStartBasis::basic) {
5460          lastws->setArtifStatus(currentNumberCuts+numberRowsAtContinuous_,
5461                                 status);
5462          addedCuts_[currentNumberCuts++]=addedCuts_[j-numberRowsAtContinuous_];
5463        }
5464      }
5465      if (currentNumberCuts<numberAdded) {
5466        printf("deleting %d rows\n",numberAdded-currentNumberCuts);
5467        lastws->resize(currentNumberCuts+numberRowsAtContinuous_,
5468                       lastws->getNumStructural());
5469        currentNumberCuts_=currentNumberCuts;
5470      }
5471#else
5472      int nDelete=0;
5473      for (int j=numberRowsAtContinuous_;
5474           j<numberAdded+numberRowsAtContinuous_;j++) {
5475        CoinWarmStartBasis::Status status = lastws->getArtifStatus(j);
5476        if (status==CoinWarmStartBasis::basic)
5477          nDelete++;
5478      }
5479      if (nDelete)
5480        printf("depth %d can delete %d\n",saveNode-1,nDelete);
5481#endif
5482    }
5483#endif
5484  }
5485  if (0) {
5486    int numberDebugValues=18;
5487    double * debugValues = new double[numberDebugValues];
5488    CoinZeroN(debugValues,numberDebugValues);
5489    debugValues[1]=6.0;
5490    debugValues[3]=60.0;
5491    debugValues[4]=6.0;
5492    debugValues[6]=60.0;
5493    debugValues[7]=16.0;
5494    debugValues[9]=70.0;
5495    debugValues[10]=7.0;
5496    debugValues[12]=70.0;
5497    debugValues[13]=12.0;
5498    debugValues[15]=75.0;
5499    int nBad=0;
5500    for (int j=0;j<numberColumns;j++) {
5501      if (integerInfo_[j]) {
5502        if(solver_->getColLower()[j]>debugValues[j]||
5503           solver_->getColUpper()[j]<debugValues[j]) {
5504          printf("** (%g) ** ",debugValues[j]);
5505          nBad++;
5506        }
5507        printf("%d bounds %g %g\n",j,solver_->getColLower()[j],solver_->getColUpper()[j]);
5508      }
5509    }
5510    if (nBad)
5511      printf("%d BAD\n",nBad);
5512    else
5513      printf("OKAY\n");
5514    delete [] debugValues;
5515  }
5516}
5517
5518/*
5519  adjustCuts might be a better name: If the node is feasible, we sift through
5520  the cuts collected by addCuts1, add the ones that are tight and omit the
5521  ones that are loose. If the node is infeasible, we just adjust the
5522  reference counts to reflect that we're about to prune this node and its
5523  descendants.
5524*/
5525int CbcModel::addCuts (CbcNode *node, CoinWarmStartBasis *&lastws,bool canFix)
5526{
5527/*
5528  addCuts1 performs step 1 of restoring the subproblem at this node; see the
5529  comments there.
5530*/
5531  addCuts1(node,lastws);
5532  int i;
5533  int numberColumns = getNumCols();
5534  if (solver_->getNumRows()>maximumRows_) {
5535    maximumRows_ = solver_->getNumRows();
5536    workingBasis_.resize(maximumRows_,numberColumns);
5537  }
5538  CbcNodeInfo * nodeInfo = node->nodeInfo();
5539  double cutoff = getCutoff() ;
5540  int currentNumberCuts=currentNumberCuts_;
5541  if (canFix) {
5542    bool feasible=true;
5543    const double *lower = solver_->getColLower() ;
5544    const double *upper = solver_->getColUpper() ;
5545    double * newLower = analyzeResults_;
5546    double * objLower = newLower+numberIntegers_;
5547    double * newUpper = objLower+numberIntegers_;
5548    double * objUpper = newUpper+numberIntegers_;
5549    int n=0;
5550    for (i=0;i<numberIntegers_;i++) {
5551      int iColumn = integerVariable_[i];
5552      bool changed=false;
5553      double lo = 0.0;
5554      double up = 0.0;
5555      if (objLower[i]>cutoff) {
5556        lo = lower[iColumn];
5557        up = upper[iColumn];
5558        if (lo<newLower[i]) {
5559          lo = newLower[i];
5560          solver_->setColLower(iColumn,lo);
5561          changed=true;
5562          n++;
5563        }
5564        if (objUpper[i]>cutoff) {
5565          if (up>newUpper[i]) {
5566            up = newUpper[i];
5567            solver_->setColUpper(iColumn,up);
5568            changed=true;
5569            n++;
5570          }
5571        }
5572      } else if (objUpper[i]>cutoff) {
5573        lo = lower[iColumn];
5574        up = upper[iColumn];
5575        if (up>newUpper[i]) {
5576          up = newUpper[i];
5577          solver_->setColUpper(iColumn,up);
5578          changed=true;
5579          n++;
5580        }
5581      }
5582      if (changed&&lo>up) {
5583        feasible=false;
5584        break;
5585      }
5586    }
5587    if (!feasible) {
5588      printf("analysis says node infeas\n");
5589      cutoff=-COIN_DBL_MAX;
5590    }
5591  }
5592/*
5593  If the node can't be fathomed by bound, reinstall tight cuts in the
5594  constraint system. Even if there are no cuts, we'll want to set the
5595  reconstructed basis in the solver.
5596*/
5597  if (node->objectiveValue() < cutoff||numberThreads_)
5598  { 
5599    //#   define CBC_CHECK_BASIS
5600#   ifdef CBC_CHECK_BASIS
5601    printf("addCuts: expanded basis; rows %d+%d\n",
5602           numberRowsAtContinuous_,currentNumberCuts);
5603    lastws->print();
5604#   endif
5605/*
5606  Adjust the basis and constraint system so that we retain only active cuts.
5607  There are three steps:
5608    1) Scan the basis. Sort the cuts into effective cuts to be kept and
5609       loose cuts to be dropped.
5610    2) Drop the loose cuts and resize the basis to fit.
5611    3) Install the tight cuts in the constraint system (applyRowCuts) and
5612       and install the basis (setWarmStart).
5613  Use of compressRows conveys we're compressing the basis and not just
5614  tweaking the artificialStatus_ array.
5615*/
5616    if (currentNumberCuts > 0) {
5617      int numberToAdd = 0;
5618      const OsiRowCut **addCuts;
5619      int numberToDrop = 0 ;
5620      int *cutsToDrop ;
5621      addCuts = new const OsiRowCut* [currentNumberCuts];
5622      cutsToDrop = new int[currentNumberCuts] ;
5623      assert (currentNumberCuts+numberRowsAtContinuous_<=lastws->getNumArtificial());
5624      for (i=0;i<currentNumberCuts;i++) {
5625        CoinWarmStartBasis::Status status = 
5626          lastws->getArtifStatus(i+numberRowsAtContinuous_);
5627        if (addedCuts_[i] &&
5628            (status != CoinWarmStartBasis::basic ||
5629             (addedCuts_[i]->effectiveness()>1.0e10&&
5630              !addedCuts_[i]->canDropCut(solver_,i+numberRowsAtContinuous_)))) {
5631#         ifdef CHECK_CUT_COUNTS
5632          printf("Using cut %d %x as row %d\n",i,addedCuts_[i],
5633                 numberRowsAtContinuous_+numberToAdd);
5634#         endif
5635          addCuts[numberToAdd++] = addedCuts_[i];
5636        } else {
5637#         ifdef CHECK_CUT_COUNTS
5638          printf("Dropping cut %d %x\n",i,addedCuts_[i]);
5639#         endif
5640          addedCuts_[i]=NULL;
5641          cutsToDrop[numberToDrop++] = numberRowsAtContinuous_+i ;
5642        }
5643      }
5644      int numberRowsNow=numberRowsAtContinuous_+numberToAdd;
5645      lastws->compressRows(numberToDrop,cutsToDrop) ;
5646      lastws->resize(numberRowsNow,numberColumns);
5647      solver_->applyRowCuts(numberToAdd,addCuts);
5648#     ifdef CBC_CHECK_BASIS
5649      printf("addCuts: stripped basis; rows %d + %d\n",
5650             numberRowsAtContinuous_,numberToAdd);
5651      lastws->print();
5652#     endif
5653      //for (i=0;i<numberToAdd;i++)
5654      //delete addCuts[i];
5655      delete [] addCuts;
5656      delete [] cutsToDrop ;
5657    }
5658/*
5659  Set the basis in the solver.
5660*/
5661    solver_->setWarmStart(lastws);
5662#if 0
5663    if ((numberNodes_%printFrequency_)==0) {
5664      printf("Objective %g, depth %d, unsatisfied %d\n",
5665             node->objectiveValue(),
5666             node->depth(),node->numberUnsatisfied());
5667    }
5668#endif
5669/*
5670  Clean up and we're out of here.
5671*/
5672    numberNodes_++;
5673    return 0;
5674  } 
5675/*
5676  This node has been fathomed by bound as we try to revive it out of the live
5677  set. Adjust the cut reference counts to reflect that we no longer need to
5678  explore the remaining branch arms, hence they will no longer reference any
5679  cuts. Cuts whose reference count falls to zero are deleted. 
5680*/
5681  else
5682  { int i;
5683    if (currentNumberCuts) {
5684#ifndef CBC_DETERMINISTIC_THREAD
5685      lockThread();
5686#endif
5687      int numberLeft = nodeInfo->numberBranchesLeft();
5688      for (i = 0 ; i < currentNumberCuts ; i++)
5689        { if (addedCuts_[i])
5690          { if (!addedCuts_[i]->decrement(numberLeft))
5691            { delete addedCuts_[i];
5692            addedCuts_[i] = NULL; } } }
5693#ifndef CBC_DETERMINISTIC_THREAD
5694      unlockThread();
5695#endif
5696    }
5697    return 1 ; }
5698}
5699
5700
5701/*
5702  Perform reduced cost fixing on integer variables.
5703
5704  The variables in question are already nonbasic at bound. We're just nailing
5705  down the current situation.
5706*/
5707#if 0
5708int CbcModel::reducedCostFix ()
5709
5710{
5711  if(!solverCharacteristics_->reducedCostsAccurate())
5712    return 0; //NLP
5713  double cutoff = getCutoff() ;
5714  double direction = solver_->getObjSense() ;
5715  double gap = cutoff - solver_->getObjValue()*direction ;
5716  double tolerance;
5717  solver_->getDblParam(OsiDualTolerance,tolerance) ;
5718  if (gap<=0.0)
5719    gap = tolerance; //return 0;
5720  gap += 100.0*tolerance;
5721  double integerTolerance = getDblParam(CbcIntegerTolerance) ;
5722
5723  const double *lower = solver_->getColLower() ;
5724  const double *upper = solver_->getColUpper() ;
5725  const double *solution = solver_->getColSolution() ;
5726  const double *reducedCost = solver_->getReducedCost() ;
5727
5728  int numberFixed = 0 ;
5729
5730# ifdef COIN_HAS_CLP
5731  OsiClpSolverInterface * clpSolver
5732    = dynamic_cast<OsiClpSolverInterface *> (solver_);
5733  ClpSimplex * clpSimplex=NULL;
5734  if (clpSolver)
5735    clpSimplex = clpSolver->getModelPtr();
5736# endif
5737  for (int i = 0 ; i < numberIntegers_ ; i++) {
5738    int iColumn = integerVariable_[i] ;
5739    double djValue = direction*reducedCost[iColumn] ;
5740    if (upper[iColumn]-lower[iColumn] > integerTolerance) {
5741      if (solution[iColumn] < lower[iColumn]+integerTolerance && djValue > gap) {
5742#ifdef COIN_HAS_CLP
5743        // may just have been fixed before
5744        if (clpSimplex) {
5745          if (clpSimplex->getColumnStatus(iColumn)==ClpSimplex::basic) {
5746#ifdef COIN_DEVELOP
5747            printf("DJfix %d has status of %d, dj of %g gap %g, bounds %g %g\n",
5748                   iColumn,clpSimplex->getColumnStatus(iColumn),
5749                   djValue,gap,lower[iColumn],upper[iColumn]);
5750#endif
5751          } else {         
5752            assert(clpSimplex->getColumnStatus(iColumn)==ClpSimplex::atLowerBound||
5753                   clpSimplex->getColumnStatus(iColumn)==ClpSimplex::isFixed);
5754          }
5755        }
5756#endif
5757        solver_->setColUpper(iColumn,lower[iColumn]) ;
5758        numberFixed++ ;
5759      } else if (solution[iColumn] > upper[iColumn]-integerTolerance && -djValue > gap) {
5760#ifdef COIN_HAS_CLP
5761        // may just have been fixed before
5762        if (clpSimplex) {
5763          if (clpSimplex->getColumnStatus(iColumn)==ClpSimplex::basic) {
5764#ifdef COIN_DEVELOP
5765            printf("DJfix %d has status of %d, dj of %g gap %g, bounds %g %g\n",
5766                   iColumn,clpSimplex->getColumnStatus(iColumn),
5767                   djValue,gap,lower[iColumn],upper[iColumn]);
5768#endif
5769          } else {         
5770            assert(clpSimplex->getColumnStatus(iColumn)==ClpSimplex::atUpperBound||
5771                   clpSimplex->getColumnStatus(iColumn)==ClpSimplex::isFixed);
5772          }
5773        }
5774#endif
5775        solver_->setColLower(iColumn,upper[iColumn]) ;
5776        numberFixed++ ;
5777      }
5778    }
5779  }
5780  numberDJFixed_ += numberFixed;
5781  return numberFixed; }
5782#else
5783int CbcModel::reducedCostFix ()
5784
5785{
5786  if(!solverCharacteristics_->reducedCostsAccurate())
5787    return 0; //NLP
5788  double cutoff = getCutoff() ;
5789  double direction = solver_->getObjSense() ;
5790  double gap = cutoff - solver_->getObjValue()*direction ;
5791  double tolerance;
5792  solver_->getDblParam(OsiDualTolerance,tolerance) ;
5793  if (gap<=0.0)
5794    gap = tolerance; //return 0;
5795  gap += 100.0*tolerance;
5796  double integerTolerance = getDblParam(CbcIntegerTolerance) ;
5797
5798  const double *lower = solver_->getColLower() ;
5799  const double *upper = solver_->getColUpper() ;
5800  const double *solution = solver_->getColSolution() ;
5801  const double *reducedCost = solver_->getReducedCost() ;
5802
5803  int numberFixed = 0 ;
5804  int numberTightened = 0 ;
5805
5806# ifdef COIN_HAS_CLP
5807  OsiClpSolverInterface * clpSolver
5808    = dynamic_cast<OsiClpSolverInterface *> (solver_);
5809  ClpSimplex * clpSimplex=NULL;
5810  if (clpSolver)
5811    clpSimplex = clpSolver->getModelPtr();
5812# endif
5813  for (int i = 0 ; i < numberIntegers_ ; i++) {
5814    int iColumn = integerVariable_[i] ;
5815    double djValue = direction*reducedCost[iColumn] ;
5816    double boundGap=upper[iColumn]-lower[iColumn];
5817    if (boundGap > integerTolerance) {
5818      if (solution[iColumn] < lower[iColumn]+integerTolerance
5819          && djValue*boundGap > gap) {
5820#ifdef COIN_HAS_CLP
5821        // may just have been fixed before
5822        if (clpSimplex) {
5823          if (clpSimplex->getColumnStatus(iColumn)==ClpSimplex::basic) {
5824#ifdef COIN_DEVELOP
5825            printf("DJfix %d has status of %d, dj of %g gap %g, bounds %g %g\n",
5826                   iColumn,clpSimplex->getColumnStatus(iColumn),
5827                   djValue,gap,lower[iColumn],upper[iColumn]);
5828#endif
5829          } else {         
5830            assert(clpSimplex->getColumnStatus(iColumn)==ClpSimplex::atLowerBound||
5831                   clpSimplex->getColumnStatus(iColumn)==ClpSimplex::isFixed);
5832          }
5833        }
5834#endif
5835        double newBound=lower[iColumn];
5836        if (boundGap>1.99) {
5837          boundGap = gap/djValue+1.0e-4*boundGap;
5838          newBound=lower[iColumn]+floor(boundGap);
5839          numberTightened++;
5840          //if (newBound)
5841          //printf("tighter - gap %g dj %g newBound %g\n",
5842          //   gap,djValue,newBound);
5843        }
5844        solver_->setColUpper(iColumn,newBound) ;
5845        numberFixed++ ;
5846      } else if (solution[iColumn] > upper[iColumn]-integerTolerance && -djValue > boundGap*gap) {
5847#ifdef COIN_HAS_CLP
5848        // may just have been fixed before
5849        if (clpSimplex) {
5850          if (clpSimplex->getColumnStatus(iColumn)==ClpSimplex::basic) {
5851#ifdef COIN_DEVELOP
5852            printf("DJfix %d has status of %d, dj of %g gap %g, bounds %g %g\n",
5853                   iColumn,clpSimplex->getColumnStatus(iColumn),
5854                   djValue,gap,lower[iColumn],upper[iColumn]);
5855#endif
5856          } else {         
5857            assert(clpSimplex->getColumnStatus(iColumn)==ClpSimplex::atUpperBound||
5858                   clpSimplex->getColumnStatus(iColumn)==ClpSimplex::isFixed);
5859          }
5860        }
5861#endif
5862        double newBound=upper[iColumn];
5863        if (boundGap>1.99) {
5864          boundGap = -gap/djValue+1.0e-4*boundGap;
5865          newBound=upper[iColumn]-floor(boundGap);
5866          //if (newBound)
5867          //printf("tighter - gap %g dj %g newBound %g\n",
5868          //   gap,djValue,newBound);
5869          numberTightened++;
5870        }
5871        solver_->setColLower(iColumn,newBound) ;
5872        numberFixed++ ;
5873      }
5874    }
5875  }
5876  numberDJFixed_ += numberFixed-numberTightened;
5877#ifdef COIN_DEVELOP
5878  //if (numberTightened)
5879  //printf("%d tightened, %d others fixed\n",numberTightened,
5880  //   numberFixed-numberTightened);
5881#endif
5882  return numberFixed; }
5883#endif
5884// Collect coding to replace whichGenerator
5885void
5886CbcModel::resizeWhichGenerator(int numberNow, int numberAfter)
5887{
5888  if (numberAfter > maximumWhich_) {
5889    maximumWhich_ = CoinMax(maximumWhich_*2+100,numberAfter) ;
5890    int * temp = new int[2*maximumWhich_] ;
5891    memcpy(temp,whichGenerator_,numberNow*sizeof(int)) ;
5892    delete [] whichGenerator_ ;
5893    whichGenerator_ = temp ;
5894    memset(whichGenerator_+numberNow,0,(maximumWhich_-numberNow)*sizeof(int));
5895  }
5896}
5897
5898/** Solve the model using cuts
5899
5900  This version takes off redundant cuts from node.
5901  Returns true if feasible.
5902
5903  \todo
5904  Why do I need to resolve the problem? What has been done between the last
5905  relaxation and calling solveWithCuts?
5906
5907  If numberTries == 0 then user did not want any cuts.
5908*/
5909
5910bool 
5911CbcModel::solveWithCuts (OsiCuts &cuts, int numberTries, CbcNode *node)
5912/*