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

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

final changes before cleaning

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