source: stable/2.7/Cbc/src/CbcModel.cpp @ 1833

Last change on this file since 1833 was 1833, checked in by stefan, 6 years ago

reset numberSavedSolutions_ to 0 after deleting saved solutions

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