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 benderscut_feasalt.c 26 * @brief Alternative feasibility cuts for Benders' decomposition 27 * @author Stephen J. Maher 28 */ 29 30 /*---+----1----+----2----+----3----+----4----+----5----+----6----+----7----+----8----+----9----+----0----+----1----+----2*/ 31 32 #include <assert.h> 33 #include <string.h> 34 35 #include "scip/pub_expr.h" 36 #include "scip/benderscut_feasalt.h" 37 #include "scip/benderscut_opt.h" 38 #include "scip/cons_linear.h" 39 #include "scip/pub_benderscut.h" 40 #include "scip/pub_benders.h" 41 #include "scip/pub_lp.h" 42 #include "scip/pub_message.h" 43 #include "scip/pub_misc.h" 44 #include "scip/pub_misc_linear.h" 45 #include "scip/pub_nlp.h" 46 #include "scip/pub_var.h" 47 #include "scip/scip_benders.h" 48 #include "scip/scip_cons.h" 49 #include "scip/scip_general.h" 50 #include "scip/scip_lp.h" 51 #include "scip/scip_mem.h" 52 #include "scip/scip_message.h" 53 #include "scip/scip_nlp.h" 54 #include "scip/scip_nlpi.h" 55 #include "scip/scip_numerics.h" 56 #include "scip/scip_param.h" 57 #include "scip/scip_prob.h" 58 #include "scip/scip_solvingstats.h" 59 #include "scip/scip_timing.h" 60 #include "scip/scip_var.h" 61 62 #define BENDERSCUT_NAME "feasalt" 63 #define BENDERSCUT_DESC "Alternative feasibility cuts for Benders' decomposition" 64 #define BENDERSCUT_PRIORITY 10001 65 #define BENDERSCUT_LPCUT TRUE 66 67 #define SCIP_DEFAULT_DISPLAYFREQ 20 68 #define SLACKVAR_NAME "##bendersslackvar" /** the name for the Benders' slack variables added to each 69 * constraints in the subproblems */ 70 71 struct SCIP_BenderscutData 72 { 73 SCIP_NLPI* nlpi; /**< nlpi used to create the nlpi problem */ 74 SCIP_NLPIPROBLEM* nlpiprob; /**< nlpi problem representing the convex NLP relaxation */ 75 SCIP_HASHMAP* var2idx; /**< mapping the variable to the index in the NLPI problem */ 76 SCIP_HASHMAP* row2idx; /**< mapping the rows to the index in the NLPI problem */ 77 SCIP_VAR** nlpivars; /**< the variables in the NLPI problem */ 78 SCIP_NLROW** nlpirows; /**< the rows in the NLPI problem */ 79 int nlpinvars; /**< the number of variables in the NPLI problem */ 80 int nlpinrows; /**< the number of rows in the NLPI problem */ 81 int nlpinslackvars; /**< the number of slack variables in the NLPI problem */ 82 int nlpiprobsubprob; /**< the index of the subproblem that the nonlinear problem belongs to */ 83 84 SCIP_Real* slackvarlbs; /**< an array of zeros for the slack variable lower bounds*/ 85 SCIP_Real* slackvarubs; /**< an array of infinity for the slack variable upper bounds*/ 86 int* slackvarinds; /**< array of indices for the slack variables */ 87 }; 88 89 /* 90 * Local methods 91 */ 92 93 /** frees the non linear problem */ 94 static 95 SCIP_RETCODE freeNonlinearProblem( 96 SCIP* masterprob, /**< the SCIP instance of the master problem */ 97 SCIP* subproblem, /**< the SCIP instance of the pricing problem */ 98 SCIP_BENDERSCUT* benderscut /**< the Benders' decomposition structure */ 99 ) 100 { 101 SCIP_BENDERSCUTDATA* benderscutdata; 102 103 assert(masterprob != NULL); 104 assert(subproblem != NULL); 105 assert(benderscut != NULL); 106 107 benderscutdata = SCIPbenderscutGetData(benderscut); 108 assert(benderscutdata != NULL); 109 110 if( benderscutdata->nlpiprob != NULL ) 111 { 112 assert(benderscutdata->nlpi != NULL); 113 114 SCIPfreeBlockMemoryArray(masterprob, &benderscutdata->slackvarinds, benderscutdata->nlpinvars); 115 SCIPfreeBlockMemoryArray(masterprob, &benderscutdata->slackvarubs, benderscutdata->nlpinvars); 116 SCIPfreeBlockMemoryArray(masterprob, &benderscutdata->slackvarlbs, benderscutdata->nlpinvars); 117 SCIPfreeBlockMemoryArray(masterprob, &benderscutdata->nlpirows, benderscutdata->nlpinrows); 118 SCIPfreeBlockMemoryArray(masterprob, &benderscutdata->nlpivars, benderscutdata->nlpinvars); 119 SCIPhashmapFree(&benderscutdata->row2idx); 120 SCIPhashmapFree(&benderscutdata->var2idx); 121 122 SCIP_CALL( SCIPfreeNlpiProblem(subproblem, benderscutdata->nlpi, &benderscutdata->nlpiprob) ); 123 124 benderscutdata->nlpinslackvars = 0; 125 benderscutdata->nlpinrows = 0; 126 benderscutdata->nlpinvars = 0; 127 128 benderscutdata->nlpi = NULL; 129 } 130 131 return SCIP_OKAY; 132 } 133 134 /** solves the auxiliary feasibility subproblem. 135 * 136 * @note: the variable fixings need to be setup before calling this function 137 */ 138 static 139 SCIP_RETCODE solveFeasibilityNonlinearSubproblem( 140 SCIP* scip, /**< SCIP data structure */ 141 SCIP_BENDERSCUTDATA* benderscutdata, /**< Benders' cut data */ 142 SCIP_Bool* success /**< returns whether solving the feasibility problem was successful */ 143 ) 144 { 145 SCIP_NLPSOLSTAT nlpsolstat; 146 147 assert(scip != NULL); 148 assert(benderscutdata != NULL); 149 150 (*success) = TRUE; 151 152 SCIP_CALL( SCIPsolveNlpi(scip, benderscutdata->nlpi, benderscutdata->nlpiprob, .iterlimit = 3000) ); /*lint !e666*/ 153 SCIPdebugMsg(scip, "NLP solstat = %d\n", SCIPgetNlpiSolstat(scip, benderscutdata->nlpi, benderscutdata->nlpiprob)); 154 155 nlpsolstat = SCIPgetNlpiSolstat(scip, benderscutdata->nlpi, benderscutdata->nlpiprob); 156 157 /* if the feasibility NLP is not feasible, then it is not possible to generate a Benders' cut. This is also an error, 158 * since the NLP should always be feasible. In debug mode, an ABORT will be thrown. 159 */ 160 if( nlpsolstat > SCIP_NLPSOLSTAT_FEASIBLE ) 161 (*success) = FALSE; 162 163 return SCIP_OKAY; 164 } 165 166 /** builds the non-linear problem to resolve to generate a cut for the infeasible subproblem */ 167 static 168 SCIP_RETCODE createAuxiliaryNonlinearSubproblem( 169 SCIP* masterprob, /**< the SCIP instance of the master problem */ 170 SCIP* subproblem, /**< the SCIP instance of the pricing problem */ 171 SCIP_BENDERSCUT* benderscut /**< the benders' decomposition cut method */ 172 ) 173 { 174 SCIP_BENDERSCUTDATA* benderscutdata; 175 SCIP_Real* obj; 176 int i; 177 178 assert(masterprob != NULL); 179 180 benderscutdata = SCIPbenderscutGetData(benderscut); 181 assert(benderscutdata != NULL); 182 183 /* first freeing the non-linear problem if it exists */ 184 SCIP_CALL( freeNonlinearProblem(masterprob, subproblem, benderscut) ); 185 186 assert(benderscutdata->nlpi == NULL); 187 assert(benderscutdata->nlpiprob == NULL); 188 189 benderscutdata->nlpinvars = SCIPgetNVars(subproblem); 190 benderscutdata->nlpinrows = SCIPgetNNLPNlRows(subproblem); 191 benderscutdata->nlpi = SCIPgetNlpis(subproblem)[0]; 192 assert(benderscutdata->nlpi != NULL); 193 194 SCIP_CALL( SCIPhashmapCreate(&benderscutdata->var2idx, SCIPblkmem(masterprob), benderscutdata->nlpinvars) ); 195 SCIP_CALL( SCIPhashmapCreate(&benderscutdata->row2idx, SCIPblkmem(masterprob), benderscutdata->nlpinrows) ); 196 197 SCIP_CALL( SCIPduplicateBlockMemoryArray(masterprob, &benderscutdata->nlpivars, SCIPgetVars(subproblem), 198 benderscutdata->nlpinvars) ); /*lint !e666*/ 199 SCIP_CALL( SCIPduplicateBlockMemoryArray(masterprob, &benderscutdata->nlpirows, SCIPgetNLPNlRows(subproblem), 200 benderscutdata->nlpinrows) ); /*lint !e666*/ 201 202 SCIP_CALL( SCIPcreateNlpiProblemFromNlRows(subproblem, benderscutdata->nlpi, &benderscutdata->nlpiprob, "benders-feascutalt-nlp", 203 SCIPgetNLPNlRows(subproblem), benderscutdata->nlpinrows, benderscutdata->var2idx, benderscutdata->row2idx, NULL, SCIPinfinity(subproblem), FALSE, 204 FALSE) ); 205 206 /* storing the slack variable bounds and indices */ 207 SCIP_CALL( SCIPallocBufferArray(masterprob, &obj, benderscutdata->nlpinvars) ); 208 209 SCIP_CALL( SCIPallocBlockMemoryArray(masterprob, &benderscutdata->slackvarlbs, benderscutdata->nlpinvars) ); 210 SCIP_CALL( SCIPallocBlockMemoryArray(masterprob, &benderscutdata->slackvarubs, benderscutdata->nlpinvars) ); 211 SCIP_CALL( SCIPallocBlockMemoryArray(masterprob, &benderscutdata->slackvarinds, benderscutdata->nlpinvars) ); 212 benderscutdata->nlpinslackvars = 0; 213 for( i = 0; i < benderscutdata->nlpinvars; i++ ) 214 { 215 if( strstr(SCIPvarGetName(benderscutdata->nlpivars[i]), SLACKVAR_NAME) ) 216 { 217 benderscutdata->slackvarlbs[benderscutdata->nlpinslackvars] = 0.0; 218 benderscutdata->slackvarubs[benderscutdata->nlpinslackvars] = SCIPinfinity(subproblem); 219 benderscutdata->slackvarinds[benderscutdata->nlpinslackvars] = SCIPhashmapGetImageInt(benderscutdata->var2idx, 220 (void*)benderscutdata->nlpivars[i]); 221 222 obj[benderscutdata->nlpinslackvars] = 1.0; 223 224 benderscutdata->nlpinslackvars++; 225 } 226 } 227 228 /* setting the objective function */ 229 SCIP_CALL( SCIPsetNlpiObjective(subproblem, benderscutdata->nlpi, benderscutdata->nlpiprob, benderscutdata->nlpinslackvars, 230 benderscutdata->slackvarinds, obj, NULL, 0.0) ); 231 232 /* unfixing the slack variables */ 233 SCIP_CALL( SCIPchgNlpiVarBounds(subproblem, benderscutdata->nlpi, benderscutdata->nlpiprob, benderscutdata->nlpinslackvars, 234 benderscutdata->slackvarinds, benderscutdata->slackvarlbs, benderscutdata->slackvarubs) ); 235 236 SCIPfreeBufferArray(masterprob, &obj); 237 238 return SCIP_OKAY; 239 } 240 241 /** updates the non-linear problem that is resolved to generate a cut for the infeasible subproblem */ 242 static 243 SCIP_RETCODE updateAuxiliaryNonlinearSubproblem( 244 SCIP* subproblem, /**< the SCIP instance of the pricing problem */ 245 SCIP_BENDERSCUT* benderscut /**< the benders' decomposition cut method */ 246 ) 247 { 248 SCIP_BENDERSCUTDATA* benderscutdata; 249 250 assert(subproblem != NULL); 251 assert(benderscut != NULL); 252 253 benderscutdata = SCIPbenderscutGetData(benderscut); 254 assert(benderscutdata != NULL); 255 assert(benderscutdata->nlpi != NULL); 256 assert(benderscutdata->nlpiprob != NULL); 257 assert(benderscutdata->var2idx != NULL); 258 assert(benderscutdata->row2idx != NULL); 259 260 /* setting the variable bounds to that from the current subproblem */ 261 SCIP_CALL( SCIPupdateNlpiProblem(subproblem, benderscutdata->nlpi, benderscutdata->nlpiprob, benderscutdata->var2idx, 262 benderscutdata->nlpivars, benderscutdata->nlpinvars, SCIPinfinity(subproblem)) ); 263 264 /* unfixing the slack variables */ 265 SCIP_CALL( SCIPchgNlpiVarBounds(subproblem, benderscutdata->nlpi, benderscutdata->nlpiprob, benderscutdata->nlpinslackvars, 266 benderscutdata->slackvarinds, benderscutdata->slackvarlbs, benderscutdata->slackvarubs) ); 267 268 return SCIP_OKAY; 269 } 270 271 /** generates and applies Benders' cuts */ 272 static 273 SCIP_RETCODE generateAndApplyBendersCuts( 274 SCIP* masterprob, /**< the SCIP instance of the master problem */ 275 SCIP* subproblem, /**< the SCIP instance of the pricing problem */ 276 SCIP_BENDERS* benders, /**< the benders' decomposition */ 277 SCIP_BENDERSCUT* benderscut, /**< the benders' decomposition cut method */ 278 SCIP_SOL* sol, /**< primal CIP solution */ 279 int probnumber, /**< the number of the pricing problem */ 280 SCIP_BENDERSENFOTYPE type, /**< the enforcement type calling this function */ 281 SCIP_RESULT* result /**< the result from solving the subproblems */ 282 ) 283 { 284 SCIP_BENDERSCUTDATA* benderscutdata; 285 SCIP_Real* primalvals; 286 SCIP_Real* consdualvals; 287 SCIP_Real* varlbdualvals; 288 SCIP_Real* varubdualvals; 289 SCIP_Real obj; 290 char cutname[SCIP_MAXSTRLEN]; 291 SCIP_Bool success; 292 #ifdef SCIP_EVENMOREDEBUG 293 int i; 294 #endif 295 296 assert(masterprob != NULL); 297 assert(subproblem != NULL); 298 assert(benders != NULL); 299 assert(result != NULL); 300 301 benderscutdata = SCIPbenderscutGetData(benderscut); 302 assert(benderscutdata != NULL); 303 304 /* creating or updating the NLPI problem */ 305 if( benderscutdata->nlpiprob == NULL || benderscutdata->nlpiprobsubprob != probnumber ) 306 { 307 SCIP_CALL( createAuxiliaryNonlinearSubproblem(masterprob, subproblem, benderscut) ); 308 benderscutdata->nlpiprobsubprob = probnumber; 309 } 310 else 311 { 312 SCIP_CALL( updateAuxiliaryNonlinearSubproblem(subproblem, benderscut) ); 313 } 314 315 /* solving the NLPI problem to get the minimum infeasible solution */ 316 SCIP_CALL( solveFeasibilityNonlinearSubproblem(subproblem, benderscutdata, &success) ); 317 318 if( !success ) 319 { 320 (*result) = SCIP_DIDNOTFIND; 321 SCIPdebugMsg(masterprob, "Error in generating Benders' feasibility cut for problem %d. " 322 "The feasibility subproblem failed to solve with a feasible solution.\n", probnumber); 323 return SCIP_OKAY; 324 } 325 326 /* getting the solution from the NLPI problem */ 327 SCIP_CALL( SCIPgetNlpiSolution(subproblem, benderscutdata->nlpi, benderscutdata->nlpiprob, &primalvals, &consdualvals, 328 &varlbdualvals, &varubdualvals, &obj) ); 329 330 #ifdef SCIP_EVENMOREDEBUG 331 SCIPdebugMsg(masterprob, "NLP Feasibility problem solution.\n"); 332 SCIPdebugMsg(masterprob, "Objective: %g.\n", obj); 333 for( i = 0; i < benderscutdata->nlpinvars; i++ ) 334 { 335 int varindex; 336 SCIP_Real solval; 337 if( SCIPhashmapExists(benderscutdata->var2idx, benderscutdata->nlpivars[i]) ) 338 { 339 varindex = SCIPhashmapGetImageInt(benderscutdata->var2idx, benderscutdata->nlpivars[i]); 340 solval = primalvals[varindex]; 341 342 if( !SCIPisZero(masterprob, solval) ) 343 { 344 SCIPdebugMsg(masterprob, "%s (obj: %g): %20g\n", SCIPvarGetName(benderscutdata->nlpivars[i]), 345 SCIPvarGetObj(benderscutdata->nlpivars[i]), solval); 346 } 347 } 348 } 349 #endif 350 351 /* setting the name of the generated cut */ 352 (void) SCIPsnprintf(cutname, SCIP_MAXSTRLEN, "altfeasibilitycut_%d_%" SCIP_LONGINT_FORMAT, probnumber, 353 SCIPbenderscutGetNFound(benderscut) ); 354 355 /* generating a Benders' decomposition cut using the classical optimality cut methods */ 356 SCIP_CALL( SCIPgenerateAndApplyBendersOptCut(masterprob, subproblem, benders, benderscut, 357 sol, probnumber, cutname, obj, primalvals, consdualvals, varlbdualvals, varubdualvals, benderscutdata->row2idx, 358 benderscutdata->var2idx, type, FALSE, TRUE, result) ); 359 360 if( (*result) == SCIP_CONSADDED ) 361 { 362 if( SCIPisInfinity(masterprob, -SCIPgetDualbound(masterprob)) 363 && SCIPbenderscutGetNFound(benderscut) % SCIP_DEFAULT_DISPLAYFREQ == 0 ) 364 { 365 if( SCIPgetStage(masterprob) == SCIP_STAGE_SOLVING ) 366 { 367 SCIP_CALL( SCIPprintDisplayLine(masterprob, NULL, SCIP_VERBLEVEL_NORMAL, TRUE) ); 368 SCIPverbMessage(masterprob, SCIP_VERBLEVEL_NORMAL, NULL, 369 "Benders' Decomposition: Master problem LP is infeasible. Added %" SCIP_LONGINT_FORMAT " feasibility cuts.\n", 370 SCIPbenderscutGetNFound(benderscut)); 371 } 372 } 373 SCIPdebugMsg(masterprob, "Constraint <%s> has been added to the master problem.\n", cutname); 374 } 375 376 return SCIP_OKAY; 377 } 378 379 /* 380 * Callback methods of Benders' decomposition cuts 381 */ 382 383 /** deinitialization method of Benders' decomposition cuts (called before transformed problem is freed) */ 384 static 385 SCIP_DECL_BENDERSCUTEXIT(benderscutExitFeasalt) 386 { /*lint --e{715}*/ 387 assert( benderscut != NULL ); 388 assert( strcmp(SCIPbenderscutGetName(benderscut), BENDERSCUT_NAME) == 0 ); 389 390 return SCIP_OKAY; 391 } 392 393 /** destructor of the Benders' decomposition cut to free user data (called when SCIP is exiting) */ 394 static 395 SCIP_DECL_BENDERSCUTFREE(benderscutFreeFeasalt) 396 { /*lint --e{715}*/ 397 SCIP_BENDERSCUTDATA* benderscutdata; 398 399 assert(scip != NULL); 400 assert(benderscut != NULL); 401 402 benderscutdata = SCIPbenderscutGetData(benderscut); 403 assert(benderscutdata != NULL); 404 405 SCIPfreeBlockMemory(scip, &benderscutdata); 406 407 return SCIP_OKAY; 408 } 409 410 /** execution method of Benders' decomposition cuts */ 411 static 412 SCIP_DECL_BENDERSCUTEXEC(benderscutExecFeasalt) 413 { /*lint --e{715}*/ 414 SCIP* subproblem; 415 SCIP_Bool nlprelaxation; 416 417 assert(scip != NULL); 418 assert(benders != NULL); 419 assert(benderscut != NULL); 420 assert(result != NULL); 421 assert(probnumber >= 0 && probnumber < SCIPbendersGetNSubproblems(benders)); 422 423 subproblem = SCIPbendersSubproblem(benders, probnumber); 424 425 /* setting a flag to indicate whether the NLP relaxation should be used to generate cuts */ 426 nlprelaxation = SCIPisNLPConstructed(subproblem) && SCIPgetNNlpis(subproblem) 427 && SCIPbendersGetSubproblemType(benders, probnumber) <= SCIP_BENDERSSUBTYPE_CONVEXDIS; 428 429 /* only generate feasibility cuts if the subproblem LP or NLP is infeasible, 430 * since we use the farkas proof from the LP or the dual solution of the NLP to construct the feasibility cut 431 */ 432 if( SCIPgetStage(subproblem) == SCIP_STAGE_SOLVING && 433 (nlprelaxation && (SCIPgetNLPSolstat(subproblem) == SCIP_NLPSOLSTAT_LOCINFEASIBLE || SCIPgetNLPSolstat(subproblem) == SCIP_NLPSOLSTAT_GLOBINFEASIBLE)) ) 434 { 435 /* generating a cut for a given subproblem */ 436 SCIP_CALL( generateAndApplyBendersCuts(scip, subproblem, benders, benderscut, sol, probnumber, type, result) ); 437 438 /* TODO this was in benderscutExitFeasalt, but freeNonlinearProblem now needs subproblem, which didn't seem to be easily available there */ 439 /* freeing the non-linear problem information */ 440 SCIP_CALL( freeNonlinearProblem(scip, subproblem, benderscut) ); 441 } 442 443 return SCIP_OKAY; 444 } 445 446 447 /* 448 * Benders' decomposition cuts specific interface methods 449 */ 450 451 /** creates the Alternative Feasibility Benders' decomposition cuts and includes it in SCIP */ 452 SCIP_RETCODE SCIPincludeBenderscutFeasalt( 453 SCIP* scip, /**< SCIP data structure */ 454 SCIP_BENDERS* benders /**< Benders' decomposition */ 455 ) 456 { 457 SCIP_BENDERSCUT* benderscut; 458 SCIP_BENDERSCUTDATA* benderscutdata; 459 460 assert(benders != NULL); 461 462 benderscut = NULL; 463 464 SCIP_CALL( SCIPallocBlockMemory(scip, &benderscutdata) ); 465 BMSclearMemory(benderscutdata); 466 benderscutdata->nlpiprobsubprob = -1; 467 468 /* include Benders' decomposition cuts */ 469 SCIP_CALL( SCIPincludeBenderscutBasic(scip, benders, &benderscut, BENDERSCUT_NAME, BENDERSCUT_DESC, 470 BENDERSCUT_PRIORITY, BENDERSCUT_LPCUT, benderscutExecFeasalt, benderscutdata) ); 471 472 /* set non fundamental callbacks via setter functions */ 473 SCIP_CALL( SCIPsetBenderscutFree(scip, benderscut, benderscutFreeFeasalt) ); 474 SCIP_CALL( SCIPsetBenderscutExit(scip, benderscut, benderscutExitFeasalt) ); 475 476 assert(benderscut != NULL); 477 478 return SCIP_OKAY; 479 } 480