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 prop_obbt.c
17 * @ingroup DEFPLUGINS_PROP
18 * @brief optimization-based bound tightening propagator
19 * @author Stefan Weltge
20 * @author Benjamin Mueller
21 */
22
23 /**@todo if bound tightenings of other propagators are the reason for lpsolstat != SCIP_LPSOLSTAT_OPTIMAL, resolve LP */
24 /**@todo only run more than once in root node if primal bound improved or many cuts were added to the LP */
25 /**@todo filter bounds of a variable already if SCIPisLbBetter()/SCIPisUbBetter() would return FALSE */
26 /**@todo improve warmstarting of LP solving */
27 /**@todo include bound value (finite/infinite) into getScore() function */
28 /**@todo use unbounded ray in filtering */
29 /**@todo do we want to run if the LP is unbounded, maybe for infinite variable bounds? */
30 /**@todo add first filter round in direction of objective function */
31 /**@todo implement conflict resolving callback by calling public method of genvbounds propagator, since the reason are
32 * exactly the variable bounds with nonnegative reduced costs stored in the right-hand side of the generated
33 * generalized variable bound (however, this only makes sense if we run locally)
34 */
35
36 /*---+----1----+----2----+----3----+----4----+----5----+----6----+----7----+----8----+----9----+----0----+----1----+----2*/
37
38 #include "blockmemshell/memory.h"
39 #include "nlpi/pub_expr.h"
40 #include "scip/cons_abspower.h"
41 #include "scip/cons_bivariate.h"
42 #include "scip/cons_nonlinear.h"
43 #include "scip/cons_quadratic.h"
44 #include "scip/intervalarith.h"
45 #include "scip/prop_genvbounds.h"
46 #include "scip/prop_obbt.h"
47 #include "scip/pub_cons.h"
48 #include "scip/pub_lp.h"
49 #include "scip/pub_message.h"
50 #include "scip/pub_misc.h"
51 #include "scip/pub_misc_sort.h"
52 #include "scip/pub_nlp.h"
53 #include "scip/pub_prop.h"
54 #include "scip/pub_tree.h"
55 #include "scip/pub_var.h"
56 #include "scip/scip_cons.h"
57 #include "scip/scip_copy.h"
58 #include "scip/scip_cut.h"
59 #include "scip/scip_general.h"
60 #include "scip/scip_lp.h"
61 #include "scip/scip_mem.h"
62 #include "scip/scip_message.h"
63 #include "scip/scip_nlp.h"
64 #include "scip/scip_numerics.h"
65 #include "scip/scip_param.h"
66 #include "scip/scip_prob.h"
67 #include "scip/scip_probing.h"
68 #include "scip/scip_prop.h"
69 #include "scip/scip_randnumgen.h"
70 #include "scip/scip_solvingstats.h"
71 #include "scip/scip_tree.h"
72 #include "scip/scip_var.h"
73 #include <string.h>
74
75 #define PROP_NAME "obbt"
76 #define PROP_DESC "optimization-based bound tightening propagator"
77 #define PROP_TIMING SCIP_PROPTIMING_AFTERLPLOOP
78 #define PROP_PRIORITY -1000000 /**< propagator priority */
79 #define PROP_FREQ 0 /**< propagator frequency */
80 #define PROP_DELAY TRUE /**< should propagation method be delayed, if other propagators
81 * found reductions? */
82
83 #define DEFAULT_CREATE_GENVBOUNDS TRUE /**< should obbt try to provide genvbounds if possible? */
84 #define DEFAULT_FILTERING_NORM TRUE /**< should coefficients in filtering be normalized w.r.t. the
85 * domains sizes? */
86 #define DEFAULT_APPLY_FILTERROUNDS FALSE /**< try to filter bounds in so-called filter rounds by solving
87 * auxiliary LPs? */
88 #define DEFAULT_APPLY_TRIVIALFITLERING TRUE /**< should obbt try to use the LP solution to filter some bounds? */
89 #define DEFAULT_GENVBDSDURINGFILTER TRUE /**< try to genrate genvbounds during trivial and aggressive filtering? */
90 #define DEFAULT_DUALFEASTOL 1e-9 /**< feasibility tolerance for reduced costs used in obbt; this value
91 * is used if SCIP's dual feastol is greater */
92 #define DEFAULT_CONDITIONLIMIT -1.0 /**< maximum condition limit used in LP solver (-1.0: no limit) */
93 #define DEFAULT_BOUNDSTREPS 0.001 /**< minimal relative improve for strengthening bounds */
94 #define DEFAULT_FILTERING_MIN 2 /**< minimal number of filtered bounds to apply another filter
95 * round */
96 #define DEFAULT_ITLIMITFACTOR 10.0 /**< multiple of root node LP iterations used as total LP iteration
97 * limit for obbt (<= 0: no limit ) */
98 #define DEFAULT_MINITLIMIT 5000L /**< minimum LP iteration limit */
99 #define DEFAULT_ONLYNONCONVEXVARS FALSE /**< only apply obbt on non-convex variables */
100 #define DEFAULT_TIGHTINTBOUNDSPROBING TRUE /**< should bounds of integral variables be tightened during
101 * the probing mode? */
102 #define DEFAULT_TIGHTCONTBOUNDSPROBING FALSE /**< should bounds of continuous variables be tightened during
103 * the probing mode? */
104 #define DEFAULT_ORDERINGALGO 1 /**< which type of ordering algorithm should we use?
105 * (0: no, 1: greedy, 2: greedy reverse) */
106 #define OBBT_SCOREBASE 5 /**< base that is used to calculate a bounds score value */
107 #define GENVBOUND_PROP_NAME "genvbounds"
108 #define INTERVALINFTY 1E+43 /**< value for infinity in interval operations */
109
110 #define DEFAULT_SEPARATESOL FALSE /**< should the obbt LP solution be separated? note that that by
111 * separating solution OBBT will apply all bound tightenings
112 * immediatly */
113 #define DEFAULT_SEPAMINITER 0 /**< minimum number of iteration spend to separate an obbt LP solution */
114 #define DEFAULT_SEPAMAXITER 10 /**< maximum number of iteration spend to separate an obbt LP solution */
115 #define DEFAULT_GENVBDSDURINGSEPA TRUE /**< try to create genvbounds during separation process? */
116 #define DEFAULT_PROPAGATEFREQ 0 /**< trigger a propagation round after that many bound tightenings
117 * (0: no propagation) */
118 #define DEFAULT_CREATE_BILININEQS TRUE /**< solve auxiliary LPs in order to find valid inequalities for bilinear terms? */
119 #define DEFAULT_ITLIMITFAC_BILININEQS 3.0 /**< multiple of OBBT LP limit used as total LP iteration limit for solving bilinear inequality LPs (< 0 for no limit) */
120 #define DEFAULT_MINNONCONVEXITY 1e-1 /**< minimum nonconvexity for choosing a bilinear term */
121 #define DEFAULT_RANDSEED 149 /**< initial random seed */
122
123
124 /** translate from one value of infinity to another
125 *
126 * if val is >= infty1, then give infty2, else give val
127 */
128 #define infty2infty(infty1, infty2, val) ((val) >= (infty1) ? (infty2) : (val))
129
130 /*
131 * Data structures
132 */
133
134 /** bound data */
135 struct Bound
136 {
137 SCIP_VAR* var; /**< variable */
138 SCIP_Real newval; /**< stores a probably tighter value for this bound */
139 SCIP_BOUNDTYPE boundtype; /**< type of bound */
140 unsigned int score; /**< score value that is used to group bounds */
141 unsigned int filtered:1; /**< thrown out during pre-filtering step */
142 unsigned int found:1; /**< stores whether a probably tighter value for this bound was found */
143 unsigned int done:1; /**< has this bound been processed already? */
144 unsigned int nonconvex:1; /**< is this bound affecting a nonconvex term? */
145 int index; /**< unique index */
146 };
147 typedef struct Bound BOUND;
148
149 /* all possible corners of a rectangular domain */
150 enum Corner
151 {
152 LEFTBOTTOM = 1,
153 RIGHTBOTTOM = 2,
154 RIGHTTOP = 4,
155 LEFTTOP = 8,
156 FILTERED = 15
157 };
158 typedef enum Corner CORNER;
159
160 /** bilinear bound data */
161 struct BilinBound
162 {
163 SCIP_VAR* x; /**< first variable */
164 SCIP_VAR* y; /**< second variable */
165 int filtered; /**< corners that could be thrown out during pre-filtering step */
166 unsigned int done:1; /**< has this bilinear term been processed already? */
167 int nunderest; /**< number of constraints that require to underestimate the bilinear term */
168 int noverest; /**< number of constraints that require to overestimate the bilinear term */
169 int index; /**< index of the bilinear term in the quadratic constraint handler */
170 SCIP_Real score; /**< score value that is used to group bilinear term bounds */
171 };
172 typedef struct BilinBound BILINBOUND;
173
174 /** propagator data */
175 struct SCIP_PropData
176 {
177 BOUND** bounds; /**< array of interesting bounds */
178 BILINBOUND** bilinbounds; /**< array of interesting bilinear bounds */
179 SCIP_ROW* cutoffrow; /**< pointer to current objective cutoff row */
180 SCIP_PROP* genvboundprop; /**< pointer to genvbound propagator */
181 SCIP_RANDNUMGEN* randnumgen; /**< random number generator */
182 SCIP_Longint lastnode; /**< number of last node where obbt was performed */
183 SCIP_Longint npropagatedomreds; /**< number of domain reductions found during propagation */
184 SCIP_Longint minitlimit; /**< minimum LP iteration limit */
185 SCIP_Longint itlimitbilin; /**< total LP iterations limit for solving bilinear inequality LPs */
186 SCIP_Longint itusedbilin; /**< total LP iterations used for solving bilinear inequality LPs */
187 SCIP_Real dualfeastol; /**< feasibility tolerance for reduced costs used in obbt; this value is
188 * used if SCIP's dual feastol is greater */
189 SCIP_Real conditionlimit; /**< maximum condition limit used in LP solver (-1.0: no limit) */
190 SCIP_Real boundstreps; /**< minimal relative improve for strengthening bounds */
191 SCIP_Real itlimitfactor; /**< LP iteration limit for obbt will be this factor times total LP
192 * iterations in root node */
193 SCIP_Real itlimitfactorbilin; /**< multiple of OBBT LP limit used as total LP iteration limit for solving bilinear inequality LPs (< 0 for no limit) */
194 SCIP_Real minnonconvexity; /**< lower bound on minimum absolute value of nonconvex eigenvalues for a bilinear term */
195 SCIP_Bool applyfilterrounds; /**< apply filter rounds? */
196 SCIP_Bool applytrivialfilter; /**< should obbt try to use the LP solution to filter some bounds? */
197 SCIP_Bool genvbdsduringfilter;/**< should we try to generate genvbounds during trivial and aggressive
198 * filtering? */
199 SCIP_Bool genvbdsduringsepa; /**< try to create genvbounds during separation process? */
200 SCIP_Bool creategenvbounds; /**< should obbt try to provide genvbounds if possible? */
201 SCIP_Bool normalize; /**< should coefficients in filtering be normalized w.r.t. the domains
202 * sizes? */
203 SCIP_Bool onlynonconvexvars; /**< only apply obbt on non-convex variables */
204 SCIP_Bool tightintboundsprobing; /**< should bounds of integral variables be tightened during
205 * the probing mode? */
206 SCIP_Bool tightcontboundsprobing;/**< should bounds of continuous variables be tightened during
207 * the probing mode? */
208 SCIP_Bool separatesol; /**< should the obbt LP solution be separated? note that that by
209 * separating solution OBBT will apply all bound tightenings
210 * immediatly */
211 SCIP_Bool createbilinineqs; /**< solve auxiliary LPs in order to find valid inequalities for bilinear terms? */
212 int orderingalgo; /**< which type of ordering algorithm should we use?
213 * (0: no, 1: greedy, 2: greedy reverse) */
214 int nbounds; /**< length of interesting bounds array */
215 int nbilinbounds; /**< length of interesting bilinear bounds array */
216 int boundssize; /**< size of bounds array */
217 int nminfilter; /**< minimal number of filtered bounds to apply another filter round */
218 int lastidx; /**< index to store the last undone and unfiltered bound */
219 int lastbilinidx; /**< index to store the last undone and unfiltered bilinear bound */
220 int sepaminiter; /**< minimum number of iteration spend to separate an obbt LP solution */
221 int sepamaxiter; /**< maximum number of iteration spend to separate an obbt LP solution */
222 int propagatefreq; /**< trigger a propagation round after that many bound tightenings
223 * (0: no propagation) */
224 int propagatecounter; /**< number of bound tightenings since the last propagation round */
225 };
226
227
228 /*
229 * Local methods
230 */
231
232 /** solves the LP and handles errors */
233 static
solveLP(SCIP * scip,int itlimit,SCIP_Bool * error,SCIP_Bool * optimal)234 SCIP_RETCODE solveLP(
235 SCIP* scip, /**< SCIP data structure */
236 int itlimit, /**< maximal number of LP iterations to perform, or -1 for no limit */
237 SCIP_Bool* error, /**< pointer to store whether an unresolved LP error occurred */
238 SCIP_Bool* optimal /**< was the LP solved to optimalilty? */
239 )
240 {
241 SCIP_LPSOLSTAT lpsolstat;
242 SCIP_RETCODE retcode;
243
244 assert(scip != NULL);
245 assert(itlimit == -1 || itlimit >= 0);
246 assert(error != NULL);
247 assert(optimal != NULL);
248
249 *optimal = FALSE;
250 *error = FALSE;
251
252 retcode = SCIPsolveProbingLP(scip, itlimit, error, NULL);
253
254 lpsolstat = SCIPgetLPSolstat(scip);
255
256 /* an error should not kill the overall solving process */
257 if( retcode != SCIP_OKAY )
258 {
259 SCIPwarningMessage(scip, " error while solving LP in obbt propagator; LP solve terminated with code <%d>\n", retcode);
260 SCIPwarningMessage(scip, " this does not affect the remaining solution procedure --> continue\n");
261
262 *error = TRUE;
263
264 return SCIP_OKAY;
265 }
266
267 if( lpsolstat == SCIP_LPSOLSTAT_OPTIMAL )
268 {
269 assert(!*error);
270 *optimal = TRUE;
271 }
272 #ifdef SCIP_DEBUG
273 else
274 {
275 switch( lpsolstat )
276 {
277 case SCIP_LPSOLSTAT_ITERLIMIT:
278 SCIPdebugMsg(scip, " reached lp iteration limit\n");
279 break;
280 case SCIP_LPSOLSTAT_TIMELIMIT:
281 SCIPdebugMsg(scip, " reached time limit while solving lp\n");
282 break;
283 case SCIP_LPSOLSTAT_UNBOUNDEDRAY:
284 SCIPdebugMsg(scip, " lp was unbounded\n");
285 break;
286 case SCIP_LPSOLSTAT_NOTSOLVED:
287 SCIPdebugMsg(scip, " lp was not solved\n");
288 break;
289 case SCIP_LPSOLSTAT_ERROR:
290 SCIPdebugMsg(scip, " an error occured during solving lp\n");
291 break;
292 case SCIP_LPSOLSTAT_INFEASIBLE:
293 case SCIP_LPSOLSTAT_OBJLIMIT:
294 case SCIP_LPSOLSTAT_OPTIMAL: /* should not appear because it is handled earlier */
295 default:
296 SCIPdebugMsg(scip, " received an unexpected solstat during solving lp: %d\n", lpsolstat);
297 }
298 }
299 #endif
300
301 return SCIP_OKAY;
302 }
303
304 /** adds the objective cutoff to the LP; must be in probing mode */
305 static
addObjCutoff(SCIP * scip,SCIP_PROPDATA * propdata)306 SCIP_RETCODE addObjCutoff(
307 SCIP* scip, /**< SCIP data structure */
308 SCIP_PROPDATA* propdata /**< data of the obbt propagator */
309 )
310 {
311 SCIP_ROW* row;
312 SCIP_VAR** vars;
313 char rowname[SCIP_MAXSTRLEN];
314
315 int nvars;
316 int i;
317
318 assert(scip != NULL);
319 assert(SCIPinProbing(scip));
320 assert(propdata != NULL);
321 assert(propdata->cutoffrow == NULL);
322
323 if( SCIPisInfinity(scip, SCIPgetCutoffbound(scip)) )
324 {
325 SCIPdebugMsg(scip, "no objective cutoff since there is no cutoff bound\n");
326 return SCIP_OKAY;
327 }
328
329 SCIPdebugMsg(scip, "create objective cutoff and add it to the LP\n");
330
331 /* get variables data */
332 SCIP_CALL( SCIPgetVarsData(scip, &vars, &nvars, NULL, NULL, NULL, NULL) );
333
334 /* create objective cutoff row; set local flag to FALSE since primal cutoff is globally valid */
335 (void) SCIPsnprintf(rowname, SCIP_MAXSTRLEN, "obbt_objcutoff");
336 SCIP_CALL( SCIPcreateEmptyRowUnspec(scip, &row, rowname, -SCIPinfinity(scip), SCIPgetCutoffbound(scip), FALSE, FALSE, FALSE) );
337 SCIP_CALL( SCIPcacheRowExtensions(scip, row) );
338
339 for( i = 0; i < nvars; i++ )
340 {
341 SCIP_CALL( SCIPaddVarToRow(scip, row, vars[i], SCIPvarGetObj(vars[i])) );
342 }
343 SCIP_CALL( SCIPflushRowExtensions(scip, row) );
344
345 /* add row to the LP */
346 SCIP_CALL( SCIPaddRowProbing(scip, row) );
347
348 propdata->cutoffrow = row;
349 assert(SCIProwIsInLP(propdata->cutoffrow));
350
351 return SCIP_OKAY;
352 }
353
354 /** determines, whether a variable is already locally fixed */
355 static
varIsFixedLocal(SCIP * scip,SCIP_VAR * var)356 SCIP_Bool varIsFixedLocal(
357 SCIP* scip, /**< SCIP data structure */
358 SCIP_VAR* var /**< variable to check */
359 )
360 {
361 return SCIPisFeasEQ(scip, SCIPvarGetLbLocal(var), SCIPvarGetUbLocal(var));
362 }
363
364 /** sets objective to minimize or maximize a single variable */
365 static
setObjProbing(SCIP * scip,SCIP_PROPDATA * propdata,BOUND * bound,SCIP_Real coef)366 SCIP_RETCODE setObjProbing(
367 SCIP* scip,
368 SCIP_PROPDATA* propdata,
369 BOUND* bound,
370 SCIP_Real coef
371 )
372 {
373 #ifdef SCIP_DEBUG
374 SCIP_VAR** vars;
375 int nvars;
376 int counter;
377 int i;
378 #endif
379
380 assert( scip != NULL );
381 assert( propdata != NULL );
382 assert( bound != NULL );
383
384 /* set the objective for bound->var */
385 if( bound->boundtype == SCIP_BOUNDTYPE_LOWER )
386 {
387 SCIP_CALL( SCIPchgVarObjProbing(scip, bound->var, coef) );
388 }
389 else
390 {
391 SCIP_CALL( SCIPchgVarObjProbing(scip, bound->var, -coef) );
392 }
393
394 #ifdef SCIP_DEBUG
395 vars = SCIPgetVars(scip);
396 nvars = SCIPgetNVars(scip);
397 counter = 0;
398
399 for( i = 0; i < nvars; ++i )
400 {
401 if( SCIPgetVarObjProbing(scip, vars[i]) != 0.0 )
402 ++counter;
403 }
404
405 assert((counter == 0 && coef == 0.0) || (counter == 1 && coef != 0.0));
406 #endif
407
408 return SCIP_OKAY;
409 }
410
411 /** determines whether variable should be included in the right-hand side of the generalized variable bound */
412 static
includeVarGenVBound(SCIP * scip,SCIP_VAR * var)413 SCIP_Bool includeVarGenVBound(
414 SCIP* scip, /**< SCIP data structure */
415 SCIP_VAR* var /**< variable to check */
416 )
417 {
418 SCIP_Real redcost;
419
420 assert(scip != NULL);
421 assert(var != NULL);
422
423 if( SCIPvarGetStatus(var) != SCIP_VARSTATUS_COLUMN )
424 return FALSE;
425
426 redcost = SCIPgetVarRedcost(scip, var);
427 assert(redcost != SCIP_INVALID); /*lint !e777 */
428
429 if( redcost == SCIP_INVALID ) /*lint !e777 */
430 return FALSE;
431
432 if( redcost < SCIPdualfeastol(scip) && redcost > -SCIPdualfeastol(scip) )
433 return FALSE;
434
435 return TRUE;
436 }
437
438 /** returns number of LP iterations left (-1: no limit ) */
439 static
getIterationsLeft(SCIP * scip,SCIP_Longint nolditerations,SCIP_Longint itlimit)440 int getIterationsLeft(
441 SCIP* scip, /**< SCIP data structure */
442 SCIP_Longint nolditerations, /**< iterations count at the beginning of the corresponding function */
443 SCIP_Longint itlimit /**< LP iteration limit (-1: no limit) */
444 )
445 {
446 SCIP_Longint itsleft;
447
448 assert(scip != NULL);
449 assert(nolditerations >= 0);
450 assert(itlimit == -1 || itlimit >= 0);
451
452 if( itlimit == -1 )
453 {
454 SCIPdebugMsg(scip, "iterations left: unlimited\n");
455 return -1;
456 }
457 else
458 {
459 itsleft = itlimit - ( SCIPgetNLPIterations(scip) - nolditerations );
460 itsleft = MAX(itsleft, 0);
461 itsleft = MIN(itsleft, INT_MAX);
462
463 SCIPdebugMsg(scip, "iterations left: %d\n", (int) itsleft);
464 return (int) itsleft;
465 }
466 }
467
468 /** returns the objective coefficient for a variable's bound that will be chosen during filtering */
469 static
getFilterCoef(SCIP * scip,SCIP_PROPDATA * propdata,SCIP_VAR * var,SCIP_BOUNDTYPE boundtype)470 SCIP_Real getFilterCoef(
471 SCIP* scip, /**< SCIP data structure */
472 SCIP_PROPDATA* propdata, /**< data of the obbt propagator */
473 SCIP_VAR* var, /**< variable */
474 SCIP_BOUNDTYPE boundtype /**< boundtype to be filtered? */
475 )
476 {
477 SCIP_Real lb;
478 SCIP_Real ub;
479
480 assert(scip != NULL);
481 assert(propdata != NULL);
482 assert(var != NULL);
483
484 lb = SCIPvarGetLbLocal(var);
485 ub = SCIPvarGetUbLocal(var);
486
487 /* this function should not be called for fixed variables */
488 assert(!varIsFixedLocal(scip, var));
489
490 /* infinite bounds will not be reached */
491 if( boundtype == SCIP_BOUNDTYPE_LOWER && SCIPisInfinity(scip, -lb) )
492 return 0.0;
493 if( boundtype == SCIP_BOUNDTYPE_UPPER && SCIPisInfinity(scip, ub) )
494 return 0.0;
495
496 if( propdata->normalize )
497 {
498 /* if the length of the domain is too large then the coefficient should be set to +/- 1.0 */
499 if( boundtype == SCIP_BOUNDTYPE_LOWER && SCIPisInfinity(scip, ub) )
500 return 1.0;
501 if( boundtype == SCIP_BOUNDTYPE_UPPER && SCIPisInfinity(scip, -lb) )
502 return -1.0;
503
504 /* otherwise the coefficient is +/- 1.0 / ( ub - lb ) */
505 return boundtype == SCIP_BOUNDTYPE_LOWER ? 1.0 / (ub - lb) : -1.0 / (ub - lb);
506 }
507 else
508 {
509 return boundtype == SCIP_BOUNDTYPE_LOWER ? 1.0 : -1.0;
510 }
511 }
512
513 /** creates a genvbound if the dual LP solution provides such information
514 *
515 * Consider the problem
516 *
517 * min { +/- x_i : obj * x <= z, lb <= Ax <= ub, l <= x <= u },
518 *
519 * where z is the current cutoff bound. Let (mu, nu, gamma, alpha, beta) >= 0 be the optimal solution of the dual of
520 * problem (P), where the variables correspond to the primal inequalities in the following way:
521 *
522 * Ax >= lb <-> mu
523 * -Ax >= -ub <-> nu
524 * -obj * x >= -z <-> gamma
525 * x >= l <-> alpha
526 * -x >= -u <-> beta
527 *
528 * Fixing these multipliers, by weak duality, we obtain the inequality
529 *
530 * +/- x_i >= lb*mu - ub*nu - z*gamma + l*alpha - u*beta
531 *
532 * that holds for all primal feasible points x with objective value at least z. Setting
533 *
534 * c = lb*mu - ub*nu, redcost_k = alpha_k - beta_k
535 *
536 * we obtain the inequality
537 *
538 * +/- x_i >= sum ( redcost_k * x_k ) + (-gamma) * cutoff_bound + c,
539 *
540 * that holds for all primal feasible points with objective value at least cutoff_bound. Therefore, the latter
541 * inequality can be added as a generalized variable bound.
542 */
543 static
createGenVBound(SCIP * scip,SCIP_PROPDATA * propdata,BOUND * bound,SCIP_Bool * found)544 SCIP_RETCODE createGenVBound(
545 SCIP* scip, /**< SCIP data structure */
546 SCIP_PROPDATA* propdata, /**< data of the obbt propagator */
547 BOUND* bound, /**< bound of x_i */
548 SCIP_Bool* found /**< pointer to store if we have found a non-trivial genvbound */
549 )
550 {
551 assert(scip != NULL);
552 assert(bound != NULL);
553 assert(propdata != NULL);
554 assert(propdata->genvboundprop != NULL);
555 assert(found != NULL);
556
557 *found = FALSE;
558
559 /* make sure we are in probing mode having an optimal LP solution */
560 assert(SCIPinProbing(scip));
561
562 assert(SCIPgetLPSolstat(scip) == SCIP_LPSOLSTAT_OPTIMAL);
563
564 /* only genvbounds created in the root node are globally valid
565 *
566 * note: depth changes to one if we use the probing mode to solve the obbt LPs
567 */
568 assert(SCIPgetDepth(scip) == 0 || (SCIPinProbing(scip) && SCIPgetDepth(scip) == 1));
569
570 SCIPdebugMsg(scip, " try to create a genvbound for <%s>...\n", SCIPvarGetName(bound->var));
571
572 /* a genvbound with a multiplier for x_i would not help us */
573 if( SCIPisZero(scip, SCIPgetVarRedcost(scip, bound->var)) )
574 {
575 SCIP_VAR** vars; /* global variables array */
576 SCIP_VAR** genvboundvars; /* genvbound variables array */
577
578 SCIP_VAR* xi; /* variable x_i */
579
580 SCIP_Real* genvboundcoefs; /* genvbound coefficients array */
581
582 SCIP_Real gamma_dual; /* dual multiplier of objective cutoff */
583
584 int k; /* variable for indexing global variables array */
585 int ncoefs; /* number of nonzero coefficients in genvbound */
586 int nvars; /* number of global variables */
587
588 /* set x_i */
589 xi = bound->var;
590
591 /* get variable data */
592 SCIP_CALL( SCIPgetVarsData(scip, &vars, &nvars, NULL, NULL, NULL, NULL) );
593
594 /* count nonzero coefficients in genvbound */
595 ncoefs = 0;
596 for( k = 0; k < nvars; k++ )
597 {
598 if( includeVarGenVBound(scip, vars[k]) )
599 {
600 assert(vars[k] != xi);
601 ncoefs++;
602 }
603 }
604
605 /* get dual multiplier for the objective cutoff (set to zero if there is no) */
606 if( propdata->cutoffrow == NULL )
607 {
608 gamma_dual = 0.0;
609 }
610 else
611 {
612 assert(!SCIPisInfinity(scip, SCIPgetCutoffbound(scip)));
613
614 /* note that the objective cutoff is of the form
615 * -inf <= obj * x <= cutoff_bound
616 * but we want the positive dual multiplier!
617 */
618 gamma_dual = -SCIProwGetDualsol(propdata->cutoffrow);
619
620 /* we need to treat gamma to be exactly 0 if it is below the dual feasibility tolerance, see #2914 */
621 if( EPSZ(gamma_dual, SCIPdualfeastol(scip)) )
622 gamma_dual = 0.0;
623 }
624
625 /* we need at least one nonzero coefficient or a nonzero dual multiplier for the objective cutoff */
626 if( ncoefs > 0 || gamma_dual != 0.0 )
627 {
628 SCIP_Bool addgenvbound; /* if everything is fine with the redcosts and the bounds, add the genvbound */
629 SCIP_Real c; /* helper variable to calculate constant term in genvbound */
630 int idx; /* variable for indexing genvbound's coefficients array */
631
632 /* add the bound if the bool is still TRUE after the loop */
633 addgenvbound = TRUE;
634
635 /* there should be no coefficient for x_i */
636 assert(SCIPisZero(scip, SCIPgetVarRedcost(scip, xi)));
637
638 /* allocate memory for storing the genvbounds right-hand side variables and coefficients */
639 SCIP_CALL( SCIPallocBufferArray(scip, &(genvboundvars), ncoefs) );
640 SCIP_CALL( SCIPallocBufferArray(scip, &(genvboundcoefs), ncoefs) );
641
642 /* set c = lb*mu - ub*nu - z*gamma + l*alpha - u*beta */
643 c = SCIPgetLPObjval(scip);
644
645 /* subtract ( - z * gamma ) from c */
646 c += SCIPgetCutoffbound(scip) * gamma_dual;
647
648 /* subtract ( l*alpha - u*beta ) from c and set the coefficients of the variables */
649 idx = 0;
650 for( k = 0; k < nvars; k++ )
651 {
652 SCIP_VAR* xk;
653
654 xk = vars[k];
655
656 if( includeVarGenVBound(scip, xk) )
657 {
658 SCIP_Real redcost;
659
660 redcost = SCIPgetVarRedcost(scip, xk);
661
662 assert(redcost != SCIP_INVALID); /*lint !e777 */
663 assert(xk != xi);
664
665 /* in this case dont add a genvbound */
666 if( ( (redcost > SCIPdualfeastol(scip)) && SCIPisInfinity(scip, -SCIPvarGetLbLocal(xk)) ) ||
667 ( (redcost < -SCIPdualfeastol(scip)) && SCIPisInfinity(scip, SCIPvarGetUbLocal(xk)) ) )
668 {
669 addgenvbound = FALSE;
670 break;
671 }
672
673 /* store coefficients */
674 assert(idx < ncoefs);
675 genvboundvars[idx] = xk;
676 genvboundcoefs[idx] = redcost;
677 idx++;
678
679 /* if redcost > 0, then redcost = alpha_k, otherwise redcost = - beta_k */
680 assert(redcost <= 0 || !SCIPisInfinity(scip, -SCIPvarGetLbLocal(xk)));
681 assert(redcost >= 0 || !SCIPisInfinity(scip, SCIPvarGetUbLocal(xk)));
682 c -= redcost > 0 ? redcost * SCIPvarGetLbLocal(xk) : redcost * SCIPvarGetUbLocal(xk);
683 }
684 }
685
686 assert(!addgenvbound || idx == ncoefs);
687
688 /* add genvbound */
689 if( addgenvbound && !SCIPisInfinity(scip, -c) )
690 {
691 #ifndef NDEBUG
692 /* check whether the activity of the LVB in the optimal solution of the LP is equal to the LP objective value */
693 SCIP_Real activity = c - gamma_dual * SCIPgetCutoffbound(scip);
694
695 for( k = 0; k < ncoefs; ++k )
696 activity += genvboundcoefs[k] * SCIPvarGetLPSol(genvboundvars[k]);
697
698 SCIPdebugMsg(scip, "LVB activity = %g lpobj = %g\n", activity, SCIPgetLPObjval(scip));
699 assert(EPSZ(SCIPrelDiff(activity, SCIPgetLPObjval(scip)), 10.0 * SCIPdualfeastol(scip)));
700 #endif
701
702 SCIPdebugMsg(scip, " adding genvbound\n");
703 SCIP_CALL( SCIPgenVBoundAdd(scip, propdata->genvboundprop, genvboundvars, xi, genvboundcoefs, ncoefs,
704 gamma_dual < SCIPdualfeastol(scip) ? 0.0 : -gamma_dual, c, bound->boundtype) );
705 *found = TRUE;
706 }
707
708 /* free arrays */
709 SCIPfreeBufferArray(scip, &genvboundcoefs);
710 SCIPfreeBufferArray(scip, &genvboundvars);
711 }
712 else
713 {
714 SCIPdebugMsg(scip, " trivial genvbound, skipping\n");
715 }
716 }
717 else
718 {
719 SCIPdebugMsg(scip, " found multiplier for <%s>: %g, skipping\n",
720 SCIPvarGetName(bound->var), SCIPgetVarRedcost(scip, bound->var));
721 }
722
723 return SCIP_OKAY;
724 }
725
726 /** exchange a bound which has been processed and updates the last undone and unfiltered bound index
727 * NOTE: this method has to be called after filtering or processing a bound
728 */
729 static
exchangeBounds(SCIP_PROPDATA * propdata,int i)730 void exchangeBounds(
731 SCIP_PROPDATA* propdata, /**< propagator data */
732 int i /**< bound that was filtered or processed */
733 )
734 {
735 assert(i >= 0 && i < propdata->nbounds);
736 assert(propdata->lastidx >= 0 && propdata->lastidx < propdata->nbounds);
737
738 /* exchange the bounds */
739 if( propdata->lastidx != i )
740 {
741 BOUND* tmp;
742
743 tmp = propdata->bounds[i];
744 propdata->bounds[i] = propdata->bounds[propdata->lastidx];
745 propdata->bounds[propdata->lastidx] = tmp;
746 }
747
748 propdata->lastidx -= 1;
749 }
750
751 /** helper function to return a corner of the domain of two variables */
752 static
getCorner(SCIP_VAR * x,SCIP_VAR * y,CORNER corner,SCIP_Real * px,SCIP_Real * py)753 void getCorner(
754 SCIP_VAR* x, /**< first variable */
755 SCIP_VAR* y, /**< second variable */
756 CORNER corner, /**< corner */
757 SCIP_Real* px, /**< buffer to store point for x */
758 SCIP_Real* py /**< buffer to store point for y */
759 )
760 {
761 assert(x != NULL);
762 assert(y != NULL);
763 assert(px != NULL);
764 assert(py != NULL);
765
766 switch( corner )
767 {
768 case LEFTBOTTOM:
769 *px = SCIPvarGetLbGlobal(x);
770 *py = SCIPvarGetLbGlobal(y);
771 break;
772 case RIGHTBOTTOM:
773 *px = SCIPvarGetUbGlobal(x);
774 *py = SCIPvarGetLbGlobal(y);
775 break;
776 case LEFTTOP:
777 *px = SCIPvarGetLbGlobal(x);
778 *py = SCIPvarGetUbGlobal(y);
779 break;
780 case RIGHTTOP:
781 *px = SCIPvarGetUbGlobal(x);
782 *py = SCIPvarGetUbGlobal(y);
783 break;
784 case FILTERED:
785 SCIPABORT();
786 }
787 }
788
789 /** helper function to return the two end points of a diagonal */
790 static
getCorners(SCIP_VAR * x,SCIP_VAR * y,CORNER corner,SCIP_Real * xs,SCIP_Real * ys,SCIP_Real * xt,SCIP_Real * yt)791 void getCorners(
792 SCIP_VAR* x, /**< first variable */
793 SCIP_VAR* y, /**< second variable */
794 CORNER corner, /**< corner */
795 SCIP_Real* xs, /**< buffer to store start point for x */
796 SCIP_Real* ys, /**< buffer to store start point for y */
797 SCIP_Real* xt, /**< buffer to store end point for x */
798 SCIP_Real* yt /**< buffer to store end point for y */
799 )
800 {
801 assert(x != NULL);
802 assert(y != NULL);
803 assert(xs != NULL);
804 assert(ys != NULL);
805 assert(xt != NULL);
806 assert(yt != NULL);
807
808 /* get end point */
809 getCorner(x,y, corner, xt, yt);
810
811 /* get start point */
812 switch( corner )
813 {
814 case LEFTBOTTOM:
815 getCorner(x,y, RIGHTTOP, xs, ys);
816 break;
817 case RIGHTBOTTOM:
818 getCorner(x,y, LEFTTOP, xs, ys);
819 break;
820 case LEFTTOP:
821 getCorner(x,y, RIGHTBOTTOM, xs, ys);
822 break;
823 case RIGHTTOP:
824 getCorner(x,y, LEFTBOTTOM, xs, ys);
825 break;
826 case FILTERED:
827 SCIPABORT();
828 }
829 }
830
831 /** trying to filter some bounds using the existing LP solution */
832 static
filterExistingLP(SCIP * scip,SCIP_PROPDATA * propdata,int * nfiltered,BOUND * currbound)833 SCIP_RETCODE filterExistingLP(
834 SCIP* scip, /**< original SCIP data structure */
835 SCIP_PROPDATA* propdata, /**< data of the obbt propagator */
836 int* nfiltered, /**< how many bounds were filtered this round? */
837 BOUND* currbound /**< bound for which OBBT LP was solved (Note: might be NULL) */
838 )
839 {
840 int i;
841
842 assert(scip != NULL);
843 assert(propdata != NULL);
844 assert(nfiltered != NULL);
845
846 *nfiltered = 0;
847
848 /* only apply filtering if an LP solution is at hand */
849 if( SCIPgetLPSolstat(scip) != SCIP_LPSOLSTAT_OPTIMAL )
850 {
851 SCIPdebugMsg(scip, "can't filter using existing lp solution since it was not solved to optimality\n");
852 return SCIP_OKAY;
853 }
854
855 /* check if a bound is tight */
856 for( i = propdata->nbounds - 1; i >= 0; --i )
857 {
858 BOUND* bound; /* shortcut for current bound */
859
860 SCIP_Real solval; /* the variables value in the current solution */
861 SCIP_Real boundval; /* current local bound for the variable */
862
863 bound = propdata->bounds[i];
864 if( bound->filtered || bound->done )
865 continue;
866
867 boundval = bound->boundtype == SCIP_BOUNDTYPE_UPPER ?
868 SCIPvarGetUbLocal(bound->var) : SCIPvarGetLbLocal(bound->var);
869 solval = SCIPvarGetLPSol(bound->var);
870
871 /* bound is tight; since this holds for all fixed variables, those are filtered here automatically; if the lp solution
872 * is infinity, then also the bound is tight */
873 if( (bound->boundtype == SCIP_BOUNDTYPE_UPPER &&
874 (SCIPisInfinity(scip, solval) || SCIPisFeasGE(scip, solval, boundval)))
875 || (bound->boundtype == SCIP_BOUNDTYPE_LOWER &&
876 (SCIPisInfinity(scip, -solval) || SCIPisFeasLE(scip, solval, boundval))) )
877 {
878 SCIP_BASESTAT basestat;
879
880 /* mark bound as filtered */
881 bound->filtered = TRUE;
882 SCIPdebugMsg(scip, "trivial filtered var: %s boundval=%e solval=%e\n", SCIPvarGetName(bound->var), boundval, solval);
883
884 /* get the basis status of the variable */
885 basestat = SCIPcolGetBasisStatus(SCIPvarGetCol(bound->var));
886
887 /* solve corresponding OBBT LP and try to generate a nontrivial genvbound */
888 if( propdata->genvbdsduringfilter && currbound != NULL && basestat == SCIP_BASESTAT_BASIC )
889 {
890 #ifndef NDEBUG
891 int j;
892 #endif
893 SCIP_Bool optimal;
894 SCIP_Bool error;
895
896 /* set objective coefficient of the bound */
897 SCIP_CALL( SCIPchgVarObjProbing(scip, currbound->var, 0.0) );
898 SCIP_CALL( setObjProbing(scip, propdata, bound, 1.0) );
899
900 #ifndef NDEBUG
901 for( j = 0; j < SCIPgetNVars(scip); ++j )
902 {
903 SCIP_VAR* var;
904
905 var = SCIPgetVars(scip)[j];
906 assert(var != NULL);
907 assert(SCIPisZero(scip, SCIPgetVarObjProbing(scip, var)) || var == bound->var);
908 }
909 #endif
910
911 /* solve the OBBT LP */
912 SCIP_CALL( solveLP(scip, -1, &error, &optimal) );
913
914 /* try to generate a genvbound if we have solved the OBBT LP */
915 if( optimal && propdata->genvboundprop != NULL
916 && (SCIPgetDepth(scip) == 0 || (SCIPinProbing(scip) && SCIPgetDepth(scip) == 1)) )
917 {
918 SCIP_Bool found;
919
920 assert(!error);
921 SCIP_CALL( createGenVBound(scip, propdata, bound, &found) );
922
923 SCIPdebugMsg(scip, "found genvbound during trivial filtering? %u\n", found);
924 } /*lint !e438*/
925
926 /* restore objective function */
927 SCIP_CALL( setObjProbing(scip, propdata, bound, 0.0) );
928 SCIP_CALL( setObjProbing(scip, propdata, currbound, 1.0) );
929 }
930
931 /* exchange bound i with propdata->bounds[propdata->lastidx] */
932 if( propdata->lastidx >= 0 )
933 exchangeBounds(propdata, i);
934
935 /* increase number of filtered variables */
936 (*nfiltered)++;
937 }
938 }
939
940 /* try to filter bilinear bounds */
941 for( i = propdata->lastbilinidx; i < propdata->nbilinbounds; ++i )
942 {
943 CORNER corners[4] = {LEFTTOP, LEFTBOTTOM, RIGHTTOP, RIGHTBOTTOM};
944 BILINBOUND* bilinbound = propdata->bilinbounds[i];
945 SCIP_Real solx;
946 SCIP_Real soly;
947 SCIPdebug(int oldfiltered;)
948 int j;
949
950 /* skip processed and filtered bounds */
951 if( bilinbound->done || bilinbound->filtered == FILTERED ) /*lint !e641*/
952 continue;
953
954 SCIPdebug(oldfiltered = bilinbound->filtered;)
955 solx = SCIPvarGetLPSol(bilinbound->x);
956 soly = SCIPvarGetLPSol(bilinbound->y);
957
958 /* check cases of unbounded solution values */
959 if( SCIPisInfinity(scip, solx) )
960 bilinbound->filtered = bilinbound->filtered | RIGHTTOP | RIGHTBOTTOM; /*lint !e641*/
961 else if( SCIPisInfinity(scip, -solx) )
962 bilinbound->filtered = bilinbound->filtered | LEFTTOP | LEFTBOTTOM; /*lint !e641*/
963
964 if( SCIPisInfinity(scip, soly) )
965 bilinbound->filtered = bilinbound->filtered | RIGHTTOP | LEFTTOP; /*lint !e641*/
966 else if( SCIPisInfinity(scip, -soly) )
967 bilinbound->filtered = bilinbound->filtered | RIGHTBOTTOM | LEFTBOTTOM; /*lint !e641*/
968
969 /* check all corners */
970 for( j = 0; j < 4; ++j )
971 {
972 SCIP_Real xt = SCIP_INVALID;
973 SCIP_Real yt = SCIP_INVALID;
974
975 getCorner(bilinbound->x, bilinbound->y, corners[j], &xt, &yt);
976
977 if( (SCIPisInfinity(scip, REALABS(solx)) || SCIPisFeasEQ(scip, xt, solx))
978 && (SCIPisInfinity(scip, REALABS(soly)) || SCIPisFeasEQ(scip, yt, soly)) )
979 bilinbound->filtered = bilinbound->filtered | corners[j]; /*lint !e641*/
980 }
981
982 #ifdef SCIP_DEBUG
983 if( oldfiltered != bilinbound->filtered )
984 {
985 SCIP_VAR* x = bilinbound->x;
986 SCIP_VAR* y = bilinbound->y;
987 SCIPdebugMessage("filtered corners %d for (%s,%s) = (%g,%g) in [%g,%g]x[%g,%g]\n",
988 bilinbound->filtered - oldfiltered, SCIPvarGetName(x), SCIPvarGetName(y), solx, soly,
989 SCIPvarGetLbGlobal(x), SCIPvarGetUbGlobal(x), SCIPvarGetLbGlobal(y), SCIPvarGetUbGlobal(y));
990 }
991 #endif
992 }
993
994 return SCIP_OKAY;
995 }
996
997 /** enforces one round of filtering */
998 static
filterRound(SCIP * scip,SCIP_PROPDATA * propdata,int itlimit,int * nfiltered,SCIP_Real * objcoefs,int * objcoefsinds,int nobjcoefs)999 SCIP_RETCODE filterRound(
1000 SCIP* scip, /**< SCIP data structure */
1001 SCIP_PROPDATA* propdata, /**< data of the obbt propagator */
1002 int itlimit, /**< LP iteration limit (-1: no limit) */
1003 int* nfiltered, /**< how many bounds were filtered this round */
1004 SCIP_Real* objcoefs, /**< array to store the nontrivial objective coefficients */
1005 int* objcoefsinds, /**< array to store bound indices for which their corresponding variables
1006 * has a nontrivial objective coefficient */
1007 int nobjcoefs /**< number of nontrivial objective coefficients */
1008 )
1009 {
1010 SCIP_VAR** vars; /* array of the problems variables */
1011 SCIP_Bool error;
1012 SCIP_Bool optimal;
1013
1014 int nvars; /* number of the problems variables */
1015 int i;
1016
1017 assert(scip != NULL);
1018 assert(SCIPinProbing(scip));
1019 assert(propdata != NULL);
1020 assert(itlimit == -1 || itlimit >= 0);
1021 assert(nfiltered != NULL);
1022 assert(objcoefs != NULL);
1023 assert(objcoefsinds != NULL);
1024 assert(nobjcoefs >= 0);
1025
1026 *nfiltered = 0;
1027
1028 /* get variable data */
1029 SCIP_CALL( SCIPgetVarsData(scip, &vars, &nvars, NULL, NULL, NULL, NULL) );
1030
1031 /* solve LP */
1032 SCIP_CALL( solveLP(scip, itlimit, &error, &optimal) );
1033
1034 if( !optimal )
1035 {
1036 SCIPdebugMsg(scip, "skipping filter round since the LP was not solved to optimality\n");
1037 return SCIP_OKAY;
1038 }
1039
1040 assert(!error);
1041
1042 /* check if a bound is tight */
1043 for( i = 0; i < propdata->nbounds; i++ )
1044 {
1045 BOUND* bound; /* shortcut for current bound */
1046
1047 SCIP_Real solval; /* the variables value in the current solution */
1048 SCIP_Real boundval; /* current local bound for the variable */
1049
1050 bound = propdata->bounds[i];
1051
1052 /* if bound is filtered it was handled already before */
1053 if( bound->filtered )
1054 continue;
1055
1056 boundval = bound->boundtype == SCIP_BOUNDTYPE_UPPER ?
1057 SCIPvarGetUbLocal(bound->var) : SCIPvarGetLbLocal(bound->var);
1058 solval = SCIPvarGetLPSol(bound->var);
1059
1060 /* bound is tight */
1061 if( (bound->boundtype == SCIP_BOUNDTYPE_UPPER && SCIPisFeasGE(scip, solval, boundval))
1062 || (bound->boundtype == SCIP_BOUNDTYPE_LOWER && SCIPisFeasLE(scip, solval, boundval)) )
1063 {
1064 SCIP_Real objcoef;
1065 SCIP_BASESTAT basestat;
1066
1067 /* mark bound as filtered */
1068 bound->filtered = TRUE;
1069
1070 /* get the basis status of the variable */
1071 basestat = SCIPcolGetBasisStatus(SCIPvarGetCol(bound->var));
1072
1073 /* increase number of filtered variables */
1074 (*nfiltered)++;
1075
1076 /* solve corresponding OBBT LP and try to generate a nontrivial genvbound */
1077 if( propdata->genvbdsduringfilter && basestat == SCIP_BASESTAT_BASIC )
1078 {
1079 int j;
1080
1081 /* set all objective coefficients to zero */
1082 for( j = 0; j < nobjcoefs; ++j )
1083 {
1084 BOUND* filterbound;
1085
1086 filterbound = propdata->bounds[ objcoefsinds[j] ];
1087 assert(filterbound != NULL);
1088
1089 SCIP_CALL( SCIPchgVarObjProbing(scip, filterbound->var, 0.0) );
1090 }
1091
1092 #ifndef NDEBUG
1093 for( j = 0; j < nvars; ++j )
1094 assert(SCIPisZero(scip, SCIPgetVarObjProbing(scip, vars[j])));
1095 #endif
1096
1097 /* set objective coefficient of the bound */
1098 SCIP_CALL( setObjProbing(scip, propdata, bound, 1.0) );
1099
1100 /* solve the OBBT LP */
1101 SCIP_CALL( solveLP(scip, -1, &error, &optimal) );
1102
1103 /* try to generate a genvbound if we have solved the OBBT LP */
1104 if( optimal && propdata->genvboundprop != NULL
1105 && (SCIPgetDepth(scip) == 0 || (SCIPinProbing(scip) && SCIPgetDepth(scip) == 1)) )
1106 {
1107 SCIP_Bool found;
1108
1109 assert(!error);
1110 SCIP_CALL( createGenVBound(scip, propdata, bound, &found) );
1111 SCIPdebugMsg(scip, "found genvbound during aggressive filtering? %u\n", found);
1112 } /*lint !e438*/
1113
1114 /* restore objective function */
1115 for( j = 0; j < nobjcoefs; ++j )
1116 {
1117 BOUND* filterbound;
1118
1119 filterbound = propdata->bounds[ objcoefsinds[j] ];
1120 assert(filterbound != NULL);
1121
1122 /* NOTE: only restore coefficients of nonfiltered bounds */
1123 if( !filterbound->filtered )
1124 {
1125 assert(!SCIPisZero(scip, objcoefs[j]));
1126 SCIP_CALL( SCIPchgVarObjProbing(scip, propdata->bounds[ objcoefsinds[j] ]->var, objcoefs[j]) );
1127 }
1128 }
1129 }
1130
1131 /* get the corresponding variable's objective coefficient */
1132 objcoef = SCIPgetVarObjProbing(scip, bound->var);
1133
1134 /* change objective coefficient if it was set up for this bound */
1135 if( (bound->boundtype == SCIP_BOUNDTYPE_UPPER && SCIPisNegative(scip, objcoef))
1136 || (bound->boundtype == SCIP_BOUNDTYPE_LOWER && SCIPisPositive(scip, objcoef)) )
1137 {
1138 SCIP_CALL( SCIPchgVarObjProbing(scip, bound->var, 0.0) );
1139 }
1140 }
1141 }
1142
1143 return SCIP_OKAY;
1144 }
1145
1146 /** filter some bounds that are not improvable by solving auxiliary LPs */
1147 static
filterBounds(SCIP * scip,SCIP_PROPDATA * propdata,SCIP_Longint itlimit)1148 SCIP_RETCODE filterBounds(
1149 SCIP* scip, /**< SCIP data structure */
1150 SCIP_PROPDATA* propdata, /**< data of the obbt propagator */
1151 SCIP_Longint itlimit /**< LP iteration limit (-1: no limit) */
1152 )
1153 {
1154 SCIP_VAR** vars;
1155 SCIP_Longint nolditerations;
1156 SCIP_Real* objcoefs; /* array to store the nontrivial objective coefficients */
1157 int* objcoefsinds; /* array to store bound indices for which the corresponding variable
1158 * has a nontrivial objective coefficient */
1159 int nobjcoefs; /* number of nontrivial objective coefficients */
1160 int nleftiterations;
1161 int i;
1162 int nfiltered;
1163 int ntotalfiltered;
1164 int nvars;
1165
1166 assert(scip != NULL);
1167 assert(SCIPinProbing(scip));
1168 assert(propdata != NULL);
1169 assert(itlimit == -1 || itlimit >= 0);
1170
1171 ntotalfiltered = 0;
1172 nolditerations = SCIPgetNLPIterations(scip);
1173 nleftiterations = getIterationsLeft(scip, nolditerations, itlimit);
1174
1175 /* get variable data */
1176 SCIP_CALL( SCIPgetVarsData(scip, &vars, &nvars, NULL, NULL, NULL, NULL) );
1177
1178 SCIPdebugMsg(scip, "start filter rounds\n");
1179
1180 SCIP_CALL( SCIPallocBufferArray(scip, &objcoefs, propdata->nbounds) );
1181 SCIP_CALL( SCIPallocBufferArray(scip, &objcoefsinds, propdata->nbounds) );
1182 nobjcoefs = 0;
1183
1184 /*
1185 * 1.) Try first to filter lower bounds of interesting variables, whose bounds are not already filtered
1186 */
1187
1188 for( i = 0; i < nvars; i++ )
1189 {
1190 SCIP_CALL( SCIPchgVarObjProbing(scip, vars[i], 0.0) );
1191 }
1192
1193 for( i = 0; i < propdata->nbounds; i++ )
1194 {
1195 if( propdata->bounds[i]->boundtype == SCIP_BOUNDTYPE_LOWER && !propdata->bounds[i]->filtered
1196 && !propdata->bounds[i]->done )
1197 {
1198 SCIP_Real objcoef;
1199
1200 objcoef = getFilterCoef(scip, propdata, propdata->bounds[i]->var, SCIP_BOUNDTYPE_LOWER);
1201
1202 if( !SCIPisZero(scip, objcoef) )
1203 {
1204 SCIP_CALL( SCIPchgVarObjProbing(scip, propdata->bounds[i]->var, objcoef) );
1205
1206 /* store nontrivial objective coefficients */
1207 objcoefs[nobjcoefs] = objcoef;
1208 objcoefsinds[nobjcoefs] = i;
1209 ++nobjcoefs;
1210 }
1211 }
1212 }
1213
1214 do
1215 {
1216 SCIPdebugMsg(scip, "doing a lower bounds round\n");
1217 SCIP_CALL( filterRound(scip, propdata, nleftiterations, &nfiltered, objcoefs, objcoefsinds, nobjcoefs) );
1218 ntotalfiltered += nfiltered;
1219 SCIPdebugMsg(scip, "filtered %d more bounds in lower bounds round\n", nfiltered);
1220
1221 /* update iterations left */
1222 nleftiterations = getIterationsLeft(scip, nolditerations, itlimit);
1223 }
1224 while( nfiltered >= propdata->nminfilter && ( nleftiterations == -1 || nleftiterations > 0 ) );
1225
1226 /*
1227 * 2.) Now try to filter the remaining upper bounds of interesting variables, whose bounds are not already filtered
1228 */
1229
1230 /* set all objective coefficients to zero */
1231 for( i = 0; i < nobjcoefs; i++ )
1232 {
1233 BOUND* bound;
1234
1235 assert(objcoefsinds[i] >= 0 && objcoefsinds[i] < propdata->nbounds);
1236 bound = propdata->bounds[ objcoefsinds[i] ];
1237 assert(bound != NULL);
1238 SCIP_CALL( SCIPchgVarObjProbing(scip, bound->var, 0.0) );
1239 }
1240
1241 /* reset number of nontrivial objective coefficients */
1242 nobjcoefs = 0;
1243
1244 #ifndef NDEBUG
1245 for( i = 0; i < nvars; ++i )
1246 assert(SCIPisZero(scip, SCIPgetVarObjProbing(scip, vars[i])));
1247 #endif
1248
1249 for( i = 0; i < propdata->nbounds; i++ )
1250 {
1251 if( propdata->bounds[i]->boundtype == SCIP_BOUNDTYPE_UPPER && !propdata->bounds[i]->filtered )
1252 {
1253 SCIP_Real objcoef;
1254
1255 objcoef = getFilterCoef(scip, propdata, propdata->bounds[i]->var, SCIP_BOUNDTYPE_UPPER);
1256
1257 if( !SCIPisZero(scip, objcoef) )
1258 {
1259 SCIP_CALL( SCIPchgVarObjProbing(scip, propdata->bounds[i]->var, objcoef) );
1260
1261 /* store nontrivial objective coefficients */
1262 objcoefs[nobjcoefs] = objcoef;
1263 objcoefsinds[nobjcoefs] = i;
1264 ++nobjcoefs;
1265 }
1266 }
1267 }
1268
1269 do
1270 {
1271 SCIPdebugMsg(scip, "doing an upper bounds round\n");
1272 SCIP_CALL( filterRound(scip, propdata, nleftiterations, &nfiltered, objcoefs, objcoefsinds, nobjcoefs) );
1273 SCIPdebugMsg(scip, "filtered %d more bounds in upper bounds round\n", nfiltered);
1274 ntotalfiltered += nfiltered;
1275 /* update iterations left */
1276 nleftiterations = getIterationsLeft(scip, nolditerations, itlimit);
1277 }
1278 while( nfiltered >= propdata->nminfilter && ( nleftiterations == -1 || nleftiterations > 0 ) );
1279
1280 SCIPdebugMsg(scip, "filtered %d this round\n", ntotalfiltered);
1281
1282 /* free array */
1283 SCIPfreeBufferArray(scip, &objcoefsinds);
1284 SCIPfreeBufferArray(scip, &objcoefs);
1285
1286 return SCIP_OKAY;
1287 }
1288
1289 /** applies possible bound changes that were found */
1290 static
applyBoundChgs(SCIP * scip,SCIP_PROPDATA * propdata,SCIP_RESULT * result)1291 SCIP_RETCODE applyBoundChgs(
1292 SCIP* scip, /**< SCIP data structure */
1293 SCIP_PROPDATA* propdata, /**< data of the obbt propagator */
1294 SCIP_RESULT* result /**< result pointer */
1295 )
1296 {
1297 #ifdef SCIP_DEBUG
1298 int ntightened; /* stores the number of successful bound changes */
1299 #endif
1300 int i;
1301
1302 assert(scip != NULL);
1303 assert(!SCIPinProbing(scip));
1304 assert(propdata != NULL);
1305 assert(result != NULL);
1306 assert(*result == SCIP_DIDNOTFIND);
1307
1308 SCIPdebug( ntightened = 0 );
1309
1310 for( i = 0; i < propdata->nbounds; i++ )
1311 {
1312 BOUND* bound; /* shortcut to the current bound */
1313 SCIP_Bool infeas; /* stores wether a tightening approach forced an infeasibilty */
1314 SCIP_Bool tightened; /* stores wether a tightening approach was successful */
1315
1316 bound = propdata->bounds[i];
1317
1318 if( bound->found )
1319 {
1320 SCIPdebug( double oldbound = (bound->boundtype == SCIP_BOUNDTYPE_LOWER)
1321 ? SCIPvarGetLbLocal(bound->var)
1322 : SCIPvarGetUbLocal(bound->var) );
1323
1324 if( bound->boundtype == SCIP_BOUNDTYPE_LOWER )
1325 {
1326 SCIP_CALL( SCIPtightenVarLb(scip, bound->var, bound->newval, FALSE, &infeas, &tightened) );
1327 }
1328 else
1329 {
1330 SCIP_CALL( SCIPtightenVarUb(scip, bound->var, bound->newval, FALSE, &infeas, &tightened) );
1331 }
1332
1333 /* handle information about the success */
1334 if( infeas )
1335 {
1336 *result = SCIP_CUTOFF;
1337 SCIPdebugMsg(scip, "cut off\n");
1338 break;
1339 }
1340
1341 if( tightened )
1342 {
1343 SCIPdebug( SCIPdebugMsg(scip, "tightended: %s old: %e new: %e\n" , SCIPvarGetName(bound->var), oldbound,
1344 bound->newval) );
1345 *result = SCIP_REDUCEDDOM;
1346 SCIPdebug( ntightened++ );
1347 }
1348 }
1349 }
1350
1351 SCIPdebug( SCIPdebugMsg(scip, "tightened bounds: %d\n", ntightened) );
1352
1353 return SCIP_OKAY;
1354 }
1355
1356 /** tries to tighten a bound in probing mode */
1357 static
tightenBoundProbing(SCIP * scip,BOUND * bound,SCIP_Real newval,SCIP_Bool * tightened)1358 SCIP_RETCODE tightenBoundProbing(
1359 SCIP* scip, /**< SCIP data structure */
1360 BOUND* bound, /**< bound that could be tightened */
1361 SCIP_Real newval, /**< new bound value */
1362 SCIP_Bool* tightened /**< was tightening successful? */
1363 )
1364 {
1365 SCIP_Real lb;
1366 SCIP_Real ub;
1367
1368 assert(scip != NULL);
1369 assert(SCIPinProbing(scip));
1370 assert(bound != NULL);
1371 assert(tightened != NULL);
1372
1373 *tightened = FALSE;
1374
1375 /* get old bounds */
1376 lb = SCIPvarGetLbLocal(bound->var);
1377 ub = SCIPvarGetUbLocal(bound->var);
1378
1379 if( bound->boundtype == SCIP_BOUNDTYPE_LOWER )
1380 {
1381 /* round bounds new value if variable is integral */
1382 if( SCIPvarIsIntegral(bound->var) )
1383 newval = SCIPceil(scip, newval);
1384
1385 /* ensure that we give consistent bounds to the LP solver */
1386 if( newval > ub )
1387 newval = ub;
1388
1389 /* tighten if really better */
1390 if( SCIPisLbBetter(scip, newval, lb, ub) )
1391 {
1392 SCIP_CALL( SCIPchgVarLbProbing(scip, bound->var, newval) );
1393 *tightened = TRUE;
1394 }
1395 }
1396 else
1397 {
1398 /* round bounds new value if variable is integral */
1399 if( SCIPvarIsIntegral(bound->var) )
1400 newval = SCIPfloor(scip, newval);
1401
1402 /* ensure that we give consistent bounds to the LP solver */
1403 if( newval < lb )
1404 newval = lb;
1405
1406 /* tighten if really better */
1407 if( SCIPisUbBetter(scip, newval, lb, ub) )
1408 {
1409 SCIP_CALL( SCIPchgVarUbProbing(scip, bound->var, newval) );
1410 *tightened = TRUE;
1411 }
1412 }
1413
1414 return SCIP_OKAY;
1415 }
1416
1417 /** comparison method for two bounds w.r.t. their scores */
1418 static
SCIP_DECL_SORTPTRCOMP(compBoundsScore)1419 SCIP_DECL_SORTPTRCOMP(compBoundsScore)
1420 {
1421 BOUND* bound1 = (BOUND*) elem1;
1422 BOUND* bound2 = (BOUND*) elem2;
1423
1424 return bound1->score == bound2->score ? 0 : ( bound1->score > bound2->score ? 1 : -1 );
1425 }
1426
1427 /** comparison method for two bilinear term bounds w.r.t. their scores */
1428 static
SCIP_DECL_SORTPTRCOMP(compBilinboundsScore)1429 SCIP_DECL_SORTPTRCOMP(compBilinboundsScore)
1430 {
1431 BILINBOUND* bound1 = (BILINBOUND*) elem1;
1432 BILINBOUND* bound2 = (BILINBOUND*) elem2;
1433
1434 return bound1->score == bound2->score ? 0 : ( bound1->score > bound2->score ? 1 : -1 ); /*lint !e777*/
1435 }
1436
1437 /** comparison method for two bounds w.r.t. their boundtype */
1438 static
SCIP_DECL_SORTPTRCOMP(compBoundsBoundtype)1439 SCIP_DECL_SORTPTRCOMP(compBoundsBoundtype)
1440 {
1441 int diff;
1442 BOUND* bound1 = (BOUND*) elem1;
1443 BOUND* bound2 = (BOUND*) elem2;
1444
1445 /* prioritize undone bounds */
1446 diff = (!bound1->done ? 1 : 0) - (!bound2->done ? 1 : 0);
1447 if( diff != 0 )
1448 return diff;
1449
1450 /* prioritize unfiltered bounds */
1451 diff = (!bound1->filtered ? 1 : 0) - (!bound2->filtered ? 1 : 0);
1452 if( diff != 0 )
1453 return diff;
1454
1455 diff = (bound1->boundtype == SCIP_BOUNDTYPE_LOWER ? 1 : 0) - (bound2->boundtype == SCIP_BOUNDTYPE_LOWER ? 1 : 0);
1456
1457 if( diff == 0 )
1458 return (bound1->score == bound2->score) ? 0 : (bound1->score > bound2->score ? 1 : -1);
1459 else
1460 return diff;
1461 }
1462
1463 /** sort the propdata->bounds array with their distance or their boundtype key */
1464 static
sortBounds(SCIP * scip,SCIP_PROPDATA * propdata)1465 SCIP_RETCODE sortBounds(
1466 SCIP* scip, /**< SCIP data structure */
1467 SCIP_PROPDATA* propdata /**< propagator data */
1468 )
1469 {
1470 assert(scip != NULL);
1471 assert(propdata != NULL);
1472
1473 SCIPdebugMsg(scip, "sort bounds\n");
1474 SCIPsortDownPtr((void**) propdata->bounds, compBoundsBoundtype, propdata->nbounds);
1475
1476 return SCIP_OKAY;
1477 }
1478
1479 /** evaluates a bound for the current LP solution */
1480 static
evalBound(SCIP * scip,BOUND * bound)1481 SCIP_Real evalBound(
1482 SCIP* scip,
1483 BOUND* bound
1484 )
1485 {
1486 assert(scip != NULL);
1487 assert(bound != NULL);
1488
1489 if( bound->boundtype == SCIP_BOUNDTYPE_LOWER )
1490 return REALABS( SCIPvarGetLPSol(bound->var) - SCIPvarGetLbLocal(bound->var) );
1491 else
1492 return REALABS( SCIPvarGetUbLocal(bound->var) - SCIPvarGetLPSol(bound->var) );
1493 }
1494
1495 /** returns the index of the next undone and unfiltered bound with the smallest distance */
1496 static
nextBound(SCIP * scip,SCIP_PROPDATA * propdata,SCIP_Bool convexphase)1497 int nextBound(
1498 SCIP* scip, /**< SCIP data structure */
1499 SCIP_PROPDATA* propdata, /**< data of the obbt propagator */
1500 SCIP_Bool convexphase /**< consider only convex variables? */
1501 )
1502 {
1503 SCIP_Real bestval;
1504 int bestidx;
1505 int k;
1506
1507 assert(scip != NULL);
1508 assert(propdata != NULL);
1509
1510 bestidx = -1;
1511 bestval = SCIPinfinity(scip);
1512
1513 for( k = 0; k <= propdata->lastidx; ++k )
1514 {
1515 BOUND* tmpbound;
1516 tmpbound = propdata->bounds[k];
1517
1518 assert(tmpbound != NULL);
1519
1520 if( !tmpbound->filtered && !tmpbound->done && (tmpbound->nonconvex == !convexphase) )
1521 {
1522 SCIP_Real boundval;
1523
1524 /* return the next bound which is not done or unfiltered yet */
1525 if( propdata->orderingalgo == 0 )
1526 return k;
1527
1528 boundval = evalBound(scip, tmpbound);
1529
1530 /* negate boundval if we use the reverse greedy algorithm */
1531 boundval = (propdata->orderingalgo == 2) ? -1.0 * boundval : boundval;
1532
1533 if( bestidx == -1 || boundval < bestval )
1534 {
1535 bestidx = k;
1536 bestval = boundval;
1537 }
1538 }
1539 }
1540
1541 return bestidx; /*lint !e438*/
1542 }
1543
1544 /** try to separate the solution of the last OBBT LP in order to learn better variable bounds; we apply additional
1545 * separation rounds as long as the routine finds better bounds; because of dual degeneracy we apply a minimum number of
1546 * separation rounds
1547 */
1548 static
applySeparation(SCIP * scip,SCIP_PROPDATA * propdata,BOUND * currbound,SCIP_Longint * nleftiterations,SCIP_Bool * success)1549 SCIP_RETCODE applySeparation(
1550 SCIP* scip, /**< SCIP data structure */
1551 SCIP_PROPDATA* propdata, /**< data of the obbt propagator */
1552 BOUND* currbound, /**< current bound */
1553 SCIP_Longint* nleftiterations, /**< number of left iterations (-1 for no limit) */
1554 SCIP_Bool* success /**< pointer to store if we have found a better bound */
1555 )
1556 {
1557 SCIP_Bool inroot;
1558 int i;
1559
1560 assert(nleftiterations != NULL);
1561 assert(success != NULL);
1562 assert(SCIPinProbing(scip));
1563
1564 *success = FALSE;
1565
1566 /* check if we are originally in the root node */
1567 inroot = SCIPgetDepth(scip) == 1;
1568
1569 for( i = 0; i <= propdata->sepamaxiter; ++i )
1570 {
1571 SCIPdebug( SCIP_Longint nlpiter; )
1572 SCIP_Real oldval;
1573 SCIP_Bool cutoff;
1574 SCIP_Bool delayed;
1575 SCIP_Bool error;
1576 SCIP_Bool optimal;
1577 SCIP_Bool tightened;
1578
1579 oldval = SCIPvarGetLPSol(currbound->var);
1580
1581 /* find and store cuts to separate the current LP solution */
1582 SCIP_CALL( SCIPseparateSol(scip, NULL, inroot, TRUE, FALSE, &delayed, &cutoff) );
1583 SCIPdebugMsg(scip, "applySeparation() - ncuts = %d\n", SCIPgetNCuts(scip));
1584
1585 /* leave if we did not found any cut */
1586 if( SCIPgetNCuts(scip) == 0 )
1587 break;
1588
1589 /* apply cuts and resolve LP */
1590 SCIP_CALL( SCIPapplyCutsProbing(scip, &cutoff) );
1591 assert(SCIPgetNCuts(scip) == 0);
1592 SCIPdebug( nlpiter = SCIPgetNLPIterations(scip); )
1593 SCIP_CALL( solveLP(scip, (int) *nleftiterations, &error, &optimal) );
1594 SCIPdebug( nlpiter = SCIPgetNLPIterations(scip) - nlpiter; )
1595 SCIPdebug( SCIPdebugMsg(scip, "applySeparation() - optimal=%u error=%u lpiter=%" SCIP_LONGINT_FORMAT "\n", optimal, error, nlpiter); )
1596 SCIPdebugMsg(scip, "oldval = %e newval = %e\n", oldval, SCIPvarGetLPSol(currbound->var));
1597
1598 /* leave if we did not solve the LP to optimality or an error occured */
1599 if( error || !optimal )
1600 break;
1601
1602 /* try to generate a genvbound */
1603 if( inroot && propdata->genvboundprop != NULL && propdata->genvbdsduringsepa )
1604 {
1605 SCIP_Bool found;
1606 SCIP_CALL( createGenVBound(scip, propdata, currbound, &found) );
1607 } /*lint !e438*/
1608
1609 /* try to tight the variable bound */
1610 tightened = FALSE;
1611 if( !SCIPisEQ(scip, oldval, SCIPvarGetLPSol(currbound->var)) )
1612 {
1613 SCIP_CALL( tightenBoundProbing(scip, currbound, SCIPvarGetLPSol(currbound->var), &tightened) );
1614 SCIPdebugMsg(scip, "apply separation - tightened=%u oldval=%e newval=%e\n", tightened, oldval,
1615 SCIPvarGetLPSol(currbound->var));
1616
1617 *success |= tightened;
1618 }
1619
1620 /* leave the separation if we did not tighten the bound and proceed at least propdata->sepaminiter iterations */
1621 if( !tightened && i >= propdata->sepaminiter )
1622 break;
1623 }
1624
1625 return SCIP_OKAY;
1626 }
1627
1628 /** finds new variable bounds until no iterations left or all bounds have been checked */
1629 static
findNewBounds(SCIP * scip,SCIP_PROPDATA * propdata,SCIP_Longint * nleftiterations,SCIP_Bool convexphase)1630 SCIP_RETCODE findNewBounds(
1631 SCIP* scip, /**< SCIP data structure */
1632 SCIP_PROPDATA* propdata, /**< data of the obbt propagator */
1633 SCIP_Longint* nleftiterations, /**< pointer to store the number of left iterations */
1634 SCIP_Bool convexphase /**< consider only convex variables? */
1635 )
1636 {
1637 SCIP_Longint nolditerations;
1638 SCIP_Bool iterationsleft;
1639 BOUND* currbound;
1640 SCIP_Longint itlimit;
1641 int nextboundidx;
1642
1643 assert(scip != NULL);
1644 assert(propdata != NULL);
1645 assert(nleftiterations != NULL);
1646
1647 /* update the number of left iterations */
1648 nolditerations = SCIPgetNLPIterations(scip);
1649 itlimit = *nleftiterations;
1650 assert(*nleftiterations == getIterationsLeft(scip, nolditerations, itlimit));
1651 iterationsleft = (*nleftiterations == -1) || (*nleftiterations > 0);
1652
1653 /* To improve the performance we sort the bound in such a way that the undone and
1654 * unfiltered bounds are at the end of propdata->bounds. We calculate and update
1655 * the position of the last unfiltered and undone bound in propdata->lastidx
1656 */
1657 if( !convexphase )
1658 {
1659 /* sort bounds */
1660 SCIP_CALL( sortBounds(scip, propdata) );
1661
1662 /* if the first bound is filtered or done then there is no bound left */
1663 if( propdata->bounds[0]->done || propdata->bounds[0]->filtered )
1664 {
1665 SCIPdebugMsg(scip, "no unprocessed/unfiltered bound left\n");
1666 return SCIP_OKAY;
1667 }
1668
1669 /* compute the last undone and unfiltered node */
1670 propdata->lastidx = 0;
1671 while( propdata->lastidx < propdata->nbounds - 1 && !propdata->bounds[propdata->lastidx]->done &&
1672 !propdata->bounds[propdata->lastidx]->filtered )
1673 ++propdata->lastidx;
1674
1675 SCIPdebugMsg(scip, "lastidx = %d\n", propdata->lastidx);
1676 }
1677
1678 /* find the first unprocessed bound */
1679 nextboundidx = nextBound(scip, propdata, convexphase);
1680
1681 /* skip if there is no bound left */
1682 if( nextboundidx == -1 )
1683 {
1684 SCIPdebugMsg(scip, "no unprocessed/unfiltered bound left\n");
1685 return SCIP_OKAY;
1686 }
1687
1688 currbound = propdata->bounds[nextboundidx];
1689 assert(!currbound->done && !currbound->filtered);
1690
1691 /* main loop */
1692 while( iterationsleft && !SCIPisStopped(scip) )
1693 {
1694 SCIP_Bool optimal;
1695 SCIP_Bool error;
1696 int nfiltered;
1697
1698 assert(currbound != NULL);
1699 assert(currbound->done == FALSE);
1700 assert(currbound->filtered == FALSE);
1701
1702 /* do not visit currbound more than once */
1703 currbound->done = TRUE;
1704 exchangeBounds(propdata, nextboundidx);
1705
1706 /* set objective for curr */
1707 SCIP_CALL( setObjProbing(scip, propdata, currbound, 1.0) );
1708
1709 SCIPdebugMsg(scip, "before solving Boundtype: %d , LB: %e , UB: %e\n",
1710 currbound->boundtype == SCIP_BOUNDTYPE_LOWER, SCIPvarGetLbLocal(currbound->var),
1711 SCIPvarGetUbLocal(currbound->var) );
1712 SCIPdebugMsg(scip, "before solving var <%s>, LP value: %f\n",
1713 SCIPvarGetName(currbound->var), SCIPvarGetLPSol(currbound->var));
1714
1715 SCIPdebugMsg(scip, "probing iterations before solve: %lld \n", SCIPgetNLPIterations(scip));
1716
1717 /* now solve the LP */
1718 SCIP_CALL( solveLP(scip, (int) *nleftiterations, &error, &optimal) );
1719
1720 SCIPdebugMsg(scip, "probing iterations after solve: %lld \n", SCIPgetNLPIterations(scip));
1721 SCIPdebugMsg(scip, "OPT: %u ERROR: %u\n" , optimal, error);
1722 SCIPdebugMsg(scip, "after solving Boundtype: %d , LB: %e , UB: %e\n",
1723 currbound->boundtype == SCIP_BOUNDTYPE_LOWER, SCIPvarGetLbLocal(currbound->var),
1724 SCIPvarGetUbLocal(currbound->var) );
1725 SCIPdebugMsg(scip, "after solving var <%s>, LP value: %f\n",
1726 SCIPvarGetName(currbound->var), SCIPvarGetLPSol(currbound->var));
1727
1728 /* update nleftiterations */
1729 *nleftiterations = getIterationsLeft(scip, nolditerations, itlimit);
1730 iterationsleft = (*nleftiterations == -1) || (*nleftiterations > 0);
1731
1732 if( error )
1733 {
1734 SCIPdebugMsg(scip, "ERROR during LP solving\n");
1735
1736 /* set the objective of currbound to zero to null the whole objective; otherwise the objective is wrong when
1737 * we call findNewBounds() for the convex phase
1738 */
1739 SCIP_CALL( SCIPchgVarObjProbing(scip, currbound->var, 0.0) );
1740
1741 return SCIP_OKAY;
1742 }
1743
1744 if( optimal )
1745 {
1746 SCIP_Bool success;
1747
1748 currbound->newval = SCIPvarGetLPSol(currbound->var);
1749 currbound->found = TRUE;
1750
1751 /* in root node we may want to create a genvbound (independent of tightening success) */
1752 if( (SCIPgetDepth(scip) == 0 || (SCIPinProbing(scip) && SCIPgetDepth(scip) == 1))
1753 && propdata->genvboundprop != NULL )
1754 {
1755 SCIP_Bool found;
1756
1757 SCIP_CALL( createGenVBound(scip, propdata, currbound, &found) );
1758 } /*lint !e438*/
1759
1760 /* try to tighten bound in probing mode */
1761 success = FALSE;
1762 if( propdata->tightintboundsprobing && SCIPvarIsIntegral(currbound->var) )
1763 {
1764 SCIPdebugMsg(scip, "tightening bound %s = %e bounds: [%e, %e]\n", SCIPvarGetName(currbound->var),
1765 currbound->newval, SCIPvarGetLbLocal(currbound->var), SCIPvarGetUbLocal(currbound->var) );
1766 SCIP_CALL( tightenBoundProbing(scip, currbound, currbound->newval, &success) );
1767 SCIPdebugMsg(scip, "tightening bound %s\n", success ? "successful" : "not successful");
1768 }
1769 else if( propdata->tightcontboundsprobing && !SCIPvarIsIntegral(currbound->var) )
1770 {
1771 SCIPdebugMsg(scip, "tightening bound %s = %e bounds: [%e, %e]\n", SCIPvarGetName(currbound->var),
1772 currbound->newval, SCIPvarGetLbLocal(currbound->var), SCIPvarGetUbLocal(currbound->var) );
1773 SCIP_CALL( tightenBoundProbing(scip, currbound, currbound->newval, &success) );
1774 SCIPdebugMsg(scip, "tightening bound %s\n", success ? "successful" : "not successful");
1775 }
1776
1777 /* separate current OBBT LP solution */
1778 if( iterationsleft && propdata->separatesol )
1779 {
1780 SCIP_CALL( applySeparation(scip, propdata, currbound, nleftiterations, &success) );
1781
1782 /* remember best solution value after solving additional separations LPs */
1783 if( success )
1784 {
1785 #ifndef NDEBUG
1786 SCIP_Real newval = SCIPvarGetLPSol(currbound->var);
1787
1788 /* round new bound if the variable is integral */
1789 if( SCIPvarIsIntegral(currbound->var) )
1790 newval = currbound->boundtype == SCIP_BOUNDTYPE_LOWER ?
1791 SCIPceil(scip, newval) : SCIPfloor(scip, newval);
1792
1793 assert((currbound->boundtype == SCIP_BOUNDTYPE_LOWER &&
1794 SCIPisGT(scip, newval, currbound->newval))
1795 || (currbound->boundtype == SCIP_BOUNDTYPE_UPPER &&
1796 SCIPisLT(scip, newval, currbound->newval)));
1797 #endif
1798
1799 currbound->newval = SCIPvarGetLPSol(currbound->var);
1800 }
1801 }
1802
1803 /* filter bound candidates by using the current LP solution */
1804 if( propdata->applytrivialfilter )
1805 {
1806 SCIP_CALL( filterExistingLP(scip, propdata, &nfiltered, currbound) );
1807 SCIPdebugMsg(scip, "filtered %d bounds via inspecting present LP solution\n", nfiltered);
1808 }
1809
1810 propdata->propagatecounter += success ? 1 : 0;
1811
1812 /* propagate if we have found enough bound tightenings */
1813 if( propdata->propagatefreq != 0 && propdata->propagatecounter >= propdata->propagatefreq )
1814 {
1815 SCIP_Longint ndomredsfound;
1816 SCIP_Bool cutoff;
1817
1818 SCIP_CALL( SCIPpropagateProbing(scip, 0, &cutoff, &ndomredsfound) );
1819 SCIPdebugMsg(scip, "propagation - cutoff %u ndomreds %" SCIP_LONGINT_FORMAT "\n", cutoff, ndomredsfound);
1820
1821 propdata->npropagatedomreds += ndomredsfound;
1822 propdata->propagatecounter = 0;
1823 }
1824 }
1825
1826 /* set objective to zero */
1827 SCIP_CALL( setObjProbing(scip, propdata, currbound, 0.0) );
1828
1829 /* find the first unprocessed bound */
1830 nextboundidx = nextBound(scip, propdata, convexphase);
1831
1832 /* check if there is no unprocessed and unfiltered node left */
1833 if( nextboundidx == -1 )
1834 {
1835 SCIPdebugMsg(scip, "NO unvisited/unfiltered bound left!\n");
1836 break;
1837 }
1838
1839 currbound = propdata->bounds[nextboundidx];
1840 assert(!currbound->done && !currbound->filtered);
1841 }
1842
1843 if( iterationsleft )
1844 {
1845 SCIPdebugMsg(scip, "still iterations left: %" SCIP_LONGINT_FORMAT "\n", *nleftiterations);
1846 }
1847 else
1848 {
1849 SCIPdebugMsg(scip, "no iterations left\n");
1850 }
1851
1852 return SCIP_OKAY;
1853 }
1854
1855
1856 /** main function of obbt */
1857 static
applyObbt(SCIP * scip,SCIP_PROPDATA * propdata,SCIP_Longint itlimit,SCIP_RESULT * result)1858 SCIP_RETCODE applyObbt(
1859 SCIP* scip, /**< SCIP data structure */
1860 SCIP_PROPDATA* propdata, /**< data of the obbt propagator */
1861 SCIP_Longint itlimit, /**< LP iteration limit (-1: no limit) */
1862 SCIP_RESULT* result /**< result pointer */
1863 )
1864 {
1865 SCIP_VAR** vars;
1866 SCIP_Real* oldlbs;
1867 SCIP_Real* oldubs;
1868 SCIP_Longint lastnpropagatedomreds;
1869 SCIP_Longint nleftiterations;
1870 SCIP_Real oldconditionlimit;
1871 SCIP_Real oldboundstreps;
1872 SCIP_Real olddualfeastol;
1873 SCIP_Bool hasconditionlimit;
1874 SCIP_Bool continuenode;
1875 SCIP_Bool boundleft;
1876 int oldpolishing;
1877 int nfiltered;
1878 int nvars;
1879 int i;
1880
1881 assert(scip != NULL);
1882 assert(propdata != NULL);
1883 assert(itlimit == -1 || itlimit >= 0);
1884
1885 SCIPdebugMsg(scip, "apply obbt\n");
1886
1887 oldlbs = NULL;
1888 oldubs = NULL;
1889 lastnpropagatedomreds = propdata->npropagatedomreds;
1890 nleftiterations = itlimit;
1891 continuenode = SCIPnodeGetNumber(SCIPgetCurrentNode(scip)) == propdata->lastnode;
1892 propdata->lastidx = -1;
1893 boundleft = FALSE;
1894 *result = SCIP_DIDNOTFIND;
1895
1896 /* store old variable bounds if we use propagation during obbt */
1897 if( propdata->propagatefreq > 0 )
1898 {
1899 SCIP_CALL( SCIPallocBufferArray(scip, &oldlbs, propdata->nbounds) );
1900 SCIP_CALL( SCIPallocBufferArray(scip, &oldubs, propdata->nbounds) );
1901 }
1902
1903 /* reset bound data structure flags; fixed variables are marked as filtered */
1904 for( i = 0; i < propdata->nbounds; i++ )
1905 {
1906 BOUND* bound = propdata->bounds[i];
1907 bound->found = FALSE;
1908
1909 /* store old variable bounds */
1910 if( oldlbs != NULL && oldubs != NULL )
1911 {
1912 oldlbs[bound->index] = SCIPvarGetLbLocal(bound->var);
1913 oldubs[bound->index] = SCIPvarGetUbLocal(bound->var);
1914 }
1915
1916 /* reset 'done' and 'filtered' flag in a new B&B node */
1917 if( !continuenode )
1918 {
1919 bound->done = FALSE;
1920 bound->filtered = FALSE;
1921 }
1922
1923 /* mark fixed variables as filtered */
1924 bound->filtered |= varIsFixedLocal(scip, bound->var);
1925
1926 /* check for an unprocessed bound */
1927 if( !bound->filtered && !bound->done )
1928 boundleft = TRUE;
1929 }
1930
1931 /* no bound left to check */
1932 if( !boundleft )
1933 goto TERMINATE;
1934
1935 /* filter variables via inspecting present LP solution */
1936 if( propdata->applytrivialfilter && !continuenode )
1937 {
1938 SCIP_CALL( filterExistingLP(scip, propdata, &nfiltered, NULL) );
1939 SCIPdebugMsg(scip, "filtered %d bounds via inspecting present LP solution\n", nfiltered);
1940 }
1941
1942 /* store old dualfeasibletol */
1943 olddualfeastol = SCIPdualfeastol(scip);
1944
1945 /* start probing */
1946 SCIP_CALL( SCIPstartProbing(scip) );
1947 SCIPdebugMsg(scip, "start probing\n");
1948
1949 /* tighten dual feastol */
1950 if( propdata->dualfeastol < olddualfeastol )
1951 {
1952 SCIP_CALL( SCIPchgDualfeastol(scip, propdata->dualfeastol) );
1953 }
1954
1955 /* tighten condition limit */
1956 hasconditionlimit = (SCIPgetRealParam(scip, "lp/conditionlimit", &oldconditionlimit) == SCIP_OKAY);
1957 if( !hasconditionlimit )
1958 {
1959 SCIPwarningMessage(scip, "obbt propagator could not set condition limit in LP solver - running without\n");
1960 }
1961 else if( propdata->conditionlimit > 0.0 && (oldconditionlimit < 0.0 || propdata->conditionlimit < oldconditionlimit) )
1962 {
1963 SCIP_CALL( SCIPsetRealParam(scip, "lp/conditionlimit", propdata->conditionlimit) );
1964 }
1965
1966 /* tighten relative bound improvement limit */
1967 SCIP_CALL( SCIPgetRealParam(scip, "numerics/boundstreps", &oldboundstreps) );
1968 if( !SCIPisEQ(scip, oldboundstreps, propdata->boundstreps) )
1969 {
1970 SCIP_CALL( SCIPsetRealParam(scip, "numerics/boundstreps", propdata->boundstreps) );
1971 }
1972
1973 /* add objective cutoff */
1974 SCIP_CALL( addObjCutoff(scip, propdata) );
1975
1976 /* deactivate LP polishing */
1977 SCIP_CALL( SCIPgetIntParam(scip, "lp/solutionpolishing", &oldpolishing) );
1978 SCIP_CALL( SCIPsetIntParam(scip, "lp/solutionpolishing", 0) );
1979
1980 /* apply filtering */
1981 if( propdata->applyfilterrounds )
1982 {
1983 SCIP_CALL( filterBounds(scip, propdata, nleftiterations) );
1984 }
1985
1986 /* set objective coefficients to zero */
1987 vars = SCIPgetVars(scip);
1988 nvars = SCIPgetNVars(scip);
1989 for( i = 0; i < nvars; ++i )
1990 {
1991 /* note that it is not possible to change the objective of non-column variables during probing; we have to take
1992 * care of the objective contribution of loose variables in createGenVBound()
1993 */
1994 if( SCIPvarGetObj(vars[i]) != 0.0 && SCIPvarGetStatus(vars[i]) == SCIP_VARSTATUS_COLUMN )
1995 {
1996 SCIP_CALL( SCIPchgVarObjProbing(scip, vars[i], 0.0) );
1997 }
1998 }
1999
2000 /* find new bounds for the variables */
2001 SCIP_CALL( findNewBounds(scip, propdata, &nleftiterations, FALSE) );
2002
2003 if( nleftiterations > 0 || itlimit < 0 )
2004 {
2005 SCIP_CALL( findNewBounds(scip, propdata, &nleftiterations, TRUE) );
2006 }
2007
2008 /* reset dual feastol and condition limit */
2009 SCIP_CALL( SCIPchgDualfeastol(scip, olddualfeastol) );
2010 if( hasconditionlimit )
2011 {
2012 SCIP_CALL( SCIPsetRealParam(scip, "lp/conditionlimit", oldconditionlimit) );
2013 }
2014
2015 /* update bound->newval if we have learned additional bound tightenings during SCIPpropagateProbing() */
2016 if( oldlbs != NULL && oldubs != NULL && propdata->npropagatedomreds - lastnpropagatedomreds > 0 )
2017 {
2018 assert(propdata->propagatefreq > 0);
2019 for( i = 0; i < propdata->nbounds; ++i )
2020 {
2021 BOUND* bound = propdata->bounds[i];
2022
2023 /* it might be the case that a bound found by the additional propagation is better than the bound found after solving an OBBT
2024 * LP
2025 */
2026 if( bound->found )
2027 {
2028 if( bound->boundtype == SCIP_BOUNDTYPE_LOWER )
2029 bound->newval = MAX(bound->newval, SCIPvarGetLbLocal(bound->var)); /*lint !e666*/
2030 else
2031 bound->newval = MIN(bound->newval, SCIPvarGetUbLocal(bound->var)); /*lint !e666*/
2032 }
2033 else
2034 {
2035 SCIP_Real oldlb;
2036 SCIP_Real oldub;
2037
2038 oldlb = oldlbs[bound->index];
2039 oldub = oldubs[bound->index];
2040
2041 if( bound->boundtype == SCIP_BOUNDTYPE_LOWER && SCIPisLbBetter(scip, SCIPvarGetLbLocal(bound->var), oldlb, oldub) )
2042 {
2043 SCIPdebugMsg(scip, "tighter lower bound due to propagation: %d - %e -> %e\n", i, oldlb, SCIPvarGetLbLocal(bound->var));
2044 bound->newval = SCIPvarGetLbLocal(bound->var);
2045 bound->found = TRUE;
2046 }
2047
2048 if( bound->boundtype == SCIP_BOUNDTYPE_UPPER && SCIPisUbBetter(scip, SCIPvarGetUbLocal(bound->var), oldlb, oldub) )
2049 {
2050 SCIPdebugMsg(scip, "tighter upper bound due to propagation: %d - %e -> %e\n", i, oldub, SCIPvarGetUbLocal(bound->var));
2051 bound->newval = SCIPvarGetUbLocal(bound->var);
2052 bound->found = TRUE;
2053 }
2054 }
2055 }
2056 }
2057
2058 /* reset relative bound improvement limit */
2059 SCIP_CALL( SCIPsetRealParam(scip, "numerics/boundstreps", oldboundstreps) );
2060
2061 /* reset original LP polishing setting */
2062 SCIP_CALL( SCIPsetIntParam(scip, "lp/solutionpolishing", oldpolishing) );
2063
2064 /* end probing */
2065 SCIP_CALL( SCIPendProbing(scip) );
2066 SCIPdebugMsg(scip, "end probing!\n");
2067
2068 /* release cutoff row if there is one */
2069 if( propdata->cutoffrow != NULL )
2070 {
2071 assert(!SCIProwIsInLP(propdata->cutoffrow));
2072 SCIP_CALL( SCIPreleaseRow(scip, &(propdata->cutoffrow)) );
2073 }
2074
2075 /* apply buffered bound changes */
2076 SCIP_CALL( applyBoundChgs(scip, propdata, result) );
2077
2078 TERMINATE:
2079 SCIPfreeBufferArrayNull(scip, &oldubs);
2080 SCIPfreeBufferArrayNull(scip, &oldlbs);
2081
2082 return SCIP_OKAY;
2083 }
2084
2085 /** computes a valid inequality from the current LP relaxation for a bilinear term xy only involving x and y; the
2086 * inequality is found by optimizing along the line connecting the points (xs,ys) and (xt,yt) over the currently given
2087 * linear relaxation of the problem; this optimization problem is an LP
2088 *
2089 * max lambda
2090 * s.t. Ax <= b
2091 * (x,y) = (xs,ys) + lambda ((xt,yt) - (xs,ys))
2092 * lambda in [0,1]
2093 *
2094 * which is equivalent to
2095 *
2096 * max x
2097 * s.t. (1) Ax <= b
2098 * (2) (x - xs) / (xt - xs) = (y - ys) / (yt - ys)
2099 *
2100 * Let x* be the optimal primal and (mu,theta) be the optimal dual solution of this LP. The KKT conditions imply that
2101 * the aggregation of the linear constraints mu*Ax <= mu*b can be written as
2102 *
2103 * x * (1 - theta / (xt - xs)) + y * theta / (yt - ys) = mu * Ax <= mu * b
2104 *
2105 * <=> alpha * x + beta * y <= mu * b = alpha * (x*) + beta * (y*)
2106 *
2107 * which is a valid inequality in the (x,y)-space; in order to avoid numerical difficulties when (xs,ys) is too close
2108 * to (xt,yt), we scale constraint (2) by min{ max{1,|xt-xs|,|yt-ys|}, 100 } beforehand
2109 */
2110 static
solveBilinearLP(SCIP * scip,SCIP_VAR * x,SCIP_VAR * y,SCIP_Real xs,SCIP_Real ys,SCIP_Real xt,SCIP_Real yt,SCIP_Real * xcoef,SCIP_Real * ycoef,SCIP_Real * constant,SCIP_Longint iterlim)2111 SCIP_RETCODE solveBilinearLP(
2112 SCIP* scip, /**< SCIP data structure */
2113 SCIP_VAR* x, /**< first variable */
2114 SCIP_VAR* y, /**< second variable */
2115 SCIP_Real xs, /**< x-coordinate of the first point */
2116 SCIP_Real ys, /**< y-coordinate of the first point */
2117 SCIP_Real xt, /**< x-coordinate of the second point */
2118 SCIP_Real yt, /**< y-coordinate of the second point */
2119 SCIP_Real* xcoef, /**< pointer to store the coefficient of x */
2120 SCIP_Real* ycoef, /**< pointer to store the coefficient of y */
2121 SCIP_Real* constant, /**< pointer to store the constant */
2122 SCIP_Longint iterlim /**< iteration limit (-1: for no limit) */
2123 )
2124 {
2125 SCIP_ROW* row;
2126 SCIP_Real signx;
2127 SCIP_Real scale;
2128 SCIP_Real side;
2129 SCIP_Bool lperror;
2130
2131 assert(xcoef != NULL);
2132 assert(ycoef != NULL);
2133 assert(constant != NULL);
2134 assert(SCIPinProbing(scip));
2135
2136 *xcoef = SCIP_INVALID;
2137 *ycoef = SCIP_INVALID;
2138 *constant= SCIP_INVALID;
2139
2140 SCIPdebugMsg(scip, " solve bilinear LP for (%s,%s) from (%g,%g) to (%g,%g)\n", SCIPvarGetName(x), SCIPvarGetName(y), xs,
2141 ys, xt, yt);
2142
2143 /* skip computations if (xs,ys) and (xt,yt) are too close to each other or contain too large values */
2144 if( SCIPisFeasEQ(scip, xs, xt) || SCIPisFeasEQ(scip, ys, yt)
2145 || SCIPisHugeValue(scip, REALABS(xs)) || SCIPisHugeValue(scip, REALABS(xt))
2146 || SCIPisHugeValue(scip, REALABS(ys)) || SCIPisHugeValue(scip, REALABS(yt)) )
2147 {
2148 SCIPdebugMsg(scip, " -> skip: bounds are too close/large\n");
2149 return SCIP_OKAY;
2150 }
2151
2152 /* compute scaler for the additional linear constraint */
2153 scale = MIN(MAX3(1.0, REALABS(xt-xs), REALABS(yt-ys)), 100.0); /*lint !e666*/
2154
2155 /* set objective function */
2156 signx = (xs > xt) ? 1.0 : -1.0;
2157 SCIP_CALL( SCIPchgVarObjProbing(scip, x, signx) );
2158
2159 /* create new probing node to remove the added LP row afterwards */
2160 SCIP_CALL( SCIPnewProbingNode(scip) );
2161
2162 /* create row */
2163 side = scale * (xs/(xt-xs) - ys/(yt-ys));
2164 SCIP_CALL( SCIPcreateEmptyRowUnspec(scip, &row, "bilinrow", side, side, FALSE, FALSE, TRUE) );
2165 SCIP_CALL( SCIPaddVarToRow(scip, row, x, scale/(xt-xs)) );
2166 SCIP_CALL( SCIPaddVarToRow(scip, row, y, -scale/(yt-ys)) );
2167 SCIP_CALL( SCIPaddRowProbing(scip, row) );
2168
2169 /* solve probing LP */
2170 #ifdef NDEBUG
2171 {
2172 SCIP_RETCODE retstat;
2173 retstat = SCIPsolveProbingLP(scip, iterlim, &lperror, NULL);
2174 if( retstat != SCIP_OKAY )
2175 {
2176 SCIPwarningMessage(scip, "Error while solving LP in quadratic constraint handler; LP solve terminated with" \
2177 "code <%d>\n", retstat);
2178 }
2179 }
2180 #else
2181 SCIP_CALL( SCIPsolveProbingLP(scip, (int)iterlim, &lperror, NULL) ); /*lint !e712*/
2182 #endif
2183
2184 SCIPdebugMsg(scip, " solved probing LP -> lperror=%u lpstat=%d\n", lperror, SCIPgetLPSolstat(scip));
2185
2186 /* collect dual and primal solution entries */
2187 if( !lperror && SCIPgetLPSolstat(scip) == SCIP_LPSOLSTAT_OPTIMAL )
2188 {
2189 SCIP_Real xval = SCIPvarGetLPSol(x);
2190 SCIP_Real yval = SCIPvarGetLPSol(y);
2191 SCIP_Real mu = -SCIProwGetDualsol(row);
2192
2193 SCIPdebugMsg(scip, " primal=(%g,%g) dual=%g\n", xval, yval, mu);
2194
2195 /* xcoef x + ycoef y <= constant */
2196 *xcoef = -signx - (mu * scale) / (xt - xs);
2197 *ycoef = (mu * scale) / (yt - ys);
2198 *constant = (*xcoef) * xval + (*ycoef) * yval;
2199
2200 /* xcoef x <= -ycoef y + constant */
2201 *ycoef = -(*ycoef);
2202
2203 /* inequality is only useful when both coefficients are different from zero; normalize inequality if possible */
2204 if( !SCIPisFeasZero(scip, *xcoef) && !SCIPisFeasZero(scip, *ycoef) )
2205 {
2206 SCIP_Real val = REALABS(*xcoef);
2207 *xcoef /= val;
2208 *ycoef /= val;
2209 *constant /= val;
2210 }
2211 else
2212 {
2213 *xcoef = SCIP_INVALID;
2214 *ycoef = SCIP_INVALID;
2215 *constant = SCIP_INVALID;
2216 }
2217 }
2218
2219 /* release row and backtrack probing node */
2220 SCIP_CALL( SCIPreleaseRow(scip, &row) );
2221 SCIP_CALL( SCIPbacktrackProbing(scip, 0) );
2222
2223 /* reset objective function */
2224 SCIP_CALL( SCIPchgVarObjProbing(scip, x, 0.0) );
2225
2226 return SCIP_OKAY;
2227 }
2228
2229 /* applies obbt for finding valid inequalities for bilinear terms; function works as follows:
2230 *
2231 * 1. start probing mode
2232 * 2. add objective cutoff (if possible)
2233 * 3. set objective function to zero
2234 * 4. set feasibility, optimality, and relative bound improvement tolerances of SCIP
2235 * 5. main loop
2236 * 6. restore old tolerances
2237 *
2238 */
2239 static
applyObbtBilinear(SCIP * scip,SCIP_PROPDATA * propdata,SCIP_Longint itlimit,SCIP_RESULT * result)2240 SCIP_RETCODE applyObbtBilinear(
2241 SCIP* scip, /**< SCIP data structure */
2242 SCIP_PROPDATA* propdata, /**< data of the obbt propagator */
2243 SCIP_Longint itlimit, /**< LP iteration limit (-1: no limit) */
2244 SCIP_RESULT* result /**< result pointer */
2245 )
2246 {
2247 SCIP_VAR** vars;
2248 SCIP_Real oldfeastol;
2249 SCIP_Bool lperror;
2250 SCIP_Longint nolditerations;
2251 SCIP_Longint nleftiterations;
2252 int nvars;
2253 int i;
2254
2255 assert(scip != NULL);
2256 assert(propdata != NULL);
2257 assert(itlimit == -1 || itlimit >= 0);
2258 assert(result != NULL);
2259
2260 if( propdata->nbilinbounds <= 0 || SCIPgetDepth(scip) != 0 || propdata->lastbilinidx >= propdata->nbilinbounds )
2261 return SCIP_OKAY;
2262
2263 SCIPdebugMsg(scip, "call applyObbtBilinear starting from %d\n", propdata->lastbilinidx);
2264
2265 vars = SCIPgetVars(scip);
2266 nvars = SCIPgetNVars(scip);
2267
2268 nolditerations = SCIPgetNLPIterations(scip);
2269 nleftiterations = getIterationsLeft(scip, nolditerations, itlimit);
2270 SCIPdebugMsg(scip, "iteration limit: %lld\n", nleftiterations);
2271
2272 /* 1. start probing */
2273 SCIP_CALL( SCIPstartProbing(scip) );
2274
2275 /* 2. add objective cutoff */
2276 SCIP_CALL( addObjCutoff(scip, propdata) );
2277
2278 /* 3. set objective function to zero */
2279 for( i = 0; i < nvars; ++i )
2280 {
2281 SCIP_CALL( SCIPchgVarObjProbing(scip, vars[i], 0.0) );
2282 }
2283
2284 /* 4. tighten LP feasibility tolerance to be at most feastol/10.0 */
2285 oldfeastol = SCIPchgRelaxfeastol(scip, SCIPfeastol(scip) / 10.0);
2286
2287 /* we need to solve the probing LP before creating new probing nodes in solveBilinearLP() */
2288 SCIP_CALL( SCIPsolveProbingLP(scip, (int)nleftiterations, &lperror, NULL) );
2289
2290 if( lperror )
2291 goto TERMINATE;
2292
2293 /* 5. main loop */
2294 for( i = propdata->lastbilinidx; i < propdata->nbilinbounds
2295 && (nleftiterations > 0 || nleftiterations == -1)
2296 && (propdata->itlimitbilin < 0 || propdata->itlimitbilin > propdata->itusedbilin )
2297 && !SCIPisStopped(scip); ++i )
2298 {
2299 CORNER corners[4] = {LEFTBOTTOM, LEFTTOP, RIGHTTOP, RIGHTBOTTOM};
2300 BILINBOUND* bilinbound;
2301 int k;
2302
2303 bilinbound = propdata->bilinbounds[i];
2304 assert(bilinbound != NULL);
2305
2306 SCIPdebugMsg(scip, "process %d: %s %s done=%u filtered=%d nunderest=%d noverest=%d\n", i,
2307 SCIPvarGetName(bilinbound->x), SCIPvarGetName(bilinbound->y), bilinbound->done, bilinbound->filtered,
2308 bilinbound->nunderest, bilinbound->noverest);
2309
2310 /* we already solved LPs for this bilinear term */
2311 if( bilinbound->done || bilinbound->filtered == (int)FILTERED )
2312 continue;
2313
2314 /* iterate through all corners
2315 *
2316 * 0: (xs,ys)=(ubx,lby) (xt,yt)=(lbx,uby) -> underestimate
2317 * 1: (xs,ys)=(ubx,uby) (xt,yt)=(lbx,lby) -> overestimate
2318 * 2: (xs,ys)=(lbx,uby) (xt,yt)=(ubx,lby) -> underestimate
2319 * 3: (xs,ys)=(lbx,lby) (xt,yt)=(ubx,uby) -> overestimate
2320 */
2321 for( k = 0; k < 4; ++k )
2322 {
2323 CORNER corner = corners[k];
2324 SCIP_Real xcoef;
2325 SCIP_Real ycoef;
2326 SCIP_Real constant;
2327 SCIP_Real xs = SCIP_INVALID;
2328 SCIP_Real ys = SCIP_INVALID;
2329 SCIP_Real xt = SCIP_INVALID;
2330 SCIP_Real yt = SCIP_INVALID;
2331
2332 /* skip corners that lead to an under- or overestimate that is not needed */
2333 if( ((corner == LEFTTOP || corner == RIGHTBOTTOM) && bilinbound->nunderest == 0)
2334 || ((corner == LEFTBOTTOM || corner == RIGHTTOP) && bilinbound->noverest == 0) )
2335 continue;
2336
2337 /* check whether corner has been filtered already */
2338 if( (bilinbound->filtered & corner) != 0 ) /*lint !e641*/
2339 continue;
2340
2341 /* get corners (xs,ys) and (xt,yt) */
2342 getCorners(bilinbound->x, bilinbound->y, corner, &xs, &ys, &xt, &yt);
2343
2344 /* skip target corner points with too large values */
2345 if( SCIPisHugeValue(scip, REALABS(xt)) || SCIPisHugeValue(scip, REALABS(yt)) )
2346 continue;
2347
2348 /* compute inequality */
2349 propdata->itusedbilin -= SCIPgetNLPIterations(scip);
2350 SCIP_CALL( solveBilinearLP(scip, bilinbound->x, bilinbound->y, xs, ys, xt, yt, &xcoef, &ycoef, &constant, -1L) );
2351 propdata->itusedbilin += SCIPgetNLPIterations(scip);
2352
2353 /* update number of LP iterations */
2354 nleftiterations = getIterationsLeft(scip, nolditerations, itlimit);
2355 SCIPdebugMsg(scip, "LP iterations left: %lld\n", nleftiterations);
2356
2357 /* add inequality to quadratic constraint handler if it separates (xt,yt) */
2358 if( !SCIPisHugeValue(scip, xcoef) && !SCIPisFeasZero(scip, xcoef)
2359 && REALABS(ycoef) < 1e+3 && REALABS(ycoef) > 1e-3 /* TODO add a parameter for this */
2360 && SCIPisFeasGT(scip, (xcoef*xt - ycoef*yt - constant) / SQRT(SQR(xcoef) + SQR(ycoef) + SQR(constant)), 1e-2) )
2361 {
2362 SCIP_Bool success;
2363
2364 SCIP_CALL( SCIPaddBilinearIneqQuadratic(scip, bilinbound->x, bilinbound->y, bilinbound->index, xcoef,
2365 ycoef, constant, &success) );
2366
2367 /* check whether the inequality has been accepted by the quadratic constraint handler */
2368 if( success )
2369 {
2370 *result = SCIP_REDUCEDDOM;
2371 SCIPdebugMsg(scip, " found %g x <= %g y + %g with violation %g\n", xcoef, ycoef, constant,
2372 (xcoef*xt - ycoef*yt - constant) / SQRT(SQR(xcoef) + SQR(ycoef) + SQR(constant)));
2373 }
2374 }
2375 }
2376
2377 /* mark the bound as processed */
2378 bilinbound->done = TRUE;
2379 }
2380
2381 /* remember last unprocessed bilinear term */
2382 propdata->lastbilinidx = i;
2383
2384 TERMINATE:
2385 /* end probing */
2386 SCIP_CALL( SCIPendProbing(scip) );
2387
2388 /* release cutoff row if there is one */
2389 if( propdata->cutoffrow != NULL )
2390 {
2391 assert(!SCIProwIsInLP(propdata->cutoffrow));
2392 SCIP_CALL( SCIPreleaseRow(scip, &(propdata->cutoffrow)) );
2393 }
2394
2395 /* 6. restore old tolerance */
2396 (void) SCIPchgRelaxfeastol(scip, oldfeastol);
2397
2398 return SCIP_OKAY;
2399 }
2400
2401 /** computes the score of a bound */
2402 static
getScore(SCIP * scip,BOUND * bound,int nlcount,int maxnlcount)2403 unsigned int getScore(
2404 SCIP* scip, /**< SCIP data structure */
2405 BOUND* bound, /**< pointer of bound */
2406 int nlcount, /**< number of nonlinear constraints containing the bounds variable */
2407 int maxnlcount /**< maximal number of nonlinear constraints a variable appears in */
2408 )
2409 {
2410 unsigned int score; /* score to be computed */
2411
2412 assert(scip != NULL);
2413 assert(bound != NULL);
2414 assert(nlcount >= 0);
2415 assert(maxnlcount >= nlcount);
2416
2417 /* score = ( nlcount * ( BASE - 1 ) / maxnlcount ) * BASE^2 + vartype * BASE + boundtype */
2418 score = (unsigned int) ( nlcount > 0 ? (OBBT_SCOREBASE * nlcount * ( OBBT_SCOREBASE - 1 )) / maxnlcount : 0 ); /*lint !e414*/
2419 switch( SCIPvarGetType(bound->var) )
2420 {
2421 case SCIP_VARTYPE_INTEGER:
2422 score += 1;
2423 break;
2424 case SCIP_VARTYPE_IMPLINT:
2425 score += 2;
2426 break;
2427 case SCIP_VARTYPE_CONTINUOUS:
2428 score += 3;
2429 break;
2430 case SCIP_VARTYPE_BINARY:
2431 score += 4;
2432 break;
2433 default:
2434 break;
2435 }
2436
2437 score *= OBBT_SCOREBASE;
2438 if( bound->boundtype == SCIP_BOUNDTYPE_UPPER )
2439 score += 1;
2440
2441 return score;
2442 }
2443
2444 /** computes the score of a bilinear term bound */
2445 static
getScoreBilinBound(SCIP * scip,SCIP_RANDNUMGEN * randnumgen,BILINBOUND * bilinbound,int nbilinterms)2446 SCIP_Real getScoreBilinBound(
2447 SCIP* scip, /**< SCIP data structure */
2448 SCIP_RANDNUMGEN* randnumgen, /**< random number generator */
2449 BILINBOUND* bilinbound, /**< bilinear term bound */
2450 int nbilinterms /**< maximal number of bilinear terms in all quadratic constraints */
2451 )
2452 {
2453 SCIP_Real lbx = SCIPvarGetLbLocal(bilinbound->x);
2454 SCIP_Real ubx = SCIPvarGetUbLocal(bilinbound->x);
2455 SCIP_Real lby = SCIPvarGetLbLocal(bilinbound->y);
2456 SCIP_Real uby = SCIPvarGetUbLocal(bilinbound->y);
2457 SCIP_Real score;
2458
2459 assert(scip != NULL);
2460 assert(randnumgen != NULL);
2461 assert(bilinbound != NULL);
2462
2463 /* consider how often a bilinear term is present in the problem */
2464 score = (bilinbound->noverest + bilinbound->nunderest) / (SCIP_Real)nbilinterms;
2465
2466 /* penalize small variable domains; TODO tune the factor in the logarithm, maybe add a parameter for it */
2467 if( ubx - lbx < 0.5 )
2468 score += log(2.0*(ubx-lbx) + SCIPepsilon(scip));
2469 if( uby - lby < 0.5 )
2470 score += log(2.0*(uby-lby) + SCIPepsilon(scip));
2471
2472 /* consider interiority of variables in the LP solution */
2473 if( SCIPgetLPSolstat(scip) == SCIP_LPSOLSTAT_OPTIMAL )
2474 {
2475 SCIP_Real solx = SCIPvarGetLPSol(bilinbound->x);
2476 SCIP_Real soly = SCIPvarGetLPSol(bilinbound->y);
2477 SCIP_Real interiorityx = MIN(solx-lbx, ubx-solx) / MAX(ubx-lbx, SCIPepsilon(scip)); /*lint !e666*/
2478 SCIP_Real interiorityy = MIN(soly-lby, uby-soly) / MAX(uby-lby, SCIPepsilon(scip)); /*lint !e666*/
2479
2480 score += interiorityx + interiorityy;
2481 }
2482
2483 /* randomize score */
2484 score *= 1.0 + SCIPrandomGetReal(randnumgen, -SCIPepsilon(scip), SCIPepsilon(scip));
2485
2486 return score;
2487 }
2488
2489 /** count the variables which appear in non-convex term of nlrow */
2490 static
countNLRowVarsNonConvexity(SCIP * scip,int * nlcounts,SCIP_NLROW * nlrow)2491 SCIP_RETCODE countNLRowVarsNonConvexity(
2492 SCIP* scip, /**< SCIP data structure */
2493 int* nlcounts, /**< store the number each variable appears in a
2494 * non-convex term */
2495 SCIP_NLROW* nlrow /**< nonlinear row */
2496 )
2497 {
2498 int t;
2499 int nexprtreevars;
2500 SCIP_VAR** exprtreevars;
2501 SCIP_EXPRTREE* exprtree;
2502
2503 assert(scip != NULL);
2504 assert(nlcounts != NULL);
2505 assert(nlrow != NULL);
2506
2507 /* go through all quadratic terms */
2508 for( t = SCIPnlrowGetNQuadElems(nlrow) - 1; t >= 0; --t )
2509 {
2510 SCIP_QUADELEM* quadelem;
2511 SCIP_VAR* bilinvar1;
2512 SCIP_VAR* bilinvar2;
2513
2514 /* get quadratic term */
2515 quadelem = &SCIPnlrowGetQuadElems(nlrow)[t];
2516
2517 /* get involved variables */
2518 bilinvar1 = SCIPnlrowGetQuadVars(nlrow)[quadelem->idx1];
2519 bilinvar2 = SCIPnlrowGetQuadVars(nlrow)[quadelem->idx2];
2520
2521 assert(bilinvar1 != NULL);
2522 assert(bilinvar2 != NULL);
2523
2524 /* we have a non-convex square term */
2525 if( bilinvar1 == bilinvar2 && !(quadelem->coef >= 0 ? SCIPisInfinity(scip, -SCIPnlrowGetLhs(nlrow)) : SCIPisInfinity(scip, SCIPnlrowGetRhs(nlrow))) )
2526 {
2527 ++nlcounts[SCIPvarGetProbindex(bilinvar1)];
2528 ++nlcounts[SCIPvarGetProbindex(bilinvar2)];
2529 }
2530
2531 /* bilinear terms are in general non-convex */
2532 if( bilinvar1 != bilinvar2 )
2533 {
2534 ++nlcounts[SCIPvarGetProbindex(bilinvar1)];
2535 ++nlcounts[SCIPvarGetProbindex(bilinvar2)];
2536 }
2537 }
2538
2539 exprtree = SCIPnlrowGetExprtree(nlrow);
2540 if( exprtree != NULL )
2541 {
2542 nexprtreevars = SCIPexprtreeGetNVars(exprtree);
2543 exprtreevars = SCIPexprtreeGetVars(exprtree);
2544
2545 /* assume that the expression tree represents a non-convex constraint */
2546 for( t = 0; t < nexprtreevars; ++t)
2547 {
2548 SCIP_VAR* var;
2549 var = exprtreevars[t];
2550 assert(var != NULL);
2551
2552 ++nlcounts[SCIPvarGetProbindex(var)];
2553 }
2554 }
2555
2556 return SCIP_OKAY;
2557 }
2558
2559 /** count how often each variable appears in a non-convex term */
2560 static
getNLPVarsNonConvexity(SCIP * scip,int * nlcounts)2561 SCIP_RETCODE getNLPVarsNonConvexity(
2562 SCIP* scip, /**< SCIP data structure */
2563 int* nlcounts /**< store the number each variable appears in a
2564 * non-convex term */
2565 )
2566 {
2567 SCIP_CONSHDLR* conshdlr;
2568 SCIP_CONS** conss;
2569 int nvars;
2570 int nconss;
2571 int i;
2572
2573 assert(scip != NULL);
2574 assert(nlcounts != NULL);
2575
2576 nvars = SCIPgetNVars(scip);
2577 BMSclearMemoryArray(nlcounts, nvars);
2578
2579 /* quadratic constraint handler */
2580 conshdlr = SCIPfindConshdlr(scip, "quadratic");
2581 if( conshdlr != NULL )
2582 {
2583 nconss = SCIPconshdlrGetNActiveConss(conshdlr);
2584 conss = SCIPconshdlrGetConss(conshdlr);
2585
2586 SCIPdebugMsg(scip, "nconss(quadratic) = %d\n", nconss);
2587
2588 for( i = 0; i < nconss; ++i )
2589 {
2590 SCIP_Bool isnonconvex;
2591
2592 isnonconvex = (!SCIPisConvexQuadratic(scip, conss[i]) && !SCIPisInfinity(scip, SCIPgetRhsQuadratic(scip, conss[i])))
2593 || (!SCIPisConcaveQuadratic(scip, conss[i]) && !SCIPisInfinity(scip, -SCIPgetLhsQuadratic(scip, conss[i])));
2594
2595 /* only check the nlrow if the constraint is not convex */
2596 if( isnonconvex )
2597 {
2598 SCIP_NLROW* nlrow;
2599 SCIP_CALL( SCIPgetNlRowQuadratic(scip, conss[i], &nlrow) );
2600 assert(nlrow != NULL);
2601
2602 SCIP_CALL( countNLRowVarsNonConvexity(scip, nlcounts, nlrow) );
2603 }
2604 }
2605 }
2606
2607 /* nonlinear constraint handler */
2608 conshdlr = SCIPfindConshdlr(scip, "nonlinear");
2609 if( conshdlr != NULL )
2610 {
2611 nconss = SCIPconshdlrGetNActiveConss(conshdlr);
2612 conss = SCIPconshdlrGetConss(conshdlr);
2613
2614 SCIPdebugMsg(scip, "nconss(nonlinear) = %d\n", nconss);
2615
2616 for( i = 0; i < nconss; ++i )
2617 {
2618 SCIP_EXPRCURV curvature;
2619 SCIP_Bool isnonconvex;
2620
2621 SCIP_CALL( SCIPgetCurvatureNonlinear(scip, conss[i], TRUE, &curvature) );
2622
2623 isnonconvex = (curvature != SCIP_EXPRCURV_CONVEX && !SCIPisInfinity(scip, SCIPgetRhsNonlinear(scip, conss[i])))
2624 || (curvature != SCIP_EXPRCURV_CONCAVE && !SCIPisInfinity(scip, -SCIPgetLhsNonlinear(scip, conss[i])));
2625
2626 /* only check the nlrow if the constraint is not convex */
2627 if( isnonconvex )
2628 {
2629 SCIP_NLROW* nlrow;
2630 SCIP_CALL( SCIPgetNlRowNonlinear(scip, conss[i], &nlrow) );
2631 assert(nlrow != NULL);
2632
2633 SCIP_CALL( countNLRowVarsNonConvexity(scip, nlcounts, nlrow) );
2634 }
2635 }
2636 }
2637
2638 /* bivariate constraint handler */
2639 conshdlr = SCIPfindConshdlr(scip, "bivariate");
2640 if( conshdlr != NULL )
2641 {
2642 nconss = SCIPconshdlrGetNActiveConss(conshdlr);
2643 conss = SCIPconshdlrGetConss(conshdlr);
2644
2645 SCIPdebugMsg(scip, "nconss(bivariate) = %d\n", nconss);
2646
2647 for( i = 0; i < nconss; ++i )
2648 {
2649 SCIP_EXPRCURV curvature;
2650 SCIP_INTERVAL* varbounds;
2651 SCIP_EXPRTREE* exprtree;
2652 int j;
2653
2654 exprtree = SCIPgetExprtreeBivariate(scip, conss[i]);
2655 if( exprtree != NULL )
2656 {
2657 SCIP_Bool isnonconvex;
2658
2659 SCIP_CALL( SCIPallocBufferArray(scip, &varbounds, SCIPexprtreeGetNVars(exprtree)) );
2660 for( j = 0; j < SCIPexprtreeGetNVars(exprtree); ++j )
2661 {
2662 SCIP_VAR* var;
2663 var = SCIPexprtreeGetVars(exprtree)[j];
2664
2665 SCIPintervalSetBounds(&varbounds[j],
2666 -infty2infty(SCIPinfinity(scip), INTERVALINFTY, -MIN(SCIPvarGetLbGlobal(var), SCIPvarGetUbGlobal(var))), /*lint !e666*/
2667 +infty2infty(SCIPinfinity(scip), INTERVALINFTY, MAX(SCIPvarGetLbGlobal(var), SCIPvarGetUbGlobal(var))) ); /*lint !e666*/
2668 }
2669
2670 SCIP_CALL( SCIPexprtreeCheckCurvature(exprtree, SCIPinfinity(scip), varbounds, &curvature, NULL) );
2671
2672 isnonconvex = (curvature != SCIP_EXPRCURV_CONVEX && !SCIPisInfinity(scip, SCIPgetRhsBivariate(scip, conss[i])))
2673 || (curvature != SCIP_EXPRCURV_CONCAVE && !SCIPisInfinity(scip, -SCIPgetLhsBivariate(scip, conss[i])));
2674
2675 /* increase counter for all variables in the expression tree if the constraint is non-convex */
2676 if( isnonconvex )
2677 {
2678 for( j = 0; j < SCIPexprtreeGetNVars(exprtree); ++j )
2679 {
2680 SCIP_VAR* var;
2681 var = SCIPexprtreeGetVars(exprtree)[j];
2682
2683 ++nlcounts[SCIPvarGetProbindex(var)];
2684 }
2685 }
2686 SCIPfreeBufferArray(scip, &varbounds);
2687 }
2688 }
2689 }
2690
2691 /* abspower constraint handler */
2692 conshdlr = SCIPfindConshdlr(scip, "abspower");
2693 if( conshdlr != NULL )
2694 {
2695 nconss = SCIPconshdlrGetNActiveConss(conshdlr);
2696 conss = SCIPconshdlrGetConss(conshdlr);
2697
2698 SCIPdebugMsg(scip, "nconss(abspower) = %d\n", nconss);
2699
2700 for( i = 0; i < nconss; ++i )
2701 {
2702 /* constraint is non-convex in general */
2703 SCIP_NLROW* nlrow;
2704 SCIP_CALL( SCIPgetNlRowAbspower(scip, conss[i], &nlrow) );
2705 assert(nlrow != NULL);
2706
2707 SCIP_CALL( countNLRowVarsNonConvexity(scip, nlcounts, nlrow) );
2708 }
2709 }
2710
2711 return SCIP_OKAY;
2712 }
2713
2714
2715 /** determines whether a variable is interesting */
2716 static
varIsInteresting(SCIP * scip,SCIP_VAR * var,int nlcount)2717 SCIP_Bool varIsInteresting(
2718 SCIP* scip, /**< SCIP data structure */
2719 SCIP_VAR* var, /**< variable to check */
2720 int nlcount /**< number of nonlinear constraints containing the variable
2721 * or number of non-convex terms containing the variable
2722 * (depends on propdata->onlynonconvexvars) */
2723 )
2724 {
2725 assert(SCIPgetDepth(scip) == 0);
2726
2727 return !SCIPvarIsBinary(var) && SCIPvarGetStatus(var) == SCIP_VARSTATUS_COLUMN && nlcount > 0
2728 && !varIsFixedLocal(scip, var);
2729 }
2730
2731 /** initializes interesting bounds */
2732 static
initBounds(SCIP * scip,SCIP_PROPDATA * propdata)2733 SCIP_RETCODE initBounds(
2734 SCIP* scip, /**< SCIP data structure */
2735 SCIP_PROPDATA* propdata /**< data of the obbt propagator */
2736 )
2737 {
2738 SCIP_VAR** vars; /* array of the problems variables */
2739 int* nlcount; /* array that stores in how many nonlinearities each variable appears */
2740 int* nccount; /* array that stores in how many nonconvexities each variable appears */
2741
2742 int bdidx; /* bound index inside propdata->bounds */
2743 int maxnlcount; /* maximal number of nonlinear constraints a variable appears in */
2744 int nvars; /* number of the problems variables */
2745 int i;
2746
2747 assert(scip != NULL);
2748 assert(propdata != NULL);
2749 assert(SCIPisNLPConstructed(scip));
2750
2751 SCIPdebugMsg(scip, "initialize bounds\n");
2752
2753 /* get variable data */
2754 SCIP_CALL( SCIPgetVarsData(scip, &vars, &nvars, NULL, NULL, NULL, NULL) );
2755
2756 /* count nonlinearities */
2757 assert(SCIPgetNNLPVars(scip) == nvars);
2758
2759 SCIP_CALL( SCIPallocBufferArray(scip, &nlcount, nvars) );
2760 SCIP_CALL( SCIPallocBufferArray(scip, &nccount, nvars) );
2761
2762 SCIP_CALL( SCIPgetNLPVarsNonlinearity(scip, nlcount) );
2763 SCIP_CALL( getNLPVarsNonConvexity(scip, nccount) );
2764
2765 maxnlcount = 0;
2766 for( i = 0; i < nvars; i++ )
2767 {
2768 if( maxnlcount < nlcount[i] )
2769 maxnlcount = nlcount[i];
2770 }
2771
2772 /* allocate interesting bounds array */
2773 propdata->boundssize = 2 * nvars;
2774 SCIP_CALL( SCIPallocBlockMemoryArray(scip, &(propdata->bounds), 2 * nvars) );
2775
2776 /* get all interesting variables and their bounds */
2777 bdidx = 0;
2778 for( i = 0; i < nvars; i++ )
2779 {
2780 if( varIsInteresting(scip, vars[i], (propdata->onlynonconvexvars ? nccount[i] : nlcount[i])) )
2781 {
2782 BOUND** bdaddress;
2783
2784 /* create lower bound */
2785 bdaddress = &(propdata->bounds[bdidx]);
2786 SCIP_CALL( SCIPallocBlockMemory(scip, bdaddress) );
2787 propdata->bounds[bdidx]->boundtype = SCIP_BOUNDTYPE_LOWER;
2788 propdata->bounds[bdidx]->var = vars[i];
2789 propdata->bounds[bdidx]->found = FALSE;
2790 propdata->bounds[bdidx]->filtered = FALSE;
2791 propdata->bounds[bdidx]->newval = 0.0;
2792 propdata->bounds[bdidx]->score = getScore(scip, propdata->bounds[bdidx], nlcount[i], maxnlcount);
2793 propdata->bounds[bdidx]->done = FALSE;
2794 propdata->bounds[bdidx]->nonconvex = (nccount[i] > 0);
2795 propdata->bounds[bdidx]->index = bdidx;
2796 bdidx++;
2797
2798 /* create upper bound */
2799 bdaddress = &(propdata->bounds[bdidx]);
2800 SCIP_CALL( SCIPallocBlockMemory(scip, bdaddress) );
2801 propdata->bounds[bdidx]->boundtype = SCIP_BOUNDTYPE_UPPER;
2802 propdata->bounds[bdidx]->var = vars[i];
2803 propdata->bounds[bdidx]->found = FALSE;
2804 propdata->bounds[bdidx]->filtered = FALSE;
2805 propdata->bounds[bdidx]->newval = 0.0;
2806 propdata->bounds[bdidx]->score = getScore(scip, propdata->bounds[bdidx], nlcount[i], maxnlcount);
2807 propdata->bounds[bdidx]->done = FALSE;
2808 propdata->bounds[bdidx]->nonconvex = (nccount[i] > 0);
2809 propdata->bounds[bdidx]->index = bdidx;
2810 bdidx++;
2811 }
2812 }
2813
2814 /* set number of interesting bounds */
2815 propdata->nbounds = bdidx;
2816
2817 /* collect all bilinear terms from quadratic constraint handler */
2818 if( propdata->nbounds > 0 && SCIPgetNAllBilinearTermsQuadratic(scip) > 0 && propdata->createbilinineqs )
2819 {
2820 SCIP_VAR** x;
2821 SCIP_VAR** y;
2822 SCIP_Real* maxnonconvexity;
2823 int* nunderest;
2824 int* noverest;
2825 int nbilins;
2826 int bilinidx;
2827 int nbilinterms;
2828
2829 nbilins = SCIPgetNAllBilinearTermsQuadratic(scip);
2830 bilinidx = 0;
2831 nbilinterms = 0;
2832
2833 /* allocate memory */
2834 SCIP_CALL( SCIPallocBufferArray(scip, &x, nbilins) );
2835 SCIP_CALL( SCIPallocBufferArray(scip, &y, nbilins) );
2836 SCIP_CALL( SCIPallocBufferArray(scip, &nunderest, nbilins) );
2837 SCIP_CALL( SCIPallocBufferArray(scip, &noverest, nbilins) );
2838 SCIP_CALL( SCIPallocBufferArray(scip, &maxnonconvexity, nbilins) );
2839
2840 /* get data for bilinear terms */
2841 SCIP_CALL( SCIPgetAllBilinearTermsQuadratic(scip, x, y, &nbilins, nunderest, noverest, maxnonconvexity) );
2842
2843 /* count the number of interesting bilinear terms */
2844 propdata->nbilinbounds = 0;
2845 for( i = 0; i < nbilins; ++i )
2846 if( nunderest[i] + noverest[i] > 0 && propdata->minnonconvexity <= maxnonconvexity[i]
2847 && varIsInteresting(scip, x[i], 1) && varIsInteresting(scip, y[i], 1) )
2848 ++(propdata->nbilinbounds);
2849
2850 if( propdata->nbilinbounds == 0 )
2851 goto TERMINATE;
2852
2853 SCIP_CALL( SCIPallocBlockMemoryArray(scip, &(propdata->bilinbounds), propdata->nbilinbounds) );
2854 BMSclearMemoryArray(propdata->bilinbounds, propdata->nbilinbounds);
2855
2856 for( i = 0; i < nbilins; ++i )
2857 {
2858 if( nunderest[i] + noverest[i] > 0 && propdata->minnonconvexity <= maxnonconvexity[i]
2859 && varIsInteresting(scip, x[i], 1) && varIsInteresting(scip, y[i], 1) )
2860 {
2861 SCIP_CALL( SCIPallocBlockMemory(scip, &propdata->bilinbounds[bilinidx]) ); /*lint !e866*/
2862 BMSclearMemory(propdata->bilinbounds[bilinidx]); /*lint !e866*/
2863
2864 propdata->bilinbounds[bilinidx]->x = x[i];
2865 propdata->bilinbounds[bilinidx]->y = y[i];
2866 propdata->bilinbounds[bilinidx]->nunderest = nunderest[i];
2867 propdata->bilinbounds[bilinidx]->noverest = noverest[i];
2868 propdata->bilinbounds[bilinidx]->index = i;
2869 ++bilinidx;
2870
2871 /* count how often bilinear terms appear in quadratic constraints */
2872 nbilinterms += nunderest[i] + noverest[i];
2873 }
2874 }
2875 assert(propdata->nbilinbounds == bilinidx);
2876
2877 /* compute scores for each term */
2878 for( i = 0; i < propdata->nbilinbounds; ++i )
2879 {
2880 propdata->bilinbounds[i]->score = getScoreBilinBound(scip, propdata->randnumgen, propdata->bilinbounds[i],
2881 nbilinterms);
2882 SCIPdebugMsg(scip, "score of %i = %g\n", i, propdata->bilinbounds[i]->score);
2883 }
2884
2885 /* sort bounds according to decreasing score */
2886 if( propdata->nbilinbounds > 1 )
2887 {
2888 SCIPsortDownPtr((void**) propdata->bilinbounds, compBilinboundsScore, propdata->nbilinbounds);
2889 }
2890
2891 TERMINATE:
2892 /* free memory */
2893 SCIPfreeBufferArray(scip, &maxnonconvexity);
2894 SCIPfreeBufferArray(scip, &noverest);
2895 SCIPfreeBufferArray(scip, &nunderest);
2896 SCIPfreeBufferArray(scip, &y);
2897 SCIPfreeBufferArray(scip, &x);
2898 }
2899
2900 /* free memory for buffering nonlinearities */
2901 assert(nlcount != NULL);
2902 assert(nccount != NULL);
2903 SCIPfreeBufferArray(scip, &nccount);
2904 SCIPfreeBufferArray(scip, &nlcount);
2905
2906 /* propdata->bounds array if empty */
2907 if( propdata->nbounds <= 0 )
2908 {
2909 assert(propdata->nbounds == 0);
2910 assert(propdata->boundssize >= 0 );
2911 SCIPfreeBlockMemoryArray(scip, &(propdata->bounds), propdata->boundssize);
2912 }
2913
2914 SCIPdebugMsg(scip, "problem has %d/%d interesting bounds\n", propdata->nbounds, 2 * nvars);
2915
2916 if( propdata->nbounds > 0 )
2917 {
2918 /* sort bounds according to decreasing score; although this initial order will be overruled by the distance
2919 * criterion later, gives a more well-defined starting situation for OBBT and might help to reduce solver
2920 * variability
2921 */
2922 SCIPsortDownPtr((void**) propdata->bounds, compBoundsScore, propdata->nbounds);
2923 }
2924
2925 return SCIP_OKAY;
2926 }
2927
2928 /*
2929 * Callback methods of propagator
2930 */
2931
2932 /** copy method for propagator plugins (called when SCIP copies plugins)
2933 *
2934 * @note The UG framework assumes that all default plug-ins of SCIP implement a copy callback. We check
2935 * SCIPgetSubscipDepth() in PROPEXEC to prevent the propagator to run in a sub-SCIP.
2936 */
2937 static
SCIP_DECL_PROPCOPY(propCopyObbt)2938 SCIP_DECL_PROPCOPY(propCopyObbt)
2939 { /*lint --e{715}*/
2940 assert(scip != NULL);
2941 assert(prop != NULL);
2942 assert(strcmp(SCIPpropGetName(prop), PROP_NAME) == 0);
2943
2944 /* call inclusion method of constraint handler */
2945 SCIP_CALL( SCIPincludePropObbt(scip) );
2946
2947 return SCIP_OKAY;
2948 }
2949
2950 /** solving process initialization method of propagator (called when branch and bound process is about to begin) */
2951 static
SCIP_DECL_PROPINITSOL(propInitsolObbt)2952 SCIP_DECL_PROPINITSOL(propInitsolObbt)
2953 { /*lint --e{715}*/
2954 SCIP_PROPDATA* propdata;
2955
2956 assert(scip != NULL);
2957 assert(prop != NULL);
2958 assert(strcmp(SCIPpropGetName(prop), PROP_NAME) == 0);
2959
2960 /* get propagator data */
2961 propdata = SCIPpropGetData(prop);
2962 assert(propdata != NULL);
2963
2964 propdata->bounds = NULL;
2965 propdata->nbounds = -1;
2966 propdata->boundssize = 0;
2967 propdata->cutoffrow = NULL;
2968 propdata->lastnode = -1;
2969
2970 /* if genvbounds propagator is not available, we cannot create genvbounds */
2971 propdata->genvboundprop = propdata->creategenvbounds ? SCIPfindProp(scip, GENVBOUND_PROP_NAME) : NULL;
2972
2973 SCIPdebugMsg(scip, "creating genvbounds: %s\n", propdata->genvboundprop != NULL ? "true" : "false");
2974
2975 /* create random number generator */
2976 SCIP_CALL( SCIPcreateRandom(scip, &propdata->randnumgen, DEFAULT_RANDSEED, TRUE) );
2977
2978 return SCIP_OKAY;
2979 }
2980
2981 /** execution method of propagator */
2982 static
SCIP_DECL_PROPEXEC(propExecObbt)2983 SCIP_DECL_PROPEXEC(propExecObbt)
2984 { /*lint --e{715}*/
2985 SCIP_PROPDATA* propdata;
2986 SCIP_Longint itlimit;
2987
2988 assert(scip != NULL);
2989 assert(prop != NULL);
2990 assert(strcmp(SCIPpropGetName(prop), PROP_NAME) == 0);
2991
2992 *result = SCIP_DIDNOTRUN;
2993
2994 /* do not run in: presolving, repropagation, probing mode, if no objective propagation is allowed */
2995 if( SCIPgetStage(scip) != SCIP_STAGE_SOLVING || SCIPinRepropagation(scip) || SCIPinProbing(scip) || !SCIPallowWeakDualReds(scip) )
2996 return SCIP_OKAY;
2997
2998 /* do not run propagator in a sub-SCIP */
2999 if( SCIPgetSubscipDepth(scip) > 0 )
3000 return SCIP_OKAY;
3001
3002 /* only run for nonlinear problems, i.e., if NLP is constructed */
3003 if( !SCIPisNLPConstructed(scip) )
3004 {
3005 SCIPdebugMsg(scip, "NLP not constructed, skipping obbt\n");
3006 return SCIP_OKAY;
3007 }
3008
3009 /* only run if LP all columns are in the LP, i.e., the LP is a relaxation; e.g., do not run if pricers are active
3010 * since pricing is not performed in probing mode
3011 */
3012 if( !SCIPallColsInLP(scip) )
3013 {
3014 SCIPdebugMsg(scip, "not all columns in LP, skipping obbt\n");
3015 return SCIP_OKAY;
3016 }
3017
3018 if( !SCIPallowWeakDualReds(scip) )
3019 return SCIP_OKAY;
3020
3021 /* get propagator data */
3022 propdata = SCIPpropGetData(prop);
3023 assert(propdata != NULL);
3024
3025 /* ensure that bounds are initialized */
3026 if( propdata->nbounds == -1 )
3027 {
3028 /* bounds must be initialized at root node */
3029 if( SCIPgetDepth(scip) == 0 )
3030 {
3031 SCIP_CALL( initBounds(scip, propdata) );
3032 }
3033 else
3034 {
3035 assert(!SCIPinProbing(scip));
3036 return SCIP_OKAY;
3037 }
3038 }
3039 assert(propdata->nbounds >= 0);
3040
3041 /* do not run if there are no interesting bounds */
3042 /**@todo disable */
3043 if( propdata->nbounds <= 0 )
3044 {
3045 SCIPdebugMsg(scip, "there are no interesting bounds\n");
3046 return SCIP_OKAY;
3047 }
3048
3049 /* only run once in a node != root */
3050 if( SCIPgetDepth(scip) > 0 && SCIPnodeGetNumber(SCIPgetCurrentNode(scip)) == propdata->lastnode )
3051 {
3052 return SCIP_OKAY;
3053 }
3054
3055 SCIPdebugMsg(scip, "applying obbt for problem <%s> at depth %d\n", SCIPgetProbName(scip), SCIPgetDepth(scip));
3056
3057 /* without an optimal LP solution we don't want to run; this may be because propagators with higher priority have
3058 * already found reductions or numerical troubles occured during LP solving
3059 */
3060 if( SCIPgetLPSolstat(scip) != SCIP_LPSOLSTAT_OPTIMAL && SCIPgetLPSolstat(scip) != SCIP_LPSOLSTAT_UNBOUNDEDRAY )
3061 {
3062 SCIPdebugMsg(scip, "aborting since no optimal LP solution is at hand\n");
3063 return SCIP_OKAY;
3064 }
3065
3066 /* compute iteration limit */
3067 if( propdata->itlimitfactor > 0.0 )
3068 itlimit = (SCIP_Longint) MAX(propdata->itlimitfactor * SCIPgetNRootLPIterations(scip),
3069 propdata->minitlimit); /*lint !e666*/
3070 else
3071 itlimit = -1;
3072
3073 /* apply obbt */
3074 SCIP_CALL( applyObbt(scip, propdata, itlimit, result) );
3075 assert(*result != SCIP_DIDNOTRUN);
3076
3077 /* compute globally inequalities for bilinear terms */
3078 if( propdata->createbilinineqs )
3079 {
3080 /* set LP iteration limit */
3081 if( propdata->itlimitbilin == 0L )
3082 {
3083 /* no iteration limit if itlimit < 0 or itlimitfactorbilin < 0 */
3084 propdata->itlimitbilin = (itlimit < 0 || propdata->itlimitfactorbilin < 0)
3085 ? -1L : (SCIP_Longint)(itlimit * propdata->itlimitfactorbilin);
3086 }
3087
3088 SCIP_CALL( applyObbtBilinear(scip, propdata, itlimit, result) );
3089 }
3090
3091 /* set current node as last node */
3092 propdata->lastnode = SCIPnodeGetNumber(SCIPgetCurrentNode(scip));
3093
3094 return SCIP_OKAY;
3095 }
3096
3097 /** propagation conflict resolving method of propagator */
3098 static
SCIP_DECL_PROPRESPROP(propRespropObbt)3099 SCIP_DECL_PROPRESPROP(propRespropObbt)
3100 { /*lint --e{715}*/
3101 *result = SCIP_DIDNOTFIND;
3102
3103 return SCIP_OKAY;
3104 }
3105
3106 /** solving process deinitialization method of propagator (called before branch and bound process data is freed) */
3107 static
SCIP_DECL_PROPEXITSOL(propExitsolObbt)3108 SCIP_DECL_PROPEXITSOL(propExitsolObbt)
3109 { /*lint --e{715}*/
3110 SCIP_PROPDATA* propdata;
3111 int i;
3112
3113 assert(scip != NULL);
3114 assert(prop != NULL);
3115 assert(strcmp(SCIPpropGetName(prop), PROP_NAME) == 0);
3116
3117 /* get propagator data */
3118 propdata = SCIPpropGetData(prop);
3119 assert(propdata != NULL);
3120
3121 /* free random number generator */
3122 SCIPfreeRandom(scip, &propdata->randnumgen);
3123 propdata->randnumgen = NULL;
3124
3125 /* free bilinear bounds */
3126 if( propdata->nbilinbounds > 0 )
3127 {
3128 for( i = propdata->nbilinbounds - 1; i >= 0; --i )
3129 {
3130 SCIPfreeBlockMemory(scip, &propdata->bilinbounds[i]); /*lint !e866*/
3131 }
3132 SCIPfreeBlockMemoryArray(scip, &propdata->bilinbounds, propdata->nbilinbounds);
3133 propdata->nbilinbounds = 0;
3134 }
3135
3136 /* free memory allocated for the bounds */
3137 if( propdata->nbounds > 0 )
3138 {
3139 /* free bounds */
3140 for( i = propdata->nbounds - 1; i >= 0; i-- )
3141 {
3142 SCIPfreeBlockMemory(scip, &(propdata->bounds[i])); /*lint !e866*/
3143 }
3144 SCIPfreeBlockMemoryArray(scip, &(propdata->bounds), propdata->boundssize);
3145 }
3146
3147 /* reset variables */
3148 propdata->lastidx = -1;
3149 propdata->lastbilinidx = 0;
3150 propdata->propagatecounter = 0;
3151 propdata->npropagatedomreds = 0;
3152 propdata->nbounds = -1;
3153 propdata->itlimitbilin = 0;
3154 propdata->itusedbilin = 0;
3155
3156 return SCIP_OKAY;
3157 }
3158
3159 /** destructor of propagator to free user data (called when SCIP is exiting) */
3160 static
SCIP_DECL_PROPFREE(propFreeObbt)3161 SCIP_DECL_PROPFREE(propFreeObbt)
3162 { /*lint --e{715}*/
3163 SCIP_PROPDATA* propdata;
3164
3165 assert(strcmp(SCIPpropGetName(prop), PROP_NAME) == 0);
3166
3167 /* free propagator data */
3168 propdata = SCIPpropGetData(prop);
3169 assert(propdata != NULL);
3170
3171 SCIPfreeBlockMemory(scip, &propdata);
3172
3173 SCIPpropSetData(prop, NULL);
3174
3175 return SCIP_OKAY;
3176 }
3177
3178
3179 /*
3180 * propagator specific interface methods
3181 */
3182
3183 /** creates the obbt propagator and includes it in SCIP */
SCIPincludePropObbt(SCIP * scip)3184 SCIP_RETCODE SCIPincludePropObbt(
3185 SCIP* scip /**< SCIP data structure */
3186 )
3187 {
3188 SCIP_PROPDATA* propdata;
3189 SCIP_PROP* prop;
3190
3191 /* create obbt propagator data */
3192 SCIP_CALL( SCIPallocBlockMemory(scip, &propdata) );
3193 BMSclearMemory(propdata);
3194
3195 /* initialize variables with a non-zero default value */
3196 propdata->lastidx = -1;
3197
3198 /* include propagator */
3199 SCIP_CALL( SCIPincludePropBasic(scip, &prop, PROP_NAME, PROP_DESC, PROP_PRIORITY, PROP_FREQ, PROP_DELAY, PROP_TIMING,
3200 propExecObbt, propdata) );
3201
3202 SCIP_CALL( SCIPsetPropCopy(scip, prop, propCopyObbt) );
3203 SCIP_CALL( SCIPsetPropFree(scip, prop, propFreeObbt) );
3204 SCIP_CALL( SCIPsetPropExitsol(scip, prop, propExitsolObbt) );
3205 SCIP_CALL( SCIPsetPropInitsol(scip, prop, propInitsolObbt) );
3206 SCIP_CALL( SCIPsetPropResprop(scip, prop, propRespropObbt) );
3207
3208 SCIP_CALL( SCIPaddBoolParam(scip, "propagating/" PROP_NAME "/creategenvbounds",
3209 "should obbt try to provide genvbounds if possible?",
3210 &propdata->creategenvbounds, TRUE, DEFAULT_CREATE_GENVBOUNDS, NULL, NULL) );
3211
3212 SCIP_CALL( SCIPaddBoolParam(scip, "propagating/" PROP_NAME "/normalize",
3213 "should coefficients in filtering be normalized w.r.t. the domains sizes?",
3214 &propdata->normalize, TRUE, DEFAULT_FILTERING_NORM, NULL, NULL) );
3215
3216 SCIP_CALL( SCIPaddBoolParam(scip, "propagating/" PROP_NAME "/applyfilterrounds",
3217 "try to filter bounds in so-called filter rounds by solving auxiliary LPs?",
3218 &propdata->applyfilterrounds, TRUE, DEFAULT_APPLY_FILTERROUNDS, NULL, NULL) );
3219
3220 SCIP_CALL( SCIPaddBoolParam(scip, "propagating/" PROP_NAME "/applytrivialfilter",
3221 "try to filter bounds with the LP solution after each solve?",
3222 &propdata->applytrivialfilter, TRUE, DEFAULT_APPLY_TRIVIALFITLERING, NULL, NULL) );
3223
3224 SCIP_CALL( SCIPaddBoolParam(scip, "propagating/" PROP_NAME "/genvbdsduringfilter",
3225 "should we try to generate genvbounds during trivial and aggressive filtering?",
3226 &propdata->genvbdsduringfilter, TRUE, DEFAULT_GENVBDSDURINGFILTER, NULL, NULL) );
3227
3228 SCIP_CALL( SCIPaddBoolParam(scip, "propagating/" PROP_NAME "/genvbdsduringsepa",
3229 "try to create genvbounds during separation process?",
3230 &propdata->genvbdsduringsepa, TRUE, DEFAULT_GENVBDSDURINGSEPA, NULL, NULL) );
3231
3232 SCIP_CALL( SCIPaddIntParam(scip, "propagating/" PROP_NAME "/minfilter",
3233 "minimal number of filtered bounds to apply another filter round",
3234 &propdata->nminfilter, TRUE, DEFAULT_FILTERING_MIN, 1, INT_MAX, NULL, NULL) );
3235
3236 SCIP_CALL( SCIPaddRealParam(scip, "propagating/" PROP_NAME "/itlimitfactor",
3237 "multiple of root node LP iterations used as total LP iteration limit for obbt (<= 0: no limit )",
3238 &propdata->itlimitfactor, FALSE, DEFAULT_ITLIMITFACTOR, SCIP_REAL_MIN, SCIP_REAL_MAX, NULL, NULL) );
3239
3240 SCIP_CALL( SCIPaddRealParam(scip, "propagating/" PROP_NAME "/itlimitfactorbilin",
3241 "multiple of OBBT LP limit used as total LP iteration limit for solving bilinear inequality LPs (< 0 for no limit)",
3242 &propdata->itlimitfactorbilin, FALSE, DEFAULT_ITLIMITFAC_BILININEQS, SCIP_REAL_MIN, SCIP_REAL_MAX, NULL, NULL) );
3243
3244 SCIP_CALL( SCIPaddRealParam(scip, "propagating/" PROP_NAME "/minnonconvexity",
3245 "minimum absolute value of nonconvex eigenvalues for a bilinear term",
3246 &propdata->minnonconvexity, FALSE, DEFAULT_MINNONCONVEXITY, 0.0, SCIP_REAL_MAX, NULL, NULL) );
3247
3248 SCIP_CALL( SCIPaddLongintParam(scip, "propagating/" PROP_NAME "/minitlimit",
3249 "minimum LP iteration limit",
3250 &propdata->minitlimit, FALSE, DEFAULT_MINITLIMIT, 0L, SCIP_LONGINT_MAX, NULL, NULL) );
3251
3252 SCIP_CALL( SCIPaddRealParam(scip, "propagating/" PROP_NAME "/dualfeastol",
3253 "feasibility tolerance for reduced costs used in obbt; this value is used if SCIP's dual feastol is greater",
3254 &propdata->dualfeastol, FALSE, DEFAULT_DUALFEASTOL, 0.0, SCIP_REAL_MAX, NULL, NULL) );
3255
3256 SCIP_CALL( SCIPaddRealParam(scip, "propagating/" PROP_NAME "/conditionlimit",
3257 "maximum condition limit used in LP solver (-1.0: no limit)",
3258 &propdata->conditionlimit, FALSE, DEFAULT_CONDITIONLIMIT, -1.0, SCIP_REAL_MAX, NULL, NULL) );
3259
3260 SCIP_CALL( SCIPaddRealParam(scip, "propagating/" PROP_NAME "/boundstreps",
3261 "minimal relative improve for strengthening bounds",
3262 &propdata->boundstreps, FALSE, DEFAULT_BOUNDSTREPS, 0.0, 1.0, NULL, NULL) );
3263
3264 SCIP_CALL( SCIPaddBoolParam(scip, "propagating/" PROP_NAME "/onlynonconvexvars",
3265 "only apply obbt on non-convex variables",
3266 &propdata->onlynonconvexvars, TRUE, DEFAULT_ONLYNONCONVEXVARS, NULL, NULL) );
3267
3268 SCIP_CALL( SCIPaddBoolParam(scip, "propagating/" PROP_NAME "/tightintboundsprobing",
3269 "should integral bounds be tightened during the probing mode?",
3270 &propdata->tightintboundsprobing, TRUE, DEFAULT_TIGHTINTBOUNDSPROBING, NULL, NULL) );
3271
3272 SCIP_CALL( SCIPaddBoolParam(scip, "propagating/" PROP_NAME "/tightcontboundsprobing",
3273 "should continuous bounds be tightened during the probing mode?",
3274 &propdata->tightcontboundsprobing, TRUE, DEFAULT_TIGHTCONTBOUNDSPROBING, NULL, NULL) );
3275
3276 SCIP_CALL( SCIPaddBoolParam(scip, "propagating/" PROP_NAME "/createbilinineqs",
3277 "solve auxiliary LPs in order to find valid inequalities for bilinear terms?",
3278 &propdata->createbilinineqs, TRUE, DEFAULT_CREATE_BILININEQS, NULL, NULL) );
3279
3280 SCIP_CALL( SCIPaddIntParam(scip, "propagating/" PROP_NAME "/orderingalgo",
3281 "select the type of ordering algorithm which should be used (0: no special ordering, 1: greedy, 2: greedy reverse)",
3282 &propdata->orderingalgo, TRUE, DEFAULT_ORDERINGALGO, 0, 2, NULL, NULL) );
3283
3284 SCIP_CALL( SCIPaddBoolParam(scip, "propagating/" PROP_NAME "/separatesol",
3285 "should the obbt LP solution be separated?",
3286 &propdata->separatesol, TRUE, DEFAULT_SEPARATESOL, NULL, NULL) );
3287
3288 SCIP_CALL( SCIPaddIntParam(scip, "propagating/" PROP_NAME "/sepaminiter",
3289 "minimum number of iteration spend to separate an obbt LP solution",
3290 &propdata->sepaminiter, TRUE, DEFAULT_SEPAMINITER, 0, INT_MAX, NULL, NULL) );
3291
3292 SCIP_CALL( SCIPaddIntParam(scip, "propagating/" PROP_NAME "/sepamaxiter",
3293 "maximum number of iteration spend to separate an obbt LP solution",
3294 &propdata->sepamaxiter, TRUE, DEFAULT_SEPAMAXITER, 0, INT_MAX, NULL, NULL) );
3295
3296 SCIP_CALL( SCIPaddIntParam(scip, "propagating/" PROP_NAME "/propagatefreq",
3297 "trigger a propagation round after that many bound tightenings (0: no propagation)",
3298 &propdata->propagatefreq, TRUE, DEFAULT_PROPAGATEFREQ, 0, INT_MAX, NULL, NULL) );
3299
3300 return SCIP_OKAY;
3301 }
3302