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