1 | // Copyright (C) 2003, International Business Machines |
---|
2 | // Corporation and others. All Rights Reserved. |
---|
3 | #if defined(_MSC_VER) |
---|
4 | // Turn off compiler warning about long names |
---|
5 | # pragma warning(disable:4786) |
---|
6 | #endif |
---|
7 | #include "CbcConfig.h" |
---|
8 | #include <cassert> |
---|
9 | #include <cstdlib> |
---|
10 | #include <cmath> |
---|
11 | #include <cfloat> |
---|
12 | |
---|
13 | #ifdef COIN_HAS_CLP |
---|
14 | #include "OsiClpSolverInterface.hpp" |
---|
15 | #else |
---|
16 | #include "OsiSolverInterface.hpp" |
---|
17 | #endif |
---|
18 | #include "CbcModel.hpp" |
---|
19 | #include "CbcMessage.hpp" |
---|
20 | #include "CbcCutGenerator.hpp" |
---|
21 | #include "CbcBranchDynamic.hpp" |
---|
22 | #include "CglProbing.hpp" |
---|
23 | #include "CoinTime.hpp" |
---|
24 | |
---|
25 | // Default Constructor |
---|
26 | CbcCutGenerator::CbcCutGenerator () |
---|
27 | : model_(NULL), |
---|
28 | generator_(NULL), |
---|
29 | whenCutGenerator_(-1), |
---|
30 | whenCutGeneratorInSub_(-100), |
---|
31 | switchOffIfLessThan_(0), |
---|
32 | depthCutGenerator_(-1), |
---|
33 | depthCutGeneratorInSub_(-1), |
---|
34 | generatorName_(NULL), |
---|
35 | normal_(true), |
---|
36 | atSolution_(false), |
---|
37 | whenInfeasible_(false), |
---|
38 | mustCallAgain_(false), |
---|
39 | switchedOff_(false), |
---|
40 | timing_(false), |
---|
41 | timeInCutGenerator_(0.0), |
---|
42 | numberTimes_(0), |
---|
43 | numberCuts_(0), |
---|
44 | numberColumnCuts_(0), |
---|
45 | numberCutsActive_(0), |
---|
46 | numberCutsAtRoot_(0), |
---|
47 | numberActiveCutsAtRoot_(0) |
---|
48 | { |
---|
49 | } |
---|
50 | // Normal constructor |
---|
51 | CbcCutGenerator::CbcCutGenerator(CbcModel * model,CglCutGenerator * generator, |
---|
52 | int howOften, const char * name, |
---|
53 | bool normal, bool atSolution, |
---|
54 | bool infeasible, int howOftenInSub, |
---|
55 | int whatDepth, int whatDepthInSub, |
---|
56 | int switchOffIfLessThan) |
---|
57 | : |
---|
58 | depthCutGenerator_(whatDepth), |
---|
59 | depthCutGeneratorInSub_(whatDepthInSub), |
---|
60 | mustCallAgain_(false), |
---|
61 | switchedOff_(false), |
---|
62 | timing_(false), |
---|
63 | timeInCutGenerator_(0.0), |
---|
64 | numberTimes_(0), |
---|
65 | numberCuts_(0), |
---|
66 | numberColumnCuts_(0), |
---|
67 | numberCutsActive_(0), |
---|
68 | numberCutsAtRoot_(0), |
---|
69 | numberActiveCutsAtRoot_(0) |
---|
70 | { |
---|
71 | model_ = model; |
---|
72 | generator_=generator->clone(); |
---|
73 | generator_->refreshSolver(model_->solver()); |
---|
74 | whenCutGenerator_=howOften; |
---|
75 | whenCutGeneratorInSub_ = howOftenInSub; |
---|
76 | switchOffIfLessThan_=switchOffIfLessThan; |
---|
77 | if (name) |
---|
78 | generatorName_=strdup(name); |
---|
79 | else |
---|
80 | generatorName_ = strdup("Unknown"); |
---|
81 | normal_=normal; |
---|
82 | atSolution_=atSolution; |
---|
83 | whenInfeasible_=infeasible; |
---|
84 | } |
---|
85 | |
---|
86 | // Copy constructor |
---|
87 | CbcCutGenerator::CbcCutGenerator ( const CbcCutGenerator & rhs) |
---|
88 | { |
---|
89 | model_ = rhs.model_; |
---|
90 | generator_=rhs.generator_->clone(); |
---|
91 | //generator_->refreshSolver(model_->solver()); |
---|
92 | whenCutGenerator_=rhs.whenCutGenerator_; |
---|
93 | whenCutGeneratorInSub_ = rhs.whenCutGeneratorInSub_; |
---|
94 | switchOffIfLessThan_ = rhs.switchOffIfLessThan_; |
---|
95 | depthCutGenerator_=rhs.depthCutGenerator_; |
---|
96 | depthCutGeneratorInSub_ = rhs.depthCutGeneratorInSub_; |
---|
97 | generatorName_=strdup(rhs.generatorName_); |
---|
98 | normal_=rhs.normal_; |
---|
99 | atSolution_=rhs.atSolution_; |
---|
100 | whenInfeasible_=rhs.whenInfeasible_; |
---|
101 | mustCallAgain_ = rhs.mustCallAgain_; |
---|
102 | switchedOff_ = rhs.switchedOff_; |
---|
103 | timing_ = rhs.timing_; |
---|
104 | timeInCutGenerator_ = rhs.timeInCutGenerator_; |
---|
105 | numberTimes_ = rhs.numberTimes_; |
---|
106 | numberCuts_ = rhs.numberCuts_; |
---|
107 | numberColumnCuts_ = rhs.numberColumnCuts_; |
---|
108 | numberCutsActive_ = rhs.numberCutsActive_; |
---|
109 | numberCutsAtRoot_ = rhs.numberCutsAtRoot_; |
---|
110 | numberActiveCutsAtRoot_ = rhs.numberActiveCutsAtRoot_; |
---|
111 | } |
---|
112 | |
---|
113 | // Assignment operator |
---|
114 | CbcCutGenerator & |
---|
115 | CbcCutGenerator::operator=( const CbcCutGenerator& rhs) |
---|
116 | { |
---|
117 | if (this!=&rhs) { |
---|
118 | delete generator_; |
---|
119 | free(generatorName_); |
---|
120 | model_ = rhs.model_; |
---|
121 | generator_=rhs.generator_->clone(); |
---|
122 | generator_->refreshSolver(model_->solver()); |
---|
123 | whenCutGenerator_=rhs.whenCutGenerator_; |
---|
124 | whenCutGeneratorInSub_ = rhs.whenCutGeneratorInSub_; |
---|
125 | switchOffIfLessThan_ = rhs.switchOffIfLessThan_; |
---|
126 | depthCutGenerator_=rhs.depthCutGenerator_; |
---|
127 | depthCutGeneratorInSub_ = rhs.depthCutGeneratorInSub_; |
---|
128 | generatorName_=strdup(rhs.generatorName_); |
---|
129 | normal_=rhs.normal_; |
---|
130 | atSolution_=rhs.atSolution_; |
---|
131 | whenInfeasible_=rhs.whenInfeasible_; |
---|
132 | mustCallAgain_ = rhs.mustCallAgain_; |
---|
133 | switchedOff_ = rhs.switchedOff_; |
---|
134 | timing_ = rhs.timing_; |
---|
135 | timeInCutGenerator_ = rhs.timeInCutGenerator_; |
---|
136 | numberTimes_ = rhs.numberTimes_; |
---|
137 | numberCuts_ = rhs.numberCuts_; |
---|
138 | numberColumnCuts_ = rhs.numberColumnCuts_; |
---|
139 | numberCutsActive_ = rhs.numberCutsActive_; |
---|
140 | numberCutsAtRoot_ = rhs.numberCutsAtRoot_; |
---|
141 | numberActiveCutsAtRoot_ = rhs.numberActiveCutsAtRoot_; |
---|
142 | } |
---|
143 | return *this; |
---|
144 | } |
---|
145 | |
---|
146 | // Destructor |
---|
147 | CbcCutGenerator::~CbcCutGenerator () |
---|
148 | { |
---|
149 | free(generatorName_); |
---|
150 | delete generator_; |
---|
151 | } |
---|
152 | |
---|
153 | /* This is used to refresh any inforamtion. |
---|
154 | It also refreshes the solver in the cut generator |
---|
155 | in case generator wants to do some work |
---|
156 | */ |
---|
157 | void |
---|
158 | CbcCutGenerator::refreshModel(CbcModel * model) |
---|
159 | { |
---|
160 | model_=model; |
---|
161 | generator_->refreshSolver(model_->solver()); |
---|
162 | } |
---|
163 | /* Generate cuts for the model data contained in si. |
---|
164 | The generated cuts are inserted into and returned in the |
---|
165 | collection of cuts cs. |
---|
166 | */ |
---|
167 | bool |
---|
168 | CbcCutGenerator::generateCuts( OsiCuts & cs , int fullScan, OsiSolverInterface * solver, CbcNode * node) |
---|
169 | { |
---|
170 | int howOften = whenCutGenerator_; |
---|
171 | if (howOften==-100) |
---|
172 | return false; |
---|
173 | if (howOften>0) |
---|
174 | howOften = howOften % 1000000; |
---|
175 | else |
---|
176 | howOften=1; |
---|
177 | if (!howOften) |
---|
178 | howOften=1; |
---|
179 | bool returnCode=false; |
---|
180 | //OsiSolverInterface * solver = model_->solver(); |
---|
181 | int depth; |
---|
182 | if (node) |
---|
183 | depth=node->depth(); |
---|
184 | else |
---|
185 | depth=0; |
---|
186 | int pass=model_->getCurrentPassNumber()-1; |
---|
187 | bool doThis=(model_->getNodeCount()%howOften)==0; |
---|
188 | CoinThreadRandom * randomNumberGenerator=NULL; |
---|
189 | #ifdef COIN_HAS_CLP |
---|
190 | { |
---|
191 | OsiClpSolverInterface * clpSolver |
---|
192 | = dynamic_cast<OsiClpSolverInterface *> (solver); |
---|
193 | if (clpSolver) |
---|
194 | randomNumberGenerator = clpSolver->getModelPtr()->randomNumberGenerator(); |
---|
195 | } |
---|
196 | #endif |
---|
197 | if (depthCutGenerator_>0) { |
---|
198 | doThis = (depth % depthCutGenerator_) ==0; |
---|
199 | if (depth<depthCutGenerator_) |
---|
200 | doThis=true; // and also at top of tree |
---|
201 | } |
---|
202 | // But turn off if 100 |
---|
203 | if (howOften==100) |
---|
204 | doThis=false; |
---|
205 | // Switch off if special setting |
---|
206 | if (whenCutGeneratorInSub_==-200) { |
---|
207 | fullScan=0; |
---|
208 | doThis=false; |
---|
209 | } |
---|
210 | if (fullScan||doThis) { |
---|
211 | double time1=0.0; |
---|
212 | if (timing_) |
---|
213 | time1 = CoinCpuTime(); |
---|
214 | //#define CBC_DEBUG |
---|
215 | int numberRowCutsBefore = cs.sizeRowCuts() ; |
---|
216 | int numberColumnCutsBefore = cs.sizeColCuts() ; |
---|
217 | #if 0 |
---|
218 | int cutsBefore = cs.sizeCuts(); |
---|
219 | #endif |
---|
220 | CglTreeInfo info; |
---|
221 | info.level = depth; |
---|
222 | info.pass = pass; |
---|
223 | info.formulation_rows = model_->numberRowsAtContinuous(); |
---|
224 | info.inTree = node!=NULL; |
---|
225 | info.randomNumberGenerator=randomNumberGenerator; |
---|
226 | incrementNumberTimesEntered(); |
---|
227 | CglProbing* generator = |
---|
228 | dynamic_cast<CglProbing*>(generator_); |
---|
229 | if (!generator) { |
---|
230 | // Pass across model information in case it could be useful |
---|
231 | //void * saveData = solver->getApplicationData(); |
---|
232 | //solver->setApplicationData(model_); |
---|
233 | generator_->generateCuts(*solver,cs,info); |
---|
234 | //solver->setApplicationData(saveData); |
---|
235 | } else { |
---|
236 | // Probing - return tight column bounds |
---|
237 | CglTreeProbingInfo * info2 = model_->probingInfo(); |
---|
238 | if (info2&&!depth) { |
---|
239 | info2->level = depth; |
---|
240 | info2->pass = pass; |
---|
241 | info2->formulation_rows = model_->numberRowsAtContinuous(); |
---|
242 | info2->inTree = node!=NULL; |
---|
243 | info2->randomNumberGenerator=randomNumberGenerator; |
---|
244 | generator->generateCutsAndModify(*solver,cs,info2); |
---|
245 | } else { |
---|
246 | if ((numberTimes_==200||(numberTimes_>200&&(numberTimes_%2000)==0)) |
---|
247 | &&!model_->parentModel()&&false) { |
---|
248 | // in tree, maxStack, maxProbe |
---|
249 | int test[]= { |
---|
250 | 100123, |
---|
251 | 199999, |
---|
252 | 200123, |
---|
253 | 299999, |
---|
254 | 500123, |
---|
255 | 599999, |
---|
256 | 1000123, |
---|
257 | 1099999, |
---|
258 | 2000123, |
---|
259 | 2099999}; |
---|
260 | int n = (int) (sizeof(test)/sizeof(int)); |
---|
261 | int saveStack = generator->getMaxLook(); |
---|
262 | int saveNumber = generator->getMaxProbe(); |
---|
263 | int kr1=0; |
---|
264 | int kc1=0; |
---|
265 | int bestStackTree=-1; |
---|
266 | int bestNumberTree=-1; |
---|
267 | for (int i=0;i<n;i++) { |
---|
268 | OsiCuts cs2 = cs; |
---|
269 | int stack = test[i]/100000; |
---|
270 | int number = test[i] - 100000*stack; |
---|
271 | generator->setMaxLook(stack); |
---|
272 | generator->setMaxProbe(number); |
---|
273 | generator_->generateCuts(*solver,cs2,info); |
---|
274 | int numberRowCuts = cs2.sizeRowCuts()-numberRowCutsBefore ; |
---|
275 | int numberColumnCuts= cs2.sizeColCuts()-numberColumnCutsBefore ; |
---|
276 | if (numberRowCuts<kr1||numberColumnCuts<kc1) |
---|
277 | printf("Odd "); |
---|
278 | if (numberRowCuts>kr1||numberColumnCuts>kc1) { |
---|
279 | printf("*** "); |
---|
280 | kr1=numberRowCuts; |
---|
281 | kc1=numberColumnCuts; |
---|
282 | bestStackTree=stack; |
---|
283 | bestNumberTree=number; |
---|
284 | } |
---|
285 | printf("maxStack %d number %d gives %d row cuts and %d column cuts\n", |
---|
286 | stack,number,numberRowCuts,numberColumnCuts); |
---|
287 | } |
---|
288 | generator->setMaxLook(saveStack); |
---|
289 | generator->setMaxProbe(saveNumber); |
---|
290 | if (bestStackTree>0) { |
---|
291 | generator->setMaxLook(bestStackTree); |
---|
292 | generator->setMaxProbe(bestNumberTree); |
---|
293 | printf("RRNumber %d -> %d, stack %d -> %d\n", |
---|
294 | saveNumber,bestNumberTree,saveStack,bestStackTree); |
---|
295 | } else { |
---|
296 | // no good |
---|
297 | generator->setMaxLook(1); |
---|
298 | printf("RRSwitching off number %d -> %d, stack %d -> %d\n", |
---|
299 | saveNumber,saveNumber,saveStack,1); |
---|
300 | } |
---|
301 | } |
---|
302 | if (generator->getMaxLook()>0) |
---|
303 | generator->generateCutsAndModify(*solver,cs,&info); |
---|
304 | } |
---|
305 | const double * tightLower = generator->tightLower(); |
---|
306 | const double * lower = solver->getColLower(); |
---|
307 | const double * tightUpper = generator->tightUpper(); |
---|
308 | const double * upper = solver->getColUpper(); |
---|
309 | const double * solution = solver->getColSolution(); |
---|
310 | int j; |
---|
311 | int numberColumns = solver->getNumCols(); |
---|
312 | double primalTolerance = 1.0e-8; |
---|
313 | const char * tightenBounds = generator->tightenBounds(); |
---|
314 | if ((model_->getThreadMode()&2)==0) { |
---|
315 | for (j=0;j<numberColumns;j++) { |
---|
316 | if (solver->isInteger(j)) { |
---|
317 | if (tightUpper[j]<upper[j]) { |
---|
318 | double nearest = floor(tightUpper[j]+0.5); |
---|
319 | //assert (fabs(tightUpper[j]-nearest)<1.0e-5); may be infeasible |
---|
320 | solver->setColUpper(j,nearest); |
---|
321 | if (nearest<solution[j]-primalTolerance) |
---|
322 | returnCode=true; |
---|
323 | } |
---|
324 | if (tightLower[j]>lower[j]) { |
---|
325 | double nearest = floor(tightLower[j]+0.5); |
---|
326 | //assert (fabs(tightLower[j]-nearest)<1.0e-5); may be infeasible |
---|
327 | solver->setColLower(j,nearest); |
---|
328 | if (nearest>solution[j]+primalTolerance) |
---|
329 | returnCode=true; |
---|
330 | } |
---|
331 | } else { |
---|
332 | if (upper[j]>lower[j]) { |
---|
333 | if (tightUpper[j]==tightLower[j]) { |
---|
334 | // fix |
---|
335 | solver->setColLower(j,tightLower[j]); |
---|
336 | solver->setColUpper(j,tightUpper[j]); |
---|
337 | if (tightLower[j]>solution[j]+primalTolerance|| |
---|
338 | tightUpper[j]<solution[j]-primalTolerance) |
---|
339 | returnCode=true; |
---|
340 | } else if (tightenBounds&&tightenBounds[j]) { |
---|
341 | solver->setColLower(j,CoinMax(tightLower[j],lower[j])); |
---|
342 | solver->setColUpper(j,CoinMin(tightUpper[j],upper[j])); |
---|
343 | if (tightLower[j]>solution[j]+primalTolerance|| |
---|
344 | tightUpper[j]<solution[j]-primalTolerance) |
---|
345 | returnCode=true; |
---|
346 | } |
---|
347 | } |
---|
348 | } |
---|
349 | } |
---|
350 | } else { |
---|
351 | CoinPackedVector lbs; |
---|
352 | CoinPackedVector ubs; |
---|
353 | int numberChanged=0; |
---|
354 | bool ifCut=false; |
---|
355 | for (j=0;j<numberColumns;j++) { |
---|
356 | if (solver->isInteger(j)) { |
---|
357 | if (tightUpper[j]<upper[j]) { |
---|
358 | double nearest = floor(tightUpper[j]+0.5); |
---|
359 | //assert (fabs(tightUpper[j]-nearest)<1.0e-5); may be infeasible |
---|
360 | ubs.insert(j,nearest); |
---|
361 | numberChanged++; |
---|
362 | if (nearest<solution[j]-primalTolerance) |
---|
363 | ifCut=true; |
---|
364 | } |
---|
365 | if (tightLower[j]>lower[j]) { |
---|
366 | double nearest = floor(tightLower[j]+0.5); |
---|
367 | //assert (fabs(tightLower[j]-nearest)<1.0e-5); may be infeasible |
---|
368 | lbs.insert(j,nearest); |
---|
369 | numberChanged++; |
---|
370 | if (nearest>solution[j]+primalTolerance) |
---|
371 | ifCut=true; |
---|
372 | } |
---|
373 | } else { |
---|
374 | if (upper[j]>lower[j]) { |
---|
375 | if (tightUpper[j]==tightLower[j]) { |
---|
376 | // fix |
---|
377 | lbs.insert(j,tightLower[j]); |
---|
378 | ubs.insert(j,tightUpper[j]); |
---|
379 | if (tightLower[j]>solution[j]+primalTolerance|| |
---|
380 | tightUpper[j]<solution[j]-primalTolerance) |
---|
381 | ifCut=true; |
---|
382 | } else if (tightenBounds&&tightenBounds[j]) { |
---|
383 | lbs.insert(j,CoinMax(tightLower[j],lower[j])); |
---|
384 | ubs.insert(j,CoinMin(tightUpper[j],upper[j])); |
---|
385 | if (tightLower[j]>solution[j]+primalTolerance|| |
---|
386 | tightUpper[j]<solution[j]-primalTolerance) |
---|
387 | ifCut=true; |
---|
388 | } |
---|
389 | } |
---|
390 | } |
---|
391 | } |
---|
392 | if (numberChanged) { |
---|
393 | OsiColCut cc; |
---|
394 | cc.setUbs(ubs); |
---|
395 | cc.setLbs(lbs); |
---|
396 | if (ifCut) { |
---|
397 | cc.setEffectiveness(100.0); |
---|
398 | } else { |
---|
399 | cc.setEffectiveness(1.0e-5); |
---|
400 | } |
---|
401 | cs.insert(cc); |
---|
402 | } |
---|
403 | } |
---|
404 | //if (!solver->basisIsAvailable()) |
---|
405 | //returnCode=true; |
---|
406 | #if 0 |
---|
407 | // Pass across info to pseudocosts |
---|
408 | char * mark = new char[numberColumns]; |
---|
409 | memset(mark,0,numberColumns); |
---|
410 | int nLook = generator->numberThisTime(); |
---|
411 | const int * lookedAt = generator->lookedAt(); |
---|
412 | const int * fixedDown = generator->fixedDown(); |
---|
413 | const int * fixedUp = generator->fixedUp(); |
---|
414 | for (j=0;j<nLook;j++) |
---|
415 | mark[lookedAt[j]]=1; |
---|
416 | int numberObjects = model_->numberObjects(); |
---|
417 | for (int i=0;i<numberObjects;i++) { |
---|
418 | CbcSimpleIntegerDynamicPseudoCost * obj1 = |
---|
419 | dynamic_cast <CbcSimpleIntegerDynamicPseudoCost *>(model_->modifiableObject(i)) ; |
---|
420 | if (obj1) { |
---|
421 | int iColumn = obj1->columnNumber(); |
---|
422 | if (mark[iColumn]) |
---|
423 | obj1->setProbingInformation(fixedDown[iColumn],fixedUp[iColumn]); |
---|
424 | } |
---|
425 | } |
---|
426 | delete [] mark; |
---|
427 | #endif |
---|
428 | } |
---|
429 | CbcCutModifier * modifier = model_->cutModifier(); |
---|
430 | if (modifier) { |
---|
431 | int numberRowCutsAfter = cs.sizeRowCuts() ; |
---|
432 | int k ; |
---|
433 | int nOdd=0; |
---|
434 | //const OsiSolverInterface * solver = model_->solver(); |
---|
435 | for (k = numberRowCutsAfter-1;k>=numberRowCutsBefore;k--) { |
---|
436 | OsiRowCut & thisCut = cs.rowCut(k) ; |
---|
437 | int returnCode = modifier->modify(solver,thisCut); |
---|
438 | if (returnCode) { |
---|
439 | nOdd++; |
---|
440 | if (returnCode==3) |
---|
441 | cs.eraseRowCut(k); |
---|
442 | } |
---|
443 | } |
---|
444 | if (nOdd) |
---|
445 | printf("Cut generator %s produced %d cuts of which %d were modified\n", |
---|
446 | generatorName_,numberRowCutsAfter-numberRowCutsBefore,nOdd); |
---|
447 | } |
---|
448 | { |
---|
449 | // make all row cuts without test for duplicate |
---|
450 | int numberRowCutsAfter = cs.sizeRowCuts() ; |
---|
451 | int k ; |
---|
452 | for (k = numberRowCutsBefore;k<numberRowCutsAfter;k++) { |
---|
453 | OsiRowCut * thisCut = cs.rowCutPtr(k) ; |
---|
454 | thisCut->mutableRow().setTestForDuplicateIndex(false); |
---|
455 | } |
---|
456 | } |
---|
457 | { |
---|
458 | int numberRowCutsAfter = cs.sizeRowCuts() ; |
---|
459 | if (numberRowCutsBefore<numberRowCutsAfter) { |
---|
460 | #if 0 |
---|
461 | printf("generator %s generated %d row cuts\n", |
---|
462 | generatorName_,numberRowCutsAfter-numberRowCutsBefore); |
---|
463 | #endif |
---|
464 | numberCuts_ += numberRowCutsAfter-numberRowCutsBefore; |
---|
465 | } |
---|
466 | int numberColumnCutsAfter = cs.sizeColCuts() ; |
---|
467 | if (numberColumnCutsBefore<numberColumnCutsAfter) { |
---|
468 | #if 0 |
---|
469 | printf("generator %s generated %d column cuts\n", |
---|
470 | generatorName_,numberColumnCutsAfter-numberColumnCutsBefore); |
---|
471 | #endif |
---|
472 | numberColumnCuts_ += numberColumnCutsAfter-numberColumnCutsBefore; |
---|
473 | } |
---|
474 | } |
---|
475 | { |
---|
476 | int numberRowCutsAfter = cs.sizeRowCuts() ; |
---|
477 | int k ; |
---|
478 | int nEls=0; |
---|
479 | int nCuts= numberRowCutsAfter-numberRowCutsBefore; |
---|
480 | // Remove NULL cuts! |
---|
481 | int nNull=0; |
---|
482 | const double * solution = solver->getColSolution(); |
---|
483 | bool feasible=true; |
---|
484 | for (k = numberRowCutsAfter-1;k>=numberRowCutsBefore;k--) { |
---|
485 | const OsiRowCut * thisCut = cs.rowCutPtr(k) ; |
---|
486 | double sum=0.0; |
---|
487 | if (thisCut->lb()<=thisCut->ub()) { |
---|
488 | int n=thisCut->row().getNumElements(); |
---|
489 | const int * column = thisCut->row().getIndices(); |
---|
490 | const double * element = thisCut->row().getElements(); |
---|
491 | if (n<=0) { |
---|
492 | // infeasible cut - give up |
---|
493 | feasible=false; |
---|
494 | break; |
---|
495 | } |
---|
496 | nEls+= n; |
---|
497 | for (int i=0;i<n;i++) { |
---|
498 | double value = element[i]; |
---|
499 | sum += value*solution[column[i]]; |
---|
500 | } |
---|
501 | if (sum>thisCut->ub()) { |
---|
502 | sum= sum-thisCut->ub(); |
---|
503 | } else if (sum<thisCut->lb()) { |
---|
504 | sum= thisCut->lb()-sum; |
---|
505 | } else { |
---|
506 | sum=0.0; |
---|
507 | cs.eraseRowCut(k); |
---|
508 | nNull++; |
---|
509 | } |
---|
510 | } |
---|
511 | } |
---|
512 | //if (nNull) |
---|
513 | //printf("%s has %d cuts and %d elements - %d null!\n",generatorName_, |
---|
514 | // nCuts,nEls,nNull); |
---|
515 | numberRowCutsAfter = cs.sizeRowCuts() ; |
---|
516 | nCuts= numberRowCutsAfter-numberRowCutsBefore; |
---|
517 | nEls=0; |
---|
518 | for (k = numberRowCutsBefore;k<numberRowCutsAfter;k++) { |
---|
519 | const OsiRowCut * thisCut = cs.rowCutPtr(k) ; |
---|
520 | int n=thisCut->row().getNumElements(); |
---|
521 | nEls+= n; |
---|
522 | } |
---|
523 | //printf("%s has %d cuts and %d elements\n",generatorName_, |
---|
524 | // nCuts,nEls); |
---|
525 | int nElsNow = solver->getMatrixByCol()->getNumElements(); |
---|
526 | int nAdd = model_->parentModel() ? 200 : 10000; |
---|
527 | int numberColumns = solver->getNumCols(); |
---|
528 | int nAdd2 = model_->parentModel() ? 2*numberColumns : 5*numberColumns; |
---|
529 | if (/*nEls>CoinMax(nAdd2,nElsNow/8+nAdd)*/nCuts&&feasible) { |
---|
530 | //printf("need to remove cuts\n"); |
---|
531 | // just add most effective |
---|
532 | int nReasonable = CoinMax(nAdd2,nElsNow/8+nAdd); |
---|
533 | int nDelete = nEls - nReasonable; |
---|
534 | |
---|
535 | nElsNow = nEls; |
---|
536 | double * sort = new double [nCuts]; |
---|
537 | int * which = new int [nCuts]; |
---|
538 | // For parallel cuts |
---|
539 | double * element2 = new double [numberColumns]; |
---|
540 | CoinZeroN(element2,numberColumns); |
---|
541 | for (k = numberRowCutsBefore;k<numberRowCutsAfter;k++) { |
---|
542 | const OsiRowCut * thisCut = cs.rowCutPtr(k) ; |
---|
543 | double sum=0.0; |
---|
544 | if (thisCut->lb()<=thisCut->ub()) { |
---|
545 | int n=thisCut->row().getNumElements(); |
---|
546 | const int * column = thisCut->row().getIndices(); |
---|
547 | const double * element = thisCut->row().getElements(); |
---|
548 | assert (n); |
---|
549 | double norm=0.0; |
---|
550 | for (int i=0;i<n;i++) { |
---|
551 | double value = element[i]; |
---|
552 | sum += value*solution[column[i]]; |
---|
553 | norm += value*value; |
---|
554 | } |
---|
555 | if (sum>thisCut->ub()) { |
---|
556 | sum= sum-thisCut->ub(); |
---|
557 | } else if (sum<thisCut->lb()) { |
---|
558 | sum= thisCut->lb()-sum; |
---|
559 | } else { |
---|
560 | sum=0.0; |
---|
561 | } |
---|
562 | // normalize |
---|
563 | sum /= sqrt(norm); |
---|
564 | // adjust for length |
---|
565 | //sum /= sqrt((double) n); |
---|
566 | // randomize |
---|
567 | //double randomNumber = |
---|
568 | //model_->randomNumberGenerator()->randomDouble(); |
---|
569 | //sum *= (0.5+randomNumber); |
---|
570 | } else { |
---|
571 | // keep |
---|
572 | sum=COIN_DBL_MAX; |
---|
573 | } |
---|
574 | sort[k-numberRowCutsBefore]=sum; |
---|
575 | which[k-numberRowCutsBefore]=k; |
---|
576 | } |
---|
577 | CoinSort_2(sort,sort+nCuts,which); |
---|
578 | // Now see which ones are too similar |
---|
579 | int nParallel=0; |
---|
580 | for (k = 0;k<nCuts;k++) { |
---|
581 | int j=which[k]; |
---|
582 | const OsiRowCut * thisCut = cs.rowCutPtr(j) ; |
---|
583 | if (thisCut->lb()>thisCut->ub()) |
---|
584 | break; // cut is infeasible |
---|
585 | int n=thisCut->row().getNumElements(); |
---|
586 | const int * column = thisCut->row().getIndices(); |
---|
587 | const double * element = thisCut->row().getElements(); |
---|
588 | assert (n); |
---|
589 | double norm=0.0; |
---|
590 | double lb = thisCut->lb(); |
---|
591 | double ub = thisCut->ub(); |
---|
592 | for (int i=0;i<n;i++) { |
---|
593 | double value = element[i]; |
---|
594 | element2[column[i]]=value; |
---|
595 | norm += value*value; |
---|
596 | } |
---|
597 | int kkk = CoinMin(nCuts,k+5); |
---|
598 | for (int kk=k+1;kk<kkk;kk++) { |
---|
599 | int jj=which[kk]; |
---|
600 | const OsiRowCut * thisCut2 = cs.rowCutPtr(jj) ; |
---|
601 | if (thisCut2->lb()>thisCut2->ub()) |
---|
602 | break; // cut is infeasible |
---|
603 | int nB=thisCut2->row().getNumElements(); |
---|
604 | const int * columnB = thisCut2->row().getIndices(); |
---|
605 | const double * elementB = thisCut2->row().getElements(); |
---|
606 | assert (nB); |
---|
607 | double normB=0.0; |
---|
608 | double product=0.0; |
---|
609 | for (int i=0;i<nB;i++) { |
---|
610 | double value = elementB[i]; |
---|
611 | normB += value*value; |
---|
612 | product += value*element2[columnB[i]]; |
---|
613 | } |
---|
614 | if (product>0.0&&product*product>0.99999*norm*normB) { |
---|
615 | bool parallel=true; |
---|
616 | double lbB = thisCut2->lb(); |
---|
617 | double ubB = thisCut2->ub(); |
---|
618 | if ((lb<-1.0e20&&lbB>-1.0e20)|| |
---|
619 | (lbB<-1.0e20&&lb>-1.0e20)) |
---|
620 | parallel = false; |
---|
621 | double tolerance; |
---|
622 | tolerance = CoinMax(fabs(lb),fabs(lbB))+1.0e-6; |
---|
623 | if (fabs(lb-lbB)>tolerance) |
---|
624 | parallel=false; |
---|
625 | if ((ub>1.0e20&&ubB<1.0e20)|| |
---|
626 | (ubB>1.0e20&&ub<1.0e20)) |
---|
627 | parallel = false; |
---|
628 | tolerance = CoinMax(fabs(ub),fabs(ubB))+1.0e-6; |
---|
629 | if (fabs(ub-ubB)>tolerance) |
---|
630 | parallel=false; |
---|
631 | if (parallel) { |
---|
632 | nParallel++; |
---|
633 | sort[k]=0.0; |
---|
634 | break; |
---|
635 | } |
---|
636 | } |
---|
637 | } |
---|
638 | for (int i=0;i<n;i++) { |
---|
639 | element2[column[i]]=0.0; |
---|
640 | } |
---|
641 | } |
---|
642 | delete [] element2; |
---|
643 | CoinSort_2(sort,sort+nCuts,which); |
---|
644 | k=0; |
---|
645 | while (nDelete>0||!sort[k]) { |
---|
646 | int iCut=which[k]; |
---|
647 | const OsiRowCut * thisCut = cs.rowCutPtr(iCut) ; |
---|
648 | int n=thisCut->row().getNumElements(); |
---|
649 | nDelete-=n; |
---|
650 | k++; |
---|
651 | if (k>=nCuts) |
---|
652 | break; |
---|
653 | } |
---|
654 | std::sort(which,which+k); |
---|
655 | k--; |
---|
656 | for (;k>=0;k--) { |
---|
657 | cs.eraseRowCut(which[k]); |
---|
658 | } |
---|
659 | delete [] sort; |
---|
660 | delete [] which; |
---|
661 | numberRowCutsAfter = cs.sizeRowCuts() ; |
---|
662 | #if 0 //def CLP_INVESTIGATE |
---|
663 | nEls=0; |
---|
664 | int nCuts2= numberRowCutsAfter-numberRowCutsBefore; |
---|
665 | for (k = numberRowCutsBefore;k<numberRowCutsAfter;k++) { |
---|
666 | const OsiRowCut * thisCut = cs.rowCutPtr(k) ; |
---|
667 | int n=thisCut->row().getNumElements(); |
---|
668 | nEls+= n; |
---|
669 | } |
---|
670 | if (!model_->parentModel()&&nCuts!=nCuts2) |
---|
671 | printf("%s NOW has %d cuts and %d elements( down from %d cuts and %d els) - %d parallel\n", |
---|
672 | generatorName_, |
---|
673 | nCuts2,nEls,nCuts,nElsNow,nParallel); |
---|
674 | #endif |
---|
675 | } |
---|
676 | } |
---|
677 | #ifdef CBC_DEBUG |
---|
678 | { |
---|
679 | int numberRowCutsAfter = cs.sizeRowCuts() ; |
---|
680 | int k ; |
---|
681 | int nBad=0; |
---|
682 | for (k = numberRowCutsBefore;k<numberRowCutsAfter;k++) { |
---|
683 | OsiRowCut thisCut = cs.rowCut(k) ; |
---|
684 | if (thisCut.lb()>thisCut.ub()|| |
---|
685 | thisCut.lb()>1.0e8|| |
---|
686 | thisCut.ub()<-1.0e8) |
---|
687 | printf("cut from %s has bounds %g and %g!\n", |
---|
688 | generatorName_,thisCut.lb(),thisCut.ub()); |
---|
689 | if (thisCut.lb()<=thisCut.ub()) { |
---|
690 | /* check size of elements. |
---|
691 | We can allow smaller but this helps debug generators as it |
---|
692 | is unsafe to have small elements */ |
---|
693 | int n=thisCut.row().getNumElements(); |
---|
694 | const int * column = thisCut.row().getIndices(); |
---|
695 | const double * element = thisCut.row().getElements(); |
---|
696 | assert (n); |
---|
697 | for (int i=0;i<n;i++) { |
---|
698 | double value = element[i]; |
---|
699 | if (fabs(value)<=1.0e-12||fabs(value)>=1.0e20) |
---|
700 | nBad++; |
---|
701 | } |
---|
702 | } |
---|
703 | if (nBad) |
---|
704 | printf("Cut generator %s produced %d cuts of which %d had tiny or large elements\n", |
---|
705 | generatorName_,numberRowCutsAfter-numberRowCutsBefore,nBad); |
---|
706 | } |
---|
707 | } |
---|
708 | #endif |
---|
709 | if (timing_) |
---|
710 | timeInCutGenerator_ += CoinCpuTime()-time1; |
---|
711 | #if 0 |
---|
712 | // switch off if first time and no good |
---|
713 | if (node==NULL&&!pass) { |
---|
714 | if (cs.sizeCuts()-cutsBefore<CoinAbs(switchOffIfLessThan_)) { |
---|
715 | whenCutGenerator_=-99; |
---|
716 | whenCutGeneratorInSub_ = -200; |
---|
717 | } |
---|
718 | } |
---|
719 | #endif |
---|
720 | } |
---|
721 | return returnCode; |
---|
722 | } |
---|
723 | void |
---|
724 | CbcCutGenerator::setHowOften(int howOften) |
---|
725 | { |
---|
726 | |
---|
727 | if (howOften>=1000000) { |
---|
728 | // leave Probing every SCANCUTS_PROBING |
---|
729 | howOften = howOften % 1000000; |
---|
730 | CglProbing* generator = |
---|
731 | dynamic_cast<CglProbing*>(generator_); |
---|
732 | |
---|
733 | if (generator&&howOften>SCANCUTS_PROBING) |
---|
734 | howOften=SCANCUTS_PROBING+1000000; |
---|
735 | else |
---|
736 | howOften += 1000000; |
---|
737 | } |
---|
738 | whenCutGenerator_ = howOften; |
---|
739 | } |
---|
740 | void |
---|
741 | CbcCutGenerator::setWhatDepth(int value) |
---|
742 | { |
---|
743 | depthCutGenerator_ = value; |
---|
744 | } |
---|
745 | void |
---|
746 | CbcCutGenerator::setWhatDepthInSub(int value) |
---|
747 | { |
---|
748 | depthCutGeneratorInSub_ = value; |
---|
749 | } |
---|
750 | |
---|
751 | |
---|
752 | // Default Constructor |
---|
753 | CbcCutModifier::CbcCutModifier() |
---|
754 | { |
---|
755 | } |
---|
756 | |
---|
757 | |
---|
758 | // Destructor |
---|
759 | CbcCutModifier::~CbcCutModifier () |
---|
760 | { |
---|
761 | } |
---|
762 | |
---|
763 | // Copy constructor |
---|
764 | CbcCutModifier::CbcCutModifier ( const CbcCutModifier & rhs) |
---|
765 | { |
---|
766 | } |
---|
767 | |
---|
768 | // Assignment operator |
---|
769 | CbcCutModifier & |
---|
770 | CbcCutModifier::operator=( const CbcCutModifier& rhs) |
---|
771 | { |
---|
772 | if (this!=&rhs) { |
---|
773 | } |
---|
774 | return *this; |
---|
775 | } |
---|
776 | |
---|
777 | // Default Constructor |
---|
778 | CbcCutSubsetModifier::CbcCutSubsetModifier () |
---|
779 | : CbcCutModifier(), |
---|
780 | firstOdd_(COIN_INT_MAX) |
---|
781 | { |
---|
782 | } |
---|
783 | |
---|
784 | // Useful constructor |
---|
785 | CbcCutSubsetModifier::CbcCutSubsetModifier (int firstOdd) |
---|
786 | : CbcCutModifier() |
---|
787 | { |
---|
788 | firstOdd_=firstOdd; |
---|
789 | } |
---|
790 | |
---|
791 | // Copy constructor |
---|
792 | CbcCutSubsetModifier::CbcCutSubsetModifier ( const CbcCutSubsetModifier & rhs) |
---|
793 | :CbcCutModifier(rhs) |
---|
794 | { |
---|
795 | firstOdd_ = rhs.firstOdd_; |
---|
796 | } |
---|
797 | |
---|
798 | // Clone |
---|
799 | CbcCutModifier * |
---|
800 | CbcCutSubsetModifier::clone() const |
---|
801 | { |
---|
802 | return new CbcCutSubsetModifier(*this); |
---|
803 | } |
---|
804 | |
---|
805 | // Assignment operator |
---|
806 | CbcCutSubsetModifier & |
---|
807 | CbcCutSubsetModifier::operator=( const CbcCutSubsetModifier& rhs) |
---|
808 | { |
---|
809 | if (this!=&rhs) { |
---|
810 | CbcCutModifier::operator=(rhs); |
---|
811 | firstOdd_ = rhs.firstOdd_; |
---|
812 | } |
---|
813 | return *this; |
---|
814 | } |
---|
815 | |
---|
816 | // Destructor |
---|
817 | CbcCutSubsetModifier::~CbcCutSubsetModifier () |
---|
818 | { |
---|
819 | } |
---|
820 | /* Returns |
---|
821 | 0 unchanged |
---|
822 | 1 strengthened |
---|
823 | 2 weakened |
---|
824 | 3 deleted |
---|
825 | */ |
---|
826 | int |
---|
827 | CbcCutSubsetModifier::modify(const OsiSolverInterface * solver, OsiRowCut & cut) |
---|
828 | { |
---|
829 | int n=cut.row().getNumElements(); |
---|
830 | if (!n) |
---|
831 | return 0; |
---|
832 | const int * column = cut.row().getIndices(); |
---|
833 | //const double * element = cut.row().getElements(); |
---|
834 | int returnCode=0; |
---|
835 | for (int i=0;i<n;i++) { |
---|
836 | if (column[i]>=firstOdd_) { |
---|
837 | returnCode=3; |
---|
838 | break; |
---|
839 | } |
---|
840 | } |
---|
841 | if (!returnCode) { |
---|
842 | const double * element = cut.row().getElements(); |
---|
843 | printf("%g <= ",cut.lb()); |
---|
844 | for (int i=0;i<n;i++) { |
---|
845 | printf("%g*x%d ",element[i],column[i]); |
---|
846 | } |
---|
847 | printf("<= %g\n",cut.ub()); |
---|
848 | } |
---|
849 | //return 3; |
---|
850 | return returnCode; |
---|
851 | } |
---|
852 | |
---|