source: stable/2.6/Cbc/src/CbcModel.cpp @ 1627

Last change on this file since 1627 was 1627, checked in by tkr, 8 years ago

Merging r1623 from trunk

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