1 | #include <stdio.h> |
---|
2 | #include <math.h> |
---|
3 | |
---|
4 | #include "PresolveMatrix.hpp" |
---|
5 | #include "PresolveEmpty.hpp" // for DROP_COL/DROP_ROW |
---|
6 | #include "PresolveFixed.hpp" |
---|
7 | #include "PresolveSubst.hpp" |
---|
8 | #include "PresolveUseless.hpp" |
---|
9 | #include "PresolveForcing.hpp" |
---|
10 | #include "ClpMessage.hpp" |
---|
11 | |
---|
12 | #if 0 |
---|
13 | inline double max(double x, double y) |
---|
14 | { |
---|
15 | return (x < y) ? y : x; |
---|
16 | } |
---|
17 | |
---|
18 | inline double min(double x, double y) |
---|
19 | { |
---|
20 | return (x < y) ? x : y; |
---|
21 | } |
---|
22 | #endif |
---|
23 | |
---|
24 | /*static*/ void implied_bounds(const double *els, |
---|
25 | const double *clo, const double *cup, |
---|
26 | const int *hcol, |
---|
27 | CoinBigIndex krs, CoinBigIndex kre, |
---|
28 | double *maxupp, double *maxdownp, |
---|
29 | int jcol, |
---|
30 | double rlo, double rup, |
---|
31 | double *iclb, double *icub) |
---|
32 | { |
---|
33 | if (rlo<=-PRESOLVE_INF&&rup>=PRESOLVE_INF) { |
---|
34 | *iclb = -PRESOLVE_INF; |
---|
35 | *icub = PRESOLVE_INF; |
---|
36 | return; |
---|
37 | } |
---|
38 | bool posinf = false; |
---|
39 | bool neginf = false; |
---|
40 | double maxup = 0.0; |
---|
41 | double maxdown = 0.0; |
---|
42 | |
---|
43 | int jcolk = -1; |
---|
44 | |
---|
45 | // compute sum of all bounds except for jcol |
---|
46 | CoinBigIndex kk; |
---|
47 | for (kk=krs; kk<kre; kk++) { |
---|
48 | if (hcol[kk] == jcol) |
---|
49 | jcolk = kk; |
---|
50 | |
---|
51 | // swap jcol with hcol[kre-1]; |
---|
52 | // that is, consider jcol last |
---|
53 | // this assumes that jcol occurs in this row |
---|
54 | CoinBigIndex k = (hcol[kk] == jcol |
---|
55 | ? kre-1 |
---|
56 | : kk == kre-1 |
---|
57 | ? jcolk |
---|
58 | : kk); |
---|
59 | |
---|
60 | PRESOLVEASSERT(k != -1); // i.e. jcol had better be in the row |
---|
61 | |
---|
62 | int col = hcol[k]; |
---|
63 | double coeff = els[k]; |
---|
64 | double lb = clo[col]; |
---|
65 | double ub = cup[col]; |
---|
66 | |
---|
67 | // quick! compute the implied col bounds before maxup/maxdown are changed |
---|
68 | if (kk == kre-1) { |
---|
69 | PRESOLVEASSERT(fabs(coeff) > ZTOLDP); |
---|
70 | |
---|
71 | double ilb = (rlo - maxup) / coeff; |
---|
72 | bool finite_ilb = (-PRESOLVE_INF < rlo && !posinf && PRESOLVEFINITE(maxup)); |
---|
73 | |
---|
74 | double iub = (rup - maxdown) / coeff; |
---|
75 | bool finite_iub = ( rup < PRESOLVE_INF && !neginf && PRESOLVEFINITE(maxdown)); |
---|
76 | |
---|
77 | if (coeff > 0.0) { |
---|
78 | *iclb = (finite_ilb ? ilb : -PRESOLVE_INF); |
---|
79 | *icub = (finite_iub ? iub : PRESOLVE_INF); |
---|
80 | } else { |
---|
81 | *iclb = (finite_iub ? iub : -PRESOLVE_INF); |
---|
82 | *icub = (finite_ilb ? ilb : PRESOLVE_INF); |
---|
83 | } |
---|
84 | } |
---|
85 | |
---|
86 | if (coeff > 0.0) { |
---|
87 | if (PRESOLVE_INF <= ub) { |
---|
88 | posinf = true; |
---|
89 | if (neginf) |
---|
90 | break; // pointless |
---|
91 | } else |
---|
92 | maxup += ub * coeff; |
---|
93 | |
---|
94 | if (lb <= -PRESOLVE_INF) { |
---|
95 | neginf = true; |
---|
96 | if (posinf) |
---|
97 | break; // pointless |
---|
98 | } else |
---|
99 | maxdown += lb * coeff; |
---|
100 | } else { |
---|
101 | if (PRESOLVE_INF <= ub) { |
---|
102 | neginf = true; |
---|
103 | if (posinf) |
---|
104 | break; // pointless |
---|
105 | } else |
---|
106 | maxdown += ub * coeff; |
---|
107 | |
---|
108 | if (lb <= -PRESOLVE_INF) { |
---|
109 | posinf = true; |
---|
110 | if (neginf) |
---|
111 | break; // pointless |
---|
112 | } else |
---|
113 | maxup += lb * coeff; |
---|
114 | } |
---|
115 | } |
---|
116 | |
---|
117 | // If we broke from the loop, then the bounds are infinite. |
---|
118 | // However, since we put the column whose implied bounds we want |
---|
119 | // to know at the end, and it doesn't matter if its own bounds |
---|
120 | // are infinite, don't worry about breaking at the last iteration. |
---|
121 | if (kk<kre-1) { |
---|
122 | *iclb = -PRESOLVE_INF; |
---|
123 | *icub = PRESOLVE_INF; |
---|
124 | } else |
---|
125 | PRESOLVEASSERT(jcolk != -1); |
---|
126 | |
---|
127 | // store row bounds |
---|
128 | *maxupp = (posinf) ? PRESOLVE_INF : maxup; |
---|
129 | *maxdownp = (neginf) ? -PRESOLVE_INF : maxdown; |
---|
130 | } |
---|
131 | |
---|
132 | static void implied_row_bounds(const double *els, |
---|
133 | const double *clo, const double *cup, |
---|
134 | const int *hcol, |
---|
135 | CoinBigIndex krs, CoinBigIndex kre, |
---|
136 | double *maxupp, double *maxdownp) |
---|
137 | { |
---|
138 | // int jcol; |
---|
139 | double iclb, icub; |
---|
140 | |
---|
141 | implied_bounds(els, clo, cup, hcol, krs, kre, maxupp, maxdownp, |
---|
142 | hcol[krs], 0.0, 0.0, &iclb, &icub); |
---|
143 | } |
---|
144 | |
---|
145 | const char *forcing_constraint_action::name() const |
---|
146 | { |
---|
147 | return ("forcing_constraint_action"); |
---|
148 | } |
---|
149 | |
---|
150 | static bool some_col_was_fixed(const int *hcol, CoinBigIndex krs, CoinBigIndex kre, |
---|
151 | const double *clo, |
---|
152 | const double *cup) |
---|
153 | { |
---|
154 | CoinBigIndex k; |
---|
155 | for (k=krs; k<kre; k++) { |
---|
156 | int jcol = hcol[k]; |
---|
157 | if (clo[jcol] == cup[jcol]) |
---|
158 | break; |
---|
159 | } |
---|
160 | return (k<kre); |
---|
161 | } |
---|
162 | |
---|
163 | |
---|
164 | // |
---|
165 | // It may be the case that the variable bounds are such that no matter what |
---|
166 | // feasible value they take, the constraint cannot be violated; |
---|
167 | // in this case, we obviously don't need to take it into account, and |
---|
168 | // we just drop it as a USELESS constraint. |
---|
169 | // |
---|
170 | // On the other hand, it may be that the only way to satisfy a constraint |
---|
171 | // is to jam all the constraint variables to one of their bounds; |
---|
172 | // in this case, these variables are essentially fixed, which |
---|
173 | // we do with a FORCING constraint. |
---|
174 | // For now, this just tightens the bounds; subsequently the fixed variables |
---|
175 | // will be removed, then the row will be dropped. |
---|
176 | // |
---|
177 | // Since both operations use similar information (the implied row bounds), |
---|
178 | // we combine them into one presolve routine. |
---|
179 | // |
---|
180 | // I claim that these checks could be performed in parallel, |
---|
181 | // that is, the tests could be carried out for all rows in parallel, |
---|
182 | // and then the rows deleted and columns tightened afterward. |
---|
183 | // Obviously, this is true for useless rows. |
---|
184 | // The potential problem is forcing constraints, but I think |
---|
185 | // that is ok. |
---|
186 | // By doing it in parallel rather than sequentially, |
---|
187 | // we may miss transformations due to variables that were fixed |
---|
188 | // by forcing constraints, though. |
---|
189 | // |
---|
190 | // Note that both of these operations will cause problems |
---|
191 | // if the variables in question really need to exceed their bounds in order |
---|
192 | // to make the problem feasible. |
---|
193 | const PresolveAction *forcing_constraint_action::presolve(PresolveMatrix *prob, |
---|
194 | const PresolveAction *next) |
---|
195 | { |
---|
196 | double *clo = prob->clo_; |
---|
197 | double *cup = prob->cup_; |
---|
198 | |
---|
199 | const double *rowels = prob->rowels_; |
---|
200 | const int *hcol = prob->hcol_; |
---|
201 | const CoinBigIndex *mrstrt = prob->mrstrt_; |
---|
202 | const int *hinrow = prob->hinrow_; |
---|
203 | const int nrows = prob->nrows_; |
---|
204 | |
---|
205 | const double *rlo = prob->rlo_; |
---|
206 | const double *rup = prob->rup_; |
---|
207 | |
---|
208 | // const char *integerType = prob->integerType_; |
---|
209 | |
---|
210 | const double tol = ZTOLDP; |
---|
211 | const double inftol = prob->feasibilityTolerance_; |
---|
212 | const int ncols = prob->ncols_; |
---|
213 | |
---|
214 | int *fixed_cols = new int[ncols]; |
---|
215 | int nfixed_cols = 0; |
---|
216 | |
---|
217 | action *actions = new action [nrows]; |
---|
218 | int nactions = 0; |
---|
219 | |
---|
220 | int *useless_rows = new int[nrows]; |
---|
221 | int nuseless_rows = 0; |
---|
222 | |
---|
223 | int numberLook = prob->numberRowsToDo_; |
---|
224 | int iLook; |
---|
225 | int * look = prob->rowsToDo_; |
---|
226 | |
---|
227 | for (iLook=0;iLook<numberLook;iLook++) { |
---|
228 | int irow = look[iLook]; |
---|
229 | if (hinrow[irow] > 0) { |
---|
230 | CoinBigIndex krs = mrstrt[irow]; |
---|
231 | CoinBigIndex kre = krs + hinrow[irow]; |
---|
232 | |
---|
233 | double maxup, maxdown; |
---|
234 | implied_row_bounds(rowels, clo, cup, hcol, krs, kre, |
---|
235 | &maxup, &maxdown); |
---|
236 | |
---|
237 | if (maxup < PRESOLVE_INF && maxup + inftol < rlo[irow]) { |
---|
238 | /* there is an upper bound and it can't be reached */ |
---|
239 | prob->status_|= 1; |
---|
240 | prob->messageHandler()->message(CLP_PRESOLVE_ROWINFEAS, |
---|
241 | prob->messages()) |
---|
242 | <<irow |
---|
243 | <<rlo[irow] |
---|
244 | <<rup[irow] |
---|
245 | <<CoinMessageEol; |
---|
246 | break; |
---|
247 | } else if (-PRESOLVE_INF < maxdown && rup[irow] < maxdown - inftol) { |
---|
248 | /* there is a lower bound and it can't be reached */ |
---|
249 | prob->status_|= 1; |
---|
250 | prob->messageHandler()->message(CLP_PRESOLVE_ROWINFEAS, |
---|
251 | prob->messages()) |
---|
252 | <<irow |
---|
253 | <<rlo[irow] |
---|
254 | <<rup[irow] |
---|
255 | <<CoinMessageEol; |
---|
256 | break; |
---|
257 | } |
---|
258 | // ADD TOLERANCE TO THESE TESTS |
---|
259 | else if ((rlo[irow] <= -PRESOLVE_INF || |
---|
260 | (-PRESOLVE_INF < maxdown && rlo[irow] <= maxdown)) && |
---|
261 | (rup[irow] >= PRESOLVE_INF || |
---|
262 | (maxup < PRESOLVE_INF && rup[irow] >= maxup))) { |
---|
263 | |
---|
264 | // I'm not sure that these transforms don't intefere with each other |
---|
265 | // We can get it next time |
---|
266 | if (some_col_was_fixed(hcol, krs, kre, clo, cup)) { |
---|
267 | // make sure on next time |
---|
268 | prob->addRow(irow); |
---|
269 | continue; |
---|
270 | } |
---|
271 | |
---|
272 | // this constraint must always be satisfied - drop it |
---|
273 | useless_rows[nuseless_rows++] = irow; |
---|
274 | |
---|
275 | } else if ((maxup < PRESOLVE_INF && fabs(rlo[irow] - maxup) < tol) || |
---|
276 | (-PRESOLVE_INF < maxdown && fabs(rup[irow] - maxdown) < tol)) { |
---|
277 | |
---|
278 | // the lower bound can just be reached, or |
---|
279 | // the upper bound can just be reached; |
---|
280 | // called a "forcing constraint" in the paper (p. 226) |
---|
281 | const int lbound_tight = (maxup < PRESOLVE_INF && |
---|
282 | fabs(rlo[irow] - maxup) < tol); |
---|
283 | |
---|
284 | // I'm not sure that these transforms don't intefere with each other |
---|
285 | // We can get it next time |
---|
286 | if (some_col_was_fixed(hcol, krs, kre, clo, cup)) { |
---|
287 | // make sure on next time |
---|
288 | prob->addRow(irow); |
---|
289 | continue; |
---|
290 | } |
---|
291 | |
---|
292 | // out of space - this probably never happens |
---|
293 | if (nfixed_cols + (kre-krs) >= ncols) |
---|
294 | break; |
---|
295 | |
---|
296 | double *bounds = new double[hinrow[irow]]; |
---|
297 | int *rowcols = new int[hinrow[irow]]; |
---|
298 | int lk = krs; // load fix-to-down in front |
---|
299 | int uk = kre; // load fix-to-up in back |
---|
300 | CoinBigIndex k; |
---|
301 | for ( k=krs; k<kre; k++) { |
---|
302 | int jcol = hcol[k]; |
---|
303 | prob->addCol(jcol); |
---|
304 | double coeff = rowels[k]; |
---|
305 | |
---|
306 | PRESOLVEASSERT(fabs(coeff) > ZTOLDP); |
---|
307 | |
---|
308 | // one of the two contributed to maxup - set the other to that |
---|
309 | if (lbound_tight == (coeff > 0.0)) { |
---|
310 | --uk; |
---|
311 | bounds[uk-krs] = clo[jcol]; |
---|
312 | rowcols[uk-krs] = jcol; |
---|
313 | |
---|
314 | clo[jcol] = cup[jcol]; |
---|
315 | } else { |
---|
316 | bounds[lk-krs] = cup[jcol]; |
---|
317 | rowcols[lk-krs] = jcol; |
---|
318 | ++lk; |
---|
319 | |
---|
320 | cup[jcol] = clo[jcol]; |
---|
321 | } |
---|
322 | |
---|
323 | fixed_cols[nfixed_cols++] = jcol; |
---|
324 | } |
---|
325 | PRESOLVEASSERT(uk == lk); |
---|
326 | |
---|
327 | action *f = &actions[nactions]; |
---|
328 | nactions++; |
---|
329 | |
---|
330 | f->row = irow; |
---|
331 | f->nlo = lk-krs; |
---|
332 | f->nup = kre-uk; |
---|
333 | f->rowcols = rowcols; |
---|
334 | f->bounds = bounds; |
---|
335 | } |
---|
336 | } |
---|
337 | } |
---|
338 | |
---|
339 | |
---|
340 | if (nactions) { |
---|
341 | #if PRESOLVE_SUMMARY |
---|
342 | printf("NFORCED: %d\n", nactions); |
---|
343 | #endif |
---|
344 | next = new forcing_constraint_action(nactions, |
---|
345 | copyOfArray(actions,nactions), next); |
---|
346 | } |
---|
347 | deleteAction(actions,action*); |
---|
348 | if (nuseless_rows) { |
---|
349 | next = useless_constraint_action::presolve(prob, |
---|
350 | useless_rows, nuseless_rows, |
---|
351 | next); |
---|
352 | } |
---|
353 | delete[]useless_rows; |
---|
354 | |
---|
355 | if (nfixed_cols) { |
---|
356 | next = remove_fixed_action::presolve(prob, fixed_cols, nfixed_cols, next); |
---|
357 | } |
---|
358 | delete[]fixed_cols; |
---|
359 | |
---|
360 | return (next); |
---|
361 | } |
---|
362 | |
---|
363 | void forcing_constraint_action::postsolve(PostsolveMatrix *prob) const |
---|
364 | { |
---|
365 | const action *const actions = actions_; |
---|
366 | const int nactions = nactions_; |
---|
367 | |
---|
368 | const double *colels = prob->colels_; |
---|
369 | const int *hrow = prob->hrow_; |
---|
370 | const CoinBigIndex *mcstrt = prob->mcstrt_; |
---|
371 | const int *hincol = prob->hincol_; |
---|
372 | const int *link = prob->link_; |
---|
373 | |
---|
374 | // CoinBigIndex free_list = prob->free_list_; |
---|
375 | |
---|
376 | double *clo = prob->clo_; |
---|
377 | double *cup = prob->cup_; |
---|
378 | |
---|
379 | const double *sol = prob->sol_; |
---|
380 | double *rcosts = prob->rcosts_; |
---|
381 | |
---|
382 | //double *acts = prob->acts_; |
---|
383 | double *rowduals = prob->rowduals_; |
---|
384 | |
---|
385 | const double ztoldj = prob->ztoldj_; |
---|
386 | const double ztolzb = prob->ztolzb_; |
---|
387 | |
---|
388 | for (const action *f = &actions[nactions-1]; actions<=f; f--) { |
---|
389 | |
---|
390 | const int irow = f->row; |
---|
391 | const int nlo = f->nlo; |
---|
392 | const int nup = f->nup; |
---|
393 | const int ninrow = nlo + nup; |
---|
394 | const int *rowcols = f->rowcols; |
---|
395 | const double *bounds= f->bounds; |
---|
396 | int k; |
---|
397 | |
---|
398 | for (k=0; k<nlo; k++) { |
---|
399 | int jcol = rowcols[k]; |
---|
400 | cup[jcol] = bounds[k]; |
---|
401 | PRESOLVEASSERT(prob->getColumnStatus(jcol)!=PrePostsolveMatrix::basic); |
---|
402 | } |
---|
403 | |
---|
404 | for (k=nlo; k<ninrow; k++) { |
---|
405 | int jcol = rowcols[k]; |
---|
406 | clo[jcol] = bounds[k]; |
---|
407 | PRESOLVEASSERT(prob->getColumnStatus(jcol)!=PrePostsolveMatrix::basic); |
---|
408 | } |
---|
409 | |
---|
410 | PRESOLVEASSERT(prob->getRowStatus(irow)==PrePostsolveMatrix::basic); |
---|
411 | PRESOLVEASSERT(rowduals[irow] == 0.0); |
---|
412 | |
---|
413 | // this is a lazy implementation. |
---|
414 | // we tightened the col bounds, then let them be eliminated |
---|
415 | // by repeated uses of FIX_VARIABLE and a final DROP_ROW. |
---|
416 | // Therefore, by this point the row has been marked basic, |
---|
417 | // the rowdual of this row is 0.0, |
---|
418 | // and the reduced costs for the cols may or may not be ok |
---|
419 | // for the relaxed column bounds. |
---|
420 | // |
---|
421 | // find the one most out of whack and fix it. |
---|
422 | int whacked = -1; |
---|
423 | double whack = 0.0; |
---|
424 | for (k=0; k<ninrow; k++) { |
---|
425 | int jcol = rowcols[k]; |
---|
426 | CoinBigIndex kk = presolve_find_row2(irow, mcstrt[jcol], hincol[jcol], hrow, link); |
---|
427 | |
---|
428 | // choose rowdual to cancel out reduced cost |
---|
429 | double whack0 = rcosts[jcol] / colels[kk]; |
---|
430 | |
---|
431 | if (((rcosts[jcol] > ztoldj && !(fabs(sol[jcol] - clo[jcol]) <= ztolzb)) || |
---|
432 | (rcosts[jcol] < -ztoldj && !(fabs(sol[jcol] - cup[jcol]) <= ztolzb))) && |
---|
433 | fabs(whack0) > fabs(whack)) { |
---|
434 | whacked = jcol; |
---|
435 | whack = whack0; |
---|
436 | } |
---|
437 | } |
---|
438 | |
---|
439 | if (whacked != -1) { |
---|
440 | prob->setColumnStatus(whacked,PrePostsolveMatrix::basic); |
---|
441 | prob->setRowStatus(irow,PrePostsolveMatrix::atLowerBound); |
---|
442 | rowduals[irow] = whack; |
---|
443 | |
---|
444 | for (k=0; k<ninrow; k++) { |
---|
445 | int jcol = rowcols[k]; |
---|
446 | CoinBigIndex kk = presolve_find_row2(irow, mcstrt[jcol], hincol[jcol], hrow, link); |
---|
447 | |
---|
448 | rcosts[jcol] -= (rowduals[irow] * colels[kk]); |
---|
449 | } |
---|
450 | } |
---|
451 | } |
---|
452 | } |
---|
453 | |
---|
454 | |
---|
455 | |
---|
456 | #if 0 |
---|
457 | // Determine the maximum and minimum values the constraint sums |
---|
458 | // may take, given the bounds on the variables. |
---|
459 | // If there are infinite terms, record where the first one is, |
---|
460 | // and whether there is more than one. |
---|
461 | // It is possible to compute implied bounds for the (one) variable |
---|
462 | // with no bound. |
---|
463 | static void implied_bounds1(PresolveMatrix * prob, const double *rowels, |
---|
464 | const int *mrstrt, |
---|
465 | const int *hrow, |
---|
466 | const int *hinrow, |
---|
467 | const double *clo, const double *cup, |
---|
468 | const int *hcol, |
---|
469 | int ncols, |
---|
470 | const double *rlo, const double *rup, |
---|
471 | const char *integerType, |
---|
472 | int nrows, |
---|
473 | double *ilbound, double *iubound) |
---|
474 | { |
---|
475 | const double tol = prob->feasibilityTolerance_; |
---|
476 | |
---|
477 | for (int irow=0; irow<nrows; irow++) { |
---|
478 | CoinBigIndex krs = mrstrt[irow]; |
---|
479 | CoinBigIndex kre = krs + hinrow[irow]; |
---|
480 | |
---|
481 | double irlo = rlo[irow]; |
---|
482 | double irup = rup[irow]; |
---|
483 | |
---|
484 | // These are used to set column bounds below. |
---|
485 | // If there are no (positive) infinite terms, |
---|
486 | // the loop will range from krs to kre; |
---|
487 | // if there is just one, it will range over that one variable; |
---|
488 | // otherwise, it will be empty. |
---|
489 | int ub_inf_index = -1; |
---|
490 | int lb_inf_index = -1; |
---|
491 | |
---|
492 | double maxup = 0.0; |
---|
493 | double maxdown = 0.0; |
---|
494 | CoinBigIndex k; |
---|
495 | for (k=krs; k<kre; k++) { |
---|
496 | int jcol = hcol[k]; |
---|
497 | double coeff = rowels[k]; |
---|
498 | double lb = clo[jcol]; |
---|
499 | double ub = cup[jcol]; |
---|
500 | |
---|
501 | // HAVE TO DEAL WITH BOUNDS OF INTEGER VARIABLES |
---|
502 | if (coeff > 0.0) { |
---|
503 | if (PRESOLVE_INF <= ub) { |
---|
504 | if (ub_inf_index == -1) { |
---|
505 | ub_inf_index = k; |
---|
506 | } else { |
---|
507 | ub_inf_index = -2; |
---|
508 | if (lb_inf_index == -2) |
---|
509 | break; // pointless |
---|
510 | } |
---|
511 | } else |
---|
512 | maxup += ub * coeff; |
---|
513 | |
---|
514 | if (lb <= -PRESOLVE_INF) { |
---|
515 | if (lb_inf_index == -1) { |
---|
516 | lb_inf_index = k; |
---|
517 | } else { |
---|
518 | lb_inf_index = -2; |
---|
519 | if (ub_inf_index == -2) |
---|
520 | break; // pointless |
---|
521 | } |
---|
522 | } else |
---|
523 | maxdown += lb * coeff; |
---|
524 | } |
---|
525 | else { |
---|
526 | if (PRESOLVE_INF <= ub) { |
---|
527 | if (lb_inf_index == -1) { |
---|
528 | lb_inf_index = k; |
---|
529 | } else { |
---|
530 | lb_inf_index = -2; |
---|
531 | if (ub_inf_index == -2) |
---|
532 | break; // pointless |
---|
533 | } |
---|
534 | } else |
---|
535 | maxdown += ub * coeff; |
---|
536 | |
---|
537 | if (lb <= -PRESOLVE_INF) { |
---|
538 | if (ub_inf_index == -1) { |
---|
539 | ub_inf_index = k; |
---|
540 | } else { |
---|
541 | ub_inf_index = -2; |
---|
542 | if (lb_inf_index == -2) |
---|
543 | break; // pointless |
---|
544 | } |
---|
545 | } else |
---|
546 | maxup += lb * coeff; |
---|
547 | } |
---|
548 | } |
---|
549 | |
---|
550 | // ub_inf says whether the sum of the "other" ub terms is infinite |
---|
551 | // in the loop below. |
---|
552 | // In the case where we only saw one infinite term, the loop |
---|
553 | // will only cover that case, in which case the other terms |
---|
554 | // are *not* infinite. |
---|
555 | // With two or more terms, it is infinite. |
---|
556 | // If we only saw one infinite term, then |
---|
557 | if (ub_inf_index == -2) |
---|
558 | maxup = PRESOLVE_INF; |
---|
559 | |
---|
560 | if (lb_inf_index == -2) |
---|
561 | maxdown = -PRESOLVE_INF; |
---|
562 | |
---|
563 | const bool maxup_finite = PRESOLVEFINITE(maxup); |
---|
564 | const bool maxdown_finite = PRESOLVEFINITE(maxdown); |
---|
565 | |
---|
566 | if (ub_inf_index == -1 && maxup_finite && maxup + tol < rlo[irow]) { |
---|
567 | /* infeasible */ |
---|
568 | prob->status_|= 1; |
---|
569 | prob->messageHandler()->message(CLP_PRESOLVE_ROWINFEAS, |
---|
570 | prob->messages()) |
---|
571 | <<irow |
---|
572 | <<rlo[irow] |
---|
573 | <<rup[irow] |
---|
574 | <<CoinMessageEol; |
---|
575 | break; |
---|
576 | } else if (lb_inf_index == -1 && maxdown_finite && rup[irow] < maxdown - tol) { |
---|
577 | /* infeasible */ |
---|
578 | prob->status_|= 1; |
---|
579 | prob->messageHandler()->message(CLP_PRESOLVE_ROWINFEAS, |
---|
580 | prob->messages()) |
---|
581 | <<irow |
---|
582 | <<rlo[irow] |
---|
583 | <<rup[irow] |
---|
584 | <<CoinMessageEol; |
---|
585 | break; |
---|
586 | } |
---|
587 | |
---|
588 | for (k = krs; k<kre; ++k) { |
---|
589 | int jcol = hcol[k]; |
---|
590 | double coeff = rowels[k]; |
---|
591 | |
---|
592 | // SHOULD GET RID OF THIS |
---|
593 | if (fabs(coeff) > ZTOLDP && |
---|
594 | !integerType[jcol]) { |
---|
595 | double maxup1 = (ub_inf_index == -1 || ub_inf_index == k |
---|
596 | ? maxup |
---|
597 | : PRESOLVE_INF); |
---|
598 | bool maxup_finite1 = (ub_inf_index == -1 || ub_inf_index == k |
---|
599 | ? maxup_finite |
---|
600 | : false); |
---|
601 | double maxdown1 = (lb_inf_index == -1 || lb_inf_index == k |
---|
602 | ? maxdown |
---|
603 | : PRESOLVE_INF); |
---|
604 | bool maxdown_finite1 = (ub_inf_index == -1 || ub_inf_index == k |
---|
605 | ? maxdown_finite |
---|
606 | : false); |
---|
607 | |
---|
608 | double ilb = (irlo - maxup1) / coeff; |
---|
609 | bool finite_ilb = (-PRESOLVE_INF < irlo && maxup_finite1); |
---|
610 | |
---|
611 | double iub = (irup - maxdown1) / coeff; |
---|
612 | bool finite_iub = ( irup < PRESOLVE_INF && maxdown_finite1); |
---|
613 | |
---|
614 | double ilb1 = (coeff > 0.0 |
---|
615 | ? (finite_ilb ? ilb : -PRESOLVE_INF) |
---|
616 | : (finite_iub ? iub : -PRESOLVE_INF)); |
---|
617 | |
---|
618 | if (ilbound[jcol] < ilb1) { |
---|
619 | ilbound[jcol] = ilb1; |
---|
620 | if (jcol == 278001) |
---|
621 | printf("JCOL LB %g\n", ilb1); |
---|
622 | } |
---|
623 | } |
---|
624 | } |
---|
625 | |
---|
626 | for (k = krs; k<kre; ++k) { |
---|
627 | int jcol = hcol[k]; |
---|
628 | double coeff = rowels[k]; |
---|
629 | |
---|
630 | // SHOULD GET RID OF THIS |
---|
631 | if (fabs(coeff) > ZTOLDP && |
---|
632 | !integerType[jcol]) { |
---|
633 | double maxup1 = (ub_inf_index == -1 || ub_inf_index == k |
---|
634 | ? maxup |
---|
635 | : PRESOLVE_INF); |
---|
636 | bool maxup_finite1 = (ub_inf_index == -1 || ub_inf_index == k |
---|
637 | ? maxup_finite |
---|
638 | : false); |
---|
639 | double maxdown1 = (lb_inf_index == -1 || lb_inf_index == k |
---|
640 | ? maxdown |
---|
641 | : PRESOLVE_INF); |
---|
642 | bool maxdown_finite1 = (ub_inf_index == -1 || ub_inf_index == k |
---|
643 | ? maxdown_finite |
---|
644 | : false); |
---|
645 | |
---|
646 | |
---|
647 | double ilb = (irlo - maxup1) / coeff; |
---|
648 | bool finite_ilb = (-PRESOLVE_INF < irlo && maxup_finite1); |
---|
649 | |
---|
650 | double iub = (irup - maxdown1) / coeff; |
---|
651 | bool finite_iub = ( irup < PRESOLVE_INF && maxdown_finite1); |
---|
652 | |
---|
653 | double iub1 = (coeff > 0.0 |
---|
654 | ? (finite_iub ? iub : PRESOLVE_INF) |
---|
655 | : (finite_ilb ? ilb : PRESOLVE_INF)); |
---|
656 | |
---|
657 | if (iub1 < iubound[jcol]) { |
---|
658 | iubound[jcol] = iub1; |
---|
659 | if (jcol == 278001) |
---|
660 | printf("JCOL UB %g\n", iub1); |
---|
661 | } |
---|
662 | } |
---|
663 | } |
---|
664 | } |
---|
665 | } |
---|
666 | |
---|
667 | #if 0 |
---|
668 | postsolve for implied_bound |
---|
669 | { |
---|
670 | double lo0 = pa->clo; |
---|
671 | double up0 = pa->cup; |
---|
672 | int irow = pa->irow; |
---|
673 | int jcol = pa->icol; |
---|
674 | int *rowcols = pa->rowcols; |
---|
675 | int ninrow = pa->ninrow; |
---|
676 | |
---|
677 | clo[jcol] = lo0; |
---|
678 | cup[jcol] = up0; |
---|
679 | |
---|
680 | if ((colstat[jcol] & PRESOLVE_XBASIC) == 0 && |
---|
681 | fabs(lo0 - sol[jcol]) > ztolzb && |
---|
682 | fabs(up0 - sol[jcol]) > ztolzb) { |
---|
683 | |
---|
684 | // this non-basic variable is now away from its bound |
---|
685 | // it is ok just to force it to be basic |
---|
686 | // informally: if this variable is at its implied bound, |
---|
687 | // then the other variables must be at their bounds, |
---|
688 | // which means the bounds will stop them even if the aren't basic. |
---|
689 | if (rowstat[irow] & PRESOLVE_XBASIC) |
---|
690 | rowstat[irow] = 0; |
---|
691 | else { |
---|
692 | int k; |
---|
693 | for (k=0; k<ninrow; k++) { |
---|
694 | int col = rowcols[k]; |
---|
695 | if (cdone[col] && |
---|
696 | (colstat[col] & PRESOLVE_XBASIC) && |
---|
697 | ((fabs(clo[col] - sol[col]) <= ztolzb && rcosts[col] >= -ztoldj) || |
---|
698 | (fabs(cup[col] - sol[col]) <= ztolzb && rcosts[col] <= ztoldj))) |
---|
699 | break; |
---|
700 | } |
---|
701 | if (k<ninrow) { |
---|
702 | int col = rowcols[k]; |
---|
703 | // steal this basic variable |
---|
704 | #if DEBUG_PRESOLVE |
---|
705 | printf("PIVOTING ON COL: %d %d -> %d\n", irow, col, jcol); |
---|
706 | #endif |
---|
707 | colstat[col] = 0; |
---|
708 | |
---|
709 | // since all vars were at their bounds, the slack must be 0 |
---|
710 | PRESOLVEASSERT(fabs(acts[irow]) < ZTOLDP); |
---|
711 | rowstat[irow] = PRESOLVE_XBASIC; |
---|
712 | } |
---|
713 | else { |
---|
714 | // should never happen? |
---|
715 | abort(); |
---|
716 | } |
---|
717 | |
---|
718 | // get rid of any remaining basic structurals, since their rcosts |
---|
719 | // are going to become non-zero in a second. |
---|
720 | abort(); |
---|
721 | ///////////////////zero_pivot(); |
---|
722 | } |
---|
723 | |
---|
724 | double rdual_adjust; |
---|
725 | { |
---|
726 | CoinBigIndex kk = presolve_find_row(irow, mcstrt[jcol], mcstrt[jcol] + hincol[jcol], hrow); |
---|
727 | // adjust rowdual to cancel out reduced cost |
---|
728 | // should probably search for col with largest factor |
---|
729 | rdual_adjust = (rcosts[jcol] / colels[kk]); |
---|
730 | rowduals[irow] += rdual_adjust; |
---|
731 | colstat[jcol] = PRESOLVE_XBASIC; |
---|
732 | } |
---|
733 | |
---|
734 | for (k=0; k<ninrow; k++) { |
---|
735 | int jcol = rowcols[k]; |
---|
736 | CoinBigIndex kk = presolve_find_row(irow, mcstrt[jcol], mcstrt[jcol] + hincol[jcol], hrow); |
---|
737 | |
---|
738 | rcosts[jcol] -= (rdual_adjust * colels[kk]); |
---|
739 | } |
---|
740 | |
---|
741 | { |
---|
742 | int k; |
---|
743 | int badbasic = -1; |
---|
744 | |
---|
745 | // we may have just screwed up the rcost of another basic variable |
---|
746 | for (k=0; k<ninrow; k++) { |
---|
747 | int col = rowcols[k]; |
---|
748 | if (col != jcol && |
---|
749 | cdone[col] && |
---|
750 | (colstat[col] & PRESOLVE_XBASIC) && |
---|
751 | !(fabs(rcosts[col]) < ztoldj)) |
---|
752 | if (badbasic == -1) |
---|
753 | badbasic = k; |
---|
754 | else |
---|
755 | abort(); // two!! what to do??? |
---|
756 | } |
---|
757 | |
---|
758 | if (badbasic != -1) { |
---|
759 | int col = rowcols[badbasic]; |
---|
760 | |
---|
761 | if (fabs(acts[irow]) < ZTOLDP) { |
---|
762 | #if DEBUG_PRESOLVE |
---|
763 | printf("PIVOTING COL TO SLACK!: %d %d\n", irow, col); |
---|
764 | #endif |
---|
765 | colstat[col] = 0; |
---|
766 | rowstat[irow] = PRESOLVE_XBASIC; |
---|
767 | } |
---|
768 | else |
---|
769 | abort(); |
---|
770 | } |
---|
771 | } |
---|
772 | } |
---|
773 | } |
---|
774 | #endif |
---|
775 | #endif |
---|
776 | |
---|
777 | |
---|
778 | |
---|