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_stuffing.c
17  * @ingroup DEFPLUGINS_PRESOL
18  * @brief  fix singleton continuous variables
19  * @author Dieter Weninger
20  *
21  * Investigate singleton continuous variables if one can be fixed at a bound.
22  *
23  * @todo enhancement from singleton continuous variables to continuous
24  *       variables with only one lock in a common row
25  */
26 
27 /*---+----1----+----2----+----3----+----4----+----5----+----6----+----7----+----8----+----9----+----0----+----1----+----2*/
28 
29 #include "blockmemshell/memory.h"
30 #include "scip/presol_stuffing.h"
31 #include "scip/pub_matrix.h"
32 #include "scip/pub_message.h"
33 #include "scip/pub_misc_sort.h"
34 #include "scip/pub_presol.h"
35 #include "scip/pub_var.h"
36 #include "scip/scip_general.h"
37 #include "scip/scip_mem.h"
38 #include "scip/scip_message.h"
39 #include "scip/scip_nlp.h"
40 #include "scip/scip_numerics.h"
41 #include "scip/scip_presol.h"
42 #include "scip/scip_pricer.h"
43 #include "scip/scip_prob.h"
44 #include "scip/scip_probing.h"
45 #include "scip/scip_var.h"
46 #include <string.h>
47 
48 #define PRESOL_NAME            "stuffing"
49 #define PRESOL_DESC            "fix redundant singleton continuous variables"
50 #define PRESOL_PRIORITY             -100     /**< priority of the presolver (>= 0: before, < 0: after constraint handlers) */
51 #define PRESOL_MAXROUNDS               0     /**< maximal number of presolving rounds the presolver participates in (-1: no limit) */
52 #define PRESOL_TIMING           SCIP_PRESOLTIMING_EXHAUSTIVE /* timing of the presolver (fast, medium, or exhaustive) */
53 
54 /** type of fixing direction */
55 enum Fixingdirection
56 {
57    FIXATLB = -1,         /**< fix variable at lower bound */
58    NOFIX   =  0,         /**< do not fix variable */
59    FIXATUB =  1          /**< fix variable at upper bound */
60 };
61 typedef enum Fixingdirection FIXINGDIRECTION;
62 
63 /*
64  * Local methods
65  */
66 
67 /** try to fix singleton continuous variables */
68 static
singletonColumnStuffing(SCIP * scip,SCIP_MATRIX * matrix,FIXINGDIRECTION * varstofix,int * nfixings)69 SCIP_RETCODE singletonColumnStuffing(
70    SCIP*                 scip,               /**< SCIP main data structure */
71    SCIP_MATRIX*          matrix,             /**< matrix containing the constraints */
72    FIXINGDIRECTION*      varstofix,          /**< array holding fixing information */
73    int*                  nfixings            /**< number of possible fixings */
74    )
75 {
76    SCIP_Real* valpnt;
77    SCIP_Real* colratios;
78    SCIP_Real* colcoeffs;
79    SCIP_Bool* rowprocessed;
80    int* rowpnt;
81    int* rowend;
82    int* colindices;
83    int* dummy;
84    SCIP_Bool* swapped;
85    SCIP_Real upperconst;
86    SCIP_Real lowerconst;
87    SCIP_Real rhs;
88    SCIP_Bool tryfixing;
89    int idx;
90    int col;
91    int row;
92    int fillcnt;
93    int k;
94    int nrows;
95    int ncols;
96 
97    assert(scip != NULL);
98    assert(matrix != NULL);
99    assert(varstofix != NULL);
100    assert(nfixings != NULL);
101 
102    nrows = SCIPmatrixGetNRows(matrix);
103    ncols = SCIPmatrixGetNColumns(matrix);
104 
105    SCIP_CALL( SCIPallocBufferArray(scip, &colindices, ncols) );
106    SCIP_CALL( SCIPallocBufferArray(scip, &colratios, ncols) );
107    SCIP_CALL( SCIPallocBufferArray(scip, &colcoeffs, ncols) );
108    SCIP_CALL( SCIPallocBufferArray(scip, &dummy, ncols) );
109 
110    SCIP_CALL( SCIPallocBufferArray(scip, &swapped, ncols) );
111    BMSclearMemoryArray(swapped, ncols);
112 
113    SCIP_CALL( SCIPallocBufferArray(scip, &rowprocessed, nrows) );
114    BMSclearMemoryArray(rowprocessed, nrows);
115 
116    for( col = 0; col < ncols; col++ )
117    {
118       /* consider only rows with minimal one continuous singleton column */
119       if( SCIPmatrixGetColNNonzs(matrix, col) == 1
120          && SCIPvarGetType(SCIPmatrixGetVar(matrix, col)) == SCIP_VARTYPE_CONTINUOUS
121          && SCIPvarGetNLocksUpType(SCIPmatrixGetVar(matrix, col), SCIP_LOCKTYPE_MODEL) == SCIPmatrixGetColNUplocks(matrix, col)
122          && SCIPvarGetNLocksDownType(SCIPmatrixGetVar(matrix, col), SCIP_LOCKTYPE_MODEL) == SCIPmatrixGetColNDownlocks(matrix, col) )
123       {
124          row = *(SCIPmatrixGetColIdxPtr(matrix, col));
125          if( rowprocessed[row] )
126             continue;
127 
128          rowprocessed[row] = TRUE;
129 
130          /* treat all >= rows from matrix, but internally we transform to <= relation */
131          if( SCIPmatrixIsRowRhsInfinity(matrix, row) )
132          {
133             fillcnt = 0;
134             tryfixing = TRUE;
135             upperconst = 0.0;
136             lowerconst = 0.0;
137             rhs = -SCIPmatrixGetRowLhs(matrix, row);
138 
139             rowpnt = SCIPmatrixGetRowIdxPtr(matrix, row);
140             rowend = rowpnt + SCIPmatrixGetRowNNonzs(matrix, row);
141             valpnt = SCIPmatrixGetRowValPtr(matrix, row);
142 
143             for( ; (rowpnt < rowend); rowpnt++, valpnt++ )
144             {
145                SCIP_Real coef;
146                SCIP_VAR* var;
147                SCIP_Real lb;
148                SCIP_Real ub;
149 
150                coef = -(*valpnt);
151                idx = *rowpnt;
152                var = SCIPmatrixGetVar(matrix, idx);
153                lb = SCIPvarGetLbGlobal(var);
154                ub = SCIPvarGetUbGlobal(var);
155 
156                /* we need to check if this is a singleton continuous variable and
157                 * all constraints containing this variable are present inside
158                 * the mixed integer linear matrix
159                 */
160                if( SCIPmatrixGetColNNonzs(matrix, idx) == 1
161                   && (SCIPvarGetNLocksUpType(var, SCIP_LOCKTYPE_MODEL) + SCIPvarGetNLocksDownType(var, SCIP_LOCKTYPE_MODEL)) == 1
162                   && SCIPvarGetType(var) == SCIP_VARTYPE_CONTINUOUS )
163                {
164                   if( SCIPisLT(scip, SCIPvarGetObj(var), 0.0) && SCIPisGT(scip, coef, 0.0) )
165                   {
166                      /* case 1: obj < 0 and coef > 0 */
167                      if( SCIPisInfinity(scip, -lb) )
168                      {
169                         tryfixing = FALSE;
170                         break;
171                      }
172 
173                      upperconst += coef * lb;
174                      lowerconst += coef * lb;
175                      colratios[fillcnt] = SCIPvarGetObj(var) / coef;
176                      colindices[fillcnt] = idx;
177                      colcoeffs[fillcnt] = coef;
178                      fillcnt++;
179                   }
180                   else if( SCIPisGT(scip, SCIPvarGetObj(var), 0.0) && SCIPisLT(scip, coef, 0.0) )
181                   {
182                      /* case 2: obj > 0 and coef < 0 */
183                      if( SCIPisInfinity(scip, ub) )
184                      {
185                         tryfixing = FALSE;
186                         break;
187                      }
188 
189                      /* multiply column by (-1) to become case 1.
190                       * now bounds are swapped: ub := -lb, lb := -ub
191                       */
192                      swapped[idx] = TRUE;
193                      upperconst += coef * ub;
194                      lowerconst += coef * ub;
195                      colratios[fillcnt] = SCIPvarGetObj(var) / coef;
196                      colindices[fillcnt] = idx;
197                      colcoeffs[fillcnt] = -coef;
198                      fillcnt++;
199                   }
200                   else if( SCIPisGE(scip, SCIPvarGetObj(var), 0.0) && SCIPisGE(scip, coef, 0.0) )
201                   {
202                      /* case 3: obj >= 0 and coef >= 0 is handled by duality fixing.
203                       *  we only consider the lower bound for the constants
204                       */
205                      if( SCIPisInfinity(scip, -lb) )
206                      {
207                         /* maybe unbounded */
208                         tryfixing = FALSE;
209                         break;
210                      }
211 
212                      upperconst += coef * lb;
213                      lowerconst += coef * lb;
214                   }
215                   else
216                   {
217                      /* case 4: obj <= 0 and coef <= 0 is also handled by duality fixing.
218                       * we only consider the upper bound for the constants
219                       */
220                      assert(SCIPisLE(scip, SCIPvarGetObj(var), 0.0));
221                      assert(SCIPisLE(scip, coef, 0.0));
222 
223                      if( SCIPisInfinity(scip, ub) )
224                      {
225                         /* maybe unbounded */
226                         tryfixing = FALSE;
227                         break;
228                      }
229 
230                      upperconst += coef * ub;
231                      lowerconst += coef * ub;
232                   }
233                }
234                else
235                {
236                   /* consider contribution of discrete variables, non-singleton
237                    * continuous variables and variables with more than one lock
238                    */
239                   if( SCIPisInfinity(scip, -lb) || SCIPisInfinity(scip, ub) )
240                   {
241                      tryfixing = FALSE;
242                      break;
243                   }
244 
245                   if( coef > 0 )
246                   {
247                      upperconst += coef * ub;
248                      lowerconst += coef * lb;
249                   }
250                   else
251                   {
252                      upperconst += coef * lb;
253                      lowerconst += coef * ub;
254                   }
255                }
256             }
257 
258             if( tryfixing )
259             {
260                SCIPsortRealRealIntInt(colratios, colcoeffs, colindices, dummy, fillcnt);
261 
262                /* verify which singleton continuous variable can be fixed */
263                for( k = 0; k < fillcnt; k++ )
264                {
265                   SCIP_VAR* var;
266                   SCIP_Real lb;
267                   SCIP_Real ub;
268                   SCIP_Real delta;
269 
270                   idx = colindices[k];
271                   var = SCIPmatrixGetVar(matrix, idx);
272                   lb = SCIPvarGetLbGlobal(var);
273                   ub = SCIPvarGetUbGlobal(var);
274 
275                   /* stop fixing if variable bounds are not finite */
276                   if( SCIPisInfinity(scip, -lb) || SCIPisInfinity(scip, ub) )
277                      break;
278 
279                   assert(SCIPmatrixGetColNNonzs(matrix, idx) == 1);
280                   assert(SCIPvarGetType(var) == SCIP_VARTYPE_CONTINUOUS);
281                   assert((SCIPvarGetNLocksUpType(var, SCIP_LOCKTYPE_MODEL) + SCIPvarGetNLocksDownType(var, SCIP_LOCKTYPE_MODEL)) == 1);
282                   assert(colcoeffs[k] >= 0);
283 
284                   /* calculate the change in the row activities if this variable changes
285                    * its value from to its other bound
286                    */
287                   if( swapped[idx] )
288                      delta = -(lb - ub) * colcoeffs[k];
289                   else
290                      delta =  (ub - lb) * colcoeffs[k];
291 
292                   assert(delta >= 0);
293 
294                   if( SCIPisLE(scip, delta, rhs - upperconst) )
295                   {
296                      if( swapped[idx] )
297                         varstofix[idx] = FIXATLB;
298                      else
299                         varstofix[idx] = FIXATUB;
300 
301                      (*nfixings)++;
302                   }
303                   else if( SCIPisLE(scip, rhs, lowerconst) )
304                   {
305                      if( swapped[idx] )
306                         varstofix[idx] = FIXATUB;
307                      else
308                         varstofix[idx] = FIXATLB;
309 
310                      (*nfixings)++;
311                   }
312 
313                   upperconst += delta;
314                   lowerconst += delta;
315                }
316             }
317          }
318       }
319    }
320 
321    SCIPfreeBufferArray(scip, &rowprocessed);
322    SCIPfreeBufferArray(scip, &swapped);
323    SCIPfreeBufferArray(scip, &dummy);
324    SCIPfreeBufferArray(scip, &colcoeffs);
325    SCIPfreeBufferArray(scip, &colratios);
326    SCIPfreeBufferArray(scip, &colindices);
327 
328    return SCIP_OKAY;
329 }
330 
331 /*
332  * Callback methods of presolver
333  */
334 
335 /** copy method for constraint handler plugins (called when SCIP copies plugins) */
336 static
SCIP_DECL_PRESOLCOPY(presolCopyStuffing)337 SCIP_DECL_PRESOLCOPY(presolCopyStuffing)
338 {  /*lint --e{715}*/
339    assert(scip != NULL);
340    assert(presol != NULL);
341    assert(strcmp(SCIPpresolGetName(presol), PRESOL_NAME) == 0);
342 
343    /* call inclusion method of presolver */
344    SCIP_CALL( SCIPincludePresolStuffing(scip) );
345 
346    return SCIP_OKAY;
347 }
348 
349 /** execution method of presolver */
350 static
SCIP_DECL_PRESOLEXEC(presolExecStuffing)351 SCIP_DECL_PRESOLEXEC(presolExecStuffing)
352 {  /*lint --e{715}*/
353    SCIP_MATRIX* matrix;
354    SCIP_Bool initialized;
355    SCIP_Bool complete;
356    SCIP_Bool infeasible;
357 
358    assert(result != NULL);
359    *result = SCIP_DIDNOTRUN;
360 
361    if( (SCIPgetStage(scip) != SCIP_STAGE_PRESOLVING) || SCIPinProbing(scip) || SCIPisNLPEnabled(scip) )
362       return SCIP_OKAY;
363 
364    if( SCIPgetNContVars(scip) == 0 || SCIPisStopped(scip) || SCIPgetNActivePricers(scip) > 0 )
365       return SCIP_OKAY;
366 
367    if( !SCIPallowStrongDualReds(scip) )
368       return SCIP_OKAY;
369 
370    *result = SCIP_DIDNOTFIND;
371 
372    matrix = NULL;
373    SCIP_CALL( SCIPmatrixCreate(scip, &matrix, FALSE, &initialized, &complete, &infeasible,
374       naddconss, ndelconss, nchgcoefs, nchgbds, nfixedvars) );
375 
376    /* if infeasibility was detected during matrix creation, return here */
377    if( infeasible )
378    {
379       if( initialized )
380          SCIPmatrixFree(scip, &matrix);
381 
382       *result = SCIP_CUTOFF;
383       return SCIP_OKAY;
384    }
385 
386    if( initialized )
387    {
388       FIXINGDIRECTION* varstofix;
389       int nfixings;
390       int ncols;
391 
392       nfixings = 0;
393       ncols = SCIPmatrixGetNColumns(matrix);
394 
395       SCIP_CALL( SCIPallocBufferArray(scip, &varstofix, ncols) );
396       BMSclearMemoryArray(varstofix, ncols);
397 
398       SCIP_CALL( singletonColumnStuffing(scip, matrix, varstofix, &nfixings) );
399 
400       if( nfixings > 0 )
401       {
402          int v;
403          int oldnfixedvars;
404 
405          oldnfixedvars = *nfixedvars;
406 
407          /* look for fixable variables */
408          for( v = ncols - 1; v >= 0; --v )
409          {
410             SCIP_Bool fixed;
411             SCIP_VAR* var;
412 
413             var = SCIPmatrixGetVar(matrix, v);
414 
415             if( varstofix[v] == FIXATLB )
416             {
417                SCIP_Real lb;
418 
419                assert(SCIPvarGetNLocksUpType(var, SCIP_LOCKTYPE_MODEL) == SCIPmatrixGetColNUplocks(matrix, v)
420                   && SCIPvarGetNLocksDownType(var, SCIP_LOCKTYPE_MODEL) == SCIPmatrixGetColNDownlocks(matrix, v));
421 
422                lb = SCIPvarGetLbGlobal(var);
423                assert(SCIPvarGetType(var) == SCIP_VARTYPE_CONTINUOUS);
424 
425                /* avoid fixings to infinite values */
426                assert(!SCIPisInfinity(scip, -lb));
427 
428                SCIPdebugMsg(scip, "Fix variable %s at lower bound %.15g\n", SCIPvarGetName(var), lb);
429 
430                /* fix at lower bound */
431                SCIP_CALL( SCIPfixVar(scip, var, lb, &infeasible, &fixed) );
432                if( infeasible )
433                {
434                   SCIPdebugMsg(scip, " -> infeasible fixing\n");
435                   *result = SCIP_CUTOFF;
436 
437                   break;
438                }
439                assert(fixed);
440                (*nfixedvars)++;
441             }
442             else if( varstofix[v] == FIXATUB )
443             {
444                SCIP_Real ub;
445 
446                assert(SCIPvarGetNLocksUpType(var, SCIP_LOCKTYPE_MODEL) == SCIPmatrixGetColNUplocks(matrix, v)
447                   && SCIPvarGetNLocksDownType(var, SCIP_LOCKTYPE_MODEL) == SCIPmatrixGetColNDownlocks(matrix, v));
448 
449                ub = SCIPvarGetUbGlobal(var);
450                assert(SCIPvarGetType(var) == SCIP_VARTYPE_CONTINUOUS);
451 
452                /* avoid fixings to infinite values */
453                assert(!SCIPisInfinity(scip, ub));
454 
455                SCIPdebugMsg(scip, "Fix variable %s at upper bound %.15g\n", SCIPvarGetName(var), ub);
456 
457                /* fix at upper bound */
458                SCIP_CALL( SCIPfixVar(scip, var, ub, &infeasible, &fixed) );
459                if( infeasible )
460                {
461                   SCIPdebugMsg(scip, " -> infeasible fixing\n");
462                   *result = SCIP_CUTOFF;
463 
464                   break;
465                }
466                assert(fixed);
467 
468                (*nfixedvars)++;
469             }
470          }
471 
472          if( *result != SCIP_CUTOFF && *nfixedvars > oldnfixedvars )
473             *result = SCIP_SUCCESS;
474       }
475 
476       SCIPfreeBufferArray(scip, &varstofix);
477    }
478 
479    SCIPmatrixFree(scip, &matrix);
480 
481    return SCIP_OKAY;
482 }
483 
484 /*
485  * presolver specific interface methods
486  */
487 
488 /** creates the stuffing presolver and includes it in SCIP */
SCIPincludePresolStuffing(SCIP * scip)489 SCIP_RETCODE SCIPincludePresolStuffing(
490    SCIP*                 scip                /**< SCIP data structure */
491    )
492 {
493    SCIP_PRESOL* presol;
494 
495    /* include presolver */
496    SCIP_CALL( SCIPincludePresolBasic(scip, &presol, PRESOL_NAME, PRESOL_DESC, PRESOL_PRIORITY, PRESOL_MAXROUNDS,
497          PRESOL_TIMING, presolExecStuffing, NULL) );
498    SCIP_CALL( SCIPsetPresolCopy(scip, presol, presolCopyStuffing) );
499 
500    return SCIP_OKAY;
501 }
502 
503 /*lint --e{749}*/
504