1 /* * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * */
2 /* */
3 /* This file is part of the program and library */
4 /* SCIP --- Solving Constraint Integer Programs */
5 /* */
6 /* Copyright (C) 2002-2021 Konrad-Zuse-Zentrum */
7 /* fuer Informationstechnik Berlin */
8 /* */
9 /* SCIP is distributed under the terms of the ZIB Academic License. */
10 /* */
11 /* You should have received a copy of the ZIB Academic License */
12 /* along with SCIP; see the file COPYING. If not visit scipopt.org. */
13 /* */
14 /* * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * */
15
16 /**@file presol_domcol.c
17 * @ingroup DEFPLUGINS_PRESOL
18 * @brief dominated column presolver
19 * @author Dieter Weninger
20 * @author Gerald Gamrath
21 *
22 * This presolver looks for dominance relations between variable pairs.
23 * From a dominance relation and certain bound/clique-constellations
24 * variable fixings mostly at the lower bound of the dominated variable can be derived.
25 * Additionally it is possible to improve bounds by predictive bound strengthening.
26 *
27 * @todo Also run on general CIPs, if the number of locks of the investigated variables comes only from (upgraded)
28 * linear constraints.
29 *
30 * @todo Instead of choosing variables for comparison out of one row, we should try to use 'hashvalues' for columns that
31 * indicate in which constraint type (<=, >=, or ranged row / ==) they are existing. Then sort the variables (and
32 * corresponding data) after the ranged row/equation hashvalue and only try to derive dominance on variables with
33 * the same hashvalue on ranged row/equation and fitting hashvalues for the other constraint types.
34 * @todo run on incomplete matrices; in order to do so, check at the time when dominance is detected that the locks are
35 * consistent; probably, it is sufficient to check one lock direction for each of the two variables
36 *
37 */
38
39 /*---+----1----+----2----+----3----+----4----+----5----+----6----+----7----+----8----+----9----+----0----+----1----+----2*/
40
41 #include "blockmemshell/memory.h"
42 #include "scip/presol_domcol.h"
43 #include "scip/pub_matrix.h"
44 #include "scip/pub_message.h"
45 #include "scip/pub_misc_sort.h"
46 #include "scip/pub_presol.h"
47 #include "scip/pub_var.h"
48 #include "scip/scip_general.h"
49 #include "scip/scip_mem.h"
50 #include "scip/scip_message.h"
51 #include "scip/scip_nlp.h"
52 #include "scip/scip_numerics.h"
53 #include "scip/scip_param.h"
54 #include "scip/scip_presol.h"
55 #include "scip/scip_pricer.h"
56 #include "scip/scip_prob.h"
57 #include "scip/scip_probing.h"
58 #include "scip/scip_var.h"
59 #include <string.h>
60
61 #define PRESOL_NAME "domcol"
62 #define PRESOL_DESC "dominated column presolver"
63 #define PRESOL_PRIORITY -1000 /**< priority of the presolver (>= 0: before, < 0: after constraint handlers) */
64 #define PRESOL_MAXROUNDS -1 /**< maximal number of presolving rounds the presolver participates in (-1: no limit) */
65 #define PRESOL_TIMING SCIP_PRESOLTIMING_EXHAUSTIVE /* timing of the presolver (fast, medium, or exhaustive) */
66
67 #define DEFAULT_NUMMINPAIRS 1024 /**< minimal number of pair comparisons */
68 #define DEFAULT_NUMMAXPAIRS 1048576 /**< maximal number of pair comparisons */
69
70 #define DEFAULT_PREDBNDSTR FALSE /**< should predictive bound strengthening be applied? */
71 #define DEFAULT_CONTINUOUS_RED TRUE /**< should reductions for continuous variables be carried out? */
72
73
74
75 /*
76 * Data structures
77 */
78
79 /** control parameters */
80 struct SCIP_PresolData
81 {
82 int numminpairs; /**< minimal number of pair comparisons */
83 int nummaxpairs; /**< maximal number of pair comparisons */
84 int numcurrentpairs; /**< current number of pair comparisons */
85 SCIP_Bool predbndstr; /**< flag indicating if predictive bound strengthening should be applied */
86 SCIP_Bool continuousred; /**< flag indicating if reductions for continuous variables should be performed */
87 };
88
89 /** type of fixing direction */
90 enum Fixingdirection
91 {
92 FIXATLB = -1, /**< fix variable at lower bound */
93 NOFIX = 0, /**< do not fix variable */
94 FIXATUB = 1 /**< fix variable at upper bound */
95 };
96 typedef enum Fixingdirection FIXINGDIRECTION;
97
98
99 /*
100 * Local methods
101 */
102 #ifdef SCIP_DEBUG
103 /** print a row from the constraint matrix */
104 static
printRow(SCIP * scip,SCIP_MATRIX * matrix,int row)105 void printRow(
106 SCIP* scip, /**< SCIP main data structure */
107 SCIP_MATRIX* matrix, /**< matrix containing the constraints */
108 int row /**< row index for printing */
109 )
110 {
111 int* rowpnt;
112 int* rowend;
113 int c;
114 SCIP_Real val;
115 SCIP_Real* valpnt;
116 char relation;
117 SCIP_VAR* var;
118
119 relation='-';
120 if( !SCIPmatrixIsRowRhsInfinity(matrix, row) &&
121 SCIPisEQ(scip, SCIPmatrixGetRowLhs(matrix, row), SCIPmatrixGetRowRhs(matrix, row)) )
122 {
123 relation='e';
124 }
125 else if( !SCIPmatrixIsRowRhsInfinity(matrix, row) &&
126 !SCIPisEQ(scip, SCIPmatrixGetRowLhs(matrix, row), SCIPmatrixGetRowRhs(matrix, row)) )
127 {
128 relation='r';
129 }
130 else
131 {
132 relation='g';
133 }
134
135 rowpnt = SCIPmatrixGetRowIdxPtr(matrix, row);
136 rowend = rowpnt + SCIPmatrixGetRowNNonzs(matrix, row);
137 valpnt = SCIPmatrixGetRowValPtr(matrix, row);
138
139 SCIPdebugMsgPrint(scip, "\n(L:%g) [%c] %g <=",
140 (SCIPmatrixGetRowNMinActPosInf(matrix, row) + SCIPmatrixGetRowNMinActNegInf(matrix,row) > 0) ?
141 -SCIPinfinity(scip) :
142 SCIPmatrixGetRowMinActivity(matrix, row), relation, SCIPmatrixGetRowLhs(matrix, row));
143 for(; (rowpnt < rowend); rowpnt++, valpnt++)
144 {
145 c = *rowpnt;
146 val = *valpnt;
147 var = SCIPmatrixGetVar(matrix, c);
148 SCIPdebugMsgPrint(scip, " %g{%s[idx:%d][bnd:%g,%g]}",
149 val, SCIPvarGetName(var), c, SCIPvarGetLbGlobal(var), SCIPvarGetUbGlobal(var));
150 }
151 SCIPdebugMsgPrint(scip, " <= %g (U:%g)", (SCIPmatrixGetRowNMaxActPosInf(matrix, row) + SCIPmatrixGetRowNMaxActNegInf(matrix, row) > 0) ?
152 SCIPinfinity(scip) :
153 SCIPmatrixGetRowRhs(matrix, row), SCIPmatrixGetRowMaxActivity(matrix , row));
154 }
155
156 /** print all rows from the constraint matrix containing a variable */
157 static
printRowsOfCol(SCIP * scip,SCIP_MATRIX * matrix,int col)158 SCIP_RETCODE printRowsOfCol(
159 SCIP* scip, /**< SCIP main data structure */
160 SCIP_MATRIX* matrix, /**< matrix containing the constraints */
161 int col /**< column index for printing */
162 )
163 {
164 int numrows;
165 int* rows;
166 int i;
167 int* colpnt;
168 int* colend;
169
170 numrows = SCIPmatrixGetColNNonzs(matrix, col);
171
172 SCIP_CALL( SCIPallocBufferArray(scip, &rows, numrows) );
173
174 colpnt = SCIPmatrixGetColIdxPtr(matrix, col);
175 colend = colpnt + SCIPmatrixGetColNNonzs(matrix, col);
176 for( i = 0; (colpnt < colend); colpnt++, i++ )
177 {
178 rows[i] = *colpnt;
179 }
180
181 SCIPdebugMsgPrint(scip, "\n-------");
182 SCIPdebugMsgPrint(scip, "\ncol %d number rows: %d",col,numrows);
183 for( i = 0; i < numrows; i++ )
184 {
185 printRow(scip, matrix, rows[i]);
186 }
187 SCIPdebugMsgPrint(scip, "\n-------");
188
189 SCIPfreeBufferArray(scip, &rows);
190
191 return SCIP_OKAY;
192 }
193
194 /** print information about a dominance relation */
195 static
printDomRelInfo(SCIP * scip,SCIP_MATRIX * matrix,SCIP_VAR * dominatingvar,int dominatingidx,SCIP_VAR * dominatedvar,int dominatedidx,SCIP_Real dominatingub,SCIP_Real dominatingwclb)196 SCIP_RETCODE printDomRelInfo(
197 SCIP* scip, /**< SCIP main data structure */
198 SCIP_MATRIX* matrix, /**< matrix containing the constraints */
199 SCIP_VAR* dominatingvar, /**< dominating variable */
200 int dominatingidx, /**< index of dominating variable */
201 SCIP_VAR* dominatedvar, /**< dominated variable */
202 int dominatedidx, /**< index of dominated variable */
203 SCIP_Real dominatingub, /**< predicted upper bound of dominating variable */
204 SCIP_Real dominatingwclb /**< worst case lower bound of dominating variable */
205 )
206 {
207 char type;
208
209 assert(SCIPvarGetType(dominatingvar)==SCIPvarGetType(dominatedvar));
210
211 switch(SCIPvarGetType(dominatingvar))
212 {
213 case SCIP_VARTYPE_CONTINUOUS:
214 type='C';
215 break;
216 case SCIP_VARTYPE_BINARY:
217 type='B';
218 break;
219 case SCIP_VARTYPE_INTEGER:
220 case SCIP_VARTYPE_IMPLINT:
221 type='I';
222 break;
223 default:
224 SCIPerrorMessage("unknown variable type\n");
225 SCIPABORT();
226 return SCIP_INVALIDDATA; /*lint !e527*/
227 }
228
229 SCIPdebugMsgPrint(scip, "\n\n### [%c], obj:%g->%g,\t%s[idx:%d](nrows:%d)->%s[idx:%d](nrows:%d)\twclb=%g, ub'=%g, ub=%g",
230 type, SCIPvarGetObj(dominatingvar), SCIPvarGetObj(dominatedvar),
231 SCIPvarGetName(dominatingvar), dominatingidx, SCIPmatrixGetColNNonzs(matrix, dominatingidx),
232 SCIPvarGetName(dominatedvar), dominatedidx, SCIPmatrixGetColNNonzs(matrix, dominatedidx),
233 dominatingwclb, dominatingub, SCIPvarGetUbGlobal(dominatingvar));
234
235 SCIP_CALL( printRowsOfCol(scip, matrix, dominatingidx) );
236
237 return SCIP_OKAY;
238 }
239 #endif
240
241
242 /** get minimum/maximum residual activity for the specified variable and with another variable set to its upper bound */
243 static
getActivityResidualsUpperBound(SCIP * scip,SCIP_MATRIX * matrix,int row,int col,SCIP_Real coef,int upperboundcol,SCIP_Real upperboundcoef,SCIP_Real * minresactivity,SCIP_Real * maxresactivity,SCIP_Bool * success)244 void getActivityResidualsUpperBound(
245 SCIP* scip,
246 SCIP_MATRIX* matrix,
247 int row,
248 int col,
249 SCIP_Real coef,
250 int upperboundcol,
251 SCIP_Real upperboundcoef,
252 SCIP_Real* minresactivity,
253 SCIP_Real* maxresactivity,
254 SCIP_Bool* success
255 )
256 {
257 SCIP_VAR* var;
258 SCIP_VAR* upperboundvar;
259 SCIP_Real lb;
260 SCIP_Real ub;
261 SCIP_Real lbupperboundvar;
262 SCIP_Real ubupperboundvar;
263 SCIP_Real tmpmaxact;
264 SCIP_Real tmpminact;
265 int maxactinf;
266 int minactinf;
267
268 assert(scip != NULL);
269 assert(matrix != NULL);
270 assert(row < SCIPmatrixGetNRows(matrix));
271 assert(minresactivity != NULL);
272 assert(maxresactivity != NULL);
273
274 var = SCIPmatrixGetVar(matrix, col);
275 upperboundvar = SCIPmatrixGetVar(matrix, upperboundcol);
276 assert(var != NULL);
277 assert(upperboundvar != NULL);
278
279 lb = SCIPvarGetLbGlobal(var);
280 ub = SCIPvarGetUbGlobal(var);
281
282 maxactinf = SCIPmatrixGetRowNMaxActPosInf(matrix, row) + SCIPmatrixGetRowNMaxActNegInf(matrix, row);
283 minactinf = SCIPmatrixGetRowNMinActPosInf(matrix, row) + SCIPmatrixGetRowNMinActNegInf(matrix, row);
284
285 assert(!SCIPisInfinity(scip, lb));
286 assert(!SCIPisInfinity(scip, -ub));
287
288 lbupperboundvar = SCIPvarGetLbGlobal(upperboundvar);
289 ubupperboundvar = SCIPvarGetUbGlobal(upperboundvar);
290
291 assert(!SCIPisInfinity(scip, lbupperboundvar));
292 assert(!SCIPisInfinity(scip, -ubupperboundvar));
293
294 tmpminact = SCIPmatrixGetRowMinActivity(matrix, row);
295 tmpmaxact = SCIPmatrixGetRowMaxActivity(matrix, row);
296
297 /* first, adjust min and max activity such that upperboundvar is always set to its upper bound */
298
299 /* abort if the upper bound is infty */
300 if( SCIPisInfinity(scip, ubupperboundvar) )
301 {
302 *success = FALSE;
303 return;
304 }
305 else
306 {
307 /* coef > 0 --> the lower bound is used for the minactivity --> update the minactivity */
308 if( upperboundcoef > 0.0 )
309 {
310 if( SCIPisInfinity(scip, -lbupperboundvar) )
311 {
312 assert(minactinf > 0);
313 minactinf--;
314 tmpminact += (upperboundcoef * ubupperboundvar);
315 }
316 else
317 {
318 tmpminact = tmpminact - (upperboundcoef * lbupperboundvar) + (upperboundcoef * ubupperboundvar);
319 }
320 }
321 /* coef < 0 --> the lower bound is used for the maxactivity --> update the maxactivity */
322 else
323 {
324 if( SCIPisInfinity(scip, -lbupperboundvar) )
325 {
326 assert(maxactinf > 0);
327 maxactinf--;
328 tmpmaxact += (upperboundcoef * ubupperboundvar);
329 }
330 else
331 {
332 tmpmaxact = tmpmaxact - (upperboundcoef * lbupperboundvar) + (upperboundcoef * ubupperboundvar);
333 }
334 }
335 }
336
337 *success = TRUE;
338
339 /* the coefficient is positive, so the upper bound contributed to the maxactivity and the lower bound to the minactivity */
340 if( coef >= 0.0 )
341 {
342 if( SCIPisInfinity(scip, ub) )
343 {
344 assert(maxactinf >= 1);
345 if( maxactinf == 1 )
346 {
347 *maxresactivity = tmpmaxact;
348 }
349 else
350 *maxresactivity = SCIPinfinity(scip);
351 }
352 else
353 {
354 if( maxactinf > 0 )
355 *maxresactivity = SCIPinfinity(scip);
356 else
357 *maxresactivity = tmpmaxact - coef * ub;
358 }
359
360 if( SCIPisInfinity(scip, -lb) )
361 {
362 assert(minactinf >= 1);
363 if( minactinf == 1 )
364 {
365 *minresactivity = tmpminact;
366 }
367 else
368 *minresactivity = -SCIPinfinity(scip);
369 }
370 else
371 {
372 if( minactinf > 0 )
373 *minresactivity = -SCIPinfinity(scip);
374 else
375 *minresactivity = tmpminact - coef * lb;
376 }
377 }
378 /* the coefficient is negative, so the lower bound contributed to the maxactivity and the upper bound to the minactivity */
379 else
380 {
381 if( SCIPisInfinity(scip, -lb) )
382 {
383 assert(maxactinf >= 1);
384 if( maxactinf == 1 )
385 {
386 *maxresactivity = tmpmaxact;
387 }
388 else
389 *maxresactivity = SCIPinfinity(scip);
390 }
391 else
392 {
393 if( maxactinf > 0 )
394 *maxresactivity = SCIPinfinity(scip);
395 else
396 *maxresactivity = tmpmaxact - coef * lb;
397 }
398
399 if( SCIPisInfinity(scip, ub) )
400 {
401 assert(minactinf >= 1);
402 if( minactinf == 1 )
403 {
404 *minresactivity = tmpminact;
405 }
406 else
407 *minresactivity = -SCIPinfinity(scip);
408 }
409 else
410 {
411 if( minactinf > 0 )
412 *minresactivity = -SCIPinfinity(scip);
413 else
414 *minresactivity = tmpminact - coef * ub;
415 }
416 }
417 }
418
419 /** get minimum/maximum residual activity for the specified variable and with another variable set to its lower bound */
420 static
getActivityResidualsLowerBound(SCIP * scip,SCIP_MATRIX * matrix,int row,int col,SCIP_Real coef,int lowerboundcol,SCIP_Real lowerboundcoef,SCIP_Real * minresactivity,SCIP_Real * maxresactivity,SCIP_Bool * success)421 void getActivityResidualsLowerBound(
422 SCIP* scip, /**< SCIP main data structure */
423 SCIP_MATRIX* matrix, /**< matrix containing the constraints */
424 int row, /**< row index */
425 int col, /**< column index */
426 SCIP_Real coef, /**< coefficient of the column in this row */
427 int lowerboundcol, /**< column index of variable to set to its lower bound */
428 SCIP_Real lowerboundcoef, /**< coefficient in this row of the column to be set to its lower bound */
429 SCIP_Real* minresactivity, /**< minimum residual activity of this row */
430 SCIP_Real* maxresactivity, /**< maximum residual activity of this row */
431 SCIP_Bool* success /**< pointer to store whether the computation was successful */
432 )
433 {
434 SCIP_VAR* var;
435 SCIP_VAR* lowerboundvar;
436 SCIP_Real lb;
437 SCIP_Real ub;
438 SCIP_Real lblowerboundvar;
439 SCIP_Real ublowerboundvar;
440 SCIP_Real tmpmaxact;
441 SCIP_Real tmpminact;
442 int maxactinf;
443 int minactinf;
444
445 assert(scip != NULL);
446 assert(matrix != NULL);
447 assert(row < SCIPmatrixGetNRows(matrix));
448 assert(minresactivity != NULL);
449 assert(maxresactivity != NULL);
450
451 var = SCIPmatrixGetVar(matrix, col);;
452 lowerboundvar = SCIPmatrixGetVar(matrix, lowerboundcol);
453 assert(var != NULL);
454 assert(lowerboundvar != NULL);
455
456 lb = SCIPvarGetLbGlobal(var);
457 ub = SCIPvarGetUbGlobal(var);
458
459 maxactinf = SCIPmatrixGetRowNMaxActPosInf(matrix, row) + SCIPmatrixGetRowNMaxActNegInf(matrix, row);
460 minactinf = SCIPmatrixGetRowNMinActPosInf(matrix, row) + SCIPmatrixGetRowNMinActNegInf(matrix, row);
461
462 assert(!SCIPisInfinity(scip, lb));
463 assert(!SCIPisInfinity(scip, -ub));
464
465 lblowerboundvar = SCIPvarGetLbGlobal(lowerboundvar);
466 ublowerboundvar = SCIPvarGetUbGlobal(lowerboundvar);
467
468 assert(!SCIPisInfinity(scip, lblowerboundvar));
469 assert(!SCIPisInfinity(scip, -ublowerboundvar));
470
471 tmpminact = SCIPmatrixGetRowMinActivity(matrix, row);
472 tmpmaxact = SCIPmatrixGetRowMaxActivity(matrix, row);
473
474 /* first, adjust min and max activity such that lowerboundvar is always set to its lower bound */
475
476 /* abort if the lower bound is -infty */
477 if( SCIPisInfinity(scip, -lblowerboundvar) )
478 {
479 *success = FALSE;
480 return;
481 }
482 else
483 {
484 /* coef > 0 --> the upper bound is used for the maxactivity --> update the maxactivity */
485 if( lowerboundcoef > 0.0 )
486 {
487 if( SCIPisInfinity(scip, ublowerboundvar) )
488 {
489 assert(maxactinf > 0);
490 maxactinf--;
491 tmpmaxact += (lowerboundcoef * lblowerboundvar);
492 }
493 else
494 {
495 tmpmaxact = tmpmaxact - (lowerboundcoef * ublowerboundvar) + (lowerboundcoef * lblowerboundvar);
496 }
497 }
498 /* coef < 0 --> the upper bound is used for the minactivity --> update the minactivity */
499 else
500 {
501 if( SCIPisInfinity(scip, ublowerboundvar) )
502 {
503 assert(minactinf > 0);
504 minactinf--;
505 tmpminact += (lowerboundcoef * lblowerboundvar);
506 }
507 else
508 {
509 tmpminact = tmpminact - (lowerboundcoef * ublowerboundvar) + (lowerboundcoef * lblowerboundvar);
510 }
511 }
512 }
513
514 *success = TRUE;
515
516 /* the coefficient is positive, so the upper bound contributed to the maxactivity and the lower bound to the minactivity */
517 if( coef >= 0.0 )
518 {
519 if( SCIPisInfinity(scip, ub) )
520 {
521 assert(maxactinf >= 1);
522 if( maxactinf == 1 )
523 {
524 *maxresactivity = tmpmaxact;
525 }
526 else
527 *maxresactivity = SCIPinfinity(scip);
528 }
529 else
530 {
531 if( maxactinf > 0 )
532 *maxresactivity = SCIPinfinity(scip);
533 else
534 *maxresactivity = tmpmaxact - coef * ub;
535 }
536
537 if( SCIPisInfinity(scip, -lb) )
538 {
539 assert(minactinf >= 1);
540 if( minactinf == 1 )
541 {
542 *minresactivity = tmpminact;
543 }
544 else
545 *minresactivity = -SCIPinfinity(scip);
546 }
547 else
548 {
549 if( minactinf > 0 )
550 *minresactivity = -SCIPinfinity(scip);
551 else
552 *minresactivity = tmpminact - coef * lb;
553 }
554 }
555 /* the coefficient is negative, so the lower bound contributed to the maxactivity and the upper bound to the minactivity */
556 else
557 {
558 if( SCIPisInfinity(scip, -lb) )
559 {
560 assert(maxactinf >= 1);
561 if( maxactinf == 1 )
562 {
563 *maxresactivity = tmpmaxact;
564 }
565 else
566 *maxresactivity = SCIPinfinity(scip);
567 }
568 else
569 {
570 if( maxactinf > 0 )
571 *maxresactivity = SCIPinfinity(scip);
572 else
573 *maxresactivity = tmpmaxact - coef * lb;
574 }
575
576 if( SCIPisInfinity(scip, ub) )
577 {
578 assert(minactinf >= 1);
579 if( minactinf == 1 )
580 {
581 *minresactivity = tmpminact;
582 }
583 else
584 *minresactivity = -SCIPinfinity(scip);
585 }
586 else
587 {
588 if( minactinf > 0 )
589 *minresactivity = -SCIPinfinity(scip);
590 else
591 *minresactivity = tmpminact - coef * ub;
592 }
593 }
594 }
595
596 /** Calculate bounds of the dominated variable by rowbound analysis.
597 * We use a special kind of predictive rowbound analysis by first setting the dominating variable to its upper bound.
598 */
599 static
calcVarBoundsDominated(SCIP * scip,SCIP_MATRIX * matrix,int row,int coldominating,SCIP_Real valdominating,int coldominated,SCIP_Real valdominated,SCIP_Bool * ubcalculated,SCIP_Real * calculatedub,SCIP_Bool * wclbcalculated,SCIP_Real * calculatedwclb,SCIP_Bool * lbcalculated,SCIP_Real * calculatedlb,SCIP_Bool * wcubcalculated,SCIP_Real * calculatedwcub)600 SCIP_RETCODE calcVarBoundsDominated(
601 SCIP* scip, /**< SCIP main data structure */
602 SCIP_MATRIX* matrix, /**< matrix containing the constraints */
603 int row, /**< current row index */
604 int coldominating, /**< column index of dominating variable */
605 SCIP_Real valdominating, /**< row coefficient of dominating variable */
606 int coldominated, /**< column index of dominated variable */
607 SCIP_Real valdominated, /**< row coefficient of dominated variable */
608 SCIP_Bool* ubcalculated, /**< was a upper bound calculated? */
609 SCIP_Real* calculatedub, /**< predicted upper bound */
610 SCIP_Bool* wclbcalculated, /**< was a lower worst case bound calculated? */
611 SCIP_Real* calculatedwclb, /**< predicted worst case lower bound */
612 SCIP_Bool* lbcalculated, /**< was a lower bound calculated? */
613 SCIP_Real* calculatedlb, /**< predicted lower bound */
614 SCIP_Bool* wcubcalculated, /**< was a worst case upper bound calculated? */
615 SCIP_Real* calculatedwcub /**< calculated worst case upper bound */
616 )
617 {
618 SCIP_Real minresactivity;
619 SCIP_Real maxresactivity;
620 SCIP_Real lhs;
621 SCIP_Real rhs;
622 SCIP_Bool success;
623
624 assert(scip != NULL);
625 assert(matrix != NULL);
626 assert(0 <= row && row < SCIPmatrixGetNRows(matrix));
627 assert(0 <= coldominating && coldominating < SCIPmatrixGetNColumns(matrix));
628 assert(0 <= coldominated && coldominated < SCIPmatrixGetNColumns(matrix));
629
630 assert(ubcalculated != NULL);
631 assert(calculatedub != NULL);
632 assert(wclbcalculated != NULL);
633 assert(calculatedwclb != NULL);
634 assert(lbcalculated != NULL);
635 assert(calculatedlb != NULL);
636 assert(wcubcalculated != NULL);
637 assert(calculatedwcub != NULL);
638
639 assert(!SCIPisZero(scip, valdominated));
640 assert(SCIPmatrixGetVar(matrix, coldominated) != NULL);
641
642 *ubcalculated = FALSE;
643 *wclbcalculated = FALSE;
644 *lbcalculated = FALSE;
645 *wcubcalculated = FALSE;
646
647 /* no rowbound analysis for multiaggregated variables, which should not exist, because the matrix only consists of
648 * active variables
649 */
650 assert(SCIPvarGetStatus(SCIPmatrixGetVar(matrix, coldominating)) != SCIP_VARSTATUS_MULTAGGR);
651 assert(SCIPvarGetStatus(SCIPmatrixGetVar(matrix, coldominated)) != SCIP_VARSTATUS_MULTAGGR);
652
653 lhs = SCIPmatrixGetRowLhs(matrix, row);
654 rhs = SCIPmatrixGetRowRhs(matrix, row);
655 assert(!SCIPisInfinity(scip, lhs));
656 assert(!SCIPisInfinity(scip, -rhs));
657
658 /* get minimum/maximum activity of this row without the dominated variable */
659 getActivityResidualsUpperBound(scip, matrix, row, coldominated, valdominated, coldominating, valdominating,
660 &minresactivity, &maxresactivity, &success);
661
662 if( !success )
663 return SCIP_OKAY;
664
665 assert(!SCIPisInfinity(scip, minresactivity));
666 assert(!SCIPisInfinity(scip, -maxresactivity));
667
668 *calculatedub = SCIPinfinity(scip);
669 *calculatedwclb = -SCIPinfinity(scip);
670 *calculatedlb = -SCIPinfinity(scip);
671 *calculatedwcub = SCIPinfinity(scip);
672
673 /* predictive rowbound analysis */
674
675 if( valdominated > 0.0 )
676 {
677 /* lower bound calculation */
678 if( !SCIPisInfinity(scip, maxresactivity) )
679 {
680 *calculatedlb = (lhs - maxresactivity)/valdominated;
681 *lbcalculated = TRUE;
682 }
683
684 /* worst case calculation of lower bound */
685 if( !SCIPisInfinity(scip, -minresactivity) )
686 {
687 *calculatedwclb = (lhs - minresactivity)/valdominated;
688 *wclbcalculated = TRUE;
689 }
690 else
691 {
692 /* worst case lower bound is infinity */
693 *calculatedwclb = SCIPinfinity(scip);
694 *wclbcalculated = TRUE;
695 }
696
697 /* we can only calculate upper bounds, of the right hand side is finite */
698 if( !SCIPmatrixIsRowRhsInfinity(matrix, row) )
699 {
700 /* upper bound calculation */
701 if( !SCIPisInfinity(scip, -minresactivity) )
702 {
703 *calculatedub = (rhs - minresactivity)/valdominated;
704 *ubcalculated = TRUE;
705 }
706
707 /* worst case calculation of upper bound */
708 if( !SCIPisInfinity(scip, maxresactivity) )
709 {
710 *calculatedwcub = (rhs - maxresactivity)/valdominated;
711 *wcubcalculated = TRUE;
712 }
713 else
714 {
715 /* worst case upper bound is -infinity */
716 *calculatedwcub = -SCIPinfinity(scip);
717 *wcubcalculated = TRUE;
718 }
719 }
720 }
721 else
722 {
723 /* upper bound calculation */
724 if( !SCIPisInfinity(scip, maxresactivity) )
725 {
726 *calculatedub = (lhs - maxresactivity)/valdominated;
727 *ubcalculated = TRUE;
728 }
729
730 /* worst case calculation of upper bound */
731 if( !SCIPisInfinity(scip, -minresactivity) )
732 {
733 *calculatedwcub = (lhs - minresactivity)/valdominated;
734 *wcubcalculated = TRUE;
735 }
736 else
737 {
738 /* worst case upper bound is -infinity */
739 *calculatedwcub = -SCIPinfinity(scip);
740 *wcubcalculated = TRUE;
741 }
742
743 /* we can only calculate lower bounds, of the right hand side is finite */
744 if( !SCIPmatrixIsRowRhsInfinity(matrix, row) )
745 {
746 /* lower bound calculation */
747 if( !SCIPisInfinity(scip, -minresactivity) )
748 {
749 *calculatedlb = (rhs - minresactivity)/valdominated;
750 *lbcalculated = TRUE;
751 }
752
753 /* worst case calculation of lower bound */
754 if( !SCIPisInfinity(scip, maxresactivity) )
755 {
756 *calculatedwclb = (rhs - maxresactivity)/valdominated;
757 *wclbcalculated = TRUE;
758 }
759 else
760 {
761 /* worst case lower bound is infinity */
762 *calculatedwclb = SCIPinfinity(scip);
763 *wclbcalculated = TRUE;
764 }
765 }
766 }
767
768 return SCIP_OKAY;
769 }
770
771 /** Calculate bounds of the dominating variable by rowbound analysis.
772 * We use a special kind of predictive rowbound analysis by first setting the dominated variable to its lower bound.
773 */
774 static
calcVarBoundsDominating(SCIP * scip,SCIP_MATRIX * matrix,int row,int coldominating,SCIP_Real valdominating,int coldominated,SCIP_Real valdominated,SCIP_Bool * ubcalculated,SCIP_Real * calculatedub,SCIP_Bool * wclbcalculated,SCIP_Real * calculatedwclb,SCIP_Bool * lbcalculated,SCIP_Real * calculatedlb,SCIP_Bool * wcubcalculated,SCIP_Real * calculatedwcub)775 SCIP_RETCODE calcVarBoundsDominating(
776 SCIP* scip, /**< SCIP main data structure */
777 SCIP_MATRIX* matrix, /**< matrix containing the constraints */
778 int row, /**< current row index */
779 int coldominating, /**< column index of dominating variable */
780 SCIP_Real valdominating, /**< row coefficient of dominating variable */
781 int coldominated, /**< column index of dominated variable */
782 SCIP_Real valdominated, /**< row coefficient of dominated variable */
783 SCIP_Bool* ubcalculated, /**< was a upper bound calculated? */
784 SCIP_Real* calculatedub, /**< predicted upper bound */
785 SCIP_Bool* wclbcalculated, /**< was a lower worst case bound calculated? */
786 SCIP_Real* calculatedwclb, /**< predicted worst case lower bound */
787 SCIP_Bool* lbcalculated, /**< was a lower bound calculated? */
788 SCIP_Real* calculatedlb, /**< predicted lower bound */
789 SCIP_Bool* wcubcalculated, /**< was a worst case upper bound calculated? */
790 SCIP_Real* calculatedwcub /**< calculated worst case upper bound */
791 )
792 {
793 SCIP_Real minresactivity;
794 SCIP_Real maxresactivity;
795 SCIP_Real lhs;
796 SCIP_Real rhs;
797 SCIP_Bool success;
798
799 assert(scip != NULL);
800 assert(matrix != NULL);
801 assert(0 <= row && row < SCIPmatrixGetNRows(matrix) );
802 assert(0 <= coldominating && coldominating < SCIPmatrixGetNColumns(matrix));
803 assert(0 <= coldominated && coldominated < SCIPmatrixGetNColumns(matrix));
804
805 assert(ubcalculated != NULL);
806 assert(calculatedub != NULL);
807 assert(wclbcalculated != NULL);
808 assert(calculatedwclb != NULL);
809 assert(lbcalculated != NULL);
810 assert(calculatedlb != NULL);
811 assert(wcubcalculated != NULL);
812 assert(calculatedwcub != NULL);
813
814 assert(!SCIPisZero(scip, valdominating));
815 assert(SCIPmatrixGetVar(matrix, coldominating) != NULL);
816
817 *ubcalculated = FALSE;
818 *wclbcalculated = FALSE;
819 *lbcalculated = FALSE;
820 *wcubcalculated = FALSE;
821
822 /* no rowbound analysis for multiaggregated variables, which should not exist, because the matrix only consists of
823 * active variables
824 */
825 assert(SCIPvarGetStatus(SCIPmatrixGetVar(matrix, coldominating)) != SCIP_VARSTATUS_MULTAGGR);
826 assert(SCIPvarGetStatus(SCIPmatrixGetVar(matrix, coldominated)) != SCIP_VARSTATUS_MULTAGGR);
827
828 lhs = SCIPmatrixGetRowLhs(matrix, row);
829 rhs = SCIPmatrixGetRowRhs(matrix, row);
830 assert(!SCIPisInfinity(scip, lhs));
831 assert(!SCIPisInfinity(scip, -rhs));
832
833 /* get minimum/maximum activity of this row without the dominating variable */
834 getActivityResidualsLowerBound(scip, matrix, row, coldominating, valdominating, coldominated, valdominated,
835 &minresactivity, &maxresactivity, &success);
836
837 if( !success )
838 return SCIP_OKAY;
839
840 assert(!SCIPisInfinity(scip, minresactivity));
841 assert(!SCIPisInfinity(scip, -maxresactivity));
842
843 *calculatedub = SCIPinfinity(scip);
844 *calculatedwclb = -SCIPinfinity(scip);
845 *calculatedlb = -SCIPinfinity(scip);
846 *calculatedwcub = SCIPinfinity(scip);
847
848 /* predictive rowbound analysis */
849
850 if( valdominating > 0.0 )
851 {
852 /* lower bound calculation */
853 if( !SCIPisInfinity(scip, maxresactivity) )
854 {
855 *calculatedlb = (lhs - maxresactivity)/valdominating;
856 *lbcalculated = TRUE;
857 }
858
859 /* worst case calculation of lower bound */
860 if( !SCIPisInfinity(scip, -minresactivity) )
861 {
862 *calculatedwclb = (lhs - minresactivity)/valdominating;
863 *wclbcalculated = TRUE;
864 }
865 else
866 {
867 /* worst case lower bound is infinity */
868 *calculatedwclb = SCIPinfinity(scip);
869 *wclbcalculated = TRUE;
870 }
871
872 /* we can only calculate upper bounds, of the right hand side is finite */
873 if( !SCIPmatrixIsRowRhsInfinity(matrix, row) )
874 {
875 /* upper bound calculation */
876 if( !SCIPisInfinity(scip, -minresactivity) )
877 {
878 *calculatedub = (rhs - minresactivity)/valdominating;
879 *ubcalculated = TRUE;
880 }
881
882 /* worst case calculation of upper bound */
883 if( !SCIPisInfinity(scip, maxresactivity) )
884 {
885 *calculatedwcub = (rhs - maxresactivity)/valdominating;
886 *wcubcalculated = TRUE;
887 }
888 else
889 {
890 /* worst case upper bound is -infinity */
891 *calculatedwcub = -SCIPinfinity(scip);
892 *wcubcalculated = TRUE;
893 }
894 }
895 }
896 else
897 {
898 /* upper bound calculation */
899 if( !SCIPisInfinity(scip, maxresactivity) )
900 {
901 *calculatedub = (lhs - maxresactivity)/valdominating;
902 *ubcalculated = TRUE;
903 }
904
905 /* worst case calculation of upper bound */
906 if( !SCIPisInfinity(scip, -minresactivity) )
907 {
908 *calculatedwcub = (lhs - minresactivity)/valdominating;
909 *wcubcalculated = TRUE;
910 }
911 else
912 {
913 /* worst case upper bound is -infinity */
914 *calculatedwcub = -SCIPinfinity(scip);
915 *wcubcalculated = TRUE;
916 }
917
918 /* we can only calculate lower bounds, of the right hand side is finite */
919 if( !SCIPmatrixIsRowRhsInfinity(matrix, row) )
920 {
921 /* lower bound calculation */
922 if( !SCIPisInfinity(scip, -minresactivity) )
923 {
924 *calculatedlb = (rhs - minresactivity)/valdominating;
925 *lbcalculated = TRUE;
926 }
927
928 /* worst case calculation of lower bound */
929 if( !SCIPisInfinity(scip, maxresactivity) )
930 {
931 *calculatedwclb = (rhs - maxresactivity)/valdominating;
932 *wclbcalculated = TRUE;
933 }
934 else
935 {
936 /* worst case lower bound is infinity */
937 *calculatedwclb = SCIPinfinity(scip);
938 *wclbcalculated = TRUE;
939 }
940 }
941 }
942
943 return SCIP_OKAY;
944 }
945
946 /** try to find new variable bounds and update them when they are better then the old bounds */
947 static
updateBounds(SCIP * scip,SCIP_MATRIX * matrix,int row,int col1,SCIP_Real val1,int col2,SCIP_Real val2,SCIP_Bool predictdominating,SCIP_Real * upperbound,SCIP_Real * wclowerbound,SCIP_Real * lowerbound,SCIP_Real * wcupperbound)948 SCIP_RETCODE updateBounds(
949 SCIP* scip, /**< SCIP main data structure */
950 SCIP_MATRIX* matrix, /**< matrix containing the constraints */
951 int row, /**< row index */
952 int col1, /**< dominating variable index */
953 SCIP_Real val1, /**< dominating variable coefficient */
954 int col2, /**< dominated variable index */
955 SCIP_Real val2, /**< dominated variable coefficient */
956 SCIP_Bool predictdominating, /**< flag indicating if bounds of dominating or dominated variable are predicted */
957 SCIP_Real* upperbound, /**< predicted upper bound */
958 SCIP_Real* wclowerbound, /**< predicted worst case lower bound */
959 SCIP_Real* lowerbound, /**< predicted lower bound */
960 SCIP_Real* wcupperbound /**< predicted worst case upper bound */
961 )
962 {
963 SCIP_Bool ubcalculated;
964 SCIP_Bool wclbcalculated;
965 SCIP_Bool lbcalculated;
966 SCIP_Bool wcubcalculated;
967 SCIP_Real newub;
968 SCIP_Real newwclb;
969 SCIP_Real newlb;
970 SCIP_Real newwcub;
971
972 assert(scip != NULL);
973 assert(matrix != NULL);
974 assert(upperbound != NULL);
975 assert(wclowerbound != NULL);
976 assert(lowerbound != NULL);
977 assert(wcupperbound != NULL);
978
979 newub = SCIPinfinity(scip);
980 newlb = -SCIPinfinity(scip);
981 newwcub = newub;
982 newwclb = newlb;
983
984 if( predictdominating )
985 {
986 /* do predictive rowbound analysis for the dominating variable */
987 SCIP_CALL( calcVarBoundsDominating(scip, matrix, row, col1, val1, col2, val2,
988 &ubcalculated, &newub, &wclbcalculated, &newwclb,
989 &lbcalculated, &newlb, &wcubcalculated, &newwcub) );
990 }
991 else
992 {
993 /* do predictive rowbound analysis for the dominated variable */
994 SCIP_CALL( calcVarBoundsDominated(scip, matrix, row, col1, val1, col2, val2,
995 &ubcalculated, &newub, &wclbcalculated, &newwclb,
996 &lbcalculated, &newlb, &wcubcalculated, &newwcub) );
997 }
998
999 /* update bounds in case if they are better */
1000 if( ubcalculated )
1001 {
1002 if( newub < *upperbound )
1003 *upperbound = newub;
1004 }
1005 if( wclbcalculated )
1006 {
1007 if( newwclb > *wclowerbound )
1008 *wclowerbound = newwclb;
1009 }
1010 if( lbcalculated )
1011 {
1012 if( newlb > *lowerbound )
1013 *lowerbound = newlb;
1014 }
1015 if( wcubcalculated )
1016 {
1017 if( newwcub < *wcupperbound )
1018 *wcupperbound = newwcub;
1019 }
1020
1021 return SCIP_OKAY;
1022 }
1023
1024 /** detect parallel columns by using the algorithm of Bixby and Wagner
1025 * see paper: "A note on Detecting Simple Redundancies in Linear Systems", June 1986
1026 */
1027 static
detectParallelCols(SCIP * scip,SCIP_MATRIX * matrix,int * pclass,SCIP_Bool * varineq)1028 SCIP_RETCODE detectParallelCols(
1029 SCIP* scip, /**< SCIP main data structure */
1030 SCIP_MATRIX* matrix, /**< matrix containing the constraints */
1031 int* pclass, /**< parallel column classes */
1032 SCIP_Bool* varineq /**< indicating if variable is within an equation */
1033 )
1034 {
1035 SCIP_Real* valpnt;
1036 SCIP_Real* values;
1037 SCIP_Real* scale;
1038 int* classsizes;
1039 int* pcset;
1040 int* rowpnt;
1041 int* rowend;
1042 int* colindices;
1043 int* pcs;
1044 SCIP_Real startval;
1045 SCIP_Real aij;
1046 int startpc;
1047 int startk;
1048 int startt;
1049 int pcsetfill;
1050 int colidx;
1051 int k;
1052 int t;
1053 int m;
1054 int i;
1055 int r;
1056 int newpclass;
1057 int pc;
1058 int nrows;
1059 int ncols;
1060
1061 assert(scip != NULL);
1062 assert(matrix != NULL);
1063 assert(pclass != NULL);
1064 assert(varineq != NULL);
1065
1066 nrows = SCIPmatrixGetNRows(matrix);
1067 ncols = SCIPmatrixGetNColumns(matrix);
1068
1069 SCIP_CALL( SCIPallocBufferArray(scip, &classsizes, ncols) );
1070 SCIP_CALL( SCIPallocBufferArray(scip, &scale, ncols) );
1071 SCIP_CALL( SCIPallocBufferArray(scip, &pcset, ncols) );
1072 SCIP_CALL( SCIPallocBufferArray(scip, &values, ncols) );
1073 SCIP_CALL( SCIPallocBufferArray(scip, &colindices, ncols) );
1074 SCIP_CALL( SCIPallocBufferArray(scip, &pcs, ncols) );
1075
1076 BMSclearMemoryArray(scale, ncols);
1077 BMSclearMemoryArray(pclass, ncols);
1078 BMSclearMemoryArray(classsizes, ncols);
1079
1080 classsizes[0] = ncols;
1081 pcsetfill = 0;
1082 for( t = 1; t < ncols; ++t )
1083 pcset[pcsetfill++] = t;
1084
1085 /* loop over all rows */
1086 for( r = 0; r < nrows; ++r )
1087 {
1088 /* we consider only non-empty equations or ranged rows */
1089 if( !SCIPmatrixIsRowRhsInfinity(matrix, r) && SCIPmatrixGetRowNNonzs(matrix, r) > 0 )
1090 {
1091 rowpnt = SCIPmatrixGetRowIdxPtr(matrix, r);
1092 rowend = rowpnt + SCIPmatrixGetRowNNonzs(matrix, r);
1093 valpnt = SCIPmatrixGetRowValPtr(matrix, r);
1094
1095 i = 0;
1096 for( ; (rowpnt < rowend); rowpnt++, valpnt++ )
1097 {
1098 aij = *valpnt;
1099 colidx = *rowpnt;
1100
1101 /* remember variable was part of an equation or ranged row */
1102 varineq[colidx] = TRUE;
1103
1104 if( scale[colidx] == 0.0 )
1105 scale[colidx] = aij;
1106 assert(scale[colidx] != 0.0);
1107
1108 colindices[i] = colidx;
1109 values[i] = aij / scale[colidx];
1110 pc = pclass[colidx];
1111 assert(pc < ncols);
1112
1113 /* update class sizes and pclass set */
1114 assert(classsizes[pc] > 0);
1115 classsizes[pc]--;
1116 if( classsizes[pc] == 0 )
1117 {
1118 assert(pcsetfill < ncols);
1119 pcset[pcsetfill++] = pc;
1120 }
1121 pcs[i] = pc;
1122
1123 i++;
1124 }
1125
1126 assert(i > 0);
1127
1128 /* sort on the pclass values */
1129 if( i > 1 )
1130 {
1131 SCIPsortIntIntReal(pcs, colindices, values, i);
1132 }
1133
1134 k = 0;
1135 while( TRUE ) /*lint !e716*/
1136 {
1137 assert(k < i);
1138 startpc = pcs[k];
1139 startk = k;
1140
1141 /* find pclass-sets */
1142 while( k < i && pcs[k] == startpc )
1143 k++;
1144
1145 /* sort on the A values which have equal pclass values */
1146 if( k - startk > 1 )
1147 SCIPsortRealInt(&(values[startk]), &(colindices[startk]), k - startk);
1148
1149 t = 0;
1150 while( TRUE ) /*lint !e716*/
1151 {
1152 assert(startk + t < i);
1153 startval = values[startk + t];
1154 startt = t;
1155
1156 /* find A-sets */
1157 while( t < k - startk && SCIPisEQ(scip, startval, values[startk + t]) )
1158 t++;
1159
1160 /* get new pclass */
1161 newpclass = pcset[0];
1162 assert(pcsetfill > 0);
1163 pcset[0] = pcset[--pcsetfill];
1164
1165 /* renumbering */
1166 for( m = startk + startt; m < startk + t; m++ )
1167 {
1168 assert(m < i);
1169 assert(colindices[m] < ncols);
1170 assert(newpclass < ncols);
1171
1172 pclass[colindices[m]] = newpclass;
1173 classsizes[newpclass]++;
1174 }
1175
1176 if( t == k - startk )
1177 break;
1178 }
1179
1180 if( k == SCIPmatrixGetRowNNonzs(matrix, r) )
1181 break;
1182 }
1183 }
1184 }
1185
1186 SCIPfreeBufferArray(scip, &pcs);
1187 SCIPfreeBufferArray(scip, &colindices);
1188 SCIPfreeBufferArray(scip, &values);
1189 SCIPfreeBufferArray(scip, &pcset);
1190 SCIPfreeBufferArray(scip, &scale);
1191 SCIPfreeBufferArray(scip, &classsizes);
1192
1193 return SCIP_OKAY;
1194 }
1195
1196
1197 /** try to improve variable bounds by predictive bound strengthening */
1198 static
predBndStr(SCIP * scip,SCIP_VAR * dominatingvar,int dominatingidx,SCIP_Real dominatingub,SCIP_Real dominatinglb,SCIP_Real dominatingwcub,SCIP_VAR * dominatedvar,int dominatedidx,SCIP_Real dominatedub,SCIP_Real dominatedwclb,SCIP_Real dominatedlb,FIXINGDIRECTION * varstofix,int * nchgbds)1199 SCIP_RETCODE predBndStr(
1200 SCIP* scip, /**< SCIP main data structure */
1201 SCIP_VAR* dominatingvar, /**< dominating variable */
1202 int dominatingidx, /**< column index of the dominating variable */
1203 SCIP_Real dominatingub, /**< predicted upper bound of the dominating variable */
1204 SCIP_Real dominatinglb, /**< predicted lower bound of the dominating variable */
1205 SCIP_Real dominatingwcub, /**< predicted worst case upper bound of the dominating variable */
1206 SCIP_VAR* dominatedvar, /**< dominated variable */
1207 int dominatedidx, /**< column index of the dominated variable */
1208 SCIP_Real dominatedub, /**< predicted upper bound of the dominated variable */
1209 SCIP_Real dominatedwclb, /**< predicted worst case lower bound of the dominated variable */
1210 SCIP_Real dominatedlb, /**< predicted lower bound of the dominated variable */
1211 FIXINGDIRECTION* varstofix, /**< array holding fixing information */
1212 int* nchgbds /**< count number of bound changes */
1213 )
1214 {
1215 /* we compare only variables from the same type */
1216 if( !(SCIPvarGetType(dominatingvar) == SCIPvarGetType(dominatedvar) ||
1217 SCIPvarIsBinary(dominatingvar) == SCIPvarIsBinary(dominatedvar) ||
1218 (SCIPvarGetType(dominatingvar) == SCIP_VARTYPE_INTEGER && SCIPvarGetType(dominatedvar) == SCIP_VARTYPE_IMPLINT) ||
1219 (SCIPvarGetType(dominatedvar) == SCIP_VARTYPE_INTEGER && SCIPvarGetType(dominatingvar) == SCIP_VARTYPE_IMPLINT)) )
1220 {
1221 return SCIP_OKAY;
1222 }
1223
1224 if( varstofix[dominatingidx] == NOFIX )
1225 {
1226 /* assume x dominates y (x->y). we get this bound from a positive coefficient
1227 * of the dominating variable. because x->y the dominated variable y has
1228 * a positive coefficient too. thus y contributes to the minactivity with its
1229 * lower bound. but this case is considered within predictive bound analysis.
1230 * thus the dominating upper bound is a common upper bound.
1231 */
1232 if( !SCIPisInfinity(scip, dominatingub) )
1233 {
1234 SCIP_Real newub;
1235 SCIP_Real oldub;
1236 SCIP_Real lb;
1237
1238 newub = dominatingub;
1239 oldub = SCIPvarGetUbGlobal(dominatingvar);
1240 lb = SCIPvarGetLbGlobal(dominatingvar);
1241
1242 /* if( SCIPvarGetType(dominatingvar) != SCIP_VARTYPE_CONTINUOUS )
1243 newub = SCIPceil(scip, newub); */
1244
1245 if( SCIPisLE(scip, lb, newub) && SCIPisLT(scip, newub, oldub) )
1246 {
1247 SCIPdebugMsg(scip, "[ub]\tupper bound for dominating variable <%s> changed: [%.17f,%.17f] -> [%.17f,%.17f]\n",
1248 SCIPvarGetName(dominatingvar), lb, oldub, lb, newub);
1249 SCIP_CALL( SCIPchgVarUb(scip, dominatingvar, newub) );
1250 (*nchgbds)++;
1251 }
1252 }
1253
1254 /* assume x dominates y (x->y). we get this lower bound of the dominating variable
1255 * from a negative coefficient within a <= relation. if y has a positive coefficient
1256 * we get a common lower bound of x from predictive bound analysis. if y has a
1257 * negative coefficient we get an improved lower bound of x because the minactivity
1258 * is greater. for discrete variables we have to round down the lower bound.
1259 */
1260 if( !SCIPisInfinity(scip, -dominatinglb) )
1261 {
1262 SCIP_Real newlb;
1263 SCIP_Real oldlb;
1264 SCIP_Real ub;
1265
1266 newlb = dominatinglb;
1267 oldlb = SCIPvarGetLbGlobal(dominatingvar);
1268 ub = SCIPvarGetUbGlobal(dominatingvar);
1269
1270 if( SCIPvarGetType(dominatingvar) != SCIP_VARTYPE_CONTINUOUS )
1271 newlb = SCIPfloor(scip, newlb);
1272
1273 if( SCIPisLT(scip, oldlb, newlb) && SCIPisLE(scip, newlb, ub) )
1274 {
1275 SCIPdebugMsg(scip, "[lb]\tlower bound of dominating variable <%s> changed: [%.17f,%.17f] -> [%.17f,%.17f]\n",
1276 SCIPvarGetName(dominatingvar), oldlb, ub, newlb, ub);
1277 SCIP_CALL( SCIPchgVarLb(scip, dominatingvar, newlb) );
1278 (*nchgbds)++;
1279 }
1280 }
1281
1282 /* assume x dominates y (x->y). we get this bound from a positive coefficient
1283 * of x within a <= relation. from x->y it follows, that y has a positive
1284 * coefficient in this row too. the worst case upper bound of x is an estimation
1285 * of how great x can be in every case. if the objective coefficient of x is
1286 * negative we get thus a lower bound of x. for discrete variables we have
1287 * to round down the lower bound before.
1288 */
1289 if( !SCIPisInfinity(scip, dominatingwcub) && SCIPisNegative(scip, SCIPvarGetObj(dominatingvar)))
1290 {
1291 SCIP_Real newlb;
1292 SCIP_Real oldlb;
1293 SCIP_Real ub;
1294
1295 newlb = dominatingwcub;
1296 oldlb = SCIPvarGetLbGlobal(dominatingvar);
1297 ub = SCIPvarGetUbGlobal(dominatingvar);
1298
1299 if( SCIPvarGetType(dominatingvar) != SCIP_VARTYPE_CONTINUOUS )
1300 newlb = SCIPfloor(scip, newlb);
1301
1302 if( SCIPisLT(scip, oldlb, newlb) && SCIPisLE(scip, newlb, ub) )
1303 {
1304 SCIPdebugMsg(scip, "[wcub]\tlower bound of dominating variable <%s> changed: [%.17f,%.17f] -> [%.17f,%.17f]\n",
1305 SCIPvarGetName(dominatingvar), oldlb, ub, newlb, ub);
1306 SCIP_CALL( SCIPchgVarLb(scip, dominatingvar, newlb) );
1307 (*nchgbds)++;
1308 }
1309 }
1310 }
1311
1312 if( varstofix[dominatedidx] == NOFIX )
1313 {
1314 /* assume x dominates y (x->y). we get this bound for a positive coefficient of y
1315 * within a <= relation. if x has a negative coefficient we get a common upper
1316 * bound of y. if x has a positive coefficient we get an improved upper bound
1317 * of y because the minactivity is greater.
1318 */
1319 if( !SCIPisInfinity(scip, dominatedub) )
1320 {
1321 SCIP_Real newub;
1322 SCIP_Real oldub;
1323 SCIP_Real lb;
1324
1325 newub = dominatedub;
1326 oldub = SCIPvarGetUbGlobal(dominatedvar);
1327 lb = SCIPvarGetLbGlobal(dominatedvar);
1328
1329 if( SCIPisLE(scip, lb, newub) && SCIPisLT(scip, newub, oldub) )
1330 {
1331 SCIPdebugMsg(scip, "[ub]\tupper bound of dominated variable <%s> changed: [%.17f,%.17f] -> [%.17f,%.17f]\n",
1332 SCIPvarGetName(dominatedvar), lb, oldub, lb, newub);
1333 SCIP_CALL( SCIPchgVarUb(scip, dominatedvar, newub) );
1334 (*nchgbds)++;
1335 }
1336 }
1337
1338 /* assume x dominates y (x->y). we get this bound only from a negative
1339 * coefficient of y within a <= relation. because of x->y then x has a negative
1340 * coefficient too. the worst case lower bound is an estimation what property
1341 * the dominated variable must have if the dominating variable is at its upper bound.
1342 * to get an upper bound of the dominated variable we need to consider a positive
1343 * objective coefficient. for discrete variables we have to round up the upper bound.
1344 */
1345 if( !SCIPisInfinity(scip, -dominatedwclb) && SCIPisPositive(scip, SCIPvarGetObj(dominatedvar)) )
1346 {
1347 SCIP_Real newub;
1348 SCIP_Real oldub;
1349 SCIP_Real lb;
1350
1351 newub = dominatedwclb;
1352 oldub = SCIPvarGetUbGlobal(dominatedvar);
1353 lb = SCIPvarGetLbGlobal(dominatedvar);
1354
1355 if( SCIPvarGetType(dominatedvar) != SCIP_VARTYPE_CONTINUOUS )
1356 newub = SCIPceil(scip, newub);
1357
1358 if( SCIPisLE(scip, lb, newub) && SCIPisLT(scip, newub, oldub) )
1359 {
1360 SCIPdebugMsg(scip, "[wclb]\tupper bound of dominated variable <%s> changed: [%.17f,%.17f] -> [%.17f,%.17f]\n",
1361 SCIPvarGetName(dominatedvar), lb, oldub, lb, newub);
1362 SCIP_CALL( SCIPchgVarUb(scip, dominatedvar, newub) );
1363 (*nchgbds)++;
1364 }
1365 }
1366
1367 /* assume x dominates y (x->y). we get a lower bound of y from a negative
1368 * coefficient within a <= relation. but if x->y then x has a negative
1369 * coefficient too and contributes with its upper bound to the minactivity.
1370 * thus in all we have a common lower bound calculation and no rounding is necessary.
1371 */
1372 if( !SCIPisInfinity(scip, -dominatedlb) )
1373 {
1374 SCIP_Real newlb;
1375 SCIP_Real oldlb;
1376 SCIP_Real ub;
1377
1378 newlb = dominatedlb;
1379 oldlb = SCIPvarGetLbGlobal(dominatedvar);
1380 ub = SCIPvarGetUbGlobal(dominatedvar);
1381
1382 if( SCIPisLT(scip, oldlb, newlb) && SCIPisLE(scip, newlb, ub) )
1383 {
1384 SCIPdebugMsg(scip, "[lb]\tlower bound of dominated variable <%s> changed: [%.17f,%.17f] -> [%.17f,%.17f]\n",
1385 SCIPvarGetName(dominatedvar), oldlb, ub, newlb, ub);
1386 SCIP_CALL( SCIPchgVarLb(scip, dominatedvar, newlb) );
1387 (*nchgbds)++;
1388 }
1389 }
1390 }
1391
1392 return SCIP_OKAY;
1393 }
1394
1395 /** try to find variable fixings */
1396 static
findFixings(SCIP * scip,SCIP_MATRIX * matrix,SCIP_VAR * dominatingvar,int dominatingidx,SCIP_Real dominatingub,SCIP_Real dominatingwclb,SCIP_Real dominatinglb,SCIP_Real dominatingwcub,SCIP_VAR * dominatedvar,int dominatedidx,FIXINGDIRECTION * varstofix,SCIP_Bool onlybinvars,SCIP_Bool onlyoneone,int * nfixings)1397 SCIP_RETCODE findFixings(
1398 SCIP* scip, /**< SCIP main data structure */
1399 SCIP_MATRIX* matrix, /**< constraint matrix structure */
1400 SCIP_VAR* dominatingvar, /**< dominating variable */
1401 int dominatingidx, /**< column index of the dominating variable */
1402 SCIP_Real dominatingub, /**< predicted upper bound of the dominating variable */
1403 SCIP_Real dominatingwclb, /**< predicted worst case lower bound of the dominating variable */
1404 SCIP_Real dominatinglb, /**< predicted lower bound of the dominating variable */
1405 SCIP_Real dominatingwcub, /**< predicted worst case upper bound of the dominating variable */
1406 SCIP_VAR* dominatedvar, /**< dominated variable */
1407 int dominatedidx, /**< column index of the dominated variable */
1408 FIXINGDIRECTION* varstofix, /**< array holding fixing information */
1409 SCIP_Bool onlybinvars, /**< flag indicating only binary variables are present */
1410 SCIP_Bool onlyoneone, /**< when onlybinvars is TRUE, flag indicates if both binary variables are in clique */
1411 int* nfixings /**< counter for possible fixings */
1412 )
1413 {
1414 /* we compare only variables from the same type */
1415 if( !(SCIPvarGetType(dominatingvar) == SCIPvarGetType(dominatedvar) ||
1416 SCIPvarIsBinary(dominatingvar) == SCIPvarIsBinary(dominatedvar) ||
1417 (SCIPvarGetType(dominatingvar) == SCIP_VARTYPE_INTEGER && SCIPvarGetType(dominatedvar) == SCIP_VARTYPE_IMPLINT) ||
1418 (SCIPvarGetType(dominatedvar) == SCIP_VARTYPE_INTEGER && SCIPvarGetType(dominatingvar) == SCIP_VARTYPE_IMPLINT)) )
1419 {
1420 return SCIP_OKAY;
1421 }
1422
1423 if( varstofix[dominatedidx] == NOFIX && SCIPmatrixGetColNNonzs(matrix, dominatingidx) == 1
1424 && SCIPmatrixGetColNNonzs(matrix, dominatedidx) == 1 )
1425 {
1426 /* We have a x->y dominance relation and only one equality constraint
1427 * where the dominating variable x with an infinity upper bound and the
1428 * dominated variable y are present. Then the dominated variable y
1429 * can be fixed at its lower bound.
1430 */
1431 int row;
1432 row = *(SCIPmatrixGetColIdxPtr(matrix, dominatedidx));
1433
1434 if( SCIPisEQ(scip, SCIPmatrixGetRowLhs(matrix, row), SCIPmatrixGetRowRhs(matrix, row)) &&
1435 SCIPisInfinity(scip, SCIPvarGetUbGlobal(dominatingvar)) )
1436 {
1437 varstofix[dominatedidx] = FIXATLB;
1438 (*nfixings)++;
1439
1440 return SCIP_OKAY;
1441 }
1442 }
1443
1444 if( varstofix[dominatedidx] == NOFIX && !SCIPisNegative(scip, SCIPvarGetObj(dominatedvar)) )
1445 {
1446 if( !SCIPisInfinity(scip, -dominatingwclb) &&
1447 SCIPisLE(scip, dominatingwclb, SCIPvarGetUbGlobal(dominatingvar)) )
1448 {
1449 /* we have a x->y dominance relation with a positive obj coefficient
1450 * of the dominated variable y. we need to secure feasibility
1451 * by testing if the predicted lower worst case bound is less equal the
1452 * current upper bound. it is possible, that the lower worst case bound
1453 * is infinity and the upper bound of the dominating variable x is
1454 * infinity too.
1455 */
1456 varstofix[dominatedidx] = FIXATLB;
1457 (*nfixings)++;
1458 }
1459 }
1460
1461 if( varstofix[dominatedidx] == NOFIX && !SCIPisInfinity(scip, dominatingub) &&
1462 SCIPisLE(scip, dominatingub, SCIPvarGetUbGlobal(dominatingvar)) )
1463 {
1464 /* we have a x->y dominance relation with an arbitrary obj coefficient
1465 * of the dominating variable x. in all cases we have to look
1466 * if the predicted upper bound of the dominating variable is great enough.
1467 * by testing, that the predicted upper bound is not infinity we avoid problems
1468 * with x->y e.g.
1469 * min -x -y
1470 * s.t. -x -y <= -1
1471 * 0<=x<=1, 0<=y<=1
1472 * where y is not at their lower bound.
1473 */
1474 varstofix[dominatedidx] = FIXATLB;
1475 (*nfixings)++;
1476 }
1477
1478 if( varstofix[dominatingidx] == NOFIX && !SCIPisPositive(scip, SCIPvarGetObj(dominatingvar)) )
1479 {
1480 /* we have a x->y dominance relation with a negative obj coefficient
1481 * of the dominating variable x. if the worst case upper bound is
1482 * greater equal than upper bound, we fix x at the upper bound
1483 */
1484 if( !SCIPisInfinity(scip, dominatingwcub) &&
1485 SCIPisGE(scip, dominatingwcub, SCIPvarGetUbGlobal(dominatingvar)) )
1486 {
1487 varstofix[dominatingidx] = FIXATUB;
1488 (*nfixings)++;
1489 }
1490 }
1491
1492 if( varstofix[dominatingidx] == NOFIX && !SCIPisInfinity(scip, -dominatinglb) &&
1493 SCIPisGE(scip, dominatinglb, SCIPvarGetUbGlobal(dominatingvar)) )
1494 {
1495 /* we have a x->y dominance relation with an arbitrary obj coefficient
1496 * of the dominating variable x. if the predicted lower bound is greater
1497 * equal than upper bound, we fix x at the upper bound.
1498 */
1499 varstofix[dominatingidx] = FIXATUB;
1500 (*nfixings)++;
1501 }
1502
1503 if( onlybinvars )
1504 {
1505 if( varstofix[dominatedidx] == NOFIX && (onlyoneone || SCIPvarsHaveCommonClique(dominatingvar, TRUE, dominatedvar, TRUE, TRUE)) )
1506 {
1507 /* We have a (1->1)-clique with dominance relation (x->y) (x dominates y).
1508 * From this dominance relation, we know (1->0) is possible and not worse than (0->1)
1509 * concerning the objective function. It follows that only (1->0) or (0->0) are possible,
1510 * but in both cases y has the value 0 => y=0.
1511 */
1512 varstofix[dominatedidx] = FIXATLB;
1513 (*nfixings)++;
1514 }
1515
1516 if( varstofix[dominatingidx] == NOFIX && SCIPvarsHaveCommonClique(dominatingvar, FALSE, dominatedvar, FALSE, TRUE) )
1517 {
1518 /* We have a (0->0)-clique with dominance relation x->y (x dominates y).
1519 * From this dominance relation, we know (1->0) is possible and not worse than (0->1)
1520 * concerning the objective function. It follows that only (1->0) or (1->1) are possible,
1521 * but in both cases x has the value 1 => x=1
1522 */
1523 varstofix[dominatingidx] = FIXATUB;
1524 (*nfixings)++;
1525 }
1526 }
1527 else
1528 assert(!onlyoneone);
1529
1530 return SCIP_OKAY;
1531 }
1532
1533 /** find dominance relation between variable pairs */
1534 static
findDominancePairs(SCIP * scip,SCIP_MATRIX * matrix,SCIP_PRESOLDATA * presoldata,int * searchcols,int searchsize,SCIP_Bool onlybinvars,FIXINGDIRECTION * varstofix,int * nfixings,SCIP_Longint * ndomrelations,int * nchgbds)1535 SCIP_RETCODE findDominancePairs(
1536 SCIP* scip, /**< SCIP main data structure */
1537 SCIP_MATRIX* matrix, /**< matrix containing the constraints */
1538 SCIP_PRESOLDATA* presoldata, /**< presolver data */
1539 int* searchcols, /**< indexes of variables for pair comparisons */
1540 int searchsize, /**< number of variables for pair comparisons */
1541 SCIP_Bool onlybinvars, /**< flag indicating searchcols contains only binary variable indexes */
1542 FIXINGDIRECTION* varstofix, /**< array holding information for later upper/lower bound fixing */
1543 int* nfixings, /**< found number of possible fixings */
1544 SCIP_Longint* ndomrelations, /**< found number of dominance relations */
1545 int* nchgbds /**< number of changed bounds */
1546 )
1547 {
1548 SCIP_Real* vals1;
1549 SCIP_Real* vals2;
1550 SCIP_Real tmpupperbounddominatingcol1;
1551 SCIP_Real tmpupperbounddominatingcol2;
1552 SCIP_Real tmpwclowerbounddominatingcol1;
1553 SCIP_Real tmpwclowerbounddominatingcol2;
1554 SCIP_Real tmplowerbounddominatingcol1;
1555 SCIP_Real tmplowerbounddominatingcol2;
1556 SCIP_Real tmpwcupperbounddominatingcol1;
1557 SCIP_Real tmpwcupperbounddominatingcol2;
1558 int* rows1;
1559 int* rows2;
1560 int nrows1;
1561 int nrows2;
1562 SCIP_Real tmpupperbounddominatedcol1;
1563 SCIP_Real tmpupperbounddominatedcol2;
1564 SCIP_Real tmpwclowerbounddominatedcol1;
1565 SCIP_Real tmpwclowerbounddominatedcol2;
1566 SCIP_Real tmplowerbounddominatedcol1;
1567 SCIP_Real tmplowerbounddominatedcol2;
1568 SCIP_Real tmpwcupperbounddominatedcol1;
1569 SCIP_Real tmpwcupperbounddominatedcol2;
1570 SCIP_Real obj1;
1571 SCIP_Bool col1domcol2;
1572 SCIP_Bool col2domcol1;
1573 SCIP_Bool onlyoneone;
1574 int cnt1;
1575 int cnt2;
1576 int col1;
1577 int col2;
1578 int r1;
1579 int r2;
1580 int paircnt;
1581 int oldnfixings;
1582
1583 assert(scip != NULL);
1584 assert(matrix != NULL);
1585 assert(presoldata != NULL);
1586 assert(searchcols != NULL);
1587 assert(varstofix != NULL);
1588 assert(nfixings != NULL);
1589 assert(ndomrelations != NULL);
1590 assert(nchgbds != NULL);
1591
1592 paircnt = 0;
1593 oldnfixings = *nfixings;
1594
1595 /* pair comparisons */
1596 for( cnt1 = 0; cnt1 < searchsize; cnt1++ )
1597 {
1598 SCIP_VAR* varcol1;
1599 SCIP_VAR* varcol2;
1600
1601 /* get index of the first variable */
1602 col1 = searchcols[cnt1];
1603
1604 if( varstofix[col1] == FIXATLB )
1605 continue;
1606
1607 varcol1 = SCIPmatrixGetVar(matrix, col1);
1608 obj1 = SCIPvarGetObj(varcol1);
1609
1610 for( cnt2 = cnt1+1; cnt2 < searchsize; cnt2++ )
1611 {
1612 /* get index of the second variable */
1613 col2 = searchcols[cnt2];
1614 varcol2 = SCIPmatrixGetVar(matrix, col2);
1615 onlyoneone = FALSE;
1616
1617 /* we always have minimize as obj sense */
1618
1619 /* column 1 dominating column 2 */
1620 col1domcol2 = (obj1 <= SCIPvarGetObj(varcol2));
1621
1622 /* column 2 dominating column 1 */
1623 col2domcol1 = (SCIPvarGetObj(varcol2) <= obj1);
1624
1625 /* search only if nothing was found yet */
1626 col1domcol2 = col1domcol2 && (varstofix[col2] == NOFIX);
1627 col2domcol1 = col2domcol1 && (varstofix[col1] == NOFIX);
1628
1629 /* we only search for a dominance relation if the lower bounds are not negative */
1630 if( !onlybinvars )
1631 {
1632 if( SCIPisLT(scip, SCIPvarGetLbGlobal(varcol1), 0.0) ||
1633 SCIPisLT(scip, SCIPvarGetLbGlobal(varcol2), 0.0) )
1634 {
1635 col1domcol2 = FALSE;
1636 col2domcol1 = FALSE;
1637 }
1638 }
1639
1640 /* pair comparison control */
1641 if( paircnt == presoldata->numcurrentpairs )
1642 {
1643 assert(*nfixings >= oldnfixings);
1644 if( *nfixings == oldnfixings )
1645 {
1646 /* not enough fixings found, decrement number of comparisons */
1647 presoldata->numcurrentpairs >>= 1; /*lint !e702*/
1648 if( presoldata->numcurrentpairs < presoldata->numminpairs )
1649 presoldata->numcurrentpairs = presoldata->numminpairs;
1650
1651 /* stop searching in this row */
1652 return SCIP_OKAY;
1653 }
1654 oldnfixings = *nfixings;
1655 paircnt = 0;
1656
1657 /* increment number of comparisons */
1658 presoldata->numcurrentpairs <<= 1; /*lint !e701*/
1659 if( presoldata->numcurrentpairs > presoldata->nummaxpairs )
1660 presoldata->numcurrentpairs = presoldata->nummaxpairs;
1661 }
1662 paircnt++;
1663
1664 if( !col1domcol2 && !col2domcol1 )
1665 continue;
1666
1667 /* get the data for both columns */
1668 vals1 = SCIPmatrixGetColValPtr(matrix, col1);
1669 rows1 = SCIPmatrixGetColIdxPtr(matrix, col1);
1670 nrows1 = SCIPmatrixGetColNNonzs(matrix, col1);
1671 r1 = 0;
1672 vals2 = SCIPmatrixGetColValPtr(matrix, col2);
1673 rows2 = SCIPmatrixGetColIdxPtr(matrix, col2);
1674 nrows2 = SCIPmatrixGetColNNonzs(matrix, col2);
1675 r2 = 0;
1676
1677 /* do we have a obj constant? */
1678 if( nrows1 == 0 || nrows2 == 0 )
1679 continue;
1680
1681 /* initialize temporary bounds of dominating variable */
1682 tmpupperbounddominatingcol1 = SCIPinfinity(scip);
1683 tmpupperbounddominatingcol2 = tmpupperbounddominatingcol1;
1684 tmpwclowerbounddominatingcol1 = -SCIPinfinity(scip);
1685 tmpwclowerbounddominatingcol2 = tmpwclowerbounddominatingcol1;
1686 tmplowerbounddominatingcol1 = -SCIPinfinity(scip);
1687 tmplowerbounddominatingcol2 = tmplowerbounddominatingcol1;
1688 tmpwcupperbounddominatingcol1 = SCIPinfinity(scip);
1689 tmpwcupperbounddominatingcol2 = tmpwcupperbounddominatingcol1;
1690
1691 /* initialize temporary bounds of dominated variable */
1692 tmpupperbounddominatedcol1 = SCIPinfinity(scip);
1693 tmpupperbounddominatedcol2 = tmpupperbounddominatedcol1;
1694 tmpwclowerbounddominatedcol1 = -SCIPinfinity(scip);
1695 tmpwclowerbounddominatedcol2 = tmpwclowerbounddominatedcol1;
1696 tmplowerbounddominatedcol1 = -SCIPinfinity(scip);
1697 tmplowerbounddominatedcol2 = tmplowerbounddominatedcol1;
1698 tmpwcupperbounddominatedcol1 = SCIPinfinity(scip);
1699 tmpwcupperbounddominatedcol2 = tmpwcupperbounddominatedcol1;
1700
1701 /* compare rows of this column pair */
1702 while( (col1domcol2 || col2domcol1) && (r1 < nrows1 || r2 < nrows2) )
1703 {
1704 assert((r1 >= nrows1-1) || (rows1[r1] < rows1[r1+1]));
1705 assert((r2 >= nrows2-1) || (rows2[r2] < rows2[r2+1]));
1706
1707 /* there is a nonredundant row containing column 1 but not column 2 */
1708 if( r1 < nrows1 && (r2 == nrows2 || rows1[r1] < rows2[r2]) )
1709 {
1710 /* dominance depends on the type of relation */
1711 if( !SCIPmatrixIsRowRhsInfinity(matrix, rows1[r1]) )
1712 {
1713 /* no dominance relation for equations or ranged rows */
1714 col2domcol1 = FALSE;
1715 col1domcol2 = FALSE;
1716 }
1717 else
1718 {
1719 /* >= relation, larger coefficients dominate smaller ones */
1720 if( vals1[r1] > 0.0 )
1721 col2domcol1 = FALSE;
1722 else
1723 col1domcol2 = FALSE;
1724 }
1725
1726 r1++;
1727 }
1728 /* there is a nonredundant row containing column 2, but not column 1 */
1729 else if( r2 < nrows2 && (r1 == nrows1 || rows1[r1] > rows2[r2]) )
1730 {
1731 /* dominance depends on the type of relation */
1732 if( !SCIPmatrixIsRowRhsInfinity(matrix, rows2[r2]) )
1733 {
1734 /* no dominance relation for equations or ranged rows */
1735 col2domcol1 = FALSE;
1736 col1domcol2 = FALSE;
1737 }
1738 else
1739 {
1740 /* >= relation, larger coefficients dominate smaller ones */
1741 if( vals2[r2] < 0.0 )
1742 col2domcol1 = FALSE;
1743 else
1744 col1domcol2 = FALSE;
1745 }
1746
1747 r2++;
1748 }
1749 /* if both columns appear in a common row, compare the coefficients */
1750 else
1751 {
1752 assert(r1 < nrows1 && r2 < nrows2);
1753 assert(rows1[r1] == rows2[r2]);
1754
1755 /* if both columns are binary variables we check if they have a common clique
1756 * and do not calculate any bounds
1757 */
1758 if( onlybinvars && !onlyoneone )
1759 {
1760 if( vals1[r1] < 0 && vals2[r2] < 0 )
1761 {
1762 if( (SCIPmatrixGetRowNMaxActPosInf(matrix, rows1[r1]) + SCIPmatrixGetRowNMaxActNegInf(matrix, rows1[r1]) == 0)
1763 && SCIPisFeasLE(scip, SCIPmatrixGetRowMaxActivity(matrix, rows1[r1]) + MAX(vals1[r1], vals2[r2]),
1764 SCIPmatrixGetRowLhs(matrix, rows1[r1])) )
1765 {
1766 onlyoneone = TRUE;
1767 }
1768 }
1769
1770 if( !onlyoneone && !SCIPmatrixIsRowRhsInfinity(matrix, rows1[r1]) )
1771 {
1772 if ( vals1[r1] > 0 && vals2[r2] > 0 )
1773 {
1774 if( (SCIPmatrixGetRowNMinActPosInf(matrix, rows1[r1]) + SCIPmatrixGetRowNMinActNegInf(matrix, rows1[r1]) == 0)
1775 && SCIPisFeasGE(scip, SCIPmatrixGetRowMinActivity(matrix, rows1[r1]) + MIN(vals1[r1], vals2[r2]),
1776 SCIPmatrixGetRowRhs(matrix, rows1[r1])) )
1777 {
1778 onlyoneone = TRUE;
1779 }
1780 }
1781 }
1782
1783 if( onlyoneone )
1784 {
1785 /* reset bounds */
1786 tmpupperbounddominatingcol1 = SCIPinfinity(scip);
1787 tmpupperbounddominatingcol2 = tmpupperbounddominatingcol1;
1788 tmpwclowerbounddominatingcol1 = -SCIPinfinity(scip);
1789 tmpwclowerbounddominatingcol2 = tmpwclowerbounddominatingcol1;
1790 tmplowerbounddominatingcol1 = -SCIPinfinity(scip);
1791 tmplowerbounddominatingcol2 = tmplowerbounddominatingcol1;
1792 tmpwcupperbounddominatingcol1 = SCIPinfinity(scip);
1793 tmpwcupperbounddominatingcol2 = tmpwcupperbounddominatingcol1;
1794
1795 tmpupperbounddominatedcol1 = SCIPinfinity(scip);
1796 tmpupperbounddominatedcol2 = tmpupperbounddominatedcol1;
1797 tmpwclowerbounddominatedcol1 = -SCIPinfinity(scip);
1798 tmpwclowerbounddominatedcol2 = tmpwclowerbounddominatedcol1;
1799 tmplowerbounddominatedcol1 = -SCIPinfinity(scip);
1800 tmplowerbounddominatedcol2 = tmplowerbounddominatedcol1;
1801 tmpwcupperbounddominatedcol1 = SCIPinfinity(scip);
1802 tmpwcupperbounddominatedcol2 = tmpwcupperbounddominatedcol1;
1803 }
1804 }
1805
1806 /* dominance depends on the type of inequality */
1807 if( !SCIPmatrixIsRowRhsInfinity(matrix, rows1[r1]) )
1808 {
1809 /* no dominance relation if coefficients differ for equations or ranged rows */
1810 if( !SCIPisEQ(scip, vals1[r1], vals2[r2]) )
1811 {
1812 col2domcol1 = FALSE;
1813 col1domcol2 = FALSE;
1814 }
1815 }
1816 else
1817 {
1818 /* >= relation, larger coefficients dominate smaller ones */
1819 if( vals1[r1] > vals2[r2] )
1820 col2domcol1 = FALSE;
1821 else if( vals1[r1] < vals2[r2] )
1822 col1domcol2 = FALSE;
1823 }
1824
1825 /* we do not use bound calulations if two binary variable are in one common clique.
1826 * for the other cases we claim the same sign for the coefficients to
1827 * achieve monotonically decreasing predictive bound functions.
1828 */
1829 if( !onlyoneone &&
1830 ((vals1[r1] < 0 && vals2[r2] < 0) || (vals1[r1] > 0 && vals2[r2] > 0)) )
1831 {
1832 if( col1domcol2 )
1833 {
1834 /* update bounds of dominating variable for column 1 */
1835 SCIP_CALL( updateBounds(scip, matrix, rows1[r1],
1836 col1, vals1[r1], col2, vals2[r2], TRUE,
1837 &tmpupperbounddominatingcol1, &tmpwclowerbounddominatingcol1,
1838 &tmplowerbounddominatingcol1, &tmpwcupperbounddominatingcol1) );
1839
1840 /* update bounds of dominated variable for column 1 */
1841 SCIP_CALL( updateBounds(scip, matrix, rows1[r1],
1842 col1, vals1[r1], col2, vals2[r2], FALSE,
1843 &tmpupperbounddominatedcol1, &tmpwclowerbounddominatedcol1,
1844 &tmplowerbounddominatedcol1, &tmpwcupperbounddominatedcol1) );
1845 }
1846
1847 if( col2domcol1 )
1848 {
1849 /* update bounds of dominating variable for column 2 */
1850 SCIP_CALL( updateBounds(scip, matrix, rows2[r2],
1851 col2, vals2[r2], col1, vals1[r1], TRUE,
1852 &tmpupperbounddominatingcol2, &tmpwclowerbounddominatingcol2,
1853 &tmplowerbounddominatingcol2, &tmpwcupperbounddominatingcol2) );
1854
1855 /* update bounds of dominated variable for column 2 */
1856 SCIP_CALL( updateBounds(scip, matrix, rows2[r2],
1857 col2, vals2[r2], col1, vals1[r1], FALSE,
1858 &tmpupperbounddominatedcol2, &tmpwclowerbounddominatedcol2,
1859 &tmplowerbounddominatedcol2, &tmpwcupperbounddominatedcol2) );
1860 }
1861 }
1862
1863 r1++;
1864 r2++;
1865 }
1866 }
1867
1868 /* a column is only dominated, if there are no more rows in which it is contained */
1869 col1domcol2 = col1domcol2 && r2 == nrows2;
1870 col2domcol1 = col2domcol1 && r1 == nrows1;
1871
1872 if( !col1domcol2 && !col2domcol1 )
1873 continue;
1874
1875 /* no dominance relation for left equations or ranged rows */
1876 while( r1 < nrows1 )
1877 {
1878 if( !SCIPmatrixIsRowRhsInfinity(matrix, rows1[r1]) )
1879 {
1880 col2domcol1 = FALSE;
1881 col1domcol2 = FALSE;
1882 break;
1883 }
1884 r1++;
1885 }
1886 if( !col1domcol2 && !col2domcol1 )
1887 continue;
1888 while( r2 < nrows2 )
1889 {
1890 if( !SCIPmatrixIsRowRhsInfinity(matrix, rows2[r2]) )
1891 {
1892 col2domcol1 = FALSE;
1893 col1domcol2 = FALSE;
1894 break;
1895 }
1896 r2++;
1897 }
1898
1899 if( col1domcol2 || col2domcol1 )
1900 (*ndomrelations)++;
1901
1902 if( col1domcol2 && col2domcol1 )
1903 {
1904 /* prefer the variable as dominating variable with the greater upper bound */
1905 if( SCIPisGE(scip, SCIPvarGetUbGlobal(varcol1), SCIPvarGetUbGlobal(varcol2)) )
1906 col2domcol1 = FALSE;
1907 else
1908 col1domcol2 = FALSE;
1909 }
1910
1911 /* use dominance relation and clique/bound-information
1912 * to find variable fixings. additionally try to strengthen
1913 * variable bounds by predictive bound strengthening.
1914 */
1915 if( col1domcol2 )
1916 {
1917 SCIP_CALL( findFixings(scip, matrix, varcol1, col1,
1918 tmpupperbounddominatingcol1, tmpwclowerbounddominatingcol1,
1919 tmplowerbounddominatingcol1, tmpwcupperbounddominatingcol1,
1920 varcol2, col2,
1921 varstofix, onlybinvars, onlyoneone, nfixings) );
1922
1923 if( presoldata->predbndstr )
1924 {
1925 SCIP_CALL( predBndStr(scip, varcol1, col1,
1926 tmpupperbounddominatingcol1,
1927 tmplowerbounddominatingcol1, tmpwcupperbounddominatingcol1,
1928 varcol2, col2,
1929 tmpupperbounddominatedcol1, tmpwclowerbounddominatedcol1,
1930 tmplowerbounddominatedcol1,
1931 varstofix, nchgbds) );
1932 }
1933 }
1934 else if( col2domcol1 )
1935 {
1936 SCIP_CALL( findFixings(scip, matrix, varcol2, col2,
1937 tmpupperbounddominatingcol2, tmpwclowerbounddominatingcol2,
1938 tmplowerbounddominatingcol2, tmpwcupperbounddominatingcol2,
1939 varcol1, col1,
1940 varstofix, onlybinvars, onlyoneone, nfixings) );
1941
1942 if( presoldata->predbndstr )
1943 {
1944 SCIP_CALL( predBndStr(scip, varcol2, col2,
1945 tmpupperbounddominatingcol2,
1946 tmplowerbounddominatingcol2, tmpwcupperbounddominatingcol2,
1947 varcol1, col1,
1948 tmpupperbounddominatedcol2, tmpwclowerbounddominatedcol2,
1949 tmplowerbounddominatedcol2,
1950 varstofix, nchgbds) );
1951 }
1952 }
1953 if( varstofix[col1] == FIXATLB )
1954 break;
1955 }
1956 }
1957
1958 return SCIP_OKAY;
1959 }
1960
1961
1962 /*
1963 * Callback methods of presolver
1964 */
1965
1966
1967 /** copy method for constraint handler plugins (called when SCIP copies plugins) */
1968 static
SCIP_DECL_PRESOLCOPY(presolCopyDomcol)1969 SCIP_DECL_PRESOLCOPY(presolCopyDomcol)
1970 { /*lint --e{715}*/
1971 assert(scip != NULL);
1972 assert(presol != NULL);
1973 assert(strcmp(SCIPpresolGetName(presol), PRESOL_NAME) == 0);
1974
1975 /* call inclusion method of presolver */
1976 SCIP_CALL( SCIPincludePresolDomcol(scip) );
1977
1978 return SCIP_OKAY;
1979 }
1980
1981 /** destructor of presolver to free user data (called when SCIP is exiting) */
1982 static
SCIP_DECL_PRESOLFREE(presolFreeDomcol)1983 SCIP_DECL_PRESOLFREE(presolFreeDomcol)
1984 { /*lint --e{715}*/
1985 SCIP_PRESOLDATA* presoldata;
1986
1987 /* free presolver data */
1988 presoldata = SCIPpresolGetData(presol);
1989 assert(presoldata != NULL);
1990
1991 SCIPfreeBlockMemory(scip, &presoldata);
1992 SCIPpresolSetData(presol, NULL);
1993
1994 return SCIP_OKAY;
1995 }
1996
1997 /** execution method of presolver */
1998 static
SCIP_DECL_PRESOLEXEC(presolExecDomcol)1999 SCIP_DECL_PRESOLEXEC(presolExecDomcol)
2000 { /*lint --e{715}*/
2001 SCIP_PRESOLDATA* presoldata;
2002 SCIP_MATRIX* matrix;
2003 SCIP_Bool initialized;
2004 SCIP_Bool complete;
2005 SCIP_Bool infeasible;
2006 int nfixings;
2007 SCIP_Longint ndomrelations;
2008 int v;
2009 int r;
2010 FIXINGDIRECTION* varstofix;
2011 SCIP_Bool* varsprocessed;
2012 int nrows;
2013 int ncols;
2014 int* rowidxsorted;
2015 int* rowsparsity;
2016 int varcount;
2017 int* consearchcols;
2018 int* intsearchcols;
2019 int* binsearchcols;
2020 int nconfill;
2021 int nintfill;
2022 int nbinfill;
2023 #ifdef SCIP_DEBUG
2024 int nconvarsfixed = 0;
2025 int nintvarsfixed = 0;
2026 int nbinvarsfixed = 0;
2027 #endif
2028 int* pclass;
2029 int* colidx;
2030 int pclassstart;
2031 int pc;
2032 SCIP_Bool* varineq;
2033
2034 assert(result != NULL);
2035 *result = SCIP_DIDNOTRUN;
2036
2037 if( !SCIPallowStrongDualReds(scip) || SCIPisStopped(scip) )
2038 return SCIP_OKAY;
2039
2040 presoldata = SCIPpresolGetData(presol);
2041 assert(presoldata != NULL);
2042
2043 /* don't run for pure LPs */
2044 if( !presoldata->continuousred && (SCIPgetNBinVars(scip) + SCIPgetNIntVars(scip) == 0) )
2045 return SCIP_OKAY;
2046
2047 *result = SCIP_DIDNOTFIND;
2048
2049 matrix = NULL;
2050 SCIP_CALL( SCIPmatrixCreate(scip, &matrix, TRUE, &initialized, &complete, &infeasible,
2051 naddconss, ndelconss, nchgcoefs, nchgbds, nfixedvars) );
2052
2053 /* if infeasibility was detected during matrix creation, return here */
2054 if( infeasible )
2055 {
2056 if( initialized )
2057 SCIPmatrixFree(scip, &matrix);
2058
2059 *result = SCIP_CUTOFF;
2060 return SCIP_OKAY;
2061 }
2062
2063 if( !initialized )
2064 return SCIP_OKAY;
2065
2066 nfixings = 0;
2067 ndomrelations = 0;
2068 nrows = SCIPmatrixGetNRows(matrix);
2069 ncols = SCIPmatrixGetNColumns(matrix);
2070 assert(SCIPgetNVars(scip) == ncols);
2071
2072 SCIP_CALL( SCIPallocBufferArray(scip, &varstofix, ncols) );
2073 BMSclearMemoryArray(varstofix, ncols);
2074
2075 SCIP_CALL( SCIPallocBufferArray(scip, &varsprocessed, ncols) );
2076
2077 /* do not process columns that do not have all locks represented in the matrix */
2078 for( v = 0; v < ncols; ++v )
2079 varsprocessed[v] = SCIPmatrixUplockConflict(matrix, v) || SCIPmatrixDownlockConflict(matrix, v);
2080
2081 SCIP_CALL( SCIPallocBufferArray(scip, &consearchcols, ncols) );
2082 SCIP_CALL( SCIPallocBufferArray(scip, &intsearchcols, ncols) );
2083 SCIP_CALL( SCIPallocBufferArray(scip, &binsearchcols, ncols) );
2084
2085 SCIP_CALL( SCIPallocBufferArray(scip, &rowidxsorted, nrows) );
2086 SCIP_CALL( SCIPallocBufferArray(scip, &rowsparsity, nrows) );
2087 for( r = 0; r < nrows; ++r )
2088 {
2089 rowidxsorted[r] = r;
2090 rowsparsity[r] = SCIPmatrixGetRowNNonzs(matrix, r);
2091 }
2092
2093 SCIP_CALL( SCIPallocBufferArray(scip, &pclass, ncols) );
2094 SCIP_CALL( SCIPallocBufferArray(scip, &colidx, ncols) );
2095 SCIP_CALL( SCIPallocBufferArray(scip, &varineq, ncols) );
2096 for( v = 0; v < ncols; v++ )
2097 {
2098 colidx[v] = v;
2099 varineq[v] = FALSE;
2100 }
2101
2102 /* init pair comparision control */
2103 presoldata->numcurrentpairs = presoldata->nummaxpairs;
2104
2105 varcount = 0;
2106
2107 /* 1.stage: search dominance relations of parallel columns
2108 * within equalities and ranged rows
2109 */
2110 if( (presoltiming & SCIP_PRESOLTIMING_EXHAUSTIVE) != 0 )
2111 {
2112 SCIP_CALL( detectParallelCols(scip, matrix, pclass, varineq) );
2113 SCIPsortIntInt(pclass, colidx, ncols);
2114
2115 pc = 0;
2116 while( pc < ncols )
2117 {
2118 int varidx;
2119
2120 varidx = 0;
2121 nconfill = 0;
2122 nintfill = 0;
2123 nbinfill = 0;
2124
2125 pclassstart = pclass[pc];
2126 while( pc < ncols && pclassstart == pclass[pc] )
2127 {
2128 SCIP_VAR* var;
2129
2130 varidx = colidx[pc];
2131 var = SCIPmatrixGetVar(matrix, varidx);
2132
2133 /* we only regard variables which were not processed yet and
2134 * are present within equalities or ranged rows
2135 */
2136 if( !varsprocessed[varidx] && varineq[varidx] )
2137 {
2138 /* we search only for dominance relations between the same variable type */
2139 if( SCIPvarGetType(var) == SCIP_VARTYPE_CONTINUOUS )
2140 {
2141 consearchcols[nconfill++] = varidx;
2142 }
2143 else if( SCIPvarIsBinary(var) )
2144 {
2145 binsearchcols[nbinfill++] = varidx;
2146 }
2147 else
2148 {
2149 assert(SCIPvarGetType(var) == SCIP_VARTYPE_INTEGER || SCIPvarGetType(var) == SCIP_VARTYPE_IMPLINT);
2150 intsearchcols[nintfill++] = varidx;
2151 }
2152 }
2153 ++pc;
2154 }
2155
2156 /* continuous variables */
2157 if( nconfill > 1 && presoldata->continuousred )
2158 {
2159 SCIP_CALL( findDominancePairs(scip, matrix, presoldata, consearchcols, nconfill, FALSE,
2160 varstofix, &nfixings, &ndomrelations, nchgbds) );
2161
2162 for( v = 0; v < nconfill; ++v )
2163 varsprocessed[consearchcols[v]] = TRUE;
2164
2165 varcount += nconfill;
2166 }
2167 else if( nconfill == 1 )
2168 {
2169 if( varineq[varidx] )
2170 varsprocessed[consearchcols[0]] = TRUE;
2171 }
2172
2173 /* integer and impl-integer variables */
2174 if( nintfill > 1 )
2175 {
2176 SCIP_CALL( findDominancePairs(scip, matrix, presoldata, intsearchcols, nintfill, FALSE,
2177 varstofix, &nfixings, &ndomrelations, nchgbds) );
2178
2179 for( v = 0; v < nintfill; ++v )
2180 varsprocessed[intsearchcols[v]] = TRUE;
2181
2182 varcount += nintfill;
2183 }
2184 else if( nintfill == 1 )
2185 {
2186 if( varineq[varidx] )
2187 varsprocessed[intsearchcols[0]] = TRUE;
2188 }
2189
2190 /* binary variables */
2191 if( nbinfill > 1 )
2192 {
2193 SCIP_CALL( findDominancePairs(scip, matrix, presoldata, binsearchcols, nbinfill, TRUE,
2194 varstofix, &nfixings, &ndomrelations, nchgbds) );
2195
2196 for( v = 0; v < nbinfill; ++v )
2197 varsprocessed[binsearchcols[v]] = TRUE;
2198
2199 varcount += nbinfill;
2200 }
2201 else if( nbinfill == 1 )
2202 {
2203 if( varineq[varidx] )
2204 varsprocessed[binsearchcols[0]] = TRUE;
2205 }
2206
2207 if( varcount >= ncols )
2208 break;
2209 }
2210 }
2211
2212 /* 2.stage: search dominance relations for the remaining columns
2213 * by increasing row-sparsity
2214 */
2215 if( (presoltiming & SCIP_PRESOLTIMING_EXHAUSTIVE) != 0 )
2216 {
2217 SCIPsortIntInt(rowsparsity, rowidxsorted, nrows);
2218
2219 for( r = 0; r < nrows; ++r )
2220 {
2221 int rowidx;
2222 int* rowpnt;
2223 int* rowend;
2224
2225 /* break if the time limit was reached; since the check is expensive,
2226 * we only check all 1000 constraints
2227 */
2228 if( (r % 1000 == 0) && SCIPisStopped(scip) )
2229 break;
2230
2231 rowidx = rowidxsorted[r];
2232 rowpnt = SCIPmatrixGetRowIdxPtr(matrix, rowidx);
2233 rowend = rowpnt + SCIPmatrixGetRowNNonzs(matrix, rowidx);
2234
2235 if( SCIPmatrixGetRowNNonzs(matrix, rowidx) == 1 )
2236 continue;
2237
2238 nconfill = 0;
2239 nintfill = 0;
2240 nbinfill = 0;
2241
2242 for( ; rowpnt < rowend; rowpnt++ )
2243 {
2244 if( !(varsprocessed[*rowpnt]) )
2245 {
2246 int varidx;
2247 SCIP_VAR* var;
2248
2249 varidx = *rowpnt;
2250 var = SCIPmatrixGetVar(matrix, varidx);
2251
2252 /* we search only for dominance relations between the same variable type */
2253 if( SCIPvarGetType(var) == SCIP_VARTYPE_CONTINUOUS )
2254 {
2255 consearchcols[nconfill++] = varidx;
2256 }
2257 else if( SCIPvarIsBinary(var) )
2258 {
2259 binsearchcols[nbinfill++] = varidx;
2260 }
2261 else
2262 {
2263 assert(SCIPvarGetType(var) == SCIP_VARTYPE_INTEGER || SCIPvarGetType(var) == SCIP_VARTYPE_IMPLINT);
2264 intsearchcols[nintfill++] = varidx;
2265 }
2266 }
2267 }
2268
2269 /* continuous variables */
2270 if( nconfill > 1 && presoldata->continuousred )
2271 {
2272 SCIP_CALL( findDominancePairs(scip, matrix, presoldata, consearchcols, nconfill, FALSE,
2273 varstofix, &nfixings, &ndomrelations, nchgbds) );
2274
2275 for( v = 0; v < nconfill; ++v )
2276 varsprocessed[consearchcols[v]] = TRUE;
2277
2278 varcount += nconfill;
2279 }
2280
2281 /* integer and impl-integer variables */
2282 if( nintfill > 1 )
2283 {
2284 SCIP_CALL( findDominancePairs(scip, matrix, presoldata, intsearchcols, nintfill, FALSE,
2285 varstofix, &nfixings, &ndomrelations, nchgbds) );
2286
2287 for( v = 0; v < nintfill; ++v )
2288 varsprocessed[intsearchcols[v]] = TRUE;
2289
2290 varcount += nintfill;
2291 }
2292
2293 /* binary variables */
2294 if( nbinfill > 1 )
2295 {
2296 SCIP_CALL( findDominancePairs(scip, matrix, presoldata, binsearchcols, nbinfill, TRUE,
2297 varstofix, &nfixings, &ndomrelations, nchgbds) );
2298
2299 for( v = 0; v < nbinfill; ++v )
2300 varsprocessed[binsearchcols[v]] = TRUE;
2301
2302 varcount += nbinfill;
2303 }
2304
2305 if( varcount >= ncols )
2306 break;
2307 }
2308 }
2309
2310 if( nfixings > 0 )
2311 {
2312 int oldnfixedvars;
2313
2314 oldnfixedvars = *nfixedvars;
2315
2316 for( v = ncols - 1; v >= 0; --v )
2317 {
2318 SCIP_Bool fixed;
2319 SCIP_VAR* var;
2320
2321 var = SCIPmatrixGetVar(matrix,v);
2322
2323 /* there should be no fixings for variables with lock conflicts,
2324 * since they where marked as processed and therefore should be skipped
2325 */
2326 assert(varstofix[v] == NOFIX || (!SCIPmatrixUplockConflict(matrix, v) && !SCIPmatrixDownlockConflict(matrix, v)) );
2327
2328 if( varstofix[v] == FIXATLB )
2329 {
2330 SCIP_Real lb;
2331
2332 lb = SCIPvarGetLbGlobal(var);
2333 assert(SCIPvarGetType(var) == SCIP_VARTYPE_CONTINUOUS || SCIPisFeasIntegral(scip, lb));
2334
2335 /* fix at lower bound */
2336 SCIP_CALL( SCIPfixVar(scip, var, lb, &infeasible, &fixed) );
2337 if( infeasible )
2338 {
2339 SCIPdebugMsg(scip, " -> infeasible fixing\n");
2340 *result = SCIP_CUTOFF;
2341
2342 break;
2343 }
2344 assert(fixed);
2345 (*nfixedvars)++;
2346
2347 #ifdef SCIP_DEBUG
2348 if( SCIPvarGetType(var) == SCIP_VARTYPE_CONTINUOUS )
2349 nconvarsfixed++;
2350 else if( SCIPvarIsBinary(var) )
2351 nbinvarsfixed++;
2352 else
2353 nintvarsfixed++;
2354 #endif
2355 }
2356 else if( varstofix[v] == FIXATUB )
2357 {
2358 SCIP_Real ub;
2359
2360 ub = SCIPvarGetUbGlobal(var);
2361 assert(SCIPvarGetType(var) == SCIP_VARTYPE_CONTINUOUS || SCIPisFeasIntegral(scip, ub));
2362
2363 /* fix at upper bound */
2364 SCIP_CALL( SCIPfixVar(scip, var, ub, &infeasible, &fixed) );
2365 if( infeasible )
2366 {
2367 SCIPdebugMsg(scip, " -> infeasible fixing\n");
2368 *result = SCIP_CUTOFF;
2369
2370 break;
2371 }
2372 assert(fixed);
2373 (*nfixedvars)++;
2374
2375 #ifdef SCIP_DEBUG
2376 if( SCIPvarGetType(var) == SCIP_VARTYPE_CONTINUOUS )
2377 nconvarsfixed++;
2378 else if( SCIPvarIsBinary(var) )
2379 nbinvarsfixed++;
2380 else
2381 nintvarsfixed++;
2382 #endif
2383 }
2384 }
2385
2386 if( *result != SCIP_CUTOFF && *nfixedvars > oldnfixedvars )
2387 *result = SCIP_SUCCESS;
2388 }
2389
2390 SCIPfreeBufferArray(scip, &varineq);
2391 SCIPfreeBufferArray(scip, &colidx);
2392 SCIPfreeBufferArray(scip, &pclass);
2393 SCIPfreeBufferArray(scip, &rowsparsity);
2394 SCIPfreeBufferArray(scip, &rowidxsorted);
2395 SCIPfreeBufferArray(scip, &binsearchcols);
2396 SCIPfreeBufferArray(scip, &intsearchcols);
2397 SCIPfreeBufferArray(scip, &consearchcols);
2398 SCIPfreeBufferArray(scip, &varsprocessed);
2399 SCIPfreeBufferArray(scip, &varstofix);
2400
2401 #ifdef SCIP_DEBUG
2402 if( (nconvarsfixed + nintvarsfixed + nbinvarsfixed) > 0 )
2403 {
2404 SCIPdebugMsg(scip, "### %d vars [%" SCIP_LONGINT_FORMAT " dom] => fixed [cont: %d, int: %d, bin: %d], %scutoff detected\n",
2405 ncols, ndomrelations, nconvarsfixed, nintvarsfixed, nbinvarsfixed, (*result != SCIP_CUTOFF) ? "no " : "");
2406 }
2407 #endif
2408
2409 SCIPmatrixFree(scip, &matrix);
2410
2411 return SCIP_OKAY;
2412 }
2413
2414 /*
2415 * presolver specific interface methods
2416 */
2417
2418 /** creates the domcol presolver and includes it in SCIP */
SCIPincludePresolDomcol(SCIP * scip)2419 SCIP_RETCODE SCIPincludePresolDomcol(
2420 SCIP* scip /**< SCIP data structure */
2421 )
2422 {
2423 SCIP_PRESOLDATA* presoldata;
2424 SCIP_PRESOL* presol;
2425
2426 /* create domcol presolver data */
2427 SCIP_CALL( SCIPallocBlockMemory(scip, &presoldata) );
2428
2429 /* include presolver */
2430 SCIP_CALL( SCIPincludePresolBasic(scip, &presol, PRESOL_NAME, PRESOL_DESC, PRESOL_PRIORITY, PRESOL_MAXROUNDS,
2431 PRESOL_TIMING, presolExecDomcol, presoldata) );
2432 SCIP_CALL( SCIPsetPresolCopy(scip, presol, presolCopyDomcol) );
2433 SCIP_CALL( SCIPsetPresolFree(scip, presol, presolFreeDomcol) );
2434
2435 SCIP_CALL( SCIPaddIntParam(scip,
2436 "presolving/domcol/numminpairs",
2437 "minimal number of pair comparisons",
2438 &presoldata->numminpairs, FALSE, DEFAULT_NUMMINPAIRS, 100, DEFAULT_NUMMAXPAIRS, NULL, NULL) );
2439
2440 SCIP_CALL( SCIPaddIntParam(scip,
2441 "presolving/domcol/nummaxpairs",
2442 "maximal number of pair comparisons",
2443 &presoldata->nummaxpairs, FALSE, DEFAULT_NUMMAXPAIRS, DEFAULT_NUMMINPAIRS, 1000000000, NULL, NULL) );
2444
2445 SCIP_CALL( SCIPaddBoolParam(scip,
2446 "presolving/domcol/predbndstr",
2447 "should predictive bound strengthening be applied?",
2448 &presoldata->predbndstr, FALSE, DEFAULT_PREDBNDSTR, NULL, NULL) );
2449
2450 SCIP_CALL( SCIPaddBoolParam(scip,
2451 "presolving/domcol/continuousred",
2452 "should reductions for continuous variables be performed?",
2453 &presoldata->continuousred, FALSE, DEFAULT_CONTINUOUS_RED, NULL, NULL) );
2454
2455 return SCIP_OKAY;
2456 }
2457