source: stable/2.9/Cbc/src/CbcNode.cpp @ 2187

Last change on this file since 2187 was 2187, checked in by stefan, 3 years ago

more sync with trunk

  • Property svn:eol-style set to native
  • Property svn:keywords set to Author Date Id Revision
File size: 249.3 KB
Line 
1/* $Id: CbcNode.cpp 2187 2015-05-05 13:04:34Z 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#include <windows.h> // for Sleep()
10#ifdef small
11#undef small
12#endif
13#else
14#include <unistd.h> // for usleep()
15#endif
16
17#include "CbcConfig.h"
18#ifdef COIN_HAS_NTY
19#include "CbcSymmetry.hpp"
20#endif
21//#define DEBUG_SOLUTION
22#ifdef DEBUG_SOLUTION
23#define COIN_DETAIL
24#endif
25#include <string>
26//#define CBC_DEBUG 1
27//#define CHECK_CUT_COUNTS
28//#define CHECK_NODE
29//#define CBC_CHECK_BASIS
30#include <cassert>
31#include <cfloat>
32#define CUTS
33#include "OsiSolverInterface.hpp"
34#include "OsiChooseVariable.hpp"
35#include "OsiAuxInfo.hpp"
36#include "OsiSolverBranch.hpp"
37#include "CoinWarmStartBasis.hpp"
38#include "CoinTime.hpp"
39#include "CbcModel.hpp"
40#include "CbcNode.hpp"
41#include "CbcStatistics.hpp"
42#include "CbcStrategy.hpp"
43#include "CbcBranchActual.hpp"
44#include "CbcBranchDynamic.hpp"
45#include "OsiRowCut.hpp"
46#include "OsiRowCutDebugger.hpp"
47#include "OsiCuts.hpp"
48#include "CbcCountRowCut.hpp"
49#include "CbcFeasibilityBase.hpp"
50#include "CbcMessage.hpp"
51#ifdef COIN_HAS_CLP
52#include "OsiClpSolverInterface.hpp"
53#include "ClpSimplexOther.hpp"
54#include "ClpSolve.hpp"
55#include "ClpDualRowSteepest.hpp"
56#include "ClpPrimalColumnPivot.hpp"
57#endif
58using namespace std;
59#include "CglCutGenerator.hpp"
60
61
62CbcNode::CbcNode() :
63        nodeInfo_(NULL),
64        objectiveValue_(1.0e100),
65        guessedObjectiveValue_(1.0e100),
66        sumInfeasibilities_(0.0),
67        branch_(NULL),
68        depth_(-1),
69        numberUnsatisfied_(0),
70        nodeNumber_(-1),
71        state_(0)
72{
73#ifdef CHECK_NODE
74    printf("CbcNode %p Constructor\n", this);
75#endif
76}
77// Print
78void
79CbcNode::print() const
80{
81    printf("number %d obj %g depth %d sumun %g nunsat %d state %d\n",
82           nodeNumber_, objectiveValue_, depth_, sumInfeasibilities_, numberUnsatisfied_, state_);
83}
84CbcNode::CbcNode(CbcModel * model,
85                 CbcNode * lastNode) :
86        nodeInfo_(NULL),
87        objectiveValue_(1.0e100),
88        guessedObjectiveValue_(1.0e100),
89        sumInfeasibilities_(0.0),
90        branch_(NULL),
91        depth_(-1),
92        numberUnsatisfied_(0),
93        nodeNumber_(-1),
94        state_(0)
95{
96#ifdef CHECK_NODE
97    printf("CbcNode %p Constructor from model\n", this);
98#endif
99    model->setObjectiveValue(this, lastNode);
100
101    if (lastNode) {
102        if (lastNode->nodeInfo_) {
103            lastNode->nodeInfo_->increment();
104        }
105    }
106    nodeNumber_ = model->getNodeCount();
107}
108
109#define CBC_NEW_CREATEINFO
110#ifdef CBC_NEW_CREATEINFO
111
112/*
113  New createInfo, with basis manipulation hidden inside mergeBasis. Allows
114  solvers to override and carry over all information from one basis to
115  another.
116*/
117
118void
119CbcNode::createInfo (CbcModel *model,
120                     CbcNode *lastNode,
121                     const CoinWarmStartBasis *lastws,
122                     const double *lastLower, const double *lastUpper,
123                     int numberOldActiveCuts, int numberNewCuts)
124
125{
126    OsiSolverInterface *solver = model->solver();
127    CbcStrategy *strategy = model->strategy();
128    /*
129      The root --- no parent. Create full basis and bounds information.
130    */
131    if (!lastNode) {
132        if (!strategy)
133            nodeInfo_ = new CbcFullNodeInfo(model, solver->getNumRows());
134        else
135            nodeInfo_ = strategy->fullNodeInfo(model, solver->getNumRows());
136    } else {
137        /*
138          Not the root. Create an edit from the parent's basis & bound information.
139          This is not quite as straightforward as it seems. We need to reintroduce
140          cuts we may have dropped out of the basis, in the correct position, because
141          this whole process is strictly positional. Start by grabbing the current
142          basis.
143        */
144        bool mustDeleteBasis;
145        const CoinWarmStartBasis *ws =
146            dynamic_cast<const CoinWarmStartBasis*>(solver->getPointerToWarmStart(mustDeleteBasis));
147        assert(ws != NULL); // make sure not volume
148        //int numberArtificials = lastws->getNumArtificial();
149        int numberColumns = solver->getNumCols();
150        int numberRowsAtContinuous = model->numberRowsAtContinuous();
151        int currentNumberCuts = model->currentNumberCuts();
152#   ifdef CBC_CHECK_BASIS
153        std::cout
154            << "Before expansion: orig " << numberRowsAtContinuous
155            << ", old " << numberOldActiveCuts
156            << ", new " << numberNewCuts
157            << ", current " << currentNumberCuts << "." << std::endl ;
158        ws->print();
159#   endif
160        /*
161          Clone the basis and resize it to hold the structural constraints, plus
162          all the cuts: old cuts, both active and inactive (currentNumberCuts),
163          and new cuts (numberNewCuts). This will become the expanded basis.
164        */
165        CoinWarmStartBasis *expanded =
166            dynamic_cast<CoinWarmStartBasis *>(ws->clone()) ;
167        int iCompact = numberRowsAtContinuous + numberOldActiveCuts + numberNewCuts ;
168        // int nPartial = numberRowsAtContinuous+currentNumberCuts;
169        int iFull = numberRowsAtContinuous + currentNumberCuts + numberNewCuts;
170        // int maxBasisLength = ((iFull+15)>>4)+((numberColumns+15)>>4);
171        // printf("l %d full %d\n",maxBasisLength,iFull);
172        expanded->resize(iFull, numberColumns);
173#   ifdef CBC_CHECK_BASIS
174        std::cout
175            << "\tFull basis " << iFull << " rows, "
176            << numberColumns << " columns; compact "
177            << iCompact << " rows." << std::endl ;
178#   endif
179        /*
180          Now flesh out the expanded basis. The clone already has the
181          correct status information for the variables and for the structural
182          (numberRowsAtContinuous) constraints. Any indices beyond nPartial must be
183          cuts created while processing this node --- they can be copied en bloc
184          into the correct position in the expanded basis. The space reserved for
185          xferRows is a gross overestimate.
186        */
187        CoinWarmStartBasis::XferVec xferRows ;
188        xferRows.reserve(iFull - numberRowsAtContinuous + 1) ;
189        if (numberNewCuts) {
190            xferRows.push_back(
191                CoinWarmStartBasis::XferEntry(iCompact - numberNewCuts,
192                                              iFull - numberNewCuts, numberNewCuts)) ;
193        }
194        /*
195          From nPartial down, record the entries we want to copy from the current
196          basis (the entries for the active cuts; non-zero in the list returned
197          by addedCuts). Fill the expanded basis with entries showing a status of
198          basic for the deactivated (loose) cuts.
199        */
200        CbcCountRowCut **cut = model->addedCuts();
201        iFull -= (numberNewCuts + 1) ;
202        iCompact -= (numberNewCuts + 1) ;
203        int runLen = 0 ;
204        CoinWarmStartBasis::XferEntry entry(-1, -1, -1) ;
205        while (iFull >= numberRowsAtContinuous) {
206            for ( ; iFull >= numberRowsAtContinuous &&
207                    cut[iFull-numberRowsAtContinuous] ; iFull--)
208                runLen++ ;
209            if (runLen) {
210                iCompact -= runLen ;
211                entry.first = iCompact + 1 ;
212                entry.second = iFull + 1 ;
213                entry.third = runLen ;
214                runLen = 0 ;
215                xferRows.push_back(entry) ;
216            }
217            for ( ; iFull >= numberRowsAtContinuous &&
218                    !cut[iFull-numberRowsAtContinuous] ; iFull--)
219                expanded->setArtifStatus(iFull, CoinWarmStartBasis::basic);
220        }
221        /*
222          Finally, call mergeBasis to copy over entries from the current basis to
223          the expanded basis. Since we cloned the expanded basis from the active basis
224          and haven't changed the number of variables, only row status entries need
225          to be copied.
226        */
227        expanded->mergeBasis(ws, &xferRows, 0) ;
228
229#ifdef CBC_CHECK_BASIS
230        std::cout << "Expanded basis:" << std::endl ;
231        expanded->print() ;
232        std::cout << "Diffing against:" << std::endl ;
233        lastws->print() ;
234#endif
235        assert (expanded->getNumArtificial() >= lastws->getNumArtificial());
236#ifdef CLP_INVESTIGATE
237        if (!expanded->fullBasis()) {
238            int iFull = numberRowsAtContinuous + currentNumberCuts + numberNewCuts;
239            printf("cont %d old %d new %d current %d full inc %d full %d\n",
240                   numberRowsAtContinuous, numberOldActiveCuts, numberNewCuts,
241                   currentNumberCuts, iFull, iFull - numberNewCuts);
242        }
243#endif
244
245        /*
246          Now that we have two bases in proper positional correspondence, creating
247          the actual diff is dead easy.
248
249          Note that we're going to compare the expanded basis here to the stripped
250          basis (lastws) produced by addCuts. It doesn't affect the correctness (the
251          diff process has no knowledge of the meaning of an entry) but it does
252          mean that we'll always generate a whack of diff entries because the expanded
253          basis is considerably larger than the stripped basis.
254        */
255        CoinWarmStartDiff *basisDiff = expanded->generateDiff(lastws) ;
256
257        /*
258          Diff the bound vectors. It's assumed the number of structural variables
259          is not changing. For branching objects that change bounds on integer
260          variables, we should see at least one bound change as a consequence
261          of applying the branch that generated this subproblem from its parent.
262          This need not hold for other types of branching objects (hyperplane
263          branches, for example).
264        */
265        const double * lower = solver->getColLower();
266        const double * upper = solver->getColUpper();
267
268        double *boundChanges = new double [2*numberColumns] ;
269        int *variables = new int [2*numberColumns] ;
270        int numberChangedBounds = 0;
271
272        int i;
273        for (i = 0; i < numberColumns; i++) {
274            if (lower[i] != lastLower[i]) {
275                variables[numberChangedBounds] = i;
276                boundChanges[numberChangedBounds++] = lower[i];
277            }
278            if (upper[i] != lastUpper[i]) {
279                variables[numberChangedBounds] = i | 0x80000000;
280                boundChanges[numberChangedBounds++] = upper[i];
281            }
282#ifdef CBC_DEBUG
283            if (lower[i] != lastLower[i]) {
284                std::cout
285                    << "lower on " << i << " changed from "
286                    << lastLower[i] << " to " << lower[i] << std::endl ;
287            }
288            if (upper[i] != lastUpper[i]) {
289                std::cout
290                    << "upper on " << i << " changed from "
291                    << lastUpper[i] << " to " << upper[i] << std::endl ;
292            }
293#endif
294        }
295#ifdef CBC_DEBUG
296        std::cout << numberChangedBounds << " changed bounds." << std::endl ;
297#endif
298        //if (lastNode->branchingObject()->boundBranch())
299        //assert (numberChangedBounds);
300        /*
301          Hand the lot over to the CbcPartialNodeInfo constructor, then clean up and
302          return.
303        */
304        if (!strategy)
305            nodeInfo_ =
306                new CbcPartialNodeInfo(lastNode->nodeInfo_, this, numberChangedBounds,
307                                       variables, boundChanges, basisDiff) ;
308        else
309            nodeInfo_ =
310                strategy->partialNodeInfo(model, lastNode->nodeInfo_, this,
311                                          numberChangedBounds, variables, boundChanges,
312                                          basisDiff) ;
313        delete basisDiff ;
314        delete [] boundChanges;
315        delete [] variables;
316        delete expanded ;
317        if  (mustDeleteBasis)
318            delete ws;
319    }
320    // Set node number
321    nodeInfo_->setNodeNumber(model->getNodeCount2());
322    state_ |= 2; // say active
323}
324
325#else   // CBC_NEW_CREATEINFO
326
327/*
328  Original createInfo, with bare manipulation of basis vectors. Fails if solver
329  maintains additional information in basis.
330*/
331
332void
333CbcNode::createInfo (CbcModel *model,
334                     CbcNode *lastNode,
335                     const CoinWarmStartBasis *lastws,
336                     const double *lastLower, const double *lastUpper,
337                     int numberOldActiveCuts, int numberNewCuts)
338{
339    OsiSolverInterface * solver = model->solver();
340    CbcStrategy * strategy = model->strategy();
341    /*
342      The root --- no parent. Create full basis and bounds information.
343    */
344    if (!lastNode) {
345        if (!strategy)
346            nodeInfo_ = new CbcFullNodeInfo(model, solver->getNumRows());
347        else
348            nodeInfo_ = strategy->fullNodeInfo(model, solver->getNumRows());
349    }
350    /*
351      Not the root. Create an edit from the parent's basis & bound information.
352      This is not quite as straightforward as it seems. We need to reintroduce
353      cuts we may have dropped out of the basis, in the correct position, because
354      this whole process is strictly positional. Start by grabbing the current
355      basis.
356    */
357    else {
358        bool mustDeleteBasis;
359        const CoinWarmStartBasis* ws =
360            dynamic_cast<const CoinWarmStartBasis*>(solver->getPointerToWarmStart(mustDeleteBasis));
361        assert(ws != NULL); // make sure not volume
362        //int numberArtificials = lastws->getNumArtificial();
363        int numberColumns = solver->getNumCols();
364
365        const double * lower = solver->getColLower();
366        const double * upper = solver->getColUpper();
367
368        int i;
369        /*
370        Create a clone and resize it to hold all the structural constraints, plus
371        all the cuts: old cuts, both active and inactive (currentNumberCuts), and
372        new cuts (numberNewCuts).
373
374        TODO: You'd think that the set of constraints (logicals) in the expanded
375        basis should match the set represented in lastws. At least, that's
376        what I thought. But at the point I first looked hard at this bit of
377        code, it turned out that lastws was the stripped basis produced at
378        the end of addCuts(), rather than the raw basis handed back by
379        addCuts1(). The expanded basis here is equivalent to the raw basis of
380        addCuts1(). I said ``whoa, that's not good, I must have introduced a
381        bug'' and went back to John's code to see where I'd gone wrong.
382        And discovered the same `error' in his code.
383
384        After a bit of thought, my conclusion is that correctness is not
385        affected by whether lastws is the stripped or raw basis. The diffs
386        have no semantics --- just a set of changes that need to be made
387        to convert lastws into expanded. I think the only effect is that we
388        store a lot more diffs (everything in expanded that's not covered by
389        the stripped basis). But I need to give this more thought. There
390        may well be some subtle error cases.
391
392        In the mean time, I've twiddled addCuts() to set lastws to the raw
393        basis. Makes me (Lou) less nervous to compare apples to apples.
394        */
395        CoinWarmStartBasis *expanded =
396            dynamic_cast<CoinWarmStartBasis *>(ws->clone()) ;
397        int numberRowsAtContinuous = model->numberRowsAtContinuous();
398        int iFull = numberRowsAtContinuous + model->currentNumberCuts() +
399                    numberNewCuts;
400        //int numberArtificialsNow = iFull;
401        //int maxBasisLength = ((iFull+15)>>4)+((numberColumns+15)>>4);
402        //printf("l %d full %d\n",maxBasisLength,iFull);
403        if (expanded)
404            expanded->resize(iFull, numberColumns);
405#ifdef CBC_CHECK_BASIS
406        printf("Before expansion: orig %d, old %d, new %d, current %d\n",
407               numberRowsAtContinuous, numberOldActiveCuts, numberNewCuts,
408               model->currentNumberCuts()) ;
409        ws->print();
410#endif
411        /*
412        Now fill in the expanded basis. Any indices beyond nPartial must
413        be cuts created while processing this node --- they can be copied directly
414        into the expanded basis. From nPartial down, pull the status of active cuts
415        from ws, interleaving with a B entry for the deactivated (loose) cuts.
416        */
417        int numberDropped = model->currentNumberCuts() - numberOldActiveCuts;
418        int iCompact = iFull - numberDropped;
419        CbcCountRowCut ** cut = model->addedCuts();
420        int nPartial = model->currentNumberCuts() + numberRowsAtContinuous;
421        iFull--;
422        for (; iFull >= nPartial; iFull--) {
423            CoinWarmStartBasis::Status status = ws->getArtifStatus(--iCompact);
424            //assert (status != CoinWarmStartBasis::basic); // may be permanent cut
425            expanded->setArtifStatus(iFull, status);
426        }
427        for (; iFull >= numberRowsAtContinuous; iFull--) {
428            if (cut[iFull-numberRowsAtContinuous]) {
429                CoinWarmStartBasis::Status status = ws->getArtifStatus(--iCompact);
430                // If no cut generator being used then we may have basic variables
431                //if (model->getMaximumCutPasses()&&
432                //  status == CoinWarmStartBasis::basic)
433                //printf("cut basic\n");
434                expanded->setArtifStatus(iFull, status);
435            } else {
436                expanded->setArtifStatus(iFull, CoinWarmStartBasis::basic);
437            }
438        }
439#ifdef CBC_CHECK_BASIS
440        printf("Expanded basis\n");
441        expanded->print() ;
442        printf("Diffing against\n") ;
443        lastws->print() ;
444#endif
445        /*
446        Now that we have two bases in proper positional correspondence, creating
447        the actual diff is dead easy.
448        */
449
450        CoinWarmStartDiff *basisDiff = expanded->generateDiff(lastws) ;
451        /*
452        Diff the bound vectors. It's assumed the number of structural variables is
453        not changing. Assuming that branching objects all involve integer variables,
454        we should see at least one bound change as a consequence of processing this
455        subproblem. Different types of branching objects could break this assertion.
456        Not true at all - we have not applied current branch - JJF.
457        */
458        double *boundChanges = new double [2*numberColumns] ;
459        int *variables = new int [2*numberColumns] ;
460        int numberChangedBounds = 0;
461        for (i = 0; i < numberColumns; i++) {
462            if (lower[i] != lastLower[i]) {
463                variables[numberChangedBounds] = i;
464                boundChanges[numberChangedBounds++] = lower[i];
465            }
466            if (upper[i] != lastUpper[i]) {
467                variables[numberChangedBounds] = i | 0x80000000;
468                boundChanges[numberChangedBounds++] = upper[i];
469            }
470#ifdef CBC_DEBUG
471            if (lower[i] != lastLower[i])
472                printf("lower on %d changed from %g to %g\n",
473                       i, lastLower[i], lower[i]);
474            if (upper[i] != lastUpper[i])
475                printf("upper on %d changed from %g to %g\n",
476                       i, lastUpper[i], upper[i]);
477#endif
478        }
479#ifdef CBC_DEBUG
480        printf("%d changed bounds\n", numberChangedBounds) ;
481#endif
482        //if (lastNode->branchingObject()->boundBranch())
483        //assert (numberChangedBounds);
484        /*
485        Hand the lot over to the CbcPartialNodeInfo constructor, then clean up and
486        return.
487        */
488        if (!strategy)
489            nodeInfo_ =
490                new CbcPartialNodeInfo(lastNode->nodeInfo_, this, numberChangedBounds,
491                                       variables, boundChanges, basisDiff) ;
492        else
493            nodeInfo_ = strategy->partialNodeInfo(model, lastNode->nodeInfo_, this, numberChangedBounds,
494                                                  variables, boundChanges, basisDiff) ;
495        delete basisDiff ;
496        delete [] boundChanges;
497        delete [] variables;
498        delete expanded ;
499        if  (mustDeleteBasis)
500            delete ws;
501    }
502    // Set node number
503    nodeInfo_->setNodeNumber(model->getNodeCount2());
504    state_ |= 2; // say active
505}
506
507#endif  // CBC_NEW_CREATEINFO
508/*
509  The routine scans through the object list of the model looking for objects
510  that indicate infeasibility. It tests each object using strong branching
511  and selects the one with the least objective degradation.  A corresponding
512  branching object is left attached to lastNode.
513
514  If strong branching is disabled, a candidate object is chosen essentially
515  at random (whatever object ends up in pos'n 0 of the candidate array).
516
517  If a branching candidate is found to be monotone, bounds are set to fix the
518  variable and the routine immediately returns (the caller is expected to
519  reoptimize).
520
521  If a branching candidate is found to result in infeasibility in both
522  directions, the routine immediately returns an indication of infeasibility.
523
524  Returns:  0   both branch directions are feasible
525  -1    branching variable is monotone
526  -2    infeasible
527
528  Original comments:
529  Here could go cuts etc etc
530  For now just fix on objective from strong branching.
531*/
532
533int CbcNode::chooseBranch (CbcModel *model, CbcNode *lastNode, int numberPassesLeft)
534
535{
536    if (lastNode)
537        depth_ = lastNode->depth_ + 1;
538    else
539        depth_ = 0;
540    delete branch_;
541    branch_ = NULL;
542    OsiSolverInterface * solver = model->solver();
543# ifdef COIN_HAS_CLP
544    OsiClpSolverInterface * osiclp = dynamic_cast< OsiClpSolverInterface*> (solver);
545    int saveClpOptions = 0;
546    if (osiclp) {
547        // for faster hot start
548        saveClpOptions = osiclp->specialOptions();
549        osiclp->setSpecialOptions(saveClpOptions | 8192);
550    }
551# else
552    OsiSolverInterface *osiclp = NULL ;
553# endif
554    double saveObjectiveValue = solver->getObjValue();
555    double objectiveValue = CoinMax(solver->getObjSense() * saveObjectiveValue, objectiveValue_);
556    const double * lower = solver->getColLower();
557    const double * upper = solver->getColUpper();
558    // See what user thinks
559    int anyAction = model->problemFeasibility()->feasible(model, 0);
560    if (anyAction) {
561        // will return -2 if infeasible , 0 if treat as integer
562        return anyAction - 1;
563    }
564    double integerTolerance =
565        model->getDblParam(CbcModel::CbcIntegerTolerance);
566    // point to useful information
567    OsiBranchingInformation usefulInfo = model->usefulInformation();
568    // and modify
569    usefulInfo.depth_ = depth_;
570    int i;
571    bool beforeSolution = model->getSolutionCount() == 0;
572    int numberStrong = model->numberStrong();
573    // switch off strong if hotstart
574    const double * hotstartSolution = model->hotstartSolution();
575    const int * hotstartPriorities = model->hotstartPriorities();
576    int numberObjects = model->numberObjects();
577    int numberColumns = model->getNumCols();
578    double * saveUpper = new double[numberColumns];
579    double * saveLower = new double[numberColumns];
580    for (i = 0; i < numberColumns; i++) {
581        saveLower[i] = lower[i];
582        saveUpper[i] = upper[i];
583    }
584
585    // Save solution in case heuristics need good solution later
586
587    double * saveSolution = new double[numberColumns];
588    memcpy(saveSolution, solver->getColSolution(), numberColumns*sizeof(double));
589    model->reserveCurrentSolution(saveSolution);
590    if (hotstartSolution) {
591        numberStrong = 0;
592        if ((model->moreSpecialOptions()&1024) != 0) {
593            int nBad = 0;
594            int nUnsat = 0;
595            int nDiff = 0;
596            for (int i = 0; i < numberObjects; i++) {
597                OsiObject * object = model->modifiableObject(i);
598                const CbcSimpleInteger * thisOne = dynamic_cast <const CbcSimpleInteger *> (object);
599                if (thisOne) {
600                    int iColumn = thisOne->columnNumber();
601                    double targetValue = hotstartSolution[iColumn];
602                    double value = saveSolution[iColumn];
603                    if (fabs(value - floor(value + 0.5)) > 1.0e-6) {
604                        nUnsat++;
605#ifdef CLP_INVESTIGATE
606                        printf("H %d is %g target %g\n", iColumn, value, targetValue);
607#endif
608                    } else if (fabs(targetValue - value) > 1.0e-6) {
609                        nDiff++;
610                    }
611                    if (targetValue < saveLower[iColumn] ||
612                            targetValue > saveUpper[iColumn]) {
613#ifdef CLP_INVESTIGATE
614                        printf("%d has target %g and current bounds %g and %g\n",
615                               iColumn, targetValue, saveLower[iColumn], saveUpper[iColumn]);
616#endif
617                        nBad++;
618                    }
619                }
620            }
621#ifdef CLP_INVESTIGATE
622            printf("Hot %d unsatisfied, %d outside limits, %d different\n",
623                   nUnsat, nBad, nDiff);
624#endif
625            if (nBad) {
626                // switch off as not possible
627                hotstartSolution = NULL;
628                model->setHotstartSolution(NULL, NULL);
629                usefulInfo.hotstartSolution_ = NULL;
630            }
631        }
632    }
633    int numberStrongDone = 0;
634    int numberUnfinished = 0;
635    int numberStrongInfeasible = 0;
636    int numberStrongIterations = 0;
637    int saveNumberStrong = numberStrong;
638    bool checkFeasibility = numberObjects > model->numberIntegers();
639    int maximumStrong = CoinMax(CoinMin(numberStrong, numberObjects), 1);
640    /*
641      Get a branching decision object. Use the default decision criteria unless
642      the user has loaded a decision method into the model.
643    */
644    CbcBranchDecision *decision = model->branchingMethod();
645    CbcDynamicPseudoCostBranchingObject * dynamicBranchingObject =
646        dynamic_cast<CbcDynamicPseudoCostBranchingObject *>(decision);
647    if (!decision || dynamicBranchingObject)
648        decision = new CbcBranchDefaultDecision();
649    decision->initialize(model);
650    CbcStrongInfo * choice = new CbcStrongInfo[maximumStrong];
651    // May go round twice if strong branching fixes all local candidates
652    bool finished = false;
653    double estimatedDegradation = 0.0;
654    while (!finished) {
655        finished = true;
656        // Some objects may compute an estimate of best solution from here
657        estimatedDegradation = 0.0;
658        //int numberIntegerInfeasibilities=0; // without odd ones
659        numberStrongDone = 0;
660        numberUnfinished = 0;
661        numberStrongInfeasible = 0;
662        numberStrongIterations = 0;
663
664        // We may go round this loop twice (only if we think we have solution)
665        for (int iPass = 0; iPass < 2; iPass++) {
666
667            // compute current state
668            //int numberObjectInfeasibilities; // just odd ones
669            //model->feasibleSolution(
670            //                      numberIntegerInfeasibilities,
671            //                      numberObjectInfeasibilities);
672            // Some objects may compute an estimate of best solution from here
673            estimatedDegradation = 0.0;
674            numberUnsatisfied_ = 0;
675            // initialize sum of "infeasibilities"
676            sumInfeasibilities_ = 0.0;
677            int bestPriority = COIN_INT_MAX;
678            /*
679              Scan for branching objects that indicate infeasibility. Choose the best
680              maximumStrong candidates, using priority as the first criteria, then
681              integer infeasibility.
682
683              The algorithm is to fill the choice array with a set of good candidates (by
684              infeasibility) with priority bestPriority.  Finding a candidate with
685              priority better (less) than bestPriority flushes the choice array. (This
686              serves as initialization when the first candidate is found.)
687
688              A new candidate is added to choices only if its infeasibility exceeds the
689              current max infeasibility (mostAway). When a candidate is added, it
690              replaces the candidate with the smallest infeasibility (tracked by
691              iSmallest).
692            */
693            int iSmallest = 0;
694            double mostAway = 1.0e-100;
695            for (i = 0 ; i < maximumStrong ; i++)
696                choice[i].possibleBranch = NULL ;
697            numberStrong = 0;
698            bool canDoOneHot = false;
699            for (i = 0; i < numberObjects; i++) {
700                OsiObject * object = model->modifiableObject(i);
701                int preferredWay;
702                double infeasibility = object->infeasibility(&usefulInfo, preferredWay);
703                int priorityLevel = object->priority();
704                if (hotstartSolution) {
705                    // we are doing hot start
706                    const CbcSimpleInteger * thisOne = dynamic_cast <const CbcSimpleInteger *> (object);
707                    if (thisOne) {
708                        int iColumn = thisOne->columnNumber();
709                        bool canDoThisHot = true;
710                        double targetValue = hotstartSolution[iColumn];
711                        if (saveUpper[iColumn] > saveLower[iColumn]) {
712                            double value = saveSolution[iColumn];
713                            // clean
714                            value = CoinMin(value,saveUpper[iColumn]);
715                            value = CoinMax(value,saveLower[iColumn]);
716                            if (hotstartPriorities)
717                                priorityLevel = hotstartPriorities[iColumn];
718                            //double originalLower = thisOne->originalLower();
719                            //double originalUpper = thisOne->originalUpper();
720                            // switch off if not possible
721                            if (targetValue >= saveLower[iColumn] && targetValue <= saveUpper[iColumn]) {
722                                /* priority outranks rest always if negative
723                                   otherwise can be downgraded if at correct level.
724                                   Infeasibility may be increased to choose 1.0 values first.
725                                   choose one near wanted value
726                                */
727                                if (fabs(value - targetValue) > integerTolerance) {
728                                    //if (infeasibility>0.01)
729                                    //infeasibility = fabs(1.0e6-fabs(value-targetValue));
730                                    //else
731                                    infeasibility = fabs(value - targetValue);
732                                    //if (targetValue==1.0)
733                                    //infeasibility += 1.0;
734                                    if (value > targetValue) {
735                                        preferredWay = -1;
736                                    } else {
737                                        preferredWay = 1;
738                                    }
739                                    priorityLevel = CoinAbs(priorityLevel);
740                                } else if (priorityLevel < 0) {
741                                    priorityLevel = CoinAbs(priorityLevel);
742                                    if (targetValue == saveLower[iColumn]) {
743                                        infeasibility = integerTolerance + 1.0e-12;
744                                        preferredWay = -1;
745                                    } else if (targetValue == saveUpper[iColumn]) {
746                                        infeasibility = integerTolerance + 1.0e-12;
747                                        preferredWay = 1;
748                                    } else {
749                                        // can't
750                                        priorityLevel += 10000000;
751                                        canDoThisHot = false;
752                                    }
753                                } else {
754                                    priorityLevel += 10000000;
755                                    canDoThisHot = false;
756                                }
757                            } else {
758                                // switch off if not possible
759                                canDoThisHot = false;
760                            }
761                            if (canDoThisHot)
762                                canDoOneHot = true;
763                        } else if (targetValue < saveLower[iColumn] || targetValue > saveUpper[iColumn]) {
764                        }
765                    } else {
766                        priorityLevel += 10000000;
767                    }
768                }
769                if (infeasibility) {
770                    // Increase estimated degradation to solution
771                    estimatedDegradation += CoinMin(object->upEstimate(), object->downEstimate());
772                    numberUnsatisfied_++;
773                    sumInfeasibilities_ += infeasibility;
774                    // Better priority? Flush choices.
775                    if (priorityLevel < bestPriority) {
776                        int j;
777                        iSmallest = 0;
778                        for (j = 0; j < maximumStrong; j++) {
779                            choice[j].upMovement = 0.0;
780                            delete choice[j].possibleBranch;
781                            choice[j].possibleBranch = NULL;
782                        }
783                        bestPriority = priorityLevel;
784                        mostAway = 1.0e-100;
785                        numberStrong = 0;
786                    } else if (priorityLevel > bestPriority) {
787                        continue;
788                    }
789                    // Check for suitability based on infeasibility.
790                    if (infeasibility > mostAway) {
791                        //add to list
792                        choice[iSmallest].upMovement = infeasibility;
793                        delete choice[iSmallest].possibleBranch;
794                        CbcObject * obj =
795                            dynamic_cast <CbcObject *>(object) ;
796                        assert (obj);
797                        choice[iSmallest].possibleBranch = obj->createCbcBranch(solver, &usefulInfo, preferredWay);
798                        numberStrong = CoinMax(numberStrong, iSmallest + 1);
799                        // Save which object it was
800                        choice[iSmallest].objectNumber = i;
801                        int j;
802                        iSmallest = -1;
803                        mostAway = 1.0e50;
804                        for (j = 0; j < maximumStrong; j++) {
805                            if (choice[j].upMovement < mostAway) {
806                                mostAway = choice[j].upMovement;
807                                iSmallest = j;
808                            }
809                        }
810                    }
811                }
812            }
813            if (!canDoOneHot && hotstartSolution) {
814                // switch off as not possible
815                hotstartSolution = NULL;
816                model->setHotstartSolution(NULL, NULL);
817                usefulInfo.hotstartSolution_ = NULL;
818            }
819            if (numberUnsatisfied_) {
820                // some infeasibilities - go to next steps
821#ifdef CLP_INVESTIGATE
822                if (hotstartSolution) {
823                    int k = choice[0].objectNumber;
824                    OsiObject * object = model->modifiableObject(k);
825                    const CbcSimpleInteger * thisOne = dynamic_cast <const CbcSimpleInteger *> (object);
826                    assert (thisOne);
827                    int iColumn = thisOne->columnNumber();
828                    double targetValue = hotstartSolution[iColumn];
829                    double value = saveSolution[iColumn];
830                    printf("Branch on %d has target %g (value %g) and current bounds %g and %g\n",
831                           iColumn, targetValue, value, saveLower[iColumn], saveUpper[iColumn]);
832                }
833#endif
834                break;
835            } else if (!iPass) {
836                // looks like a solution - get paranoid
837                bool roundAgain = false;
838                // get basis
839                CoinWarmStartBasis * ws = dynamic_cast<CoinWarmStartBasis*>(solver->getWarmStart());
840                if (!ws)
841                    break;
842                for (i = 0; i < numberColumns; i++) {
843                    double value = saveSolution[i];
844                    if (value < lower[i]) {
845                        saveSolution[i] = lower[i];
846                        roundAgain = true;
847                        ws->setStructStatus(i, CoinWarmStartBasis::atLowerBound);
848                    } else if (value > upper[i]) {
849                        saveSolution[i] = upper[i];
850                        roundAgain = true;
851                        ws->setStructStatus(i, CoinWarmStartBasis::atUpperBound);
852                    }
853                }
854                if (roundAgain && saveNumberStrong) {
855                    // restore basis
856                    solver->setWarmStart(ws);
857                    delete ws;
858                    solver->resolve();
859                    memcpy(saveSolution, solver->getColSolution(), numberColumns*sizeof(double));
860                    model->reserveCurrentSolution(saveSolution);
861                    if (!solver->isProvenOptimal()) {
862                        // infeasible
863                        anyAction = -2;
864                        break;
865                    }
866                } else {
867                    delete ws;
868                    break;
869                }
870            }
871        }
872        /* Some solvers can do the strong branching calculations faster if
873           they do them all at once.  At present only Clp does for ordinary
874           integers but I think this coding would be easy to modify
875        */
876        bool allNormal = true; // to say if we can do fast strong branching
877        // Say which one will be best
878        int bestChoice = 0;
879        double worstInfeasibility = 0.0;
880        for (i = 0; i < numberStrong; i++) {
881            choice[i].numIntInfeasUp = numberUnsatisfied_;
882            choice[i].numIntInfeasDown = numberUnsatisfied_;
883            choice[i].fix = 0; // say not fixed
884            if (!dynamic_cast <const CbcSimpleInteger *> (model->object(choice[i].objectNumber)))
885                allNormal = false; // Something odd so lets skip clever fast branching
886            if ( !model->object(choice[i].objectNumber)->boundBranch())
887                numberStrong = 0; // switch off
888            if ( choice[i].possibleBranch->numberBranches() > 2)
889                numberStrong = 0; // switch off
890            // Do best choice in case switched off
891            if (choice[i].upMovement > worstInfeasibility) {
892                worstInfeasibility = choice[i].upMovement;
893                bestChoice = i;
894            }
895        }
896        //if (!model->parentModel())
897        //solver->writeMps("query");
898        // If we have hit max time don't do strong branching
899        bool hitMaxTime = (model->getCurrentSeconds() >
900                            model->getDblParam(CbcModel::CbcMaximumSeconds));
901        // also give up if we are looping round too much
902        if (hitMaxTime || numberPassesLeft <= 0)
903            numberStrong = 0;
904        /*
905          Is strong branching enabled? If so, set up and do it. Otherwise, we'll
906          fall through to simple branching.
907
908          Setup for strong branching involves saving the current basis (for restoration
909          afterwards) and setting up for hot starts.
910        */
911        if (numberStrong && saveNumberStrong) {
912
913            bool solveAll = false; // set true to say look at all even if some fixed (experiment)
914            solveAll = true;
915            // worth trying if too many times
916            // Save basis
917            CoinWarmStart * ws = solver->getWarmStart();
918            // save limit
919            int saveLimit;
920            solver->getIntParam(OsiMaxNumIterationHotStart, saveLimit);
921            if (beforeSolution && saveLimit < 100)
922                solver->setIntParam(OsiMaxNumIterationHotStart, 100); // go to end
923#     ifdef COIN_HAS_CLP
924            /* If we are doing all strong branching in one go then we create new arrays
925               to store information.  If clp NULL then doing old way.
926               Going down -
927               outputSolution[2*i] is final solution.
928               outputStuff[2*i] is status (0 - finished, 1 infeas, other unknown
929               outputStuff[2*i+numberStrong] is number iterations
930               On entry newUpper[i] is new upper bound, on exit obj change
931               Going up -
932               outputSolution[2*i+1] is final solution.
933               outputStuff[2*i+1] is status (0 - finished, 1 infeas, other unknown
934               outputStuff[2*i+1+numberStrong] is number iterations
935            On entry newLower[i] is new lower bound, on exit obj change
936            */
937            ClpSimplex * clp = NULL;
938            double * newLower = NULL;
939            double * newUpper = NULL;
940            double ** outputSolution = NULL;
941            int * outputStuff = NULL;
942            // Go back to normal way if user wants it
943            if (osiclp && (osiclp->specialOptions()&16) != 0 && osiclp->specialOptions() > 0)
944                allNormal = false;
945            if (osiclp && !allNormal) {
946                // say do fast
947                int easy = 1;
948                osiclp->setHintParam(OsiDoInBranchAndCut, true, OsiHintDo, &easy) ;
949            }
950            if (osiclp && allNormal) {
951                clp = osiclp->getModelPtr();
952                // Clp - do a different way
953                newLower = new double[numberStrong];
954                newUpper = new double[numberStrong];
955                outputSolution = new double * [2*numberStrong];
956                outputStuff = new int [4*numberStrong];
957                int * which = new int[numberStrong];
958                int startFinishOptions;
959                int specialOptions = osiclp->specialOptions();
960                int clpOptions = clp->specialOptions();
961                int returnCode = 0;
962#define CRUNCH
963#ifdef CRUNCH
964                // Crunch down problem
965                int numberRows = clp->numberRows();
966                // Use dual region
967                double * rhs = clp->dualRowSolution();
968                int * whichRow = new int[3*numberRows];
969                int * whichColumn = new int[2*numberColumns];
970                int nBound;
971                ClpSimplex * small = static_cast<ClpSimplexOther *> (clp)->crunch(rhs, whichRow, whichColumn, nBound, true);
972                if (!small) {
973                    anyAction = -2;
974                    //printf("XXXX Inf by inspection\n");
975                    delete [] whichColumn;
976                    whichColumn = NULL;
977                    delete [] whichRow;
978                    whichRow = NULL;
979                    break;
980                } else {
981                    clp = small;
982                }
983#else
984                int saveLogLevel = clp->logLevel();
985                int saveMaxIts = clp->maximumIterations();
986#endif
987                clp->setLogLevel(0);
988                if ((specialOptions&1) == 0) {
989                    startFinishOptions = 0;
990                    clp->setSpecialOptions(clpOptions | (64 | 1024));
991                } else {
992                    startFinishOptions = 1 + 2 + 4;
993                    //startFinishOptions=1+4; // for moment re-factorize
994                    if ((specialOptions&4) == 0)
995                        clp->setSpecialOptions(clpOptions | (64 | 128 | 512 | 1024 | 4096));
996                    else
997                        clp->setSpecialOptions(clpOptions | (64 | 128 | 512 | 1024 | 2048 | 4096));
998                }
999                // User may want to clean up before strong branching
1000                if ((clp->specialOptions()&32) != 0) {
1001                    clp->primal(1);
1002                    if (clp->numberIterations())
1003                        model->messageHandler()->message(CBC_ITERATE_STRONG, *model->messagesPointer())
1004                        << clp->numberIterations()
1005                        << CoinMessageEol;
1006                }
1007                clp->setMaximumIterations(saveLimit);
1008#ifdef CRUNCH
1009                int * backColumn = whichColumn + numberColumns;
1010#endif
1011                for (i = 0; i < numberStrong; i++) {
1012                    int iObject = choice[i].objectNumber;
1013                    const OsiObject * object = model->object(iObject);
1014                    const CbcSimpleInteger * simple = static_cast <const CbcSimpleInteger *> (object);
1015                    int iSequence = simple->columnNumber();
1016                    newLower[i] = ceil(saveSolution[iSequence]);
1017                    newUpper[i] = floor(saveSolution[iSequence]);
1018#ifdef CRUNCH
1019                    iSequence = backColumn[iSequence];
1020                    assert (iSequence >= 0);
1021#endif
1022                    which[i] = iSequence;
1023                    outputSolution[2*i] = new double [numberColumns];
1024                    outputSolution[2*i+1] = new double [numberColumns];
1025                }
1026                //clp->writeMps("bad");
1027                returnCode = clp->strongBranching(numberStrong, which,
1028                                                  newLower, newUpper, outputSolution,
1029                                                  outputStuff, outputStuff + 2 * numberStrong, !solveAll, false,
1030                                                  startFinishOptions);
1031#ifndef CRUNCH
1032                clp->setSpecialOptions(clpOptions); // restore
1033                clp->setMaximumIterations(saveMaxIts);
1034                clp->setLogLevel(saveLogLevel);
1035#endif
1036                if (returnCode == -2) {
1037                    // bad factorization!!!
1038                    // Doing normal way
1039                    // Mark hot start
1040                    solver->markHotStart();
1041                    clp = NULL;
1042                } else {
1043#ifdef CRUNCH
1044                    // extract solution
1045                    //bool checkSol=true;
1046                    for (i = 0; i < numberStrong; i++) {
1047                        int iObject = choice[i].objectNumber;
1048                        const OsiObject * object = model->object(iObject);
1049                        const CbcSimpleInteger * simple = static_cast <const CbcSimpleInteger *> (object);
1050                        int iSequence = simple->columnNumber();
1051                        which[i] = iSequence;
1052                        double * sol = outputSolution[2*i];
1053                        double * sol2 = outputSolution[2*i+1];
1054                        //bool x=true;
1055                        //bool x2=true;
1056                        for (int iColumn = numberColumns - 1; iColumn >= 0; iColumn--) {
1057                            int jColumn = backColumn[iColumn];
1058                            if (jColumn >= 0) {
1059                                sol[iColumn] = sol[jColumn];
1060                                sol2[iColumn] = sol2[jColumn];
1061                            } else {
1062                                sol[iColumn] = saveSolution[iColumn];
1063                                sol2[iColumn] = saveSolution[iColumn];
1064                            }
1065                        }
1066                    }
1067#endif
1068                }
1069#ifdef CRUNCH
1070                delete [] whichColumn;
1071                delete [] whichRow;
1072                delete small;
1073#endif
1074                delete [] which;
1075            } else {
1076                // Doing normal way
1077                // Mark hot start
1078                solver->markHotStart();
1079            }
1080#     else      /* COIN_HAS_CLP */
1081
1082            OsiSolverInterface *clp = NULL ;
1083            double **outputSolution = NULL ;
1084            int *outputStuff = NULL ;
1085            double * newLower = NULL ;
1086            double * newUpper = NULL ;
1087
1088            solver->markHotStart();
1089
1090#     endif     /* COIN_HAS_CLP */
1091            /*
1092              Open a loop to do the strong branching LPs. For each candidate variable,
1093              solve an LP with the variable forced down, then up. If a direction turns
1094              out to be infeasible or monotonic (i.e., over the dual objective cutoff),
1095              force the objective change to be big (1.0e100). If we determine the problem
1096              is infeasible, or find a monotone variable, escape the loop.
1097
1098              TODO: The `restore bounds' part might be better encapsulated as an
1099            unbranch() method. Branching objects more exotic than simple integers
1100            or cliques might not restrict themselves to variable bounds.
1101
1102              TODO: Virtuous solvers invalidate the current solution (or give bogus
1103            results :-) when the bounds are changed out from under them. So we
1104            need to do all the work associated with finding a new solution before
1105            restoring the bounds.
1106            */
1107            for (i = 0 ; i < numberStrong ; i++) {
1108                double objectiveChange ;
1109                double newObjectiveValue = 1.0e100;
1110                // status is 0 finished, 1 infeasible and other
1111                int iStatus;
1112                /*
1113                  Try the down direction first. (Specify the initial branching alternative as
1114                  down with a call to way(-1). Each subsequent call to branch() performs the
1115                  specified branch and advances the branch object state to the next branch
1116                  alternative.)
1117                */
1118                if (!clp) {
1119                    choice[i].possibleBranch->way(-1) ;
1120                    choice[i].possibleBranch->branch() ;
1121                    bool feasible = true;
1122                    if (checkFeasibility) {
1123                        // check branching did not make infeasible
1124                        int iColumn;
1125                        int numberColumns = solver->getNumCols();
1126                        const double * columnLower = solver->getColLower();
1127                        const double * columnUpper = solver->getColUpper();
1128                        for (iColumn = 0; iColumn < numberColumns; iColumn++) {
1129                            if (columnLower[iColumn] > columnUpper[iColumn] + 1.0e-5)
1130                                feasible = false;
1131                        }
1132                    }
1133                    if (feasible) {
1134                        solver->solveFromHotStart() ;
1135                        numberStrongDone++;
1136                        numberStrongIterations += solver->getIterationCount();
1137                        /*
1138                        We now have an estimate of objective degradation that we can use for strong
1139                        branching. If we're over the cutoff, the variable is monotone up.
1140                        If we actually made it to optimality, check for a solution, and if we have
1141                        a good one, call setBestSolution to process it. Note that this may reduce the
1142                        cutoff, so we check again to see if we can declare this variable monotone.
1143                        */
1144                        if (solver->isProvenOptimal())
1145                            iStatus = 0; // optimal
1146                        else if (solver->isIterationLimitReached()
1147                                 && !solver->isDualObjectiveLimitReached())
1148                            iStatus = 2; // unknown
1149                        else
1150                            iStatus = 1; // infeasible
1151                        newObjectiveValue = solver->getObjSense() * solver->getObjValue();
1152                        choice[i].numItersDown = solver->getIterationCount();
1153                    } else {
1154                        iStatus = 1; // infeasible
1155                        newObjectiveValue = 1.0e100;
1156                        choice[i].numItersDown = 0;
1157                    }
1158                } else {
1159                    iStatus = outputStuff[2*i];
1160                    choice[i].numItersDown = outputStuff[2*numberStrong+2*i];
1161                    numberStrongDone++;
1162                    numberStrongIterations += choice[i].numItersDown;
1163                    newObjectiveValue = objectiveValue + newUpper[i];
1164                    solver->setColSolution(outputSolution[2*i]);
1165                }
1166                objectiveChange = CoinMax(newObjectiveValue  - objectiveValue_, 0.0);
1167                if (!iStatus) {
1168                    choice[i].finishedDown = true ;
1169                    if (newObjectiveValue >= model->getCutoff()) {
1170                        objectiveChange = 1.0e100; // say infeasible
1171                        numberStrongInfeasible++;
1172                    } else {
1173                        // See if integer solution
1174                        if (model->feasibleSolution(choice[i].numIntInfeasDown,
1175                                                    choice[i].numObjInfeasDown)
1176                                && model->problemFeasibility()->feasible(model, -1) >= 0) {
1177                            model->setBestSolution(CBC_STRONGSOL,
1178                                                   newObjectiveValue,
1179                                                   solver->getColSolution()) ;
1180                            // only needed for odd solvers
1181                            newObjectiveValue = solver->getObjSense() * solver->getObjValue();
1182                            objectiveChange = CoinMax(newObjectiveValue - objectiveValue_, 0.0) ;
1183                            model->setLastHeuristic(NULL);
1184                            model->incrementUsed(solver->getColSolution());
1185                            if (newObjectiveValue >= model->getCutoff()) {      //  *new* cutoff
1186                                objectiveChange = 1.0e100 ;
1187                                numberStrongInfeasible++;
1188                            }
1189                        }
1190                    }
1191                } else if (iStatus == 1) {
1192                    objectiveChange = 1.0e100 ;
1193                    numberStrongInfeasible++;
1194                } else {
1195                    // Can't say much as we did not finish
1196                    choice[i].finishedDown = false ;
1197                    numberUnfinished++;
1198                }
1199                choice[i].downMovement = objectiveChange ;
1200
1201                // restore bounds
1202                if (!clp) {
1203                    for (int j = 0; j < numberColumns; j++) {
1204                        if (saveLower[j] != lower[j])
1205                            solver->setColLower(j, saveLower[j]);
1206                        if (saveUpper[j] != upper[j])
1207                            solver->setColUpper(j, saveUpper[j]);
1208                    }
1209                }
1210                //printf("Down on %d, status is %d, obj %g its %d cost %g finished %d inf %d infobj %d\n",
1211                //     choice[i].objectNumber,iStatus,newObjectiveValue,choice[i].numItersDown,
1212                //     choice[i].downMovement,choice[i].finishedDown,choice[i].numIntInfeasDown,
1213                //     choice[i].numObjInfeasDown);
1214
1215                // repeat the whole exercise, forcing the variable up
1216                if (!clp) {
1217                    bool feasible = true;
1218                    // If odd branching then maybe just one possibility
1219                    if (choice[i].possibleBranch->numberBranchesLeft() > 0) {
1220                        choice[i].possibleBranch->branch();
1221                        if (checkFeasibility) {
1222                            // check branching did not make infeasible
1223                            int iColumn;
1224                            int numberColumns = solver->getNumCols();
1225                            const double * columnLower = solver->getColLower();
1226                            const double * columnUpper = solver->getColUpper();
1227                            for (iColumn = 0; iColumn < numberColumns; iColumn++) {
1228                                if (columnLower[iColumn] > columnUpper[iColumn] + 1.0e-5)
1229                                    feasible = false;
1230                            }
1231                        }
1232                    } else {
1233                        // second branch infeasible
1234                        feasible = false;
1235                    }
1236                    if (feasible) {
1237                        solver->solveFromHotStart() ;
1238                        numberStrongDone++;
1239                        numberStrongIterations += solver->getIterationCount();
1240                        /*
1241                        We now have an estimate of objective degradation that we can use for strong
1242                        branching. If we're over the cutoff, the variable is monotone up.
1243                        If we actually made it to optimality, check for a solution, and if we have
1244                        a good one, call setBestSolution to process it. Note that this may reduce the
1245                        cutoff, so we check again to see if we can declare this variable monotone.
1246                        */
1247                        if (solver->isProvenOptimal())
1248                            iStatus = 0; // optimal
1249                        else if (solver->isIterationLimitReached()
1250                                 && !solver->isDualObjectiveLimitReached())
1251                            iStatus = 2; // unknown
1252                        else
1253                            iStatus = 1; // infeasible
1254                        newObjectiveValue = solver->getObjSense() * solver->getObjValue();
1255                        choice[i].numItersUp = solver->getIterationCount();
1256                    } else {
1257                        iStatus = 1; // infeasible
1258                        newObjectiveValue = 1.0e100;
1259                        choice[i].numItersDown = 0;
1260                    }
1261                } else {
1262                    iStatus = outputStuff[2*i+1];
1263                    choice[i].numItersUp = outputStuff[2*numberStrong+2*i+1];
1264                    numberStrongDone++;
1265                    numberStrongIterations += choice[i].numItersUp;
1266                    newObjectiveValue = objectiveValue + newLower[i];
1267                    solver->setColSolution(outputSolution[2*i+1]);
1268                }
1269                objectiveChange = CoinMax(newObjectiveValue  - objectiveValue_, 0.0);
1270                if (!iStatus) {
1271                    choice[i].finishedUp = true ;
1272                    if (newObjectiveValue >= model->getCutoff()) {
1273                        objectiveChange = 1.0e100; // say infeasible
1274                        numberStrongInfeasible++;
1275                    } else {
1276                        // See if integer solution
1277                        if (model->feasibleSolution(choice[i].numIntInfeasUp,
1278                                                    choice[i].numObjInfeasUp)
1279                                && model->problemFeasibility()->feasible(model, -1) >= 0) {
1280                            model->setBestSolution(CBC_STRONGSOL,
1281                                                   newObjectiveValue,
1282                                                   solver->getColSolution()) ;
1283                            // only needed for odd solvers
1284                            newObjectiveValue = solver->getObjSense() * solver->getObjValue();
1285                            objectiveChange = CoinMax(newObjectiveValue - objectiveValue_, 0.0) ;
1286                            model->setLastHeuristic(NULL);
1287                            model->incrementUsed(solver->getColSolution());
1288                            if (newObjectiveValue >= model->getCutoff()) {      //  *new* cutoff
1289                                objectiveChange = 1.0e100 ;
1290                                numberStrongInfeasible++;
1291                            }
1292                        }
1293                    }
1294                } else if (iStatus == 1) {
1295                    objectiveChange = 1.0e100 ;
1296                    numberStrongInfeasible++;
1297                } else {
1298                    // Can't say much as we did not finish
1299                    choice[i].finishedUp = false ;
1300                    numberUnfinished++;
1301                }
1302                choice[i].upMovement = objectiveChange ;
1303
1304                // restore bounds
1305                if (!clp) {
1306                    for (int j = 0; j < numberColumns; j++) {
1307                        if (saveLower[j] != lower[j])
1308                            solver->setColLower(j, saveLower[j]);
1309                        if (saveUpper[j] != upper[j])
1310                            solver->setColUpper(j, saveUpper[j]);
1311                    }
1312                }
1313
1314                //printf("Up on %d, status is %d, obj %g its %d cost %g finished %d inf %d infobj %d\n",
1315                //     choice[i].objectNumber,iStatus,newObjectiveValue,choice[i].numItersUp,
1316                //     choice[i].upMovement,choice[i].finishedUp,choice[i].numIntInfeasUp,
1317                //     choice[i].numObjInfeasUp);
1318
1319                /*
1320                  End of evaluation for this candidate variable. Possibilities are:
1321                  * Both sides below cutoff; this variable is a candidate for branching.
1322                  * Both sides infeasible or above the objective cutoff: no further action
1323                  here. Break from the evaluation loop and assume the node will be purged
1324                  by the caller.
1325                  * One side below cutoff: Install the branch (i.e., fix the variable). Break
1326                  from the evaluation loop and assume the node will be reoptimised by the
1327                  caller.
1328                */
1329                // reset
1330                choice[i].possibleBranch->resetNumberBranchesLeft();
1331                if (choice[i].upMovement < 1.0e100) {
1332                    if (choice[i].downMovement < 1.0e100) {
1333                        // feasible - no action
1334                    } else {
1335                        // up feasible, down infeasible
1336                        anyAction = -1;
1337                        //printf("Down infeasible for choice %d sequence %d\n",i,
1338                        // model->object(choice[i].objectNumber)->columnNumber());
1339                        if (!solveAll) {
1340                            choice[i].possibleBranch->way(1);
1341                            choice[i].possibleBranch->branch();
1342                            break;
1343                        } else {
1344                            choice[i].fix = 1;
1345                        }
1346                    }
1347                } else {
1348                    if (choice[i].downMovement < 1.0e100) {
1349                        // down feasible, up infeasible
1350                        anyAction = -1;
1351                        //printf("Up infeasible for choice %d sequence %d\n",i,
1352                        // model->object(choice[i].objectNumber)->columnNumber());
1353                        if (!solveAll) {
1354                            choice[i].possibleBranch->way(-1);
1355                            choice[i].possibleBranch->branch();
1356                            break;
1357                        } else {
1358                            choice[i].fix = -1;
1359                        }
1360                    } else {
1361                        // neither side feasible
1362                        anyAction = -2;
1363                        //printf("Both infeasible for choice %d sequence %d\n",i,
1364                        // model->object(choice[i].objectNumber)->columnNumber());
1365                        break;
1366                    }
1367                }
1368                bool hitMaxTime = (model->getCurrentSeconds() >
1369                                    model->getDblParam(CbcModel::CbcMaximumSeconds));
1370                if (hitMaxTime) {
1371                    numberStrong = i + 1;
1372                    break;
1373                }
1374            }
1375            if (!clp) {
1376                // Delete the snapshot
1377                solver->unmarkHotStart();
1378            } else {
1379                delete [] newLower;
1380                delete [] newUpper;
1381                delete [] outputStuff;
1382                int i;
1383                for (i = 0; i < 2*numberStrong; i++)
1384                    delete [] outputSolution[i];
1385                delete [] outputSolution;
1386            }
1387            solver->setIntParam(OsiMaxNumIterationHotStart, saveLimit);
1388            // restore basis
1389            solver->setWarmStart(ws);
1390            // Unless infeasible we will carry on
1391            // But we could fix anyway
1392            if (anyAction == -1 && solveAll) {
1393                // apply and take off
1394                for (i = 0 ; i < numberStrong ; i++) {
1395                    if (choice[i].fix) {
1396                        choice[i].possibleBranch->way(choice[i].fix) ;
1397                        choice[i].possibleBranch->branch() ;
1398                    }
1399                }
1400                bool feasible = true;
1401                if (checkFeasibility) {
1402                    // check branching did not make infeasible
1403                    int iColumn;
1404                    int numberColumns = solver->getNumCols();
1405                    const double * columnLower = solver->getColLower();
1406                    const double * columnUpper = solver->getColUpper();
1407                    for (iColumn = 0; iColumn < numberColumns; iColumn++) {
1408                        if (columnLower[iColumn] > columnUpper[iColumn] + 1.0e-5)
1409                            feasible = false;
1410                    }
1411                }
1412                if (feasible) {
1413                    // can do quick optimality check
1414                    int easy = 2;
1415                    solver->setHintParam(OsiDoInBranchAndCut, true, OsiHintDo, &easy) ;
1416                    solver->resolve() ;
1417                    solver->setHintParam(OsiDoInBranchAndCut, true, OsiHintDo, NULL) ;
1418                    feasible = solver->isProvenOptimal();
1419                }
1420                if (feasible) {
1421                    memcpy(saveSolution, solver->getColSolution(), numberColumns*sizeof(double));
1422                    model->reserveCurrentSolution(saveSolution);
1423                    memcpy(saveLower, solver->getColLower(), numberColumns*sizeof(double));
1424                    memcpy(saveUpper, solver->getColUpper(), numberColumns*sizeof(double));
1425                    // Clean up all candidates whih are fixed
1426                    int numberLeft = 0;
1427                    for (i = 0 ; i < numberStrong ; i++) {
1428                        CbcStrongInfo thisChoice = choice[i];
1429                        choice[i].possibleBranch = NULL;
1430                        const OsiObject * object = model->object(thisChoice.objectNumber);
1431                        double infeasibility = object->checkInfeasibility(&usefulInfo);
1432                        if (!infeasibility) {
1433                            // take out
1434                            delete thisChoice.possibleBranch;
1435                        } else {
1436                            choice[numberLeft++] = thisChoice;
1437                        }
1438                    }
1439                    numberStrong = numberLeft;
1440                    for (; i < maximumStrong; i++) {
1441                        delete choice[i].possibleBranch;
1442                        choice[i].possibleBranch = NULL;
1443                    }
1444                    // If all fixed then round again
1445                    if (!numberLeft) {
1446                        finished = false;
1447                        numberStrong = 0;
1448                        saveNumberStrong = 0;
1449                        maximumStrong = 1;
1450                    } else {
1451                        anyAction = 0;
1452                    }
1453                    // If these two uncommented then different action
1454                    anyAction = -1;
1455                    finished = true;
1456                    //printf("some fixed but continuing %d left\n",numberLeft);
1457                } else {
1458                    anyAction = -2; // say infeasible
1459                }
1460            }
1461            delete ws;
1462            //int numberNodes = model->getNodeCount();
1463            // update number of strong iterations etc
1464            model->incrementStrongInfo(numberStrongDone, numberStrongIterations,
1465                                       anyAction == -2 ? 0 : numberStrongInfeasible, anyAction == -2);
1466
1467            /*
1468              anyAction >= 0 indicates that strong branching didn't produce any monotone
1469              variables. Sift through the candidates for the best one.
1470
1471              QUERY: Setting numberNodes looks to be a distributed noop. numberNodes is
1472              local to this code block. Perhaps should be numberNodes_ from model?
1473              Unclear what this calculation is doing.
1474            */
1475            if (anyAction >= 0) {
1476
1477                // get average cost per iteration and assume stopped ones
1478                // would stop after 50% more iterations at average cost??? !!! ???
1479                double averageCostPerIteration = 0.0;
1480                double totalNumberIterations = 1.0;
1481                int smallestNumberInfeasibilities = COIN_INT_MAX;
1482                for (i = 0; i < numberStrong; i++) {
1483                    totalNumberIterations += choice[i].numItersDown +
1484                                             choice[i].numItersUp ;
1485                    averageCostPerIteration += choice[i].downMovement +
1486                                               choice[i].upMovement;
1487                    smallestNumberInfeasibilities =
1488                        CoinMin(CoinMin(choice[i].numIntInfeasDown ,
1489                                        choice[i].numIntInfeasUp ),
1490                                smallestNumberInfeasibilities);
1491                }
1492                //if (smallestNumberInfeasibilities>=numberIntegerInfeasibilities)
1493                //numberNodes=1000000; // switch off search for better solution
1494                averageCostPerIteration /= totalNumberIterations;
1495                // all feasible - choose best bet
1496
1497                // New method does all at once so it can be more sophisticated
1498                // in deciding how to balance actions.
1499                // But it does need arrays
1500                double * changeUp = new double [numberStrong];
1501                int * numberInfeasibilitiesUp = new int [numberStrong];
1502                double * changeDown = new double [numberStrong];
1503                int * numberInfeasibilitiesDown = new int [numberStrong];
1504                CbcBranchingObject ** objects = new CbcBranchingObject * [ numberStrong];
1505                for (i = 0 ; i < numberStrong ; i++) {
1506                    int iColumn = choice[i].possibleBranch->variable() ;
1507                    model->messageHandler()->message(CBC_STRONG, *model->messagesPointer())
1508                    << i << iColumn
1509                    << choice[i].downMovement << choice[i].numIntInfeasDown
1510                    << choice[i].upMovement << choice[i].numIntInfeasUp
1511                    << choice[i].possibleBranch->value()
1512                    << CoinMessageEol;
1513                    changeUp[i] = choice[i].upMovement;
1514                    numberInfeasibilitiesUp[i] = choice[i].numIntInfeasUp;
1515                    changeDown[i] = choice[i].downMovement;
1516                    numberInfeasibilitiesDown[i] = choice[i].numIntInfeasDown;
1517                    objects[i] = choice[i].possibleBranch;
1518                }
1519                int whichObject = decision->bestBranch(objects, numberStrong, numberUnsatisfied_,
1520                                                       changeUp, numberInfeasibilitiesUp,
1521                                                       changeDown, numberInfeasibilitiesDown,
1522                                                       objectiveValue_);
1523                // move branching object and make sure it will not be deleted
1524                if (whichObject >= 0) {
1525                    branch_ = objects[whichObject];
1526                    if (model->messageHandler()->logLevel() > 3)
1527                        printf("Choosing column %d\n", choice[whichObject].possibleBranch->variable()) ;
1528                    choice[whichObject].possibleBranch = NULL;
1529                }
1530                delete [] changeUp;
1531                delete [] numberInfeasibilitiesUp;
1532                delete [] changeDown;
1533                delete [] numberInfeasibilitiesDown;
1534                delete [] objects;
1535            }
1536#     ifdef COIN_HAS_CLP
1537            if (osiclp && !allNormal) {
1538                // back to normal
1539                osiclp->setHintParam(OsiDoInBranchAndCut, true, OsiHintDo, NULL) ;
1540            }
1541#     endif
1542        }
1543        /*
1544          Simple branching. Probably just one, but we may have got here
1545          because of an odd branch e.g. a cut
1546        */
1547        else {
1548            // not strong
1549            // C) create branching object
1550            branch_ = choice[bestChoice].possibleBranch;
1551            choice[bestChoice].possibleBranch = NULL;
1552        }
1553    }
1554    // Set guessed solution value
1555    guessedObjectiveValue_ = objectiveValue_ + estimatedDegradation;
1556    //printf ("Node %d depth %d unsatisfied %d sum %g obj %g guess %g\n",
1557    //      model->getNodeCount(),depth_,numberUnsatisfied_,
1558    //      sumInfeasibilities_,objectiveValue_,guessedObjectiveValue_);
1559    /*
1560      Cleanup, then we're outta here.
1561    */
1562    if (!model->branchingMethod() || dynamicBranchingObject)
1563        delete decision;
1564
1565    for (i = 0; i < maximumStrong; i++)
1566        delete choice[i].possibleBranch;
1567    delete [] choice;
1568    delete [] saveLower;
1569    delete [] saveUpper;
1570
1571    // restore solution
1572    solver->setColSolution(saveSolution);
1573    delete [] saveSolution;
1574# ifdef COIN_HAS_CLP
1575    if (osiclp)
1576        osiclp->setSpecialOptions(saveClpOptions);
1577# endif
1578    return anyAction;
1579}
1580
1581/*
1582  Version for dynamic pseudo costs.
1583
1584  **** For now just return if anything odd
1585  later allow even if odd
1586
1587  The routine scans through the object list of the model looking for objects
1588  that indicate infeasibility. It tests each object using strong branching
1589  and selects the one with the least objective degradation.  A corresponding
1590  branching object is left attached to lastNode.
1591  This version gives preference in evaluation to variables which
1592  have not been evaluated many times.  It also uses numberStrong
1593  to say give up if last few tries have not changed incumbent.
1594  See Achterberg, Koch and Martin.
1595
1596  If strong branching is disabled, a candidate object is chosen essentially
1597  at random (whatever object ends up in pos'n 0 of the candidate array).
1598
1599  If a branching candidate is found to be monotone, bounds are set to fix the
1600  variable and the routine immediately returns (the caller is expected to
1601  reoptimize).
1602
1603  If a branching candidate is found to result in infeasibility in both
1604  directions, the routine immediately returns an indication of infeasibility.
1605
1606  Returns:  0   both branch directions are feasible
1607  -1    branching variable is monotone
1608  -2    infeasible
1609  -3   Use another method
1610
1611  For now just fix on objective from strong branching.
1612*/
1613
1614int CbcNode::chooseDynamicBranch (CbcModel *model, CbcNode *lastNode,
1615                                  OsiSolverBranch * & /*branches*/,
1616                                  int numberPassesLeft)
1617
1618{
1619    if (lastNode)
1620        depth_ = lastNode->depth_ + 1;
1621    else
1622        depth_ = 0;
1623    // Go to other choose if hot start
1624    if (model->hotstartSolution() &&
1625            (((model->moreSpecialOptions()&1024) == 0) || false))
1626        return -3;
1627    delete branch_;
1628    branch_ = NULL;
1629    OsiSolverInterface * solver = model->solver();
1630    //#define CHECK_DEBUGGER_PATH
1631#ifdef CHECK_DEBUGGER_PATH
1632    bool onOptimalPath=false;
1633    if ((model->specialOptions()&1) != 0) {
1634      const OsiRowCutDebugger *debugger = solver->getRowCutDebugger() ;
1635      if (debugger) {
1636        onOptimalPath = true;
1637      } 
1638    }
1639#endif
1640    // get information on solver type
1641    const OsiAuxInfo * auxInfo = solver->getAuxiliaryInfo();
1642    const OsiBabSolver * auxiliaryInfo = dynamic_cast<const OsiBabSolver *> (auxInfo);
1643    if (!auxiliaryInfo) {
1644        // use one from CbcModel
1645        auxiliaryInfo = model->solverCharacteristics();
1646    }
1647    int numberObjects = model->numberObjects();
1648    // If very odd set of objects then use older chooseBranch
1649    bool useOldWay = false;
1650    // point to useful information
1651    OsiBranchingInformation usefulInfo = model->usefulInformation();
1652    if (numberObjects > model->numberIntegers()) {
1653        for (int i = model->numberIntegers(); i < numberObjects; i++) {
1654            OsiObject * object = model->modifiableObject(i);
1655            CbcObject * obj =   dynamic_cast <CbcObject *>(object) ;
1656            if (!obj || !obj->optionalObject()) {
1657                double infeasibility = object->checkInfeasibility(&usefulInfo);
1658                if (infeasibility) {
1659                    useOldWay = true;
1660                    break;
1661                }
1662            } else {
1663              obj->initializeForBranching(model);
1664            }
1665        }
1666    }
1667    if ((model->specialOptions()&128) != 0)
1668        useOldWay = false; // allow
1669    // For now return if not simple
1670    if (useOldWay)
1671        return -3;
1672    // Modify useful info
1673    usefulInfo.depth_ = depth_;
1674    if ((model->specialOptions()&128) != 0) {
1675        // SOS - shadow prices
1676        int numberRows = solver->getNumRows();
1677        const double * pi = usefulInfo.pi_;
1678        double sumPi = 0.0;
1679        for (int i = 0; i < numberRows; i++)
1680            sumPi += fabs(pi[i]);
1681        sumPi /= static_cast<double> (numberRows);
1682        // and scale back
1683        sumPi *= 0.01;
1684        usefulInfo.defaultDual_ = sumPi; // switch on
1685        int numberColumns = solver->getNumCols();
1686        int size = CoinMax(numberColumns, 2 * numberRows);
1687        usefulInfo.usefulRegion_ = new double [size];
1688        CoinZeroN(usefulInfo.usefulRegion_, size);
1689        usefulInfo.indexRegion_ = new int [size];
1690        // pi may change
1691        usefulInfo.pi_ = CoinCopyOfArray(usefulInfo.pi_, numberRows);
1692    }
1693    assert (auxiliaryInfo);
1694    double cutoff = model->getCutoff();
1695    const double * lower = solver->getColLower();
1696    const double * upper = solver->getColUpper();
1697    // See if user thinks infeasible
1698    int anyAction = model->problemFeasibility()->feasible(model, 0);
1699    if (anyAction) {
1700        // will return -2 if infeasible , 0 if treat as integer
1701        return anyAction - 1;
1702    }
1703    int i;
1704    int saveStateOfSearch = model->stateOfSearch() % 10;
1705    int numberStrong = model->numberStrong();
1706    /* Ranging is switched off.
1707       The idea is that you can find out the effect of one iteration
1708       on each unsatisfied variable cheaply.  Then use this
1709       if you have not got much else to go on.
1710    */
1711    //#define RANGING
1712#ifdef RANGING
1713    // must have clp
1714#ifndef COIN_HAS_CLP
1715#  warning("Ranging switched off as not Clp");
1716#undef RANGING
1717#endif
1718    // Pass number
1719    int kPass = 0;
1720    int numberRows = solver->getNumRows();
1721#endif
1722    int numberColumns = model->getNumCols();
1723    double * saveUpper = new double[numberColumns];
1724    double * saveLower = new double[numberColumns];
1725    for (i = 0; i < numberColumns; i++) {
1726        saveLower[i] = lower[i];
1727        saveUpper[i] = upper[i];
1728    }
1729
1730    // Save solution in case heuristics need good solution later
1731
1732    double * saveSolution = new double[numberColumns];
1733    memcpy(saveSolution, solver->getColSolution(), numberColumns*sizeof(double));
1734    model->reserveCurrentSolution(saveSolution);
1735    const double * hotstartSolution = model->hotstartSolution();
1736    const int * hotstartPriorities = model->hotstartPriorities();
1737    double integerTolerance =
1738        model->getDblParam(CbcModel::CbcIntegerTolerance);
1739    if (hotstartSolution) {
1740        if ((model->moreSpecialOptions()&1024) != 0) {
1741            int nBad = 0;
1742            int nUnsat = 0;
1743            int nDiff = 0;
1744            for (int i = 0; i < numberObjects; i++) {
1745                OsiObject * object = model->modifiableObject(i);
1746                const CbcSimpleInteger * thisOne = dynamic_cast <const CbcSimpleInteger *> (object);
1747                if (thisOne) {
1748                    int iColumn = thisOne->columnNumber();
1749                    double targetValue = hotstartSolution[iColumn];
1750                    double value = saveSolution[iColumn];
1751                    if (fabs(value - floor(value + 0.5)) > 1.0e-6) {
1752                        nUnsat++;
1753#ifdef CLP_INVESTIGATE
1754                        printf("H %d is %g target %g\n", iColumn, value, targetValue);
1755#endif
1756                    } else if (fabs(targetValue - value) > 1.0e-6) {
1757                        nDiff++;
1758                    }
1759                    if (targetValue < saveLower[iColumn] ||
1760                            targetValue > saveUpper[iColumn]) {
1761#ifdef CLP_INVESTIGATE
1762                        printf("%d has target %g and current bounds %g and %g\n",
1763                               iColumn, targetValue, saveLower[iColumn], saveUpper[iColumn]);
1764#endif
1765                        nBad++;
1766                    }
1767                }
1768            }
1769#ifdef CLP_INVESTIGATE
1770            printf("Hot %d unsatisfied, %d outside limits, %d different\n",
1771                   nUnsat, nBad, nDiff);
1772#endif
1773            if (nBad) {
1774                // switch off as not possible
1775                hotstartSolution = NULL;
1776                model->setHotstartSolution(NULL, NULL);
1777                usefulInfo.hotstartSolution_ = NULL;
1778            }
1779        }
1780    }
1781    /*
1782      Get a branching decision object. Use the default dynamic decision criteria unless
1783      the user has loaded a decision method into the model.
1784    */
1785    CbcBranchDecision *decision = model->branchingMethod();
1786    if (!decision)
1787        decision = new CbcBranchDynamicDecision();
1788    int xMark = 0;
1789    // Get arrays to sort
1790    double * sort = new double[numberObjects];
1791    int * whichObject = new int[numberObjects];
1792#ifdef RANGING
1793    int xPen = 0;
1794    int * objectMark = new int[2*numberObjects+1];
1795#endif
1796    // Arrays with movements
1797    double * upEstimate = new double[numberObjects];
1798    double * downEstimate = new double[numberObjects];
1799    double estimatedDegradation = 0.0;
1800    int numberNodes = model->getNodeCount();
1801    int saveLogLevel = model->logLevel();
1802#ifdef JJF_ZERO
1803    if ((numberNodes % 500) == 0) {
1804        model->setLogLevel(6);
1805        // Get average up and down costs
1806        double averageUp = 0.0;
1807        double averageDown = 0.0;
1808        int numberUp = 0;
1809        int numberDown = 0;
1810        int i;
1811        for ( i = 0; i < numberObjects; i++) {
1812            OsiObject * object = model->modifiableObject(i);
1813            CbcSimpleIntegerDynamicPseudoCost * dynamicObject =
1814                dynamic_cast <CbcSimpleIntegerDynamicPseudoCost *>(object) ;
1815            assert(dynamicObject);
1816            int  numberUp2 = 0;
1817            int numberDown2 = 0;
1818            double up = 0.0;
1819            double down = 0.0;
1820            if (dynamicObject->numberTimesUp()) {
1821                numberUp++;
1822                averageUp += dynamicObject->upDynamicPseudoCost();
1823                numberUp2 += dynamicObject->numberTimesUp();
1824                up = dynamicObject->upDynamicPseudoCost();
1825            }
1826            if (dynamicObject->numberTimesDown()) {
1827                numberDown++;
1828                averageDown += dynamicObject->downDynamicPseudoCost();
1829                numberDown2 += dynamicObject->numberTimesDown();
1830                down = dynamicObject->downDynamicPseudoCost();
1831            }
1832            if (numberUp2 || numberDown2)
1833                printf("col %d - up %d times cost %g, - down %d times cost %g\n",
1834                       dynamicObject->columnNumber(), numberUp2, up, numberDown2, down);
1835        }
1836        if (numberUp)
1837            averageUp /= static_cast<double> (numberUp);
1838        else
1839            averageUp = 1.0;
1840        if (numberDown)
1841            averageDown /= static_cast<double> (numberDown);
1842        else
1843            averageDown = 1.0;
1844        printf("total - up %d vars average %g, - down %d vars average %g\n",
1845               numberUp, averageUp, numberDown, averageDown);
1846    }
1847#endif
1848    int numberBeforeTrust = model->numberBeforeTrust();
1849    // May go round twice if strong branching fixes all local candidates
1850    bool finished = false;
1851    int numberToFix = 0;
1852# ifdef COIN_HAS_CLP
1853    OsiClpSolverInterface * osiclp = dynamic_cast< OsiClpSolverInterface*> (solver);
1854    int saveClpOptions = 0;
1855    if (osiclp) {
1856        // for faster hot start
1857        saveClpOptions = osiclp->specialOptions();
1858        osiclp->setSpecialOptions(saveClpOptions | 8192);
1859    }
1860# else
1861    OsiSolverInterface *osiclp = NULL ;
1862# endif
1863    //const CglTreeProbingInfo * probingInfo = NULL; //model->probingInfo();
1864    // Old code left in with DEPRECATED_STRATEGY
1865    assert (model->searchStrategy() == -1 ||
1866            model->searchStrategy() == 1 ||
1867            model->searchStrategy() == 2);
1868#ifdef DEPRECATED_STRATEGY
1869    int saveSearchStrategy2 = model->searchStrategy();
1870#endif
1871    // Get average up and down costs
1872    {
1873        double averageUp = 0.0;
1874        double averageDown = 0.0;
1875        int numberUp = 0;
1876        int numberDown = 0;
1877        int i;
1878        for ( i = 0; i < numberObjects; i++) {
1879            OsiObject * object = model->modifiableObject(i);
1880            CbcSimpleIntegerDynamicPseudoCost * dynamicObject =
1881                dynamic_cast <CbcSimpleIntegerDynamicPseudoCost *>(object) ;
1882            if (dynamicObject) {
1883                if (dynamicObject->numberTimesUp()) {
1884                    numberUp++;
1885                    averageUp += dynamicObject->upDynamicPseudoCost();
1886                }
1887                if (dynamicObject->numberTimesDown()) {
1888                    numberDown++;
1889                    averageDown += dynamicObject->downDynamicPseudoCost();
1890                }
1891            }
1892        }
1893        if (numberUp)
1894            averageUp /= static_cast<double> (numberUp);
1895        else
1896            averageUp = 1.0;
1897        if (numberDown)
1898            averageDown /= static_cast<double> (numberDown);
1899        else
1900            averageDown = 1.0;
1901        for ( i = 0; i < numberObjects; i++) {
1902            OsiObject * object = model->modifiableObject(i);
1903            CbcSimpleIntegerDynamicPseudoCost * dynamicObject =
1904                dynamic_cast <CbcSimpleIntegerDynamicPseudoCost *>(object) ;
1905            if (dynamicObject) {
1906                if (!dynamicObject->numberTimesUp())
1907                    dynamicObject->setUpDynamicPseudoCost(averageUp);
1908                if (!dynamicObject->numberTimesDown())
1909                    dynamicObject->setDownDynamicPseudoCost(averageDown);
1910            }
1911        }
1912    }
1913    /*
1914      1 strong
1915      2 no strong
1916      3 strong just before solution
1917      4 no strong just before solution
1918      5 strong first time or before solution
1919      6 strong first time
1920    */
1921    int useShadow = model->moreSpecialOptions() & 7;
1922    if (useShadow > 2) {
1923        if (model->getSolutionCount()) {
1924            if (numberNodes || useShadow < 5) {
1925                useShadow = 0;
1926                // zap pseudo shadow prices
1927                model->pseudoShadow(-1);
1928                // and switch off
1929                model->setMoreSpecialOptions(model->moreSpecialOptions()&(~1023));
1930            } else {
1931                useShadow = 1;
1932            }
1933        } else if (useShadow < 5) {
1934            useShadow -= 2;
1935        } else {
1936            useShadow = 1;
1937        }
1938    }
1939    if (useShadow) {
1940        // pseudo shadow prices
1941        model->pseudoShadow((model->moreSpecialOptions() >> 3)&63);
1942    }
1943#ifdef DEPRECATED_STRATEGY
1944    { // in for tabbing
1945    } else if (saveSearchStrategy2 < 1999) {
1946        // pseudo shadow prices
1947        model->pseudoShadow(NULL, NULL);
1948    } else if (saveSearchStrategy2 < 2999) {
1949        // leave old ones
1950    } else if (saveSearchStrategy2 < 3999) {
1951        // pseudo shadow prices at root
1952        if (!numberNodes)
1953            model->pseudoShadow(NULL, NULL);
1954    } else {
1955        abort();
1956    }
1957    if (saveSearchStrategy2 >= 0)
1958        saveSearchStrategy2 = saveSearchStrategy2 % 1000;
1959    if (saveSearchStrategy2 == 999)
1960        saveSearchStrategy2 = -1;
1961    int saveSearchStrategy = saveSearchStrategy2 < 99 ? saveSearchStrategy2 : saveSearchStrategy2 - 100;
1962#endif //DEPRECATED_STRATEGY
1963    int numberNotTrusted = 0;
1964    int numberStrongDone = 0;
1965    int numberUnfinished = 0;
1966    int numberStrongInfeasible = 0;
1967    int numberStrongIterations = 0;
1968    int strongType=0;
1969#define DO_ALL_AT_ROOT
1970#ifdef DO_ALL_AT_ROOT
1971    int saveSatisfiedVariables=0;
1972    int saveNumberToDo=0;
1973#endif
1974    // so we can save lots of stuff
1975    CbcStrongInfo choice;
1976    CbcDynamicPseudoCostBranchingObject * choiceObject = NULL;
1977    if (model->allDynamic()) {
1978        CbcSimpleIntegerDynamicPseudoCost * object = NULL;
1979        choiceObject = new CbcDynamicPseudoCostBranchingObject(model, 0, -1, 0.5, object);
1980    }
1981    choice.possibleBranch = choiceObject;
1982    numberPassesLeft = CoinMax(numberPassesLeft, 2);
1983    /* How dogged to be in strong branching
1984       0 - default
1985       1 - go to end on first time
1986       2 - always go to end
1987     */
1988    int goToEndInStrongBranching = (model->moreSpecialOptions2()&(3*8192))>>13;
1989#ifdef COIN_HAS_NTY
1990    // 1 after, 2 strong, 3 until depth 5
1991    int orbitOption = (model->moreSpecialOptions2()&(128|256))>>7;
1992#endif
1993    //#define DEBUG_SOLUTION
1994#ifdef DEBUG_SOLUTION
1995    bool onOptimalPath=false;
1996    if ((model->specialOptions()&1) != 0) {
1997      const OsiRowCutDebugger *debugger = model->continuousSolver()->getRowCutDebugger() ;
1998      if (debugger) {
1999        const OsiRowCutDebugger *debugger2 = model->solver()->getRowCutDebugger() ;
2000        printf("On optimal in CbcNode %s\n",debugger2 ? "" : "but bad cuts");
2001        onOptimalPath=true;
2002      }
2003    }
2004#endif
2005    while (!finished) {
2006        numberPassesLeft--;
2007        finished = true;
2008        decision->initialize(model);
2009        // Some objects may compute an estimate of best solution from here
2010        estimatedDegradation = 0.0;
2011        numberToFix = 0;
2012        int numberToDo = 0;
2013        int iBestNot = -1;
2014        int iBestGot = -1;
2015        double best = 0.0;
2016        numberNotTrusted = 0;
2017        numberStrongDone = 0;
2018        numberUnfinished = 0;
2019        numberStrongInfeasible = 0;
2020        numberStrongIterations = 0;
2021#ifdef RANGING
2022        int * which = objectMark + numberObjects + 1;
2023        int neededPenalties;
2024        int optionalPenalties;
2025#endif
2026        // We may go round this loop three times (only if we think we have solution)
2027        for (int iPass = 0; iPass < 3; iPass++) {
2028
2029            // Some objects may compute an estimate of best solution from here
2030            estimatedDegradation = 0.0;
2031            numberUnsatisfied_ = 0;
2032            // initialize sum of "infeasibilities"
2033            sumInfeasibilities_ = 0.0;
2034            int bestPriority = COIN_INT_MAX;
2035#ifdef JJF_ZERO
2036            int number01 = 0;
2037            const cliqueEntry * entry = NULL;
2038            const int * toZero = NULL;
2039            const int * toOne = NULL;
2040            const int * backward = NULL;
2041            int numberUnsatisProbed = 0;
2042            int numberUnsatisNotProbed = 0; // 0-1
2043            if (probingInfo) {
2044                number01 = probingInfo->numberIntegers();
2045                entry = probingInfo->fixEntries();
2046                toZero = probingInfo->toZero();
2047                toOne = probingInfo->toOne();
2048                backward = probingInfo->backward();
2049                if (!toZero[number01] || number01 < numberObjects || true) {
2050                    // no info
2051                    probingInfo = NULL;
2052                }
2053            }
2054#endif
2055            /*
2056              Scan for branching objects that indicate infeasibility. Choose candidates
2057              using priority as the first criteria, then integer infeasibility.
2058
2059              The algorithm is to fill the array with a set of good candidates (by
2060              infeasibility) with priority bestPriority.  Finding a candidate with
2061              priority better (less) than bestPriority flushes the choice array. (This
2062              serves as initialization when the first candidate is found.)
2063
2064            */
2065            numberToDo = 0;
2066#ifdef RANGING
2067            neededPenalties = 0;
2068            optionalPenalties = numberObjects;
2069#endif
2070            iBestNot = -1;
2071            double bestNot = 0.0;
2072            iBestGot = -1;
2073            best = 0.0;
2074            /* Problem type as set by user or found by analysis.  This will be extended
2075            0 - not known
2076            1 - Set partitioning <=
2077            2 - Set partitioning ==
2078            3 - Set covering
2079            4 - all +- 1 or all +1 and odd
2080            */
2081            int problemType = model->problemType();
2082            bool canDoOneHot = false;
2083            for (i = 0; i < numberObjects; i++) {
2084                OsiObject * object = model->modifiableObject(i);
2085                CbcSimpleIntegerDynamicPseudoCost * dynamicObject =
2086                    dynamic_cast <CbcSimpleIntegerDynamicPseudoCost *>(object) ;
2087                double infeasibility = object->checkInfeasibility(&usefulInfo);
2088                int priorityLevel = object->priority();
2089                if (hotstartSolution) {
2090                    // we are doing hot start
2091                    const CbcSimpleInteger * thisOne = dynamic_cast <const CbcSimpleInteger *> (object);
2092                    if (thisOne) {
2093                        int iColumn = thisOne->columnNumber();
2094                        bool canDoThisHot = true;
2095                        double targetValue = hotstartSolution[iColumn];
2096                        if (saveUpper[iColumn] > saveLower[iColumn]) {
2097                            double value = saveSolution[iColumn];
2098                            if (hotstartPriorities)
2099                                priorityLevel = hotstartPriorities[iColumn];
2100                            //double originalLower = thisOne->originalLower();
2101                            //double originalUpper = thisOne->originalUpper();
2102                            // switch off if not possible
2103                            if (targetValue >= saveLower[iColumn] && targetValue <= saveUpper[iColumn]) {
2104                                /* priority outranks rest always if negative
2105                                   otherwise can be downgraded if at correct level.
2106                                   Infeasibility may be increased to choose 1.0 values first.
2107                                   choose one near wanted value
2108                                */
2109                                if (fabs(value - targetValue) > integerTolerance) {
2110                                    //if (infeasibility>0.01)
2111                                    //infeasibility = fabs(1.0e6-fabs(value-targetValue));
2112                                    //else
2113                                    infeasibility = fabs(value - targetValue);
2114                                    priorityLevel = CoinAbs(priorityLevel);
2115                                } else if (priorityLevel < 0) {
2116                                    priorityLevel = CoinAbs(priorityLevel);
2117                                    if (targetValue == saveLower[iColumn] ||
2118                                    targetValue == saveUpper[iColumn]) {
2119                                        infeasibility = integerTolerance + 1.0e-12;
2120                                    } else {
2121                                        // can't
2122                                        priorityLevel += 10000000;
2123                                        canDoThisHot = false;
2124                                    }
2125                                } else {
2126                                    priorityLevel += 10000000;
2127                                    canDoThisHot = false;
2128                                }
2129                            } else {
2130                                // switch off if not possible
2131                                canDoThisHot = false;
2132                            }
2133                            if (canDoThisHot)
2134                                canDoOneHot = true;
2135                        } else if (targetValue < saveLower[iColumn] || targetValue > saveUpper[iColumn]) {
2136                        }
2137                    } else {
2138                        priorityLevel += 10000000;
2139                    }
2140                }
2141#define ZERO_ONE 0
2142#define ZERO_FAKE 1.0e20;
2143#if ZERO_ONE==1
2144                // branch on 0-1 first (temp)
2145                if (fabs(saveSolution[dynamicObject->columnNumber()]) < 1.0)
2146                    priorityLevel--;
2147#endif
2148#if ZERO_ONE==2
2149                if (fabs(saveSolution[dynamicObject->columnNumber()]) < 1.0)
2150                    infeasibility *= ZERO_FAKE;
2151#endif
2152                if (infeasibility) {
2153                    int iColumn = numberColumns + i;
2154                    bool gotDown = false;
2155                    int numberThisDown = 0;
2156                    bool gotUp = false;
2157                    int numberThisUp = 0;
2158                    double downGuess = object->downEstimate();
2159                    double upGuess = object->upEstimate();
2160                    if (dynamicObject) {
2161                        // Use this object's numberBeforeTrust
2162                        int numberBeforeTrustThis = dynamicObject->numberBeforeTrust();
2163                        iColumn = dynamicObject->columnNumber();
2164                        gotDown = false;
2165                        numberThisDown = dynamicObject->numberTimesDown();
2166                        if (numberThisDown >= numberBeforeTrustThis)
2167                            gotDown = true;
2168                        gotUp = false;
2169                        numberThisUp = dynamicObject->numberTimesUp();
2170                        if (numberThisUp >= numberBeforeTrustThis)
2171                            gotUp = true;
2172                        if (!depth_ && false) {
2173                            // try closest to 0.5
2174                            double part = saveSolution[iColumn] - floor(saveSolution[iColumn]);
2175                            infeasibility = fabs(0.5 - part);
2176                        }
2177                        if (problemType > 0 && problemType < 4 && false) {
2178                            // try closest to 0.5
2179                            double part = saveSolution[iColumn] - floor(saveSolution[iColumn]);
2180                            infeasibility = 0.5 - fabs(0.5 - part);
2181                        }
2182#ifdef JJF_ZERO
2183                        if (probingInfo) {
2184                            int iSeq = backward[iColumn];
2185                            assert (iSeq >= 0);
2186                            infeasibility = 1.0 + (toZero[iSeq+1] - toZero[iSeq]) +
2187                                            5.0 * CoinMin(toOne[iSeq] - toZero[iSeq], toZero[iSeq+1] - toOne[iSeq]);
2188                            if (toZero[iSeq+1] > toZero[iSeq]) {
2189                                numberUnsatisProbed++;
2190                            } else {
2191                                numberUnsatisNotProbed++;
2192                            }
2193                        }
2194#endif
2195                    } else {
2196                        // see if SOS
2197                        CbcSOS * sosObject =
2198                            dynamic_cast <CbcSOS *>(object) ;
2199                        if (sosObject) {
2200                            gotDown = false;
2201                            numberThisDown = sosObject->numberTimesDown();
2202                            if (numberThisDown >= numberBeforeTrust)
2203                                gotDown = true;
2204                            gotUp = false;
2205                            numberThisUp = sosObject->numberTimesUp();
2206                            if (numberThisUp >= numberBeforeTrust)
2207                                gotUp = true;
2208                        } else {
2209                            gotDown = true;
2210                            numberThisDown = 999999;
2211                            downGuess = 1.0e20;
2212                            gotUp = true;
2213                            numberThisUp = 999999;
2214                            upGuess = 1.0e20;
2215                            numberPassesLeft = 0;
2216                        }
2217                    }
2218                    // Increase estimated degradation to solution
2219                    estimatedDegradation += CoinMin(downGuess, upGuess);
2220                    downEstimate[i] = downGuess;
2221                    upEstimate[i] = upGuess;
2222                    numberUnsatisfied_++;
2223                    sumInfeasibilities_ += infeasibility;
2224                    // Better priority? Flush choices.
2225                    if (priorityLevel < bestPriority) {
2226                        numberToDo = 0;
2227                        bestPriority = priorityLevel;
2228                        iBestGot = -1;
2229                        best = 0.0;
2230                        numberNotTrusted = 0;
2231#ifdef RANGING
2232                        neededPenalties = 0;
2233                        optionalPenalties = numberObjects;
2234#endif
2235                    } else if (priorityLevel > bestPriority) {
2236                        continue;
2237                    }
2238                    if (!gotUp || !gotDown)
2239                        numberNotTrusted++;
2240                    // Check for suitability based on infeasibility.
2241                    if ((gotDown && gotUp) && numberStrong > 0) {
2242                        sort[numberToDo] = -infeasibility;
2243                        if (infeasibility > best) {
2244                            best = infeasibility;
2245                            iBestGot = numberToDo;
2246                        }
2247#ifdef RANGING
2248                        if (dynamicObject) {
2249                            objectMark[--optionalPenalties] = numberToDo;
2250                            which[optionalPenalties] = iColumn;
2251                        }
2252#endif
2253                    } else {
2254#ifdef RANGING
2255                        if (dynamicObject) {
2256                            objectMark[neededPenalties] = numberToDo;
2257                            which[neededPenalties++] = iColumn;
2258                        }
2259#endif
2260                        sort[numberToDo] = -10.0 * infeasibility;
2261                        if (!(numberThisUp + numberThisDown))
2262                            sort[numberToDo] *= 100.0; // make even more likely
2263                        if (iColumn < numberColumns) {
2264                            double part = saveSolution[iColumn] - floor(saveSolution[iColumn]);
2265                            if (1.0 - fabs(part - 0.5) > bestNot) {
2266                                iBestNot = numberToDo;
2267                                bestNot = 1.0 - fabs(part - 0.5);
2268                            }
2269                        } else {
2270                            // SOS
2271                            if (-sort[numberToDo] > bestNot) {
2272                                iBestNot = numberToDo;
2273                                bestNot = -sort[numberToDo];
2274                            }
2275                        }
2276                    }
2277                    if (model->messageHandler()->logLevel() > 3) {
2278                        printf("%d (%d) down %d %g up %d %g - infeas %g - sort %g solution %g\n",
2279                               i, iColumn, numberThisDown, object->downEstimate(), numberThisUp, object->upEstimate(),
2280                               infeasibility, sort[numberToDo], saveSolution[iColumn]);
2281                    }
2282                    whichObject[numberToDo++] = i;
2283                } else {
2284                    // for debug
2285                    downEstimate[i] = -1.0;
2286                    upEstimate[i] = -1.0;
2287                }
2288            }
2289            if (numberUnsatisfied_) {
2290                //if (probingInfo&&false)
2291                //printf("nunsat %d, %d probed, %d other 0-1\n",numberUnsatisfied_,
2292                // numberUnsatisProbed,numberUnsatisNotProbed);
2293                // some infeasibilities - go to next steps
2294                if (!canDoOneHot && hotstartSolution) {
2295                    // switch off as not possible
2296                    hotstartSolution = NULL;
2297                    model->setHotstartSolution(NULL, NULL);
2298                    usefulInfo.hotstartSolution_ = NULL;
2299                }
2300                break;
2301            } else if (!iPass) {
2302                // may just need resolve
2303                model->resolve(NULL, 11, saveSolution, saveLower, saveUpper);
2304                double newObjValue = solver->getObjSense()*solver->getObjValue();
2305                objectiveValue_ = CoinMax(objectiveValue_,newObjValue);
2306                if (!solver->isProvenOptimal()) {
2307                    // infeasible
2308                    anyAction = -2;
2309                    break;
2310                }
2311                // Double check looks OK - just look at rows with all integers
2312                if (model->allDynamic()) {
2313                    double * solution = CoinCopyOfArray(saveSolution, numberColumns);
2314                    for (int i = 0; i < numberColumns; i++) {
2315                        if (model->isInteger(i))
2316                            solution[i] = floor(solution[i] + 0.5);
2317                    }
2318                    int numberRows = solver->getNumRows();
2319                    double * rowActivity = new double [numberRows];
2320                    CoinZeroN(rowActivity, numberRows);
2321                    solver->getMatrixByCol()->times(solution, rowActivity);
2322                    //const double * element = model->solver()->getMatrixByCol()->getElements();
2323                    const int * row = model->solver()->getMatrixByCol()->getIndices();
2324                    const CoinBigIndex * columnStart = model->solver()->getMatrixByCol()->getVectorStarts();
2325                    const int * columnLength = model->solver()->getMatrixByCol()->getVectorLengths();
2326                    int nFree = 0;
2327                    int nFreeNon = 0;
2328                    int nFixedNon = 0;
2329                    double mostAway = 0.0;
2330                    int whichAway = -1;
2331                    const double * columnLower = solver->getColLower();
2332                    const double * columnUpper = solver->getColUpper();
2333                    for (int i = 0; i < numberColumns; i++) {
2334                        if (!model->isInteger(i)) {
2335                            // mark rows as flexible
2336                            CoinBigIndex start = columnStart[i];
2337                            CoinBigIndex end = start + columnLength[i];
2338                            for (CoinBigIndex j = start; j < end; j++) {
2339                                int iRow = row[j];
2340                                rowActivity[iRow] = COIN_DBL_MAX;
2341                            }
2342                        } else if (columnLower[i] < columnUpper[i]) {
2343                            double solutionValue = saveSolution[i];
2344                            if (fabs(solution[i] - solutionValue) > 
2345                                integerTolerance && 
2346                                (solutionValue - columnLower[i]) > 
2347                                integerTolerance &&
2348                                (columnUpper[i] - solutionValue) > 
2349                                integerTolerance) {
2350                                nFreeNon++;
2351                                if (fabs(solution[i] - saveSolution[i]) > mostAway) {
2352                                    mostAway = fabs(solution[i] - saveSolution[i]);
2353                                    whichAway = i;
2354                                }
2355                            } else {
2356                                nFree++;
2357                            }
2358                        } else if (solution[i] != saveSolution[i]) {
2359                            nFixedNon++;
2360                        }
2361                    }
2362                    const double * lower = solver->getRowLower();
2363                    const double * upper = solver->getRowUpper();
2364                    bool satisfied = true;
2365                    for (int i = 0; i < numberRows; i++) {
2366                        double value = rowActivity[i];
2367                        if (value != COIN_DBL_MAX) {
2368                            if (value > upper[i] + 1.0e-5 || value < lower[i] - 1.0e-5) {
2369                                satisfied = false;
2370                            }
2371                        }
2372                    }
2373                    delete [] rowActivity;
2374                    if (!satisfied) {
2375#ifdef CLP_INVESTIGATE
2376                        printf("%d free ok %d free off target %d fixed off target\n",
2377                               nFree, nFreeNon, nFixedNon);
2378#endif
2379                        if (nFreeNon) {
2380                            // try branching on these
2381                            delete branch_;
2382                            for (int i = 0; i < numberObjects; i++) {
2383                                OsiObject * object = model->modifiableObject(i);
2384                                CbcSimpleIntegerDynamicPseudoCost * obj =
2385                                    dynamic_cast <CbcSimpleIntegerDynamicPseudoCost *>(object) ;
2386                                assert (obj);
2387                                int iColumn = obj->columnNumber();
2388                                if (iColumn == whichAway) {
2389                                    int preferredWay = (saveSolution[iColumn] > solution[iColumn])
2390                                                       ? -1 : +1;
2391                                    usefulInfo.integerTolerance_ = 0.0;
2392                                    branch_ = obj->createCbcBranch(solver, &usefulInfo, preferredWay);
2393                                    break;
2394                                }
2395                            }
2396                            anyAction = 0;
2397                            break;
2398                        }
2399                    }
2400                    delete [] solution;
2401                }
2402            } else if (iPass == 1) {
2403                // looks like a solution - get paranoid
2404                bool roundAgain = false;
2405                // get basis
2406                CoinWarmStartBasis * ws = dynamic_cast<CoinWarmStartBasis*>(solver->getWarmStart());
2407                if (!ws)
2408                    break;
2409                double tolerance;
2410                solver->getDblParam(OsiPrimalTolerance, tolerance);
2411                for (i = 0; i < numberColumns; i++) {
2412                    double value = saveSolution[i];
2413                    if (value < lower[i] - tolerance) {
2414                        saveSolution[i] = lower[i];
2415                        roundAgain = true;
2416                        ws->setStructStatus(i, CoinWarmStartBasis::atLowerBound);
2417                    } else if (value > upper[i] + tolerance) {
2418                        saveSolution[i] = upper[i];
2419                        roundAgain = true;
2420                        ws->setStructStatus(i, CoinWarmStartBasis::atUpperBound);
2421                    }
2422                }
2423                if (roundAgain) {
2424                    // restore basis
2425                    solver->setWarmStart(ws);
2426                    solver->setColSolution(saveSolution);
2427                    delete ws;
2428                    bool takeHint;
2429                    OsiHintStrength strength;
2430                    solver->getHintParam(OsiDoDualInResolve, takeHint, strength);
2431                    solver->setHintParam(OsiDoDualInResolve, false, OsiHintDo) ;
2432                    model->resolve(NULL, 11, saveSolution, saveLower, saveUpper);
2433                    double newObjValue = solver->getObjSense()*solver->getObjValue();
2434                    objectiveValue_ = CoinMax(objectiveValue_,newObjValue);
2435                    solver->setHintParam(OsiDoDualInResolve, takeHint, strength) ;
2436                    if (!solver->isProvenOptimal()) {
2437                        // infeasible
2438                        anyAction = -2;
2439                        break;
2440                    }
2441                } else {
2442                    delete ws;
2443                    break;
2444                }
2445            }
2446        }
2447        if (anyAction == -2) {
2448            break;
2449        }
2450        // skip if solution
2451        if (!numberUnsatisfied_)
2452            break;
2453        int skipAll = (numberNotTrusted == 0 || numberToDo == 1) ? 1 : 0;
2454        bool doneHotStart = false;
2455        //DEPRECATED_STRATEGYint searchStrategy = saveSearchStrategy>=0 ? (saveSearchStrategy%10) : -1;
2456        int searchStrategy = model->searchStrategy();
2457        // But adjust depending on ratio of iterations
2458        if (searchStrategy > 0) {
2459          if (numberBeforeTrust >= /*5*/ 10 && numberBeforeTrust <= 10) {
2460                if (searchStrategy != 2) {
2461                    assert (searchStrategy == 1);
2462                    if (depth_ > 5) {
2463                        int numberIterations = model->getIterationCount();
2464                        int numberStrongIterations = model->numberStrongIterations();
2465                        if (numberStrongIterations > numberIterations + 10000) {
2466                            searchStrategy = 2;
2467                            skipAll = 1;
2468                        } else if (numberStrongIterations*4 + 1000 < numberIterations) {
2469                            searchStrategy = 3;
2470                            skipAll = 0;
2471                        }
2472                    } else {
2473                        searchStrategy = 3;
2474                        skipAll = 0;
2475                    }
2476                }
2477            }
2478        }
2479        // worth trying if too many times
2480        // Save basis
2481        CoinWarmStart * ws = NULL;
2482        // save limit
2483        int saveLimit = 0;
2484        solver->getIntParam(OsiMaxNumIterationHotStart, saveLimit);
2485        if (!numberPassesLeft)
2486            skipAll = 1;
2487        if (!skipAll) {
2488            ws = solver->getWarmStart();
2489            int limit = 100;
2490            if (!saveStateOfSearch && saveLimit < limit && saveLimit == 100)
2491                solver->setIntParam(OsiMaxNumIterationHotStart, limit);
2492        }
2493        // Say which one will be best
2494        int whichChoice = 0;
2495        int bestChoice;
2496        if (iBestGot >= 0)
2497            bestChoice = iBestGot;
2498        else
2499            bestChoice = iBestNot;
2500        assert (bestChoice >= 0);
2501        // If we have hit max time don't do strong branching
2502        bool hitMaxTime = (model->getCurrentSeconds() >
2503                            model->getDblParam(CbcModel::CbcMaximumSeconds));
2504        // also give up if we are looping round too much
2505        if (hitMaxTime || numberPassesLeft <= 0 || useShadow == 2) {
2506            int iObject = whichObject[bestChoice];
2507            OsiObject * object = model->modifiableObject(iObject);
2508            int preferredWay;
2509            object->infeasibility(&usefulInfo, preferredWay);
2510            CbcObject * obj =
2511                dynamic_cast <CbcObject *>(object) ;
2512            assert (obj);
2513            branch_ = obj->createCbcBranch(solver, &usefulInfo, preferredWay);
2514            {
2515                CbcBranchingObject * branchObj =
2516                    dynamic_cast <CbcBranchingObject *>(branch_) ;
2517                assert (branchObj);
2518                branchObj->way(preferredWay);
2519            }
2520            delete ws;
2521            ws = NULL;
2522            break;
2523        } else {
2524            // say do fast
2525            int easy = 1;
2526            if (!skipAll)
2527                solver->setHintParam(OsiDoInBranchAndCut, true, OsiHintDo, &easy) ;
2528            int iDo;
2529#define RESET_BOUNDS
2530#ifdef RANGING
2531            bool useRanging = model->allDynamic() && !skipAll;
2532            if (useRanging) {
2533                double currentObjective = solver->getObjValue() * solver->getObjSense();
2534                double gap = cutoff - currentObjective;
2535                // relax a bit
2536                gap *= 1.0000001;
2537                gap = CoinMax(1.0e-5, gap);
2538                // off penalties if too much
2539                double needed = neededPenalties;
2540                needed *= numberRows;
2541                if (numberNodes) {
2542                    if (needed > 1.0e6) {
2543                        neededPenalties = 0;
2544                    } else if (gap < 1.0e5) {
2545                        // maybe allow some not needed
2546                        int extra = static_cast<int> ((1.0e6 - needed) / numberRows);
2547                        int nStored = numberObjects - optionalPenalties;
2548                        extra = CoinMin(extra, nStored);
2549                        for (int i = 0; i < extra; i++) {
2550                            objectMark[neededPenalties] = objectMark[optionalPenalties+i];
2551                            which[neededPenalties++] = which[optionalPenalties+i];;
2552                        }
2553                    }
2554                }
2555                if (osiclp && neededPenalties) {
2556                    assert (!doneHotStart);
2557                    xPen += neededPenalties;
2558                    which--;
2559                    which[0] = neededPenalties;
2560                    osiclp->passInRanges(which);
2561                    // Mark hot start and get ranges
2562                    if (kPass) {
2563                        // until can work out why solution can go funny
2564                        int save = osiclp->specialOptions();
2565                        osiclp->setSpecialOptions(save | 256);
2566                        solver->markHotStart();
2567#ifdef RESET_BOUNDS
2568                        memcpy(saveLower,solver->getColLower(),solver->getNumCols()*sizeof(double));
2569                        memcpy(saveUpper,solver->getColUpper(),solver->getNumCols()*sizeof(double));
2570#endif
2571                        osiclp->setSpecialOptions(save);
2572                    } else {
2573                        solver->markHotStart();
2574#ifdef RESET_BOUNDS
2575                        memcpy(saveLower,solver->getColLower(),solver->getNumCols()*sizeof(double));
2576                        memcpy(saveUpper,solver->getColUpper(),solver->getNumCols()*sizeof(double));
2577#endif
2578                    }
2579                    doneHotStart = true;
2580                    xMark++;
2581                    kPass++;
2582                    osiclp->passInRanges(NULL);
2583                    const double * downCost = osiclp->upRange();
2584                    const double * upCost = osiclp->downRange();
2585                    bool problemFeasible = true;
2586                    int numberFixed = 0;
2587                    for (int i = 0; i < neededPenalties; i++) {
2588                        int j = objectMark[i];
2589                        int iObject = whichObject[j];
2590                        OsiObject * object = model->modifiableObject(iObject);
2591                        CbcSimpleIntegerDynamicPseudoCost * dynamicObject =
2592                            dynamic_cast <CbcSimpleIntegerDynamicPseudoCost *>(object) ;
2593                        // Use this object's numberBeforeTrust
2594                        int numberBeforeTrustThis = dynamicObject->numberBeforeTrust();
2595                        int iSequence = dynamicObject->columnNumber();
2596                        double value = saveSolution[iSequence];
2597                        value -= floor(value);
2598                        double upPenalty = CoinMin(upCost[i], 1.0e110) * (1.0 - value);
2599                        double downPenalty = CoinMin(downCost[i], 1.0e110) * value;
2600                        int numberThisDown = dynamicObject->numberTimesDown();
2601                        int numberThisUp = dynamicObject->numberTimesUp();
2602                        if (!numberBeforeTrustThis) {
2603                            // override
2604                            downEstimate[iObject] = downPenalty;
2605                            upEstimate[iObject] = upPenalty;
2606                            double min1 = CoinMin(downEstimate[iObject],
2607                                                  upEstimate[iObject]);
2608                            double max1 = CoinMax(downEstimate[iObject],
2609                                                  upEstimate[iObject]);
2610                            min1 = 0.8 * min1 + 0.2 * max1;
2611                            sort[j] = - min1;
2612                        } else if (numberThisDown < numberBeforeTrustThis ||
2613                                   numberThisUp < numberBeforeTrustThis) {
2614                            double invTrust = 1.0 / static_cast<double> (numberBeforeTrustThis);
2615                            if (numberThisDown < numberBeforeTrustThis) {
2616                                double fraction = numberThisDown * invTrust;
2617                                downEstimate[iObject] = fraction * downEstimate[iObject] + (1.0 - fraction) * downPenalty;
2618                            }
2619                            if (numberThisUp < numberBeforeTrustThis) {
2620                                double fraction = numberThisUp * invTrust;
2621                                upEstimate[iObject] = fraction * upEstimate[iObject] + (1.0 - fraction) * upPenalty;
2622                            }
2623                            double min1 = CoinMin(downEstimate[iObject],
2624                                                  upEstimate[iObject]);
2625                            double max1 = CoinMax(downEstimate[iObject],
2626                                                  upEstimate[iObject]);
2627                            min1 = 0.8 * min1 + 0.2 * max1;
2628                            min1 *= 10.0;
2629                            if (!(numberThisDown + numberThisUp))
2630                                min1 *= 100.0;
2631                            sort[j] = - min1;
2632                        }
2633                        // seems unreliable
2634                        if (false&&CoinMax(downPenalty, upPenalty) > gap) {
2635                            COIN_DETAIL_PRINT(printf("gap %g object %d has down range %g, up %g\n",
2636                                                     gap, i, downPenalty, upPenalty));
2637                            //sort[j] -= 1.0e50; // make more likely to be chosen
2638                            int number;
2639                            if (downPenalty > gap) {
2640                                number = dynamicObject->numberTimesDown();
2641                                if (upPenalty > gap)
2642                                    problemFeasible = false;
2643                                CbcBranchingObject * branch = dynamicObject->createCbcBranch(solver, &usefulInfo, 1);
2644                                //branch->fix(solver,saveLower,saveUpper,1);
2645                                delete branch;
2646                            } else {
2647                                number = dynamicObject->numberTimesUp();
2648                                CbcBranchingObject * branch = dynamicObject->createCbcBranch(solver, &usefulInfo, 1);
2649                                //branch->fix(solver,saveLower,saveUpper,-1);
2650                                delete branch;
2651                            }
2652                            if (number >= numberBeforeTrustThis)
2653                              dynamicObject->setNumberBeforeTrust(CoinMin(number + 1,5*numberBeforeTrust));
2654                            numberFixed++;
2655                        }
2656                        if (!numberNodes)
2657                            COIN_DETAIL_PRINT(printf("%d pen down ps %g -> %g up ps %g -> %g\n",
2658                                                     iObject, downPenalty, downPenalty, upPenalty, upPenalty));
2659                    }
2660                    if (numberFixed && problemFeasible) {
2661                        assert(doneHotStart);
2662                        solver->unmarkHotStart();
2663                        model->resolve(NULL, 11, saveSolution, saveLower, saveUpper);
2664#ifdef CHECK_DEBUGGER_PATH
2665                        if ((model->specialOptions()&1) != 0 && onOptimalPath) {
2666                          const OsiRowCutDebugger *debugger = solver->getRowCutDebugger() ;
2667                          if (!debugger) {
2668                            printf("Strong branching down on %d went off optimal path\n",iObject);
2669                            abort();
2670                          }
2671                        }
2672#endif
2673                        double newObjValue = solver->getObjSense()*solver->getObjValue();
2674                        objectiveValue_ = CoinMax(objectiveValue_,newObjValue);
2675                        solver->markHotStart();
2676#ifdef RESET_BOUNDS
2677                        memcpy(saveLower,solver->getColLower(),solver->getNumCols()*sizeof(double));
2678                        memcpy(saveUpper,solver->getColUpper(),solver->getNumCols()*sizeof(double));
2679#endif
2680                        problemFeasible = solver->isProvenOptimal();
2681                    }
2682                    if (!problemFeasible) {
2683                      COIN_DETAIL_PRINT(fprintf(stdout, "both ways infeas on ranging - code needed\n"));
2684                        anyAction = -2;
2685                        if (!choiceObject) {
2686                            delete choice.possibleBranch;
2687                            choice.possibleBranch = NULL;
2688                        }
2689                        //printf("Both infeasible for choice %d sequence %d\n",i,
2690                        // model->object(choice.objectNumber)->columnNumber());
2691                        // Delete the snapshot
2692                        solver->unmarkHotStart();
2693                        // back to normal
2694                        solver->setHintParam(OsiDoInBranchAndCut, true, OsiHintDo, NULL) ;
2695                        // restore basis
2696                        solver->setWarmStart(ws);
2697                        doneHotStart = false;
2698                        delete ws;
2699                        ws = NULL;
2700                        break;
2701                    }
2702                }
2703            }
2704#endif          /* RANGING */
2705            {
2706                int numberIterations = model->getIterationCount();
2707                //numberIterations += (model->numberExtraIterations()>>2);
2708                const int * strongInfo = model->strongInfo();
2709                //int numberDone = strongInfo[0]-strongInfo[3];
2710                int numberFixed = strongInfo[1] - strongInfo[4];
2711                int numberInfeasible = strongInfo[2] - strongInfo[5];
2712                assert (!strongInfo[3]);
2713                assert (!strongInfo[4]);
2714                assert (!strongInfo[5]);
2715                int numberStrongIterations = model->numberStrongIterations();
2716                int numberRows = solver->getNumRows();
2717                if (numberStrongIterations > numberIterations + CoinMin(100, 10*numberRows) && depth_ >= 4 && numberNodes > 100) {
2718                    if (20*numberInfeasible + 4*numberFixed < numberNodes) {
2719                        // Say never do
2720                        if (numberBeforeTrust == 10)
2721                          skipAll = -1;
2722                    }
2723                }
2724            }
2725            // make sure best will be first
2726            if (iBestGot >= 0)
2727                sort[iBestGot] = -COIN_DBL_MAX;
2728            // Actions 0 - exit for repeat, 1 resolve and try old choice,2 exit for continue
2729            if (anyAction)
2730                numberToDo = 0; // skip as we will be trying again
2731            // Sort
2732            CoinSort_2(sort, sort + numberToDo, whichObject);
2733            // Change in objective opposite infeasible
2734            double worstFeasible = 0.0;
2735            // Just first if strong off
2736            if (!numberStrong)
2737                numberToDo = CoinMin(numberToDo, 1);
2738            if (searchStrategy == 2)
2739                numberToDo = CoinMin(numberToDo, 20);
2740            iDo = 0;
2741            int saveLimit2;
2742            solver->getIntParam(OsiMaxNumIterationHotStart, saveLimit2);
2743            int numberTest = numberNotTrusted > 0 ? numberStrong : (numberStrong + 1) / 2;
2744            if (searchStrategy == 3) {
2745                // Previously decided we need strong
2746                numberTest = numberStrong;
2747            }
2748            // Try nearly always off
2749            if (skipAll >= 0) {
2750                if (searchStrategy < 2) {
2751                    //if ((numberNodes%20)!=0) {
2752                    if ((model->specialOptions()&8) == 0) {
2753                        numberTest = 0;
2754                    }
2755                    //} else {
2756                    //numberTest=2*numberStrong;
2757                    //skipAll=0;
2758                    //}
2759                }
2760            } else {
2761                // Just take first
2762                numberTest = 1;
2763            }
2764            int testDepth = (skipAll >= 0) ? 8 : 4;
2765            if (depth_ < testDepth && numberStrong) {
2766                if (searchStrategy != 2) {
2767                    int numberRows = solver->getNumRows();
2768                    // whether to do this or not is important - think
2769                    if (numberRows < 300 || numberRows + numberColumns < 2500) {
2770                        if (depth_ < 7)
2771                            numberStrong = CoinMin(3 * numberStrong, numberToDo);
2772                        if (!depth_)
2773                            numberStrong = CoinMin(6 * numberStrong, numberToDo);
2774                    }
2775                    numberTest = numberStrong;
2776                    skipAll = 0;
2777                }
2778            }
2779            // Do at least 5 strong
2780            if (numberColumns < 1000 && (depth_ < 15 || numberNodes < 1000000))
2781                numberTest = CoinMax(numberTest, 5);
2782            if ((model->specialOptions()&8) == 0) {
2783                if (skipAll) {
2784                    numberTest = 0;
2785                }
2786            } else {
2787                // do 5 as strong is fixing
2788                numberTest = CoinMax(numberTest, 5);
2789            }
2790            // see if switched off
2791            if (skipAll < 0) {
2792                numberTest = 0;
2793            }
2794            int realMaxHotIterations = 999999;
2795            if (skipAll < 0)
2796                numberToDo = 1;
2797            strongType=0;
2798#ifdef DO_ALL_AT_ROOT
2799            if (model->strongStrategy()) {
2800              int iStrategy=model->strongStrategy();
2801              int kDepth = iStrategy/100;
2802              if (kDepth)
2803                iStrategy -= 100*kDepth;
2804              else
2805                kDepth=5;
2806              double objValue = solver->getObjSense()*solver->getObjValue();
2807              double bestPossible = model->getBestPossibleObjValue();
2808              bestPossible += 1.0e-7*(1.0+fabs(bestPossible));
2809              int jStrategy = iStrategy/10;
2810              if (jStrategy) {
2811                if ((jStrategy&1)!=0&&!depth_)
2812                  strongType=2;
2813                else if ((jStrategy&2)!=0&&depth_<=kDepth)
2814                  strongType=2;
2815                else if ((jStrategy&4)!=0&&objValue<bestPossible)
2816                  strongType=2;
2817                iStrategy-=10*jStrategy;
2818              }
2819              if (!strongType) {
2820                if ((iStrategy&1)!=0&&!depth_)
2821                  strongType=1;
2822                else if ((iStrategy&2)!=0&&depth_<=kDepth)
2823                  strongType=1;
2824                else if ((iStrategy&4)!=0&&objValue<bestPossible)
2825                  strongType=1;
2826              }
2827              saveNumberToDo=numberToDo;
2828              if (strongType==2) {
2829                // add in satisfied
2830                const int * integerVariable = model->integerVariable();
2831                int numberIntegers = model->numberIntegers();
2832                if (numberIntegers==numberObjects) {
2833                  numberToDo=0;
2834                  for (int i=0;i<numberIntegers;i++) {
2835                    int iColumn=integerVariable[i];
2836                    if (saveUpper[iColumn]>saveLower[iColumn]) {
2837                      whichObject [numberToDo++]=i;
2838                    }
2839                  }
2840                  saveSatisfiedVariables=numberToDo-saveNumberToDo;
2841                } else {
2842                  strongType=1;
2843                }
2844              }
2845              if (strongType) {
2846                numberTest = numberToDo;
2847                numberStrong=numberToDo;
2848                skipAll=0;
2849                searchStrategy=0;
2850                solver->setIntParam(OsiMaxNumIterationHotStart, 100000);
2851                //printf("Strong branching type %d\n",strongType);
2852              }
2853            }
2854#endif
2855#ifdef COIN_HAS_NTY
2856            const int * orbits = NULL;
2857#endif
2858#ifdef COIN_HAS_NTY
2859            if (orbitOption==2 /* was >1*/) {
2860              CbcSymmetry * symmetryInfo = model->symmetryInfo();
2861              CbcNodeInfo * infoX = lastNode ? lastNode->nodeInfo() : NULL;
2862              bool worthTrying = false;
2863              if (infoX) {
2864                CbcNodeInfo * info = infoX;
2865                for (int i=0;i<NTY_BAD_DEPTH;i++) {
2866                  if (!info->parent()) {
2867                    worthTrying = true;
2868                    break;
2869                  }
2870                  info = info->parent();
2871                  if (info->symmetryWorked()) {
2872                    worthTrying = true;
2873                    break;
2874                  }
2875                }
2876              } else {
2877                worthTrying=true;
2878              }
2879              if (symmetryInfo && worthTrying) {
2880                symmetryInfo->ChangeBounds(solver->getColLower(),
2881                                            solver->getColUpper(),
2882                                            solver->getNumCols(),false);
2883                symmetryInfo->Compute_Symmetry();
2884                symmetryInfo->fillOrbits();
2885                orbits = symmetryInfo->whichOrbit();
2886                int iColumn=-1;
2887                if (orbits && symmetryInfo->numberUsefulOrbits()) {
2888                  bool doBranch=true;
2889                  int numberUsefulOrbits = symmetryInfo->numberUsefulOrbits();
2890                  if (numberUsefulOrbits<2) {
2891                    assert (numberUsefulOrbits);
2892                    double largest=-1.0;
2893                    for (int i=0;i<numberColumns;i++) {
2894                      if (orbits[i]>=0) {
2895                        if (saveSolution[i]>largest) {
2896                          largest=saveSolution[i];
2897                          iColumn=i;
2898                        }
2899                      }
2900                    }
2901                  } else {
2902#if COIN_HAS_NTY2 == 1
2903                    // take largest
2904                    int iOrbit=symmetryInfo->largestOrbit(solver->getColLower(),
2905                                                          solver->getColUpper());
2906                    double largest=-1.0;
2907                    for (int i=0;i<numberColumns;i++) {
2908                      if (orbits[i]==iOrbit) {
2909                        if (saveSolution[i]>largest) {
2910                          largest=saveSolution[i];
2911                          iColumn=i;
2912                        }
2913                      }
2914                    }
2915#endif
2916                    if (orbitOption==2) {
2917                      // strong
2918                      int nDo=0;
2919                      const double * lower = solver->getColLower();
2920                      const double * upper = solver->getColUpper();
2921                      const int * integerVariable = model->integerVariable();
2922                      for (int iOrbit = 0; iOrbit < numberUsefulOrbits; iOrbit++) {
2923                        double distance=1.0;
2924                        int iColumn = -1;
2925                        int numberIntegers = model->numberIntegers();
2926                        for (int j=0;j<numberIntegers;j++) {
2927                          int i=integerVariable[j];
2928                          if (orbits[i]==iOrbit &&lower[i]==0.0&&upper[i]==1.0) {
2929                            double away = fabs(saveSolution[i]-0.5);
2930                            if (away<distance&&away<0.4999) {
2931                              distance=away;
2932                              iColumn=j;
2933                            }
2934                          }
2935                        }
2936                        if (iColumn>=0) 
2937                          whichObject[nDo++]=iColumn;
2938                      }
2939                      if (nDo)
2940                        numberToDo=nDo;
2941                      doBranch=false;
2942                    } else if (orbitOption==3) {
2943                      // subset
2944                      int nDo=0;
2945                      for (int iDo = 0; iDo < numberToDo; iDo++) {
2946                        int iObject = whichObject[iDo];
2947                        OsiObject * object = model->modifiableObject(iObject);
2948                        CbcSimpleIntegerDynamicPseudoCost * dynamicObject =
2949                          dynamic_cast <CbcSimpleIntegerDynamicPseudoCost *>(object) ;
2950                        int iColumn = dynamicObject ? dynamicObject->columnNumber() : -1;
2951                        if (iColumn<0||orbits[iColumn]>=0)
2952                          whichObject[nDo++]=whichObject[iDo];
2953                      }
2954                      assert(nDo);
2955                      //printf("nDo %d\n",nDo);
2956                      numberToDo=nDo;
2957                      doBranch=false;
2958                      /* need NULL as if two in same orbit and strong branching fixes
2959                         then we may be in trouble.
2960                         Strong option should be OK as only one in set done.
2961                       */ 
2962                      orbits=NULL;
2963                    }
2964                  }
2965                  if(doBranch) {
2966                    orbitOption=0;
2967                    branch_ = new CbcOrbitalBranchingObject(model,iColumn,1,0,NULL);
2968                    if (infoX)
2969                      infoX->setSymmetryWorked();
2970                    numberToDo=0;
2971                  }
2972                }
2973              }
2974            }
2975#endif
2976            for ( iDo = 0; iDo < numberToDo; iDo++) {
2977                int iObject = whichObject[iDo];
2978                OsiObject * object = model->modifiableObject(iObject);
2979                CbcSimpleIntegerDynamicPseudoCost * dynamicObject =
2980                    dynamic_cast <CbcSimpleIntegerDynamicPseudoCost *>(object) ;
2981                int iColumn = dynamicObject ? dynamicObject->columnNumber() : numberColumns + iObject;
2982                int preferredWay;
2983                double infeasibility = object->infeasibility(&usefulInfo, preferredWay);
2984                bool feasibleSolution=false;
2985                double predictedChange=0.0;
2986                // may have become feasible
2987                if (!infeasibility) {
2988                  if(strongType!=2||solver->getColLower()[iColumn]==solver->getColUpper()[iColumn])
2989                    continue;
2990                }
2991#ifndef NDEBUG
2992                if (iColumn < numberColumns) {
2993                    const double * solution = model->testSolution();
2994                    assert (saveSolution[iColumn] == solution[iColumn]);
2995                }
2996#endif
2997                CbcSimpleInteger * obj =
2998                    dynamic_cast <CbcSimpleInteger *>(object) ;
2999                if (obj) {
3000                    if (choiceObject) {
3001                        obj->fillCreateBranch(choiceObject, &usefulInfo, preferredWay);
3002                        choiceObject->setObject(dynamicObject);
3003                    } else {
3004                        choice.possibleBranch = obj->createCbcBranch(solver, &usefulInfo, preferredWay);
3005                    }
3006                } else {
3007                    CbcObject * obj =
3008                        dynamic_cast <CbcObject *>(object) ;
3009                    assert (obj);
3010                    choice.possibleBranch = obj->createCbcBranch(solver, &usefulInfo, preferredWay);
3011                }
3012                // Save which object it was
3013                choice.objectNumber = iObject;
3014                choice.numIntInfeasUp = numberUnsatisfied_;
3015                choice.numIntInfeasDown = numberUnsatisfied_;
3016                if (strongType!=2) {
3017                  choice.upMovement = upEstimate[iObject];
3018                  choice.downMovement = downEstimate[iObject];
3019                } else {
3020                  choice.upMovement = 0.1;
3021                  choice.downMovement = 0.1;
3022                }
3023                assert (choice.upMovement >= 0.0);
3024                assert (choice.downMovement >= 0.0);
3025                choice.fix = 0; // say not fixed
3026                // see if can skip strong branching
3027                int canSkip = choice.possibleBranch->fillStrongInfo(choice);
3028                if ((numberTest <= 0 || skipAll)) {
3029                    if (iDo > 20) {
3030                        if (!choiceObject) {
3031                            delete choice.possibleBranch;
3032                            choice.possibleBranch = NULL;
3033                        }
3034                        break; // give up anyway
3035                    }
3036                }
3037                if (model->messageHandler()->logLevel() > 3 && numberBeforeTrust && dynamicObject)
3038                    dynamicObject->print(1, choice.possibleBranch->value());
3039                if (strongType)
3040                  canSkip=0;
3041                if (skipAll < 0)
3042                    canSkip = 1;
3043                if (!canSkip) {
3044                    if (!doneHotStart) {
3045                        // Mark hot start
3046                        doneHotStart = true;
3047                        solver->markHotStart();
3048#ifdef RESET_BOUNDS
3049                        memcpy(saveLower,solver->getColLower(),solver->getNumCols()*sizeof(double));
3050                        memcpy(saveUpper,solver->getColUpper(),solver->getNumCols()*sizeof(double));
3051#endif
3052                        if (!solver->isProvenOptimal()) {
3053                          skipAll=-2;
3054                          canSkip = 1;
3055                        }
3056                        xMark++;
3057                    }
3058                }
3059                if (!canSkip) {
3060                    numberTest--;
3061                    // just do a few
3062                    if (searchStrategy == 2)
3063                        solver->setIntParam(OsiMaxNumIterationHotStart, 10);
3064                    double objectiveChange ;
3065                    double newObjectiveValue = 1.0e100;
3066                    int j;
3067#ifdef COIN_HAS_CLP
3068                    int saveMaxHotIts=0;
3069                    int saveOsiClpOptions=0;
3070                    if (osiclp && goToEndInStrongBranching) {
3071                      /* How dogged to be in strong branching
3072                         0 - default
3073                         1 - go to end on first time
3074                         2 - always go to end
3075                      */
3076                      osiclp->getIntParam(OsiMaxNumIterationHotStart, saveMaxHotIts);
3077                      saveOsiClpOptions=osiclp->specialOptions();
3078                      if (goToEndInStrongBranching==2 ||
3079                          dynamicObject->numberTimesBranched()==0) {
3080                        osiclp->setIntParam(OsiMaxNumIterationHotStart, 
3081                                            10*(osiclp->getNumRows()+numberColumns));
3082                        osiclp->setSpecialOptions(saveOsiClpOptions & (~32));
3083                      }
3084                    }
3085#endif
3086                    // status is 0 finished, 1 infeasible and other
3087                    int iStatus;
3088                    /*
3089                      Try the down direction first. (Specify the initial branching alternative as
3090                      down with a call to way(-1). Each subsequent call to branch() performs the
3091                      specified branch and advances the branch object state to the next branch
3092                      alternative.)
3093                    */
3094                    choice.possibleBranch->way(-1) ;
3095                    predictedChange = choice.possibleBranch->branch() ;
3096#ifdef COIN_HAS_NTY
3097                    if (orbits) {
3098                      // can fix all in orbit
3099                      int fixOrbit = orbits[iObject];
3100                      if (fixOrbit>=0) {
3101                        //printf("fixing all in orbit %d for column %d\n",fixOrbit,iObject);
3102                        for (int i=0;i<numberColumns;i++) {
3103                          if (orbits[i]==fixOrbit)
3104                            solver->setColUpper(i,0.0);
3105                        }
3106                      }
3107                    }
3108#endif
3109                    solver->solveFromHotStart() ;
3110                    bool needHotStartUpdate = false;
3111                    numberStrongDone++;
3112                    numberStrongIterations += solver->getIterationCount();
3113                    /*
3114                      We now have an estimate of objective degradation that we can use for strong
3115                      branching. If we're over the cutoff, the variable is monotone up.
3116                      If we actually made it to optimality, check for a solution, and if we have
3117                      a good one, call setBestSolution to process it. Note that this may reduce the
3118                      cutoff, so we check again to see if we can declare this variable monotone.
3119                    */
3120                    if (solver->isProvenOptimal())
3121                        iStatus = 0; // optimal
3122                    else if (solver->isIterationLimitReached() 
3123                             && !solver->isDualObjectiveLimitReached()) {
3124                        iStatus = 2; // unknown
3125                    } else {
3126                        iStatus = 1; // infeasible
3127#ifdef CONFLICT_CUTS
3128#undef CONFLICT_CUTS
3129                        //#define CONFLICT_CUTS 2
3130#endif
3131#ifdef CONFLICT_CUTS
3132# ifdef COIN_HAS_CLP
3133                        if (osiclp&&(model->moreSpecialOptions()&4194304)!=0) {
3134                          const CbcFullNodeInfo * topOfTree =
3135                            model->topOfTree();
3136                          if (topOfTree) {
3137#if CONFLICT_CUTS==2
3138                            OsiRowCut * cut = osiclp->smallModelCut(topOfTree->lower(),
3139                                                                    topOfTree->upper(),
3140                                                                    model->numberRowsAtContinuous(),
3141                                                                    model->whichGenerator());
3142#else
3143                            OsiRowCut * cut = osiclp->modelCut(topOfTree->lower(),
3144                                                               topOfTree->upper(),
3145                                                               model->numberRowsAtContinuous(),
3146                                                               model->whichGenerator(),0);
3147#endif
3148                            if (cut) {
3149                              if (model->messageHandler()->logLevel() > 1)
3150                                printf("Conflict cut found in strong branching (%d elements)\n",
3151                                       cut->row().getNumElements());
3152                              //cut->print();
3153                              if ((model->specialOptions()&1) != 0) {
3154                                const OsiRowCutDebugger *debugger = model->continuousSolver()->getRowCutDebugger() ;
3155                                if (debugger) {
3156                                  if (debugger->invalidCut(*cut)) {
3157                                    model->continuousSolver()->applyRowCuts(1,cut);
3158                                    model->continuousSolver()->writeMps("bad");
3159                                  }
3160                                  CoinAssert (!debugger->invalidCut(*cut));
3161                                }
3162                              }
3163                              model->makeGlobalCut(cut) ;
3164                            }
3165                          }
3166                        }
3167#endif
3168#endif
3169                    }
3170                    // say infeasible if branch says so
3171                    if (predictedChange==COIN_DBL_MAX)
3172                      iStatus=1;
3173                    if (iStatus != 2 && solver->getIterationCount() >
3174                            realMaxHotIterations)
3175                        numberUnfinished++;
3176                    newObjectiveValue = solver->getObjSense() * solver->getObjValue();
3177                    choice.numItersDown = solver->getIterationCount();
3178                    objectiveChange = CoinMax(newObjectiveValue  - objectiveValue_, 0.0);
3179                    // Update branching information if wanted
3180                    CbcBranchingObject * cbcobj = dynamic_cast<CbcBranchingObject *> (choice.possibleBranch);
3181                    if (cbcobj) {
3182                        CbcObject * object = cbcobj->object();
3183                        assert (object) ;
3184                        CbcObjectUpdateData update = object->createUpdateInformation(solver, this, cbcobj);
3185                        update.objectNumber_ = choice.objectNumber;
3186                        model->addUpdateInformation(update);
3187                    } else {
3188                        decision->updateInformation( solver, this);
3189                    }
3190                    if (!iStatus) {
3191                        choice.finishedDown = true ;
3192                        if (newObjectiveValue >= cutoff) {
3193                            objectiveChange = 1.0e100; // say infeasible
3194                            numberStrongInfeasible++;
3195                        } else {
3196#define CBCNODE_TIGHTEN_BOUNDS
3197#ifdef CBCNODE_TIGHTEN_BOUNDS
3198                            // Can we tighten bounds?
3199                            if (iColumn < numberColumns && cutoff < 1.0e20
3200                                    && objectiveChange > 1.0e-5) {
3201                                double value = saveSolution[iColumn];
3202                                double down = value - floor(value-integerTolerance);
3203                                double changePer = objectiveChange / (down + 1.0e-7);
3204                                double distance = (cutoff - objectiveValue_) / changePer;
3205                                distance += 1.0e-3;
3206                                if (distance < 5.0) {
3207                                    double newLower = ceil(value - distance);
3208                                    if (newLower > saveLower[iColumn]) {
3209                                        //printf("Could increase lower bound on %d from %g to %g\n",
3210                                        //   iColumn,saveLower[iColumn],newLower);
3211                                        saveLower[iColumn] = newLower;
3212                                        solver->setColLower(iColumn, newLower);
3213                                    }
3214                                }
3215                            }
3216#endif
3217                            // See if integer solution
3218                            feasibleSolution = 
3219                              model->feasibleSolution(choice.numIntInfeasDown,
3220                                                      choice.numObjInfeasDown);
3221                            if (feasibleSolution
3222                                    && model->problemFeasibility()->feasible(model, -1) >= 0) {
3223                                if (auxiliaryInfo->solutionAddsCuts()) {
3224                                    needHotStartUpdate = true;
3225                                    solver->unmarkHotStart();
3226                                }
3227                                model->setLogLevel(saveLogLevel);
3228                                model->setBestSolution(CBC_STRONGSOL,
3229                                                       newObjectiveValue,
3230                                                       solver->getColSolution()) ;
3231                                if (needHotStartUpdate) {
3232                                    model->resolve(NULL, 11, saveSolution, saveLower, saveUpper);
3233                                    newObjectiveValue = solver->getObjSense() * solver->getObjValue();
3234                                    objectiveValue_ = CoinMax(objectiveValue_,newObjectiveValue);
3235                                    objectiveChange = CoinMax(newObjectiveValue  - objectiveValue_, 0.0);
3236                                    model->feasibleSolution(choice.numIntInfeasDown,
3237                                                            choice.numObjInfeasDown);
3238                                }
3239                                model->setLastHeuristic(NULL);
3240                                model->incrementUsed(solver->getColSolution());
3241                                cutoff = model->getCutoff();
3242                                if (newObjectiveValue >= cutoff) {      //  *new* cutoff
3243                                    objectiveChange = 1.0e100 ;
3244                                    numberStrongInfeasible++;
3245                                }
3246                            }
3247                        }
3248                    } else if (iStatus == 1) {
3249                        objectiveChange = 1.0e100 ;
3250                        numberStrongInfeasible++;
3251                    } else {
3252                        // Can't say much as we did not finish
3253                        choice.finishedDown = false ;
3254                        numberUnfinished++;
3255                    }
3256                    choice.downMovement = objectiveChange ;
3257
3258                    // restore bounds
3259                    for ( j = 0; j < numberColumns; j++) {
3260                        if (saveLower[j] != lower[j])
3261                            solver->setColLower(j, saveLower[j]);
3262                        if (saveUpper[j] != upper[j])
3263                            solver->setColUpper(j, saveUpper[j]);
3264                    }
3265                    if (needHotStartUpdate) {
3266                        needHotStartUpdate = false;
3267                        model->resolve(NULL, 11, saveSolution, saveLower, saveUpper);
3268#ifdef CHECK_DEBUGGER_PATH
3269                        if ((model->specialOptions()&1) != 0 && onOptimalPath) {
3270                          const OsiRowCutDebugger *debugger = solver->getRowCutDebugger() ;
3271                          if (!debugger) {
3272                            printf("Strong branching down on %d went off optimal path\n",iObject);
3273                            model->solver()->writeMps("query");
3274                            abort();
3275                          }
3276                        }
3277#endif
3278                        double newObjValue = solver->getObjSense()*solver->getObjValue();
3279                        objectiveValue_ = CoinMax(objectiveValue_,newObjValue);
3280                        //we may again have an integer feasible solution
3281                        int numberIntegerInfeasibilities;
3282                        int numberObjectInfeasibilities;
3283                        if (model->feasibleSolution(
3284                                    numberIntegerInfeasibilities,
3285                                    numberObjectInfeasibilities)) {
3286#ifdef BONMIN
3287                            //In this case node has become integer feasible, let us exit the loop
3288                            std::cout << "Node has become integer feasible" << std::endl;
3289                            numberUnsatisfied_ = 0;
3290                            break;
3291#endif
3292                            double objValue = solver->getObjValue();
3293                            model->setLogLevel(saveLogLevel);
3294                            model->setBestSolution(CBC_STRONGSOL,
3295                                                   objValue,
3296                                                   solver->getColSolution()) ;
3297                            model->resolve(NULL, 11, saveSolution, saveLower, saveUpper);
3298                            double newObjValue = solver->getObjSense()*solver->getObjValue();
3299                            objectiveValue_ = CoinMax(objectiveValue_,newObjValue);
3300                            cutoff = model->getCutoff();
3301                        }
3302                        solver->markHotStart();
3303#ifdef RESET_BOUNDS
3304                        memcpy(saveLower,solver->getColLower(),solver->getNumCols()*sizeof(double));
3305                        memcpy(saveUpper,solver->getColUpper(),solver->getNumCols()*sizeof(double));
3306#endif
3307                        if (!solver->isProvenOptimal()) {
3308                          skipAll=-2;
3309                          canSkip = 1;
3310                        }
3311                        xMark++;
3312                    }
3313#if 0 //def DO_ALL_AT_ROOT
3314                    if (strongType)
3315                        printf("Down on %d, status is %d, obj %g its %d cost %g finished %d inf %d infobj %d\n",
3316                               choice.objectNumber, iStatus, newObjectiveValue, choice.numItersDown,
3317                               choice.downMovement, choice.finishedDown, choice.numIntInfeasDown,
3318                               choice.numObjInfeasDown);
3319#endif
3320
3321                    // repeat the whole exercise, forcing the variable up
3322                    predictedChange=choice.possibleBranch->branch();
3323                    solver->solveFromHotStart() ;
3324#ifdef COIN_HAS_CLP
3325                    if (osiclp && goToEndInStrongBranching) {
3326                      osiclp->setIntParam(OsiMaxNumIterationHotStart, saveMaxHotIts);
3327                      osiclp->setSpecialOptions(saveOsiClpOptions);
3328                    }
3329#endif
3330                    numberStrongDone++;
3331                    numberStrongIterations += solver->getIterationCount();
3332                    /*
3333                      We now have an estimate of objective degradation that we can use for strong
3334                      branching. If we're over the cutoff, the variable is monotone up.
3335                      If we actually made it to optimality, check for a solution, and if we have
3336                      a good one, call setBestSolution to process it. Note that this may reduce the
3337                      cutoff, so we check again to see if we can declare this variable monotone.
3338                    */
3339                    if (solver->isProvenOptimal())
3340                        iStatus = 0; // optimal
3341                    else if (solver->isIterationLimitReached()
3342                             && !solver->isDualObjectiveLimitReached()) {
3343                        iStatus = 2; // unknown
3344                    } else {
3345                        iStatus = 1; // infeasible
3346#ifdef CONFLICT_CUTS
3347# ifdef COIN_HAS_CLP
3348                        if (osiclp&&(model->moreSpecialOptions()&4194304)!=0) {
3349                          const CbcFullNodeInfo * topOfTree =
3350                            model->topOfTree();
3351                          if (topOfTree) {
3352#if CONFLICT_CUTS==2
3353                            OsiRowCut * cut = osiclp->smallModelCut(topOfTree->lower(),
3354                                                                    topOfTree->upper(),
3355                                                                    model->numberRowsAtContinuous(),
3356                                                                    model->whichGenerator());
3357#else
3358                            OsiRowCut * cut = osiclp->modelCut(topOfTree->lower(),
3359                                                               topOfTree->upper(),
3360                                                               model->numberRowsAtContinuous(),
3361                                                               model->whichGenerator(),0);
3362#endif
3363                            if (cut) {
3364                              //printf("XXXXXX found conflict cut in strong branching\n");
3365                              //cut->print();
3366                              if ((model->specialOptions()&1) != 0) {
3367                                const OsiRowCutDebugger *debugger = model->continuousSolver()->getRowCutDebugger() ;
3368                                if (debugger) {
3369                                  if (debugger->invalidCut(*cut)) {
3370                                    model->continuousSolver()->applyRowCuts(1,cut);
3371                                    model->continuousSolver()->writeMps("bad");
3372                                  }
3373                                  CoinAssert (!debugger->invalidCut(*cut));
3374                                }
3375                              }
3376                              model->makeGlobalCut(cut) ;
3377                            }
3378                          }
3379                        }
3380#endif
3381#endif
3382                    }
3383                    // say infeasible if branch says so
3384                    if (predictedChange==COIN_DBL_MAX)
3385                      iStatus=1;
3386                    if (iStatus != 2 && solver->getIterationCount() >
3387                            realMaxHotIterations)
3388                        numberUnfinished++;
3389                    newObjectiveValue = solver->getObjSense() * solver->getObjValue();
3390                    choice.numItersUp = solver->getIterationCount();
3391                    objectiveChange = CoinMax(newObjectiveValue  - objectiveValue_, 0.0);
3392                    // Update branching information if wanted
3393                    cbcobj = dynamic_cast<CbcBranchingObject *> (choice.possibleBranch);
3394                    if (cbcobj) {
3395                        CbcObject * object = cbcobj->object();
3396                        assert (object) ;
3397                        CbcObjectUpdateData update = object->createUpdateInformation(solver, this, cbcobj);
3398                        update.objectNumber_ = choice.objectNumber;
3399                        model->addUpdateInformation(update);
3400                    } else {
3401                        decision->updateInformation( solver, this);
3402                    }
3403                    if (!iStatus) {
3404                        choice.finishedUp = true ;
3405                        if (newObjectiveValue >= cutoff) {
3406                            objectiveChange = 1.0e100; // say infeasible
3407                            numberStrongInfeasible++;
3408                        } else {
3409#ifdef CBCNODE_TIGHTEN_BOUNDS
3410                            // Can we tighten bounds?
3411                            if (iColumn < numberColumns && cutoff < 1.0e20
3412                                    && objectiveChange > 1.0e-5) {
3413                                double value = saveSolution[iColumn];
3414                                double up = ceil(value+integerTolerance) - value;
3415                                double changePer = objectiveChange / (up + 1.0e-7);
3416                                double distance = (cutoff - objectiveValue_) / changePer;
3417                                distance += 1.0e-3;
3418                                if (distance < 5.0) {
3419                                    double newUpper = floor(value + distance);
3420                                    if (newUpper < saveUpper[iColumn]) {
3421                                        //printf("Could decrease upper bound on %d from %g to %g\n",
3422                                        //   iColumn,saveUpper[iColumn],newUpper);
3423                                        saveUpper[iColumn] = newUpper;
3424                                        solver->setColUpper(iColumn, newUpper);
3425                                    }
3426                                }
3427                            }
3428#endif
3429                            // See if integer solution
3430                            feasibleSolution = 
3431                              model->feasibleSolution(choice.numIntInfeasUp,
3432                                                      choice.numObjInfeasUp);
3433                            if (feasibleSolution
3434                                    && model->problemFeasibility()->feasible(model, -1) >= 0) {
3435#ifdef BONMIN
3436                                std::cout << "Node has become integer feasible" << std::endl;
3437                                numberUnsatisfied_ = 0;
3438                                break;
3439#endif
3440                                if (auxiliaryInfo->solutionAddsCuts()) {
3441                                    needHotStartUpdate = true;
3442                                    solver->unmarkHotStart();
3443                                }
3444                                model->setLogLevel(saveLogLevel);
3445                                model->setBestSolution(CBC_STRONGSOL,
3446                                                       newObjectiveValue,
3447                                                       solver->getColSolution()) ;
3448                                if (choice.finishedDown) {
3449                                  double cutoff = model->getCutoff();
3450                                  double downObj = objectiveValue_
3451                                    + choice.downMovement ;
3452                                  if (downObj >= cutoff) {     
3453                                    choice.downMovement = 1.0e100 ;
3454                                    numberStrongInfeasible++;
3455                                }
3456
3457                                }
3458                                if (needHotStartUpdate) {
3459                                    model->resolve(NULL, 11, saveSolution, saveLower, saveUpper);
3460#ifdef CHECK_DEBUGGER_PATH
3461                                    if ((model->specialOptions()&1) != 0 && onOptimalPath) {
3462                                      const OsiRowCutDebugger *debugger = solver->getRowCutDebugger() ;
3463                                      if (!debugger) {
3464                                        printf("Strong branching up on %d went off optimal path\n",iObject);
3465                                        abort();
3466                                      }
3467                                    }
3468#endif
3469                                    newObjectiveValue = solver->getObjSense() * solver->getObjValue();
3470                                    objectiveValue_ = CoinMax(objectiveValue_,newObjectiveValue);
3471                                    objectiveChange = CoinMax(newObjectiveValue  - objectiveValue_, 0.0);
3472                                    model->feasibleSolution(choice.numIntInfeasDown,
3473                                                            choice.numObjInfeasDown);
3474                                }
3475                                model->setLastHeuristic(NULL);
3476                                model->incrementUsed(solver->getColSolution());
3477                                cutoff = model->getCutoff();
3478                                if (newObjectiveValue >= cutoff) {      //  *new* cutoff
3479                                    objectiveChange = 1.0e100 ;
3480                                    numberStrongInfeasible++;
3481                                }
3482                            }
3483                        }
3484                    } else if (iStatus == 1) {
3485                        objectiveChange = 1.0e100 ;
3486                        numberStrongInfeasible++;
3487                    } else {
3488                        // Can't say much as we did not finish
3489                        choice.finishedUp = false ;
3490                        numberUnfinished++;
3491                    }
3492                    choice.upMovement = objectiveChange ;
3493
3494                    // restore bounds
3495                    for ( j = 0; j < numberColumns; j++) {
3496                        if (saveLower[j] != lower[j])
3497                            solver->setColLower(j, saveLower[j]);
3498                        if (saveUpper[j] != upper[j])
3499                            solver->setColUpper(j, saveUpper[j]);
3500                    }
3501                    if (needHotStartUpdate) {
3502                        needHotStartUpdate = false;
3503                        model->resolve(NULL, 11, saveSolution, saveLower, saveUpper);
3504#ifdef CHECK_DEBUGGER_PATH
3505                        if ((model->specialOptions()&1) != 0 && onOptimalPath) {
3506                          const OsiRowCutDebugger *debugger = solver->getRowCutDebugger() ;
3507                          if (!debugger) {
3508                            printf("Strong branching up on %d went off optimal path\n",iObject);
3509                            abort();
3510                          }
3511                        }
3512#endif
3513                        double newObjValue = solver->getObjSense()*solver->getObjValue();
3514                        objectiveValue_ = CoinMax(objectiveValue_,newObjValue);
3515                        //we may again have an integer feasible solution
3516                        int numberIntegerInfeasibilities;
3517                        int numberObjectInfeasibilities;
3518                        if (model->feasibleSolution(
3519                                    numberIntegerInfeasibilities,
3520                                    numberObjectInfeasibilities)) {
3521                            double objValue = solver->getObjValue();
3522                            model->setLogLevel(saveLogLevel);
3523                            model->setBestSolution(CBC_STRONGSOL,
3524                                                   objValue,
3525                                                   solver->getColSolution()) ;
3526                            model->resolve(NULL, 11, saveSolution, saveLower, saveUpper);
3527#ifdef CHECK_DEBUGGER_PATH
3528                            if ((model->specialOptions()&1) != 0 && onOptimalPath) {
3529                              const OsiRowCutDebugger *debugger = solver->getRowCutDebugger() ;
3530                              if (!debugger) {
3531                                printf("Strong branching up on %d went off optimal path\n",iObject);
3532                                abort();
3533                              }
3534                            }
3535#endif
3536                            double newObjValue = solver->getObjSense()*solver->getObjValue();
3537                            objectiveValue_ = CoinMax(objectiveValue_,newObjValue);
3538                            cutoff = model->getCutoff();
3539                        }
3540                        solver->markHotStart();
3541#ifdef RESET_BOUNDS
3542                        memcpy(saveLower,solver->getColLower(),solver->getNumCols()*sizeof(double));
3543                        memcpy(saveUpper,solver->getColUpper(),solver->getNumCols()*sizeof(double));
3544#endif
3545                        if (!solver->isProvenOptimal()) {
3546                          skipAll=-2;
3547                          canSkip = 1;
3548                        }
3549                        xMark++;
3550                    }
3551
3552#if 0 //def DO_ALL_AT_ROOT
3553                    if (strongType)
3554                        printf("Up on %d, status is %d, obj %g its %d cost %g finished %d inf %d infobj %d\n",
3555                               choice.objectNumber, iStatus, newObjectiveValue, choice.numItersUp,
3556                               choice.upMovement, choice.finishedUp, choice.numIntInfeasUp,
3557                               choice.numObjInfeasUp);
3558#endif
3559                }
3560
3561                solver->setIntParam(OsiMaxNumIterationHotStart, saveLimit2);
3562                /*
3563                  End of evaluation for this candidate variable. Possibilities are:
3564                  * Both sides below cutoff; this variable is a candidate for branching.
3565                  * Both sides infeasible or above the objective cutoff: no further action
3566                  here. Break from the evaluation loop and assume the node will be purged
3567                  by the caller.
3568                  * One side below cutoff: Install the branch (i.e., fix the variable). Break
3569                  from the evaluation loop and assume the node will be reoptimised by the
3570                  caller.
3571                */
3572                // reset
3573                choice.possibleBranch->resetNumberBranchesLeft();
3574                if (choice.upMovement < 1.0e100) {
3575                    if (choice.downMovement < 1.0e100) {
3576                        // In case solution coming in was odd
3577                        choice.upMovement = CoinMax(0.0, choice.upMovement);
3578                        choice.downMovement = CoinMax(0.0, choice.downMovement);
3579#if ZERO_ONE==2
3580                        // branch on 0-1 first (temp)
3581                        if (fabs(choice.possibleBranch->value()) < 1.0) {
3582                            choice.upMovement *= ZERO_FAKE;
3583                            choice.downMovement *= ZERO_FAKE;
3584                        }
3585#endif
3586                        // feasible - see which best
3587                        if (!canSkip) {
3588                            if (model->messageHandler()->logLevel() > 3)
3589                                printf("sort %g downest %g upest %g ", sort[iDo], downEstimate[iObject],
3590                                       upEstimate[iObject]);
3591                            model->messageHandler()->message(CBC_STRONG, *model->messagesPointer())
3592                            << iObject << iColumn
3593                            << choice.downMovement << choice.numIntInfeasDown
3594                            << choice.upMovement << choice.numIntInfeasUp
3595                            << choice.possibleBranch->value()
3596                            << CoinMessageEol;
3597                        }
3598                        int betterWay=0;
3599                        // If was feasible (extra strong branching) skip
3600                        if (infeasibility) {
3601                            CbcBranchingObject * branchObj =
3602                                dynamic_cast <CbcBranchingObject *>(branch_) ;
3603                            if (branch_)
3604                                assert (branchObj);
3605                            betterWay = decision->betterBranch(choice.possibleBranch,
3606                                                               branchObj,
3607                                                               choice.upMovement,
3608                                                               choice.numIntInfeasUp ,
3609                                                               choice.downMovement,
3610                                                               choice.numIntInfeasDown );
3611                        }
3612                        if (betterWay) {
3613                            // C) create branching object
3614                            if (choiceObject) {
3615                                delete branch_;
3616                                branch_ = choice.possibleBranch->clone();
3617                            } else {
3618                                delete branch_;
3619                                branch_ = choice.possibleBranch;
3620                                choice.possibleBranch = NULL;
3621                            }
3622                            {
3623                                CbcBranchingObject * branchObj =
3624                                    dynamic_cast <CbcBranchingObject *>(branch_) ;
3625                                assert (branchObj);
3626                                //branchObj->way(preferredWay);
3627                                branchObj->way(betterWay);
3628                            }
3629                            bestChoice = choice.objectNumber;
3630                            whichChoice = iDo;
3631                            if (numberStrong <= 1) {
3632                                delete ws;
3633                                ws = NULL;
3634                                break;
3635                            }
3636                        } else {
3637                            if (!choiceObject) {
3638                                delete choice.possibleBranch;
3639                                choice.possibleBranch = NULL;
3640                            }
3641                            if (iDo >= 2*numberStrong) {
3642                                delete ws;
3643                                ws = NULL;
3644                                break;
3645                            }
3646                            if (!dynamicObject || dynamicObject->numberTimesUp() > 1) {
3647                                if (iDo - whichChoice >= numberStrong) {
3648                                    if (!choiceObject) {
3649                                        delete choice.possibleBranch;
3650                                        choice.possibleBranch = NULL;
3651                                    }
3652                                    break; // give up
3653                                }
3654                            } else {
3655                                if (iDo - whichChoice >= 2*numberStrong) {
3656                                    delete ws;
3657                                    ws = NULL;
3658                                    if (!choiceObject) {
3659                                        delete choice.possibleBranch;
3660                                        choice.possibleBranch = NULL;
3661                                    }
3662                                    break; // give up
3663                                }
3664                            }
3665                        }
3666                    } else {
3667                        // up feasible, down infeasible
3668                        anyAction = -1;
3669                        worstFeasible = CoinMax(worstFeasible, choice.upMovement);
3670                        model->messageHandler()->message(CBC_STRONG, *model->messagesPointer())
3671                        << iObject << iColumn
3672                        << choice.downMovement << choice.numIntInfeasDown
3673                        << choice.upMovement << choice.numIntInfeasUp
3674                        << choice.possibleBranch->value()
3675                        << CoinMessageEol;
3676                        //printf("Down infeasible for choice %d sequence %d\n",i,
3677                        // model->object(choice.objectNumber)->columnNumber());
3678                        choice.fix = 1;
3679                        numberToFix++;
3680                        choice.possibleBranch->fix(solver, saveLower, saveUpper, 1);
3681                        if (!choiceObject) {
3682                            delete choice.possibleBranch;
3683                            choice.possibleBranch = NULL;
3684                        } else {
3685                            //choiceObject = new CbcDynamicPseudoCostBranchingObject(*choiceObject);
3686                            choice.possibleBranch = choiceObject;
3687                        }
3688                        assert(doneHotStart);
3689                        solver->unmarkHotStart();
3690                        model->resolve(NULL, 11, saveSolution, saveLower, saveUpper);
3691#ifdef CHECK_DEBUGGER_PATH
3692                        if ((model->specialOptions()&1) != 0 && onOptimalPath) {
3693                          const OsiRowCutDebugger *debugger = solver->getRowCutDebugger() ;
3694                          if (!debugger) {
3695                            printf("Strong branching down on %d went off optimal path\n",iObject);
3696                            abort();
3697                          }
3698                        }
3699#endif
3700                        double newObjValue = solver->getObjSense()*solver->getObjValue();
3701                        objectiveValue_ = CoinMax(objectiveValue_,newObjValue);
3702                        bool goneInfeasible = (!solver->isProvenOptimal()||solver->isDualObjectiveLimitReached());
3703                        solver->markHotStart();
3704#ifdef RESET_BOUNDS
3705                        memcpy(saveLower,solver->getColLower(),solver->getNumCols()*sizeof(double));
3706                        memcpy(saveUpper,solver->getColUpper(),solver->getNumCols()*sizeof(double));
3707#endif
3708                        if (!solver->isProvenOptimal()) {
3709                          skipAll=-2;
3710                          canSkip = 1;
3711                        }
3712                        xMark++;
3713                        // may be infeasible (if other way stopped on iterations)
3714                        if (goneInfeasible) {
3715                            // neither side feasible
3716                            anyAction = -2;
3717                            if (!choiceObject) {
3718                                delete choice.possibleBranch;
3719                                choice.possibleBranch = NULL;
3720                            }
3721                            //printf("Both infeasible for choice %d sequence %d\n",i,
3722                            // model->object(choice.objectNumber)->columnNumber());
3723                            delete ws;
3724                            ws = NULL;
3725                            break;
3726                        }
3727                    }
3728                } else {
3729                    if (choice.downMovement < 1.0e100) {
3730                        // down feasible, up infeasible
3731                        anyAction = -1;
3732                        worstFeasible = CoinMax(worstFeasible, choice.downMovement);
3733                        model->messageHandler()->message(CBC_STRONG, *model->messagesPointer())
3734                        << iObject << iColumn
3735                        << choice.downMovement << choice.numIntInfeasDown
3736                        << choice.upMovement << choice.numIntInfeasUp
3737                        << choice.possibleBranch->value()
3738                        << CoinMessageEol;
3739                        choice.fix = -1;
3740                        numberToFix++;
3741                        choice.possibleBranch->fix(solver, saveLower, saveUpper, -1);
3742                        if (!choiceObject) {
3743                            delete choice.possibleBranch;
3744                            choice.possibleBranch = NULL;
3745                        } else {
3746                            //choiceObject = new CbcDynamicPseudoCostBranchingObject(*choiceObject);
3747                            choice.possibleBranch = choiceObject;
3748                        }
3749                        assert(doneHotStart);
3750                        solver->unmarkHotStart();
3751                        model->resolve(NULL, 11, saveSolution, saveLower, saveUpper);
3752#ifdef CHECK_DEBUGGER_PATH
3753                        if ((model->specialOptions()&1) != 0 && onOptimalPath) {
3754                          const OsiRowCutDebugger *debugger = solver->getRowCutDebugger() ;
3755                          if (!debugger) {
3756                            printf("Strong branching down on %d went off optimal path\n",iObject);
3757                            solver->writeMps("query");
3758                            abort();
3759                          }
3760                        }
3761#endif
3762                        double newObjValue = solver->getObjSense()*solver->getObjValue();
3763                        objectiveValue_ = CoinMax(objectiveValue_,newObjValue);
3764                        bool goneInfeasible = (!solver->isProvenOptimal()||solver->isDualObjectiveLimitReached());
3765                        solver->markHotStart();
3766#ifdef RESET_BOUNDS
3767                        memcpy(saveLower,solver->getColLower(),solver->getNumCols()*sizeof(double));
3768                        memcpy(saveUpper,solver->getColUpper(),solver->getNumCols()*sizeof(double));
3769#endif
3770                        if (!solver->isProvenOptimal()) {
3771                          skipAll=-2;
3772                          canSkip = 1;
3773                        }
3774                        xMark++;
3775                        // may be infeasible (if other way stopped on iterations)
3776                        if (goneInfeasible) {
3777                            // neither side feasible
3778                            anyAction = -2;
3779                            if (!choiceObject) {
3780                                delete choice.possibleBranch;
3781                                choice.possibleBranch = NULL;
3782                            }
3783                            delete ws;
3784                            ws = NULL;
3785                            break;
3786                        }
3787                    } else {
3788                        // neither side feasible
3789                        anyAction = -2;
3790                        if (!choiceObject) {
3791                            delete choice.possibleBranch;
3792                            choice.possibleBranch = NULL;
3793                        }
3794                        delete ws;
3795                        ws = NULL;
3796                        break;
3797                    }
3798                }
3799                // Check max time
3800                hitMaxTime = (model->getCurrentSeconds() >
3801                               model->getDblParam(CbcModel::CbcMaximumSeconds));
3802                if (hitMaxTime) {
3803                    // make sure rest are fast
3804                    for ( int jDo = iDo + 1; jDo < numberToDo; jDo++) {
3805                        int iObject = whichObject[iDo];
3806                        OsiObject * object = model->modifiableObject(iObject);
3807                        CbcSimpleIntegerDynamicPseudoCost * dynamicObject =
3808                            dynamic_cast <CbcSimpleIntegerDynamicPseudoCost *>(object) ;
3809                        if (dynamicObject)
3810                            dynamicObject->setNumberBeforeTrust(0);
3811                    }
3812                    numberTest = 0;
3813                }
3814                if (!choiceObject) {
3815                    delete choice.possibleBranch;
3816                }
3817            }
3818            if (model->messageHandler()->logLevel() > 3) {
3819                if (anyAction == -2) {
3820                    printf("infeasible\n");
3821                } else if (anyAction == -1) {
3822                    printf("%d fixed AND choosing %d iDo %d iChosenWhen %d numberToDo %d\n", numberToFix, bestChoice,
3823                           iDo, whichChoice, numberToDo);
3824                } else {
3825                    int iObject = whichObject[whichChoice];
3826                    OsiObject * object = model->modifiableObject(iObject);
3827                    CbcSimpleIntegerDynamicPseudoCost * dynamicObject =
3828                        dynamic_cast <CbcSimpleIntegerDynamicPseudoCost *>(object) ;
3829                    if (dynamicObject) {
3830                        int iColumn = dynamicObject->columnNumber();
3831                        printf("choosing %d (column %d) iChosenWhen %d numberToDo %d\n", bestChoice,
3832                               iColumn, whichChoice, numberToDo);
3833                    }
3834                }
3835            }
3836            if (doneHotStart) {
3837                // Delete the snapshot
3838                solver->unmarkHotStart();
3839                // back to normal
3840                solver->setHintParam(OsiDoInBranchAndCut, true, OsiHintDo, NULL) ;
3841                // restore basis
3842                solver->setWarmStart(ws);
3843            }
3844            solver->setIntParam(OsiMaxNumIterationHotStart, saveLimit);
3845            // Unless infeasible we will carry on
3846            // But we could fix anyway
3847            if (numberToFix && !hitMaxTime) {
3848                if (anyAction != -2) {
3849                    // apply and take off
3850                    bool feasible = true;
3851                    // can do quick optimality check
3852                    int easy = 2;
3853                    solver->setHintParam(OsiDoInBranchAndCut, true, OsiHintDo, &easy) ;
3854                    model->resolve(NULL, 11, saveSolution, saveLower, saveUpper) ;
3855#ifdef CHECK_DEBUGGER_PATH
3856                    if ((model->specialOptions()&1) != 0 && onOptimalPath) {
3857                      const OsiRowCutDebugger *debugger = solver->getRowCutDebugger() ;
3858                      if (!debugger) {
3859                        printf("Strong branching went off optimal path\n");
3860                        abort();
3861                      }
3862                    }
3863#endif
3864                    double newObjValue = solver->getObjSense()*solver->getObjValue();
3865                    objectiveValue_ = CoinMax(objectiveValue_,newObjValue);
3866                    solver->setHintParam(OsiDoInBranchAndCut, true, OsiHintDo, NULL) ;
3867                    feasible = solver->isProvenOptimal();
3868                    if (feasible) {
3869                        anyAction = 0;
3870                    } else {
3871                        anyAction = -2;
3872                        finished = true;
3873                    }
3874                }
3875            }
3876            // If  fixed then round again
3877            // See if candidate still possible
3878            if (branch_) {
3879                 const OsiObject * object = model->object(bestChoice);
3880                 double infeasibility = object->checkInfeasibility(&usefulInfo);
3881                 if (!infeasibility) {
3882                   // take out
3883                   delete branch_;
3884                   branch_ = NULL;
3885                 } else {
3886                   // get preferred way
3887                   int preferredWay;
3888                   object->infeasibility(&usefulInfo, preferredWay);
3889                   CbcBranchingObject * branchObj =
3890                     dynamic_cast <CbcBranchingObject *>(branch_) ;
3891                   assert (branchObj);
3892                   branchObj->way(preferredWay);
3893#ifdef CBCNODE_TIGHTEN_BOUNDS
3894                   bool fixed = branchObj->tighten(solver);
3895                   if (fixed) {
3896                     //printf("Variable now fixed!\n");
3897                     // take out
3898                     delete branch_;
3899                     branch_ = NULL;
3900                   }
3901#endif
3902                 }
3903            }
3904            if (!branch_ && anyAction != -2 && !hitMaxTime) {
3905                finished = false;
3906            }
3907            delete ws;
3908        }
3909    }
3910    // update number of strong iterations etc
3911    model->incrementStrongInfo(numberStrongDone, numberStrongIterations,
3912                               anyAction == -2 ? 0 : numberToFix, anyAction == -2);
3913    if (model->searchStrategy() == -1) {
3914#ifndef COIN_DEVELOP
3915        if (solver->messageHandler()->logLevel() > 1)
3916#endif
3917            printf("%d strong, %d iters, %d inf, %d not finished, %d not trusted\n",
3918                   numberStrongDone, numberStrongIterations, numberStrongInfeasible, numberUnfinished,
3919                   numberNotTrusted);
3920        // decide what to do
3921        int strategy = 1;
3922        if (((numberUnfinished*4 > numberStrongDone &&
3923                numberStrongInfeasible*40 < numberStrongDone) ||
3924                numberStrongInfeasible < 0) && model->numberStrong() < 10 && model->numberBeforeTrust() <= 20 && model->numberObjects() > CoinMax(1000, solver->getNumRows())) {
3925            strategy = 2;
3926#ifdef COIN_DEVELOP
3927            //if (model->logLevel()>1)
3928            printf("going to strategy 2\n");
3929#endif
3930            // Weaken
3931            model->setNumberStrong(2);
3932            model->setNumberBeforeTrust(1);
3933            model->synchronizeNumberBeforeTrust();
3934        }
3935        if (numberNodes)
3936            strategy = 1;  // should only happen after hot start
3937        model->setSearchStrategy(strategy);
3938    } else if (numberStrongDone) {
3939        //printf("%d strongB, %d iters, %d inf, %d not finished, %d not trusted\n",
3940        //   numberStrongDone,numberStrongIterations,numberStrongInfeasible,numberUnfinished,
3941        //   numberNotTrusted);
3942    }
3943    if (model->searchStrategy() == 1 && numberNodes > 500 && numberNodes < -510) {
3944#ifndef COIN_DEVELOP
3945        if (solver->messageHandler()->logLevel() > 1)
3946#endif
3947            printf("after %d nodes - %d strong, %d iters, %d inf, %d not finished, %d not trusted\n",
3948                   numberNodes, numberStrongDone, numberStrongIterations, numberStrongInfeasible, numberUnfinished,
3949                   numberNotTrusted);
3950        // decide what to do
3951        if (numberUnfinished*10 > numberStrongDone + 1 ||
3952                !numberStrongInfeasible) {
3953          COIN_DETAIL_PRINT(printf("going to strategy 2\n"));
3954            // Weaken
3955            model->setNumberStrong(2);
3956            model->setNumberBeforeTrust(1);
3957            model->synchronizeNumberBeforeTrust();
3958            model->setSearchStrategy(2);
3959        }
3960    }
3961    if (numberUnfinished*10 < numberStrongDone &&
3962        model->numberStrongIterations()*20 < model->getIterationCount()&&
3963                                !auxiliaryInfo->solutionAddsCuts()) {
3964        //printf("increasing trust\n");
3965        model->synchronizeNumberBeforeTrust(2);
3966    }
3967
3968    // Set guessed solution value
3969    guessedObjectiveValue_ = objectiveValue_ + estimatedDegradation;
3970    int kColumn=-1;
3971    if (branch_) {
3972      CbcObject * obj = (dynamic_cast<CbcBranchingObject *>(branch_))->object(); 
3973      CbcSimpleInteger * branchObj =
3974        dynamic_cast <CbcSimpleInteger *>(obj) ;
3975      if (branchObj) {
3976        kColumn=branchObj->columnNumber();
3977      }
3978    }
3979#ifdef COIN_HAS_NTY
3980    if (orbitOption&&kColumn>=0) {
3981      CbcSymmetry * symmetryInfo = model->symmetryInfo();
3982      CbcNodeInfo * infoX = lastNode ? lastNode->nodeInfo() : NULL;
3983      bool worthTrying = false;
3984      if (infoX) {
3985        CbcNodeInfo * info = infoX;
3986        for (int i=0;i<NTY_BAD_DEPTH;i++) {
3987          if (!info->parent()) {
3988            worthTrying = true;
3989            break;
3990          }
3991          info = info->parent();
3992          if (info->symmetryWorked()) {
3993            worthTrying = true;
3994            break;
3995          }
3996        }
3997      } else {
3998        worthTrying=true;
3999      }
4000      if (orbitOption==3&&depth_>5)
4001        worthTrying=false;
4002      if (symmetryInfo && worthTrying) {
4003        if ((orbitOption&1)==1) {
4004          symmetryInfo->ChangeBounds(solver->getColLower(),
4005                                     solver->getColUpper(),
4006                                     solver->getNumCols(),false);
4007          symmetryInfo->Compute_Symmetry();
4008          symmetryInfo->fillOrbits();
4009        }
4010        const int * orbits = symmetryInfo->whichOrbit();
4011        if (orbits && orbits[kColumn]>=0) {
4012          int numberUsefulOrbits = symmetryInfo->numberUsefulOrbits();
4013          if (solver->messageHandler()->logLevel() > 1)
4014            printf("Orbital Branching on %d - way %d n %d\n",kColumn,way(),numberUsefulOrbits);
4015          if (numberUsefulOrbits<1000||orbitOption==3) {
4016            delete branch_;
4017            branch_ = new CbcOrbitalBranchingObject(model,kColumn,1,0,NULL);
4018            if (infoX)
4019              infoX->setSymmetryWorked();
4020          }
4021        }
4022      }
4023    }
4024#endif
4025    if (model->logLevel()>1)
4026    printf ("Node %d depth %d unsatisfied %d sum %g obj %g guess %g branching on %d\n",
4027            model->getNodeCount(),depth_,numberUnsatisfied_,
4028            sumInfeasibilities_,objectiveValue_,guessedObjectiveValue_,
4029            kColumn);
4030#ifdef DO_ALL_AT_ROOT
4031    if (strongType) {
4032      char general[200];
4033      if (strongType==1)
4034        sprintf(general,"Strong branching on all %d unsatisfied, %d iterations (depth %d)\n",
4035                saveNumberToDo,numberStrongIterations,depth_);
4036      else
4037        sprintf(general,"Strong branching on all %d unfixed variables (%d unsatisfied), %d iterations (depth %d)\n",
4038                saveNumberToDo+saveSatisfiedVariables,saveNumberToDo,numberStrongIterations,depth_);
4039      model->messageHandler()->message(CBC_FPUMP2,model->messages())
4040        << general << CoinMessageEol ;
4041    }
4042#endif
4043#ifdef DEBUG_SOLUTION
4044    if(onOptimalPath&&anyAction==-2) {
4045      printf("Gone off optimal path in CbcNode\n");
4046      assert(!onOptimalPath||anyAction!=-2);
4047    }
4048#endif
4049    /*
4050      Cleanup, then we're finished
4051    */
4052    if (!model->branchingMethod())
4053        delete decision;
4054
4055    delete choiceObject;
4056    delete [] sort;
4057    delete [] whichObject;
4058#ifdef RANGING
4059    delete [] objectMark;
4060#endif
4061    delete [] saveLower;
4062    delete [] saveUpper;
4063    delete [] upEstimate;
4064    delete [] downEstimate;
4065# ifdef COIN_HAS_CLP
4066    if (osiclp) {
4067        osiclp->setSpecialOptions(saveClpOptions);
4068    }
4069# endif
4070    // restore solution
4071    solver->setColSolution(saveSolution);
4072    model->reserveCurrentSolution(saveSolution);
4073    delete [] saveSolution;
4074    model->setStateOfSearch(saveStateOfSearch);
4075    model->setLogLevel(saveLogLevel);
4076    // delete extra regions
4077    if (usefulInfo.usefulRegion_) {
4078        delete [] usefulInfo.usefulRegion_;
4079        delete [] usefulInfo.indexRegion_;
4080        delete [] usefulInfo.pi_;
4081        usefulInfo.usefulRegion_ = NULL;
4082        usefulInfo.indexRegion_ = NULL;
4083        usefulInfo.pi_ = NULL;
4084    }
4085    useShadow = model->moreSpecialOptions() & 7;
4086    if ((useShadow == 5 && model->getSolutionCount()) || useShadow == 6) {
4087        // zap pseudo shadow prices
4088        model->pseudoShadow(-1);
4089        // and switch off
4090        model->setMoreSpecialOptions(model->moreSpecialOptions()&(~1023));
4091    }
4092    return anyAction;
4093}
4094// 0 is down, 1 is up
4095typedef struct {
4096  double initialValue; // initial value
4097  double upLowerBound; // Lower bound when going up
4098  double downUpperBound; // Upper bound when going down
4099  double movement[2]; // cost  (and initial away from feasible)
4100  double sumModified[2]; // Sum of integer changes
4101  int modified[2]; // Number integers changed
4102  int numIntInfeas[2]; // without odd ones
4103  int numObjInfeas[2]; // just odd ones
4104  bool finished[2]; // true if solver finished
4105  int numIters[2]; // number of iterations in solver (-1 if never solved)
4106  double * integerSolution; // output if thinks integer solution
4107# ifdef COIN_HAS_CLP
4108  ClpDualRowSteepest * steepest;
4109#endif
4110  int columnNumber; // Which column it is
4111} StrongInfo;
4112typedef struct {
4113  double integerTolerance;
4114  double * originalSolution;
4115  CoinWarmStart * ws;
4116  double * newObjective;
4117# ifdef COIN_HAS_CLP
4118  ClpDualRowSteepest * dualRowPivot;
4119  ClpPrimalColumnPivot * primalColumnPivot;
4120# endif
4121  int * back;
4122  int solveType;
4123} StrongStaticInfo;
4124typedef struct {
4125  StrongStaticInfo *staticInfo;
4126  StrongInfo * choice;
4127  OsiSolverInterface * solver;
4128  double * tempSolution;
4129  CoinWarmStart * tempBasis;
4130  int whichChoice;
4131} StrongBundle;
4132/* return 1 if possible solution (for solveType 100 if infeasible)
4133   2 set if down was infeasible
4134   4 set if up was infeasible
4135 */
4136int solveAnalyze(void * info) {
4137  StrongBundle * bundle = reinterpret_cast<StrongBundle *>(info);
4138  StrongInfo * choice = bundle->choice;
4139  StrongStaticInfo * staticInfo = bundle->staticInfo;
4140  OsiSolverInterface * solver = bundle->solver;
4141  int solveType = staticInfo->solveType;
4142  if (solveType==77) {
4143    return 0;
4144  }
4145  const double * saveSolution = staticInfo->originalSolution;
4146  int iColumn = choice->columnNumber;
4147  const int * back = staticInfo->back;
4148  double newObjectiveValue = 1.0e100;
4149  double integerTolerance = staticInfo->integerTolerance;
4150  double bestSolutionValue=COIN_DBL_MAX;
4151  int returnStatus=0;
4152  // status is 0 finished, 1 infeasible and other
4153  int iStatus;
4154  /*
4155    Try the down direction first. (Specify the initial branching alternative as
4156    down with a call to way(-1). Each subsequent call to branch() performs the
4157    specified branch and advances the branch object state to the next branch
4158    alternative.)
4159  */
4160  for (int iWay=0;iWay<2;iWay++) {
4161    if (choice->numIters[iWay]==0) {
4162      int numberColumns=solver->getNumCols();
4163      if (solveType!=100) {
4164        double saveBound;
4165        if (iWay==0) {
4166          saveBound = solver->getColUpper()[iColumn];
4167          solver->setColUpper(iColumn,choice->downUpperBound);
4168        } else {
4169          saveBound = solver->getColLower()[iColumn];
4170          solver->setColLower(iColumn,choice->upLowerBound);
4171        }
4172        if ((solveType&2)==0) {
4173          solver->solveFromHotStart() ;
4174        } else {
4175          // restore basis
4176          solver->setWarmStart(staticInfo->ws);
4177# ifdef COIN_HAS_CLP
4178          if (staticInfo->dualRowPivot) {
4179            OsiClpSolverInterface * osiclp = dynamic_cast<OsiClpSolverInterface *>(solver);
4180            ClpSimplex * simplex = osiclp->getModelPtr();
4181            simplex->setDualRowPivotAlgorithm(*staticInfo->dualRowPivot);
4182            //simplex->dualRowPivot()->saveWeights(simplex,4);
4183            simplex->setWhatsChanged(ALL_SAME_EXCEPT_COLUMN_BOUNDS);
4184            simplex->dual(0,5);
4185          } else {
4186#endif
4187          solver->resolve();
4188# ifdef COIN_HAS_CLP
4189          }
4190#endif
4191        }
4192        if (iWay==0)
4193          solver->setColUpper(iColumn,saveBound);
4194        else
4195          solver->setColLower(iColumn,saveBound);
4196        /*
4197          We now have an estimate of objective degradation that we can use for strong
4198          branching. If we're over the cutoff, the variable is monotone up.
4199          If we actually made it to optimality, check for a solution, and if we have
4200          a good one, call setBestSolution to process it. Note that this may reduce the
4201          cutoff, so we check again to see if we can declare this variable monotone.
4202        */
4203        if (solver->isProvenOptimal()) {
4204          iStatus = 0; // optimal
4205        } else if (solver->isIterationLimitReached()
4206                   && !solver->isDualObjectiveLimitReached()) {
4207          iStatus = 2; // unknown
4208        } else {
4209          iStatus = 1; // infeasible
4210        }
4211        newObjectiveValue = solver->getObjSense() * solver->getObjValue();
4212        choice->numIters[iWay] = solver->getIterationCount();
4213        // Look at interaction
4214        const double * thisSolution = solver->getColSolution();
4215        int numberModified=0;
4216        double sumModified=0.0;
4217        int numberInfeas=0;
4218        for (int i=0;i<numberColumns;i++) {
4219          if (back[i]>=0) {
4220            double value = thisSolution[i];
4221            if (iColumn!=i) {
4222              double difference = fabs(saveSolution[i]-value);
4223              if (difference>integerTolerance) {
4224                numberModified++;
4225                sumModified += difference;
4226              }
4227            }
4228            if (fabs(value-floor(value+0.5))>integerTolerance)
4229              numberInfeas++;;
4230          }
4231        }
4232        choice->numIntInfeas[iWay]=numberInfeas;
4233        choice->sumModified[iWay] = sumModified;
4234        choice->modified[iWay] = numberModified;
4235        if (!iStatus) {
4236          choice->finished[iWay] = true ;
4237          if (!numberInfeas) {
4238            returnStatus=1;
4239            if (!choice->integerSolution) {
4240              bestSolutionValue=newObjectiveValue;
4241              choice->integerSolution=CoinCopyOfArray(thisSolution,numberColumns);;
4242            } else if (bestSolutionValue>newObjectiveValue) {
4243              memcpy(choice->integerSolution,thisSolution,numberColumns*sizeof(double));
4244            }
4245          }
4246        } else if (iStatus == 1) {
4247          newObjectiveValue = 1.0e100 ;
4248        } else {
4249          // Can't say much as we did not finish
4250          choice->finished[iWay] = false ;
4251        }
4252        choice->movement[iWay] = newObjectiveValue ;
4253      } else {
4254# ifdef COIN_HAS_CLP
4255        OsiClpSolverInterface * osiclp = dynamic_cast< OsiClpSolverInterface*> (solver);
4256        ClpSimplex * simplex = osiclp ? osiclp->getModelPtr() : NULL;
4257#endif
4258        // doing continuous and general integer
4259        solver->setColSolution(staticInfo->originalSolution);
4260        solver->setWarmStart(staticInfo->ws);
4261        double saveBound;
4262        double newBound;
4263        if (iWay==0) {
4264          saveBound=solver->getColUpper()[iColumn];
4265          solver->setColUpper(iColumn,choice->downUpperBound);
4266          newBound=choice->downUpperBound;
4267        } else {
4268          saveBound=solver->getColLower()[iColumn];
4269          solver->setColLower(iColumn,choice->upLowerBound);
4270          newBound=choice->upLowerBound;
4271        }
4272# if 0 //def COIN_HAS_CLP
4273        if (simplex) {
4274          // set solution to new bound (if basic will be recomputed)
4275          simplex->primalColumnSolution()[iColumn]=newBound;
4276        }
4277#endif
4278        solver->setHintParam(OsiDoDualInResolve, true, OsiHintDo) ;
4279#define PRINT_ANALYZE  0
4280#if PRINT_ANALYZE>0
4281        osiclp->getModelPtr()->setLogLevel(1);
4282        solver->setHintParam(OsiDoReducePrint, false, OsiHintTry);
4283#endif
4284        solver->resolve();
4285        if (iWay==0) {
4286#if PRINT_ANALYZE>0
4287          printf("column %d down original %g <= %g <= %g upper now %g - result %s\n",
4288                 iColumn,solver->getColLower()[iColumn],
4289                 staticInfo->originalSolution[iColumn],saveBound,
4290                 newBound,solver->isProvenOptimal() ? "ok" : "infeas");
4291#endif
4292          solver->setColUpper(iColumn,saveBound);
4293        } else {
4294#if PRINT_ANALYZE>0
4295          printf("column %d up original %g <= %g <= %g lower now %g - result %s\n",
4296                 iColumn,saveBound,staticInfo->originalSolution[iColumn],
4297                 solver->getColUpper()[iColumn],
4298                 newBound,solver->isProvenOptimal() ? "ok" : "infeas");
4299#endif
4300          solver->setColLower(iColumn,saveBound);
4301        }
4302        choice->numIters[iWay] = solver->getIterationCount();
4303        if (solver->isProvenOptimal()) {
4304          //printf("Way %d - all way %d iterations - column %d\n",
4305          //     iWay,solver->getIterationCount(),iColumn);
4306          // can go all way
4307          choice->movement[iWay] = newBound;
4308        } else {
4309          // zero objective
4310          double offset;
4311          solver->getDblParam(OsiObjOffset,offset);
4312          solver->setDblParam(OsiObjOffset, 0.0);
4313          solver->setObjective(staticInfo->newObjective+numberColumns);
4314          if (iWay==0) {
4315            solver->setObjCoeff(iColumn,1.0);
4316          } else {
4317            solver->setObjCoeff(iColumn,-1.0);
4318          }
4319          solver->setColSolution(staticInfo->originalSolution);
4320          solver->setWarmStart(staticInfo->ws);
4321          solver->setHintParam(OsiDoDualInResolve, false, OsiHintDo) ;
4322          solver->resolve();
4323          //printf("Way %d - first solve %d iterations, second %d - column %d\n",
4324          //     iWay,choice->numIters[iWay],solver->getIterationCount(),iColumn);
4325          choice->movement[iWay] = solver->getColSolution()[iColumn];
4326          choice->numIters[iWay] += solver->getIterationCount();
4327#if PRINT_ANALYZE>0
4328          if (iWay==0) {
4329            printf("column %d down can get to %g - result %s\n",
4330                 iColumn,solver->getColSolution()[iColumn],solver->isProvenOptimal() ? "ok" : "infeas");
4331          } else {
4332            printf("column %d up can get to %g - result %s\n",
4333                 iColumn,solver->getColSolution()[iColumn],solver->isProvenOptimal() ? "ok" : "infeas");
4334          }
4335#endif
4336          // reset objective
4337          solver->setDblParam(OsiObjOffset, offset);
4338          solver->setObjective(staticInfo->newObjective);
4339          if (!solver->isProvenOptimal()) {
4340# ifdef COIN_HAS_CLP
4341            OsiClpSolverInterface * osiclp = dynamic_cast< OsiClpSolverInterface*> (solver);
4342            ClpSimplex * simplex = osiclp->getModelPtr();
4343            double sum = simplex->sumPrimalInfeasibilities();
4344            sum /= static_cast<double>(simplex->numberPrimalInfeasibilities());
4345            if (sum>1.0e-3) {
4346#endif
4347            choice->modified[0]=1;
4348            returnStatus=1;
4349            solver->writeMps("bad","mps");
4350            abort();
4351# ifdef COIN_HAS_CLP
4352            }
4353#endif
4354          }
4355        }
4356        //solver->setObjCoeff(iColumn,0.0);
4357      }
4358    }
4359  }
4360  return returnStatus;
4361}
4362#ifdef THREADS_IN_ANALYZE
4363void * cbc_parallelManager(void * stuff)
4364{
4365  CoinPthreadStuff * driver = reinterpret_cast<CoinPthreadStuff *>(stuff);
4366  int whichThread=driver->whichThread();
4367  CoinThreadInfo * threadInfo = driver->threadInfoPointer(whichThread);
4368  threadInfo->status=-1;
4369  int * which = threadInfo->stuff;
4370  pthread_barrier_wait(driver->barrierPointer());
4371#if 0
4372  int status=-1;
4373  while (status!=100)
4374    status=timedWait(driver,1000,2);
4375  pthread_cond_signal(driver->conditionPointer(1));
4376  pthread_mutex_unlock(driver->mutexPointer(1,whichThread));
4377#endif
4378  // so now mutex_ is locked
4379  int whichLocked=0;
4380  while (true) {
4381    pthread_mutex_t * mutexPointer = driver->mutexPointer(whichLocked,whichThread);
4382    // wait
4383    //printf("Child waiting for %d - status %d %d %d\n",
4384    //     whichLocked,lockedX[0],lockedX[1],lockedX[2]);
4385#ifdef DETAIL_THREAD
4386    printf("thread %d about to lock mutex %d\n",whichThread,whichLocked);
4387#endif
4388    pthread_mutex_lock (mutexPointer);
4389    whichLocked++;
4390    if (whichLocked==3)
4391      whichLocked=0;
4392    int unLock=whichLocked+1;
4393    if (unLock==3)
4394      unLock=0;
4395    //printf("child pointer %p status %d\n",threadInfo,threadInfo->status);
4396    assert(threadInfo->status>=0);
4397    if (threadInfo->status==1000)
4398      pthread_exit(NULL);
4399    int type=threadInfo->status;
4400    int & returnCode=which[0];
4401    int iPass=which[1];
4402    //CoinIndexedVector * array;
4403    //double dummy;
4404    switch(type) {
4405      // dummy
4406    case 0:
4407      break;
4408    case 1:
4409      returnCode=solveAnalyze(threadInfo->extraInfo);
4410      threadInfo->stuff[3]=0;
4411      break;
4412    case 100:
4413      // initialization
4414      break;
4415    }
4416    threadInfo->status= (type!=1) ? -1 : -2;
4417#ifdef DETAIL_THREAD
4418    printf("thread %d about to unlock mutex %d\n",whichThread,unLock);
4419#endif
4420    pthread_mutex_unlock (driver->mutexPointer(unLock,whichThread));
4421  }
4422}
4423#endif
4424int CbcNode::analyze (CbcModel *model, double * results)
4425{
4426#define COIN_DETAIL
4427  int i;
4428  int numberIterationsAllowed = model->numberAnalyzeIterations();
4429  int numberColumns = model->getNumCols();
4430  int numberRows = model->getNumRows();
4431  int numberObjects = model->numberObjects();
4432  int numberIntegers = model->numberIntegers();
4433  int numberLookIntegers=0;
4434  int highestPriority=COIN_INT_MAX;
4435  int * back = new int[numberColumns];
4436  const int * integerVariable = model->integerVariable();
4437  for (i = 0; i < numberIntegers; i++) {
4438    highestPriority = CoinMin(highestPriority,model->modifiableObject(i)->priority());
4439  }
4440  for (i = 0; i < numberColumns; i++) 
4441    back[i] = -1;
4442  for (i = 0; i < numberIntegers; i++) {
4443    int iColumn = integerVariable[i];
4444    back[iColumn] = i;
4445    if (model->modifiableObject(i)->priority()==highestPriority) {
4446      numberLookIntegers++;
4447    } else {
4448      back[iColumn] = i+numberColumns;
4449    }
4450  }
4451  /*
4452    0 - just look
4453    0 (1) bit - use to set priorities
4454    1 (2) bit - look at bounds on all variables and more iterations
4455    2 (4) bit - do threaded (if parallelMode()==1 then not repeatable if any fixed)
4456    3 (8) bit -
4457    4 (16) bit - do even if m*n>1,000,000
4458    5 (32) bit - printing time
4459    6 (64) bit - save mps file
4460  */
4461  int solveType;
4462  char general[200];
4463  if (numberIterationsAllowed>0) {
4464    solveType = 0;
4465  } else {
4466    solveType = - numberIterationsAllowed;
4467    if ((solveType&16)==0) {
4468      double size=numberRows;
4469      size*=numberLookIntegers;
4470      if (size>1000000) {
4471        if ((solveType&32)!=0)
4472          model->messageHandler()->message(CBC_GENERAL, *model->messagesPointer())
4473            << "Skipping analyze as problem too large"
4474            << CoinMessageEol;
4475        return 0;
4476      }
4477    }
4478    sprintf(general,"Analyze options %d",solveType);
4479    model->messageHandler()->message(CBC_GENERAL, *model->messagesPointer())
4480      << general
4481      << CoinMessageEol;
4482    if ((solveType&1)!=0)
4483      model->messageHandler()->message(CBC_GENERAL, *model->messagesPointer())
4484        << "Using to set priorities (probably bad idea)"
4485        << CoinMessageEol;
4486    if ((solveType&2)!=0)
4487      model->messageHandler()->message(CBC_GENERAL, *model->messagesPointer())
4488        << "Use more iterations and look at continuous/general integer variables"
4489        << CoinMessageEol;
4490    if ((solveType&4)!=0)
4491      model->messageHandler()->message(CBC_GENERAL, *model->messagesPointer())
4492        << "Use threads"
4493        << CoinMessageEol;
4494    if ((solveType&32)!=0)
4495      model->messageHandler()->message(CBC_GENERAL, *model->messagesPointer())
4496        << "32 switches on more printing, (16 bit allows large problems)"
4497        << CoinMessageEol;
4498  }
4499  OsiSolverInterface * solver = model->solver();
4500  objectiveValue_ = solver->getObjSense() * solver->getObjValue();
4501  const double * lower = solver->getColLower();
4502  const double * upper = solver->getColUpper();
4503  const double * dj = solver->getReducedCost();
4504  // What results is
4505  double * newLower = results;
4506  double * objLower = newLower + numberIntegers;
4507  double * newUpper = objLower + numberIntegers;
4508  double * objUpper = newUpper + numberIntegers;
4509  double * interAction = objUpper + numberIntegers;
4510  for (i = 0; i < numberIntegers; i++) {
4511    int iColumn = integerVariable[i];
4512    newLower[i] = lower[iColumn];
4513    objLower[i] = -COIN_DBL_MAX;
4514    newUpper[i] = upper[iColumn];
4515    objUpper[i] = -COIN_DBL_MAX;
4516    interAction[i] = 0.0;
4517  }
4518  double * objMovement=new double[2*numberIntegers];
4519  memset(objMovement,0,2*numberIntegers*sizeof(double));
4520  double * saveUpper = new double[numberColumns];
4521  double * saveLower = new double[numberColumns];
4522  // Save solution in case heuristics need good solution later
4523 
4524  double * saveSolution = new double[numberColumns];
4525  memcpy(saveSolution, solver->getColSolution(), numberColumns*sizeof(double));
4526  model->reserveCurrentSolution(saveSolution);
4527  for (i = 0; i < numberColumns; i++) {
4528    saveLower[i] = lower[i];
4529    saveUpper[i] = upper[i];
4530  }
4531  // Get arrays to sort
4532  double * sort = new double[numberObjects];
4533  int * whichObject = new int[numberObjects];
4534  int numberToFix = 0;
4535  int numberToDo = 0;
4536  double integerTolerance =
4537    model->getDblParam(CbcModel::CbcIntegerTolerance);
4538  // point to useful information
4539  OsiBranchingInformation usefulInfo = model->usefulInformation();
4540  // and modify
4541  usefulInfo.depth_ = depth_;
4542 
4543  // compute current state
4544  int numberObjectInfeasibilities; // just odd ones
4545  int numberIntegerInfeasibilities;
4546  model->feasibleSolution(
4547                          numberIntegerInfeasibilities,
4548                          numberObjectInfeasibilities);
4549  if (solveType) {
4550    if ((solveType&2)==0)
4551      numberIterationsAllowed=200*numberIntegerInfeasibilities;
4552    else
4553      numberIterationsAllowed=COIN_INT_MAX;
4554  }
4555  int saveAllowed=numberIterationsAllowed;
4556# ifdef COIN_HAS_CLP
4557  OsiClpSolverInterface * osiclp = dynamic_cast< OsiClpSolverInterface*> (solver);
4558  int saveClpOptions = 0;
4559  bool fastIterations = (model->specialOptions() & 8) != 0;
4560  if (osiclp) {
4561    saveClpOptions = osiclp->specialOptions();
4562    // for faster hot start
4563    if (fastIterations) 
4564      osiclp->setSpecialOptions(saveClpOptions | 8192);
4565    else
4566      osiclp->setSpecialOptions(saveClpOptions | 2048); // switch off crunch
4567  }
4568# else
4569  bool fastIterations = false ;
4570# endif
4571  /*
4572    Scan for branching objects that indicate infeasibility.
4573   
4574    The algorithm is to fill the array with a set of good candidates (by
4575    infeasibility).
4576   
4577  */
4578  numberToDo = 0;
4579  for (i = 0; i < numberObjects; i++) {
4580    OsiObject * object = model->modifiableObject(i);
4581    CbcSimpleIntegerDynamicPseudoCost * dynamicObject =
4582      dynamic_cast <CbcSimpleIntegerDynamicPseudoCost *>(object) ;
4583    if (!dynamicObject)
4584      continue;
4585    if (dynamicObject->priority()!=highestPriority)
4586      continue;
4587    double infeasibility = object->checkInfeasibility(&usefulInfo);
4588    int iColumn = dynamicObject->columnNumber();
4589    if (saveUpper[iColumn] == saveLower[iColumn])
4590      continue;
4591    if (infeasibility)
4592      sort[numberToDo] = -1.0e10 - infeasibility;
4593    else
4594      sort[numberToDo] = -fabs(dj[iColumn]);
4595    whichObject[numberToDo++] = i;
4596  }
4597  // Save basis
4598  CoinWarmStart * ws = solver->getWarmStart();
4599  int saveLimit;
4600  solver->getIntParam(OsiMaxNumIterationHotStart, saveLimit);
4601  int targetIterations = CoinMax(500, numberIterationsAllowed / numberObjects);
4602  if (saveLimit < targetIterations)
4603    solver->setIntParam(OsiMaxNumIterationHotStart, targetIterations);
4604  if ((solveType&2)==0) {
4605    // Mark hot start
4606    solver->markHotStart();
4607  }
4608  solver->setHintParam(OsiDoDualInResolve, true, OsiHintDo) ;
4609  // Sort
4610  CoinSort_2(sort, sort + numberToDo, whichObject);
4611  double * currentSolution = model->currentSolution();
4612  double objMin = 1.0e50;
4613  double objMax = -1.0e50;
4614  bool needResolve = false;
4615  int maxChoices=1;
4616  int currentChoice=0;
4617  int numberThreads=0;
4618  bool doAtEnd=false;
4619  if (model->parallelMode() && (solveType&4)!=0) {
4620    numberThreads=model->getNumberThreads();
4621    sprintf(general,"Using %d threads in analysis\n",numberThreads);
4622    model->messageHandler()->message(CBC_GENERAL, *model->messagesPointer())
4623      << general
4624      << CoinMessageEol;
4625    if (model->parallelMode()==1) {
4626      maxChoices=numberThreads;
4627    } else {
4628      maxChoices = numberToDo;
4629      if ((solveType&2)!=0) 
4630        maxChoices = numberColumns;
4631      doAtEnd=true;
4632    }
4633  }
4634  StrongInfo * choices = new StrongInfo[maxChoices];
4635  StrongStaticInfo staticInfo;
4636  int numberBundles = CoinMax(1,numberThreads);
4637  StrongBundle * bundles = new StrongBundle[numberBundles];
4638  /*
4639    0 - available - no need to look at results
4640    1 - not available
4641    2 - available - need to look at results
4642  */
4643#ifndef NUMBER_THREADS
4644#define NUMBER_THREADS 4
4645#endif
4646  int status[NUMBER_THREADS];
4647  memset(status,0,sizeof(status));
4648  memset(&staticInfo,0,sizeof(staticInfo));
4649  staticInfo.solveType = solveType;
4650  staticInfo.originalSolution=saveSolution;
4651  staticInfo.back=back;
4652  staticInfo.ws=ws;
4653  staticInfo.integerTolerance=integerTolerance;
4654  double time1 = model->getCurrentSeconds();
4655#define DO_STEEPEST_SERIAL 1
4656# ifdef COIN_HAS_CLP
4657  if (osiclp&&(solveType&2)!=0&&(!numberThreads||DO_STEEPEST_SERIAL)) {
4658    ClpSimplex * simplex = osiclp->getModelPtr();
4659    simplex->setLogLevel(0);
4660    simplex->dual(0,1);
4661    ClpDualRowPivot * dualRowPivot=simplex->dualRowPivot();
4662    ClpDualRowSteepest * steep = dynamic_cast<ClpDualRowSteepest *>(dualRowPivot);
4663    if (steep) {
4664      staticInfo.dualRowPivot=new ClpDualRowSteepest (*steep);
4665      staticInfo.dualRowPivot->setMode(1); // full steepest edge
4666      simplex->spareIntArray_[0]=0;
4667      simplex->spareIntArray_[1]=numberRows;
4668      staticInfo.dualRowPivot->saveWeights(simplex,7);
4669    }
4670  }
4671#endif
4672  for (int i=0;i<numberBundles;i++) {
4673    memset(bundles+i,0,sizeof(StrongBundle));
4674    bundles[i].staticInfo=&staticInfo;
4675  }
4676#if defined (THREADS_IN_ANALYZE) && defined (COIN_HAS_CLP)
4677#define USE_STRONG_THREADS
4678  CoinPthreadStuff threadInfo(numberThreads,cbc_parallelManager);
4679  int threadNeedsRefreshing[NUMBER_THREADS];
4680  for (int i=0;i<numberThreads;i++) {
4681    threadInfo.threadInfo_[i].extraInfo2 = solver->clone();
4682    threadNeedsRefreshing[i]=0;
4683  }
4684# ifdef COIN_HAS_CLP
4685  int numberSteepThreads=0;
4686  int step=numberThreads ? (numberRows+numberThreads-1)/numberThreads : 0;
4687  int first=0;
4688  for (int i=0;i<numberThreads;i++) {
4689    if (osiclp&&(solveType&2)!=0&&!DO_STEEPEST_SERIAL) {
4690      OsiSolverInterface * solver=
4691        reinterpret_cast<OsiSolverInterface *>(threadInfo.threadInfo_[i].extraInfo2);
4692      OsiClpSolverInterface * osiclp = dynamic_cast< OsiClpSolverInterface*> (solver);
4693      ClpSimplex * simplex = osiclp->getModelPtr();
4694      simplex->setLogLevel(0);
4695      simplex->dual(0,1);
4696      ClpDualRowPivot * dualRowPivot=simplex->dualRowPivot();
4697      ClpDualRowSteepest * steep = dynamic_cast<ClpDualRowSteepest *>(dualRowPivot);
4698      if (steep) {
4699        numberSteepThreads=numberThreads;
4700        ClpDualRowSteepest * dualRowPivot=new ClpDualRowSteepest (*steep);
4701        dualRowPivot->setMode(1); // full steepest edge
4702        simplex->spareIntArray_[0]=0;
4703        simplex->spareIntArray_[1]=numberRows;
4704        simplex->spareIntArray_[0]=first;
4705        simplex->spareIntArray_[1]=CoinMin(first+step,numberRows);
4706        first += step;
4707        if (i==0)
4708          staticInfo.dualRowPivot=dualRowPivot;
4709        choices[i].steepest=dualRowPivot;
4710        dualRowPivot->saveWeights(simplex,7);
4711      }
4712    }
4713  }
4714  if (numberSteepThreads&&false) {
4715    int numberDone=0;
4716    int iDo=0;
4717    staticInfo.solveType = 200;
4718    while (numberDone<numberSteepThreads) {
4719      int iThread;
4720      threadInfo.waitParallelTask(1,iThread,iDo<numberToDo);
4721      int threadStatus = 1+threadInfo.threadInfo_[iThread].status;
4722      iThread=iThread;
4723      if (threadStatus==0&&iDo<numberSteepThreads) {
4724        StrongInfo & choice = choices[iThread];
4725        StrongBundle & bundle = bundles[iThread];
4726        bundle.whichChoice=iThread;
4727        memset(&choice,0,sizeof(StrongInfo));
4728        iDo++; //started this one
4729        bundle.choice=&choice;
4730        bundle.solver = solver;
4731        bundle.solver=reinterpret_cast<OsiSolverInterface *>(threadInfo.threadInfo_[iThread].extraInfo2);
4732        threadStatus=0;
4733#ifdef DETAIL_THREAD
4734        printf("Starting steep task on thread %d\n",
4735               choice.iThread);
4736#endif
4737        threadInfo.startParallelTask(1,iThread,&bundle);
4738      }
4739      if (!threadStatus) {
4740#ifdef _MSC_VER
4741        Sleep(1);
4742#else
4743        usleep(1000);
4744#endif
4745        continue;
4746      }
4747      if (threadStatus) {
4748        numberDone++;
4749        // say available
4750        threadInfo.sayIdle(iThread);
4751      }
4752      staticInfo.solveType = solveType;
4753    }
4754    OsiSolverInterface * solver0=
4755      reinterpret_cast<OsiSolverInterface *>(threadInfo.threadInfo_[0].extraInfo2);
4756    CoinIndexedVector * savedWeights0 = staticInfo.dualRowPivot->savedWeights();
4757    int * index0 = savedWeights0->getIndices();
4758    double * weight0 = savedWeights0->denseVector();
4759    int step=(numberRows+numberSteepThreads-1)/numberSteepThreads;
4760    int first=step;
4761    //memset(weight0+first,0,(numberRows-first)*sizeof(double));
4762    for (int i=1;i<numberSteepThreads;i++) {
4763      int n=CoinMin(step,numberRows-first);
4764      CoinIndexedVector * savedWeights = choices[i].steepest->savedWeights();
4765      int * index = savedWeights->getIndices();
4766      double * weight = savedWeights->denseVector();
4767      memcpy(index0+first,index+first,n*sizeof(int));
4768      memcpy(weight0+first,weight+first,n*sizeof(double));
4769      first += step;
4770      delete choices[i].steepest;
4771      choices[i].steepest=NULL;
4772    }
4773    //for (int j=0;j<numberRows;j++)
4774    //weight0[j]=1.0;
4775  }
4776#endif
4777#endif
4778  double bestSolutionValue = model->getMinimizationObjValue();
4779  double * bestSolution = NULL;
4780  double cutoff;
4781  solver->getDblParam(OsiDualObjectiveLimit,cutoff);
4782  double maxMovement = 2.0*(cutoff-objectiveValue_)+1.0e-6;
4783  /*
4784    Now calculate the cost forcing the variable up and down.
4785  */
4786  int iDo=0;
4787  int iDone=-1;
4788  int numberDone=0;
4789  int iThread=0;
4790  int threadStatus=0;
4791  int whenPrint = (numberToDo+9)/10;
4792  while (numberDone<numberToDo) {
4793    if ((solveType&32)!=0&&numberDone==whenPrint) {
4794      whenPrint += (numberToDo+9)/10;
4795      sprintf(general,"%d variables looked at - %d changed by initial strong branching (%.2f seconds - %d iterations)",numberDone,numberToFix,model->getCurrentSeconds()-time1,saveAllowed-numberIterationsAllowed);
4796      model->messageHandler()->message(CBC_GENERAL, *model->messagesPointer())
4797        << general
4798        << CoinMessageEol;
4799    }
4800#ifdef USE_STRONG_THREADS
4801    if (numberThreads) {
4802      threadInfo.waitParallelTask(1,iThread,iDo<numberToDo);
4803      threadStatus = 1+threadInfo.threadInfo_[iThread].status;
4804      if (!doAtEnd)
4805        currentChoice=iThread;
4806      if (threadNeedsRefreshing[iThread]) { 
4807        OsiSolverInterface * solver=
4808          reinterpret_cast<OsiSolverInterface *>(threadInfo.threadInfo_[iThread].extraInfo2);
4809        if ((threadNeedsRefreshing[iThread]&1)!=0) 
4810          solver->setColLower(saveLower);
4811        if ((threadNeedsRefreshing[iThread]&2)!=0) 
4812          solver->setColUpper(saveUpper);
4813        threadNeedsRefreshing[iThread]=0;
4814      }
4815    }
4816#endif
4817    if (threadStatus==0&&iDo<numberToDo) {
4818      StrongInfo & choice = choices[currentChoice];
4819      StrongBundle & bundle = bundles[iThread];
4820      bundle.whichChoice=currentChoice;
4821      memset(&choice,0,sizeof(StrongInfo));
4822      int iObject = whichObject[iDo];
4823      iDo++; //started this one
4824      OsiObject * object = model->modifiableObject(iObject);
4825      CbcSimpleIntegerDynamicPseudoCost * dynamicObject =
4826        dynamic_cast <CbcSimpleIntegerDynamicPseudoCost *>(object) ;
4827      int iColumn = dynamicObject->columnNumber();
4828      double value = currentSolution[iColumn];
4829      double nearest = floor(value + 0.5);
4830      double lowerValue = floor(value);
4831      bool satisfied = false;
4832      if (fabs(value - nearest) <= integerTolerance ||
4833          value < saveLower[iColumn] || value > saveUpper[iColumn]) {
4834        satisfied = true;
4835        if (nearest < saveUpper[iColumn]) {
4836          lowerValue = nearest;
4837        } else {
4838          lowerValue = nearest - 1;
4839        }
4840      }
4841      double upperValue = lowerValue + 1.0;
4842      // Save which object it was
4843      choice.columnNumber = iColumn;
4844      choice.initialValue=value;
4845      choice.upLowerBound=upperValue;
4846      choice.downUpperBound=lowerValue;
4847      choice.numIntInfeas[1] = numberUnsatisfied_;
4848      choice.numIntInfeas[0] = numberUnsatisfied_;
4849      choice.movement[0] = 0.0;
4850      choice.movement[1] = 0.0;
4851      choice.numIters[0] = 0;
4852      choice.numIters[1] = 0;
4853      if (fabs(value - lowerValue) <= integerTolerance) 
4854        choice.numIters[0]=-1; // mark as not done
4855      if (fabs(value - upperValue) <= integerTolerance) 
4856        choice.numIters[1]=-1; // mark as not done
4857      bundle.choice=&choice;
4858      bundle.solver = solver;
4859#ifdef USE_STRONG_THREADS
4860      if (numberThreads) {
4861        bundle.solver=reinterpret_cast<OsiSolverInterface *>(threadInfo.threadInfo_[iThread].extraInfo2);
4862        threadStatus=0;
4863#ifdef DETAIL_THREAD
4864        printf("Starting task for column %d on thread %d\n",
4865               choice.columnNumber,iThread);
4866#endif
4867        threadInfo.startParallelTask(1,iThread,&bundle);
4868      } else {
4869#endif
4870        threadStatus=2;
4871        solveAnalyze(&bundle);
4872#ifdef USE_STRONG_THREADS
4873      }
4874#endif
4875    }
4876    if (!threadStatus) {
4877#ifdef _MSC_VER
4878        Sleep(1);
4879#else
4880        usleep(1000);
4881#endif
4882      continue;
4883    }
4884    if (threadStatus) {
4885      int whichChoice = bundles[iThread].whichChoice; 
4886      StrongInfo & choice = choices[whichChoice];
4887      int iColumn=choice.columnNumber;
4888      if (choice.integerSolution) {
4889        double * foundSolution = choice.integerSolution;
4890        solver->setColSolution(foundSolution);
4891        // See if integer solution
4892        int numberInfeas=0;
4893        int numberOddInfeas=0;
4894        if (model->feasibleSolution(numberInfeas,numberOddInfeas)
4895            && model->problemFeasibility()->feasible(model, -1) >= 0) {
4896          double newObjectiveValue;
4897          solver->getDblParam(OsiObjOffset,newObjectiveValue);
4898          newObjectiveValue=-newObjectiveValue;
4899          const double * cost = solver->getObjCoefficients();
4900          for ( int i = 0 ; i < numberColumns ; i++ )
4901            newObjectiveValue += cost[i] * foundSolution[i];
4902          if (newObjectiveValue<bestSolutionValue) {
4903            if (doAtEnd) {
4904              if (!bestSolution)
4905                bestSolution = CoinCopyOfArray(foundSolution,numberColumns);
4906              else
4907                memcpy(bestSolution,foundSolution,numberColumns*sizeof(double));
4908              bestSolutionValue=newObjectiveValue;
4909            } else {
4910              model->setBestSolution(CBC_STRONGSOL,
4911                                     newObjectiveValue,
4912                                     foundSolution) ;
4913              model->setLastHeuristic(NULL);
4914              model->incrementUsed(solver->getColSolution());
4915              bestSolutionValue = model->getMinimizationObjValue();
4916            }
4917          }
4918        }
4919        delete [] foundSolution;
4920      }
4921      for (int iWay=0;iWay<2;iWay++) {
4922        numberIterationsAllowed -= choice.numIters[iWay];
4923        choice.movement[iWay] -= objectiveValue_;
4924      }
4925      // If objective goes above certain amount we can set bound
4926      int jInt = back[iColumn];
4927      OsiObject * object = model->modifiableObject(jInt);
4928      CbcSimpleIntegerDynamicPseudoCost * dynamicObject =
4929        dynamic_cast <CbcSimpleIntegerDynamicPseudoCost *>(object) ;
4930      if (dynamicObject) {
4931        if (choice.numIters[0]>=0) {
4932          dynamicObject->addToSumDownCost(CoinMin(choice.movement[0],maxMovement));
4933          dynamicObject->addToSumDownChange(choice.initialValue-choice.downUpperBound);
4934        }
4935        if (choice.numIters[1]>=0) {
4936          dynamicObject->addToSumUpCost(CoinMin(choice.movement[1],maxMovement));
4937          dynamicObject->addToSumUpChange(choice.upLowerBound-choice.initialValue);
4938        }
4939      }
4940      newLower[jInt] = choice.upLowerBound;
4941      if (choice.finished[0])
4942        objLower[jInt] = choice.movement[0] + objectiveValue_;
4943      else
4944        objLower[jInt] = objectiveValue_;
4945      newUpper[jInt] = choice.downUpperBound;
4946      if (choice.finished[1])
4947        objUpper[jInt] = choice.movement[1] + objectiveValue_;
4948      else
4949        objUpper[jInt] = objectiveValue_;
4950      objMin = CoinMin(CoinMin(objLower[jInt], objUpper[jInt]), objMin);
4951      objMovement[2*jInt]=choice.movement[0];
4952      objMovement[2*jInt+1]=choice.movement[1];
4953      double sumModified = choice.modified[0] + choice.modified[1] +
4954        1.0e-15*(choice.sumModified[0]+choice.sumModified[1]);
4955      if (choice.numIters[0]>=0&&choice.numIters[1]>=0)
4956        sumModified *= 0.6;
4957      interAction[jInt] = sumModified;
4958      /*
4959        End of evaluation for this candidate variable. Possibilities are:
4960        * Both sides below cutoff; this variable is a candidate for branching.
4961        * Both sides infeasible or above the objective cutoff: no further action
4962        here. Break from the evaluation loop and assume the node will be purged
4963        by the caller.
4964        * One side below cutoff: Install the branch (i.e., fix the variable). Break
4965        from the evaluation loop and assume the node will be reoptimised by the
4966        caller.
4967      */
4968      threadStatus=0;
4969      currentChoice++;
4970      numberDone++;
4971#ifdef USE_STRONG_THREADS
4972      // say available
4973      if (numberThreads) {
4974        threadInfo.sayIdle(iThread);
4975      }
4976#endif
4977      if (doAtEnd)
4978        continue;
4979      if (choice.movement[1] < 1.0e100) {
4980        if (choice.movement[0] < 1.0e100) {
4981          objMax = CoinMax(CoinMax(objLower[jInt], objUpper[jInt]), objMax);
4982          // In case solution coming in was odd
4983          choice.movement[1] = CoinMax(0.0, choice.movement[1]);
4984          choice.movement[0] = CoinMax(0.0, choice.movement[0]);
4985          // feasible -
4986          model->messageHandler()->message(CBC_STRONG, *model->messagesPointer())
4987            << iColumn << iColumn
4988            << choice.movement[0] << choice.numIntInfeas[0]
4989            << choice.movement[1] << choice.numIntInfeas[1]
4990            << choice.initialValue
4991            << CoinMessageEol;
4992        } else {
4993          // up feasible, down infeasible
4994          needResolve = true;
4995          numberToFix++;
4996          saveLower[iColumn] = choice.upLowerBound;
4997          solver->setColLower(iColumn, choice.upLowerBound);
4998#ifdef USE_STRONG_THREADS
4999          for (int i=0;i<numberThreads;i++) {
5000            threadNeedsRefreshing[i] |= 1;
5001          }
5002#endif
5003        }
5004      } else {
5005        if (choice.movement[0] < 1.0e100) {
5006          // down feasible, up infeasible
5007          needResolve = true;
5008          numberToFix++;
5009          saveUpper[iColumn] = choice.downUpperBound;
5010          solver->setColUpper(iColumn, choice.downUpperBound);
5011#ifdef USE_STRONG_THREADS
5012          for (int i=0;i<numberThreads;i++) {
5013            threadNeedsRefreshing[i] |= 2;
5014          }
5015#endif
5016        } else {
5017          // neither side feasible
5018          COIN_DETAIL_PRINT(printf("Both infeasible for choice %d sequence %d\n", i,
5019                                   model->object(choice.objectNumber)->columnNumber()));
5020          //solver->writeMps("bad");
5021          numberToFix = -1;
5022          break;
5023        }
5024      }
5025      if (numberIterationsAllowed <= 0)
5026        break;
5027      if (currentChoice==maxChoices)
5028        currentChoice=0;
5029    }
5030    //printf("obj %d, col %d, down %g up %g value %g\n",iObject,iColumn,
5031    //     choice.downMovement,choice.upMovement,value);
5032  }
5033  // Do at end if deterministic
5034  if (doAtEnd) {
5035    if (bestSolution) {
5036      model->setBestSolution(CBC_STRONGSOL,
5037                             bestSolutionValue,
5038                             bestSolution) ;
5039      model->setLastHeuristic(NULL);
5040      model->incrementUsed(solver->getColSolution());
5041      delete [] bestSolution;
5042    }
5043    for (int iDo = 0; iDo < numberLookIntegers; iDo++) {
5044      StrongInfo & choice = choices[iDo];
5045      int iColumn = choice.columnNumber;
5046      int iObject = iColumn;
5047      int jInt = back[iColumn];
5048      double value = choice.initialValue;
5049      double lowerValue = choice.downUpperBound;
5050      double upperValue = choice.upLowerBound;
5051      if (choice.movement[1] < 1.0e100) {
5052        if (choice.movement[0] < 1.0e100) {
5053          objMax = CoinMax(CoinMax(objLower[jInt], objUpper[jInt]), objMax);
5054          // In case solution coming in was odd
5055          choice.movement[1] = CoinMax(0.0, choice.movement[1]);
5056          choice.movement[0] = CoinMax(0.0, choice.movement[0]);
5057          // feasible -
5058          model->messageHandler()->message(CBC_STRONG, *model->messagesPointer())
5059            << iObject << iColumn
5060            << choice.movement[0] << choice.numIntInfeas[0]
5061            << choice.movement[1] << choice.numIntInfeas[1]
5062            << value
5063            << CoinMessageEol;
5064        } else {
5065          // up feasible, down infeasible
5066          numberToFix++;
5067          saveLower[iColumn] = upperValue;
5068          solver->setColLower(iColumn, upperValue);
5069        }
5070      } else {
5071        if (choice.movement[0] < 1.0e100) {
5072          // down feasible, up infeasible
5073          needResolve = true;
5074          numberToFix++;
5075          saveUpper[iColumn] = lowerValue;
5076          solver->setColUpper(iColumn, lowerValue);
5077        } else {
5078          // neither side feasible
5079          COIN_DETAIL_PRINT(printf("Both infeasible for choice %d sequence %d\n", i,
5080                                   model->object(choice.objectNumber)->columnNumber()));
5081          //solver->writeMps("bad");
5082          numberToFix = -1;
5083          break;
5084        }
5085      }
5086    }
5087  }
5088  if (false) {
5089    const double * lower = solver->getColLower();
5090    const double * upper = solver->getColUpper();
5091    for (int i=0;i<numberColumns;i++) {
5092      if (lower[i]!=saveLower[i]||upper[i]!=saveUpper[i]) {
5093        printf("%d changed- saved %g,%g now %g,%g\n",i,saveLower[i],saveUpper[i],
5094               lower[i],upper[i]);
5095      }
5096    }
5097  }
5098  COIN_DETAIL_PRINT(printf("Best possible solution %g, can fix more if solution of %g found - looked at %d variables in %d iterations\n",
5099                           objMin, objMax, iDo, model->numberAnalyzeIterations() - numberIterationsAllowed));
5100  model->setNumberAnalyzeIterations(numberIterationsAllowed);
5101  if (numberToFix>0) {
5102    sprintf(general,"%d variable bounds modified by initial strong branching (%.2f seconds - %d iterations)",numberToFix,model->getCurrentSeconds()-time1,saveAllowed-numberIterationsAllowed);
5103  } else if (numberToFix<0) {
5104    sprintf(general,"initial strong branching found to be infeasible (%.2f seconds - %d iterations)",model->getCurrentSeconds()-time1,saveAllowed-numberIterationsAllowed);
5105  } else if ((solveType&32)!=0) {
5106    sprintf(general,"No variables fixed by initial strong branching (%.2f seconds - %d iterations)",model->getCurrentSeconds()-time1,saveAllowed-numberIterationsAllowed);
5107  } else {
5108    general[0]='\0';
5109  }
5110  if (general[0]!='\0')
5111    model->messageHandler()->message(CBC_GENERAL, *model->messagesPointer())
5112      << general
5113      << CoinMessageEol;
5114  double smallestEffect=COIN_DBL_MAX;
5115  double largestEffect=0.0;
5116  for (i = 0; i < numberIntegers; i++) {
5117    int iColumn=integerVariable[i];
5118    if (back[iColumn]>=numberColumns)
5119      continue;
5120    smallestEffect = CoinMin(smallestEffect,interAction[i]);
5121    largestEffect = CoinMax(largestEffect,interAction[i]);
5122  }
5123  double groupValue[11];
5124  int groupCounts[11]={0,0,0,0,0,0,0,0,0,0,0};
5125  groupValue[10]=largestEffect;
5126  for (int i=0;i<10;i++)
5127    groupValue[i]=smallestEffect+i*0.1*(largestEffect-smallestEffect);
5128  sprintf(general,"Looked at %d integer variables - smallest interaction %g",
5129          numberLookIntegers,smallestEffect);
5130  model->messageHandler()->message((solveType&32)==0 ? CBC_FPUMP2 : CBC_FPUMP1, 
5131                                    *model->messagesPointer())
5132    << general
5133    << CoinMessageEol;
5134  for (int i = 0; i < numberIntegers; i++) {
5135    int iColumn=integerVariable[i];
5136    if (back[iColumn]>=numberColumns)
5137      continue;
5138    double value = interAction[i];
5139    int j;
5140    for (j=0;j<11;j++) {
5141      if (value<=groupValue[j]||j==10)
5142        break;
5143    }
5144    groupCounts[j]++;
5145  }
5146  general[0]='\0';
5147  for (int i=0;i<11;i++)
5148    sprintf(general+strlen(general),"%d <= %g ",groupCounts[i],groupValue[i]);
5149  model->messageHandler()->message((solveType&32)==0 ? CBC_FPUMP2 : CBC_FPUMP1, 
5150                                    *model->messagesPointer())
5151    << general
5152    << CoinMessageEol;
5153  smallestEffect=COIN_DBL_MAX;
5154  largestEffect=0.0;
5155  int numberChanged=0;
5156  int numberZeroMoved=0;
5157  for (i = 0; i < numberIntegers; i++) {
5158    int iColumn=integerVariable[i];
5159    if (back[iColumn]>=numberColumns)
5160      continue;
5161    for (int iWay=0;iWay<2