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