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

Last change on this file since 1573 was 1573, checked in by lou, 8 years ago

Change to EPL license notice.

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