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 sepa_gomory.c 26 * @ingroup DEFPLUGINS_SEPA 27 * @brief Gomory MIR Cuts 28 * @author Tobias Achterberg 29 * @author Stefan Heinz 30 * @author Domenico Salvagnin 31 * @author Marc Pfetsch 32 * @author Leona Gottwald 33 */ 34 35 /**@todo try k-Gomory-cuts (s. Cornuejols: K-Cuts: A Variation of Gomory Mixed Integer Cuts from the LP Tableau) 36 * 37 * @todo Try cuts on the objective tableau row. 38 * 39 * @todo Also try negative basis inverse row? 40 * 41 * @todo It happens that the SCIPcalcMIR() function returns with the same cut for different calls. Check if this is a 42 * bug or do not use it for the MIP below and turn off presolving and all heuristics: 43 * 44 * Max y 45 * Subject to 46 * c1: -x + y <= 1 47 * c2: 2x + 3y <= 12 48 * c3: 3x + 2y <= 12 49 * Bounds 50 * 0 <= x 51 * 0 <= y 52 * General 53 * x 54 * y 55 * END 56 */ 57 58 /*---+----1----+----2----+----3----+----4----+----5----+----6----+----7----+----8----+----9----+----0----+----1----+----2*/ 59 60 #include "blockmemshell/memory.h" 61 #include "scip/cuts.h" 62 #include "scip/pub_lp.h" 63 #include "scip/pub_message.h" 64 #include "scip/pub_misc.h" 65 #include "scip/pub_misc_sort.h" 66 #include "scip/pub_sepa.h" 67 #include "scip/pub_var.h" 68 #include "scip/scip_branch.h" 69 #include "scip/scip_cut.h" 70 #include "scip/scip_general.h" 71 #include "scip/scip_lp.h" 72 #include "scip/scip_mem.h" 73 #include "scip/scip_message.h" 74 #include "scip/scip_numerics.h" 75 #include "scip/scip_param.h" 76 #include "scip/scip_prob.h" 77 #include "scip/scip_randnumgen.h" 78 #include "scip/scip_sepa.h" 79 #include "scip/scip_solvingstats.h" 80 #include "scip/scip_tree.h" 81 #include "scip/scip_var.h" 82 #include "scip/sepa_gomory.h" 83 #include "struct_history.h" 84 #include <string.h> 85 86 #define SEPA_NAME "gomory" 87 #define SEPA_DESC "separator for Gomory mixed-integer and strong CG cuts from LP tableau rows" 88 #define SEPA_PRIORITY -1000 89 #define SEPA_FREQ 10 90 #define SEPA_MAXBOUNDDIST 1.0 91 #define SEPA_USESSUBSCIP FALSE /**< does the separator use a secondary SCIP instance? */ 92 #define SEPA_DELAY FALSE /**< should separation method be delayed, if other separators found cuts? */ 93 94 #define DEFAULT_MAXROUNDS 5 /**< maximal number of gomory separation rounds per node (-1: unlimited) */ 95 #define DEFAULT_MAXROUNDSROOT 10 /**< maximal number of gomory separation rounds in the root node (-1: unlimited) */ 96 #define DEFAULT_MAXSEPACUTS 50 /**< maximal number of gomory cuts separated per separation round */ 97 #define DEFAULT_MAXSEPACUTSROOT 200 /**< maximal number of gomory cuts separated per separation round in root node */ 98 #define DEFAULT_MAXRANK -1 /**< maximal rank of a gomory cut that could not be scaled to integral coefficients (-1: unlimited) */ 99 #define DEFAULT_MAXRANKINTEGRAL -1 /**< maximal rank of a gomory cut that could be scaled to integral coefficients (-1: unlimited) */ 100 #define DEFAULT_DYNAMICCUTS TRUE /**< should generated cuts be removed from the LP if they are no longer tight? */ 101 #define DEFAULT_AWAY 0.01 /**< minimal integrality violation of a basis variable in order to try Gomory cut */ 102 #define DEFAULT_MAKEINTEGRAL FALSE /**< try to scale all cuts to integral coefficients */ 103 #define DEFAULT_FORCECUTS TRUE /**< if conversion to integral coefficients failed still consider the cut */ 104 #define DEFAULT_SEPARATEROWS TRUE /**< separate rows with integral slack */ 105 #define DEFAULT_DELAYEDCUTS FALSE /**< should cuts be added to the delayed cut pool? */ 106 #define DEFAULT_SIDETYPEBASIS TRUE /**< choose side types of row (lhs/rhs) based on basis information? */ 107 #define DEFAULT_TRYSTRONGCG TRUE /**< try to generate strengthened Chvatal-Gomory cuts? */ 108 #define DEFAULT_GENBOTHGOMSCG TRUE /**< should both Gomory and strong CG cuts be generated (otherwise take best) */ 109 #define DEFAULT_RANDSEED 53 /**< initial random seed */ 110 111 #define BOUNDSWITCH 0.9999 /**< threshold for bound switching - see SCIPcalcMIR() */ 112 #define POSTPROCESS TRUE /**< apply postprocessing after MIR calculation - see SCIPcalcMIR() */ 113 #define USEVBDS TRUE /**< use variable bounds - see SCIPcalcMIR() */ 114 #define FIXINTEGRALRHS FALSE /**< try to generate an integral rhs - see SCIPcalcMIR() */ 115 #define MAKECONTINTEGRAL FALSE /**< convert continuous variable to integral variables in SCIPmakeRowIntegral() */ 116 117 #define MAXAGGRLEN(nvars) (0.1*(nvars)+1000) /**< maximal length of base inequality */ 118 119 120 /** separator data */ 121 struct SCIP_SepaData 122 { 123 SCIP_RANDNUMGEN* randnumgen; /**< random number generator */ 124 SCIP_SEPA* strongcg; /**< strong CG cut separator */ 125 SCIP_SEPA* gomory; /**< gomory cut separator */ 126 SCIP_Real away; /**< minimal integrality violation of a basis variable in order to try Gomory cut */ 127 int maxrounds; /**< maximal number of gomory separation rounds per node (-1: unlimited) */ 128 int maxroundsroot; /**< maximal number of gomory separation rounds in the root node (-1: unlimited) */ 129 int maxsepacuts; /**< maximal number of gomory cuts separated per separation round */ 130 int maxsepacutsroot; /**< maximal number of gomory cuts separated per separation round in root node */ 131 int maxrank; /**< maximal rank of a gomory cut that could not be scaled to integral coefficients (-1: unlimited) */ 132 int maxrankintegral; /**< maximal rank of a gomory cut that could be scaled to integral coefficients (-1: unlimited) */ 133 int lastncutsfound; /**< total number of cuts found after last call of separator */ 134 SCIP_Bool dynamiccuts; /**< should generated cuts be removed from the LP if they are no longer tight? */ 135 SCIP_Bool makeintegral; /**< try to scale all cuts to integral coefficients */ 136 SCIP_Bool forcecuts; /**< if conversion to integral coefficients failed still consider the cut */ 137 SCIP_Bool separaterows; /**< separate rows with integral slack */ 138 SCIP_Bool delayedcuts; /**< should cuts be added to the delayed cut pool? */ 139 SCIP_Bool sidetypebasis; /**< choose side types of row (lhs/rhs) based on basis information? */ 140 SCIP_Bool trystrongcg; /**< try to generate strengthened Chvatal-Gomory cuts? */ 141 SCIP_Bool genbothgomscg; /**< should both Gomory and strong CG cuts be generated (otherwise take best) */ 142 }; 143 144 145 /** returns TRUE if the cut can be taken, otherwise FALSE if there some numerical evidences */ 146 static 147 SCIP_RETCODE evaluateCutNumerics( 148 SCIP* scip, /**< SCIP data structure */ 149 SCIP_SEPADATA* sepadata, /**< data of the separator */ 150 SCIP_ROW* cut, /**< cut to check */ 151 SCIP_Longint maxdnom, /**< maximal denominator to use for scaling */ 152 SCIP_Real maxscale, /**< maximal scaling factor */ 153 SCIP_Bool* useful /**< pointer to store if the cut is useful */ 154 ) 155 { 156 SCIP_Bool madeintegral = FALSE; 157 158 assert(useful != NULL); 159 160 *useful = FALSE; 161 162 if( sepadata->makeintegral && SCIPgetRowNumIntCols(scip, cut) == SCIProwGetNNonz(cut) ) 163 { 164 /* try to scale the cut to integral values */ 165 SCIP_CALL( SCIPmakeRowIntegral(scip, cut, -SCIPepsilon(scip), SCIPsumepsilon(scip), 166 maxdnom, maxscale, MAKECONTINTEGRAL, &madeintegral) ); 167 168 if( !madeintegral && !sepadata->forcecuts ) 169 return SCIP_OKAY; 170 171 /* in case the right hand side is plus infinity (due to scaling) the cut is useless so we are not taking it at all */ 172 if( madeintegral && SCIPisInfinity(scip, SCIProwGetRhs(cut)) ) 173 return SCIP_OKAY; 174 } 175 176 /* discard integral cut if the rank is too high */ 177 if( madeintegral && sepadata->maxrankintegral != -1 && (SCIProwGetRank(cut) > sepadata->maxrankintegral) ) 178 return SCIP_OKAY; 179 180 /* discard cut if the rank is too high */ 181 if( !madeintegral && (sepadata->maxrank != -1) && (SCIProwGetRank(cut) > sepadata->maxrank) ) 182 return SCIP_OKAY; 183 184 *useful = TRUE; 185 186 return SCIP_OKAY; 187 } 188 189 190 /** add cut */ 191 static 192 SCIP_RETCODE addCut( 193 SCIP* scip, /**< SCIP instance */ 194 SCIP_SEPADATA* sepadata, /**< separator data */ 195 SCIP_VAR** vars, /**< array of variables */ 196 int c, /**< index of basic variable (< 0 for slack variables) */ 197 SCIP_Longint maxdnom, /**< maximal denominator to use for scaling */ 198 SCIP_Real maxscale, /**< maximal scaling factor */ 199 int cutnnz, /**< number of nonzeros in cut */ 200 int* cutinds, /**< variable indices in cut */ 201 SCIP_Real* cutcoefs, /**< cut cofficients */ 202 SCIP_Real cutefficacy, /**< cut efficacy */ 203 SCIP_Real cutrhs, /**< rhs of cut */ 204 SCIP_Bool cutislocal, /**< whether cut is local */ 205 int cutrank, /**< rank of cut */ 206 SCIP_Bool strongcg, /**< whether the cut arises from the strong-CG procedure */ 207 SCIP_Bool* cutoff, /**< pointer to store whether a cutoff appeared */ 208 int* naddedcuts /**< pointer to store number of added cuts */ 209 ) 210 { 211 assert(scip != NULL); 212 assert(cutoff != NULL); 213 assert(naddedcuts != NULL); 214 215 if( cutnnz == 0 && SCIPisFeasNegative(scip, cutrhs) ) /*lint !e644*/ 216 { 217 SCIPdebugMsg(scip, " -> gomory cut detected infeasibility with cut 0 <= %g.\n", cutrhs); 218 *cutoff = TRUE; 219 return SCIP_OKAY; 220 } 221 222 /* Only take efficient cuts, except for cuts with one non-zero coefficient (= bound 223 * changes); the latter cuts will be handled internally in sepastore. */ 224 if( SCIPisEfficacious(scip, cutefficacy) || ( cutnnz == 1 && SCIPisFeasPositive(scip, cutefficacy) ) ) 225 { 226 SCIP_ROW* cut; 227 SCIP_SEPA* cutsepa; 228 char cutname[SCIP_MAXSTRLEN]; 229 int v; 230 231 /* construct cut name */ 232 if( strongcg ) 233 { 234 cutsepa = sepadata->strongcg; 235 236 if( c >= 0 ) 237 (void) SCIPsnprintf(cutname, SCIP_MAXSTRLEN, "scg%" SCIP_LONGINT_FORMAT "_x%d", SCIPgetNLPs(scip), c); 238 else 239 (void) SCIPsnprintf(cutname, SCIP_MAXSTRLEN, "scg%" SCIP_LONGINT_FORMAT "_s%d", SCIPgetNLPs(scip), -c-1); 240 } 241 else 242 { 243 cutsepa = sepadata->gomory; 244 245 if( c >= 0 ) 246 (void) SCIPsnprintf(cutname, SCIP_MAXSTRLEN, "gom%" SCIP_LONGINT_FORMAT "_x%d", SCIPgetNLPs(scip), c); 247 else 248 (void) SCIPsnprintf(cutname, SCIP_MAXSTRLEN, "gom%" SCIP_LONGINT_FORMAT "_s%d", SCIPgetNLPs(scip), -c-1); 249 } 250 251 /* create empty cut */ 252 SCIP_CALL( SCIPcreateEmptyRowSepa(scip, &cut, cutsepa, cutname, -SCIPinfinity(scip), cutrhs, 253 cutislocal, FALSE, sepadata->dynamiccuts) ); 254 255 /* set cut rank */ 256 SCIProwChgRank(cut, cutrank); /*lint !e644*/ 257 258 /* cache the row extension and only flush them if the cut gets added */ 259 SCIP_CALL( SCIPcacheRowExtensions(scip, cut) ); 260 261 /* collect all non-zero coefficients */ 262 for( v = 0; v < cutnnz; ++v ) 263 { 264 SCIP_CALL( SCIPaddVarToRow(scip, cut, vars[cutinds[v]], cutcoefs[v]) ); 265 } 266 267 /* flush all changes before adding the cut */ 268 SCIP_CALL( SCIPflushRowExtensions(scip, cut) ); 269 270 if( SCIProwGetNNonz(cut) == 0 ) 271 { 272 assert( SCIPisFeasNegative(scip, cutrhs) ); 273 SCIPdebugMsg(scip, " -> gomory cut detected infeasibility with cut 0 <= %g.\n", cutrhs); 274 *cutoff = TRUE; 275 return SCIP_OKAY; 276 } 277 else if( SCIProwGetNNonz(cut) == 1 ) 278 { 279 /* Add the bound change as cut to avoid that the LP gets modified. This would mean that the LP is not flushed 280 * and the method SCIPgetLPBInvRow() fails; SCIP internally will apply this bound change automatically. */ 281 SCIP_CALL( SCIPaddRow(scip, cut, TRUE, cutoff) ); 282 ++(*naddedcuts); 283 } 284 else 285 { 286 SCIP_Bool useful; 287 288 assert(SCIPisInfinity(scip, -SCIProwGetLhs(cut))); 289 assert(!SCIPisInfinity(scip, SCIProwGetRhs(cut))); 290 291 SCIPdebugMsg(scip, " -> %s cut <%s>: rhs=%f, eff=%f\n", strongcg ? "strong-CG" : "gomory", cutname, cutrhs, cutefficacy); 292 293 SCIP_CALL( evaluateCutNumerics(scip, sepadata, cut, maxdnom, maxscale, &useful) ); 294 295 if( useful ) 296 { 297 SCIPdebugMsg(scip, " -> found %s cut <%s>: act=%f, rhs=%f, norm=%f, eff=%f, min=%f, max=%f (range=%f)\n", 298 strongcg ? "strong-CG" : "gomory", cutname, SCIPgetRowLPActivity(scip, cut), SCIProwGetRhs(cut), 299 SCIProwGetNorm(cut), SCIPgetCutEfficacy(scip, NULL, cut), 300 SCIPgetRowMinCoef(scip, cut), SCIPgetRowMaxCoef(scip, cut), 301 SCIPgetRowMaxCoef(scip, cut)/SCIPgetRowMinCoef(scip, cut)); 302 303 if( SCIPisCutNew(scip, cut) ) 304 { 305 /* add global cuts which are not implicit bound changes to the cut pool */ 306 if( !cutislocal ) 307 { 308 if( sepadata->delayedcuts ) 309 { 310 SCIP_CALL( SCIPaddDelayedPoolCut(scip, cut) ); 311 } 312 else 313 { 314 SCIP_CALL( SCIPaddPoolCut(scip, cut) ); 315 } 316 } 317 else 318 { 319 /* local cuts we add to the sepastore */ 320 SCIP_CALL( SCIPaddRow(scip, cut, FALSE, cutoff) ); 321 } 322 323 ++(*naddedcuts); 324 } 325 } 326 } 327 /* release the row */ 328 SCIP_CALL( SCIPreleaseRow(scip, &cut) ); 329 } 330 331 return SCIP_OKAY; 332 } 333 334 /* 335 * Callback methods 336 */ 337 338 /** copy method for separator plugins (called when SCIP copies plugins) */ 339 static 340 SCIP_DECL_SEPACOPY(sepaCopyGomory) 341 { /*lint --e{715}*/ 342 assert(scip != NULL); 343 assert(sepa != NULL); 344 assert(strcmp(SCIPsepaGetName(sepa), SEPA_NAME) == 0); 345 346 /* call inclusion method of separator */ 347 SCIP_CALL( SCIPincludeSepaGomory(scip) ); 348 349 return SCIP_OKAY; 350 } 351 352 /** destructor of separator to free user data (called when SCIP is exiting) */ 353 /**! [SnippetSepaFreeGomory] */ 354 static 355 SCIP_DECL_SEPAFREE(sepaFreeGomory) 356 { /*lint --e{715}*/ 357 SCIP_SEPADATA* sepadata; 358 359 assert(strcmp(SCIPsepaGetName(sepa), SEPA_NAME) == 0); 360 361 /* free separator data */ 362 sepadata = SCIPsepaGetData(sepa); 363 assert(sepadata != NULL); 364 365 SCIPfreeBlockMemory(scip, &sepadata); 366 367 SCIPsepaSetData(sepa, NULL); 368 369 return SCIP_OKAY; 370 } 371 /**! [SnippetSepaFreeGomory] */ 372 373 /** initialization method of separator (called after problem was transformed) */ 374 static 375 SCIP_DECL_SEPAINIT(sepaInitGomory) 376 { 377 SCIP_SEPADATA* sepadata; 378 379 sepadata = SCIPsepaGetData(sepa); 380 assert(sepadata != NULL); 381 382 /* create and initialize random number generator */ 383 SCIP_CALL( SCIPcreateRandom(scip, &sepadata->randnumgen, DEFAULT_RANDSEED, TRUE) ); 384 385 return SCIP_OKAY; 386 } 387 388 /** deinitialization method of separator (called before transformed problem is freed) */ 389 static 390 SCIP_DECL_SEPAEXIT(sepaExitGomory) 391 { /*lint --e{715}*/ 392 SCIP_SEPADATA* sepadata; 393 394 sepadata = SCIPsepaGetData(sepa); 395 assert(sepadata != NULL); 396 397 SCIPfreeRandom(scip, &sepadata->randnumgen); 398 399 return SCIP_OKAY; 400 } 401 402 403 /** LP solution separation method of separator */ 404 static 405 SCIP_DECL_SEPAEXECLP(sepaExeclpGomory) 406 { /*lint --e{715}*/ 407 SCIP_SEPADATA* sepadata; 408 SCIP_VAR** vars; 409 SCIP_COL** cols; 410 SCIP_ROW** rows; 411 SCIP_AGGRROW* aggrrow; 412 SCIP_VAR* var; 413 SCIP_Real* binvrow; 414 SCIP_Real* cutcoefs; 415 SCIP_Real* basisfrac; 416 SCIP_Real* cutefficacies; 417 int* basisind; 418 int* basisperm; 419 int* inds; 420 int* cutinds; 421 int* colindsproducedcut; 422 SCIP_Real maxscale; 423 SCIP_Real minfrac; 424 SCIP_Real maxfrac; 425 SCIP_Real maxcutefficacy; 426 SCIP_Longint maxdnom; 427 SCIP_Bool cutoff; 428 SCIP_Bool separatescg; 429 SCIP_Bool separategmi; 430 int naddedcuts; 431 int nvars; 432 int ncols; 433 int nrows; 434 int ncalls; 435 int maxdepth; 436 int maxsepacuts; 437 int freq; 438 int c; 439 int i; 440 int j; 441 442 assert(sepa != NULL); 443 assert(strcmp(SCIPsepaGetName(sepa), SEPA_NAME) == 0); 444 assert(scip != NULL); 445 assert(result != NULL); 446 447 *result = SCIP_DIDNOTRUN; 448 449 sepadata = SCIPsepaGetData(sepa); 450 assert(sepadata != NULL); 451 452 ncalls = SCIPsepaGetNCallsAtNode(sepa); 453 454 minfrac = sepadata->away; 455 maxfrac = 1.0 - sepadata->away; 456 457 /* only call separator, if we are not close to terminating */ 458 if( SCIPisStopped(scip) ) 459 return SCIP_OKAY; 460 461 /* only call the gomory cut separator a given number of times at each node */ 462 if( (depth == 0 && sepadata->maxroundsroot >= 0 && ncalls >= sepadata->maxroundsroot) 463 || (depth > 0 && sepadata->maxrounds >= 0 && ncalls >= sepadata->maxrounds) ) 464 return SCIP_OKAY; 465 466 /* only call separator, if an optimal LP solution is at hand */ 467 if( SCIPgetLPSolstat(scip) != SCIP_LPSOLSTAT_OPTIMAL ) 468 return SCIP_OKAY; 469 470 /* only call separator, if the LP solution is basic */ 471 if( !SCIPisLPSolBasic(scip) ) 472 return SCIP_OKAY; 473 474 /* only call separator, if there are fractional variables */ 475 if( SCIPgetNLPBranchCands(scip) == 0 ) 476 return SCIP_OKAY; 477 478 /* check whether strong CG cuts should be separated */ 479 freq = SCIPsepaGetFreq(sepadata->strongcg); 480 if( freq > 0 ) 481 separatescg = (depth % freq == 0); 482 else 483 separatescg = (freq == depth); 484 485 /* check whether Gomory MI cuts should be separated */ 486 freq = SCIPsepaGetFreq(sepadata->gomory); 487 if( freq > 0 ) 488 separategmi = (depth % freq == 0); 489 else 490 separategmi = (freq == depth); 491 492 if( !separatescg && !separategmi ) 493 return SCIP_OKAY; 494 495 /* get variables data */ 496 SCIP_CALL( SCIPgetVarsData(scip, &vars, &nvars, NULL, NULL, NULL, NULL) ); 497 498 /* get LP data */ 499 SCIP_CALL( SCIPgetLPColsData(scip, &cols, &ncols) ); 500 SCIP_CALL( SCIPgetLPRowsData(scip, &rows, &nrows) ); 501 if( ncols == 0 || nrows == 0 ) 502 return SCIP_OKAY; 503 504 /* set the maximal denominator in rational representation of gomory cut and the maximal scale factor to 505 * scale resulting cut to integral values to avoid numerical instabilities 506 */ 507 /**@todo find better but still stable gomory cut settings: look at dcmulti, gesa3, khb0525, misc06, p2756 */ 508 maxdepth = SCIPgetMaxDepth(scip); 509 if( depth == 0 ) 510 { 511 maxdnom = 1000; 512 maxscale = 1000.0; 513 } 514 else if( depth <= maxdepth/4 ) 515 { 516 maxdnom = 1000; 517 maxscale = 1000.0; 518 } 519 else if( depth <= maxdepth/2 ) 520 { 521 maxdnom = 100; 522 maxscale = 100.0; 523 } 524 else 525 { 526 maxdnom = 10; 527 maxscale = 10.0; 528 } 529 530 /* allocate temporary memory */ 531 SCIP_CALL( SCIPallocBufferArray(scip, &cutcoefs, nvars) ); 532 SCIP_CALL( SCIPallocBufferArray(scip, &cutinds, nvars) ); 533 SCIP_CALL( SCIPallocBufferArray(scip, &basisind, nrows) ); 534 SCIP_CALL( SCIPallocBufferArray(scip, &basisperm, nrows) ); 535 SCIP_CALL( SCIPallocBufferArray(scip, &basisfrac, nrows) ); 536 SCIP_CALL( SCIPallocBufferArray(scip, &binvrow, nrows) ); 537 SCIP_CALL( SCIPallocBufferArray(scip, &inds, nrows) ); 538 SCIP_CALL( SCIPallocBufferArray(scip, &cutefficacies, nrows) ); 539 SCIP_CALL( SCIPallocBufferArray(scip, &colindsproducedcut, nrows) ); 540 SCIP_CALL( SCIPaggrRowCreate(scip, &aggrrow) ); 541 542 /* get basis indices */ 543 SCIP_CALL( SCIPgetLPBasisInd(scip, basisind) ); 544 545 for( i = 0; i < nrows; ++i ) 546 { 547 SCIP_Real frac = 0.0; 548 549 c = basisind[i]; 550 cutefficacies[i] = 0.0; 551 552 basisperm[i] = i; 553 554 colindsproducedcut[i] = -1; 555 556 if( c >= 0 ) 557 { 558 559 assert(c < ncols); 560 var = SCIPcolGetVar(cols[c]); 561 if( SCIPvarGetType(var) != SCIP_VARTYPE_CONTINUOUS ) 562 { 563 frac = SCIPfeasFrac(scip, SCIPcolGetPrimsol(cols[c])); 564 frac = MIN(frac, 1.0 - frac); 565 } 566 } 567 else if( sepadata->separaterows ) 568 { 569 SCIP_ROW* row; 570 571 assert(0 <= -c-1 && -c-1 < nrows); 572 row = rows[-c-1]; 573 if( SCIProwIsIntegral(row) && !SCIProwIsModifiable(row) ) 574 { 575 frac = SCIPfeasFrac(scip, SCIPgetRowActivity(scip, row)); 576 frac = MIN(frac, 1.0 - frac); 577 } 578 } 579 580 if( frac >= minfrac ) 581 { 582 /* slightly change fractionality to have random order for equal fractions */ 583 basisfrac[i] = frac + SCIPrandomGetReal(sepadata->randnumgen, -1e-6, 1e-6); 584 } 585 else 586 { 587 basisfrac[i] = 0.0; 588 } 589 } 590 591 /* sort basis indices by fractionality */ 592 SCIPsortDownRealInt(basisfrac, basisperm, nrows); 593 594 /* get the maximal number of cuts allowed in a separation round */ 595 if( depth == 0 ) 596 maxsepacuts = sepadata->maxsepacutsroot; 597 else 598 maxsepacuts = sepadata->maxsepacuts; 599 600 SCIPdebugMsg(scip, "searching gomory cuts: %d cols, %d rows, maxdnom=%" SCIP_LONGINT_FORMAT ", maxscale=%g, maxcuts=%d\n", 601 ncols, nrows, maxdnom, maxscale, maxsepacuts); 602 603 cutoff = FALSE; 604 naddedcuts = 0; 605 606 /* for all basic columns belonging to integer variables, try to generate a gomory cut */ 607 for( i = 0; i < nrows && naddedcuts < maxsepacuts && !SCIPisStopped(scip) && !cutoff; ++i ) 608 { 609 SCIP_Real cutrhs; 610 SCIP_Real cutefficacy = 0.0; 611 SCIP_Bool success; 612 SCIP_Bool cutislocal; 613 SCIP_Bool strongcgsuccess = FALSE; 614 int ninds = -1; 615 int cutnnz; 616 int cutrank; 617 618 if( basisfrac[i] == 0.0 ) 619 break; 620 621 j = basisperm[i]; 622 c = basisind[j]; 623 624 /* get the row of B^-1 for this basic integer variable with fractional solution value */ 625 SCIP_CALL( SCIPgetLPBInvRow(scip, j, binvrow, inds, &ninds) ); 626 627 SCIP_CALL( SCIPaggrRowSumRows(scip, aggrrow, binvrow, inds, ninds, 628 sepadata->sidetypebasis, allowlocal, 2, (int) MAXAGGRLEN(nvars), &success) ); 629 630 if( !success ) 631 continue; 632 633 /* try to create a strong CG cut out of the aggregation row */ 634 if( separatescg ) 635 { 636 SCIP_CALL( SCIPcalcStrongCG(scip, NULL, POSTPROCESS, BOUNDSWITCH, USEVBDS, allowlocal, minfrac, maxfrac, 637 1.0, aggrrow, cutcoefs, &cutrhs, cutinds, &cutnnz, &cutefficacy, &cutrank, &cutislocal, &strongcgsuccess) ); 638 639 /* if we want to generate both cuts, add cut and reset cutefficacy and strongcgsuccess */ 640 if( strongcgsuccess && sepadata->genbothgomscg ) 641 { 642 assert(allowlocal || !cutislocal); /*lint !e644*/ 643 SCIP_CALL( addCut(scip, sepadata, vars, c, maxdnom, maxscale, cutnnz, cutinds, cutcoefs, cutefficacy, cutrhs, 644 cutislocal, cutrank, TRUE, &cutoff, &naddedcuts) ); 645 if( c >= 0 ) 646 { 647 cutefficacies[i] = cutefficacy; 648 colindsproducedcut[i] = c; 649 } 650 cutefficacy = 0.0; 651 strongcgsuccess = FALSE; 652 if( cutoff ) 653 break; 654 } 655 } 656 657 /* @todo Currently we are using the SCIPcalcMIR() function to compute the coefficients of the Gomory 658 * cut. Alternatively, we could use the direct version (see thesis of Achterberg formula (8.4)) which 659 * leads to cut a of the form \sum a_i x_i \geq 1. Rumor has it that these cuts are better. 660 */ 661 662 /* try to create Gomory cut out of the aggregation row */ 663 if( separategmi ) 664 { 665 /* SCIPcalcMIR will only override the cut if its efficacy is larger than the one of the strongcg cut */ 666 SCIP_CALL( SCIPcalcMIR(scip, NULL, POSTPROCESS, BOUNDSWITCH, USEVBDS, allowlocal, FIXINTEGRALRHS, NULL, NULL, 667 minfrac, maxfrac, 1.0, aggrrow, cutcoefs, &cutrhs, cutinds, &cutnnz, &cutefficacy, &cutrank, &cutislocal, &success) ); 668 669 if( success || strongcgsuccess ) 670 { 671 assert(allowlocal || !cutislocal); /*lint !e644*/ 672 if( success ) 673 strongcgsuccess = FALSE; /* Set strongcgsuccess to FALSE, since the MIR cut has overriden the strongcg cut. */ 674 675 SCIP_CALL( addCut(scip, sepadata, vars, c, maxdnom, maxscale, cutnnz, cutinds, cutcoefs, cutefficacy, cutrhs, 676 cutislocal, cutrank, strongcgsuccess, &cutoff, &naddedcuts) ); 677 if( c >= 0 ) 678 { 679 cutefficacies[i] = cutefficacy; 680 colindsproducedcut[i] = c; 681 } 682 } 683 } 684 } 685 686 /* Add normalized efficacy GMI statistics to history */ 687 maxcutefficacy = 0.0; 688 for( i = 0; i < nrows; ++i ) 689 { 690 if( cutefficacies[i] > maxcutefficacy && colindsproducedcut[i] >= 0 ) 691 { 692 maxcutefficacy = cutefficacies[i]; 693 } 694 } 695 696 697 for( i = 0; i < nrows; ++i ) 698 { 699 if( colindsproducedcut[i] >= 0 && SCIPisEfficacious(scip, cutefficacies[i]) ) 700 { 701 assert( maxcutefficacy > 0.0 ); 702 var = SCIPcolGetVar(cols[colindsproducedcut[i]]); 703 SCIP_CALL( SCIPsetVarLastGMIScore(scip, var, cutefficacies[i] / maxcutefficacy) ); 704 SCIP_CALL( SCIPincVarGMISumScore(scip, var, cutefficacies[i] / maxcutefficacy) ); 705 } 706 } 707 708 /* free temporary memory */ 709 SCIPfreeBufferArray(scip, &inds); 710 SCIPfreeBufferArray(scip, &binvrow); 711 SCIPfreeBufferArray(scip, &basisfrac); 712 SCIPfreeBufferArray(scip, &basisperm); 713 SCIPfreeBufferArray(scip, &basisind); 714 SCIPfreeBufferArray(scip, &cutinds); 715 SCIPfreeBufferArray(scip, &cutcoefs); 716 SCIPfreeBufferArray(scip, &cutefficacies); 717 SCIPfreeBufferArray(scip, &colindsproducedcut); 718 SCIPaggrRowFree(scip, &aggrrow); 719 720 SCIPdebugMsg(scip, "end searching gomory cuts: found %d cuts\n", naddedcuts); 721 722 sepadata->lastncutsfound = SCIPgetNCutsFound(scip); 723 724 /* evaluate the result of the separation */ 725 if( cutoff ) 726 *result = SCIP_CUTOFF; 727 else if ( naddedcuts > 0 ) 728 *result = SCIP_SEPARATED; 729 else 730 *result = SCIP_DIDNOTFIND; 731 732 return SCIP_OKAY; 733 } 734 735 736 /* 737 * separator specific interface methods 738 */ 739 740 /** LP solution separation method of dummy separator */ 741 static 742 SCIP_DECL_SEPAEXECLP(sepaExeclpDummy) 743 { /*lint --e{715}*/ 744 assert( result != NULL ); 745 746 *result = SCIP_DIDNOTRUN; 747 748 return SCIP_OKAY; 749 } 750 751 /** arbitrary primal solution separation method of dummy separator */ 752 static 753 SCIP_DECL_SEPAEXECSOL(sepaExecsolDummy) 754 { /*lint --e{715}*/ 755 assert( result != NULL ); 756 757 *result = SCIP_DIDNOTRUN; 758 759 return SCIP_OKAY; 760 } 761 762 /** creates the Gomory MIR cut separator and includes it in SCIP */ 763 SCIP_RETCODE SCIPincludeSepaGomory( 764 SCIP* scip /**< SCIP data structure */ 765 ) 766 { 767 SCIP_SEPADATA* sepadata; 768 SCIP_SEPA* sepa; 769 770 /* create separator data */ 771 SCIP_CALL( SCIPallocBlockMemory(scip, &sepadata) ); 772 sepadata->lastncutsfound = 0; 773 774 /* include separator */ 775 SCIP_CALL( SCIPincludeSepaBasic(scip, &sepa, SEPA_NAME, SEPA_DESC, SEPA_PRIORITY, SEPA_FREQ, SEPA_MAXBOUNDDIST, 776 SEPA_USESSUBSCIP, SEPA_DELAY, 777 sepaExeclpGomory, NULL, 778 sepadata) ); 779 assert(sepa != NULL); 780 781 SCIP_CALL( SCIPincludeSepaBasic(scip, &sepadata->strongcg, "strongcg", "separator for strong CG cuts", -100000, SEPA_FREQ, 0.0, 782 SEPA_USESSUBSCIP, FALSE, sepaExeclpDummy, sepaExecsolDummy, NULL) ); 783 assert(sepadata->strongcg != NULL); 784 785 SCIP_CALL( SCIPincludeSepaBasic(scip, &sepadata->gomory, "gomorymi", "separator for Gomory mixed-integer cuts", -100000, SEPA_FREQ, 0.0, 786 SEPA_USESSUBSCIP, FALSE, sepaExeclpDummy, sepaExecsolDummy, NULL) ); 787 assert(sepadata->gomory != NULL); 788 789 /* set non-NULL pointers to callback methods */ 790 SCIP_CALL( SCIPsetSepaCopy(scip, sepa, sepaCopyGomory) ); 791 SCIP_CALL( SCIPsetSepaFree(scip, sepa, sepaFreeGomory) ); 792 SCIP_CALL( SCIPsetSepaInit(scip, sepa, sepaInitGomory) ); 793 SCIP_CALL( SCIPsetSepaExit(scip, sepa, sepaExitGomory) ); 794 795 /* mark main separator as a parent */ 796 SCIPsetSepaIsParentsepa(scip, sepa); 797 798 /* set pointer from child separators to main separator */ 799 SCIPsetSepaParentsepa(scip, sepadata->strongcg, sepa); 800 SCIPsetSepaParentsepa(scip, sepadata->gomory, sepa); 801 802 /* add separator parameters */ 803 SCIP_CALL( SCIPaddIntParam(scip, 804 "separating/gomory/maxrounds", 805 "maximal number of gomory separation rounds per node (-1: unlimited)", 806 &sepadata->maxrounds, FALSE, DEFAULT_MAXROUNDS, -1, INT_MAX, NULL, NULL) ); 807 SCIP_CALL( SCIPaddIntParam(scip, 808 "separating/gomory/maxroundsroot", 809 "maximal number of gomory separation rounds in the root node (-1: unlimited)", 810 &sepadata->maxroundsroot, FALSE, DEFAULT_MAXROUNDSROOT, -1, INT_MAX, NULL, NULL) ); 811 SCIP_CALL( SCIPaddIntParam(scip, 812 "separating/gomory/maxsepacuts", 813 "maximal number of gomory cuts separated per separation round", 814 &sepadata->maxsepacuts, FALSE, DEFAULT_MAXSEPACUTS, 0, INT_MAX, NULL, NULL) ); 815 SCIP_CALL( SCIPaddIntParam(scip, 816 "separating/gomory/maxsepacutsroot", 817 "maximal number of gomory cuts separated per separation round in the root node", 818 &sepadata->maxsepacutsroot, FALSE, DEFAULT_MAXSEPACUTSROOT, 0, INT_MAX, NULL, NULL) ); 819 SCIP_CALL( SCIPaddIntParam(scip, 820 "separating/gomory/maxrank", 821 "maximal rank of a gomory cut that could not be scaled to integral coefficients (-1: unlimited)", 822 &sepadata->maxrank, FALSE, DEFAULT_MAXRANK, -1, INT_MAX, NULL, NULL) ); 823 SCIP_CALL( SCIPaddIntParam(scip, 824 "separating/gomory/maxrankintegral", 825 "maximal rank of a gomory cut that could be scaled to integral coefficients (-1: unlimited)", 826 &sepadata->maxrankintegral, FALSE, DEFAULT_MAXRANKINTEGRAL, -1, INT_MAX, NULL, NULL) ); 827 SCIP_CALL( SCIPaddRealParam(scip, 828 "separating/gomory/away", 829 "minimal integrality violation of a basis variable in order to try Gomory cut", 830 &sepadata->away, FALSE, DEFAULT_AWAY, 1e-4, 0.5, NULL, NULL) ); 831 SCIP_CALL( SCIPaddBoolParam(scip, 832 "separating/gomory/dynamiccuts", 833 "should generated cuts be removed from the LP if they are no longer tight?", 834 &sepadata->dynamiccuts, FALSE, DEFAULT_DYNAMICCUTS, NULL, NULL) ); 835 SCIP_CALL( SCIPaddBoolParam(scip, 836 "separating/gomory/makeintegral", 837 "try to scale cuts to integral coefficients", 838 &sepadata->makeintegral, TRUE, DEFAULT_MAKEINTEGRAL, NULL, NULL) ); 839 SCIP_CALL( SCIPaddBoolParam(scip, 840 "separating/gomory/forcecuts", 841 "if conversion to integral coefficients failed still consider the cut", 842 &sepadata->forcecuts, TRUE, DEFAULT_FORCECUTS, NULL, NULL) ); 843 SCIP_CALL( SCIPaddBoolParam(scip, 844 "separating/gomory/separaterows", 845 "separate rows with integral slack", 846 &sepadata->separaterows, TRUE, DEFAULT_SEPARATEROWS, NULL, NULL) ); 847 SCIP_CALL( SCIPaddBoolParam(scip, 848 "separating/gomory/delayedcuts", 849 "should cuts be added to the delayed cut pool?", 850 &sepadata->delayedcuts, TRUE, DEFAULT_DELAYEDCUTS, NULL, NULL) ); 851 SCIP_CALL( SCIPaddBoolParam(scip, 852 "separating/gomory/sidetypebasis", 853 "choose side types of row (lhs/rhs) based on basis information?", 854 &sepadata->sidetypebasis, TRUE, DEFAULT_SIDETYPEBASIS, NULL, NULL) ); 855 SCIP_CALL( SCIPaddBoolParam(scip, 856 "separating/gomory/trystrongcg", 857 "try to generate strengthened Chvatal-Gomory cuts?", 858 &sepadata->trystrongcg, TRUE, DEFAULT_TRYSTRONGCG, NULL, NULL) ); 859 SCIP_CALL( SCIPaddBoolParam(scip, 860 "separating/gomory/genbothgomscg", 861 "Should both Gomory and strong CG cuts be generated (otherwise take best)?", 862 &sepadata->genbothgomscg, TRUE, DEFAULT_GENBOTHGOMSCG, NULL, NULL) ); 863 864 return SCIP_OKAY; 865 } 866