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