1 | // Copyright (C) 2002, 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 <string> |
---|
8 | //#define CBC_DEBUG 1 |
---|
9 | //#define CHECK_CUT_COUNTS |
---|
10 | #include <cassert> |
---|
11 | #include <cfloat> |
---|
12 | #define CUTS |
---|
13 | #include "OsiSolverInterface.hpp" |
---|
14 | #include "CoinWarmStartBasis.hpp" |
---|
15 | #include "CoinTime.hpp" |
---|
16 | #include "CbcModel.hpp" |
---|
17 | #include "CbcNode.hpp" |
---|
18 | #include "CbcBranchActual.hpp" |
---|
19 | #include "OsiRowCut.hpp" |
---|
20 | #include "OsiRowCutDebugger.hpp" |
---|
21 | #include "OsiCuts.hpp" |
---|
22 | #include "CbcCountRowCut.hpp" |
---|
23 | #include "CbcMessage.hpp" |
---|
24 | #include "OsiClpSolverInterface.hpp" |
---|
25 | using namespace std; |
---|
26 | #include "CglCutGenerator.hpp" |
---|
27 | // Default Constructor |
---|
28 | CbcNodeInfo::CbcNodeInfo () |
---|
29 | : |
---|
30 | numberPointingToThis_(0), |
---|
31 | parent_(NULL), |
---|
32 | owner_(NULL), |
---|
33 | numberCuts_(0), |
---|
34 | cuts_(NULL), |
---|
35 | numberRows_(0), |
---|
36 | numberBranchesLeft_(0) |
---|
37 | { |
---|
38 | #ifdef CHECK_NODE |
---|
39 | printf("CbcNodeInfo %x Constructor\n",this); |
---|
40 | #endif |
---|
41 | } |
---|
42 | // Constructor given parent |
---|
43 | CbcNodeInfo::CbcNodeInfo (CbcNodeInfo * parent) |
---|
44 | : |
---|
45 | numberPointingToThis_(2), |
---|
46 | parent_(parent), |
---|
47 | owner_(NULL), |
---|
48 | numberCuts_(0), |
---|
49 | cuts_(NULL), |
---|
50 | numberRows_(0), |
---|
51 | numberBranchesLeft_(2) |
---|
52 | { |
---|
53 | #ifdef CHECK_NODE |
---|
54 | printf("CbcNodeInfo %x Constructor from parent %x\n",this,parent_); |
---|
55 | #endif |
---|
56 | if (parent_) { |
---|
57 | numberRows_ = parent_->numberRows_+parent_->numberCuts_; |
---|
58 | } |
---|
59 | } |
---|
60 | // Copy Constructor |
---|
61 | CbcNodeInfo::CbcNodeInfo (const CbcNodeInfo & rhs) |
---|
62 | : |
---|
63 | numberPointingToThis_(rhs.numberPointingToThis_), |
---|
64 | parent_(rhs.parent_), |
---|
65 | owner_(rhs.owner_), |
---|
66 | numberCuts_(rhs.numberCuts_), |
---|
67 | cuts_(NULL), |
---|
68 | numberRows_(rhs.numberRows_), |
---|
69 | numberBranchesLeft_(rhs.numberBranchesLeft_) |
---|
70 | { |
---|
71 | #ifdef CHECK_NODE |
---|
72 | printf("CbcNodeInfo %x Copy constructor\n",this); |
---|
73 | #endif |
---|
74 | if (numberCuts_) { |
---|
75 | cuts_ = new CbcCountRowCut * [numberCuts_]; |
---|
76 | int n=0; |
---|
77 | for (int i=0;i<numberCuts_;i++) { |
---|
78 | CbcCountRowCut * thisCut = rhs.cuts_[i]; |
---|
79 | if (thisCut) { |
---|
80 | // I think this is correct - new one should take priority |
---|
81 | thisCut->setInfo(this,n); |
---|
82 | thisCut->increment(numberBranchesLeft_); |
---|
83 | cuts_[n++] = thisCut; |
---|
84 | } |
---|
85 | } |
---|
86 | numberCuts_=n; |
---|
87 | } |
---|
88 | } |
---|
89 | // Constructor given parent and owner |
---|
90 | CbcNodeInfo::CbcNodeInfo (CbcNodeInfo * parent, CbcNode * owner) |
---|
91 | : |
---|
92 | numberPointingToThis_(2), |
---|
93 | parent_(parent), |
---|
94 | owner_(owner), |
---|
95 | numberCuts_(0), |
---|
96 | cuts_(NULL), |
---|
97 | numberRows_(0), |
---|
98 | numberBranchesLeft_(2) |
---|
99 | { |
---|
100 | #ifdef CHECK_NODE |
---|
101 | printf("CbcNodeInfo %x Constructor from parent %x\n",this,parent_); |
---|
102 | #endif |
---|
103 | if (parent_) { |
---|
104 | numberRows_ = parent_->numberRows_+parent_->numberCuts_; |
---|
105 | } |
---|
106 | } |
---|
107 | |
---|
108 | /** |
---|
109 | Take care to detach from the owning CbcNode and decrement the reference |
---|
110 | count in the parent. If this is the last nodeInfo object pointing to the |
---|
111 | parent, make a recursive call to delete the parent. |
---|
112 | */ |
---|
113 | CbcNodeInfo::~CbcNodeInfo() |
---|
114 | { |
---|
115 | #ifdef CHECK_NODE |
---|
116 | printf("CbcNodeInfo %x Destructor parent %x\n",this,parent_); |
---|
117 | #endif |
---|
118 | |
---|
119 | assert(!numberPointingToThis_); |
---|
120 | // But they may be some left (max nodes?) |
---|
121 | for (int i=0;i<numberCuts_;i++) |
---|
122 | delete cuts_[i]; |
---|
123 | |
---|
124 | delete [] cuts_; |
---|
125 | if (owner_) |
---|
126 | owner_->nullNodeInfo(); |
---|
127 | if (parent_) { |
---|
128 | int numberLinks = parent_->decrement(); |
---|
129 | if (!numberLinks) delete parent_; |
---|
130 | } |
---|
131 | } |
---|
132 | |
---|
133 | |
---|
134 | //#define ALLCUTS |
---|
135 | void |
---|
136 | CbcNodeInfo::decrementCuts(int change) |
---|
137 | { |
---|
138 | int i; |
---|
139 | // get rid of all remaining if negative |
---|
140 | int changeThis; |
---|
141 | if (change<0) |
---|
142 | changeThis = numberBranchesLeft_; |
---|
143 | else |
---|
144 | changeThis = change; |
---|
145 | // decrement cut counts |
---|
146 | for (i=0;i<numberCuts_;i++) { |
---|
147 | if (cuts_[i]) { |
---|
148 | int number = cuts_[i]->decrement(changeThis); |
---|
149 | if (!number) { |
---|
150 | //printf("info %x del cut %d %x\n",this,i,cuts_[i]); |
---|
151 | delete cuts_[i]; |
---|
152 | cuts_[i]=NULL; |
---|
153 | } |
---|
154 | } |
---|
155 | } |
---|
156 | } |
---|
157 | void |
---|
158 | CbcNodeInfo::decrementParentCuts(int change) |
---|
159 | { |
---|
160 | if (parent_) { |
---|
161 | // get rid of all remaining if negative |
---|
162 | int changeThis; |
---|
163 | if (change<0) |
---|
164 | changeThis = numberBranchesLeft_; |
---|
165 | else |
---|
166 | changeThis = change; |
---|
167 | int i; |
---|
168 | // Get over-estimate of space needed for basis |
---|
169 | CoinWarmStartBasis dummy; |
---|
170 | dummy.setSize(0,numberRows_+numberCuts_); |
---|
171 | buildRowBasis(dummy); |
---|
172 | /* everything is zero (i.e. free) so we can use to see |
---|
173 | if latest basis */ |
---|
174 | CbcNodeInfo * thisInfo = parent_; |
---|
175 | while (thisInfo) |
---|
176 | thisInfo = thisInfo->buildRowBasis(dummy); |
---|
177 | // decrement cut counts |
---|
178 | thisInfo = parent_; |
---|
179 | int numberRows=numberRows_; |
---|
180 | while (thisInfo) { |
---|
181 | for (i=thisInfo->numberCuts_-1;i>=0;i--) { |
---|
182 | CoinWarmStartBasis::Status status = dummy.getArtifStatus(--numberRows); |
---|
183 | #ifdef ALLCUTS |
---|
184 | status = CoinWarmStartBasis::isFree; |
---|
185 | #endif |
---|
186 | if (thisInfo->cuts_[i]) { |
---|
187 | int number=1; |
---|
188 | if (status!=CoinWarmStartBasis::basic) { |
---|
189 | // tight - drop 1 or 2 |
---|
190 | if (change<0) |
---|
191 | number = thisInfo->cuts_[i]->decrement(changeThis); |
---|
192 | else |
---|
193 | number = thisInfo->cuts_[i]->decrement(change); |
---|
194 | } |
---|
195 | if (!number) { |
---|
196 | delete thisInfo->cuts_[i]; |
---|
197 | thisInfo->cuts_[i]=NULL; |
---|
198 | } |
---|
199 | } |
---|
200 | } |
---|
201 | thisInfo = thisInfo->parent_; |
---|
202 | } |
---|
203 | } |
---|
204 | } |
---|
205 | |
---|
206 | void |
---|
207 | CbcNodeInfo::incrementParentCuts(int change) |
---|
208 | { |
---|
209 | if (parent_) { |
---|
210 | int i; |
---|
211 | // Get over-estimate of space needed for basis |
---|
212 | CoinWarmStartBasis dummy; |
---|
213 | dummy.setSize(0,numberRows_+numberCuts_); |
---|
214 | /* everything is zero (i.e. free) so we can use to see |
---|
215 | if latest basis */ |
---|
216 | buildRowBasis(dummy); |
---|
217 | CbcNodeInfo * thisInfo = parent_; |
---|
218 | while (thisInfo) |
---|
219 | thisInfo = thisInfo->buildRowBasis(dummy); |
---|
220 | // increment cut counts |
---|
221 | thisInfo = parent_; |
---|
222 | int numberRows=numberRows_; |
---|
223 | while (thisInfo) { |
---|
224 | for (i=thisInfo->numberCuts_-1;i>=0;i--) { |
---|
225 | CoinWarmStartBasis::Status status = dummy.getArtifStatus(--numberRows); |
---|
226 | #ifdef ALLCUTS |
---|
227 | status = CoinWarmStartBasis::isFree; |
---|
228 | #endif |
---|
229 | if (thisInfo->cuts_[i]&&status!=CoinWarmStartBasis::basic) { |
---|
230 | thisInfo->cuts_[i]->increment(change); |
---|
231 | } |
---|
232 | } |
---|
233 | thisInfo = thisInfo->parent_; |
---|
234 | } |
---|
235 | } |
---|
236 | } |
---|
237 | |
---|
238 | /* |
---|
239 | Append cuts to the cuts_ array in a nodeInfo. The initial reference count |
---|
240 | is set to numberToBranchOn, which will normally be the number of arms |
---|
241 | defined for the CbcBranchingObject attached to the CbcNode that owns this |
---|
242 | CbcNodeInfo. |
---|
243 | */ |
---|
244 | void |
---|
245 | CbcNodeInfo::addCuts (OsiCuts & cuts, int numberToBranchOn, |
---|
246 | int * whichGenerator) |
---|
247 | { |
---|
248 | int numberCuts = cuts.sizeRowCuts(); |
---|
249 | if (numberCuts) { |
---|
250 | int i; |
---|
251 | if (!numberCuts_) { |
---|
252 | cuts_ = new CbcCountRowCut * [numberCuts]; |
---|
253 | } else { |
---|
254 | CbcCountRowCut ** temp = new CbcCountRowCut * [numberCuts+numberCuts_]; |
---|
255 | memcpy(temp,cuts_,numberCuts_*sizeof(CbcCountRowCut *)); |
---|
256 | delete [] cuts_; |
---|
257 | cuts_ = temp; |
---|
258 | } |
---|
259 | for (i=0;i<numberCuts;i++) { |
---|
260 | CbcCountRowCut * thisCut = new CbcCountRowCut(*cuts.rowCutPtr(i), |
---|
261 | this,numberCuts_); |
---|
262 | thisCut->increment(numberToBranchOn); |
---|
263 | cuts_[numberCuts_++] = thisCut; |
---|
264 | #ifdef CBC_DEBUG |
---|
265 | int n=thisCut->row().getNumElements(); |
---|
266 | #if CBC_DEBUG>1 |
---|
267 | printf("Cut %d has %d entries, rhs %g %g =>",i,n,thisCut->lb(), |
---|
268 | thisCut->ub()); |
---|
269 | #endif |
---|
270 | int j; |
---|
271 | #if CBC_DEBUG>1 |
---|
272 | const int * index = thisCut->row().getIndices(); |
---|
273 | #endif |
---|
274 | const double * element = thisCut->row().getElements(); |
---|
275 | for (j=0;j<n;j++) { |
---|
276 | #if CBC_DEBUG>1 |
---|
277 | printf(" (%d,%g)",index[j],element[j]); |
---|
278 | #endif |
---|
279 | assert(fabs(element[j])>1.00e-12); |
---|
280 | } |
---|
281 | printf("\n"); |
---|
282 | #endif |
---|
283 | } |
---|
284 | } |
---|
285 | } |
---|
286 | |
---|
287 | void |
---|
288 | CbcNodeInfo::addCuts(int numberCuts, CbcCountRowCut ** cut, |
---|
289 | int numberToBranchOn) |
---|
290 | { |
---|
291 | if (numberCuts) { |
---|
292 | int i; |
---|
293 | if (!numberCuts_) { |
---|
294 | cuts_ = new CbcCountRowCut * [numberCuts]; |
---|
295 | } else { |
---|
296 | CbcCountRowCut ** temp = new CbcCountRowCut * [numberCuts+numberCuts_]; |
---|
297 | memcpy(temp,cuts_,numberCuts_*sizeof(CbcCountRowCut *)); |
---|
298 | delete [] cuts_; |
---|
299 | cuts_ = temp; |
---|
300 | } |
---|
301 | for (i=0;i<numberCuts;i++) { |
---|
302 | CbcCountRowCut * thisCut = cut[i]; |
---|
303 | thisCut->setInfo(this,numberCuts_); |
---|
304 | //printf("info %x cut %d %x\n",this,i,thisCut); |
---|
305 | thisCut->increment(numberToBranchOn); |
---|
306 | cuts_[numberCuts_++] = thisCut; |
---|
307 | #ifdef CBC_DEBUG |
---|
308 | int n=thisCut->row().getNumElements(); |
---|
309 | #if CBC_DEBUG>1 |
---|
310 | printf("Cut %d has %d entries, rhs %g %g =>",i,n,thisCut->lb(), |
---|
311 | thisCut->ub()); |
---|
312 | #endif |
---|
313 | int j; |
---|
314 | #if CBC_DEBUG>1 |
---|
315 | const int * index = thisCut->row().getIndices(); |
---|
316 | #endif |
---|
317 | const double * element = thisCut->row().getElements(); |
---|
318 | for (j=0;j<n;j++) { |
---|
319 | #if CBC_DEBUG>1 |
---|
320 | printf(" (%d,%g)",index[j],element[j]); |
---|
321 | #endif |
---|
322 | assert(fabs(element[j])>1.00e-12); |
---|
323 | } |
---|
324 | printf("\n"); |
---|
325 | #endif |
---|
326 | } |
---|
327 | } |
---|
328 | } |
---|
329 | |
---|
330 | // delete cuts |
---|
331 | void |
---|
332 | CbcNodeInfo::deleteCuts(int numberToDelete, CbcCountRowCut ** cuts) |
---|
333 | { |
---|
334 | int i; |
---|
335 | int j; |
---|
336 | int last=-1; |
---|
337 | for (i=0;i<numberToDelete;i++) { |
---|
338 | CbcCountRowCut * next = cuts[i]; |
---|
339 | for (j=last+1;j<numberCuts_;j++) { |
---|
340 | if (next==cuts_[j]) |
---|
341 | break; |
---|
342 | } |
---|
343 | if (j==numberCuts_) { |
---|
344 | // start from beginning |
---|
345 | for (j=0;j<last;j++) { |
---|
346 | if (next==cuts_[j]) |
---|
347 | break; |
---|
348 | } |
---|
349 | assert(j<last); |
---|
350 | } |
---|
351 | last=j; |
---|
352 | int number = cuts_[j]->decrement(); |
---|
353 | if (!number) |
---|
354 | delete cuts_[j]; |
---|
355 | cuts_[j]=NULL; |
---|
356 | } |
---|
357 | j=0; |
---|
358 | for (i=0;i<numberCuts_;i++) { |
---|
359 | if (cuts_[i]) |
---|
360 | cuts_[j++]=cuts_[i]; |
---|
361 | } |
---|
362 | numberCuts_ = j; |
---|
363 | } |
---|
364 | |
---|
365 | // delete cuts |
---|
366 | void |
---|
367 | CbcNodeInfo::deleteCuts(int numberToDelete, int * which) |
---|
368 | { |
---|
369 | int i; |
---|
370 | for (i=0;i<numberToDelete;i++) { |
---|
371 | int iCut=which[i]; |
---|
372 | int number = cuts_[iCut]->decrement(); |
---|
373 | if (!number) |
---|
374 | delete cuts_[iCut]; |
---|
375 | cuts_[iCut]=NULL; |
---|
376 | } |
---|
377 | int j=0; |
---|
378 | for (i=0;i<numberCuts_;i++) { |
---|
379 | if (cuts_[i]) |
---|
380 | cuts_[j++]=cuts_[i]; |
---|
381 | } |
---|
382 | numberCuts_ = j; |
---|
383 | } |
---|
384 | |
---|
385 | // Really delete a cut |
---|
386 | void |
---|
387 | CbcNodeInfo::deleteCut(int whichOne) |
---|
388 | { |
---|
389 | assert(whichOne<numberCuts_); |
---|
390 | cuts_[whichOne]=NULL; |
---|
391 | } |
---|
392 | |
---|
393 | CbcFullNodeInfo::CbcFullNodeInfo() : |
---|
394 | CbcNodeInfo(), |
---|
395 | basis_(), |
---|
396 | numberIntegers_(0), |
---|
397 | lower_(NULL), |
---|
398 | upper_(NULL) |
---|
399 | { |
---|
400 | } |
---|
401 | CbcFullNodeInfo::CbcFullNodeInfo(CbcModel * model, |
---|
402 | int numberRowsAtContinuous) : |
---|
403 | CbcNodeInfo() |
---|
404 | { |
---|
405 | OsiSolverInterface * solver = model->solver(); |
---|
406 | numberRows_ = numberRowsAtContinuous; |
---|
407 | numberIntegers_ = model->numberIntegers(); |
---|
408 | int numberColumns = model->getNumCols(); |
---|
409 | lower_ = new double [numberColumns]; |
---|
410 | upper_ = new double [numberColumns]; |
---|
411 | const double * lower = solver->getColLower(); |
---|
412 | const double * upper = solver->getColUpper(); |
---|
413 | int i; |
---|
414 | |
---|
415 | for (i=0;i<numberColumns;i++) { |
---|
416 | lower_[i]=lower[i]; |
---|
417 | upper_[i]=upper[i]; |
---|
418 | } |
---|
419 | |
---|
420 | basis_ = dynamic_cast<CoinWarmStartBasis*>(solver->getWarmStart()); |
---|
421 | } |
---|
422 | |
---|
423 | CbcFullNodeInfo::CbcFullNodeInfo(const CbcFullNodeInfo & rhs) : |
---|
424 | CbcNodeInfo(rhs) |
---|
425 | { |
---|
426 | basis_= dynamic_cast<CoinWarmStartBasis *>(rhs.basis_->clone()) ; |
---|
427 | numberIntegers_=rhs.numberIntegers_; |
---|
428 | lower_=NULL; |
---|
429 | upper_=NULL; |
---|
430 | if (rhs.lower_!=NULL) { |
---|
431 | int numberColumns = basis_->getNumStructural(); |
---|
432 | lower_ = new double [numberColumns]; |
---|
433 | upper_ = new double [numberColumns]; |
---|
434 | assert (upper_!=NULL); |
---|
435 | memcpy(lower_,rhs.lower_,numberColumns*sizeof(double)); |
---|
436 | memcpy(upper_,rhs.upper_,numberColumns*sizeof(double)); |
---|
437 | } |
---|
438 | } |
---|
439 | |
---|
440 | CbcNodeInfo * |
---|
441 | CbcFullNodeInfo::clone() const |
---|
442 | { |
---|
443 | return (new CbcFullNodeInfo(*this)); |
---|
444 | } |
---|
445 | |
---|
446 | CbcFullNodeInfo::~CbcFullNodeInfo () |
---|
447 | { |
---|
448 | delete basis_ ; |
---|
449 | delete [] lower_; |
---|
450 | delete [] upper_; |
---|
451 | } |
---|
452 | |
---|
453 | /* |
---|
454 | The basis supplied as a parameter is deleted and replaced with a new basis |
---|
455 | appropriate for the node, and lower and upper bounds on variables are |
---|
456 | reset according to the stored bounds arrays. Any cuts associated with this |
---|
457 | node are added to the list in addCuts, but not actually added to the |
---|
458 | constraint system in the model. |
---|
459 | |
---|
460 | Why pass in a basis at all? The short answer is ``We need the parameter to |
---|
461 | pass out a basis, so might as well use it to pass in the size.'' |
---|
462 | |
---|
463 | A longer answer is that in practice we take a memory allocation hit up in |
---|
464 | addCuts1 (the only place applyToModel is called) when we setSize() the |
---|
465 | basis that's passed in. It's immediately tossed here in favour of a clone |
---|
466 | of the basis attached to this nodeInfo. This can probably be fixed, given |
---|
467 | a bit of thought. |
---|
468 | */ |
---|
469 | |
---|
470 | void CbcFullNodeInfo::applyToModel (CbcModel *model, |
---|
471 | CoinWarmStartBasis *&basis, |
---|
472 | CbcCountRowCut **addCuts, |
---|
473 | int ¤tNumberCuts) const |
---|
474 | |
---|
475 | { OsiSolverInterface *solver = model->solver() ; |
---|
476 | |
---|
477 | // branch - do bounds |
---|
478 | int i; |
---|
479 | int numberColumns = model->getNumCols(); |
---|
480 | for (i=0;i<numberColumns;i++) { |
---|
481 | solver->setColBounds(i,lower_[i],upper_[i]); |
---|
482 | } |
---|
483 | // move basis - but make sure size stays |
---|
484 | int numberRows = basis->getNumArtificial(); |
---|
485 | delete basis ; |
---|
486 | basis = dynamic_cast<CoinWarmStartBasis *>(basis_->clone()) ; |
---|
487 | basis->resize(numberRows,numberColumns); |
---|
488 | for (i=0;i<numberCuts_;i++) |
---|
489 | addCuts[currentNumberCuts+i]= cuts_[i]; |
---|
490 | currentNumberCuts += numberCuts_; |
---|
491 | assert(!parent_); |
---|
492 | return ; |
---|
493 | } |
---|
494 | |
---|
495 | /* Builds up row basis backwards (until original model). |
---|
496 | Returns NULL or previous one to apply . |
---|
497 | Depends on Free being 0 and impossible for cuts |
---|
498 | */ |
---|
499 | CbcNodeInfo * |
---|
500 | CbcFullNodeInfo::buildRowBasis(CoinWarmStartBasis & basis ) const |
---|
501 | { |
---|
502 | const unsigned int * saved = |
---|
503 | (const unsigned int *) basis_->getArtificialStatus(); |
---|
504 | unsigned int * now = |
---|
505 | (unsigned int *) basis.getArtificialStatus(); |
---|
506 | int number=basis_->getNumArtificial()>>4;; |
---|
507 | int i; |
---|
508 | for (i=0;i<number;i++) { |
---|
509 | if (!now[i]) |
---|
510 | now[i] = saved[i]; |
---|
511 | } |
---|
512 | return NULL; |
---|
513 | } |
---|
514 | |
---|
515 | // Default constructor |
---|
516 | CbcPartialNodeInfo::CbcPartialNodeInfo() |
---|
517 | |
---|
518 | : CbcNodeInfo(), |
---|
519 | basisDiff_(NULL), |
---|
520 | variables_(NULL), |
---|
521 | newBounds_(NULL), |
---|
522 | numberChangedBounds_(0) |
---|
523 | |
---|
524 | { /* this space intentionally left blank */ } |
---|
525 | |
---|
526 | // Constructor from current state |
---|
527 | CbcPartialNodeInfo::CbcPartialNodeInfo (CbcNodeInfo *parent, CbcNode *owner, |
---|
528 | int numberChangedBounds, |
---|
529 | const int *variables, |
---|
530 | const double *boundChanges, |
---|
531 | const CoinWarmStartDiff *basisDiff) |
---|
532 | : CbcNodeInfo(parent,owner) |
---|
533 | { |
---|
534 | basisDiff_ = basisDiff->clone() ; |
---|
535 | |
---|
536 | numberChangedBounds_ = numberChangedBounds; |
---|
537 | variables_ = new int [numberChangedBounds_]; |
---|
538 | newBounds_ = new double [numberChangedBounds_]; |
---|
539 | |
---|
540 | int i ; |
---|
541 | for (i=0;i<numberChangedBounds_;i++) { |
---|
542 | variables_[i]=variables[i]; |
---|
543 | newBounds_[i]=boundChanges[i]; |
---|
544 | } |
---|
545 | } |
---|
546 | |
---|
547 | CbcPartialNodeInfo::CbcPartialNodeInfo (const CbcPartialNodeInfo & rhs) |
---|
548 | |
---|
549 | : CbcNodeInfo(rhs.parent_) |
---|
550 | |
---|
551 | { basisDiff_ = rhs.basisDiff_->clone() ; |
---|
552 | |
---|
553 | numberChangedBounds_ = rhs.numberChangedBounds_; |
---|
554 | variables_ = new int [numberChangedBounds_]; |
---|
555 | newBounds_ = new double [numberChangedBounds_]; |
---|
556 | |
---|
557 | int i ; |
---|
558 | for (i=0;i<numberChangedBounds_;i++) { |
---|
559 | variables_[i]=rhs.variables_[i]; |
---|
560 | newBounds_[i]=rhs.newBounds_[i]; |
---|
561 | } |
---|
562 | } |
---|
563 | |
---|
564 | CbcNodeInfo * |
---|
565 | CbcPartialNodeInfo::clone() const |
---|
566 | { |
---|
567 | return (new CbcPartialNodeInfo(*this)); |
---|
568 | } |
---|
569 | |
---|
570 | |
---|
571 | CbcPartialNodeInfo::~CbcPartialNodeInfo () |
---|
572 | { |
---|
573 | delete basisDiff_ ; |
---|
574 | delete [] variables_; |
---|
575 | delete [] newBounds_; |
---|
576 | } |
---|
577 | |
---|
578 | |
---|
579 | /** |
---|
580 | The basis supplied as a parameter is incrementally modified, and lower and |
---|
581 | upper bounds on variables in the model are incrementally modified. Any |
---|
582 | cuts associated with this node are added to the list in addCuts. |
---|
583 | */ |
---|
584 | |
---|
585 | void CbcPartialNodeInfo::applyToModel (CbcModel *model, |
---|
586 | CoinWarmStartBasis *&basis, |
---|
587 | CbcCountRowCut **addCuts, |
---|
588 | int ¤tNumberCuts) const |
---|
589 | |
---|
590 | { OsiSolverInterface *solver = model->solver(); |
---|
591 | |
---|
592 | basis->applyDiff(basisDiff_) ; |
---|
593 | |
---|
594 | // branch - do bounds |
---|
595 | int i; |
---|
596 | for (i=0;i<numberChangedBounds_;i++) { |
---|
597 | int variable = variables_[i]; |
---|
598 | if ((variable&0x80000000)==0) { |
---|
599 | // lower bound changing |
---|
600 | solver->setColLower(variable,newBounds_[i]); |
---|
601 | } else { |
---|
602 | // upper bound changing |
---|
603 | solver->setColUpper(variable&0x7fffffff,newBounds_[i]); |
---|
604 | } |
---|
605 | } |
---|
606 | for (i=0;i<numberCuts_;i++) |
---|
607 | addCuts[currentNumberCuts+i]= cuts_[i]; |
---|
608 | currentNumberCuts += numberCuts_; |
---|
609 | return ; |
---|
610 | } |
---|
611 | |
---|
612 | /* Builds up row basis backwards (until original model). |
---|
613 | Returns NULL or previous one to apply . |
---|
614 | Depends on Free being 0 and impossible for cuts |
---|
615 | */ |
---|
616 | |
---|
617 | CbcNodeInfo * |
---|
618 | CbcPartialNodeInfo::buildRowBasis(CoinWarmStartBasis & basis ) const |
---|
619 | |
---|
620 | { basis.applyDiff(basisDiff_) ; |
---|
621 | |
---|
622 | return parent_ ; } |
---|
623 | |
---|
624 | |
---|
625 | CbcNode::CbcNode() : |
---|
626 | nodeInfo_(NULL), |
---|
627 | objectiveValue_(1.0e100), |
---|
628 | guessedObjectiveValue_(1.0e100), |
---|
629 | branch_(NULL), |
---|
630 | depth_(-1), |
---|
631 | numberUnsatisfied_(0), |
---|
632 | nodeNumber_(-1) |
---|
633 | { |
---|
634 | #ifdef CHECK_NODE |
---|
635 | printf("CbcNode %x Constructor\n",this); |
---|
636 | #endif |
---|
637 | } |
---|
638 | |
---|
639 | CbcNode::CbcNode(CbcModel * model, |
---|
640 | CbcNode * lastNode) : |
---|
641 | nodeInfo_(NULL), |
---|
642 | objectiveValue_(1.0e100), |
---|
643 | guessedObjectiveValue_(1.0e100), |
---|
644 | branch_(NULL), |
---|
645 | depth_(-1), |
---|
646 | numberUnsatisfied_(0), |
---|
647 | nodeNumber_(model->getNodeCount()) |
---|
648 | { |
---|
649 | #ifdef CHECK_NODE |
---|
650 | printf("CbcNode %x Constructor from model\n",this); |
---|
651 | #endif |
---|
652 | OsiSolverInterface * solver = model->solver(); |
---|
653 | objectiveValue_ = solver->getObjSense()*solver->getObjValue(); |
---|
654 | |
---|
655 | if (lastNode) |
---|
656 | lastNode->nodeInfo_->increment(); |
---|
657 | } |
---|
658 | |
---|
659 | |
---|
660 | void |
---|
661 | CbcNode::createInfo (CbcModel *model, |
---|
662 | CbcNode *lastNode, |
---|
663 | const CoinWarmStartBasis *lastws, |
---|
664 | const double *lastLower, const double *lastUpper, |
---|
665 | int numberOldActiveCuts,int numberNewCuts) |
---|
666 | { OsiSolverInterface * solver = model->solver(); |
---|
667 | /* |
---|
668 | The root --- no parent. Create full basis and bounds information. |
---|
669 | */ |
---|
670 | if (!lastNode) |
---|
671 | { nodeInfo_=new CbcFullNodeInfo(model,solver->getNumRows()); } |
---|
672 | /* |
---|
673 | Not the root. Create an edit from the parent's basis & bound information. |
---|
674 | This is not quite as straightforward as it seems. We need to reintroduce |
---|
675 | cuts we may have dropped out of the basis, in the correct position, because |
---|
676 | this whole process is strictly positional. Start by grabbing the current |
---|
677 | basis. |
---|
678 | */ |
---|
679 | else |
---|
680 | { const CoinWarmStartBasis* ws = |
---|
681 | dynamic_cast<const CoinWarmStartBasis*>(solver->getWarmStart()); |
---|
682 | assert(ws!=NULL); // make sure not volume |
---|
683 | //int numberArtificials = lastws->getNumArtificial(); |
---|
684 | int numberColumns = solver->getNumCols(); |
---|
685 | |
---|
686 | const double * lower = solver->getColLower(); |
---|
687 | const double * upper = solver->getColUpper(); |
---|
688 | |
---|
689 | int i; |
---|
690 | /* |
---|
691 | Create a clone and resize it to hold all the structural constraints, plus |
---|
692 | all the cuts: old cuts, both active and inactive (currentNumberCuts), and |
---|
693 | new cuts (numberNewCuts). |
---|
694 | |
---|
695 | TODO: You'd think that the set of constraints (logicals) in the expanded |
---|
696 | basis should match the set represented in lastws. At least, that's |
---|
697 | what I thought. But at the point I first looked hard at this bit of |
---|
698 | code, it turned out that lastws was the stripped basis produced at |
---|
699 | the end of addCuts(), rather than the raw basis handed back by |
---|
700 | addCuts1(). The expanded basis here is equivalent to the raw basis of |
---|
701 | addCuts1(). I said ``whoa, that's not good, I must have introduced a |
---|
702 | bug'' and went back to John's code to see where I'd gone wrong. |
---|
703 | And discovered the same `error' in his code. |
---|
704 | |
---|
705 | After a bit of thought, my conclusion is that correctness is not |
---|
706 | affected by whether lastws is the stripped or raw basis. The diffs |
---|
707 | have no semantics --- just a set of changes that need to be made |
---|
708 | to convert lastws into expanded. I think the only effect is that we |
---|
709 | store a lot more diffs (everything in expanded that's not covered by |
---|
710 | the stripped basis). But I need to give this more thought. There |
---|
711 | may well be some subtle error cases. |
---|
712 | |
---|
713 | In the mean time, I've twiddled addCuts() to set lastws to the raw |
---|
714 | basis. Makes me (Lou) less nervous to compare apples to apples. |
---|
715 | */ |
---|
716 | CoinWarmStartBasis *expanded = |
---|
717 | dynamic_cast<CoinWarmStartBasis *>(ws->clone()) ; |
---|
718 | int numberRowsAtContinuous = model->numberRowsAtContinuous(); |
---|
719 | int iFull = numberRowsAtContinuous+model->currentNumberCuts()+ |
---|
720 | numberNewCuts; |
---|
721 | //int numberArtificialsNow = iFull; |
---|
722 | //int maxBasisLength = ((iFull+15)>>4)+((numberColumns+15)>>4); |
---|
723 | //printf("l %d full %d\n",maxBasisLength,iFull); |
---|
724 | expanded->resize(iFull,numberColumns); |
---|
725 | #ifdef FULL_DEBUG |
---|
726 | printf("Before expansion: orig %d, old %d, new %d, current %d\n", |
---|
727 | numberRowsAtContinuous,numberOldActiveCuts,numberNewCuts, |
---|
728 | model->currentNumberCuts()) ; |
---|
729 | ws->print(); |
---|
730 | #endif |
---|
731 | /* |
---|
732 | Now fill in the expanded basis. Any indices beyond nPartial must |
---|
733 | be cuts created while processing this node --- they can be copied directly |
---|
734 | into the expanded basis. From nPartial down, pull the status of active cuts |
---|
735 | from ws, interleaving with a B entry for the deactivated (loose) cuts. |
---|
736 | */ |
---|
737 | int numberDropped = model->currentNumberCuts()-numberOldActiveCuts; |
---|
738 | int iCompact=iFull-numberDropped; |
---|
739 | CbcCountRowCut ** cut = model->addedCuts(); |
---|
740 | int nPartial = model->currentNumberCuts()+numberRowsAtContinuous; |
---|
741 | iFull--; |
---|
742 | for (;iFull>=nPartial;iFull--) { |
---|
743 | CoinWarmStartBasis::Status status = ws->getArtifStatus(--iCompact); |
---|
744 | //assert (status != CoinWarmStartBasis::basic); // may be permanent cut |
---|
745 | expanded->setArtifStatus(iFull,status); |
---|
746 | } |
---|
747 | for (;iFull>=numberRowsAtContinuous;iFull--) { |
---|
748 | if (cut[iFull-numberRowsAtContinuous]) { |
---|
749 | CoinWarmStartBasis::Status status = ws->getArtifStatus(--iCompact); |
---|
750 | // If no cut generator being used then we may have basic variables |
---|
751 | //if (model->getMaximumCutPasses()&& |
---|
752 | // status == CoinWarmStartBasis::basic) |
---|
753 | //printf("cut basic\n"); |
---|
754 | expanded->setArtifStatus(iFull,status); |
---|
755 | } else { |
---|
756 | expanded->setArtifStatus(iFull,CoinWarmStartBasis::basic); |
---|
757 | } |
---|
758 | } |
---|
759 | #ifdef FULL_DEBUG |
---|
760 | printf("Expanded basis\n"); |
---|
761 | expanded->print() ; |
---|
762 | printf("Diffing against\n") ; |
---|
763 | lastws->print() ; |
---|
764 | #endif |
---|
765 | /* |
---|
766 | Now that we have two bases in proper positional correspondence, creating |
---|
767 | the actual diff is dead easy. |
---|
768 | */ |
---|
769 | |
---|
770 | CoinWarmStartDiff *basisDiff = expanded->generateDiff(lastws) ; |
---|
771 | /* |
---|
772 | Diff the bound vectors. It's assumed the number of structural variables is |
---|
773 | not changing. Assuming that branching objects all involve integer variables, |
---|
774 | we should see at least one bound change as a consequence of processing this |
---|
775 | subproblem. Different types of branching objects could break this assertion. |
---|
776 | Yes - cuts will break |
---|
777 | */ |
---|
778 | double *boundChanges = new double [2*numberColumns] ; |
---|
779 | int *variables = new int [2*numberColumns] ; |
---|
780 | int numberChangedBounds=0; |
---|
781 | for (i=0;i<numberColumns;i++) { |
---|
782 | if (lower[i]!=lastLower[i]) { |
---|
783 | variables[numberChangedBounds]=i; |
---|
784 | boundChanges[numberChangedBounds++]=lower[i]; |
---|
785 | } |
---|
786 | if (upper[i]!=lastUpper[i]) { |
---|
787 | variables[numberChangedBounds]=i|0x80000000; |
---|
788 | boundChanges[numberChangedBounds++]=upper[i]; |
---|
789 | } |
---|
790 | #ifdef CBC_DEBUG |
---|
791 | if (lower[i]!=lastLower[i]) |
---|
792 | printf("lower on %d changed from %g to %g\n", |
---|
793 | i,lastLower[i],lower[i]); |
---|
794 | if (upper[i]!=lastUpper[i]) |
---|
795 | printf("upper on %d changed from %g to %g\n", |
---|
796 | i,lastUpper[i],upper[i]); |
---|
797 | #endif |
---|
798 | } |
---|
799 | #ifdef CBC_DEBUG |
---|
800 | printf("%d changed bounds\n",numberChangedBounds) ; |
---|
801 | #endif |
---|
802 | if (lastNode->branchingObject()->boundBranch()) |
---|
803 | assert (numberChangedBounds); |
---|
804 | /* |
---|
805 | Hand the lot over to the CbcPartialNodeInfo constructor, then clean up and |
---|
806 | return. |
---|
807 | */ |
---|
808 | nodeInfo_ = |
---|
809 | new CbcPartialNodeInfo(lastNode->nodeInfo_,this,numberChangedBounds, |
---|
810 | variables,boundChanges,basisDiff) ; |
---|
811 | delete basisDiff ; |
---|
812 | delete [] boundChanges; |
---|
813 | delete [] variables; |
---|
814 | delete expanded ; |
---|
815 | delete ws; |
---|
816 | } |
---|
817 | } |
---|
818 | |
---|
819 | /* |
---|
820 | The routine scans through the object list of the model looking for objects |
---|
821 | that indicate infeasibility. It tests each object using strong branching |
---|
822 | and selects the one with the least objective degradation. A corresponding |
---|
823 | branching object is left attached to lastNode. |
---|
824 | |
---|
825 | If strong branching is disabled, a candidate object is chosen essentially |
---|
826 | at random (whatever object ends up in pos'n 0 of the candidate array). |
---|
827 | |
---|
828 | If a branching candidate is found to be monotone, bounds are set to fix the |
---|
829 | variable and the routine immediately returns (the caller is expected to |
---|
830 | reoptimize). |
---|
831 | |
---|
832 | If a branching candidate is found to result in infeasibility in both |
---|
833 | directions, the routine immediately returns an indication of infeasibility. |
---|
834 | |
---|
835 | Returns: 0 both branch directions are feasible |
---|
836 | -1 branching variable is monotone |
---|
837 | -2 infeasible |
---|
838 | |
---|
839 | Original comments: |
---|
840 | Here could go cuts etc etc |
---|
841 | For now just fix on objective from strong branching. |
---|
842 | */ |
---|
843 | |
---|
844 | int CbcNode::chooseBranch (CbcModel *model, CbcNode *lastNode,int numberPassesLeft) |
---|
845 | |
---|
846 | { if (lastNode) |
---|
847 | depth_ = lastNode->depth_+1; |
---|
848 | else |
---|
849 | depth_ = 0; |
---|
850 | delete branch_; |
---|
851 | branch_=NULL; |
---|
852 | OsiSolverInterface * solver = model->solver(); |
---|
853 | double saveObjectiveValue = solver->getObjValue(); |
---|
854 | double objectiveValue = solver->getObjSense()*saveObjectiveValue; |
---|
855 | const double * lower = solver->getColLower(); |
---|
856 | const double * upper = solver->getColUpper(); |
---|
857 | int anyAction=0; |
---|
858 | double integerTolerance = |
---|
859 | model->getDblParam(CbcModel::CbcIntegerTolerance); |
---|
860 | int i; |
---|
861 | bool beforeSolution = model->getSolutionCount()==0; |
---|
862 | int numberStrong=model->numberStrong(); |
---|
863 | int numberObjects = model->numberObjects(); |
---|
864 | bool checkFeasibility = numberObjects>model->numberIntegers(); |
---|
865 | int maximumStrong = CoinMax(CoinMin(model->numberStrong(),numberObjects),1); |
---|
866 | int numberColumns = model->getNumCols(); |
---|
867 | double * saveUpper = new double[numberColumns]; |
---|
868 | double * saveLower = new double[numberColumns]; |
---|
869 | |
---|
870 | // Save solution in case heuristics need good solution later |
---|
871 | |
---|
872 | double * saveSolution = new double[numberColumns]; |
---|
873 | memcpy(saveSolution,solver->getColSolution(),numberColumns*sizeof(double)); |
---|
874 | |
---|
875 | /* |
---|
876 | Get a branching decision object. Use the default decision criteria unless |
---|
877 | the user has loaded a decision method into the model. |
---|
878 | */ |
---|
879 | CbcBranchDecision *decision = model->branchingMethod(); |
---|
880 | if (!decision) |
---|
881 | decision = new CbcBranchDefaultDecision(); |
---|
882 | |
---|
883 | typedef struct { |
---|
884 | CbcBranchingObject * possibleBranch; // what a branch would do |
---|
885 | double upMovement; // cost going up (and initial away from feasible) |
---|
886 | double downMovement; // cost going down |
---|
887 | int numIntInfeasUp ; // without odd ones |
---|
888 | int numObjInfeasUp ; // just odd ones |
---|
889 | bool finishedUp; // true if solver finished |
---|
890 | int numItersUp ; // number of iterations in solver |
---|
891 | int numIntInfeasDown ; // without odd ones |
---|
892 | int numObjInfeasDown ; // just odd ones |
---|
893 | bool finishedDown; // true if solver finished |
---|
894 | int numItersDown; // number of iterations in solver |
---|
895 | int objectNumber; // Which object it is |
---|
896 | } Strong; |
---|
897 | Strong * choice = new Strong[maximumStrong]; |
---|
898 | for (i=0;i<numberColumns;i++) { |
---|
899 | saveLower[i] = lower[i]; |
---|
900 | saveUpper[i] = upper[i]; |
---|
901 | } |
---|
902 | // Some objects may compute an estimate of best solution from here |
---|
903 | double estimatedDegradation=0.0; |
---|
904 | int numberIntegerInfeasibilities=0; // without odd ones |
---|
905 | |
---|
906 | // We may go round this loop twice (only if we think we have solution) |
---|
907 | for (int iPass=0;iPass<2;iPass++) { |
---|
908 | |
---|
909 | // compute current state |
---|
910 | int numberObjectInfeasibilities; // just odd ones |
---|
911 | model->feasibleSolution( |
---|
912 | numberIntegerInfeasibilities, |
---|
913 | numberObjectInfeasibilities); |
---|
914 | // If forcePriority > 0 then we want best solution |
---|
915 | const double * bestSolution = NULL; |
---|
916 | int hotstartStrategy=model->getHotstartStrategy(); |
---|
917 | if (hotstartStrategy>0) { |
---|
918 | bestSolution = model->bestSolution(); |
---|
919 | } |
---|
920 | |
---|
921 | // Some objects may compute an estimate of best solution from here |
---|
922 | estimatedDegradation=0.0; |
---|
923 | numberUnsatisfied_ = 0; |
---|
924 | int bestPriority=INT_MAX; |
---|
925 | /* |
---|
926 | Scan for branching objects that indicate infeasibility. Choose the best |
---|
927 | maximumStrong candidates, using priority as the first criteria, then |
---|
928 | integer infeasibility. |
---|
929 | |
---|
930 | The algorithm is to fill the choice array with a set of good candidates (by |
---|
931 | infeasibility) with priority bestPriority. Finding a candidate with |
---|
932 | priority better (less) than bestPriority flushes the choice array. (This |
---|
933 | serves as initialization when the first candidate is found.) |
---|
934 | |
---|
935 | A new candidate is added to choices only if its infeasibility exceeds the |
---|
936 | current max infeasibility (mostAway). When a candidate is added, it |
---|
937 | replaces the candidate with the smallest infeasibility (tracked by |
---|
938 | iSmallest). |
---|
939 | */ |
---|
940 | int iSmallest = 0; |
---|
941 | double mostAway = 1.0e-100; |
---|
942 | for (i = 0 ; i < maximumStrong ; i++) choice[i].possibleBranch = NULL ; |
---|
943 | numberStrong=0; |
---|
944 | for (i=0;i<numberObjects;i++) { |
---|
945 | const CbcObject * object = model->object(i); |
---|
946 | int preferredWay; |
---|
947 | double infeasibility = object->infeasibility(preferredWay); |
---|
948 | int priorityLevel = model->priority(i); |
---|
949 | if (bestSolution) { |
---|
950 | // we are doing hot start |
---|
951 | const CbcSimpleInteger * thisOne = dynamic_cast <const CbcSimpleInteger *> (object); |
---|
952 | if (thisOne) { |
---|
953 | int iColumn = thisOne->modelSequence(); |
---|
954 | if (saveUpper[iColumn]>saveLower[iColumn]) { |
---|
955 | double value = saveSolution[iColumn]; |
---|
956 | double targetValue = bestSolution[iColumn]; |
---|
957 | //double originalLower = thisOne->originalLower(); |
---|
958 | //double originalUpper = thisOne->originalUpper(); |
---|
959 | // switch off if not possible |
---|
960 | if (targetValue>=saveLower[iColumn]&&targetValue<=saveUpper[iColumn]) { |
---|
961 | /* priority outranks rest always if hotstartStrategy >1 |
---|
962 | otherwise can be downgraded if at correct level. |
---|
963 | Infeasibility may be increased by targetValue to choose 1.0 values first. |
---|
964 | */ |
---|
965 | if (fabs(value-targetValue)>integerTolerance) { |
---|
966 | if (value>targetValue) { |
---|
967 | infeasibility += value; |
---|
968 | preferredWay=-1; |
---|
969 | } else { |
---|
970 | infeasibility += targetValue; |
---|
971 | preferredWay=1; |
---|
972 | } |
---|
973 | } else if (hotstartStrategy>1) { |
---|
974 | if (targetValue==saveLower[iColumn]) { |
---|
975 | infeasibility += integerTolerance+1.0e-12; |
---|
976 | preferredWay=-1; |
---|
977 | } else if (targetValue==saveUpper[iColumn]) { |
---|
978 | infeasibility += integerTolerance+1.0e-12; |
---|
979 | preferredWay=1; |
---|
980 | } else { |
---|
981 | infeasibility += integerTolerance+1.0e-12; |
---|
982 | preferredWay=1; |
---|
983 | } |
---|
984 | } else { |
---|
985 | priorityLevel += 10000000; |
---|
986 | } |
---|
987 | } else { |
---|
988 | // switch off if not possible |
---|
989 | bestSolution=NULL; |
---|
990 | model->setHotstartStrategy(0); |
---|
991 | } |
---|
992 | } |
---|
993 | } |
---|
994 | } |
---|
995 | if (infeasibility) { |
---|
996 | // Increase estimated degradation to solution |
---|
997 | estimatedDegradation += CoinMin(object->upEstimate(),object->downEstimate()); |
---|
998 | numberUnsatisfied_++; |
---|
999 | // Better priority? Flush choices. |
---|
1000 | if (priorityLevel<bestPriority) { |
---|
1001 | int j; |
---|
1002 | iSmallest=0; |
---|
1003 | for (j=0;j<maximumStrong;j++) { |
---|
1004 | choice[j].upMovement=0.0; |
---|
1005 | delete choice[j].possibleBranch; |
---|
1006 | choice[j].possibleBranch=NULL; |
---|
1007 | } |
---|
1008 | bestPriority = priorityLevel; |
---|
1009 | mostAway=1.0e-100; |
---|
1010 | numberStrong=0; |
---|
1011 | } else if (priorityLevel>bestPriority) { |
---|
1012 | continue; |
---|
1013 | } |
---|
1014 | // Check for suitability based on infeasibility. |
---|
1015 | if (infeasibility>mostAway) { |
---|
1016 | //add to list |
---|
1017 | choice[iSmallest].upMovement=infeasibility; |
---|
1018 | delete choice[iSmallest].possibleBranch; |
---|
1019 | choice[iSmallest].possibleBranch=object->createBranch(preferredWay); |
---|
1020 | numberStrong = CoinMax(numberStrong,iSmallest+1); |
---|
1021 | // Save which object it was |
---|
1022 | choice[iSmallest].objectNumber=i; |
---|
1023 | int j; |
---|
1024 | iSmallest=-1; |
---|
1025 | mostAway = 1.0e50; |
---|
1026 | for (j=0;j<maximumStrong;j++) { |
---|
1027 | if (choice[j].upMovement<mostAway) { |
---|
1028 | mostAway=choice[j].upMovement; |
---|
1029 | iSmallest=j; |
---|
1030 | } |
---|
1031 | } |
---|
1032 | } |
---|
1033 | } |
---|
1034 | } |
---|
1035 | if (numberUnsatisfied_) { |
---|
1036 | // some infeasibilities - go to next steps |
---|
1037 | break; |
---|
1038 | } else if (!iPass) { |
---|
1039 | // looks like a solution - get paranoid |
---|
1040 | bool roundAgain=false; |
---|
1041 | // get basis |
---|
1042 | CoinWarmStartBasis * ws = dynamic_cast<CoinWarmStartBasis*>(solver->getWarmStart()); |
---|
1043 | if (!ws) |
---|
1044 | break; |
---|
1045 | for (i=0;i<numberColumns;i++) { |
---|
1046 | double value = saveSolution[i]; |
---|
1047 | if (value<lower[i]) { |
---|
1048 | saveSolution[i]=lower[i]; |
---|
1049 | roundAgain=true; |
---|
1050 | ws->setStructStatus(i,CoinWarmStartBasis::atLowerBound); |
---|
1051 | } else if (value>upper[i]) { |
---|
1052 | saveSolution[i]=upper[i]; |
---|
1053 | roundAgain=true; |
---|
1054 | ws->setStructStatus(i,CoinWarmStartBasis::atUpperBound); |
---|
1055 | } |
---|
1056 | } |
---|
1057 | if (roundAgain) { |
---|
1058 | // restore basis |
---|
1059 | solver->setWarmStart(ws); |
---|
1060 | delete ws; |
---|
1061 | solver->resolve(); |
---|
1062 | memcpy(saveSolution,solver->getColSolution(),numberColumns*sizeof(double)); |
---|
1063 | if (!solver->isProvenOptimal()) { |
---|
1064 | // infeasible |
---|
1065 | anyAction=-2; |
---|
1066 | break; |
---|
1067 | } |
---|
1068 | } else { |
---|
1069 | delete ws; |
---|
1070 | break; |
---|
1071 | } |
---|
1072 | } |
---|
1073 | } |
---|
1074 | /* Some solvers can do the strong branching calculations faster if |
---|
1075 | they do them all at once. At present only Clp does for ordinary |
---|
1076 | integers but I think this coding would be easy to modify |
---|
1077 | */ |
---|
1078 | bool allNormal=true; // to say if we can do fast strong branching |
---|
1079 | // Say which one will be best |
---|
1080 | int bestChoice=0; |
---|
1081 | double worstInfeasibility=0.0; |
---|
1082 | for (i=0;i<numberStrong;i++) { |
---|
1083 | choice[i].numIntInfeasUp = numberUnsatisfied_; |
---|
1084 | choice[i].numIntInfeasDown = numberUnsatisfied_; |
---|
1085 | if (!dynamic_cast <const CbcSimpleInteger *> (model->object(choice[i].objectNumber))) |
---|
1086 | allNormal=false; // Something odd so lets skip clever fast branching |
---|
1087 | if ( !model->object(choice[i].objectNumber)->boundBranch()) |
---|
1088 | numberStrong=0; // switch off |
---|
1089 | // Do best choice in case switched off |
---|
1090 | if (choice[i].upMovement>worstInfeasibility) { |
---|
1091 | worstInfeasibility=choice[i].upMovement; |
---|
1092 | bestChoice=i; |
---|
1093 | } |
---|
1094 | } |
---|
1095 | // If we have hit max time don't do strong branching |
---|
1096 | bool hitMaxTime = ( CoinCpuTime()-model->getDblParam(CbcModel::CbcStartSeconds) > |
---|
1097 | model->getDblParam(CbcModel::CbcMaximumSeconds)); |
---|
1098 | // also give up if we are looping round too much |
---|
1099 | if (hitMaxTime||numberPassesLeft<=0) |
---|
1100 | numberStrong=0; |
---|
1101 | /* |
---|
1102 | Is strong branching enabled? If so, set up and do it. Otherwise, we'll |
---|
1103 | fall through to simple branching. |
---|
1104 | |
---|
1105 | Setup for strong branching involves saving the current basis (for restoration |
---|
1106 | afterwards) and setting up for hot starts. |
---|
1107 | */ |
---|
1108 | if (numberStrong&&model->numberStrong()) { |
---|
1109 | |
---|
1110 | bool solveAll=false; // set true to say look at all even if some fixed (experiment) |
---|
1111 | // worth trying if too many times |
---|
1112 | // Save basis |
---|
1113 | CoinWarmStart * ws = solver->getWarmStart(); |
---|
1114 | // save limit |
---|
1115 | int saveLimit; |
---|
1116 | solver->getIntParam(OsiMaxNumIterationHotStart,saveLimit); |
---|
1117 | if (beforeSolution) |
---|
1118 | solver->setIntParam(OsiMaxNumIterationHotStart,10000); // go to end |
---|
1119 | |
---|
1120 | /* If we are doing all strong branching in one go then we create new arrays |
---|
1121 | to store information. If clp NULL then doing old way. |
---|
1122 | Going down - |
---|
1123 | outputSolution[2*i] is final solution. |
---|
1124 | outputStuff[2*i] is status (0 - finished, 1 infeas, other unknown |
---|
1125 | outputStuff[2*i+numberStrong] is number iterations |
---|
1126 | On entry newUpper[i] is new upper bound, on exit obj change |
---|
1127 | Going up - |
---|
1128 | outputSolution[2*i+1] is final solution. |
---|
1129 | outputStuff[2*i+1] is status (0 - finished, 1 infeas, other unknown |
---|
1130 | outputStuff[2*i+1+numberStrong] is number iterations |
---|
1131 | On entry newLower[i] is new lower bound, on exit obj change |
---|
1132 | */ |
---|
1133 | OsiClpSolverInterface * osiclp = dynamic_cast< OsiClpSolverInterface*> (solver); |
---|
1134 | ClpSimplex * clp=NULL; |
---|
1135 | double * newLower = NULL; |
---|
1136 | double * newUpper = NULL; |
---|
1137 | double ** outputSolution=NULL; |
---|
1138 | int * outputStuff=NULL; |
---|
1139 | int saveLogLevel=0; |
---|
1140 | // Go back to normal way if user wants it |
---|
1141 | if (osiclp&&(osiclp->specialOptions()&16)!=0) |
---|
1142 | allNormal=false; |
---|
1143 | if (osiclp&& allNormal) { |
---|
1144 | clp = osiclp->getModelPtr(); |
---|
1145 | saveLogLevel = clp->logLevel(); |
---|
1146 | clp->setLogLevel(0); |
---|
1147 | int saveOptions = clp->specialOptions(); |
---|
1148 | int startFinishOptions; |
---|
1149 | int specialOptions = osiclp->specialOptions(); |
---|
1150 | if((specialOptions&1)==0) { |
---|
1151 | startFinishOptions=0; |
---|
1152 | clp->setSpecialOptions(saveOptions|(64|1024)); |
---|
1153 | } else { |
---|
1154 | startFinishOptions=1+2+4; |
---|
1155 | //startFinishOptions=1+4; // for moment re-factorize |
---|
1156 | if((specialOptions&4)==0) |
---|
1157 | clp->setSpecialOptions(saveOptions|(64|128|512|1024|4096)); |
---|
1158 | else |
---|
1159 | clp->setSpecialOptions(saveOptions|(64|128|512|1024|2048|4096)); |
---|
1160 | } |
---|
1161 | // User may want to clean up before strong branching |
---|
1162 | if ((clp->specialOptions()&32)!=0) { |
---|
1163 | clp->primal(1); |
---|
1164 | if (clp->numberIterations()) |
---|
1165 | model->messageHandler()->message(CBC_ITERATE_STRONG,model->messages()) |
---|
1166 | << clp->numberIterations() |
---|
1167 | <<CoinMessageEol; |
---|
1168 | } |
---|
1169 | int saveMaxIts = clp->maximumIterations(); |
---|
1170 | clp->setMaximumIterations(saveLimit); |
---|
1171 | // Clp - do a different way |
---|
1172 | newLower = new double[numberStrong]; |
---|
1173 | newUpper = new double[numberStrong]; |
---|
1174 | outputSolution = new double * [2*numberStrong]; |
---|
1175 | outputStuff = new int [4*numberStrong]; |
---|
1176 | int * which = new int[numberStrong]; |
---|
1177 | for (i=0;i<numberStrong;i++) { |
---|
1178 | int iObject = choice[i].objectNumber; |
---|
1179 | const CbcObject * object = model->object(iObject); |
---|
1180 | const CbcSimpleInteger * simple = dynamic_cast <const CbcSimpleInteger *> (object); |
---|
1181 | int iSequence = simple->modelSequence(); |
---|
1182 | newLower[i]= ceil(saveSolution[iSequence]); |
---|
1183 | newUpper[i]= floor(saveSolution[iSequence]); |
---|
1184 | which[i]=iSequence; |
---|
1185 | outputSolution[2*i]= new double [numberColumns]; |
---|
1186 | outputSolution[2*i+1]= new double [numberColumns]; |
---|
1187 | } |
---|
1188 | int returnCode=clp->strongBranching(numberStrong,which, |
---|
1189 | newLower, newUpper,outputSolution, |
---|
1190 | outputStuff,outputStuff+2*numberStrong,!solveAll,false, |
---|
1191 | startFinishOptions); |
---|
1192 | clp->setSpecialOptions(saveOptions); // restore |
---|
1193 | clp->setMaximumIterations(saveMaxIts); |
---|
1194 | if (returnCode==-2) { |
---|
1195 | // bad factorization!!! |
---|
1196 | // Doing normal way |
---|
1197 | // Mark hot start |
---|
1198 | solver->markHotStart(); |
---|
1199 | clp->setLogLevel(saveLogLevel); |
---|
1200 | clp = NULL; |
---|
1201 | } |
---|
1202 | delete [] which; |
---|
1203 | } else { |
---|
1204 | // Doing normal way |
---|
1205 | // Mark hot start |
---|
1206 | solver->markHotStart(); |
---|
1207 | } |
---|
1208 | /* |
---|
1209 | Open a loop to do the strong branching LPs. For each candidate variable, |
---|
1210 | solve an LP with the variable forced down, then up. If a direction turns |
---|
1211 | out to be infeasible or monotonic (i.e., over the dual objective cutoff), |
---|
1212 | force the objective change to be big (1.0e100). If we determine the problem |
---|
1213 | is infeasible, or find a monotone variable, escape the loop. |
---|
1214 | |
---|
1215 | TODO: The `restore bounds' part might be better encapsulated as an |
---|
1216 | unbranch() method. Branching objects more exotic than simple integers |
---|
1217 | or cliques might not restrict themselves to variable bounds. |
---|
1218 | |
---|
1219 | TODO: Virtuous solvers invalidate the current solution (or give bogus |
---|
1220 | results :-) when the bounds are changed out from under them. So we |
---|
1221 | need to do all the work associated with finding a new solution before |
---|
1222 | restoring the bounds. |
---|
1223 | */ |
---|
1224 | for (i = 0 ; i < numberStrong ; i++) |
---|
1225 | { double objectiveChange ; |
---|
1226 | double newObjectiveValue=1.0e100; |
---|
1227 | // status is 0 finished, 1 infeasible and other |
---|
1228 | int iStatus; |
---|
1229 | /* |
---|
1230 | Try the down direction first. (Specify the initial branching alternative as |
---|
1231 | down with a call to way(-1). Each subsequent call to branch() performs the |
---|
1232 | specified branch and advances the branch object state to the next branch |
---|
1233 | alternative.) |
---|
1234 | */ |
---|
1235 | if (!clp) { |
---|
1236 | choice[i].possibleBranch->way(-1) ; |
---|
1237 | choice[i].possibleBranch->branch() ; |
---|
1238 | bool feasible=true; |
---|
1239 | if (checkFeasibility) { |
---|
1240 | // check branching did not make infeasible |
---|
1241 | int iColumn; |
---|
1242 | int numberColumns = solver->getNumCols(); |
---|
1243 | const double * columnLower = solver->getColLower(); |
---|
1244 | const double * columnUpper = solver->getColUpper(); |
---|
1245 | for (iColumn= 0;iColumn<numberColumns;iColumn++) { |
---|
1246 | if (columnLower[iColumn]>columnUpper[iColumn]+1.0e-5) |
---|
1247 | feasible=false; |
---|
1248 | } |
---|
1249 | } |
---|
1250 | if (feasible) { |
---|
1251 | solver->solveFromHotStart() ; |
---|
1252 | /* |
---|
1253 | We now have an estimate of objective degradation that we can use for strong |
---|
1254 | branching. If we're over the cutoff, the variable is monotone up. |
---|
1255 | If we actually made it to optimality, check for a solution, and if we have |
---|
1256 | a good one, call setBestSolution to process it. Note that this may reduce the |
---|
1257 | cutoff, so we check again to see if we can declare this variable monotone. |
---|
1258 | */ |
---|
1259 | if (solver->isProvenOptimal()) |
---|
1260 | iStatus=0; // optimal |
---|
1261 | else if (solver->isIterationLimitReached() |
---|
1262 | &&!solver->isDualObjectiveLimitReached()) |
---|
1263 | iStatus=2; // unknown |
---|
1264 | else |
---|
1265 | iStatus=1; // infeasible |
---|
1266 | newObjectiveValue = solver->getObjSense()*solver->getObjValue(); |
---|
1267 | choice[i].numItersDown = solver->getIterationCount(); |
---|
1268 | } else { |
---|
1269 | iStatus=1; // infeasible |
---|
1270 | newObjectiveValue = 1.0e100; |
---|
1271 | choice[i].numItersDown = 0; |
---|
1272 | } |
---|
1273 | objectiveChange = newObjectiveValue-objectiveValue ; |
---|
1274 | } else { |
---|
1275 | iStatus = outputStuff[2*i]; |
---|
1276 | choice[i].numItersDown = outputStuff[2*numberStrong+2*i]; |
---|
1277 | newObjectiveValue = objectiveValue+newUpper[i]; |
---|
1278 | solver->setColSolution(outputSolution[2*i]); |
---|
1279 | } |
---|
1280 | objectiveChange = newObjectiveValue - objectiveValue; |
---|
1281 | if (!iStatus) { |
---|
1282 | choice[i].finishedDown = true ; |
---|
1283 | if (newObjectiveValue>=model->getCutoff()) { |
---|
1284 | objectiveChange = 1.0e100; // say infeasible |
---|
1285 | } else { |
---|
1286 | // See if integer solution |
---|
1287 | if (model->feasibleSolution(choice[i].numIntInfeasDown, |
---|
1288 | choice[i].numObjInfeasDown)) |
---|
1289 | { model->setBestSolution(CBC_STRONGSOL, |
---|
1290 | newObjectiveValue, |
---|
1291 | solver->getColSolution()) ; |
---|
1292 | if (newObjectiveValue >= model->getCutoff()) // *new* cutoff |
---|
1293 | objectiveChange = 1.0e100 ; |
---|
1294 | } |
---|
1295 | } |
---|
1296 | } else if (iStatus==1) { |
---|
1297 | objectiveChange = 1.0e100 ; |
---|
1298 | } else { |
---|
1299 | // Can't say much as we did not finish |
---|
1300 | choice[i].finishedDown = false ; |
---|
1301 | } |
---|
1302 | choice[i].downMovement = objectiveChange ; |
---|
1303 | |
---|
1304 | // restore bounds |
---|
1305 | if (!clp) |
---|
1306 | { for (int j=0;j<numberColumns;j++) { |
---|
1307 | if (saveLower[j] != lower[j]) |
---|
1308 | solver->setColLower(j,saveLower[j]); |
---|
1309 | if (saveUpper[j] != upper[j]) |
---|
1310 | solver->setColUpper(j,saveUpper[j]); |
---|
1311 | } |
---|
1312 | } |
---|
1313 | //printf("Down on %d, status is %d, obj %g its %d cost %g finished %d inf %d infobj %d\n", |
---|
1314 | // choice[i].objectNumber,iStatus,newObjectiveValue,choice[i].numItersDown, |
---|
1315 | // choice[i].downMovement,choice[i].finishedDown,choice[i].numIntInfeasDown, |
---|
1316 | // choice[i].numObjInfeasDown); |
---|
1317 | |
---|
1318 | // repeat the whole exercise, forcing the variable up |
---|
1319 | if (!clp) { |
---|
1320 | choice[i].possibleBranch->branch(); |
---|
1321 | bool feasible=true; |
---|
1322 | if (checkFeasibility) { |
---|
1323 | // check branching did not make infeasible |
---|
1324 | int iColumn; |
---|
1325 | int numberColumns = solver->getNumCols(); |
---|
1326 | const double * columnLower = solver->getColLower(); |
---|
1327 | const double * columnUpper = solver->getColUpper(); |
---|
1328 | for (iColumn= 0;iColumn<numberColumns;iColumn++) { |
---|
1329 | if (columnLower[iColumn]>columnUpper[iColumn]+1.0e-5) |
---|
1330 | feasible=false; |
---|
1331 | } |
---|
1332 | } |
---|
1333 | if (feasible) { |
---|
1334 | solver->solveFromHotStart() ; |
---|
1335 | /* |
---|
1336 | We now have an estimate of objective degradation that we can use for strong |
---|
1337 | branching. If we're over the cutoff, the variable is monotone up. |
---|
1338 | If we actually made it to optimality, check for a solution, and if we have |
---|
1339 | a good one, call setBestSolution to process it. Note that this may reduce the |
---|
1340 | cutoff, so we check again to see if we can declare this variable monotone. |
---|
1341 | */ |
---|
1342 | if (solver->isProvenOptimal()) |
---|
1343 | iStatus=0; // optimal |
---|
1344 | else if (solver->isIterationLimitReached() |
---|
1345 | &&!solver->isDualObjectiveLimitReached()) |
---|
1346 | iStatus=2; // unknown |
---|
1347 | else |
---|
1348 | iStatus=1; // infeasible |
---|
1349 | newObjectiveValue = solver->getObjSense()*solver->getObjValue(); |
---|
1350 | choice[i].numItersUp = solver->getIterationCount(); |
---|
1351 | } else { |
---|
1352 | iStatus=1; // infeasible |
---|
1353 | newObjectiveValue = 1.0e100; |
---|
1354 | choice[i].numItersDown = 0; |
---|
1355 | } |
---|
1356 | objectiveChange = newObjectiveValue-objectiveValue ; |
---|
1357 | } else { |
---|
1358 | iStatus = outputStuff[2*i+1]; |
---|
1359 | choice[i].numItersUp = outputStuff[2*numberStrong+2*i+1]; |
---|
1360 | newObjectiveValue = objectiveValue+newLower[i]; |
---|
1361 | solver->setColSolution(outputSolution[2*i+1]); |
---|
1362 | } |
---|
1363 | objectiveChange = newObjectiveValue - objectiveValue; |
---|
1364 | if (!iStatus) { |
---|
1365 | choice[i].finishedUp = true ; |
---|
1366 | if (newObjectiveValue>=model->getCutoff()) { |
---|
1367 | objectiveChange = 1.0e100; // say infeasible |
---|
1368 | } else { |
---|
1369 | // See if integer solution |
---|
1370 | if (model->feasibleSolution(choice[i].numIntInfeasUp, |
---|
1371 | choice[i].numObjInfeasUp)) |
---|
1372 | { model->setBestSolution(CBC_STRONGSOL, |
---|
1373 | newObjectiveValue, |
---|
1374 | solver->getColSolution()) ; |
---|
1375 | if (newObjectiveValue >= model->getCutoff()) // *new* cutoff |
---|
1376 | objectiveChange = 1.0e100 ; |
---|
1377 | } |
---|
1378 | } |
---|
1379 | } else if (iStatus==1) { |
---|
1380 | objectiveChange = 1.0e100 ; |
---|
1381 | } else { |
---|
1382 | // Can't say much as we did not finish |
---|
1383 | choice[i].finishedUp = false ; |
---|
1384 | } |
---|
1385 | choice[i].upMovement = objectiveChange ; |
---|
1386 | |
---|
1387 | // restore bounds |
---|
1388 | if (!clp) |
---|
1389 | { for (int j=0;j<numberColumns;j++) { |
---|
1390 | if (saveLower[j] != lower[j]) |
---|
1391 | solver->setColLower(j,saveLower[j]); |
---|
1392 | if (saveUpper[j] != upper[j]) |
---|
1393 | solver->setColUpper(j,saveUpper[j]); |
---|
1394 | } |
---|
1395 | } |
---|
1396 | |
---|
1397 | //printf("Up on %d, status is %d, obj %g its %d cost %g finished %d inf %d infobj %d\n", |
---|
1398 | // choice[i].objectNumber,iStatus,newObjectiveValue,choice[i].numItersUp, |
---|
1399 | // choice[i].upMovement,choice[i].finishedUp,choice[i].numIntInfeasUp, |
---|
1400 | // choice[i].numObjInfeasUp); |
---|
1401 | |
---|
1402 | /* |
---|
1403 | End of evaluation for this candidate variable. Possibilities are: |
---|
1404 | * Both sides below cutoff; this variable is a candidate for branching. |
---|
1405 | * Both sides infeasible or above the objective cutoff: no further action |
---|
1406 | here. Break from the evaluation loop and assume the node will be purged |
---|
1407 | by the caller. |
---|
1408 | * One side below cutoff: Install the branch (i.e., fix the variable). Break |
---|
1409 | from the evaluation loop and assume the node will be reoptimised by the |
---|
1410 | caller. |
---|
1411 | */ |
---|
1412 | if (choice[i].upMovement<1.0e100) { |
---|
1413 | if(choice[i].downMovement<1.0e100) { |
---|
1414 | // feasible - no action |
---|
1415 | } else { |
---|
1416 | // up feasible, down infeasible |
---|
1417 | anyAction=-1; |
---|
1418 | if (!solveAll) { |
---|
1419 | choice[i].possibleBranch->way(1); |
---|
1420 | choice[i].possibleBranch->branch(); |
---|
1421 | break; |
---|
1422 | } |
---|
1423 | } |
---|
1424 | } else { |
---|
1425 | if(choice[i].downMovement<1.0e100) { |
---|
1426 | // down feasible, up infeasible |
---|
1427 | anyAction=-1; |
---|
1428 | if (!solveAll) { |
---|
1429 | choice[i].possibleBranch->way(-1); |
---|
1430 | choice[i].possibleBranch->branch(); |
---|
1431 | break; |
---|
1432 | } |
---|
1433 | } else { |
---|
1434 | // neither side feasible |
---|
1435 | anyAction=-2; |
---|
1436 | break; |
---|
1437 | } |
---|
1438 | } |
---|
1439 | } |
---|
1440 | |
---|
1441 | if (!clp) { |
---|
1442 | // Delete the snapshot |
---|
1443 | solver->unmarkHotStart(); |
---|
1444 | } else { |
---|
1445 | clp->setLogLevel(saveLogLevel); |
---|
1446 | delete [] newLower; |
---|
1447 | delete [] newUpper; |
---|
1448 | delete [] outputStuff; |
---|
1449 | int i; |
---|
1450 | for (i=0;i<2*numberStrong;i++) |
---|
1451 | delete [] outputSolution[i]; |
---|
1452 | delete [] outputSolution; |
---|
1453 | } |
---|
1454 | solver->setIntParam(OsiMaxNumIterationHotStart,saveLimit); |
---|
1455 | // restore basis |
---|
1456 | solver->setWarmStart(ws); |
---|
1457 | delete ws; |
---|
1458 | |
---|
1459 | /* |
---|
1460 | anyAction >= 0 indicates that strong branching didn't produce any monotone |
---|
1461 | variables. Sift through the candidates for the best one. |
---|
1462 | |
---|
1463 | QUERY: Setting numberNodes looks to be a distributed noop. numberNodes is |
---|
1464 | local to this code block. Perhaps should be numberNodes_ from model? |
---|
1465 | Unclear what this calculation is doing. |
---|
1466 | */ |
---|
1467 | if (anyAction>=0) { |
---|
1468 | |
---|
1469 | int numberNodes = model->getNodeCount(); |
---|
1470 | // get average cost per iteration and assume stopped ones |
---|
1471 | // would stop after 50% more iterations at average cost??? !!! ??? |
---|
1472 | double averageCostPerIteration=0.0; |
---|
1473 | double totalNumberIterations=1.0; |
---|
1474 | int smallestNumberInfeasibilities=INT_MAX; |
---|
1475 | for (i=0;i<numberStrong;i++) { |
---|
1476 | totalNumberIterations += choice[i].numItersDown + |
---|
1477 | choice[i].numItersUp ; |
---|
1478 | averageCostPerIteration += choice[i].downMovement + |
---|
1479 | choice[i].upMovement; |
---|
1480 | smallestNumberInfeasibilities= |
---|
1481 | CoinMin(CoinMin(choice[i].numIntInfeasDown , |
---|
1482 | choice[i].numIntInfeasUp ), |
---|
1483 | smallestNumberInfeasibilities); |
---|
1484 | } |
---|
1485 | if (smallestNumberInfeasibilities>=numberIntegerInfeasibilities) |
---|
1486 | numberNodes=1000000; // switch off search for better solution |
---|
1487 | numberNodes=1000000; // switch off anyway |
---|
1488 | averageCostPerIteration /= totalNumberIterations; |
---|
1489 | // all feasible - choose best bet |
---|
1490 | |
---|
1491 | #if 0 |
---|
1492 | for (i = 0 ; i < numberStrong ; i++) |
---|
1493 | { int iColumn = |
---|
1494 | model->integerVariable()[choice[i].possibleBranch->variable()] ; |
---|
1495 | model->messageHandler()->message(CBC_STRONG,model->messages()) |
---|
1496 | << i << iColumn |
---|
1497 | <<choice[i].downMovement<<choice[i].numIntInfeasDown |
---|
1498 | <<choice[i].upMovement<<choice[i].numIntInfeasUp |
---|
1499 | <<choice[i].possibleBranch->value() |
---|
1500 | <<CoinMessageEol; |
---|
1501 | int betterWay = decision->betterBranch(choice[i].possibleBranch, |
---|
1502 | branch_, |
---|
1503 | choice[i].upMovement, |
---|
1504 | choice[i].numIntInfeasUp , |
---|
1505 | choice[i].downMovement, |
---|
1506 | choice[i].numIntInfeasDown ); |
---|
1507 | if (betterWay) { |
---|
1508 | delete branch_; |
---|
1509 | // C) create branching object |
---|
1510 | branch_ = choice[i].possibleBranch; |
---|
1511 | choice[i].possibleBranch=NULL; |
---|
1512 | branch_->way(betterWay); |
---|
1513 | } |
---|
1514 | } |
---|
1515 | #else |
---|
1516 | // New method does all at once so it can be more sophisticated |
---|
1517 | // in deciding how to balance actions. |
---|
1518 | // But it does need arrays |
---|
1519 | double * changeUp = new double [numberStrong]; |
---|
1520 | int * numberInfeasibilitiesUp = new int [numberStrong]; |
---|
1521 | double * changeDown = new double [numberStrong]; |
---|
1522 | int * numberInfeasibilitiesDown = new int [numberStrong]; |
---|
1523 | CbcBranchingObject ** objects = new CbcBranchingObject * [ numberStrong]; |
---|
1524 | for (i = 0 ; i < numberStrong ; i++) { |
---|
1525 | int iColumn = choice[i].possibleBranch->variable() ; |
---|
1526 | model->messageHandler()->message(CBC_STRONG,model->messages()) |
---|
1527 | << i << iColumn |
---|
1528 | <<choice[i].downMovement<<choice[i].numIntInfeasDown |
---|
1529 | <<choice[i].upMovement<<choice[i].numIntInfeasUp |
---|
1530 | <<choice[i].possibleBranch->value() |
---|
1531 | <<CoinMessageEol; |
---|
1532 | changeUp[i]=choice[i].upMovement; |
---|
1533 | numberInfeasibilitiesUp[i] = choice[i].numIntInfeasUp; |
---|
1534 | changeDown[i]=choice[i].downMovement; |
---|
1535 | numberInfeasibilitiesDown[i] = choice[i].numIntInfeasDown; |
---|
1536 | objects[i] = choice[i].possibleBranch; |
---|
1537 | } |
---|
1538 | int whichObject = decision->bestBranch(objects,numberStrong,numberUnsatisfied_, |
---|
1539 | changeUp,numberInfeasibilitiesUp, |
---|
1540 | changeDown,numberInfeasibilitiesDown, |
---|
1541 | objectiveValue); |
---|
1542 | // move branching object and make sure it will not be deleted |
---|
1543 | if (whichObject>=0) { |
---|
1544 | branch_ = objects[whichObject]; |
---|
1545 | choice[whichObject].possibleBranch=NULL; |
---|
1546 | } |
---|
1547 | delete [] changeUp; |
---|
1548 | delete [] numberInfeasibilitiesUp; |
---|
1549 | delete [] changeDown; |
---|
1550 | delete [] numberInfeasibilitiesDown; |
---|
1551 | delete [] objects; |
---|
1552 | #endif |
---|
1553 | } |
---|
1554 | } |
---|
1555 | /* |
---|
1556 | Simple branching. Probably just one, but we may have got here |
---|
1557 | because of an odd branch e.g. a cut |
---|
1558 | */ |
---|
1559 | else { |
---|
1560 | // not strong |
---|
1561 | // C) create branching object |
---|
1562 | branch_ = choice[bestChoice].possibleBranch; |
---|
1563 | choice[bestChoice].possibleBranch=NULL; |
---|
1564 | } |
---|
1565 | // Set guessed solution value |
---|
1566 | guessedObjectiveValue_ = objectiveValue_+estimatedDegradation; |
---|
1567 | /* |
---|
1568 | Cleanup, then we're outta here. |
---|
1569 | */ |
---|
1570 | if (!model->branchingMethod()) |
---|
1571 | delete decision; |
---|
1572 | |
---|
1573 | for (i=0;i<maximumStrong;i++) |
---|
1574 | delete choice[i].possibleBranch; |
---|
1575 | delete [] choice; |
---|
1576 | delete [] saveLower; |
---|
1577 | delete [] saveUpper; |
---|
1578 | |
---|
1579 | // restore solution |
---|
1580 | solver->setColSolution(saveSolution); |
---|
1581 | delete [] saveSolution; |
---|
1582 | return anyAction; |
---|
1583 | } |
---|
1584 | |
---|
1585 | |
---|
1586 | CbcNode::CbcNode(const CbcNode & rhs) |
---|
1587 | { |
---|
1588 | #ifdef CHECK_NODE |
---|
1589 | printf("CbcNode %x Constructor from rhs %x\n",this,&rhs); |
---|
1590 | #endif |
---|
1591 | if (rhs.nodeInfo_) |
---|
1592 | nodeInfo_ = rhs.nodeInfo_->clone(); |
---|
1593 | else |
---|
1594 | nodeInfo_=NULL; |
---|
1595 | objectiveValue_=rhs.objectiveValue_; |
---|
1596 | guessedObjectiveValue_ = rhs.guessedObjectiveValue_; |
---|
1597 | if (rhs.branch_) |
---|
1598 | branch_=rhs.branch_->clone(); |
---|
1599 | else |
---|
1600 | branch_=NULL; |
---|
1601 | depth_ = rhs.depth_; |
---|
1602 | numberUnsatisfied_ = rhs.numberUnsatisfied_; |
---|
1603 | nodeNumber_=rhs.nodeNumber_; |
---|
1604 | } |
---|
1605 | |
---|
1606 | CbcNode & |
---|
1607 | CbcNode::operator=(const CbcNode & rhs) |
---|
1608 | { |
---|
1609 | if (this != &rhs) { |
---|
1610 | delete nodeInfo_; |
---|
1611 | if (nodeInfo_) |
---|
1612 | nodeInfo_ = rhs.nodeInfo_->clone(); |
---|
1613 | else |
---|
1614 | nodeInfo_ = NULL; |
---|
1615 | objectiveValue_=rhs.objectiveValue_; |
---|
1616 | guessedObjectiveValue_ = rhs.guessedObjectiveValue_; |
---|
1617 | if (rhs.branch_) |
---|
1618 | branch_=rhs.branch_->clone(); |
---|
1619 | else |
---|
1620 | branch_=NULL, |
---|
1621 | depth_ = rhs.depth_; |
---|
1622 | numberUnsatisfied_ = rhs.numberUnsatisfied_; |
---|
1623 | nodeNumber_=rhs.nodeNumber_; |
---|
1624 | } |
---|
1625 | return *this; |
---|
1626 | } |
---|
1627 | |
---|
1628 | |
---|
1629 | CbcNode::~CbcNode () |
---|
1630 | { |
---|
1631 | #ifdef CHECK_NODE |
---|
1632 | if (nodeInfo_) |
---|
1633 | printf("CbcNode %x Destructor nodeInfo %x (%d)\n", |
---|
1634 | this,nodeInfo_,nodeInfo_->numberPointingToThis()); |
---|
1635 | else |
---|
1636 | printf("CbcNode %x Destructor nodeInfo %x (?)\n", |
---|
1637 | this,nodeInfo_); |
---|
1638 | #endif |
---|
1639 | if (nodeInfo_) { |
---|
1640 | nodeInfo_->nullOwner(); |
---|
1641 | int numberToDelete=nodeInfo_->numberBranchesLeft(); |
---|
1642 | // CbcNodeInfo * parent = nodeInfo_->parent(); |
---|
1643 | if (nodeInfo_->decrement(numberToDelete)==0) { |
---|
1644 | delete nodeInfo_; |
---|
1645 | } else { |
---|
1646 | //printf("node %x nodeinfo %x parent %x\n",this,nodeInfo_,parent); |
---|
1647 | // anyway decrement parent |
---|
1648 | //if (parent) |
---|
1649 | ///parent->decrement(1); |
---|
1650 | } |
---|
1651 | } |
---|
1652 | delete branch_; |
---|
1653 | } |
---|
1654 | // Decrement active cut counts |
---|
1655 | void |
---|
1656 | CbcNode::decrementCuts(int change) |
---|
1657 | { |
---|
1658 | if(nodeInfo_) { |
---|
1659 | nodeInfo_->decrementCuts(change); |
---|
1660 | } |
---|
1661 | } |
---|
1662 | void |
---|
1663 | CbcNode::decrementParentCuts(int change) |
---|
1664 | { |
---|
1665 | if(nodeInfo_) { |
---|
1666 | nodeInfo_->decrementParentCuts(change); |
---|
1667 | } |
---|
1668 | } |
---|
1669 | |
---|
1670 | /* |
---|
1671 | Initialize reference counts (numberPointingToThis, numberBranchesLeft_) |
---|
1672 | in the attached nodeInfo_. |
---|
1673 | */ |
---|
1674 | void |
---|
1675 | CbcNode::initializeInfo() |
---|
1676 | { |
---|
1677 | assert(nodeInfo_ && branch_) ; |
---|
1678 | nodeInfo_->initializeInfo(branch_->numberBranches()); |
---|
1679 | } |
---|
1680 | // Nulls out node info |
---|
1681 | void |
---|
1682 | CbcNode::nullNodeInfo() |
---|
1683 | { |
---|
1684 | nodeInfo_=NULL; |
---|
1685 | } |
---|
1686 | |
---|
1687 | int |
---|
1688 | CbcNode::branch() |
---|
1689 | { |
---|
1690 | double changeInGuessed=branch_->branch(true); |
---|
1691 | guessedObjectiveValue_+= changeInGuessed; |
---|
1692 | return nodeInfo_->branchedOn(); |
---|
1693 | } |
---|