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 heur_rounding.c 26 * @ingroup DEFPLUGINS_HEUR 27 * @brief LP rounding heuristic that tries to recover from intermediate infeasibilities 28 * @author Tobias Achterberg 29 */ 30 31 /*---+----1----+----2----+----3----+----4----+----5----+----6----+----7----+----8----+----9----+----0----+----1----+----2*/ 32 33 #include "blockmemshell/memory.h" 34 #include "scip/heur_rounding.h" 35 #include "scip/pub_heur.h" 36 #include "scip/pub_lp.h" 37 #include "scip/pub_message.h" 38 #include "scip/pub_var.h" 39 #include "scip/scip_branch.h" 40 #include "scip/scip_heur.h" 41 #include "scip/scip_lp.h" 42 #include "scip/scip_mem.h" 43 #include "scip/scip_message.h" 44 #include "scip/scip_numerics.h" 45 #include "scip/scip_param.h" 46 #include "scip/scip_sol.h" 47 #include "scip/scip_solvingstats.h" 48 #include <string.h> 49 50 #define HEUR_NAME "rounding" 51 #define HEUR_DESC "LP rounding heuristic with infeasibility recovering" 52 #define HEUR_DISPCHAR SCIP_HEURDISPCHAR_ROUNDING 53 #define HEUR_PRIORITY -1000 54 #define HEUR_FREQ 1 55 #define HEUR_FREQOFS 0 56 #define HEUR_MAXDEPTH -1 57 #define HEUR_TIMING SCIP_HEURTIMING_DURINGLPLOOP 58 #define HEUR_USESSUBSCIP FALSE /**< does the heuristic use a secondary SCIP instance? */ 59 60 #define DEFAULT_SUCCESSFACTOR 100 /**< number of calls per found solution that are considered as standard success, 61 * a higher factor causes the heuristic to be called more often 62 */ 63 #define DEFAULT_ONCEPERNODE FALSE /**< should the heuristic only be called once per node? */ 64 65 /* locally defined heuristic data */ 66 struct SCIP_HeurData 67 { 68 SCIP_SOL* sol; /**< working solution */ 69 SCIP_Longint lastlp; /**< last LP number where the heuristic was applied */ 70 int successfactor; /**< number of calls per found solution that are considered as standard success, 71 * a higher factor causes the heuristic to be called more often 72 */ 73 SCIP_Bool oncepernode; /**< should the heuristic only be called once per node? */ 74 }; 75 76 77 /* 78 * local methods 79 */ 80 81 /** update row violation arrays after a row's activity value changed */ 82 static 83 void updateViolations( 84 SCIP* scip, /**< SCIP data structure */ 85 SCIP_ROW* row, /**< LP row */ 86 SCIP_ROW** violrows, /**< array with currently violated rows */ 87 int* violrowpos, /**< position of LP rows in violrows array */ 88 int* nviolrows, /**< pointer to the number of currently violated rows */ 89 SCIP_Real oldactivity, /**< old activity value of LP row */ 90 SCIP_Real newactivity /**< new activity value of LP row */ 91 ) 92 { 93 SCIP_Real lhs; 94 SCIP_Real rhs; 95 SCIP_Bool oldviol; 96 SCIP_Bool newviol; 97 98 assert(violrows != NULL); 99 assert(violrowpos != NULL); 100 assert(nviolrows != NULL); 101 102 lhs = SCIProwGetLhs(row); 103 rhs = SCIProwGetRhs(row); 104 oldviol = (SCIPisFeasLT(scip, oldactivity, lhs) || SCIPisFeasGT(scip, oldactivity, rhs)); 105 newviol = (SCIPisFeasLT(scip, newactivity, lhs) || SCIPisFeasGT(scip, newactivity, rhs)); 106 if( oldviol != newviol ) 107 { 108 int rowpos; 109 int violpos; 110 111 rowpos = SCIProwGetLPPos(row); 112 assert(rowpos >= 0); 113 114 if( oldviol ) 115 { 116 /* the row violation was repaired: remove row from violrows array, decrease violation count */ 117 violpos = violrowpos[rowpos]; 118 assert(0 <= violpos && violpos < *nviolrows); 119 assert(violrows[violpos] == row); 120 violrowpos[rowpos] = -1; 121 if( violpos != *nviolrows-1 ) 122 { 123 violrows[violpos] = violrows[*nviolrows-1]; 124 violrowpos[SCIProwGetLPPos(violrows[violpos])] = violpos; 125 } 126 (*nviolrows)--; 127 } 128 else 129 { 130 /* the row is now violated: add row to violrows array, increase violation count */ 131 assert(violrowpos[rowpos] == -1); 132 violrows[*nviolrows] = row; 133 violrowpos[rowpos] = *nviolrows; 134 (*nviolrows)++; 135 } 136 } 137 } 138 139 /** update row activities after a variable's solution value changed */ 140 static 141 SCIP_RETCODE updateActivities( 142 SCIP* scip, /**< SCIP data structure */ 143 SCIP_Real* activities, /**< LP row activities */ 144 SCIP_ROW** violrows, /**< array with currently violated rows */ 145 int* violrowpos, /**< position of LP rows in violrows array */ 146 int* nviolrows, /**< pointer to the number of currently violated rows */ 147 int nlprows, /**< number of rows in current LP */ 148 SCIP_VAR* var, /**< variable that has been changed */ 149 SCIP_Real oldsolval, /**< old solution value of variable */ 150 SCIP_Real newsolval /**< new solution value of variable */ 151 ) 152 { 153 SCIP_COL* col; 154 SCIP_ROW** colrows; 155 SCIP_Real* colvals; 156 SCIP_Real delta; 157 int ncolrows; 158 int r; 159 160 assert(activities != NULL); 161 assert(nviolrows != NULL); 162 assert(0 <= *nviolrows && *nviolrows <= nlprows); 163 164 delta = newsolval - oldsolval; 165 col = SCIPvarGetCol(var); 166 colrows = SCIPcolGetRows(col); 167 colvals = SCIPcolGetVals(col); 168 ncolrows = SCIPcolGetNLPNonz(col); 169 assert(ncolrows == 0 || (colrows != NULL && colvals != NULL)); 170 171 for( r = 0; r < ncolrows; ++r ) 172 { 173 SCIP_ROW* row; 174 int rowpos; 175 176 row = colrows[r]; 177 rowpos = SCIProwGetLPPos(row); 178 assert(-1 <= rowpos && rowpos < nlprows); 179 180 if( rowpos >= 0 && !SCIProwIsLocal(row) ) 181 { 182 SCIP_Real oldactivity; 183 SCIP_Real newactivity; 184 185 assert(SCIProwIsInLP(row)); 186 187 /* update row activity */ 188 oldactivity = activities[rowpos]; 189 if( !SCIPisInfinity(scip, -oldactivity) && !SCIPisInfinity(scip, oldactivity) ) 190 { 191 newactivity = oldactivity + delta * colvals[r]; 192 if( SCIPisInfinity(scip, newactivity) ) 193 newactivity = SCIPinfinity(scip); 194 else if( SCIPisInfinity(scip, -newactivity) ) 195 newactivity = -SCIPinfinity(scip); 196 activities[rowpos] = newactivity; 197 198 /* update row violation arrays */ 199 updateViolations(scip, row, violrows, violrowpos, nviolrows, oldactivity, newactivity); 200 } 201 } 202 } 203 204 return SCIP_OKAY; 205 } 206 207 /** returns a variable, that pushes activity of the row in the given direction with minimal negative impact on other rows; 208 * if variables have equal impact, chooses the one with best objective value improvement in corresponding direction; 209 * rounding in a direction is forbidden, if this forces the objective value over the upper bound 210 */ 211 static 212 SCIP_RETCODE selectRounding( 213 SCIP* scip, /**< SCIP data structure */ 214 SCIP_SOL* sol, /**< primal solution */ 215 SCIP_Real minobj, /**< minimal objective value possible after rounding remaining fractional vars */ 216 SCIP_ROW* row, /**< LP row */ 217 int direction, /**< should the activity be increased (+1) or decreased (-1)? */ 218 SCIP_VAR** roundvar, /**< pointer to store the rounding variable, returns NULL if impossible */ 219 SCIP_Real* oldsolval, /**< pointer to store old (fractional) solution value of rounding variable */ 220 SCIP_Real* newsolval /**< pointer to store new (rounded) solution value of rounding variable */ 221 ) 222 { 223 SCIP_COL* col; 224 SCIP_VAR* var; 225 SCIP_Real val; 226 SCIP_COL** rowcols; 227 SCIP_Real* rowvals; 228 SCIP_Real solval; 229 SCIP_Real roundval; 230 SCIP_Real obj; 231 SCIP_Real deltaobj; 232 SCIP_Real bestdeltaobj; 233 SCIP_VARTYPE vartype; 234 int nrowcols; 235 int nlocks; 236 int minnlocks; 237 int c; 238 239 assert(direction == +1 || direction == -1); 240 assert(roundvar != NULL); 241 assert(oldsolval != NULL); 242 assert(newsolval != NULL); 243 244 /* get row entries */ 245 rowcols = SCIProwGetCols(row); 246 rowvals = SCIProwGetVals(row); 247 nrowcols = SCIProwGetNLPNonz(row); 248 249 /* select rounding variable */ 250 minnlocks = INT_MAX; 251 bestdeltaobj = SCIPinfinity(scip); 252 *roundvar = NULL; 253 for( c = 0; c < nrowcols; ++c ) 254 { 255 col = rowcols[c]; 256 var = SCIPcolGetVar(col); 257 258 vartype = SCIPvarGetType(var); 259 if( vartype == SCIP_VARTYPE_BINARY || vartype == SCIP_VARTYPE_INTEGER ) 260 { 261 solval = SCIPgetSolVal(scip, sol, var); 262 263 if( !SCIPisFeasIntegral(scip, solval) ) 264 { 265 val = rowvals[c]; 266 obj = SCIPvarGetObj(var); 267 268 if( direction * val < 0.0 ) 269 { 270 /* rounding down */ 271 nlocks = SCIPvarGetNLocksDownType(var, SCIP_LOCKTYPE_MODEL); 272 if( nlocks <= minnlocks ) 273 { 274 roundval = SCIPfeasFloor(scip, solval); 275 deltaobj = obj * (roundval - solval); 276 if( (nlocks < minnlocks || deltaobj < bestdeltaobj) && minobj - obj < SCIPgetCutoffbound(scip) ) 277 { 278 minnlocks = nlocks; 279 bestdeltaobj = deltaobj; 280 *roundvar = var; 281 *oldsolval = solval; 282 *newsolval = roundval; 283 } 284 } 285 } 286 else 287 { 288 /* rounding up */ 289 assert(direction * val > 0.0); 290 nlocks = SCIPvarGetNLocksUpType(var, SCIP_LOCKTYPE_MODEL); 291 if( nlocks <= minnlocks ) 292 { 293 roundval = SCIPfeasCeil(scip, solval); 294 deltaobj = obj * (roundval - solval); 295 if( (nlocks < minnlocks || deltaobj < bestdeltaobj) && minobj + obj < SCIPgetCutoffbound(scip) ) 296 { 297 minnlocks = nlocks; 298 bestdeltaobj = deltaobj; 299 *roundvar = var; 300 *oldsolval = solval; 301 *newsolval = roundval; 302 } 303 } 304 } 305 } 306 } 307 } 308 309 return SCIP_OKAY; 310 } 311 312 /** returns a variable, that increases the activity of the row */ 313 static 314 SCIP_RETCODE selectIncreaseRounding( 315 SCIP* scip, /**< SCIP data structure */ 316 SCIP_SOL* sol, /**< primal solution */ 317 SCIP_Real minobj, /**< minimal objective value possible after rounding remaining fractional vars */ 318 SCIP_ROW* row, /**< LP row */ 319 SCIP_VAR** roundvar, /**< pointer to store the rounding variable, returns NULL if impossible */ 320 SCIP_Real* oldsolval, /**< old (fractional) solution value of rounding variable */ 321 SCIP_Real* newsolval /**< new (rounded) solution value of rounding variable */ 322 ) 323 { 324 return selectRounding(scip, sol, minobj, row, +1, roundvar, oldsolval, newsolval); 325 } 326 327 /** returns a variable, that decreases the activity of the row */ 328 static 329 SCIP_RETCODE selectDecreaseRounding( 330 SCIP* scip, /**< SCIP data structure */ 331 SCIP_SOL* sol, /**< primal solution */ 332 SCIP_Real minobj, /**< minimal objective value possible after rounding remaining fractional vars */ 333 SCIP_ROW* row, /**< LP row */ 334 SCIP_VAR** roundvar, /**< pointer to store the rounding variable, returns NULL if impossible */ 335 SCIP_Real* oldsolval, /**< old (fractional) solution value of rounding variable */ 336 SCIP_Real* newsolval /**< new (rounded) solution value of rounding variable */ 337 ) 338 { 339 return selectRounding(scip, sol, minobj, row, -1, roundvar, oldsolval, newsolval); 340 } 341 342 /** returns a fractional variable, that has most impact on rows in opposite direction, i.e. that is most crucial to 343 * fix in the other direction; 344 * if variables have equal impact, chooses the one with best objective value improvement in corresponding direction; 345 * rounding in a direction is forbidden, if this forces the objective value over the upper bound 346 */ 347 static 348 SCIP_RETCODE selectEssentialRounding( 349 SCIP* scip, /**< SCIP data structure */ 350 SCIP_SOL* sol, /**< primal solution */ 351 SCIP_Real minobj, /**< minimal objective value possible after rounding remaining fractional vars */ 352 SCIP_VAR** lpcands, /**< fractional variables in LP */ 353 int nlpcands, /**< number of fractional variables in LP */ 354 SCIP_VAR** roundvar, /**< pointer to store the rounding variable, returns NULL if impossible */ 355 SCIP_Real* oldsolval, /**< old (fractional) solution value of rounding variable */ 356 SCIP_Real* newsolval /**< new (rounded) solution value of rounding variable */ 357 ) 358 { 359 SCIP_VAR* var; 360 SCIP_Real solval; 361 SCIP_Real roundval; 362 SCIP_Real obj; 363 SCIP_Real deltaobj; 364 SCIP_Real bestdeltaobj; 365 int maxnlocks; 366 int nlocks; 367 int v; 368 369 assert(roundvar != NULL); 370 assert(oldsolval != NULL); 371 assert(newsolval != NULL); 372 373 /* select rounding variable */ 374 maxnlocks = -1; 375 bestdeltaobj = SCIPinfinity(scip); 376 *roundvar = NULL; 377 for( v = 0; v < nlpcands; ++v ) 378 { 379 var = lpcands[v]; 380 assert(SCIPvarGetType(var) == SCIP_VARTYPE_BINARY || SCIPvarGetType(var) == SCIP_VARTYPE_INTEGER); 381 382 solval = SCIPgetSolVal(scip, sol, var); 383 if( !SCIPisFeasIntegral(scip, solval) ) 384 { 385 obj = SCIPvarGetObj(var); 386 387 /* rounding down */ 388 nlocks = SCIPvarGetNLocksUpType(var, SCIP_LOCKTYPE_MODEL); 389 if( nlocks >= maxnlocks ) 390 { 391 roundval = SCIPfeasFloor(scip, solval); 392 deltaobj = obj * (roundval - solval); 393 if( (nlocks > maxnlocks || deltaobj < bestdeltaobj) && minobj - obj < SCIPgetCutoffbound(scip) ) 394 { 395 maxnlocks = nlocks; 396 bestdeltaobj = deltaobj; 397 *roundvar = var; 398 *oldsolval = solval; 399 *newsolval = roundval; 400 } 401 } 402 403 /* rounding up */ 404 nlocks = SCIPvarGetNLocksDownType(var, SCIP_LOCKTYPE_MODEL); 405 if( nlocks >= maxnlocks ) 406 { 407 roundval = SCIPfeasCeil(scip, solval); 408 deltaobj = obj * (roundval - solval); 409 if( (nlocks > maxnlocks || deltaobj < bestdeltaobj) && minobj + obj < SCIPgetCutoffbound(scip) ) 410 { 411 maxnlocks = nlocks; 412 bestdeltaobj = deltaobj; 413 *roundvar = var; 414 *oldsolval = solval; 415 *newsolval = roundval; 416 } 417 } 418 } 419 } 420 421 return SCIP_OKAY; 422 } 423 424 425 /* 426 * Callback methods 427 */ 428 429 /** copy method for primal heuristic plugins (called when SCIP copies plugins) */ 430 static 431 SCIP_DECL_HEURCOPY(heurCopyRounding) 432 { /*lint --e{715}*/ 433 assert(scip != NULL); 434 assert(heur != NULL); 435 assert(strcmp(SCIPheurGetName(heur), HEUR_NAME) == 0); 436 437 /* call inclusion method of primal heuristic */ 438 SCIP_CALL( SCIPincludeHeurRounding(scip) ); 439 440 return SCIP_OKAY; 441 } 442 443 /** destructor of primal heuristic to free user data (called when SCIP is exiting) */ 444 static 445 SCIP_DECL_HEURFREE(heurFreeRounding) /*lint --e{715}*/ 446 { /*lint --e{715}*/ 447 SCIP_HEURDATA* heurdata; 448 449 assert(heur != NULL); 450 assert(strcmp(SCIPheurGetName(heur), HEUR_NAME) == 0); 451 assert(scip != NULL); 452 453 /* free heuristic data */ 454 heurdata = SCIPheurGetData(heur); 455 assert(heurdata != NULL); 456 SCIPfreeBlockMemory(scip, &heurdata); 457 SCIPheurSetData(heur, NULL); 458 459 return SCIP_OKAY; 460 } 461 462 /** initialization method of primal heuristic (called after problem was transformed) */ 463 static 464 SCIP_DECL_HEURINIT(heurInitRounding) /*lint --e{715}*/ 465 { /*lint --e{715}*/ 466 SCIP_HEURDATA* heurdata; 467 468 assert(strcmp(SCIPheurGetName(heur), HEUR_NAME) == 0); 469 470 heurdata = SCIPheurGetData(heur); 471 assert(heurdata != NULL); 472 473 /* create heuristic data */ 474 SCIP_CALL( SCIPcreateSol(scip, &heurdata->sol, heur) ); 475 heurdata->lastlp = -1; 476 477 return SCIP_OKAY; 478 } 479 480 /** deinitialization method of primal heuristic (called before transformed problem is freed) */ 481 static 482 SCIP_DECL_HEUREXIT(heurExitRounding) /*lint --e{715}*/ 483 { /*lint --e{715}*/ 484 SCIP_HEURDATA* heurdata; 485 486 assert(strcmp(SCIPheurGetName(heur), HEUR_NAME) == 0); 487 488 /* free heuristic data */ 489 heurdata = SCIPheurGetData(heur); 490 assert(heurdata != NULL); 491 SCIP_CALL( SCIPfreeSol(scip, &heurdata->sol) ); 492 493 return SCIP_OKAY; 494 } 495 496 /** solving process initialization method of primal heuristic (called when branch and bound process is about to begin) */ 497 static 498 SCIP_DECL_HEURINITSOL(heurInitsolRounding) 499 { 500 SCIP_HEURDATA* heurdata; 501 502 assert(strcmp(SCIPheurGetName(heur), HEUR_NAME) == 0); 503 504 heurdata = SCIPheurGetData(heur); 505 assert(heurdata != NULL); 506 heurdata->lastlp = -1; 507 508 /* change the heuristic's timingmask, if nit should be called only once per node */ 509 if( heurdata->oncepernode ) 510 SCIPheurSetTimingmask(heur, SCIP_HEURTIMING_AFTERLPNODE); 511 512 return SCIP_OKAY; 513 } 514 515 /** solving process deinitialization method of primal heuristic (called before branch and bound process data is freed) */ 516 static 517 SCIP_DECL_HEUREXITSOL(heurExitsolRounding) 518 { 519 /* reset the timing mask to its default value */ 520 SCIPheurSetTimingmask(heur, HEUR_TIMING); 521 522 return SCIP_OKAY; 523 } 524 525 526 /** execution method of primal heuristic */ 527 static 528 SCIP_DECL_HEUREXEC(heurExecRounding) /*lint --e{715}*/ 529 { /*lint --e{715}*/ 530 SCIP_HEURDATA* heurdata; 531 SCIP_SOL* sol; 532 SCIP_VAR** lpcands; 533 SCIP_Real* lpcandssol; 534 SCIP_ROW** lprows; 535 SCIP_Real* activities; 536 SCIP_ROW** violrows; 537 int* violrowpos; 538 SCIP_Real obj; 539 SCIP_Real bestroundval; 540 SCIP_Real minobj; 541 int nlpcands; 542 int nlprows; 543 int nfrac; 544 int nviolrows; 545 int c; 546 int r; 547 SCIP_Longint nlps; 548 SCIP_Longint ncalls; 549 SCIP_Longint nsolsfound; 550 SCIP_Longint nnodes; 551 552 assert(strcmp(SCIPheurGetName(heur), HEUR_NAME) == 0); 553 assert(scip != NULL); 554 assert(result != NULL); 555 assert(SCIPhasCurrentNodeLP(scip)); 556 557 *result = SCIP_DIDNOTRUN; 558 559 /* only call heuristic, if an optimal LP solution is at hand */ 560 if( SCIPgetLPSolstat(scip) != SCIP_LPSOLSTAT_OPTIMAL ) 561 return SCIP_OKAY; 562 563 /* only call heuristic, if the LP objective value is smaller than the cutoff bound */ 564 if( SCIPisGE(scip, SCIPgetLPObjval(scip), SCIPgetCutoffbound(scip)) ) 565 return SCIP_OKAY; 566 567 /* get heuristic data */ 568 heurdata = SCIPheurGetData(heur); 569 assert(heurdata != NULL); 570 571 /* don't call heuristic, if we have already processed the current LP solution */ 572 nlps = SCIPgetNLPs(scip); 573 if( nlps == heurdata->lastlp ) 574 return SCIP_OKAY; 575 heurdata->lastlp = nlps; 576 577 /* don't call heuristic, if it was not successful enough in the past */ 578 ncalls = SCIPheurGetNCalls(heur); 579 nsolsfound = 10*SCIPheurGetNBestSolsFound(heur) + SCIPheurGetNSolsFound(heur); 580 nnodes = SCIPgetNNodes(scip); 581 if( nnodes % ((ncalls/heurdata->successfactor)/(nsolsfound+1)+1) != 0 ) 582 return SCIP_OKAY; 583 584 /* get fractional variables, that should be integral */ 585 SCIP_CALL( SCIPgetLPBranchCands(scip, &lpcands, &lpcandssol, NULL, &nlpcands, NULL, NULL) ); 586 nfrac = nlpcands; 587 588 /* only call heuristic, if LP solution is fractional */ 589 if( nfrac == 0 ) 590 return SCIP_OKAY; 591 592 *result = SCIP_DIDNOTFIND; 593 594 /* get LP rows */ 595 SCIP_CALL( SCIPgetLPRowsData(scip, &lprows, &nlprows) ); 596 597 SCIPdebugMsg(scip, "executing rounding heuristic: %d LP rows, %d fractionals\n", nlprows, nfrac); 598 599 /* get memory for activities, violated rows, and row violation positions */ 600 SCIP_CALL( SCIPallocBufferArray(scip, &activities, nlprows) ); 601 SCIP_CALL( SCIPallocBufferArray(scip, &violrows, nlprows) ); 602 SCIP_CALL( SCIPallocBufferArray(scip, &violrowpos, nlprows) ); 603 604 /* get the activities for all globally valid rows; 605 * the rows should be feasible, but due to numerical inaccuracies in the LP solver, they can be violated 606 */ 607 nviolrows = 0; 608 for( r = 0; r < nlprows; ++r ) 609 { 610 SCIP_ROW* row; 611 612 row = lprows[r]; 613 assert(SCIProwGetLPPos(row) == r); 614 615 if( !SCIProwIsLocal(row) ) 616 { 617 activities[r] = SCIPgetRowActivity(scip, row); 618 if( SCIPisFeasLT(scip, activities[r], SCIProwGetLhs(row)) 619 || SCIPisFeasGT(scip, activities[r], SCIProwGetRhs(row)) ) 620 { 621 violrows[nviolrows] = row; 622 violrowpos[r] = nviolrows; 623 nviolrows++; 624 } 625 else 626 violrowpos[r] = -1; 627 } 628 } 629 630 /* get the working solution from heuristic's local data */ 631 sol = heurdata->sol; 632 assert(sol != NULL); 633 634 /* copy the current LP solution to the working solution */ 635 SCIP_CALL( SCIPlinkLPSol(scip, sol) ); 636 637 /* calculate the minimal objective value possible after rounding fractional variables */ 638 minobj = SCIPgetSolTransObj(scip, sol); 639 assert(minobj < SCIPgetCutoffbound(scip)); 640 for( c = 0; c < nlpcands; ++c ) 641 { 642 obj = SCIPvarGetObj(lpcands[c]); 643 bestroundval = obj > 0.0 ? SCIPfeasFloor(scip, lpcandssol[c]) : SCIPfeasCeil(scip, lpcandssol[c]); 644 minobj += obj * (bestroundval - lpcandssol[c]); 645 } 646 647 /* try to round remaining variables in order to become/stay feasible */ 648 while( nfrac > 0 ) 649 { 650 SCIP_VAR* roundvar; 651 SCIP_Real oldsolval; 652 SCIP_Real newsolval; 653 654 SCIPdebugMsg(scip, "rounding heuristic: nfrac=%d, nviolrows=%d, obj=%g (best possible obj: %g)\n", 655 nfrac, nviolrows, SCIPgetSolOrigObj(scip, sol), SCIPretransformObj(scip, minobj)); 656 657 /* minobj < SCIPgetCutoffbound(scip) should be true, otherwise the rounding variable selection 658 * should have returned NULL. Due to possible cancellation we use SCIPisLE. */ 659 assert( SCIPisLE(scip, minobj, SCIPgetCutoffbound(scip)) ); 660 661 /* choose next variable to process: 662 * - if a violated row exists, round a variable decreasing the violation, that has least impact on other rows 663 * - otherwise, round a variable, that has strongest devastating impact on rows in opposite direction 664 */ 665 if( nviolrows > 0 ) 666 { 667 SCIP_ROW* row; 668 int rowpos; 669 670 row = violrows[nviolrows-1]; 671 rowpos = SCIProwGetLPPos(row); 672 assert(0 <= rowpos && rowpos < nlprows); 673 assert(violrowpos[rowpos] == nviolrows-1); 674 675 SCIPdebugMsg(scip, "rounding heuristic: try to fix violated row <%s>: %g <= %g <= %g\n", 676 SCIProwGetName(row), SCIProwGetLhs(row), activities[rowpos], SCIProwGetRhs(row)); 677 if( SCIPisFeasLT(scip, activities[rowpos], SCIProwGetLhs(row)) ) 678 { 679 /* lhs is violated: select a variable rounding, that increases the activity */ 680 SCIP_CALL( selectIncreaseRounding(scip, sol, minobj, row, &roundvar, &oldsolval, &newsolval) ); 681 } 682 else 683 { 684 assert(SCIPisFeasGT(scip, activities[rowpos], SCIProwGetRhs(row))); 685 /* rhs is violated: select a variable rounding, that decreases the activity */ 686 SCIP_CALL( selectDecreaseRounding(scip, sol, minobj, row, &roundvar, &oldsolval, &newsolval) ); 687 } 688 } 689 else 690 { 691 SCIPdebugMsg(scip, "rounding heuristic: search rounding variable and try to stay feasible\n"); 692 SCIP_CALL( selectEssentialRounding(scip, sol, minobj, lpcands, nlpcands, &roundvar, &oldsolval, &newsolval) ); 693 } 694 695 /* check, whether rounding was possible */ 696 if( roundvar == NULL ) 697 { 698 SCIPdebugMsg(scip, "rounding heuristic: -> didn't find a rounding variable\n"); 699 break; 700 } 701 702 SCIPdebugMsg(scip, "rounding heuristic: -> round var <%s>, oldval=%g, newval=%g, obj=%g\n", 703 SCIPvarGetName(roundvar), oldsolval, newsolval, SCIPvarGetObj(roundvar)); 704 705 /* update row activities of globally valid rows */ 706 SCIP_CALL( updateActivities(scip, activities, violrows, violrowpos, &nviolrows, nlprows, 707 roundvar, oldsolval, newsolval) ); 708 709 /* store new solution value and decrease fractionality counter */ 710 SCIP_CALL( SCIPsetSolVal(scip, sol, roundvar, newsolval) ); 711 nfrac--; 712 713 /* update minimal objective value possible after rounding remaining variables */ 714 obj = SCIPvarGetObj(roundvar); 715 if( obj > 0.0 && newsolval > oldsolval ) 716 minobj += obj; 717 else if( obj < 0.0 && newsolval < oldsolval ) 718 minobj -= obj; 719 720 SCIPdebugMsg(scip, "rounding heuristic: -> nfrac=%d, nviolrows=%d, obj=%g (best possible obj: %g)\n", 721 nfrac, nviolrows, SCIPgetSolOrigObj(scip, sol), SCIPretransformObj(scip, minobj)); 722 } 723 724 /* check, if the new solution is feasible */ 725 if( nfrac == 0 && nviolrows == 0 ) 726 { 727 SCIP_Bool stored; 728 729 /* check solution for feasibility, and add it to solution store if possible 730 * neither integrality nor feasibility of LP rows has to be checked, because this is already 731 * done in the rounding heuristic itself; however, be better check feasibility of LP rows, 732 * because of numerical problems with activity updating 733 */ 734 SCIP_CALL( SCIPtrySol(scip, sol, FALSE, FALSE, FALSE, FALSE, TRUE, &stored) ); 735 736 if( stored ) 737 { 738 #ifdef SCIP_DEBUG 739 SCIPdebugMsg(scip, "found feasible rounded solution:\n"); 740 SCIP_CALL( SCIPprintSol(scip, sol, NULL, FALSE) ); 741 #endif 742 *result = SCIP_FOUNDSOL; 743 } 744 } 745 746 /* free memory buffers */ 747 SCIPfreeBufferArray(scip, &violrowpos); 748 SCIPfreeBufferArray(scip, &violrows); 749 SCIPfreeBufferArray(scip, &activities); 750 751 return SCIP_OKAY; 752 } 753 754 755 /* 756 * heuristic specific interface methods 757 */ 758 759 /** creates the rounding heuristic with infeasibility recovering and includes it in SCIP */ 760 SCIP_RETCODE SCIPincludeHeurRounding( 761 SCIP* scip /**< SCIP data structure */ 762 ) 763 { 764 SCIP_HEURDATA* heurdata; 765 SCIP_HEUR* heur; 766 767 /* create Rounding primal heuristic data */ 768 SCIP_CALL( SCIPallocBlockMemory(scip, &heurdata) ); 769 770 /* include primal heuristic */ 771 SCIP_CALL( SCIPincludeHeurBasic(scip, &heur, 772 HEUR_NAME, HEUR_DESC, HEUR_DISPCHAR, HEUR_PRIORITY, HEUR_FREQ, HEUR_FREQOFS, 773 HEUR_MAXDEPTH, HEUR_TIMING, HEUR_USESSUBSCIP, heurExecRounding, heurdata) ); 774 775 assert(heur != NULL); 776 777 /* set non-NULL pointers to callback methods */ 778 SCIP_CALL( SCIPsetHeurCopy(scip, heur, heurCopyRounding) ); 779 SCIP_CALL( SCIPsetHeurFree(scip, heur, heurFreeRounding) ); 780 SCIP_CALL( SCIPsetHeurInit(scip, heur, heurInitRounding) ); 781 SCIP_CALL( SCIPsetHeurExit(scip, heur, heurExitRounding) ); 782 SCIP_CALL( SCIPsetHeurInitsol(scip, heur, heurInitsolRounding) ); 783 SCIP_CALL( SCIPsetHeurExitsol(scip, heur, heurExitsolRounding) ); 784 785 /* add rounding primal heuristic parameters */ 786 SCIP_CALL( SCIPaddIntParam(scip, "heuristics/" HEUR_NAME "/successfactor", 787 "number of calls per found solution that are considered as standard success, a higher factor causes the heuristic to be called more often", 788 &heurdata->successfactor, TRUE, DEFAULT_SUCCESSFACTOR, -1, INT_MAX, NULL, NULL) ); 789 790 SCIP_CALL( SCIPaddBoolParam(scip, "heuristics/" HEUR_NAME "/oncepernode", 791 "should the heuristic only be called once per node?", 792 &heurdata->oncepernode, TRUE, DEFAULT_ONCEPERNODE, NULL, NULL) ); 793 794 return SCIP_OKAY; 795 } 796