1 | /* $Id: ClpDynamicExampleMatrix.cpp 1753 2011-06-19 16:27:26Z stefan $ */ |
2 | // Copyright (C) 2004, International Business Machines |
3 | // Corporation and others. All Rights Reserved. |
4 | // This code is licensed under the terms of the Eclipse Public License (EPL). |
5 | |
6 | #include <cstdio> |
7 | |
8 | #include "CoinPragma.hpp" |
9 | #include "CoinIndexedVector.hpp" |
10 | #include "CoinHelperFunctions.hpp" |
11 | |
12 | #include "ClpSimplex.hpp" |
13 | #include "ClpFactorization.hpp" |
14 | #include "ClpQuadraticObjective.hpp" |
15 | #include "ClpNonLinearCost.hpp" |
16 | // at end to get min/max! |
17 | #include "ClpDynamicExampleMatrix.hpp" |
18 | #include "ClpMessage.hpp" |
19 | //#define CLP_DEBUG |
20 | //#define CLP_DEBUG_PRINT |
21 | //############################################################################# |
22 | // Constructors / Destructor / Assignment |
23 | //############################################################################# |
24 | |
25 | //------------------------------------------------------------------- |
26 | // Default Constructor |
27 | //------------------------------------------------------------------- |
28 | ClpDynamicExampleMatrix::ClpDynamicExampleMatrix () |
29 | : ClpDynamicMatrix(), |
30 | numberColumns_(0), |
31 | startColumnGen_(NULL), |
32 | rowGen_(NULL), |
33 | elementGen_(NULL), |
34 | costGen_(NULL), |
35 | fullStartGen_(NULL), |
36 | dynamicStatusGen_(NULL), |
37 | idGen_(NULL), |
38 | columnLowerGen_(NULL), |
39 | columnUpperGen_(NULL) |
40 | { |
41 | setType(25); |
42 | } |
43 | |
44 | //------------------------------------------------------------------- |
45 | // Copy constructor |
46 | //------------------------------------------------------------------- |
47 | ClpDynamicExampleMatrix::ClpDynamicExampleMatrix (const ClpDynamicExampleMatrix & rhs) |
48 | : ClpDynamicMatrix(rhs) |
49 | { |
50 | numberColumns_ = rhs.numberColumns_; |
51 | startColumnGen_ = ClpCopyOfArray(rhs.startColumnGen_, numberColumns_ + 1); |
52 | CoinBigIndex numberElements = startColumnGen_[numberColumns_]; |
53 | rowGen_ = ClpCopyOfArray(rhs.rowGen_, numberElements); |
54 | elementGen_ = ClpCopyOfArray(rhs.elementGen_, numberElements); |
55 | costGen_ = ClpCopyOfArray(rhs.costGen_, numberColumns_); |
56 | fullStartGen_ = ClpCopyOfArray(rhs.fullStartGen_, numberSets_ + 1); |
57 | dynamicStatusGen_ = ClpCopyOfArray(rhs.dynamicStatusGen_, numberColumns_); |
58 | idGen_ = ClpCopyOfArray(rhs.idGen_, maximumGubColumns_); |
59 | columnLowerGen_ = ClpCopyOfArray(rhs.columnLowerGen_, numberColumns_); |
60 | columnUpperGen_ = ClpCopyOfArray(rhs.columnUpperGen_, numberColumns_); |
61 | } |
62 | |
63 | /* This is the real constructor*/ |
64 | ClpDynamicExampleMatrix::ClpDynamicExampleMatrix(ClpSimplex * model, int numberSets, |
65 | int numberGubColumns, const int * starts, |
66 | const double * lower, const double * upper, |
67 | const CoinBigIndex * startColumn, const int * row, |
68 | const double * element, const double * cost, |
69 | const double * columnLower, const double * columnUpper, |
70 | const unsigned char * status, |
71 | const unsigned char * dynamicStatus, |
72 | int numberIds, const int *ids) |
73 | : ClpDynamicMatrix(model, numberSets, 0, NULL, lower, upper, NULL, NULL, NULL, NULL, NULL, NULL, |
74 | NULL, NULL) |
75 | { |
76 | setType(25); |
77 | numberColumns_ = numberGubColumns; |
78 | // start with safe values - then experiment |
79 | maximumGubColumns_ = numberColumns_; |
80 | maximumElements_ = startColumn[numberColumns_]; |
81 | // delete odd stuff created by ClpDynamicMatrix constructor |
82 | delete [] startSet_; |
83 | startSet_ = new int [numberSets_]; |
84 | delete [] next_; |
85 | next_ = new int [maximumGubColumns_]; |
86 | delete [] row_; |
87 | delete [] element_; |
88 | delete [] startColumn_; |
89 | delete [] cost_; |
90 | delete [] columnLower_; |
91 | delete [] columnUpper_; |
92 | delete [] dynamicStatus_; |
93 | delete [] status_; |
94 | delete [] id_; |
95 | // and size correctly |
96 | row_ = new int [maximumElements_]; |
97 | element_ = new double [maximumElements_]; |
98 | startColumn_ = new CoinBigIndex [maximumGubColumns_+1]; |
99 | // say no columns yet |
100 | numberGubColumns_ = 0; |
101 | startColumn_[0] = 0; |
102 | cost_ = new double[maximumGubColumns_]; |
103 | dynamicStatus_ = new unsigned char [maximumGubColumns_]; |
104 | memset(dynamicStatus_, 0, maximumGubColumns_); |
105 | id_ = new int[maximumGubColumns_]; |
106 | if (columnLower) |
107 | columnLower_ = new double[maximumGubColumns_]; |
108 | else |
109 | columnLower_ = NULL; |
110 | if (columnUpper) |
111 | columnUpper_ = new double[maximumGubColumns_]; |
112 | else |
113 | columnUpper_ = NULL; |
114 | // space for ids |
115 | idGen_ = new int [maximumGubColumns_]; |
116 | int iSet; |
117 | for (iSet = 0; iSet < numberSets_; iSet++) |
118 | startSet_[iSet] = -1; |
119 | // This starts code specific to this storage method |
120 | CoinBigIndex i; |
121 | fullStartGen_ = ClpCopyOfArray(starts, numberSets_ + 1); |
122 | startColumnGen_ = ClpCopyOfArray(startColumn, numberColumns_ + 1); |
123 | CoinBigIndex numberElements = startColumnGen_[numberColumns_]; |
124 | rowGen_ = ClpCopyOfArray(row, numberElements); |
125 | elementGen_ = new double[numberElements]; |
126 | for (i = 0; i < numberElements; i++) |
127 | elementGen_[i] = element[i]; |
128 | costGen_ = new double[numberColumns_]; |
129 | for (i = 0; i < numberColumns_; i++) { |
130 | costGen_[i] = cost[i]; |
131 | // I don't think I need sorted but ... |
132 | CoinSort_2(rowGen_ + startColumnGen_[i], rowGen_ + startColumnGen_[i+1], elementGen_ + startColumnGen_[i]); |
133 | } |
134 | if (columnLower) { |
135 | columnLowerGen_ = new double[numberColumns_]; |
136 | for (i = 0; i < numberColumns_; i++) { |
137 | columnLowerGen_[i] = columnLower[i]; |
138 | if (columnLowerGen_[i]) { |
139 | printf("Non-zero lower bounds not allowed - subtract from model\n" ); |
140 | abort(); |
141 | } |
142 | } |
143 | } else { |
144 | columnLowerGen_ = NULL; |
145 | } |
146 | if (columnUpper) { |
147 | columnUpperGen_ = new double[numberColumns_]; |
148 | for (i = 0; i < numberColumns_; i++) |
149 | columnUpperGen_[i] = columnUpper[i]; |
150 | } else { |
151 | columnUpperGen_ = NULL; |
152 | } |
153 | // end specific coding |
154 | if (columnUpper_) { |
155 | // set all upper bounds so we have enough space |
156 | double * columnUpper = model->columnUpper(); |
157 | for(i = firstDynamic_; i < lastDynamic_; i++) |
158 | columnUpper[i] = 1.0e10; |
159 | } |
160 | if (status) { |
161 | status_ = ClpCopyOfArray(status, numberSets_); |
162 | assert (dynamicStatus); |
163 | CoinMemcpyN(dynamicStatus, numberIds, dynamicStatus_); |
164 | assert (numberIds); |
165 | } else { |
166 | assert (!numberIds); |
167 | status_ = new unsigned char [numberSets_]; |
168 | memset(status_, 0, numberSets_); |
169 | for (i = 0; i < numberSets_; i++) { |
170 | // make slack key |
171 | setStatus(i, ClpSimplex::basic); |
172 | } |
173 | } |
174 | dynamicStatusGen_ = new unsigned char [numberColumns_]; |
175 | memset(dynamicStatusGen_, 0, numberColumns_); // for clarity |
176 | for (i = 0; i < numberColumns_; i++) |
177 | setDynamicStatusGen(i, atLowerBound); |
178 | // Populate with enough columns |
179 | if (!numberIds) { |
180 | // This could be made more sophisticated |
181 | for (iSet = 0; iSet < numberSets_; iSet++) { |
182 | int sequence = fullStartGen_[iSet]; |
183 | CoinBigIndex start = startColumnGen_[sequence]; |
184 | addColumn(startColumnGen_[sequence+1] - start, |
185 | rowGen_ + start, |
186 | elementGen_ + start, |
187 | costGen_[sequence], |
188 | columnLowerGen_ ? columnLowerGen_[sequence] : 0, |
189 | columnUpperGen_ ? columnUpperGen_[sequence] : 1.0e30, |
190 | iSet, getDynamicStatusGen(sequence)); |
191 | idGen_[iSet] = sequence; // say which one in |
192 | setDynamicStatusGen(sequence, inSmall); |
193 | } |
194 | } else { |
195 | // put back old ones |
196 | int * set = new int[numberColumns_]; |
197 | for (iSet = 0; iSet < numberSets_; iSet++) { |
198 | for (CoinBigIndex j = fullStartGen_[iSet]; j < fullStartGen_[iSet+1]; j++) |
199 | set[j] = iSet; |
200 | } |
201 | for (int i = 0; i < numberIds; i++) { |
202 | int sequence = ids[i]; |
203 | CoinBigIndex start = startColumnGen_[sequence]; |
204 | addColumn(startColumnGen_[sequence+1] - start, |
205 | rowGen_ + start, |
206 | elementGen_ + start, |
207 | costGen_[sequence], |
208 | columnLowerGen_ ? columnLowerGen_[sequence] : 0, |
209 | columnUpperGen_ ? columnUpperGen_[sequence] : 1.0e30, |
210 | set[sequence], getDynamicStatus(i)); |
211 | idGen_[iSet] = sequence; // say which one in |
212 | setDynamicStatusGen(sequence, inSmall); |
213 | } |
214 | delete [] set; |
215 | } |
216 | if (!status) { |
217 | gubCrash(); |
218 | } else { |
219 | initialProblem(); |
220 | } |
221 | } |
222 | // This constructor just takes over ownership |
223 | ClpDynamicExampleMatrix::ClpDynamicExampleMatrix(ClpSimplex * model, int numberSets, |
224 | int numberGubColumns, int * starts, |
225 | const double * lower, const double * upper, |
226 | int * startColumn, int * row, |
227 | double * element, double * cost, |
228 | double * columnLower, double * columnUpper, |
229 | const unsigned char * status, |
230 | const unsigned char * dynamicStatus, |
231 | int numberIds, const int *ids) |
232 | : ClpDynamicMatrix(model, numberSets, 0, NULL, lower, upper, NULL, NULL, NULL, NULL, NULL, NULL, |
233 | NULL, NULL) |
234 | { |
235 | setType(25); |
236 | numberColumns_ = numberGubColumns; |
237 | // start with safe values - then experiment |
238 | maximumGubColumns_ = numberColumns_; |
239 | maximumElements_ = startColumn[numberColumns_]; |
240 | // delete odd stuff created by ClpDynamicMatrix constructor |
241 | delete [] startSet_; |
242 | startSet_ = new int [numberSets_]; |
243 | delete [] next_; |
244 | next_ = new int [maximumGubColumns_]; |
245 | delete [] row_; |
246 | delete [] element_; |
247 | delete [] startColumn_; |
248 | delete [] cost_; |
249 | delete [] columnLower_; |
250 | delete [] columnUpper_; |
251 | delete [] dynamicStatus_; |
252 | delete [] status_; |
253 | delete [] id_; |
254 | // and size correctly |
255 | row_ = new int [maximumElements_]; |
256 | element_ = new double [maximumElements_]; |
257 | startColumn_ = new CoinBigIndex [maximumGubColumns_+1]; |
258 | // say no columns yet |
259 | numberGubColumns_ = 0; |
260 | startColumn_[0] = 0; |
261 | cost_ = new double[maximumGubColumns_]; |
262 | dynamicStatus_ = new unsigned char [maximumGubColumns_]; |
263 | memset(dynamicStatus_, 0, maximumGubColumns_); |
264 | id_ = new int[maximumGubColumns_]; |
265 | if (columnLower) |
266 | columnLower_ = new double[maximumGubColumns_]; |
267 | else |
268 | columnLower_ = NULL; |
269 | if (columnUpper) |
270 | columnUpper_ = new double[maximumGubColumns_]; |
271 | else |
272 | columnUpper_ = NULL; |
273 | // space for ids |
274 | idGen_ = new int [maximumGubColumns_]; |
275 | int iSet; |
276 | for (iSet = 0; iSet < numberSets_; iSet++) |
277 | startSet_[iSet] = -1; |
278 | // This starts code specific to this storage method |
279 | CoinBigIndex i; |
280 | fullStartGen_ = starts; |
281 | startColumnGen_ = startColumn; |
282 | rowGen_ = row; |
283 | elementGen_ = element; |
284 | costGen_ = cost; |
285 | for (i = 0; i < numberColumns_; i++) { |
286 | // I don't think I need sorted but ... |
287 | CoinSort_2(rowGen_ + startColumnGen_[i], rowGen_ + startColumnGen_[i+1], elementGen_ + startColumnGen_[i]); |
288 | } |
289 | if (columnLower) { |
290 | columnLowerGen_ = columnLower; |
291 | for (i = 0; i < numberColumns_; i++) { |
292 | if (columnLowerGen_[i]) { |
293 | printf("Non-zero lower bounds not allowed - subtract from model\n" ); |
294 | abort(); |
295 | } |
296 | } |
297 | } else { |
298 | columnLowerGen_ = NULL; |
299 | } |
300 | if (columnUpper) { |
301 | columnUpperGen_ = columnUpper; |
302 | } else { |
303 | columnUpperGen_ = NULL; |
304 | } |
305 | // end specific coding |
306 | if (columnUpper_) { |
307 | // set all upper bounds so we have enough space |
308 | double * columnUpper = model->columnUpper(); |
309 | for(i = firstDynamic_; i < lastDynamic_; i++) |
310 | columnUpper[i] = 1.0e10; |
311 | } |
312 | if (status) { |
313 | status_ = ClpCopyOfArray(status, numberSets_); |
314 | assert (dynamicStatus); |
315 | CoinMemcpyN(dynamicStatus, numberIds, dynamicStatus_); |
316 | assert (numberIds); |
317 | } else { |
318 | assert (!numberIds); |
319 | status_ = new unsigned char [numberSets_]; |
320 | memset(status_, 0, numberSets_); |
321 | for (i = 0; i < numberSets_; i++) { |
322 | // make slack key |
323 | setStatus(i, ClpSimplex::basic); |
324 | } |
325 | } |
326 | dynamicStatusGen_ = new unsigned char [numberColumns_]; |
327 | memset(dynamicStatusGen_, 0, numberColumns_); // for clarity |
328 | for (i = 0; i < numberColumns_; i++) |
329 | setDynamicStatusGen(i, atLowerBound); |
330 | // Populate with enough columns |
331 | if (!numberIds) { |
332 | // This could be made more sophisticated |
333 | for (iSet = 0; iSet < numberSets_; iSet++) { |
334 | int sequence = fullStartGen_[iSet]; |
335 | CoinBigIndex start = startColumnGen_[sequence]; |
336 | addColumn(startColumnGen_[sequence+1] - start, |
337 | rowGen_ + start, |
338 | elementGen_ + start, |
339 | costGen_[sequence], |
340 | columnLowerGen_ ? columnLowerGen_[sequence] : 0, |
341 | columnUpperGen_ ? columnUpperGen_[sequence] : 1.0e30, |
342 | iSet, getDynamicStatusGen(sequence)); |
343 | idGen_[iSet] = sequence; // say which one in |
344 | setDynamicStatusGen(sequence, inSmall); |
345 | } |
346 | } else { |
347 | // put back old ones |
348 | int * set = new int[numberColumns_]; |
349 | for (iSet = 0; iSet < numberSets_; iSet++) { |
350 | for (CoinBigIndex j = fullStartGen_[iSet]; j < fullStartGen_[iSet+1]; j++) |
351 | set[j] = iSet; |
352 | } |
353 | for (int i = 0; i < numberIds; i++) { |
354 | int sequence = ids[i]; |
355 | int iSet = set[sequence]; |
356 | CoinBigIndex start = startColumnGen_[sequence]; |
357 | addColumn(startColumnGen_[sequence+1] - start, |
358 | rowGen_ + start, |
359 | elementGen_ + start, |
360 | costGen_[sequence], |
361 | columnLowerGen_ ? columnLowerGen_[sequence] : 0, |
362 | columnUpperGen_ ? columnUpperGen_[sequence] : 1.0e30, |
363 | iSet, getDynamicStatus(i)); |
364 | idGen_[i] = sequence; // say which one in |
365 | setDynamicStatusGen(sequence, inSmall); |
366 | } |
367 | delete [] set; |
368 | } |
369 | if (!status) { |
370 | gubCrash(); |
371 | } else { |
372 | initialProblem(); |
373 | } |
374 | } |
375 | //------------------------------------------------------------------- |
376 | // Destructor |
377 | //------------------------------------------------------------------- |
378 | ClpDynamicExampleMatrix::~ClpDynamicExampleMatrix () |
379 | { |
380 | delete [] startColumnGen_; |
381 | delete [] rowGen_; |
382 | delete [] elementGen_; |
383 | delete [] costGen_; |
384 | delete [] fullStartGen_; |
385 | delete [] dynamicStatusGen_; |
386 | delete [] idGen_; |
387 | delete [] columnLowerGen_; |
388 | delete [] columnUpperGen_; |
389 | } |
390 | |
391 | //---------------------------------------------------------------- |
392 | // Assignment operator |
393 | //------------------------------------------------------------------- |
394 | ClpDynamicExampleMatrix & |
395 | ClpDynamicExampleMatrix::operator=(const ClpDynamicExampleMatrix& rhs) |
396 | { |
397 | if (this != &rhs) { |
398 | ClpDynamicMatrix::operator=(rhs); |
399 | numberColumns_ = rhs.numberColumns_; |
400 | delete [] startColumnGen_; |
401 | delete [] rowGen_; |
402 | delete [] elementGen_; |
403 | delete [] costGen_; |
404 | delete [] fullStartGen_; |
405 | delete [] dynamicStatusGen_; |
406 | delete [] idGen_; |
407 | delete [] columnLowerGen_; |
408 | delete [] columnUpperGen_; |
409 | startColumnGen_ = ClpCopyOfArray(rhs.startColumnGen_, numberColumns_ + 1); |
410 | CoinBigIndex numberElements = startColumnGen_[numberColumns_]; |
411 | rowGen_ = ClpCopyOfArray(rhs.rowGen_, numberElements); |
412 | elementGen_ = ClpCopyOfArray(rhs.elementGen_, numberElements); |
413 | costGen_ = ClpCopyOfArray(rhs.costGen_, numberColumns_); |
414 | fullStartGen_ = ClpCopyOfArray(rhs.fullStartGen_, numberSets_ + 1); |
415 | dynamicStatusGen_ = ClpCopyOfArray(rhs.dynamicStatusGen_, numberColumns_); |
416 | idGen_ = ClpCopyOfArray(rhs.idGen_, maximumGubColumns_); |
417 | columnLowerGen_ = ClpCopyOfArray(rhs.columnLowerGen_, numberColumns_); |
418 | columnUpperGen_ = ClpCopyOfArray(rhs.columnUpperGen_, numberColumns_); |
419 | } |
420 | return *this; |
421 | } |
422 | //------------------------------------------------------------------- |
423 | // Clone |
424 | //------------------------------------------------------------------- |
425 | ClpMatrixBase * ClpDynamicExampleMatrix::clone() const |
426 | { |
427 | return new ClpDynamicExampleMatrix(*this); |
428 | } |
429 | // Partial pricing |
430 | void |
431 | ClpDynamicExampleMatrix::partialPricing(ClpSimplex * model, double startFraction, double endFraction, |
432 | int & bestSequence, int & numberWanted) |
433 | { |
434 | numberWanted = currentWanted_; |
435 | assert(!model->rowScale()); |
436 | if (!numberSets_) { |
437 | // no gub |
438 | ClpPackedMatrix::partialPricing(model, startFraction, endFraction, bestSequence, numberWanted); |
439 | } else { |
440 | // and do some proportion of full set |
441 | int startG2 = static_cast<int> (startFraction * numberSets_); |
442 | int endG2 = static_cast<int> (endFraction * numberSets_ + 0.1); |
443 | endG2 = CoinMin(endG2, numberSets_); |
444 | //printf("gub price - set start %d end %d\n", |
445 | // startG2,endG2); |
446 | double tolerance = model->currentDualTolerance(); |
447 | double * reducedCost = model->djRegion(); |
448 | const double * duals = model->dualRowSolution(); |
449 | double bestDj; |
450 | int numberRows = model->numberRows(); |
451 | int slackOffset = lastDynamic_ + numberRows; |
452 | int structuralOffset = slackOffset + numberSets_; |
453 | int structuralOffset2 = structuralOffset + maximumGubColumns_; |
454 | // If nothing found yet can go all the way to end |
455 | int endAll = endG2; |
456 | if (bestSequence < 0 && !startG2) |
457 | endAll = numberSets_; |
458 | if (bestSequence >= 0) { |
459 | if (bestSequence != savedBestSequence_) |
460 | bestDj = fabs(reducedCost[bestSequence]); // dj from slacks or permanent |
461 | else |
462 | bestDj = savedBestDj_; |
463 | } else { |
464 | bestDj = tolerance; |
465 | } |
466 | int saveSequence = bestSequence; |
467 | double djMod = 0.0; |
468 | double bestDjMod = 0.0; |
469 | //printf("iteration %d start %d end %d - wanted %d\n",model->numberIterations(), |
470 | // startG2,endG2,numberWanted); |
471 | int bestSet = -1; |
472 | int minSet = minimumObjectsScan_ < 0 ? 5 : minimumObjectsScan_; |
473 | int minNeg = minimumGoodReducedCosts_ < 0 ? 5 : minimumGoodReducedCosts_; |
474 | for (int iSet = startG2; iSet < endAll; iSet++) { |
475 | if (numberWanted + minNeg < originalWanted_ && iSet > startG2 + minSet) { |
476 | // give up |
477 | numberWanted = 0; |
478 | break; |
479 | } else if (iSet == endG2 && bestSequence >= 0) { |
480 | break; |
481 | } |
482 | int gubRow = toIndex_[iSet]; |
483 | if (gubRow >= 0) { |
484 | djMod = duals[gubRow+numberStaticRows_]; // have I got sign right? |
485 | } else { |
486 | int iBasic = keyVariable_[iSet]; |
487 | if (iBasic >= numberColumns_) { |
488 | djMod = 0.0; // set not in |
489 | } else { |
490 | // get dj without |
491 | djMod = 0.0; |
492 | for (CoinBigIndex j = startColumn_[iBasic]; |
493 | j < startColumn_[iBasic+1]; j++) { |
494 | int jRow = row_[j]; |
495 | djMod -= duals[jRow] * element_[j]; |
496 | } |
497 | djMod += cost_[iBasic]; |
498 | // See if gub slack possible - dj is djMod |
499 | if (getStatus(iSet) == ClpSimplex::atLowerBound) { |
500 | double value = -djMod; |
501 | if (value > tolerance) { |
502 | numberWanted--; |
503 | if (value > bestDj) { |
504 | // check flagged variable and correct dj |
505 | if (!flagged(iSet)) { |
506 | bestDj = value; |
507 | bestSequence = slackOffset + iSet; |
508 | bestDjMod = djMod; |
509 | bestSet = iSet; |
510 | } else { |
511 | // just to make sure we don't exit before got something |
512 | numberWanted++; |
513 | abort(); |
514 | } |
515 | } |
516 | } |
517 | } else if (getStatus(iSet) == ClpSimplex::atUpperBound) { |
518 | double value = djMod; |
519 | if (value > tolerance) { |
520 | numberWanted--; |
521 | if (value > bestDj) { |
522 | // check flagged variable and correct dj |
523 | if (!flagged(iSet)) { |
524 | bestDj = value; |
525 | bestSequence = slackOffset + iSet; |
526 | bestDjMod = djMod; |
527 | bestSet = iSet; |
528 | } else { |
529 | // just to make sure we don't exit before got something |
530 | numberWanted++; |
531 | abort(); |
532 | } |
533 | } |
534 | } |
535 | } |
536 | } |
537 | } |
538 | // do ones in small |
539 | int iSequence = startSet_[iSet]; |
540 | while (iSequence >= 0) { |
541 | DynamicStatus status = getDynamicStatus(iSequence); |
542 | if (status == atLowerBound || status == atUpperBound) { |
543 | double value = cost_[iSequence] - djMod; |
544 | for (CoinBigIndex j = startColumn_[iSequence]; |
545 | j < startColumn_[iSequence+1]; j++) { |
546 | int jRow = row_[j]; |
547 | value -= duals[jRow] * element_[j]; |
548 | } |
549 | // change sign if at lower bound |
550 | if (status == atLowerBound) |
551 | value = -value; |
552 | if (value > tolerance) { |
553 | numberWanted--; |
554 | if (value > bestDj) { |
555 | // check flagged variable and correct dj |
556 | if (!flagged(iSequence)) { |
557 | bestDj = value; |
558 | bestSequence = structuralOffset + iSequence; |
559 | bestDjMod = djMod; |
560 | bestSet = iSet; |
561 | } else { |
562 | // just to make sure we don't exit before got something |
563 | numberWanted++; |
564 | } |
565 | } |
566 | } |
567 | } |
568 | iSequence = next_[iSequence]; //onto next in set |
569 | } |
570 | // and now get best by column generation |
571 | // If no upper bounds we may not need status test |
572 | for (iSequence = fullStartGen_[iSet]; iSequence < fullStartGen_[iSet+1]; iSequence++) { |
573 | DynamicStatus status = getDynamicStatusGen(iSequence); |
574 | assert (status != atUpperBound && status != soloKey); |
575 | if (status == atLowerBound) { |
576 | double value = costGen_[iSequence] - djMod; |
577 | for (CoinBigIndex j = startColumnGen_[iSequence]; |
578 | j < startColumnGen_[iSequence+1]; j++) { |
579 | int jRow = rowGen_[j]; |
580 | value -= duals[jRow] * elementGen_[j]; |
581 | } |
582 | // change sign as at lower bound |
583 | value = -value; |
584 | if (value > tolerance) { |
585 | numberWanted--; |
586 | if (value > bestDj) { |
587 | // check flagged variable and correct dj |
588 | if (!flaggedGen(iSequence)) { |
589 | bestDj = value; |
590 | bestSequence = structuralOffset2 + iSequence; |
591 | bestDjMod = djMod; |
592 | bestSet = iSet; |
593 | } else { |
594 | // just to make sure we don't exit before got something |
595 | numberWanted++; |
596 | } |
597 | } |
598 | } |
599 | } |
600 | } |
601 | if (numberWanted <= 0) { |
602 | numberWanted = 0; |
603 | break; |
604 | } |
605 | } |
606 | if (bestSequence != saveSequence) { |
607 | savedBestGubDual_ = bestDjMod; |
608 | savedBestDj_ = bestDj; |
609 | savedBestSequence_ = bestSequence; |
610 | savedBestSet_ = bestSet; |
611 | } |
612 | // Do packed part before gub |
613 | // always??? |
614 | // Resize so just do to gub |
615 | numberActiveColumns_ = firstDynamic_; |
616 | int saveMinNeg = minimumGoodReducedCosts_; |
617 | if (bestSequence >= 0) |
618 | minimumGoodReducedCosts_ = -2; |
619 | currentWanted_ = numberWanted; |
620 | ClpPackedMatrix::partialPricing(model, startFraction, endFraction, bestSequence, numberWanted); |
621 | numberActiveColumns_ = matrix_->getNumCols(); |
622 | minimumGoodReducedCosts_ = saveMinNeg; |
623 | // See if may be finished |
624 | if (!startG2 && bestSequence < 0) |
625 | infeasibilityWeight_ = model_->infeasibilityCost(); |
626 | else if (bestSequence >= 0) |
627 | infeasibilityWeight_ = -1.0; |
628 | currentWanted_ = numberWanted; |
629 | } |
630 | } |
631 | /* Creates a variable. This is called after partial pricing and may modify matrix. |
632 | May update bestSequence. |
633 | */ |
634 | void |
635 | ClpDynamicExampleMatrix::createVariable(ClpSimplex * model, int & bestSequence) |
636 | { |
637 | int numberRows = model->numberRows(); |
638 | int slackOffset = lastDynamic_ + numberRows; |
639 | int structuralOffset = slackOffset + numberSets_; |
640 | int bestSequence2 = savedBestSequence_ - structuralOffset; |
641 | if (bestSequence2 >= 0) { |
642 | // See if needs new |
643 | if (bestSequence2 >= maximumGubColumns_) { |
644 | bestSequence2 -= maximumGubColumns_; |
645 | int sequence = addColumn(startColumnGen_[bestSequence2+1] - startColumnGen_[bestSequence2], |
646 | rowGen_ + startColumnGen_[bestSequence2], |
647 | elementGen_ + startColumnGen_[bestSequence2], |
648 | costGen_[bestSequence2], |
649 | columnLowerGen_ ? columnLowerGen_[bestSequence2] : 0, |
650 | columnUpperGen_ ? columnUpperGen_[bestSequence2] : 1.0e30, |
651 | savedBestSet_, getDynamicStatusGen(bestSequence2)); |
652 | savedBestSequence_ = structuralOffset + sequence; |
653 | idGen_[sequence] = bestSequence2; |
654 | setDynamicStatusGen(bestSequence2, inSmall); |
655 | } |
656 | } |
657 | ClpDynamicMatrix::createVariable(model, bestSequence/*, bestSequence2*/); |
658 | // clear for next iteration |
659 | savedBestSequence_ = -1; |
660 | } |
661 | /* If addColumn forces compression then this allows descendant to know what to do. |
662 | If >=0 then entry stayed in, if -1 then entry went out to lower bound.of zero. |
663 | Entries at upper bound (really nonzero) never go out (at present). |
664 | */ |
665 | void |
666 | ClpDynamicExampleMatrix::packDown(const int * in, int numberToPack) |
667 | { |
668 | int put = 0; |
669 | for (int i = 0; i < numberToPack; i++) { |
670 | int id = idGen_[i]; |
671 | if (in[i] >= 0) { |
672 | // stays in |
673 | assert (put == in[i]); // true for now |
674 | idGen_[put++] = id; |
675 | } else { |
676 | // out to lower bound |
677 | setDynamicStatusGen(id, atLowerBound); |
678 | } |
679 | } |
680 | assert (put == numberGubColumns_); |
681 | } |
682 | |