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

Last change on this file since 1623 was 1623, checked in by forrest, 9 years ago

allow for elapsed time on unix

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