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

Last change on this file since 1521 was 1521, checked in by lou, 9 years ago

Changes to externals to make it easier to generate a new stable. Catches a few
changes in comments down in src.

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