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

Last change on this file since 1255 was 1255, checked in by forrest, 10 years ago

fix cut over counting and a compiler warning

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