source: branches/sandbox/Cbc/src/CbcModel.cpp @ 1338

Last change on this file since 1338 was 1338, checked in by ladanyi, 10 years ago

temporarily exclude broken code for parallel cbc

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