source: trunk/Cbc/src/CbcGenBaB.cpp @ 1464

Last change on this file since 1464 was 1464, checked in by stefan, 9 years ago

merge split branch into trunk; fix some examples

File size: 30.2 KB
Line 
1/*
2  Copyright (C) 2007, Lou Hafer, International Business Machines Corporation
3  and others.  All Rights Reserved.
4
5  This file is part of cbc-generic.
6*/
7
8#include <iostream>
9
10#include "CoinTime.hpp"
11
12#include "OsiSolverInterface.hpp"
13#include "OsiChooseVariable.hpp"
14
15#include "CglPreProcess.hpp"
16
17#include "CbcModel.hpp"
18#include "CbcCutGenerator.hpp"
19#include "CbcBranchActual.hpp"
20#include "CbcStrategy.hpp"
21
22#include "CbcGenCtlBlk.hpp"
23#include "CbcGenParam.hpp"
24#include "CbcGenCbcParam.hpp"
25
26#define CBC_TRACK_SOLVERS 1
27// #define COIN_CBC_VERBOSITY 5
28
29/*
30  The support functions for the main branch-and-cut action routine.
31*/
32
33namespace {
34
35char svnid[] = "$Id: CbcGenBaB.cpp 1173 2009-06-04 09:44:10Z forrest $" ;
36
37/*
38  A hack to fix variables based on reduced cost prior to branch-and-cut. Note
39  that we're *not* looking at the integrality gap here. Given the reduced costs
40  of the root relaxation, we're simply placing a bet that variables with really
41  unfavourable reduced costs that are at their most favourable bound in the
42  root relaxation will never move from that bound.
43
44  For the standard OsiSolverInterface, this requires a bit of effort as the
45  solution and bounds arrays are const and the functions to change them have
46  incompatible interfaces.
47*/
48
49void reducedCostHack (OsiSolverInterface *osi, double threshold)
50
51{
52    int numCols = osi->getNumCols() ;
53    int i ;
54    const double *lower = osi->getColLower() ;
55    const double *upper = osi->getColUpper() ;
56    const double *solution = osi->getColSolution() ;
57    const double *dj = osi->getReducedCost() ;
58    /*
59      First task: scan the columns looking for variables that are at their
60      favourable bound and have reduced cost that exceeds the threshold. Remember
61      the column index and the value.
62    */
63    double *chgBnds = new double [numCols] ;
64    int *chgCols = new int [numCols] ;
65
66    int numFixed = 0 ;
67    for (i = 0 ; i < numCols ; i++) {
68        if (osi->isInteger(i)) {
69            double value = solution[i] ;
70            if (value < lower[i] + 1.0e-5 && dj[i] > threshold) {
71                chgCols[numFixed] = i ;
72                chgBnds[numFixed] = lower[i] ;
73                numFixed++ ;
74            } else if (value > upper[i] - 1.0e-5 && dj[i] < -threshold) {
75                chgCols[numFixed] = i ;
76                chgBnds[numFixed] = upper[i] ;
77                numFixed++ ;
78            }
79        }
80    }
81    /*
82      Second task: For variables that we want to fix, we need to:
83        * Prepare an array with the new lower and upper bounds for variables that
84          will be fixed. setColSetBounds requires an array with column indices and
85          an array with new values for both bounds.
86        * Set the correct value in a copy of the current solution. setColSolution
87          requires a complete solution.
88    */
89    if (numFixed > 0) {
90        double *newSoln = CoinCopyOfArray(solution, numCols) ;
91        double *newBnds = new double [2*numFixed] ;
92        double *bndPtr = &newBnds[0] ;
93        for (i = 0 ; i < numFixed ; i++) {
94            int j = chgCols[i] ;
95            double val = chgBnds[i] ;
96            *bndPtr++ = val ;
97            *bndPtr++ = val ;
98            newSoln[j] = val ;
99        }
100        osi->setColSetBounds(&chgCols[0], &chgCols[numFixed], &newBnds[0]) ;
101        osi->setColSolution(&newSoln[0]) ;
102
103        std::cout
104            << "Reduced cost fixing prior to B&C: " << numFixed
105            << " columns fixed." << std::endl ;
106
107        delete[] newSoln ;
108        delete[] newBnds ;
109    }
110
111    delete[] chgBnds ;
112    delete[] chgCols ;
113
114    return ;
115}
116
117/*
118  Helper routine to solve a continuous relaxation and print something
119  intelligent when the result is other than optimal. Returns true if the
120  result is optimal, false otherwise.
121*/
122
123bool solveRelaxation (CbcModel *model)
124
125{
126    OsiSolverInterface *osi = model->solver() ;
127
128    model->initialSolve() ;
129
130    if (!(osi->isProvenOptimal())) {
131        bool reason = false ;
132        if (osi->isProvenPrimalInfeasible()) {
133            std::cout
134                << "Continuous relaxation is primal infeasible." << std::endl ;
135            reason = true ;
136        }
137        if (osi->isProvenDualInfeasible()) {
138            std::cout
139                << "Continuous relaxation is dual infeasible." << std::endl ;
140            reason = true ;
141        }
142        if (osi->isIterationLimitReached()) {
143            std::cout
144                << "Continuous solver reached iteration limit." << std::endl ;
145            reason = true ;
146        }
147        if (osi->isAbandoned()) {
148            std::cout
149                << "Continuous solver abandoned the problem." << std::endl ;
150            reason = true ;
151        }
152        if (reason == false) {
153            std::cout
154                << "Continuous solver failed for unknown reason." << std::endl ;
155        }
156        return (false) ;
157    }
158
159    return (true) ;
160}
161
162
163/*
164  Helper routine to establish a priority vector.
165*/
166
167void setupPriorities (CbcModel *model, CbcGenCtlBlk::BPControl how)
168
169{
170    int numCols = model->getNumCols() ;
171    int *sort = new int[numCols] ;
172    double *dsort = new double[numCols] ;
173    int *priority = new int[numCols] ;
174    const double *objective = model->getObjCoefficients() ;
175    int iColumn ;
176    int n = 0 ;
177    bool priorityOK = true ;
178
179    for (iColumn = 0 ; iColumn < numCols ; iColumn++) {
180        if (model->isInteger(iColumn)) {
181            sort[n] = n ;
182            if (how == CbcGenCtlBlk::BPCost) {
183                dsort[n++] = -objective[iColumn] ;
184            } else if (how == CbcGenCtlBlk::BPOrder) {
185                dsort[n++] = iColumn ;
186            } else {
187                std::cerr
188                    << "setupPriorities: Unrecognised priority specification."
189                    << std::endl ;
190                priorityOK = false ;
191            }
192        }
193    }
194
195    if (priorityOK) {
196        CoinSort_2(dsort, dsort + n, sort) ;
197
198        int level = 0 ;
199        double last = -1.0e100 ;
200        for (int i = 0 ; i < n ; i++) {
201            int iPut = sort[i] ;
202            if (dsort[i] != last) {
203                level++ ;
204                last = dsort[i] ;
205            }
206            priority[iPut] = level ;
207        }
208
209        model->passInPriorities(priority, false) ;
210    }
211
212    delete [] priority ;
213    delete [] sort ;
214    delete [] dsort ;
215
216    return ;
217}
218
219
220/*
221  Helper routine to install a batch of heuristics. Each call to getXXXHeuristic
222  will return a pointer to the heuristic object in gen iff the heuristic is
223  enabled.
224*/
225
226void installHeuristics (CbcGenCtlBlk *ctlBlk, CbcModel *model)
227
228{
229    CbcGenCtlBlk::CGControl action ;
230    CbcHeuristic *gen ;
231    CbcTreeLocal *localTree ;
232    /*
233      FPump goes first because it only works before there's a solution.
234    */
235    action = ctlBlk->getFPump(gen, model) ;
236    if (action != CbcGenCtlBlk::CGOff) {
237        model->addHeuristic(gen, "FPump") ;
238    }
239    action = ctlBlk->getRounding(gen, model) ;
240    if (action != CbcGenCtlBlk::CGOff) {
241        model->addHeuristic(gen, "Rounding") ;
242    }
243    action = ctlBlk->getCombine(gen, model) ;
244    if (action != CbcGenCtlBlk::CGOff) {
245        model->addHeuristic(gen, "Combine") ;
246    }
247    action = ctlBlk->getGreedyCover(gen, model) ;
248    if (action != CbcGenCtlBlk::CGOff) {
249        model->addHeuristic(gen, "GCov") ;
250    }
251    action = ctlBlk->getGreedyEquality(gen, model) ;
252    if (action != CbcGenCtlBlk::CGOff) {
253        model->addHeuristic(gen, "GEq") ;
254    }
255    /*
256      This one's a bit different. We acquire the local tree and install it in the
257      model.
258    */
259    action = ctlBlk->getTreeLocal(localTree, model) ;
260    if (action != CbcGenCtlBlk::CGOff) {
261        model->passInTreeHandler(*localTree) ;
262    }
263
264    return ;
265}
266
267
268/*
269  Helper routine to install cut generators.
270
271  I need to install the new lift-and-project generator (LandP). Also need to
272  figure out stored cuts.
273*/
274
275void installCutGenerators (CbcGenCtlBlk *ctlBlk, CbcModel *model)
276
277{
278    int switches[20] ;
279    int genCnt = 0 ;
280    CbcGenCtlBlk::CGControl action ;
281    CglCutGenerator *gen ;
282
283    /*
284      The magic numbers for the howOften parameter that determines how often the
285      generator is invoked. -100 is disabled, -99 is root only, -98 will stay
286      active only so long as it generates cuts that improve the objective. A value
287      1 <= k <= 90 means the generator will be called every k nodes. If k is
288      negative, then it can be switched off if unproductive. If k is positive,
289      it'll carry on regardless.
290    */
291    int howOften[CbcGenCtlBlk::CGMarker] ;
292    howOften[CbcGenCtlBlk::CGOff] = -100 ;
293    howOften[CbcGenCtlBlk::CGOn] = -1 ;
294    howOften[CbcGenCtlBlk::CGRoot] = -99 ;
295    howOften[CbcGenCtlBlk::CGIfMove] = -98 ;
296    howOften[CbcGenCtlBlk::CGForceOn] = 1 ;
297    howOften[CbcGenCtlBlk::CGForceBut] = 1 ;
298
299    /*
300      A negative value for rowCuts means that the specified actions happen only at
301      the root.
302    */
303    action = ctlBlk->getProbing(gen) ;
304    if (action != CbcGenCtlBlk::CGOff) {
305        if (action == CbcGenCtlBlk::CGForceBut) {
306            CglProbing *probingGen = dynamic_cast<CglProbing *>(gen) ;
307            probingGen->setRowCuts(-3) ;
308        }
309        model->addCutGenerator(gen, howOften[action], "Probing") ;
310        switches[genCnt++] = 0 ;
311    }
312    action = ctlBlk->getGomory(gen) ;
313    if (action != CbcGenCtlBlk::CGOff) {
314        model->addCutGenerator(gen, howOften[action], "Gomory") ;
315        switches[genCnt++] = -1 ;
316    }
317    action = ctlBlk->getKnapsack(gen) ;
318    if (action != CbcGenCtlBlk::CGOff) {
319        model->addCutGenerator(gen, howOften[action], "Knapsack") ;
320        switches[genCnt++] = 0 ;
321    }
322    action = ctlBlk->getRedSplit(gen) ;
323    if (action != CbcGenCtlBlk::CGOff) {
324        model->addCutGenerator(gen, howOften[action], "RedSplit") ;
325        switches[genCnt++] = 1 ;
326    }
327    action = ctlBlk->getClique(gen) ;
328    if (action != CbcGenCtlBlk::CGOff) {
329        model->addCutGenerator(gen, howOften[action], "Clique") ;
330        switches[genCnt++] = 0 ;
331    }
332    action = ctlBlk->getMir(gen) ;
333    if (action != CbcGenCtlBlk::CGOff) {
334        model->addCutGenerator(gen, howOften[action], "MIR2") ;
335        switches[genCnt++] = -1 ;
336    }
337    action = ctlBlk->getFlow(gen) ;
338    if (action != CbcGenCtlBlk::CGOff) {
339        model->addCutGenerator(gen, howOften[action], "Flow") ;
340        switches[genCnt++] = 1 ;
341    }
342    action = ctlBlk->getTwomir(gen) ;
343    if (action != CbcGenCtlBlk::CGOff) {
344        model->addCutGenerator(gen, howOften[action], "2-MIR") ;
345        switches[genCnt++] = 1 ;
346    }
347    /*
348      Set control parameters on cut generators. cutDepth says `use this generator
349      when (depth in tree) mod cutDepth == 0'. setSwitchOffIfLessThan says `switch
350      this generator off if the number of cuts at the root is less than the given
351      value'. Sort of. I need to document the magic numbers for howOften , etc.
352    */
353    genCnt = model->numberCutGenerators() ;
354    int iGen ;
355    for (iGen = 0 ; iGen < genCnt ; iGen++) {
356        CbcCutGenerator *generator = model->cutGenerator(iGen) ;
357        int howOften = generator->howOften() ;
358        if (howOften == -98 || howOften == -99) {
359            generator->setSwitchOffIfLessThan(switches[iGen]) ;
360        }
361        generator->setTiming(true) ;
362        int cutDepth = ctlBlk->getCutDepth() ;
363        if (cutDepth >= 0) {
364            generator->setWhatDepth(cutDepth) ;
365        }
366    }
367    /*
368      Now some additional control parameters that affect cut generation activity.
369
370      Minimum drop is the minimum objective degradation required to continue with
371      cut passes.  We want at least .05 unless the objective is tiny, in which
372      case we'll drop down to a floor of .0001.
373    */
374    {
375        double objFrac = fabs(model->getMinimizationObjValue()) * .001 + .0001 ;
376        double minDrop = CoinMin(.05, objFrac) ;
377        model->setMinimumDrop(minDrop) ;
378    }
379    /*
380      Set the maximum number of rounds of cut generation at the root and at nodes
381      in the tree. If the value is positive, cut generation will terminate early
382      if the objective degradation doesn't meet the minimum drop requirement. If
383      the value is negatie, minimum drop is not considered.
384
385      At the root, for small problems, push for 100 passes (really we're betting
386      that we'll stop because no cuts were generated). For medium size problems,
387      the same but say we can quit if we're not achieving the minimum drop.  For
388      big problems, cut the number of rounds to 20.  The user may have expressed
389      an opinion; if so, it's already set.
390
391      Once we're in the tree, aim for one pass per activation.
392    */
393    if (ctlBlk->setByUser_[CbcCbcParam::CUTPASS] == false) {
394        int numCols = model->getNumCols() ;
395        if (numCols < 500)
396            model->setMaximumCutPassesAtRoot(-100) ;
397        else if (numCols < 5000)
398            model->setMaximumCutPassesAtRoot(100) ;
399        else
400            model->setMaximumCutPassesAtRoot(20) ;
401    }
402
403    model->setMaximumCutPasses(1) ;
404
405    return ;
406}
407
408/*
409  Install `objects' (integers, SOS sets, etc.) in the OSI. Cribbed from
410  CoinSolve 061216 and subjected to moderate rewriting. A substantial amount
411  of code that's only relevant for AMPL has been deleted. We're only supporting
412  OsiObjects in cbc-generic.
413*/
414
415void setupObjects (OsiSolverInterface *osi,
416                   bool didIPP, CglPreProcess *ippObj)
417
418{
419    int numInts = osi->getNumIntegers() ;
420    int numObjs = osi->numberObjects() ;
421    /*
422      Does this OSI have defined objects already? If not, we'd best define the
423      basic integer objects.
424    */
425    if (numInts == 0 || numObjs == 0) {
426        osi->findIntegers(false) ;
427        numInts = osi->getNumIntegers() ;
428        numObjs = osi->numberObjects() ;
429    }
430    /*
431      If we did preprocessing and discovered SOS sets, create SOS objects and
432      install them in the OSI. The priority of SOS objects is set so that larger
433      sets have higher (lower numeric value) priority. The priority of the
434      original objects is reset to be lower than the priority of any SOS object.
435      Since the SOS objects are copied into the OSI, we need to delete our
436      originals once they've been installed.
437
438      It's not clear to me that this is the right thing to do, particularly if
439      the OSI comes equipped with complex objects.  -- lh, 061216 --
440    */
441    if (didIPP && ippObj->numberSOS()) {
442        OsiObject **oldObjs = osi->objects() ;
443        int numCols = osi->getNumCols() ;
444
445        for (int iObj = 0 ; iObj < numObjs ; iObj++) {
446            oldObjs[iObj]->setPriority(numCols + 1) ;
447        }
448
449        int numSOS = ippObj->numberSOS() ;
450        OsiObject **sosObjs = new OsiObject *[numSOS] ;
451        const int *starts = ippObj->startSOS() ;
452        const int *which = ippObj->whichSOS() ;
453        const int *type = ippObj->typeSOS() ;
454        const double *weight = ippObj->weightSOS() ;
455        int iSOS ;
456        for (iSOS = 0 ; iSOS < numSOS ; iSOS++) {
457            int iStart = starts[iSOS] ;
458            int sosLen = starts[iSOS+1] - iStart ;
459            sosObjs[iSOS] =
460                new OsiSOS(osi, sosLen, which + iStart, weight + iStart, type[iSOS]) ;
461            sosObjs[iSOS]->setPriority(numCols - sosLen) ;
462        }
463        osi->addObjects(numSOS, sosObjs) ;
464
465        for (iSOS = 0 ; iSOS < numSOS ; iSOS++)
466            delete sosObjs[iSOS] ;
467        delete [] sosObjs ;
468    }
469
470    return ;
471}
472
473} // end local namespace
474
475
476namespace CbcGenParamUtils {
477
478/*
479  Run branch-and-cut.
480*/
481
482int doBaCParam (CoinParam *param)
483
484{
485    assert (param != 0) ;
486    CbcGenParam *genParam = dynamic_cast<CbcGenParam *>(param) ;
487    assert (genParam != 0) ;
488    CbcGenCtlBlk *ctlBlk = genParam->obj() ;
489    assert (ctlBlk != 0) ;
490    CbcModel *model = ctlBlk->model_ ;
491    assert (model != 0) ;
492    /*
493      Setup to return nonfatal/fatal error (1/-1) by default.
494    */
495    int retval ;
496    if (CoinParamUtils::isInteractive()) {
497        retval = 1 ;
498    } else {
499        retval = -1 ;
500    }
501    ctlBlk->setBaBStatus(CbcGenCtlBlk::BACAbandon, CbcGenCtlBlk::BACmInvalid,
502                         CbcGenCtlBlk::BACwNotStarted, false, 0) ;
503    /*
504      We ain't gonna do squat without a good model.
505    */
506    if (!ctlBlk->goodModel_) {
507        std::cout << "** Current model not valid!" << std::endl ;
508        return (retval) ;
509    }
510    /*
511      Start the clock ticking.
512    */
513    double time1 = CoinCpuTime() ;
514    double time2 ;
515    /*
516      Create a clone of the model which we can modify with impunity. Extract
517      the underlying solver for convenient access.
518    */
519    CbcModel babModel(*model) ;
520    OsiSolverInterface *babSolver = babModel.solver() ;
521    assert (babSolver != 0) ;
522# if CBC_TRACK_SOLVERS > 0
523    std::cout
524        << "doBaCParam: initial babSolver is "
525        << std::hex << babSolver << std::dec
526        << ", log level " << babSolver->messageHandler()->logLevel()
527        << "." << std::endl ;
528# endif
529    /*
530      Solve the root relaxation. Bail unless it solves to optimality.
531    */
532    if (!solveRelaxation(&babModel)) {
533        ctlBlk->setBaBStatus(&babModel, CbcGenCtlBlk::BACwBareRoot) ;
534        return (0) ;
535    }
536# if COIN_CBC_VERBOSITY > 0
537    std::cout
538        << "doBaCParam: initial relaxation z = "
539        << babSolver->getObjValue() << "." << std::endl ;
540# endif
541    /*
542      Are we up for fixing variables based on reduced cost alone?
543    */
544    if (ctlBlk->djFix_.action_ == true) {
545        reducedCostHack(babSolver, ctlBlk->djFix_.threshold_) ;
546    }
547    /*
548      Time to consider preprocessing. We'll do a bit of setup before getting to
549      the meat of the issue.
550
551      preIppSolver will hold a clone of the unpreprocessed constraint system.
552      We'll need it when we postprocess. ippSolver holds the preprocessed
553      constraint system.  Again, we clone it and give the clone to babModel for
554      B&C. Presumably we need an unmodified copy of the preprocessed system to
555      do postprocessing, but the copy itself is hidden inside the preprocess
556      object.
557    */
558    OsiSolverInterface *preIppSolver = 0 ;
559    CglPreProcess ippObj ;
560    bool didIPP = false ;
561
562    int numberChanged = 0 ;
563    int numberOriginalColumns = babSolver->getNumCols() ;
564    CbcGenCtlBlk::IPPControl ippAction = ctlBlk->getIPPAction() ;
565
566    if (!(ippAction == CbcGenCtlBlk::IPPOff ||
567            ippAction == CbcGenCtlBlk::IPPStrategy)) {
568        double timeLeft = babModel.getMaximumSeconds() ;
569        preIppSolver = babSolver->clone() ;
570        OsiSolverInterface *ippSolver ;
571#   if CBC_TRACK_SOLVERS > 0
572        std::cout
573            << "doBaCParam: clone made prior to IPP is "
574            << std::hex << preIppSolver << std::dec
575            << ", log level " << preIppSolver->messageHandler()->logLevel()
576            << "." << std::endl ;
577#   endif
578
579        preIppSolver->setHintParam(OsiDoInBranchAndCut, true, OsiHintDo) ;
580        ippObj.messageHandler()->setLogLevel(babModel.logLevel()) ;
581
582        CglProbing probingGen ;
583        probingGen.setUsingObjective(true) ;
584        probingGen.setMaxPass(3) ;
585        probingGen.setMaxProbeRoot(preIppSolver->getNumCols()) ;
586        probingGen.setMaxElements(100) ;
587        probingGen.setMaxLookRoot(50) ;
588        probingGen.setRowCuts(3) ;
589        ippObj.addCutGenerator(&probingGen) ;
590        /*
591          For preProcessNonDefault, the 2nd parameter controls the conversion of
592          clique and SOS constraints. 0 does nothing, -1 converts <= to ==, and
593          2 and 3 form SOS sets under strict and not-so-strict conditions,
594          respectively.
595        */
596        int convert = 0 ;
597        if (ippAction == CbcGenCtlBlk::IPPEqual) {
598            convert = -1 ;
599        } else if (ippAction == CbcGenCtlBlk::IPPEqualAll) {
600            convert = -2 ;
601        } else if (ippAction == CbcGenCtlBlk::IPPSOS) {
602            convert = 2 ;
603        } else if (ippAction == CbcGenCtlBlk::IPPTrySOS) {
604            convert = 3 ;
605        }
606
607        ippSolver = ippObj.preProcessNonDefault(*preIppSolver, convert, 10) ;
608#   if CBC_TRACK_SOLVERS > 0
609        std::cout
610            << "doBaCParam: solver returned from IPP is "
611            << std::hex << ippSolver << std::dec ;
612        if (ippSolver) {
613            std::cout
614                << ", log level " << ippSolver->messageHandler()->logLevel() ;
615        }
616        std::cout << "." << std::endl ;
617#   endif
618        /*
619          ippSolver == 0 is success of a sort --- integer preprocess has found the
620          problem to be infeasible or unbounded. Need to think about how to indicate
621          status.
622        */
623        if (!ippSolver) {
624            std::cout
625                << "Integer preprocess says infeasible or unbounded" << std::endl ;
626            delete preIppSolver ;
627            ctlBlk->setBaBStatus(&babModel, CbcGenCtlBlk::BACwIPP) ;
628            return (0) ;
629        }
630#   if COIN_CBC_VERBOSITY > 0
631        else {
632            std::cout
633                << "After integer preprocessing, model has "
634                << ippSolver->getNumRows()
635                << " rows, " << ippSolver->getNumCols() << " columns, and "
636                << ippSolver->getNumElements() << " elements." << std::endl ;
637        }
638#   endif
639
640        preIppSolver->setHintParam(OsiDoInBranchAndCut, false, OsiHintDo) ;
641        ippSolver->setHintParam(OsiDoInBranchAndCut, false, OsiHintDo) ;
642
643        if (ippAction == CbcGenCtlBlk::IPPSave) {
644            ippSolver->writeMps("presolved", "mps", 1.0) ;
645            std::cout
646                << "Integer preprocessed model written to `presolved.mps' "
647                << "as minimisation problem." << std::endl ;
648        }
649
650        OsiSolverInterface *osiTmp = ippSolver->clone() ;
651        babModel.assignSolver(osiTmp) ;
652        babSolver = babModel.solver() ;
653#   if CBC_TRACK_SOLVERS > 0
654        std::cout
655            << "doBaCParam: clone of IPP solver passed to babModel is "
656            << std::hex << babSolver << std::dec
657            << ", log level " << babSolver->messageHandler()->logLevel()
658            << "." << std::endl ;
659#   endif
660        if (!solveRelaxation(&babModel)) {
661            delete preIppSolver ;
662            ctlBlk->setBaBStatus(&babModel, CbcGenCtlBlk::BACwIPPRelax) ;
663            return (0) ;
664        }
665#   if COIN_CBC_VERBOSITY > 0
666        std::cout
667            << "doBaCParam: presolved relaxation z = "
668            << babSolver->getObjValue() << "." << std::endl ;
669#   endif
670        babModel.setMaximumSeconds(timeLeft - (CoinCpuTime() - time1)) ;
671        didIPP = true ;
672    }
673    /*
674      At this point, babModel and babSolver hold the constraint system we'll use
675      for B&C (either the original system or the preprocessed system) and we have
676      a solution to the lp relaxation.
677
678      If we're using the COSTSTRATEGY option, set up priorities here and pass
679      them to the babModel.
680    */
681    if (ctlBlk->priorityAction_ != CbcGenCtlBlk::BPOff) {
682        setupPriorities(&babModel, ctlBlk->priorityAction_) ;
683    }
684    /*
685      Install heuristics and cutting planes.
686    */
687    installHeuristics(ctlBlk, &babModel) ;
688    installCutGenerators(ctlBlk, &babModel) ;
689    /*
690      Set up status print frequency for babModel.
691    */
692    if (babModel.getNumCols() > 2000 || babModel.getNumRows() > 1500 ||
693            babModel.messageHandler()->logLevel() > 1)
694        babModel.setPrintFrequency(100) ;
695    /*
696      If we've read in a known good solution for debugging, activate the row cut
697      debugger.
698    */
699    if (ctlBlk->debugSol_.values_) {
700        if (ctlBlk->debugSol_.numCols_ == babModel.getNumCols()) {
701            babSolver->activateRowCutDebugger(ctlBlk->debugSol_.values_) ;
702        } else {
703            std::cout
704                << "doBaCParam: debug file has incorrect number of columns."
705                << std::endl ;
706        }
707    }
708    /*
709      Set ratio-based integrality gap, if specified by user.
710    */
711    if (ctlBlk->setByUser_[CbcCbcParam::GAPRATIO] == true) {
712        double obj = babSolver->getObjValue() ;
713        double gapRatio = babModel.getDblParam(CbcModel::CbcAllowableFractionGap) ;
714        double gap = gapRatio * (1.0e-5 + fabs(obj)) ;
715        babModel.setAllowableGap(gap) ;
716        std::cout
717            << "doBaCParam: Continuous objective = " << obj
718            << ", so allowable gap set to " << gap << std::endl ;
719    }
720    /*
721      A bit of mystery code. As best I can figure, setSpecialOptions(2) suppresses
722      the removal of warm start information when checkSolution runs an lp to check
723      a solution. John's comment, ``probably faster to use a basis to get integer
724      solutions'' makes some sense in this context. Didn't try to track down
725      moreMipOptions just yet.
726    */
727    babModel.setSpecialOptions(babModel.specialOptions() | 2) ;
728    /*
729      { int ndx = whichParam(MOREMIPOPTIONS,numberParameters,parameters) ;
730        int moreMipOptions = parameters[ndx].intValue() ;
731        if (moreMipOptions >= 0)
732        { printf("more mip options %d\n",moreMipOptions);
733          babModel.setSearchStrategy(moreMipOptions); } }
734    */
735    /*
736      Begin the final run-up to branch-and-cut.
737
738      Make sure that objects are set up in the solver. It's possible that whoever
739      loaded the model into the solver also set up objects. But it's also
740      entirely likely that none exist to this point (and interesting to note that
741      IPP doesn't need to know anything about objects).
742    */
743    setupObjects(babSolver, didIPP, &ippObj) ;
744    /*
745      Set the branching method. We can't do this until we establish objects,
746      because the constructor will set up arrays based on the number of objects,
747      and there's no provision to set this information after creation. Arguably not
748      good --- it'd be nice to set this in the prototype model that's cloned for
749      this routine. In CoinSolve, shadowPriceMode is handled with the TESTOSI
750      option.
751    */
752    OsiChooseStrong strong(babSolver) ;
753    strong.setNumberBeforeTrusted(babModel.numberBeforeTrust()) ;
754    strong.setNumberStrong(babModel.numberStrong()) ;
755    strong.setShadowPriceMode(ctlBlk->chooseStrong_.shadowPriceMode_) ;
756    CbcBranchDefaultDecision decision ;
757    decision.setChooseMethod(strong) ;
758    babModel.setBranchingMethod(decision) ;
759    /*
760      Here I've deleted a huge block of code that deals with external priorities,
761      branch direction, pseudocosts, and solution. (PRIORITYIN) Also a block of
762      code that generates C++ code.
763    */
764    /*
765      Set up strategy for branch-and-cut. Note that the integer code supplied to
766      setupPreProcessing is *not* compatible with the IPPAction enum. But at least
767      it's documented. See desiredPreProcess_ in CbcStrategyDefault. `1' is
768      accidentally equivalent to IPPOn.
769    */
770
771    if (ippAction == CbcGenCtlBlk::IPPStrategy) {
772        CbcStrategyDefault strategy(true, 5, 5) ;
773        strategy.setupPreProcessing(1) ;
774        babModel.setStrategy(strategy) ;
775    }
776    /*
777      Yes! At long last, we're ready for the big call. Do branch and cut. In
778      general, the solver used to return the solution will not be the solver we
779      passed in, so reset babSolver here.
780    */
781    int statistics = (ctlBlk->printOpt_ > 0) ? ctlBlk->printOpt_ : 0 ;
782# if CBC_TRACK_SOLVERS > 0
783    std::cout
784        << "doBaCParam: solver at call to branchAndBound is "
785        << std::hex << babModel.solver() << std::dec
786        << ", log level " << babModel.solver()->messageHandler()->logLevel()
787        << "." << std::endl ;
788# endif
789    babModel.branchAndBound(statistics) ;
790    babSolver = babModel.solver() ;
791# if CBC_TRACK_SOLVERS > 0
792    std::cout
793        << "doBaCParam: solver at return from branchAndBound is "
794        << std::hex << babModel.solver() << std::dec
795        << ", log level " << babModel.solver()->messageHandler()->logLevel()
796        << "." << std::endl ;
797# endif
798    /*
799      Write out solution to preprocessed model.
800    */
801    if (ctlBlk->debugCreate_ == "createAfterPre" &&
802            babModel.bestSolution()) {
803        CbcGenParamUtils::saveSolution(babSolver, "debug.file") ;
804    }
805    /*
806      Print some information about branch-and-cut.
807    */
808# if COIN_CBC_VERBOSITY > 0
809    std::cout
810        << "Cuts at root node changed objective from "
811        << babModel.getContinuousObjective()
812        << " to " << babModel.rootObjectiveAfterCuts() << std::endl ;
813
814    for (int iGen = 0 ; iGen < babModel.numberCutGenerators() ; iGen++) {
815        CbcCutGenerator *generator = babModel.cutGenerator(iGen) ;
816        std::cout
817            << generator->cutGeneratorName() << " was tried "
818            << generator->numberTimesEntered() << " times and created "
819            << generator->numberCutsInTotal() << " cuts of which "
820            << generator->numberCutsActive()
821            << " were active after adding rounds of cuts" ;
822        if (generator->timing()) {
823            std::cout << " ( " << generator->timeInCutGenerator() << " seconds)" ;
824        }
825        std::cout << std::endl ;
826    }
827# endif
828
829    time2 = CoinCpuTime();
830    ctlBlk->totalTime_ += time2 - time1;
831    /*
832      If we performed integer preprocessing, time to back it out.
833    */
834    if (ippAction != CbcGenCtlBlk::IPPOff) {
835#   if CBC_TRACK_SOLVERS > 0
836        std::cout
837            << "doBaCParam: solver passed to IPP postprocess is "
838            << std::hex << babSolver << std::dec << "." << std::endl ;
839#   endif
840        ippObj.postProcess(*babSolver);
841        babModel.assignSolver(preIppSolver) ;
842        babSolver = babModel.solver() ;
843#   if CBC_TRACK_SOLVERS > 0
844        std::cout
845            << "doBaCParam: solver in babModel after IPP postprocess is "
846            << std::hex << babSolver << std::dec << "." << std::endl ;
847#   endif
848    }
849    /*
850      Write out postprocessed solution to debug file, if requested.
851    */
852    if (ctlBlk->debugCreate_ == "create" && babModel.bestSolution()) {
853        CbcGenParamUtils::saveSolution(babSolver, "debug.file") ;
854    }
855    /*
856      If we have a good solution, detach the solver with the answer. Fill in the
857      rest of the status information for the benefit of the wider world.
858    */
859    bool keepAnswerSolver = false ;
860    OsiSolverInterface *answerSolver = 0 ;
861    if (babModel.bestSolution()) {
862        babModel.setModelOwnsSolver(false) ;
863        keepAnswerSolver = true ;
864        answerSolver = babSolver ;
865    }
866    ctlBlk->setBaBStatus(&babModel, CbcGenCtlBlk::BACwBAC,
867                         keepAnswerSolver, answerSolver) ;
868    /*
869      And one last bit of information & statistics.
870    */
871    ctlBlk->printBaBStatus() ;
872    std::cout << "    " ;
873    if (keepAnswerSolver) {
874        std::cout
875            << "objective " << babModel.getObjValue() << "; " ;
876    }
877    std::cout
878        << babModel.getNodeCount() << " nodes and "
879        << babModel.getIterationCount() << " iterations - took "
880        << time2 - time1 << " seconds" << std::endl ;
881
882    return (0) ;
883}
884
885} // end namespace CbcGenParamutils
886
Note: See TracBrowser for help on using the repository browser.