source: trunk/Cbc/src/CbcNode.cpp @ 2313

Last change on this file since 2313 was 2313, checked in by forrest, 23 months ago

memory leaks

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