Scippy

SCIP

Solving Constraint Integer Programs

benders.c
Go to the documentation of this file.
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 benders.c
26  * @ingroup OTHER_CFILES
27  * @brief methods for Benders' decomposition
28  * @author Stephen J. Maher
29  */
30 
31 /*---+----1----+----2----+----3----+----4----+----5----+----6----+----7----+----8----+----9----+----0----+----1----+----2*/
32 
33 #include <assert.h>
34 #include <string.h>
35 
36 #include "scip/def.h"
37 #include "scip/set.h"
38 #include "scip/clock.h"
39 #include "scip/dcmp.h"
40 #include "scip/paramset.h"
41 #include "scip/lp.h"
42 #include "scip/prob.h"
43 #include "scip/pricestore.h"
44 #include "scip/scip.h"
45 #include "scip/scipdefplugins.h"
46 #include "scip/benders.h"
47 #include "scip/pub_message.h"
48 #include "scip/pub_misc.h"
49 #include "scip/cons_linear.h"
50 #include "scip/cons_nonlinear.h"
51 
52 #include "scip/struct_benders.h"
53 #include "scip/struct_benderscut.h"
54 
55 #include "scip/benderscut.h"
56 
57 /* Defaults for parameters */
58 #define SCIP_DEFAULT_TRANSFERCUTS FALSE /** should Benders' cuts generated in LNS heuristics be transferred to the main SCIP instance? */
59 #define SCIP_DEFAULT_CUTSASCONSS TRUE /** should the transferred cuts be added as constraints? */
60 #define SCIP_DEFAULT_LNSCHECK TRUE /** should the Benders' decomposition be used in LNS heuristics */
61 #define SCIP_DEFAULT_LNSMAXDEPTH -1 /** maximum depth at which the LNS check is performed */
62 #define SCIP_DEFAULT_LNSMAXCALLS 10 /** the maximum number of Benders' decomposition calls in LNS heuristics */
63 #define SCIP_DEFAULT_LNSMAXCALLSROOT 0 /** the maximum number of root node Benders' decomposition calls in LNS heuristics */
64 #define SCIP_DEFAULT_SUBPROBFRAC 1.0 /** fraction of subproblems that are solved in each iteration */
65 #define SCIP_DEFAULT_UPDATEAUXVARBOUND FALSE /** should the auxiliary variable lower bound be updated by solving the subproblem */
66 #define SCIP_DEFAULT_AUXVARSIMPLINT FALSE /** set the auxiliary variables as implint if the subproblem objective is integer */
67 #define SCIP_DEFAULT_CUTCHECK TRUE /** should cuts be generated during the checking of solutions? */
68 #define SCIP_DEFAULT_STRENGTHENMULT 0.5 /** the convex combination multiplier for the cut strengthening */
69 #define SCIP_DEFAULT_NOIMPROVELIMIT 5 /** the maximum number of cut strengthening without improvement */
70 #define SCIP_DEFAULT_STRENGTHENPERTURB 1e-06 /** the amount by which the cut strengthening solution is perturbed */
71 #define SCIP_DEFAULT_STRENGTHENENABLED FALSE /** enable the core point cut strengthening approach */
72 #define SCIP_DEFAULT_STRENGTHENINTPOINT 'r' /** where should the strengthening interior point be sourced from ('l'p relaxation, 'f'irst solution, 'i'ncumbent solution, 'r'elative interior point, vector of 'o'nes, vector of 'z'eros) */
73 #define SCIP_DEFAULT_NUMTHREADS 1 /** the number of parallel threads to use when solving the subproblems */
74 #define SCIP_DEFAULT_EXECFEASPHASE FALSE /** should a feasibility phase be executed during the root node processing */
75 #define SCIP_DEFAULT_SLACKVARCOEF 1e+6 /** the initial objective coefficient of the slack variables in the subproblem */
76 #define SCIP_DEFAULT_MAXSLACKVARCOEF 1e+9 /** the maximal objective coefficient of the slack variables in the subproblem */
77 #define SCIP_DEFAULT_CHECKCONSCONVEXITY TRUE /** should the constraints of the subproblem be checked for convexity? */
78 
79 #define BENDERS_MAXPSEUDOSOLS 5 /** the maximum number of pseudo solutions checked before suggesting
80  * merge candidates */
81 
82 #define BENDERS_ARRAYSIZE 1000 /**< the initial size of the added constraints/cuts arrays */
83 
84 #define AUXILIARYVAR_NAME "##bendersauxiliaryvar" /** the name for the Benders' auxiliary variables in the master problem */
85 #define SLACKVAR_NAME "##bendersslackvar" /** the name for the Benders' slack variables added to each
86  * constraints in the subproblems */
87 #define NLINEARCONSHDLRS 5
88 
89 /* event handler properties */
90 #define NODEFOCUS_EVENTHDLR_NAME "bendersnodefocus"
91 #define NODEFOCUS_EVENTHDLR_DESC "node focus event handler for Benders' decomposition"
92 
93 #define MIPNODEFOCUS_EVENTHDLR_NAME "bendersmipsolvenodefocus"
94 #define MIPNODEFOCUS_EVENTHDLR_DESC "node focus event handler for the MIP solve method for Benders' decomposition"
95 
96 #define UPPERBOUND_EVENTHDLR_NAME "bendersupperbound"
97 #define UPPERBOUND_EVENTHDLR_DESC "found solution event handler to terminate subproblem solve for a given upper bound"
98 
99 #define NODESOLVED_EVENTHDLR_NAME "bendersnodesolved"
100 #define NODESOLVED_EVENTHDLR_DESC "node solved event handler for the Benders' integer cuts"
101 
102 
103 /** event handler data */
104 struct SCIP_EventhdlrData
105 {
106  int filterpos; /**< the event filter entry */
107  int numruns; /**< the number of times that the problem has been solved */
108  SCIP_Real upperbound; /**< an upper bound for the problem */
109  SCIP_Bool solvecip; /**< is the event called from a MIP subproblem solve*/
110 };
111 
112 
113 /* ---------------- Local methods for event handlers ---------------- */
114 
115 /** initialises the members of the eventhandler data */
116 static
118  SCIP* scip, /**< the SCIP data structure */
119  SCIP_EVENTHDLRDATA* eventhdlrdata /**< the event handler data */
120  )
121 {
122  assert(scip != NULL);
123  assert(eventhdlrdata != NULL);
124 
125  eventhdlrdata->filterpos = -1;
126  eventhdlrdata->numruns = 0;
127  eventhdlrdata->upperbound = -SCIPinfinity(scip);
128  eventhdlrdata->solvecip = FALSE;
129 
130  return SCIP_OKAY;
131 }
132 
133 /** initsol method for the event handlers */
134 static
136  SCIP* scip, /**< the SCIP data structure */
137  SCIP_EVENTHDLR* eventhdlr, /**< the event handlers data structure */
138  SCIP_EVENTTYPE eventtype /**< event type mask to select events to catch */
139  )
140 {
141  SCIP_EVENTHDLRDATA* eventhdlrdata;
142 
143  assert(scip != NULL);
144  assert(eventhdlr != NULL);
145 
146  eventhdlrdata = SCIPeventhdlrGetData(eventhdlr);
147 
148  SCIP_CALL(SCIPcatchEvent(scip, eventtype, eventhdlr, NULL, &eventhdlrdata->filterpos));
149 
150  return SCIP_OKAY;
151 }
152 
153 /** the exit sol method for the event handlers */
154 static
156  SCIP* scip, /**< the SCIP data structure */
157  SCIP_EVENTHDLR* eventhdlr, /**< the event handlers data structure */
158  SCIP_EVENTTYPE eventtype /**< event type mask to select events to catch */
159  )
160 {
161  SCIP_EVENTHDLRDATA* eventhdlrdata;
162 
163  assert(scip != NULL);
164  assert(eventhdlr != NULL);
165 
166  eventhdlrdata = SCIPeventhdlrGetData(eventhdlr);
167 
168  if( eventhdlrdata->filterpos >= 0 )
169  {
170  SCIP_CALL(SCIPdropEvent(scip, eventtype, eventhdlr, NULL, eventhdlrdata->filterpos));
171  eventhdlrdata->filterpos = -1;
172  }
173 
174  return SCIP_OKAY;
175 }
176 
177 /** the exit method for the event handlers */
178 static
180  SCIP* scip, /**< the SCIP data structure */
181  SCIP_EVENTHDLR* eventhdlr /**< the event handlers data structure */
182  )
183 {
184  SCIP_EVENTHDLRDATA* eventhdlrdata;
185 
186  assert(scip != NULL);
187  assert(eventhdlr != NULL);
188 
189  eventhdlrdata = SCIPeventhdlrGetData(eventhdlr);
190 
191  /* reinitialise the event handler data */
192  SCIP_CALL( initEventhandlerData(scip, eventhdlrdata) );
193 
194  return SCIP_OKAY;
195 }
196 
197 /** free method for the event handler */
198 static
200  SCIP* scip, /**< the SCIP data structure */
201  SCIP_EVENTHDLR* eventhdlr /**< the event handlers data structure */
202  )
203 {
204  SCIP_EVENTHDLRDATA* eventhdlrdata;
205 
206  assert(scip != NULL);
207  assert(eventhdlr != NULL);
208 
209  eventhdlrdata = SCIPeventhdlrGetData(eventhdlr);
210  assert(eventhdlrdata != NULL);
211 
212  SCIPfreeBlockMemory(scip, &eventhdlrdata);
213 
214  SCIPeventhdlrSetData(eventhdlr, NULL);
215 
216  return SCIP_OKAY;
217 }
218 
219 
220 
221 /* ---------------- Callback methods of node focus event handler ---------------- */
222 
223 /** exec the event handler */
224 static
225 SCIP_DECL_EVENTEXEC(eventExecBendersNodefocus)
226 { /*lint --e{715}*/
227  SCIP_EVENTHDLRDATA* eventhdlrdata;
228 
229  assert(scip != NULL);
230  assert(eventhdlr != NULL);
231  assert(strcmp(SCIPeventhdlrGetName(eventhdlr), NODEFOCUS_EVENTHDLR_NAME) == 0);
232 
233  eventhdlrdata = SCIPeventhdlrGetData(eventhdlr);
234 
235  /* sending an interrupt solve signal to return the control back to the Benders' decomposition plugin.
236  * This will ensure the SCIP stage is SCIP_STAGE_SOLVING, allowing the use of probing mode. */
238 
239  SCIP_CALL(SCIPdropEvent(scip, SCIP_EVENTTYPE_NODEFOCUSED, eventhdlr, NULL, eventhdlrdata->filterpos));
240  eventhdlrdata->filterpos = -1;
241 
242  return SCIP_OKAY;
243 }
244 
245 /** solving process initialization method of event handler (called when branch and bound process is about to begin) */
246 static
247 SCIP_DECL_EVENTINITSOL(eventInitsolBendersNodefocus)
248 {
249  assert(scip != NULL);
250  assert(eventhdlr != NULL);
251  assert(strcmp(SCIPeventhdlrGetName(eventhdlr), NODEFOCUS_EVENTHDLR_NAME) == 0);
252 
254 
255  return SCIP_OKAY;
256 }
257 
258 /** solving process deinitialization method of event handler (called before branch and bound process data is freed) */
259 static
260 SCIP_DECL_EVENTEXITSOL(eventExitsolBendersNodefocus)
261 {
262  assert(scip != NULL);
263  assert(eventhdlr != NULL);
264  assert(strcmp(SCIPeventhdlrGetName(eventhdlr), NODEFOCUS_EVENTHDLR_NAME) == 0);
265 
267 
268  return SCIP_OKAY;
269 }
270 
271 /** deinitialization method of event handler (called before transformed problem is freed) */
272 static
273 SCIP_DECL_EVENTEXIT(eventExitBendersNodefocus)
274 {
275  assert(scip != NULL);
276  assert(eventhdlr != NULL);
277  assert(strcmp(SCIPeventhdlrGetName(eventhdlr), NODEFOCUS_EVENTHDLR_NAME) == 0);
278 
279  SCIP_CALL( exitEventhandler(scip, eventhdlr) );
280 
281  return SCIP_OKAY;
282 }
283 
284 /** deinitialization method of event handler (called before transformed problem is freed) */
285 static
286 SCIP_DECL_EVENTFREE(eventFreeBendersNodefocus)
287 {
288  assert(scip != NULL);
289  assert(eventhdlr != NULL);
290  assert(strcmp(SCIPeventhdlrGetName(eventhdlr), NODEFOCUS_EVENTHDLR_NAME) == 0);
291 
292  SCIP_CALL( freeEventhandler(scip, eventhdlr) );
293 
294  return SCIP_OKAY;
295 }
296 
297 
298 /* ---------------- Callback methods of MIP solve node focus event handler ---------------- */
299 
300 /** exec the event handler */
301 static
302 SCIP_DECL_EVENTEXEC(eventExecBendersMipnodefocus)
303 { /*lint --e{715}*/
304  SCIP_EVENTHDLRDATA* eventhdlrdata;
305 
306  assert(scip != NULL);
307  assert(eventhdlr != NULL);
308  assert(strcmp(SCIPeventhdlrGetName(eventhdlr), MIPNODEFOCUS_EVENTHDLR_NAME) == 0);
309 
310  eventhdlrdata = SCIPeventhdlrGetData(eventhdlr);
311 
312  /* interrupting the solve so that the control is returned back to the Benders' core. */
313  if( eventhdlrdata->numruns == 0 && !eventhdlrdata->solvecip )
314  {
316  }
317 
318  SCIP_CALL(SCIPdropEvent(scip, SCIP_EVENTTYPE_NODEFOCUSED, eventhdlr, NULL, eventhdlrdata->filterpos));
319  eventhdlrdata->filterpos = -1;
320 
321  eventhdlrdata->numruns++;
322 
323  return SCIP_OKAY;
324 }
325 
326 /** solving process initialization method of event handler (called when branch and bound process is about to begin) */
327 static
328 SCIP_DECL_EVENTINITSOL(eventInitsolBendersMipnodefocus)
329 {
330  assert(scip != NULL);
331  assert(eventhdlr != NULL);
332  assert(strcmp(SCIPeventhdlrGetName(eventhdlr), MIPNODEFOCUS_EVENTHDLR_NAME) == 0);
333 
335 
336  return SCIP_OKAY;
337 }
338 
339 /** solving process deinitialization method of event handler (called before branch and bound process data is freed) */
340 static
341 SCIP_DECL_EVENTEXITSOL(eventExitsolBendersMipnodefocus)
342 {
343  assert(scip != NULL);
344  assert(eventhdlr != NULL);
345  assert(strcmp(SCIPeventhdlrGetName(eventhdlr), MIPNODEFOCUS_EVENTHDLR_NAME) == 0);
346 
348 
349  return SCIP_OKAY;
350 }
351 
352 /** deinitialization method of event handler (called before transformed problem is freed) */
353 static
354 SCIP_DECL_EVENTEXIT(eventExitBendersMipnodefocus)
355 {
356  assert(scip != NULL);
357  assert(eventhdlr != NULL);
358  assert(strcmp(SCIPeventhdlrGetName(eventhdlr), MIPNODEFOCUS_EVENTHDLR_NAME) == 0);
359 
360  SCIP_CALL( exitEventhandler(scip, eventhdlr) );
361 
362  return SCIP_OKAY;
363 }
364 
365 /** deinitialization method of event handler (called before transformed problem is freed) */
366 static
367 SCIP_DECL_EVENTFREE(eventFreeBendersMipnodefocus)
368 {
369  assert(scip != NULL);
370  assert(eventhdlr != NULL);
371  assert(strcmp(SCIPeventhdlrGetName(eventhdlr), MIPNODEFOCUS_EVENTHDLR_NAME) == 0);
372 
373  SCIP_CALL( freeEventhandler(scip, eventhdlr) );
374 
375  return SCIP_OKAY;
376 }
377 
378 /* ---------------- Callback methods of solution found event handler ---------------- */
379 
380 /** exec the event handler */
381 static
382 SCIP_DECL_EVENTEXEC(eventExecBendersUpperbound)
383 { /*lint --e{715}*/
384  SCIP_EVENTHDLRDATA* eventhdlrdata;
385  SCIP_SOL* bestsol;
386 
387  assert(scip != NULL);
388  assert(eventhdlr != NULL);
389  assert(strcmp(SCIPeventhdlrGetName(eventhdlr), UPPERBOUND_EVENTHDLR_NAME) == 0);
390 
391  eventhdlrdata = SCIPeventhdlrGetData(eventhdlr);
392  assert(eventhdlrdata != NULL);
393 
394  bestsol = SCIPgetBestSol(scip);
395 
396  if( SCIPisLT(scip, SCIPgetSolOrigObj(scip, bestsol)*(int)SCIPgetObjsense(scip), eventhdlrdata->upperbound) )
397  {
399  }
400 
401  return SCIP_OKAY;
402 }
403 
404 /** solving process initialization method of event handler (called when branch and bound process is about to begin) */
405 static
406 SCIP_DECL_EVENTINITSOL(eventInitsolBendersUpperbound)
407 {
408  assert(scip != NULL);
409  assert(eventhdlr != NULL);
410  assert(strcmp(SCIPeventhdlrGetName(eventhdlr), UPPERBOUND_EVENTHDLR_NAME) == 0);
411 
413 
414  return SCIP_OKAY;
415 }
416 
417 /** solving process deinitialization method of event handler (called before branch and bound process data is freed) */
418 static
419 SCIP_DECL_EVENTEXITSOL(eventExitsolBendersUpperbound)
420 {
421  assert(scip != NULL);
422  assert(eventhdlr != NULL);
423  assert(strcmp(SCIPeventhdlrGetName(eventhdlr), UPPERBOUND_EVENTHDLR_NAME) == 0);
424 
426 
427  return SCIP_OKAY;
428 }
429 
430 /** deinitialization method of event handler (called before transformed problem is freed) */
431 static
432 SCIP_DECL_EVENTEXIT(eventExitBendersUpperbound)
433 {
434  assert(scip != NULL);
435  assert(eventhdlr != NULL);
436  assert(strcmp(SCIPeventhdlrGetName(eventhdlr), UPPERBOUND_EVENTHDLR_NAME) == 0);
437 
438  SCIP_CALL( exitEventhandler(scip, eventhdlr) );
439 
440  return SCIP_OKAY;
441 }
442 
443 /** deinitialization method of event handler (called before transformed problem is freed) */
444 static
445 SCIP_DECL_EVENTFREE(eventFreeBendersUpperbound)
446 {
447  assert(scip != NULL);
448  assert(eventhdlr != NULL);
449  assert(strcmp(SCIPeventhdlrGetName(eventhdlr), UPPERBOUND_EVENTHDLR_NAME) == 0);
450 
451  SCIP_CALL( freeEventhandler(scip, eventhdlr) );
452 
453  return SCIP_OKAY;
454 }
455 
456 /** updates the upper bound in the event handler data */
457 static
459  SCIP_BENDERS* benders, /**< Benders' decomposition */
460  int probnumber, /**< the subproblem number */
461  SCIP_Real upperbound /**< the upper bound value */
462  )
463 {
464  SCIP_EVENTHDLR* eventhdlr;
465  SCIP_EVENTHDLRDATA* eventhdlrdata;
466 
467  assert(benders != NULL);
468  assert(probnumber >= 0 && probnumber < benders->nsubproblems);
469 
470  eventhdlr = SCIPfindEventhdlr(SCIPbendersSubproblem(benders, probnumber), UPPERBOUND_EVENTHDLR_NAME);
471  assert(eventhdlr != NULL);
472 
473  eventhdlrdata = SCIPeventhdlrGetData(eventhdlr);
474  assert(eventhdlrdata != NULL);
475 
476  eventhdlrdata->upperbound = upperbound;
477 
478  return SCIP_OKAY;
479 }
480 
481 /* ---------------- Callback methods of the node solved event handler ---------------- */
482 
483 /** Updates the cut constant of the Benders' cuts data.
484  * This function solves the master problem with only the auxiliary variables in the objective function.
485  */
486 static
488  SCIP* masterprob, /**< the SCIP instance of the master problem */
489  SCIP_BENDERS* benders /**< Benders' decomposition */
490  )
491 {
492  SCIP_VAR** vars;
493  int nvars;
494  int nsubproblems;
495  int i;
496  SCIP_Bool lperror;
497  SCIP_Bool cutoff;
498 
499  assert(masterprob != NULL);
500  assert(benders != NULL);
501 
502  /* don't run in probing or in repropagation */
503  if( SCIPinProbing(masterprob) || SCIPinRepropagation(masterprob) || SCIPinDive(masterprob) )
504  return SCIP_OKAY;
505 
506  nsubproblems = SCIPbendersGetNSubproblems(benders);
507 
508  SCIP_CALL( SCIPstartProbing(masterprob) );
509 
510  /* change the master problem variables to 0 */
511  nvars = SCIPgetNVars(masterprob);
512  vars = SCIPgetVars(masterprob);
513 
514  /* setting the objective function coefficient to 0 for all variables */
515  for( i = 0; i < nvars; i++ )
516  {
517  if( SCIPvarGetStatus(vars[i]) == SCIP_VARSTATUS_COLUMN )
518  {
519  SCIP_CALL( SCIPchgVarObjProbing(masterprob, vars[i], 0.0) );
520  }
521  }
522 
523  /* solving an LP for all subproblems to find the lower bound */
524  for( i = 0; i < nsubproblems; i++)
525  {
526  SCIP_VAR* auxiliaryvar;
527 
528  auxiliaryvar = SCIPbendersGetAuxiliaryVar(benders, i);
529 
530  if( SCIPvarGetStatus(auxiliaryvar) != SCIP_VARSTATUS_COLUMN )
531  continue;
532 
533  SCIP_CALL( SCIPchgVarObjProbing(masterprob, auxiliaryvar, 1.0) );
534 
535  /* solving the probing LP to get a lower bound on the auxiliary variables */
536  SCIP_CALL( SCIPsolveProbingLP(masterprob, -1, &lperror, &cutoff) );
537 
538  if( !SCIPisInfinity(masterprob, -SCIPgetSolTransObj(masterprob, NULL)) )
540 
541  SCIPdebugMsg(masterprob, "Cut constant for subproblem %d: %g\n", i,
543 
544  SCIP_CALL( SCIPchgVarObjProbing(masterprob, auxiliaryvar, 0.0) );
545  }
546 
547  SCIP_CALL( SCIPendProbing(masterprob) );
548 
549  return SCIP_OKAY;
550 }
551 
552 /** exec the event handler */
553 static
554 SCIP_DECL_EVENTEXEC(eventExecBendersNodesolved)
555 { /*lint --e{715}*/
556  SCIP_BENDERS* benders;
557 
558  assert(scip != NULL);
559  assert(eventhdlr != NULL);
560  assert(strcmp(SCIPeventhdlrGetName(eventhdlr), NODESOLVED_EVENTHDLR_NAME) == 0);
561 
562  benders = (SCIP_BENDERS*)SCIPeventhdlrGetData(eventhdlr); /*lint !e826*/
563 
564  if( SCIPbendersGetNSubproblems(benders) > 0
566  {
568  }
569 
571 
572  return SCIP_OKAY;
573 }
574 
575 /** solving process initialization method of event handler (called when branch and bound process is about to begin) */
576 static
577 SCIP_DECL_EVENTINITSOL(eventInitsolBendersNodesolved)
578 {
579  SCIP_BENDERS* benders;
580 
581  assert(scip != NULL);
582  assert(eventhdlr != NULL);
583  assert(strcmp(SCIPeventhdlrGetName(eventhdlr), NODESOLVED_EVENTHDLR_NAME) == 0);
584 
585  /* getting the Benders' decomposition data structure */
586  benders = (SCIP_BENDERS*)SCIPeventhdlrGetData(eventhdlr); /*lint !e826*/
587 
588  /* The event is only caught if there is an active Benders' decomposition, the integer subproblem are solved and
589  * the Benders' decomposition has not been copied in thread safe mode
590  */
592  && !benders->threadsafe )
593  {
595  }
596 
597  return SCIP_OKAY;
598 }
599 
600 
601 /* ---------------- Methods for the parallelisation of Benders' decomposition ---------------- */
602 
603 /** comparison method for sorting the subproblems.
604  * The subproblem that has been called the least is prioritised
605  */
606 static
607 SCIP_DECL_SORTPTRCOMP(benderssubcompdefault)
608 {
609  SCIP_SUBPROBLEMSOLVESTAT* solvestat1;
610  SCIP_SUBPROBLEMSOLVESTAT* solvestat2;
611 
612  assert(elem1 != NULL);
613  assert(elem2 != NULL);
614 
615  solvestat1 = (SCIP_SUBPROBLEMSOLVESTAT*)elem1;
616  solvestat2 = (SCIP_SUBPROBLEMSOLVESTAT*)elem2;
617 
618  /* prefer subproblems with fewer calls, using the index as tie breaker */
619  if( MAX(solvestat1->ncalls, solvestat2->ncalls) == 0 )
620  return solvestat1->idx - solvestat2->idx;
621  else if( solvestat1->ncalls != solvestat2->ncalls )
622  return solvestat1->ncalls - solvestat2->ncalls;
623  else
624  {
625  /* prefer the harder problem (with more average iterations) */
626  int avgiterdiff = (int)solvestat2->avgiter - (int)solvestat1->avgiter;
627 
628  if( avgiterdiff != 0 )
629  return avgiterdiff;
630 
631  return solvestat1->idx - solvestat2->idx;
632  }
633 
634 /* the code below does not give a total order of the elements */
635 #ifdef SCIP_DISABLED_CODE
636  if( solvestat1->ncalls == 0 )
637  if( solvestat2->ncalls == 0 )
638  if( solvestat1->idx < solvestat2->idx )
639  return -1;
640  else
641  return 1;
642  else
643  return -1;
644  else if( solvestat2->ncalls == 0 )
645  return 1;
646  else
647  {
648  if( solvestat1->ncalls < solvestat2->ncalls )
649  return -1;
650  else if( solvestat2->ncalls < solvestat1->ncalls )
651  return 1;
652  else
653  {
654  /* we want to execute the hard subproblems first */
655  if( solvestat1->avgiter > solvestat2->avgiter )
656  return 1;
657  else
658  return -1;
659  }
660  }
661 #endif
662 }
663 
664 /* Local methods */
665 
666 /** A workaround for GCG. This is a temp vardata that is set for the auxiliary variables */
667 struct SCIP_VarData
668 {
669  int vartype; /**< the variable type. In GCG this indicates whether the variable is a
670  * master problem or subproblem variable. */
671 };
672 
673 /** adds the auxiliary variables to the Benders' decomposition master problem */
674 static
676  SCIP* scip, /**< SCIP data structure */
677  SCIP_BENDERS* benders /**< Benders' decomposition structure */
678  )
679 {
680  SCIP_BENDERS* topbenders; /* the highest priority Benders' decomposition */
681  SCIP_VAR* auxiliaryvar;
682  SCIP_VARDATA* vardata;
683  char varname[SCIP_MAXSTRLEN]; /* the name of the auxiliary variable */
684  SCIP_Bool shareauxvars;
685  int i;
686 
687  /* this is a workaround for GCG. GCG expects that the variable has vardata when added. So a dummy vardata is created */
688  SCIP_CALL( SCIPallocBlockMemory(scip, &vardata) );
689  vardata->vartype = -1;
690 
691  /* getting the highest priority Benders' decomposition */
692  topbenders = SCIPgetBenders(scip)[0];
693 
694  /* if the current Benders is the highest priority Benders, then we need to create the auxiliary variables.
695  * Otherwise, if the shareauxvars flag is set, then the auxiliary variables from the highest priority Benders' are
696  * stored with this Benders. */
697  shareauxvars = FALSE;
698  if( topbenders != benders && SCIPbendersShareAuxVars(benders) )
699  shareauxvars = TRUE;
700 
701  for( i = 0; i < SCIPbendersGetNSubproblems(benders); i++ )
702  {
703  /* if the auxiliary variables are shared, then a pointer to the variable is retrieved from topbenders,
704  * otherwise the auxiliaryvariable is created. */
705  if( shareauxvars )
706  {
707  auxiliaryvar = SCIPbendersGetAuxiliaryVar(topbenders, i);
708 
709  SCIP_CALL( SCIPcaptureVar(scip, auxiliaryvar) );
710  }
711  else
712  {
713  SCIP_VARTYPE vartype;
714 
715  /* set the variable type of the auxiliary variables to implicit integer if the objective function of the
716  * subproblem is guaranteed to be integer. This behaviour is controlled through a user parameter.
717  * NOTE: It is only possible to determine if the objective function is integral if the subproblem is defined as
718  * a SCIP instance, i.e. not NULL.
719  */
720  if( benders->auxvarsimplint && SCIPbendersSubproblem(benders, i) != NULL
721  && SCIPisObjIntegral(SCIPbendersSubproblem(benders, i)) )
722  vartype = SCIP_VARTYPE_IMPLINT;
723  else
724  vartype = SCIP_VARTYPE_CONTINUOUS;
725 
726  (void) SCIPsnprintf(varname, SCIP_MAXSTRLEN, "%s_%d_%s", AUXILIARYVAR_NAME, i, SCIPbendersGetName(benders) );
727  SCIP_CALL( SCIPcreateVarBasic(scip, &auxiliaryvar, varname, benders->subproblowerbound[i], SCIPinfinity(scip),
728  1.0, vartype) );
729 
730  SCIPvarSetData(auxiliaryvar, vardata);
731 
732  SCIP_CALL( SCIPaddVar(scip, auxiliaryvar) );
733 
734  /* adding the down lock for the Benders' decomposition constraint handler */
735  SCIP_CALL( SCIPaddVarLocksType(scip, auxiliaryvar, SCIP_LOCKTYPE_MODEL, 1, 0) );
736  }
737 
738  benders->auxiliaryvars[i] = auxiliaryvar;
739  }
740 
741  SCIPfreeBlockMemory(scip, &vardata);
742 
743  return SCIP_OKAY;
744 }
745 
746 /** assigns the copied auxiliary variables in the target SCIP to the target Benders' decomposition data */
747 static
749  SCIP* scip, /**< SCIP data structure, the target scip */
750  SCIP_BENDERS* benders /**< Benders' decomposition */
751  )
752 {
753  SCIP_BENDERS* topbenders; /* the highest priority Benders' decomposition */
754  SCIP_VAR* targetvar;
755  SCIP_VARDATA* vardata;
756  char varname[SCIP_MAXSTRLEN]; /* the name of the auxiliary variable */
757  SCIP_Bool shareauxvars;
758  int subscipdepth;
759  int i;
760  int j;
761 
762  assert(scip != NULL);
763  assert(benders != NULL);
764 
765  /* this is a workaround for GCG. GCG expects that the variable has vardata when added. So a dummy vardata is created */
766  SCIP_CALL( SCIPallocBlockMemory(scip, &vardata) );
767  vardata->vartype = -1;
768 
769  /* getting the highest priority Benders' decomposition */
770  topbenders = SCIPgetBenders(scip)[0];
771 
772  /* if the auxiliary variable are shared, then the variable name will have a suffix of the highest priority Benders'
773  * name. So the shareauxvars flag indicates how to search for the auxiliary variables */
774  shareauxvars = FALSE;
775  if( topbenders != benders && SCIPbendersShareAuxVars(benders) )
776  shareauxvars = TRUE;
777 
778  subscipdepth = SCIPgetSubscipDepth(scip);
779 
780  for( i = 0; i < SCIPbendersGetNSubproblems(benders); i++ )
781  {
782  char prefix[SCIP_MAXSTRLEN];
783  char tmpprefix[SCIP_MAXSTRLEN];
784  int len = 1;
785 
786  j = 0;
787  targetvar = NULL;
788 
789  /* the prefix for the variable names is required for UG, since we don't know how many copies have been made. To
790  * find the target variable, we start with an empty prefix. Then t_ is prepended until the target variable is
791  * found
792  */
793  prefix[0] = '\0';
794  while( targetvar == NULL && j <= subscipdepth )
795  {
796  if( shareauxvars )
797  (void) SCIPsnprintf(varname, SCIP_MAXSTRLEN, "%s%s_%d_%s", prefix, AUXILIARYVAR_NAME, i, SCIPbendersGetName(topbenders));
798  else
799  (void) SCIPsnprintf(varname, SCIP_MAXSTRLEN, "%s%s_%d_%s", prefix, AUXILIARYVAR_NAME, i, SCIPbendersGetName(benders));
800 
801  /* finding the variable in the copied problem that has the same name as the auxiliary variable */
802  targetvar = SCIPfindVar(scip, varname);
803 
804  (void) SCIPsnprintf(tmpprefix, len, "t_%s", prefix);
805  len += 2;
806  (void) strncpy(prefix, tmpprefix, len); /*lint !e732*/
807 
808  j++;
809  }
810 
811  if( targetvar != NULL )
812  {
813  SCIPvarSetData(targetvar, vardata);
814 
815  benders->auxiliaryvars[i] = SCIPvarGetTransVar(targetvar);
816 
817  SCIP_CALL( SCIPcaptureVar(scip, benders->auxiliaryvars[i]) );
818  }
819  else
820  {
821  SCIPABORT();
822  }
823  }
824 
825  SCIPfreeBlockMemory(scip, &vardata);
826 
827  return SCIP_OKAY;
828 }
829 
830 /** sets the subproblem objective value array to -infinity */
831 static
833  SCIP_BENDERS* benders, /**< the Benders' decomposition structure */
834  SCIP_SET* set /**< global SCIP settings */
835  )
836 {
837  SCIP* subproblem;
838  SCIP_Real inf;
839  int nsubproblems;
840  int i;
841 
842  assert(benders != NULL);
843 
844  nsubproblems = SCIPbendersGetNSubproblems(benders);
845 
846  for( i = 0; i < nsubproblems; i++ )
847  {
848  subproblem = SCIPbendersSubproblem(benders, i);
849  if( subproblem != NULL )
850  inf = SCIPinfinity(subproblem);
851  else
852  inf = SCIPsetInfinity(set);
853 
854  SCIPbendersSetSubproblemObjval(benders, i, inf);
855  }
856 }
857 
858 /** compares two Benders' decompositions w. r. to their priority */
859 SCIP_DECL_SORTPTRCOMP(SCIPbendersComp)
860 { /*lint --e{715}*/
861  return ((SCIP_BENDERS*)elem2)->priority - ((SCIP_BENDERS*)elem1)->priority;
862 }
863 
864 /** comparison method for sorting Benders' decompositions w.r.t. to their name */
865 SCIP_DECL_SORTPTRCOMP(SCIPbendersCompName)
866 {
867  return strcmp(SCIPbendersGetName((SCIP_BENDERS*)elem1), SCIPbendersGetName((SCIP_BENDERS*)elem2));
868 }
869 
870 /** method to call, when the priority of a Benders' decomposition was changed */
871 static
872 SCIP_DECL_PARAMCHGD(paramChgdBendersPriority)
873 { /*lint --e{715}*/
874  SCIP_PARAMDATA* paramdata;
875 
876  paramdata = SCIPparamGetData(param);
877  assert(paramdata != NULL);
878 
879  /* use SCIPsetBendersPriority() to mark the Benders' decompositions as unsorted */
880  SCIPsetBendersPriority(scip, (SCIP_BENDERS*)paramdata, SCIPparamGetInt(param)); /*lint !e740*/
881 
882  return SCIP_OKAY;
883 }
884 
885 /** creates a variable mapping between the master problem variables of the source scip and the sub scip */
886 static
888  SCIP_BENDERS* benders, /**< Benders' decomposition of the target SCIP instance */
889  SCIP_SET* sourceset, /**< global SCIP settings from the source SCIP */
890  SCIP_HASHMAP* varmap /**< a hashmap to store the mapping of source variables corresponding
891  * target variables; must not be NULL */
892  )
893 {
894  SCIP_VAR** vars;
895  SCIP_VAR* targetvar;
896  int nvars;
897  int i;
898 
899  assert(benders != NULL);
900  assert(sourceset != NULL);
901  assert(benders->iscopy);
902  assert(benders->mastervarsmap == NULL);
903 
904  /* getting the master problem variable data */
905  vars = SCIPgetVars(sourceset->scip);
906  nvars = SCIPgetNVars(sourceset->scip);
907 
908  /* creating the hashmap for the mapping between the master variable of the target and source scip */
909  SCIP_CALL( SCIPhashmapCreate(&benders->mastervarsmap, SCIPblkmem(sourceset->scip), nvars) );
910 
911  for( i = 0; i < nvars; i++ )
912  {
913  /* getting the variable pointer for the target SCIP variables. The variable mapping returns the target SCIP
914  * variable for a given source SCIP variable. */
915  targetvar = (SCIP_VAR*) SCIPhashmapGetImage(varmap, vars[i]);
916  if( targetvar != NULL )
917  {
918  SCIP_CALL( SCIPhashmapInsert(benders->mastervarsmap, targetvar, vars[i]) );
919  SCIP_CALL( SCIPcaptureVar(sourceset->scip, vars[i]) );
920  }
921  }
922 
923  return SCIP_OKAY;
924 }
925 
926 /** copies the given Benders' decomposition to a new SCIP */
928  SCIP_BENDERS* benders, /**< Benders' decomposition */
929  SCIP_SET* sourceset, /**< SCIP_SET of SCIP to copy from */
930  SCIP_SET* targetset, /**< SCIP_SET of SCIP to copy to */
931  SCIP_HASHMAP* varmap, /**< a hashmap to store the mapping of source variables corresponding
932  * target variables; if NULL, then the transfer of cuts is not possible */
933  SCIP_Bool threadsafe, /**< must the Benders' decomposition copy be thread safe */
934  SCIP_Bool* valid /**< was the copying process valid? */
935  )
936 {
937  SCIP_BENDERS* targetbenders; /* the copy of the Benders' decomposition struct in the target set */
938  int i;
939 
940  assert(benders != NULL);
941  assert(targetset != NULL);
942  assert(valid != NULL);
943  assert(targetset->scip != NULL);
944 
945  (*valid) = FALSE;
946 
947  if( benders->benderscopy != NULL && targetset->benders_copybenders && SCIPbendersIsActive(benders) )
948  {
949  SCIPsetDebugMsg(targetset, "including Benders' decomposition %s in subscip %p\n", SCIPbendersGetName(benders), (void*)targetset->scip);
950  SCIP_CALL( benders->benderscopy(targetset->scip, benders, threadsafe) );
951 
952  /* copying the Benders' cuts */
953  targetbenders = SCIPsetFindBenders(targetset, SCIPbendersGetName(benders));
954 
955  /* storing the pointer to the source scip instance */
956  targetbenders->sourcescip = sourceset->scip;
957 
958  /* the flag is set to indicate that the Benders' decomposition is a copy */
959  targetbenders->iscopy = TRUE;
960 
961  /* storing whether the lnscheck should be performed */
962  targetbenders->lnscheck = benders->lnscheck;
963  targetbenders->lnsmaxdepth = benders->lnsmaxdepth;
964  targetbenders->lnsmaxcalls = benders->lnsmaxcalls;
965  targetbenders->lnsmaxcallsroot = benders->lnsmaxcallsroot;
966 
967  /* storing whether the Benders' copy required thread safety */
968  targetbenders->threadsafe = threadsafe;
969 
970  /* calling the copy method for the Benders' cuts */
972  for( i = 0; i < benders->nbenderscuts; i++ )
973  {
974  SCIP_CALL( SCIPbenderscutCopyInclude(targetbenders, benders->benderscuts[i], targetset) );
975  }
976 
977  /* When the Benders' decomposition is copied then a variable mapping between the master problem variables is
978  * required. This variable mapping is used to transfer the cuts generated in the target SCIP to the source SCIP.
979  * The variable map is stored in the target Benders' decomposition. This will be freed when the sub-SCIP is freed.
980  */
981  if( varmap != NULL )
982  {
983  SCIP_CALL( createMasterVarMapping(targetbenders, sourceset, varmap) );
984  }
985 
986  assert((varmap != NULL && targetbenders->mastervarsmap != NULL)
987  || (varmap == NULL && targetbenders->mastervarsmap == NULL));
988  }
989 
990  /* if the Benders' decomposition is active, then copy is not valid. */
991  (*valid) = !SCIPbendersIsActive(benders);
992 
993  return SCIP_OKAY;
994 }
995 
996 /** internal method for creating a Benders' decomposition structure */
997 static
999  SCIP_BENDERS** benders, /**< pointer to Benders' decomposition data structure */
1000  SCIP_SET* set, /**< global SCIP settings */
1001  SCIP_MESSAGEHDLR* messagehdlr, /**< message handler */
1002  BMS_BLKMEM* blkmem, /**< block memory for parameter settings */
1003  const char* name, /**< name of Benders' decomposition */
1004  const char* desc, /**< description of Benders' decomposition */
1005  int priority, /**< priority of the Benders' decomposition */
1006  SCIP_Bool cutlp, /**< should Benders' cuts be generated for LP solutions */
1007  SCIP_Bool cutpseudo, /**< should Benders' cuts be generated for pseudo solutions */
1008  SCIP_Bool cutrelax, /**< should Benders' cuts be generated for relaxation solutions */
1009  SCIP_Bool shareauxvars, /**< should this Benders' use the highest priority Benders aux vars */
1010  SCIP_DECL_BENDERSCOPY ((*benderscopy)), /**< copy method of Benders' decomposition or NULL if you don't want to copy your plugin into sub-SCIPs */
1011  SCIP_DECL_BENDERSFREE ((*bendersfree)), /**< destructor of Benders' decomposition */
1012  SCIP_DECL_BENDERSINIT ((*bendersinit)), /**< initialize Benders' decomposition */
1013  SCIP_DECL_BENDERSEXIT ((*bendersexit)), /**< deinitialize Benders' decomposition */
1014  SCIP_DECL_BENDERSINITPRE((*bendersinitpre)),/**< presolving initialization method for Benders' decomposition */
1015  SCIP_DECL_BENDERSEXITPRE((*bendersexitpre)),/**< presolving deinitialization method for Benders' decomposition */
1016  SCIP_DECL_BENDERSINITSOL((*bendersinitsol)),/**< solving process initialization method of Benders' decomposition */
1017  SCIP_DECL_BENDERSEXITSOL((*bendersexitsol)),/**< solving process deinitialization method of Benders' decomposition */
1018  SCIP_DECL_BENDERSGETVAR((*bendersgetvar)),/**< returns the master variable for a given subproblem variable */
1019  SCIP_DECL_BENDERSCREATESUB((*benderscreatesub)),/**< creates a Benders' decomposition subproblem */
1020  SCIP_DECL_BENDERSPRESUBSOLVE((*benderspresubsolve)),/**< called prior to the subproblem solving loop */
1021  SCIP_DECL_BENDERSSOLVESUBCONVEX((*benderssolvesubconvex)),/**< the solving method for convex Benders' decomposition subproblems */
1022  SCIP_DECL_BENDERSSOLVESUB((*benderssolvesub)),/**< the solving method for the Benders' decomposition subproblems */
1023  SCIP_DECL_BENDERSPOSTSOLVE((*benderspostsolve)),/**< called after the subproblems are solved. */
1024  SCIP_DECL_BENDERSFREESUB((*bendersfreesub)),/**< the freeing method for the Benders' decomposition subproblems */
1025  SCIP_BENDERSDATA* bendersdata /**< Benders' decomposition data */
1026  )
1027 {
1028  char paramname[SCIP_MAXSTRLEN];
1029  char paramdesc[SCIP_MAXSTRLEN];
1030 
1031  assert(benders != NULL);
1032  assert(name != NULL);
1033  assert(desc != NULL);
1034 
1035  /* Checking whether the benderssolvesub and the bendersfreesub are both implemented or both are not implemented */
1036  if( (benderssolvesubconvex == NULL && benderssolvesub == NULL && bendersfreesub != NULL)
1037  || ((benderssolvesubconvex != NULL || benderssolvesub != NULL) && bendersfreesub == NULL) )
1038  {
1039  SCIPerrorMessage("Benders' decomposition <%s> requires that if bendersFreesub%s is implemented, then at least "
1040  "one of bendersSolvesubconvex%s or bendersSolvesub%s are implemented.\n", name, name, name, name);
1041  return SCIP_INVALIDCALL;
1042  }
1043 
1044  SCIP_ALLOC( BMSallocMemory(benders) );
1045  BMSclearMemory(*benders);
1046  SCIP_ALLOC( BMSduplicateMemoryArray(&(*benders)->name, name, strlen(name)+1) );
1047  SCIP_ALLOC( BMSduplicateMemoryArray(&(*benders)->desc, desc, strlen(desc)+1) );
1048  (*benders)->priority = priority;
1049  (*benders)->cutlp = cutlp;
1050  (*benders)->cutpseudo = cutpseudo;
1051  (*benders)->cutrelax = cutrelax;
1052  (*benders)->shareauxvars = shareauxvars;
1053  (*benders)->benderscopy = benderscopy;
1054  (*benders)->bendersfree = bendersfree;
1055  (*benders)->bendersinit = bendersinit;
1056  (*benders)->bendersexit = bendersexit;
1057  (*benders)->bendersinitpre = bendersinitpre;
1058  (*benders)->bendersexitpre = bendersexitpre;
1059  (*benders)->bendersinitsol = bendersinitsol;
1060  (*benders)->bendersexitsol = bendersexitsol;
1061  (*benders)->bendersgetvar = bendersgetvar;
1062  (*benders)->benderscreatesub = benderscreatesub;
1063  (*benders)->benderspresubsolve = benderspresubsolve;
1064  (*benders)->benderssolvesubconvex = benderssolvesubconvex;
1065  (*benders)->benderssolvesub = benderssolvesub;
1066  (*benders)->benderspostsolve = benderspostsolve;
1067  (*benders)->bendersfreesub = bendersfreesub;
1068  (*benders)->bendersdata = bendersdata;
1069  SCIP_CALL( SCIPclockCreate(&(*benders)->setuptime, SCIP_CLOCKTYPE_DEFAULT) );
1070  SCIP_CALL( SCIPclockCreate(&(*benders)->bendersclock, SCIP_CLOCKTYPE_DEFAULT) );
1071 
1072  /* add parameters */
1073  (void) SCIPsnprintf(paramname, SCIP_MAXSTRLEN, "benders/%s/priority", name);
1074  (void) SCIPsnprintf(paramdesc, SCIP_MAXSTRLEN, "priority of Benders' decomposition <%s>", name);
1075  SCIP_CALL( SCIPsetAddIntParam(set, messagehdlr, blkmem, paramname, paramdesc,
1076  &(*benders)->priority, FALSE, priority, INT_MIN/4, INT_MAX/4,
1077  paramChgdBendersPriority, (SCIP_PARAMDATA*)(*benders)) ); /*lint !e740*/
1078 
1079  (void) SCIPsnprintf(paramname, SCIP_MAXSTRLEN, "benders/%s/cutlp", name);
1080  SCIP_CALL( SCIPsetAddBoolParam(set, messagehdlr, blkmem, paramname,
1081  "should Benders' cuts be generated for LP solutions?", &(*benders)->cutlp, FALSE, cutlp, NULL, NULL) ); /*lint !e740*/
1082 
1083  (void) SCIPsnprintf(paramname, SCIP_MAXSTRLEN, "benders/%s/cutpseudo", name);
1084  SCIP_CALL( SCIPsetAddBoolParam(set, messagehdlr, blkmem, paramname,
1085  "should Benders' cuts be generated for pseudo solutions?", &(*benders)->cutpseudo, FALSE, cutpseudo, NULL, NULL) ); /*lint !e740*/
1086 
1087  (void) SCIPsnprintf(paramname, SCIP_MAXSTRLEN, "benders/%s/cutrelax", name);
1088  SCIP_CALL( SCIPsetAddBoolParam(set, messagehdlr, blkmem, paramname,
1089  "should Benders' cuts be generated for relaxation solutions?", &(*benders)->cutrelax, FALSE, cutrelax, NULL, NULL) ); /*lint !e740*/
1090 
1091  /* These parameters are left for the user to decide in a settings file. This departs from the usual SCIP convention
1092  * where the settings available at the creation of the plugin can be set in the function call.
1093  */
1094  (void) SCIPsnprintf(paramname, SCIP_MAXSTRLEN, "benders/%s/transfercuts", name);
1095  SCIP_CALL( SCIPsetAddBoolParam(set, messagehdlr, blkmem, paramname,
1096  "should Benders' cuts from LNS heuristics be transferred to the main SCIP instance?", &(*benders)->transfercuts,
1097  FALSE, SCIP_DEFAULT_TRANSFERCUTS, NULL, NULL) ); /*lint !e740*/
1098 
1099  (void) SCIPsnprintf(paramname, SCIP_MAXSTRLEN, "benders/%s/lnscheck", name);
1100  SCIP_CALL( SCIPsetAddBoolParam(set, messagehdlr, blkmem, paramname,
1101  "should Benders' decomposition be used in LNS heurisics?", &(*benders)->lnscheck, FALSE, SCIP_DEFAULT_LNSCHECK,
1102  NULL, NULL) ); /*lint !e740*/
1103 
1104  (void) SCIPsnprintf(paramname, SCIP_MAXSTRLEN, "benders/%s/lnsmaxdepth", name);
1105  SCIP_CALL( SCIPsetAddIntParam(set, messagehdlr, blkmem, paramname,
1106  "maximum depth at which the LNS check is performed (-1: no limit)", &(*benders)->lnsmaxdepth, TRUE,
1108 
1109  (void) SCIPsnprintf(paramname, SCIP_MAXSTRLEN, "benders/%s/lnsmaxcalls", name);
1110  SCIP_CALL( SCIPsetAddIntParam(set, messagehdlr, blkmem, paramname,
1111  "the maximum number of Benders' decomposition calls in LNS heuristics (-1: no limit)", &(*benders)->lnsmaxcalls,
1112  TRUE, SCIP_DEFAULT_LNSMAXCALLS, -1, INT_MAX, NULL, NULL) );
1113 
1114  (void) SCIPsnprintf(paramname, SCIP_MAXSTRLEN, "benders/%s/lnsmaxcallsroot", name);
1115  SCIP_CALL( SCIPsetAddIntParam(set, messagehdlr, blkmem, paramname,
1116  "the maximum number of root node Benders' decomposition calls in LNS heuristics (-1: no limit)",
1117  &(*benders)->lnsmaxcallsroot, TRUE, SCIP_DEFAULT_LNSMAXCALLSROOT, -1, INT_MAX, NULL, NULL) );
1118 
1119  (void) SCIPsnprintf(paramname, SCIP_MAXSTRLEN, "benders/%s/cutsasconss", name);
1120  SCIP_CALL( SCIPsetAddBoolParam(set, messagehdlr, blkmem, paramname,
1121  "should the transferred cuts be added as constraints?", &(*benders)->cutsasconss, FALSE,
1122  SCIP_DEFAULT_CUTSASCONSS, NULL, NULL) ); /*lint !e740*/
1123 
1124  (void) SCIPsnprintf(paramname, SCIP_MAXSTRLEN, "benders/%s/subprobfrac", name);
1125  SCIP_CALL( SCIPsetAddRealParam(set, messagehdlr, blkmem, paramname,
1126  "fraction of subproblems that are solved in each iteration", &(*benders)->subprobfrac, FALSE,
1127  SCIP_DEFAULT_SUBPROBFRAC, 0.0, 1.0, NULL, NULL) ); /*lint !e740*/
1128 
1129  (void) SCIPsnprintf(paramname, SCIP_MAXSTRLEN, "benders/%s/updateauxvarbound", name);
1130  SCIP_CALL( SCIPsetAddBoolParam(set, messagehdlr, blkmem, paramname,
1131  "should the auxiliary variable bound be updated by solving the subproblem?", &(*benders)->updateauxvarbound,
1132  FALSE, SCIP_DEFAULT_UPDATEAUXVARBOUND, NULL, NULL) ); /*lint !e740*/
1133 
1134  (void) SCIPsnprintf(paramname, SCIP_MAXSTRLEN, "benders/%s/auxvarsimplint", name);
1135  SCIP_CALL( SCIPsetAddBoolParam(set, messagehdlr, blkmem, paramname,
1136  "if the subproblem objective is integer, then define the auxiliary variables as implicit integers?",
1137  &(*benders)->auxvarsimplint, FALSE, SCIP_DEFAULT_AUXVARSIMPLINT, NULL, NULL) ); /*lint !e740*/
1138 
1139  (void) SCIPsnprintf(paramname, SCIP_MAXSTRLEN, "benders/%s/cutcheck", name);
1140  SCIP_CALL( SCIPsetAddBoolParam(set, messagehdlr, blkmem, paramname,
1141  "should Benders' cuts be generated while checking solutions?",
1142  &(*benders)->cutcheck, FALSE, SCIP_DEFAULT_CUTCHECK, NULL, NULL) ); /*lint !e740*/
1143 
1144  (void) SCIPsnprintf(paramname, SCIP_MAXSTRLEN, "benders/%s/cutstrengthenmult", name);
1145  SCIP_CALL( SCIPsetAddRealParam(set, messagehdlr, blkmem, paramname,
1146  "the convex combination multiplier for the cut strengthening", &(*benders)->convexmult, FALSE,
1147  SCIP_DEFAULT_STRENGTHENMULT, 0.0, 1.0, NULL, NULL) ); /*lint !e740*/
1148 
1149  (void) SCIPsnprintf(paramname, SCIP_MAXSTRLEN, "benders/%s/noimprovelimit", name);
1150  SCIP_CALL( SCIPsetAddIntParam(set, messagehdlr, blkmem, paramname,
1151  "the maximum number of cut strengthening without improvement", &(*benders)->noimprovelimit, TRUE,
1152  SCIP_DEFAULT_NOIMPROVELIMIT, 0, INT_MAX, NULL, NULL) );
1153 
1154  (void) SCIPsnprintf(paramname, SCIP_MAXSTRLEN, "benders/%s/corepointperturb", name);
1155  SCIP_CALL( SCIPsetAddRealParam(set, messagehdlr, blkmem, paramname,
1156  "the constant use to perturb the cut strengthening core point", &(*benders)->perturbeps, FALSE,
1157  SCIP_DEFAULT_STRENGTHENPERTURB, 0.0, 1.0, NULL, NULL) ); /*lint !e740*/
1158 
1159  (void) SCIPsnprintf(paramname, SCIP_MAXSTRLEN, "benders/%s/cutstrengthenenabled", name);
1160  SCIP_CALL( SCIPsetAddBoolParam(set, messagehdlr, blkmem, paramname,
1161  "should the core point cut strengthening be employed (only applied to fractional solutions or continuous subproblems)?",
1162  &(*benders)->strengthenenabled, FALSE, SCIP_DEFAULT_STRENGTHENENABLED, NULL, NULL) ); /*lint !e740*/
1163 
1164  (void) SCIPsnprintf(paramname, SCIP_MAXSTRLEN, "benders/%s/cutstrengthenintpoint", name);
1165  SCIP_CALL( SCIPsetAddCharParam(set, messagehdlr, blkmem, paramname,
1166  "where should the strengthening interior point be sourced from ('l'p relaxation, 'f'irst solution, 'i'ncumbent solution, 'r'elative interior point, vector of 'o'nes, vector of 'z'eros)",
1167  &(*benders)->strengthenintpoint, FALSE, SCIP_DEFAULT_STRENGTHENINTPOINT, "lfiroz", NULL, NULL) ); /*lint !e740*/
1168 
1169  (void) SCIPsnprintf(paramname, SCIP_MAXSTRLEN, "benders/%s/numthreads", name);
1170  SCIP_CALL( SCIPsetAddIntParam(set, messagehdlr, blkmem, paramname,
1171  "the number of threads to use when solving the subproblems", &(*benders)->numthreads, TRUE,
1172  SCIP_DEFAULT_NUMTHREADS, 1, INT_MAX, NULL, NULL) );
1173 
1174  (void) SCIPsnprintf(paramname, SCIP_MAXSTRLEN, "benders/%s/execfeasphase", name);
1175  SCIP_CALL( SCIPsetAddBoolParam(set, messagehdlr, blkmem, paramname,
1176  "should a feasibility phase be executed during the root node, i.e. adding slack variables to constraints to ensure feasibility",
1177  &(*benders)->execfeasphase, FALSE, SCIP_DEFAULT_EXECFEASPHASE, NULL, NULL) ); /*lint !e740*/
1178 
1179  (void) SCIPsnprintf(paramname, SCIP_MAXSTRLEN, "benders/%s/slackvarcoef", name);
1180  SCIP_CALL( SCIPsetAddRealParam(set, messagehdlr, blkmem, paramname,
1181  "the initial objective coefficient of the slack variables in the subproblem", &(*benders)->slackvarcoef, FALSE,
1182  SCIP_DEFAULT_SLACKVARCOEF, 0.0, SCIPsetInfinity(set), NULL, NULL) ); /*lint !e740*/
1183 
1184  (void) SCIPsnprintf(paramname, SCIP_MAXSTRLEN, "benders/%s/maxslackvarcoef", name);
1185  SCIP_CALL( SCIPsetAddRealParam(set, messagehdlr, blkmem, paramname,
1186  "the maximal objective coefficient of the slack variables in the subproblem", &(*benders)->maxslackvarcoef, FALSE,
1187  SCIP_DEFAULT_MAXSLACKVARCOEF, 0.0, SCIPsetInfinity(set), NULL, NULL) ); /*lint !e740*/
1188 
1189  (void) SCIPsnprintf(paramname, SCIP_MAXSTRLEN, "benders/%s/checkconsconvexity", name);
1190  SCIP_CALL( SCIPsetAddBoolParam(set, messagehdlr, blkmem, paramname,
1191  "should the constraints of the subproblems be checked for convexity?", &(*benders)->checkconsconvexity, FALSE,
1192  SCIP_DEFAULT_CHECKCONSCONVEXITY, NULL, NULL) ); /*lint !e740*/
1193 
1194  return SCIP_OKAY;
1195 }
1196 
1197 /** creates a Benders' decomposition structure
1198  *
1199  * To use the Benders' decomposition for solving a problem, it first has to be activated with a call to SCIPactivateBenders().
1200  */
1202  SCIP_BENDERS** benders, /**< pointer to Benders' decomposition data structure */
1203  SCIP_SET* set, /**< global SCIP settings */
1204  SCIP_MESSAGEHDLR* messagehdlr, /**< message handler */
1205  BMS_BLKMEM* blkmem, /**< block memory for parameter settings */
1206  const char* name, /**< name of Benders' decomposition */
1207  const char* desc, /**< description of Benders' decomposition */
1208  int priority, /**< priority of the Benders' decomposition */
1209  SCIP_Bool cutlp, /**< should Benders' cuts be generated for LP solutions */
1210  SCIP_Bool cutpseudo, /**< should Benders' cuts be generated for pseudo solutions */
1211  SCIP_Bool cutrelax, /**< should Benders' cuts be generated for relaxation solutions */
1212  SCIP_Bool shareauxvars, /**< should this Benders' use the highest priority Benders aux vars */
1213  SCIP_DECL_BENDERSCOPY ((*benderscopy)), /**< copy method of Benders' decomposition or NULL if you don't want to copy your plugin into sub-SCIPs */
1214  SCIP_DECL_BENDERSFREE ((*bendersfree)), /**< destructor of Benders' decomposition */
1215  SCIP_DECL_BENDERSINIT ((*bendersinit)), /**< initialize Benders' decomposition */
1216  SCIP_DECL_BENDERSEXIT ((*bendersexit)), /**< deinitialize Benders' decomposition */
1217  SCIP_DECL_BENDERSINITPRE((*bendersinitpre)),/**< presolving initialization method for Benders' decomposition */
1218  SCIP_DECL_BENDERSEXITPRE((*bendersexitpre)),/**< presolving deinitialization method for Benders' decomposition */
1219  SCIP_DECL_BENDERSINITSOL((*bendersinitsol)),/**< solving process initialization method of Benders' decomposition */
1220  SCIP_DECL_BENDERSEXITSOL((*bendersexitsol)),/**< solving process deinitialization method of Benders' decomposition */
1221  SCIP_DECL_BENDERSGETVAR((*bendersgetvar)),/**< returns the master variable for a given subproblem variable */
1222  SCIP_DECL_BENDERSCREATESUB((*benderscreatesub)),/**< creates a Benders' decomposition subproblem */
1223  SCIP_DECL_BENDERSPRESUBSOLVE((*benderspresubsolve)),/**< called prior to the subproblem solving loop */
1224  SCIP_DECL_BENDERSSOLVESUBCONVEX((*benderssolvesubconvex)),/**< the solving method for convex Benders' decomposition subproblems */
1225  SCIP_DECL_BENDERSSOLVESUB((*benderssolvesub)),/**< the solving method for the Benders' decomposition subproblems */
1226  SCIP_DECL_BENDERSPOSTSOLVE((*benderspostsolve)),/**< called after the subproblems are solved. */
1227  SCIP_DECL_BENDERSFREESUB((*bendersfreesub)),/**< the freeing method for the Benders' decomposition subproblems */
1228  SCIP_BENDERSDATA* bendersdata /**< Benders' decomposition data */
1229  )
1230 {
1231  assert(benders != NULL);
1232  assert(name != NULL);
1233  assert(desc != NULL);
1234 
1235  SCIP_CALL_FINALLY( doBendersCreate(benders, set, messagehdlr, blkmem, name, desc, priority, cutlp, cutpseudo,
1236  cutrelax, shareauxvars, benderscopy, bendersfree, bendersinit, bendersexit, bendersinitpre, bendersexitpre,
1237  bendersinitsol, bendersexitsol, bendersgetvar, benderscreatesub, benderspresubsolve, benderssolvesubconvex,
1238  benderssolvesub, benderspostsolve, bendersfreesub, bendersdata), (void) SCIPbendersFree(benders, set) );
1239 
1240  return SCIP_OKAY;
1241 }
1242 
1243 
1244 /** releases the variables that have been captured in the hashmap */
1245 static
1247  SCIP* scip, /**< the SCIP data structure */
1248  SCIP_BENDERS* benders /**< Benders' decomposition */
1249  )
1250 {
1251  int nentries;
1252  int i;
1253 
1254  assert(scip != NULL);
1255  assert(benders != NULL);
1256 
1257  assert(benders->mastervarsmap != NULL);
1258 
1259  nentries = SCIPhashmapGetNEntries(benders->mastervarsmap);
1260 
1261  for( i = 0; i < nentries; ++i )
1262  {
1263  SCIP_HASHMAPENTRY* entry;
1264  entry = SCIPhashmapGetEntry(benders->mastervarsmap, i);
1265 
1266  if( entry != NULL )
1267  {
1268  SCIP_VAR* var;
1269  var = (SCIP_VAR*) SCIPhashmapEntryGetImage(entry);
1270 
1271  SCIP_CALL( SCIPreleaseVar(scip, &var) );
1272  }
1273  }
1274 
1275  return SCIP_OKAY;
1276 }
1277 
1278 
1279 /** calls destructor and frees memory of Benders' decomposition */
1281  SCIP_BENDERS** benders, /**< pointer to Benders' decomposition data structure */
1282  SCIP_SET* set /**< global SCIP settings */
1283  )
1284 {
1285  int i;
1286 
1287  assert(benders != NULL);
1288  assert(*benders != NULL);
1289  assert(!(*benders)->initialized);
1290  assert(set != NULL);
1291 
1292  /* call destructor of Benders' decomposition */
1293  if( (*benders)->bendersfree != NULL )
1294  {
1295  SCIP_CALL( (*benders)->bendersfree(set->scip, *benders) );
1296  }
1297 
1298  /* if the Benders' decomposition is a copy and a varmap has been passed to SCIP_BENDERS, then the variable map
1299  * between the source and the target SCIP needs to be freed.
1300  */
1301  if( (*benders)->iscopy && (*benders)->mastervarsmap != NULL )
1302  {
1303  SCIP_CALL( releaseVarMappingHashmapVars((*benders)->sourcescip, (*benders)) );
1304  SCIPhashmapFree(&(*benders)->mastervarsmap);
1305  }
1306 
1307  /* freeing the Benders' cuts */
1308  for( i = 0; i < (*benders)->nbenderscuts; i++ )
1309  {
1310  SCIP_CALL( SCIPbenderscutFree(&((*benders)->benderscuts[i]), set) );
1311  }
1312  BMSfreeMemoryArrayNull(&(*benders)->benderscuts);
1313 
1314  SCIPclockFree(&(*benders)->bendersclock);
1315  SCIPclockFree(&(*benders)->setuptime);
1316  BMSfreeMemoryArray(&(*benders)->name);
1317  BMSfreeMemoryArray(&(*benders)->desc);
1318  BMSfreeMemory(benders);
1319 
1320  return SCIP_OKAY;
1321 }
1322 
1323 /* adds a slack variable to the given constraint */
1324 static
1326  SCIP* scip, /**< the SCIP data structure */
1327  SCIP_BENDERS* benders, /**< Benders' decomposition */
1328  SCIP_CONS* cons, /**< constraint to which the slack variable(s) is added to */
1329  SCIP_CONSHDLR** linearconshdlrs, /**< an array storing the linear constraint handlers */
1330  SCIP_CONSHDLR* nlconshdlr, /**< pointer to the nonlinear constraint handler */
1331  int nlinearconshdlrs /**< the number of linear constraint handlers */
1332  )
1333 {
1334  SCIP_CONSHDLR* conshdlr;
1335  SCIP_VAR* var;
1336  SCIP_Real rhs;
1337  SCIP_Real lhs;
1338  SCIP_Real objcoef;
1339  int i;
1340  SCIP_Bool linearcons;
1341  SCIP_Bool success;
1342  char name[SCIP_MAXSTRLEN];
1343 
1344  conshdlr = SCIPconsGetHdlr(cons);
1345 
1346  /* assume that the constraint is not linear, then we check whether it is linear */
1347  linearcons = FALSE;
1348 
1349  /* checking whether the constraint is a linear constraint. If so, we add a coefficient to the constraint */
1350  for( i = 0; i < nlinearconshdlrs; ++i )
1351  {
1352  if( conshdlr == linearconshdlrs[i] )
1353  {
1354  linearcons = TRUE;
1355  break;
1356  }
1357  }
1358 
1359  if( !linearcons && conshdlr != nlconshdlr )
1360  {
1361  SCIPwarningMessage(scip, "The subproblem includes constraint <%s>. "
1362  "This is not supported and the slack variable will not be added to the constraint. Feasibility cuts may be invalid.\n",
1363  SCIPconshdlrGetName(conshdlr));
1364  }
1365 
1366  if( linearcons )
1367  {
1368  rhs = SCIPconsGetRhs(scip, cons, &success);
1369  assert(success);
1370  lhs = SCIPconsGetLhs(scip, cons, &success);
1371  assert(success);
1372  }
1373  else
1374  {
1375  rhs = SCIPgetRhsNonlinear(cons);
1376  lhs = SCIPgetLhsNonlinear(cons);
1377  }
1378 
1379  /* getting the objective coefficient for the slack variables */
1380  objcoef = benders->slackvarcoef;
1381 
1382  /* if the right hand side is finite, then we need to add a slack variable with a negative coefficient */
1383  if( !SCIPisInfinity(scip, rhs) )
1384  {
1385  (void) SCIPsnprintf(name, SCIP_MAXSTRLEN, "%s_%s_neg", SLACKVAR_NAME, SCIPconsGetName(cons) );
1386 
1387  SCIP_CALL( SCIPcreateVarBasic(scip, &var, name, 0.0, SCIPinfinity(scip), objcoef, SCIP_VARTYPE_CONTINUOUS) );
1388 
1389  /* adding the slack variable to the subproblem */
1390  SCIP_CALL( SCIPaddVar(scip, var) );
1391 
1392  /* adds the slack variable to the constraint */
1393  if( linearcons )
1394  {
1395  SCIP_CALL( SCIPconsAddCoef(scip, cons, var, -1.0) );
1396  }
1397  else
1398  {
1399  SCIP_CALL( SCIPaddLinearVarNonlinear(scip, cons, var, -1.0) );
1400  }
1401 
1402  /* releasing the variable */
1403  SCIP_CALL( SCIPreleaseVar(scip, &var) );
1404  }
1405 
1406  /* if the left hand side if finite, then we need to add a slack variable with a positive coefficient */
1407  if( !SCIPisInfinity(scip, -lhs) )
1408  {
1409  (void) SCIPsnprintf(name, SCIP_MAXSTRLEN, "%s_%s_pos", SLACKVAR_NAME, SCIPconsGetName(cons) );
1410 
1411  SCIP_CALL( SCIPcreateVarBasic(scip, &var, name, 0.0, SCIPinfinity(scip), objcoef, SCIP_VARTYPE_CONTINUOUS) );
1412 
1413  /* adding the slack variable to the subproblem */
1414  SCIP_CALL( SCIPaddVar(scip, var) );
1415 
1416  /* adds the slack variable to the constraint */
1417  if( linearcons )
1418  {
1419  SCIP_CALL( SCIPconsAddCoef(scip, cons, var, 1.0) );
1420  }
1421  else
1422  {
1423  SCIP_CALL( SCIPaddLinearVarNonlinear(scip, cons, var, 1.0) );
1424  }
1425 
1426  /* releasing the variable */
1427  SCIP_CALL( SCIPreleaseVar(scip, &var) );
1428  }
1429 
1430  return SCIP_OKAY;
1431 }
1432 
1433 /** adds the slack variables to each of the constraints for the generation of feasibility cuts for the given non-linear
1434  * subproblem
1435  */
1436 static
1438  SCIP_BENDERS* benders, /**< Benders' decomposition */
1439  SCIP_SET* set, /**< global SCIP settings */
1440  int probnumber /**< the subproblem number */
1441  )
1442 {
1443  SCIP* subproblem;
1444  SCIP_CONSHDLR* linearconshdlrs[NLINEARCONSHDLRS];
1445  SCIP_CONSHDLR* nlconshdlr;
1446  SCIP_CONS* cons;
1447  int i;
1448 
1449  assert(benders != NULL);
1450  assert(set != NULL);
1451  assert(probnumber >= 0 && probnumber < SCIPbendersGetNSubproblems(benders));
1452 
1453  subproblem = SCIPbendersSubproblem(benders, probnumber);
1454 
1455  /* get pointers to linear constraints handlers, so can avoid string comparisons */
1456  linearconshdlrs[0] = SCIPfindConshdlr(subproblem, "knapsack");
1457  linearconshdlrs[1] = SCIPfindConshdlr(subproblem, "linear");
1458  linearconshdlrs[2] = SCIPfindConshdlr(subproblem, "logicor");
1459  linearconshdlrs[3] = SCIPfindConshdlr(subproblem, "setppc");
1460  linearconshdlrs[4] = SCIPfindConshdlr(subproblem, "varbound");
1461 
1462  nlconshdlr = SCIPfindConshdlr(subproblem, "nonlinear");
1463 
1464  for( i = 0; i < SCIPgetNOrigConss(subproblem); ++i )
1465  {
1466  cons = SCIPgetOrigConss(subproblem)[i];
1467 
1468  /* adding the slack variables to the constraint */
1469  SCIP_CALL( addSlackVars(subproblem, benders, cons, linearconshdlrs, nlconshdlr, NLINEARCONSHDLRS) );
1470  }
1471 
1472  return SCIP_OKAY;
1473 }
1474 
1475 /** initialises a MIP subproblem by putting the problem into SCIP_STAGE_SOLVING. This is achieved by calling SCIPsolve
1476  * and then interrupting the solve in a node focus event handler.
1477  * The LP subproblem is also initialised using this method; however, a different event handler is added. This event
1478  * handler will put the LP subproblem into probing mode.
1479  * The MIP solving function is called to initialise the subproblem because this function calls SCIPsolve with the
1480  * appropriate parameter settings for Benders' decomposition.
1481  */
1482 static
1484  SCIP_BENDERS* benders, /**< Benders' decomposition */
1485  SCIP_SET* set, /**< global SCIP settings */
1486  int probnumber, /**< the subproblem number */
1487  SCIP_Bool* success /**< was the initialisation process successful */
1488  )
1489 {
1490  SCIP* subproblem;
1491  SCIP_STATUS solvestatus;
1492  SCIP_Bool cutoff;
1493 
1494  assert(benders != NULL);
1495  assert(probnumber >= 0 && probnumber < SCIPbendersGetNSubproblems(benders));
1496  assert(success != NULL);
1497 
1498  (*success) = FALSE;
1499 
1500  subproblem = SCIPbendersSubproblem(benders, probnumber);
1501  assert(subproblem != NULL);
1502 
1503  /* Getting the problem into the right SCIP stage for solving */
1504  SCIP_CALL( SCIPbendersSolveSubproblemCIP(set->scip, benders, probnumber, &solvestatus, FALSE) );
1505 
1506  /* Constructing the LP that can be solved in later iterations */
1507  if( solvestatus != SCIP_STATUS_BESTSOLLIMIT && solvestatus != SCIP_STATUS_TIMELIMIT
1508  && solvestatus != SCIP_STATUS_MEMLIMIT )
1509  {
1510  assert(SCIPgetStage(subproblem) == SCIP_STAGE_SOLVING);
1511 
1512  SCIP_CALL( SCIPconstructLP(subproblem, &cutoff) );
1513  (*success) = TRUE;
1514  }
1515 
1516  return SCIP_OKAY;
1517 }
1518 
1519 
1520 /** initialises an LP subproblem by putting the problem into probing mode. The probing mode is invoked in a node focus
1521  * event handler. This event handler is added just prior to calling the initialise subproblem function.
1522  */
1523 static
1525  SCIP_BENDERS* benders, /**< Benders' decomposition */
1526  SCIP_SET* set, /**< global SCIP settings */
1527  int probnumber /**< the subproblem number */
1528  )
1529 {
1530  SCIP* subproblem;
1531  SCIP_EVENTHDLR* eventhdlr;
1532  SCIP_EVENTHDLRDATA* eventhdlrdata;
1533  SCIP_Bool success;
1534 
1535  assert(benders != NULL);
1536  assert(probnumber >= 0 && probnumber < SCIPbendersGetNSubproblems(benders));
1537 
1538  subproblem = SCIPbendersSubproblem(benders, probnumber);
1539  assert(subproblem != NULL);
1540 
1541  /* include event handler into SCIP */
1542  SCIP_CALL( SCIPallocBlockMemory(subproblem, &eventhdlrdata) );
1543 
1544  SCIP_CALL( initEventhandlerData(subproblem, eventhdlrdata) );
1545 
1547  eventExecBendersNodefocus, eventhdlrdata) );
1548  SCIP_CALL( SCIPsetEventhdlrInitsol(subproblem, eventhdlr, eventInitsolBendersNodefocus) );
1549  SCIP_CALL( SCIPsetEventhdlrExitsol(subproblem, eventhdlr, eventExitsolBendersNodefocus) );
1550  SCIP_CALL( SCIPsetEventhdlrExit(subproblem, eventhdlr, eventExitBendersNodefocus) );
1551  SCIP_CALL( SCIPsetEventhdlrFree(subproblem, eventhdlr, eventFreeBendersNodefocus) );
1552  assert(eventhdlr != NULL);
1553 
1554  /* calling an initial solve to put the problem into probing mode */
1555  SCIP_CALL( initialiseSubproblem(benders, set, probnumber, &success) );
1556 
1557  return SCIP_OKAY; /*lint !e438*/
1558 }
1559 
1560 /** checks whether the convex relaxation of the subproblem is sufficient to solve the original problem to optimality
1561  *
1562  * We check whether we can conclude that the CIP is actually an LP or a convex NLP.
1563  * To do this, we check that all variables are of continuous type and that every constraint is either handled by known
1564  * linear constraint handler (knapsack, linear, logicor, setppc, varbound) or the nonlinear constraint handler.
1565  * In the latter case, we also check whether the nonlinear constraint is convex.
1566  * Further, nonlinear constraints are only considered if an NLP solver interface is available, i.e., and NLP could
1567  * be solved.
1568  * If constraints are present that cannot be identified as linear or convex nonlinear, then we assume that the
1569  * problem is not convex, thus solving its LP or NLP relaxation will not be sufficient.
1570  */
1571 static
1573  SCIP_BENDERS* benders, /**< Benders' decomposition */
1574  SCIP_SET* set, /**< global SCIP settings */
1575  int probnumber /**< the subproblem number, or -1 for the master problem */
1576  )
1577 {
1578  SCIP* subproblem;
1579  SCIP_CONSHDLR* conshdlr;
1580  SCIP_CONS* cons;
1581  SCIP_HASHMAP* assumevarfixed;
1582  SCIP_VAR** vars;
1583  int nvars;
1584  int nbinvars;
1585  int nintvars;
1586  int nimplintvars;
1587  int i;
1588  int j;
1589  SCIP_Bool convexcons;
1590  SCIP_Bool discretevar;
1591  SCIP_Bool isnonlinear;
1592  SCIP_CONSHDLR* linearconshdlrs[NLINEARCONSHDLRS];
1593  SCIP_CONSHDLR* nlconshdlr = NULL;
1594 
1595  assert(benders != NULL);
1596  assert(set != NULL);
1597  assert(probnumber >= -1 && probnumber < SCIPbendersGetNSubproblems(benders));
1598 
1599  assumevarfixed = NULL;
1600  if( probnumber >= 0 )
1601  subproblem = SCIPbendersSubproblem(benders, probnumber);
1602  else
1603  subproblem = set->scip;
1604 
1605  assert(subproblem != NULL);
1606 
1607  convexcons = FALSE;
1608  discretevar = FALSE;
1609  isnonlinear = FALSE;
1610 
1611  /* getting the number of integer and binary variables to determine the problem type */
1612  SCIP_CALL( SCIPgetVarsData(subproblem, &vars, &nvars, &nbinvars, &nintvars, &nimplintvars, NULL) );
1613 
1614  /* if there are any binary, integer or implicit integer variables, then the subproblems is marked as non-convex */
1615  if( nbinvars != 0 || nintvars != 0 || nimplintvars != 0 )
1616  {
1617  discretevar = TRUE;
1618  }
1619 
1620  /* get pointers to linear constraints handlers, so can avoid string comparisons */
1621  linearconshdlrs[0] = SCIPfindConshdlr(subproblem, "knapsack");
1622  linearconshdlrs[1] = SCIPfindConshdlr(subproblem, "linear");
1623  linearconshdlrs[2] = SCIPfindConshdlr(subproblem, "logicor");
1624  linearconshdlrs[3] = SCIPfindConshdlr(subproblem, "setppc");
1625  linearconshdlrs[4] = SCIPfindConshdlr(subproblem, "varbound");
1626 
1627  /* Get pointer to the nonlinear constraint handler if we also have an NLP solver to solve NLPs.
1628  * If there is no NLP solver, but there are (convex) nonlinear constraints, then the LP relaxation of subproblems
1629  * will (currently) not be sufficient to solve subproblems to optimality. Thus, we also take the presence of convex
1630  * nonlinear constraints as signal for having to solve the CIP eventually, thus, by abuse of notation,
1631  * return not-convex here. In summary, we do not need to have a special look onto non-linear constraints
1632  * if no NLP solver is present, and can treat them as any other constraint that is not of linear type.
1633  */
1634  if( SCIPgetNNlpis(subproblem) > 0 )
1635  {
1636  nlconshdlr = SCIPfindConshdlr(subproblem, "nonlinear");
1637  }
1638 
1639  /* if the nonlinear constraint handler exists, then we create a hashmap of variables that can be assumed to be fixed.
1640  * These variables correspond to the copies of the master variables in the subproblem
1641  */
1642  if( probnumber >= 0 && nlconshdlr != NULL )
1643  {
1644  SCIP_VAR* mappedvar;
1645 
1646  SCIP_CALL( SCIPhashmapCreate(&assumevarfixed, SCIPblkmem(set->scip), SCIPgetNVars(subproblem)) );
1647 
1648  /* finding the subproblem variables that correspond to master variables */
1649  for( i = 0; i < nvars; i++ )
1650  {
1651  /* getting the corresponding master problem variable for the given variable */
1652  SCIP_CALL( SCIPbendersGetVar(benders, set, vars[i], &mappedvar, -1) );
1653 
1654  /* if the mapped variable is not NULL, then it must be stored as a possible fixed variable */
1655  if( mappedvar != NULL )
1656  {
1657  SCIP_CALL( SCIPhashmapInsert(assumevarfixed, vars[i], vars[i]) );
1658  }
1659  }
1660  }
1661 
1662  for( i = 0; i < SCIPgetNOrigConss(subproblem); ++i )
1663  {
1664  cons = SCIPgetOrigConss(subproblem)[i];
1665  conshdlr = SCIPconsGetHdlr(cons);
1666 
1667  for( j = 0; j < NLINEARCONSHDLRS; ++j )
1668  if( conshdlr == linearconshdlrs[j] )
1669  break;
1670 
1671  /* if linear constraint, then we are good */
1672  if( j < NLINEARCONSHDLRS )
1673  {
1674 #ifdef SCIP_MOREDEBUG
1675  SCIPdebugMsg(subproblem, "subproblem <%s>: constraint <%s> is linear\n", SCIPgetProbName(subproblem), SCIPconsGetName(cons));
1676 #endif
1677  continue;
1678  }
1679 
1680  /* if cons_nonlinear (and nlconshdlr != NULL), then check whether convex */
1681  if( conshdlr == nlconshdlr )
1682  {
1683  SCIP_Bool isconvex;
1684  SCIP_EXPRCURV curv;
1685  SCIP_Bool havelhs;
1686  SCIP_Bool haverhs;
1687 
1688  isnonlinear = TRUE;
1689 
1690  havelhs = !SCIPisInfinity(subproblem, -SCIPgetLhsNonlinear(cons));
1691  haverhs = !SCIPisInfinity(subproblem, SCIPgetRhsNonlinear(cons));
1692  if( havelhs && haverhs )
1693  {
1694  isconvex = FALSE;
1695  }
1696  else
1697  {
1698  /* look at curvature stored in cons, though at this stage this will be unknown a.a. */
1699  curv = SCIPgetCurvatureNonlinear(cons);
1700  isconvex = ((!havelhs || (curv & SCIP_EXPRCURV_CONCAVE) == SCIP_EXPRCURV_CONCAVE)) &&
1701  ((!haverhs || (curv & SCIP_EXPRCURV_CONVEX) == SCIP_EXPRCURV_CONVEX));
1702 
1703  if( !isconvex )
1704  {
1705  /* if not found convex, compute curvature via nlhdlr_convex and decide again */
1706 
1707  /* make sure activities are up to date. SCIPhasExprCurvature currently assumes that this is already the case */
1708  SCIP_CALL( SCIPevalExprActivity(subproblem, SCIPgetExprNonlinear(cons)) );
1709 
1710  SCIP_CALL( SCIPhasExprCurvature(subproblem, SCIPgetExprNonlinear(cons), havelhs ? SCIP_EXPRCURV_CONCAVE : SCIP_EXPRCURV_CONVEX, &isconvex, assumevarfixed) );
1711  }
1712  }
1713 
1714  if( isconvex )
1715  {
1716 #ifdef SCIP_MOREDEBUG
1717  SCIPdebugMsg(subproblem, "subproblem <%s>: nonlinear constraint <%s> is convex\n", SCIPgetProbName(subproblem), SCIPconsGetName(cons));
1718 #endif
1719  continue;
1720  }
1721  else
1722  {
1723 #ifdef SCIP_MOREDEBUG
1724  SCIPdebugMsg(subproblem, "subproblem <%s>: nonlinear constraint <%s> not convex\n", SCIPgetProbName(subproblem), SCIPconsGetName(cons));
1725 #endif
1726  goto TERMINATE;
1727  }
1728  }
1729 
1730 #ifdef SCIP_MOREDEBUG
1731  SCIPdebugMsg(subproblem, "subproblem <%s>: potentially nonconvex constraint <%s>\n", SCIPgetProbName(subproblem), SCIPconsGetName(cons));
1732 #endif
1733  goto TERMINATE;
1734  }
1735 
1736  /* if we made it until here, then all constraints are known and convex */
1737  convexcons = TRUE;
1738 
1739 TERMINATE:
1740  /* setting the flag for the convexity of the subproblem. If convexity doesn't need to be checked, then it is assumed
1741  * that the subproblems are convex. However, if there are discrete variables, then the problem must be set as
1742  * non-convex. The discrete master variables will be changed to continuous, but this will happen at the first call to
1743  * SCIPbendersSetupSubproblem
1744  */
1745  if( probnumber >= 0 )
1746  {
1747  convexcons = convexcons || !benders->checkconsconvexity;
1748 
1749  if( convexcons && !discretevar )
1751  else if( convexcons && discretevar )
1753  else if( !convexcons && !discretevar )
1755  else if( !convexcons && discretevar )
1757  else
1758  SCIPABORT();
1759  }
1760 
1761  /* setting the non-linear subproblem flag */
1762  if( probnumber >= 0 )
1763  SCIPbendersSetSubproblemIsNonlinear(benders, probnumber, isnonlinear);
1764  else
1765  SCIPbendersSetMasterIsNonlinear(benders, isnonlinear);
1766 
1767  if( probnumber >= 0 )
1768  {
1769  SCIPsetDebugMsg(set, "subproblem <%s> has been found to be of type %d\n", SCIPgetProbName(subproblem),
1770  SCIPbendersGetSubproblemType(benders, probnumber));
1771  }
1772 
1773  /* releasing the fixed variable hashmap */
1774  if( assumevarfixed != NULL )
1775  SCIPhashmapFree(&assumevarfixed);
1776 
1777  return SCIP_OKAY;
1778 }
1779 
1780 /** creates the subproblems and registers it with the Benders' decomposition struct */
1781 static
1783  SCIP_BENDERS* benders, /**< Benders' decomposition */
1784  SCIP_SET* set /**< global SCIP settings */
1785  )
1786 {
1787  SCIP* subproblem;
1788  SCIP_EVENTHDLR* eventhdlr;
1789  SCIP_VAR* mastervar;
1790  SCIP_VAR** vars;
1791  int nvars;
1792  int nsubproblems;
1793  int i;
1794  int j;
1795 
1796  assert(benders != NULL);
1797  assert(set != NULL);
1798 
1799  /* if the subproblems have already been created, then they will not be created again. This is the case if the
1800  * transformed problem has been freed and then retransformed. The subproblems should only be created when the problem
1801  * is first transformed. */
1802  if( benders->subprobscreated )
1803  return SCIP_OKAY;
1804 
1805  nsubproblems = SCIPbendersGetNSubproblems(benders);
1806 
1807  /* creating all subproblems */
1808  for( i = 0; i < nsubproblems; i++ )
1809  {
1810  /* calling the create subproblem call back method */
1811  SCIP_CALL( benders->benderscreatesub(set->scip, benders, i) );
1812 
1813  subproblem = SCIPbendersSubproblem(benders, i);
1814 
1815  /* the subproblem SCIP instance could be set to NULL. This is because user defined subproblem solving methods
1816  * could be used that don't solve a SCIP instance. Thus, the following setup of the subproblem SCIP instance is
1817  * not required.
1818  *
1819  * NOTE: since the subproblems are supplied as NULL pointers, the internal convexity check can not be performed.
1820  * The user needs to explicitly specify the subproblem type.
1821  */
1822  if( subproblem != NULL )
1823  {
1824  /* setting global limits for the subproblems. This overwrites the limits set by the user */
1825  SCIP_CALL( SCIPsetIntParam(subproblem, "limits/maxorigsol", 0) );
1826 
1827  /* getting the number of integer and binary variables to determine the problem type */
1828  SCIP_CALL( SCIPgetVarsData(subproblem, &vars, &nvars, NULL, NULL, NULL, NULL) );
1829 
1830  /* The objective function coefficients of the master problem are set to zero. This is necessary for the Benders'
1831  * decomposition algorithm, since the cut methods and the objective function check assumes that the objective
1832  * coefficients of the master problem variables are zero.
1833  *
1834  * This only occurs if the Benders' decomposition is not a copy. It is assumed that the correct objective
1835  * coefficients are given during the first subproblem creation.
1836  *
1837  * If the subproblems were copied, then the master variables will be checked to ensure that they have a zero
1838  * objective value.
1839  */
1840  if( !benders->iscopy || benders->threadsafe )
1841  {
1842  SCIP_Bool objchanged = FALSE;
1843 
1844  assert(SCIPgetStage(subproblem) == SCIP_STAGE_PROBLEM);
1845  for( j = 0; j < nvars; j++ )
1846  {
1847  /* retrieving the master problem variable */
1848  SCIP_CALL( SCIPbendersGetVar(benders, set, vars[j], &mastervar, -1) );
1849 
1850  /* if mastervar is not NULL, then the subproblem variable has a corresponding master problem variable */
1851  if( mastervar != NULL && SCIPvarGetObj(vars[j]) != 0.0 )
1852  {
1853  SCIPverbMessage(subproblem, SCIP_VERBLEVEL_FULL, NULL, "Benders' decomposition: Changing the objective "
1854  "coefficient of copy of master problem variable <%s> in subproblem %d to zero.\n",
1855  SCIPvarGetName(mastervar), i);
1856  /* changing the subproblem variable objective coefficient to zero */
1857  SCIP_CALL( SCIPchgVarObj(subproblem, vars[j], 0.0) );
1858 
1859  objchanged = TRUE;
1860  }
1861  }
1862 
1863  if( objchanged )
1864  {
1865  SCIPverbMessage(subproblem, SCIP_VERBLEVEL_HIGH, NULL, "Benders' decomposition: Objective coefficients of "
1866  "copied of master problem variables has been changed to zero.\n");
1867  }
1868  }
1869 
1870  /* changing all of the master problem variable to continuous. */
1871  SCIP_CALL( SCIPbendersChgMastervarsToCont(benders, set, i) );
1872 
1873  /* checking the convexity of the subproblem. The convexity of the subproblem indicates whether the convex
1874  * relaxation is a valid relaxation for the problem
1875  */
1876  SCIP_CALL( checkSubproblemConvexity(benders, set, i) );
1877 
1878  /* if the problem is convex and has nonlinear constraints, then slack variables must be added to each of the
1879  * constraints
1880  */
1881  if( benders->execfeasphase ||
1883  && SCIPbendersSubproblemIsNonlinear(benders, i)) )
1884  {
1885  /* the slack variables are only added to the subproblem once. If the initialisation methods are called from a
1886  * copy, then the slack variables are not re-added. Alternatively, if the copy must be threadsafe, then the
1887  * subproblems are created from scratch again, so the slack variables need to be added.
1888  */
1889  if( !benders->iscopy || benders->threadsafe )
1890  {
1891  SCIP_CALL( addSlackVarsToConstraints(benders, set, i) );
1892  }
1893 
1894  /* setting the flag to indicate that slack variables have been added to the subproblem constraints. This is only
1895  * set if the slack variables have been added at the request of the user.
1896  */
1897  if( benders->execfeasphase )
1898  benders->feasibilityphase = TRUE;
1899  }
1900 
1901  /* after checking the subproblem for convexity, if the subproblem has convex constraints and continuous variables,
1902  * then the problem is entered into probing mode. Otherwise, it is initialised as a CIP
1903  */
1905  {
1906  /* if the user has not implemented a solve subproblem callback, then the subproblem solves are performed
1907  * internally. To be more efficient the subproblem is put into probing mode. */
1908  if( benders->benderssolvesubconvex == NULL && benders->benderssolvesub == NULL
1909  && SCIPgetStage(subproblem) <= SCIP_STAGE_PROBLEM )
1910  {
1911  SCIP_CALL( initialiseLPSubproblem(benders, set, i) );
1912  }
1913  }
1914  else
1915  {
1916  SCIP_EVENTHDLRDATA* eventhdlrdata_mipnodefocus;
1917  SCIP_EVENTHDLRDATA* eventhdlrdata_upperbound;
1918 
1919  /* because the subproblems could be reused in the copy, the event handler is not created again. If the
1920  * threadsafe is TRUE, then it is assumed that the subproblems are not reused.
1921  * NOTE: This currently works with the benders_default implementation. It may not be very general. */
1922  if( benders->benderssolvesubconvex == NULL && benders->benderssolvesub == NULL
1923  && (!benders->iscopy || benders->threadsafe) )
1924  {
1925  SCIP_CALL( SCIPallocBlockMemory(subproblem, &eventhdlrdata_mipnodefocus) );
1926  SCIP_CALL( SCIPallocBlockMemory(subproblem, &eventhdlrdata_upperbound) );
1927 
1928  SCIP_CALL( initEventhandlerData(subproblem, eventhdlrdata_mipnodefocus) );
1929  SCIP_CALL( initEventhandlerData(subproblem, eventhdlrdata_upperbound) );
1930 
1931  /* include the first LP solved event handler into the subproblem */
1933  MIPNODEFOCUS_EVENTHDLR_DESC, eventExecBendersMipnodefocus, eventhdlrdata_mipnodefocus) );
1934  SCIP_CALL( SCIPsetEventhdlrInitsol(subproblem, eventhdlr, eventInitsolBendersMipnodefocus) );
1935  SCIP_CALL( SCIPsetEventhdlrExitsol(subproblem, eventhdlr, eventExitsolBendersMipnodefocus) );
1936  SCIP_CALL( SCIPsetEventhdlrExit(subproblem, eventhdlr, eventExitBendersMipnodefocus) );
1937  SCIP_CALL( SCIPsetEventhdlrFree(subproblem, eventhdlr, eventFreeBendersMipnodefocus) );
1938  assert(eventhdlr != NULL);
1939 
1940  /* include the upper bound interrupt event handler into the subproblem */
1942  UPPERBOUND_EVENTHDLR_DESC, eventExecBendersUpperbound, eventhdlrdata_upperbound) );
1943  SCIP_CALL( SCIPsetEventhdlrInitsol(subproblem, eventhdlr, eventInitsolBendersUpperbound) );
1944  SCIP_CALL( SCIPsetEventhdlrExitsol(subproblem, eventhdlr, eventExitsolBendersUpperbound) );
1945  SCIP_CALL( SCIPsetEventhdlrExit(subproblem, eventhdlr, eventExitBendersUpperbound) );
1946  SCIP_CALL( SCIPsetEventhdlrFree(subproblem, eventhdlr, eventFreeBendersUpperbound) );
1947  assert(eventhdlr != NULL);
1948  }
1949  }
1950  }
1951  else
1952  {
1953  /* a user must specify the subproblem type if they are not supplying a SCIP instance. */
1955  {
1956  SCIPerrorMessage("If the subproblem is set to NULL, then the subproblem type must be specified.\n");
1957  SCIPerrorMessage("In the subproblem creation callback, call SCIPbendersSetSubproblemType with the appropriate problem type.\n");
1958 
1959  return SCIP_ERROR;
1960  }
1961  }
1962  }
1963 
1964  /* checking the convexity of the master problem. This information is useful for the cut generation methods, such as
1965  * non-good and integer cuts
1966  */
1967  SCIP_CALL( checkSubproblemConvexity(benders, set, -1) );
1968 
1969  benders->subprobscreated = TRUE;
1970 
1971  return SCIP_OKAY;
1972 }
1973 
1974 
1975 /** initializes Benders' decomposition */
1977  SCIP_BENDERS* benders, /**< Benders' decomposition */
1978  SCIP_SET* set /**< global SCIP settings */
1979  )
1980 {
1981  int i;
1982 
1983  assert(benders != NULL);
1984  assert(set != NULL);
1985 
1986  if( benders->initialized )
1987  {
1988  SCIPerrorMessage("Benders' decomposition <%s> already initialized\n", benders->name);
1989  return SCIP_INVALIDCALL;
1990  }
1991 
1992  if( set->misc_resetstat )
1993  {
1994  SCIPclockReset(benders->setuptime);
1995  SCIPclockReset(benders->bendersclock);
1996 
1997  benders->ncalls = 0;
1998  benders->ncutsfound = 0;
1999  benders->ntransferred = 0;
2000  }
2001 
2002  /* start timing */
2003  SCIPclockStart(benders->setuptime, set);
2004 
2005  if( benders->bendersinit != NULL )
2006  {
2007  SCIP_CALL( benders->bendersinit(set->scip, benders) );
2008  }
2009 
2010  benders->initialized = TRUE;
2011 
2012  /* if the Benders' decomposition is a copy, then the auxiliary variables already exist. So they are registered with
2013  * the Benders' decomposition struct during the init stage. If the Benders' decomposition is not a copy, then the
2014  * auxiliary variables need to be created, which occurs in the initpre stage
2015  */
2016  if( benders->iscopy )
2017  {
2018  /* the copied auxiliary variables must be assigned to the target Benders' decomposition */
2019  SCIP_CALL( assignAuxiliaryVariables(set->scip, benders) );
2020  }
2021 
2022  /* creates the subproblems and sets up the probing mode for LP subproblems. This function calls the benderscreatesub
2023  * callback. */
2024  SCIP_CALL( createSubproblems(benders, set) );
2025 
2026  /* storing the solution tolerance set by the SCIP parameters */
2027  SCIP_CALL( SCIPsetGetRealParam(set, "benders/solutiontol", &benders->solutiontol) );
2028 
2029  /* allocating memory for the stored constraints array */
2030  if( benders->storedcutssize == 0 )
2031  {
2033  benders->storedcutssize = BENDERS_ARRAYSIZE;
2034  benders->nstoredcuts = 0;
2035  }
2036 
2037  /* initialising the Benders' cuts */
2038  SCIPbendersSortBenderscuts(benders);
2039  for( i = 0; i < benders->nbenderscuts; i++ )
2040  {
2041  SCIP_CALL( SCIPbenderscutInit(benders->benderscuts[i], set) );
2042  }
2043 
2044  /* stop timing */
2045  SCIPclockStop(benders->setuptime, set);
2046 
2047  return SCIP_OKAY;
2048 }
2049 
2050 
2051 /** Transfers Benders' cuts that were generated while solving a sub-SCIP to the original SCIP instance. This involves
2052  * creating a constraint/cut that is equivalent to the generated cut in the sub-SCIP. This new constraint/cut is then
2053  * added to the original SCIP instance.
2054  */
2055 static
2057  SCIP* sourcescip, /**< the source SCIP from when the Benders' decomposition was copied */
2058  SCIP_BENDERS* benders, /**< the Benders' decomposition structure of the sub SCIP */
2059  SCIP_VAR** vars, /**< the variables from the source constraint */
2060  SCIP_Real* vals, /**< the coefficients of the variables in the source constriant */
2061  SCIP_Real lhs, /**< the LHS of the source constraint */
2062  SCIP_Real rhs, /**< the RHS of the source constraint */
2063  int nvars /**< the number of variables in the source constraint */
2064  )
2065 {
2066  SCIP_BENDERS* sourcebenders; /* the Benders' decomposition of the source SCIP */
2067  SCIP_CONSHDLR* consbenders; /* a helper variable for the Benders' decomposition constraint handler */
2068  SCIP_CONS* transfercons = NULL; /* the constraint that is generated to transfer the constraints/cuts */
2069  SCIP_ROW* transfercut = NULL; /* the cut that is generated to transfer the constraints/cuts */
2070  SCIP_VAR* sourcevar; /* the source variable that will be added to the transferred cut */
2071  SCIP_VAR* origvar;
2072  SCIP_Real scalar;
2073  SCIP_Real constant;
2074  char cutname[SCIP_MAXSTRLEN]; /* the name of the transferred cut */
2075  int i;
2076  SCIP_Bool fail;
2077 
2078  assert(sourcescip != NULL);
2079  assert(benders != NULL);
2080  assert(vars != NULL);
2081  assert(vals != NULL);
2082 
2083  /* retrieving the source Benders' decomposition structure */
2084  sourcebenders = SCIPfindBenders(sourcescip, SCIPbendersGetName(benders));
2085 
2086  /* retrieving the Benders' decomposition constraint handler */
2087  consbenders = SCIPfindConshdlr(sourcescip, "benders");
2088 
2089  /* setting the name of the transferred cut */
2090  (void) SCIPsnprintf(cutname, SCIP_MAXSTRLEN, "transferredcut_%d",
2091  SCIPbendersGetNTransferredCuts(sourcebenders) );
2092 
2093  /* TODO: It could be more efficient to pass an updated vars array with the vals array to the
2094  * SCIPcreateConsBasicLinear/SCIPcreateEmptyRowConshdlr. This should be implemented to improve the performance of the
2095  * Large Neighbourhood Benders Search.
2096  */
2097 
2098  /* creating an empty row/constraint for the transferred cut */
2099  if( sourcebenders->cutsasconss )
2100  {
2101  SCIP_CALL( SCIPcreateConsBasicLinear(sourcescip, &transfercons, cutname, 0, NULL, NULL, lhs, rhs) );
2102  SCIP_CALL( SCIPsetConsRemovable(sourcescip, transfercons, TRUE) );
2103  }
2104  else
2105  {
2106  SCIP_CALL( SCIPcreateEmptyRowConshdlr(sourcescip, &transfercut, consbenders, cutname, lhs, rhs, FALSE,
2107  FALSE, TRUE) );
2108  }
2109 
2110  fail = FALSE;
2111  for( i = 0; i < nvars; i++ )
2112  {
2113  /* getting the original variable for the transformed variable */
2114  origvar = vars[i];
2115  scalar = 1.0;
2116  constant = 0.0;
2117  SCIP_CALL( SCIPvarGetOrigvarSum(&origvar, &scalar, &constant) );
2118 
2119  /* getting the source var from the hash map */
2120  sourcevar = (SCIP_VAR*) SCIPhashmapGetImage(benders->mastervarsmap, origvar);
2121 
2122  /* if the source variable is not found, then the mapping in incomplete. So the constraint can not be
2123  * transferred. */
2124  if( sourcevar == NULL )
2125  {
2126  fail = TRUE;
2127  break;
2128  }
2129 
2130  if( sourcebenders->cutsasconss )
2131  {
2132  assert( transfercons != NULL );
2133  SCIP_CALL( SCIPaddCoefLinear(sourcescip, transfercons, sourcevar, vals[i]) ); /*lint !e644*/
2134  }
2135  else
2136  {
2137  assert( transfercut != NULL );
2138  SCIP_CALL( SCIPaddVarToRow(sourcescip, transfercut, sourcevar, vals[i]) ); /*lint !e644*/
2139  }
2140  }
2141 
2142  /* if all of the source variables were found to generate the cut */
2143  if( !fail )
2144  {
2145  if( sourcebenders->cutsasconss )
2146  {
2147  SCIP_CALL( SCIPaddCons(sourcescip, transfercons) );
2148  }
2149  else
2150  {
2151  SCIP_CALL( SCIPaddPoolCut(sourcescip, transfercut) );
2152  }
2153 
2154  sourcebenders->ntransferred++;
2155  }
2156 
2157  /* release the row/constraint */
2158  if( sourcebenders->cutsasconss )
2159  {
2160  /* only release if the creation of the constraint failed. */
2161  SCIP_CALL( SCIPreleaseCons(sourcescip, &transfercons) );
2162  }
2163  else
2164  {
2165  SCIP_CALL( SCIPreleaseRow(sourcescip, &transfercut) );
2166  }
2167 
2168  return SCIP_OKAY;
2169 }
2170 
2171 
2172 /** transfers the cuts generated in a subscip to the source scip */
2173 static
2175  SCIP* sourcescip, /**< the source SCIP from when the Benders' decomposition was copied */
2176  SCIP* subscip, /**< the sub SCIP where the Benders' cuts were generated */
2177  SCIP_BENDERS* benders /**< the Benders' decomposition structure of the sub SCIP */
2178  )
2179 {
2180  SCIP_BENDERS* sourcebenders; /* the Benders' decomposition of the source SCIP */
2181  SCIP_VAR** vars; /* the variables of the added constraint/row */
2182  SCIP_Real* vals; /* the values of the added constraint/row */
2183  SCIP_Real lhs; /* the LHS of the added constraint/row */
2184  SCIP_Real rhs; /* the RHS of the added constraint/row */
2185  int naddedcuts;
2186  int nvars;
2187  int i;
2188 
2189  assert(subscip != NULL);
2190  assert(benders != NULL);
2191 
2192  /* retrieving the source Benders' decomposition structure */
2193  sourcebenders = SCIPfindBenders(sourcescip, SCIPbendersGetName(benders));
2194 
2195  /* exit if the cuts should not be transferred from the sub SCIP to the source SCIP. */
2196  if( !sourcebenders->transfercuts || benders->mastervarsmap == NULL )
2197  return SCIP_OKAY;
2198 
2199  /* retrieving the number of stored Benders' cuts */
2200  naddedcuts = SCIPbendersGetNStoredCuts(benders);
2201 
2202  /* looping over all added cuts to construct the cut for the source scip */
2203  for( i = 0; i < naddedcuts; i++ )
2204  {
2205  /* collecting the variable information from the constraint */
2206  SCIP_CALL( SCIPbendersGetStoredCutData(benders, i, &vars, &vals, &lhs, &rhs, &nvars) );
2207 
2208  if( nvars > 0 )
2209  {
2210  /* create and add the cut to be transferred from the sub SCIP to the source SCIP */
2211  SCIP_CALL( createAndAddTransferredCut(sourcescip, benders, vars, vals, lhs, rhs, nvars) );
2212  }
2213  }
2214 
2215  return SCIP_OKAY;
2216 }
2217 
2218 
2219 /** calls exit method of Benders' decomposition */
2221  SCIP_BENDERS* benders, /**< Benders' decomposition */
2222  SCIP_SET* set /**< global SCIP settings */
2223  )
2224 {
2225  int nsubproblems;
2226  int i;
2227 
2228  assert(benders != NULL);
2229  assert(set != NULL);
2230 
2231  if( !benders->initialized )
2232  {
2233  SCIPerrorMessage("Benders' decomposition <%s> not initialized\n", benders->name);
2234  return SCIP_INVALIDCALL;
2235  }
2236 
2237  /* start timing */
2238  SCIPclockStart(benders->setuptime, set);
2239 
2240  if( benders->bendersexit != NULL )
2241  {
2242  SCIP_CALL( benders->bendersexit(set->scip, benders) );
2243  }
2244 
2245  /* if the Benders' decomposition is a copy, then is a variable mapping was provided, then the generated cuts will
2246  * be transferred to the source scip
2247  */
2248  if( benders->iscopy && benders->mastervarsmap != NULL )
2249  {
2250  SCIP_CALL( transferBendersCuts(benders->sourcescip, set->scip, benders) );
2251  }
2252 
2253  /* releasing the stored constraints */
2254  for( i = benders->nstoredcuts - 1; i >= 0; i-- )
2255  {
2256  SCIPfreeBlockMemoryArray(set->scip, &benders->storedcuts[i]->vals, benders->storedcuts[i]->nvars);
2257  SCIPfreeBlockMemoryArray(set->scip, &benders->storedcuts[i]->vars, benders->storedcuts[i]->nvars);
2258  SCIPfreeBlockMemory(set->scip, &benders->storedcuts[i]); /*lint !e866*/
2259  }
2260 
2261  BMSfreeBlockMemoryArray(SCIPblkmem(set->scip), &benders->storedcuts, benders->storedcutssize);
2262  benders->storedcutssize = 0;
2263  benders->nstoredcuts = 0;
2264 
2265  /* releasing all of the auxiliary variables */
2266  nsubproblems = SCIPbendersGetNSubproblems(benders);
2267  for( i = 0; i < nsubproblems; i++ )
2268  {
2269  /* it is possible that the master problem is not solved. As such, the auxiliary variables will not be created. So
2270  * we don't need to release the variables
2271  */
2272  if( benders->auxiliaryvars[i] != NULL )
2273  {
2274  /* we need to remove the locks from the auxiliary variables. This will be called always for the highest priority
2275  * Benders' plugin and others if the auxiliary variables are not shared
2276  */
2277  if( !benders->iscopy && SCIPvarGetNLocksDown(benders->auxiliaryvars[i]) > 0 )
2278  SCIP_CALL( SCIPaddVarLocksType(set->scip, benders->auxiliaryvars[i], SCIP_LOCKTYPE_MODEL, -1, 0) );
2279 
2280  SCIP_CALL( SCIPreleaseVar(set->scip, &benders->auxiliaryvars[i]) );
2281  }
2282  }
2283 
2284  /* if a corepoint has been used for cut strengthening, then this needs to be freed */
2285  if( benders->corepoint != NULL )
2286  {
2287  SCIP_CALL( SCIPfreeSol(set->scip, &benders->corepoint) );
2288  }
2289 
2290  /* calling the exit method for the Benders' cuts */
2291  SCIPbendersSortBenderscuts(benders);
2292  for( i = 0; i < benders->nbenderscuts; i++ )
2293  {
2294  SCIP_CALL( SCIPbenderscutExit(benders->benderscuts[i], set) );
2295  }
2296 
2297  benders->initialized = FALSE;
2298 
2299  /* stop timing */
2300  SCIPclockStop(benders->setuptime, set);
2301 
2302  return SCIP_OKAY;
2303 }
2304 
2305 /** Checks whether a subproblem is independent. */
2306 static
2308  SCIP* scip, /**< the SCIP data structure */
2309  SCIP_BENDERS* benders /**< Benders' decomposition */
2310  )
2311 {
2312  SCIP_VAR** vars;
2313  int nvars;
2314  int nsubproblems;
2315  int i;
2316  int j;
2317 
2318  assert(scip != NULL);
2319  assert(benders != NULL);
2320 
2321  /* retrieving the master problem variables */
2322  SCIP_CALL( SCIPgetVarsData(scip, &vars, &nvars, NULL, NULL, NULL, NULL) );
2323 
2324  nsubproblems = SCIPbendersGetNSubproblems(benders);
2325 
2326  /* looping over all subproblems to check whether there exists at least one master problem variable */
2327  for( i = 0; i < nsubproblems; i++ )
2328  {
2329  SCIP_Bool independent = FALSE;
2330 
2331  /* if there are user defined solving or freeing functions, then it is not possible to declare the independence of
2332  * the subproblems.
2333  */
2334  if( benders->benderssolvesubconvex == NULL && benders->benderssolvesub == NULL
2335  && benders->bendersfreesub == NULL )
2336  {
2337  independent = TRUE;
2338 
2339  for( j = 0; j < nvars; j++ )
2340  {
2341  SCIP_VAR* subprobvar;
2342 
2343  /* getting the subproblem problem variable corresponding to the master problem variable */
2344  SCIP_CALL( SCIPgetBendersSubproblemVar(scip, benders, vars[j], &subprobvar, i) );
2345 
2346  /* if the subporblem variable is not NULL, then the subproblem depends on the master problem */
2347  if( subprobvar != NULL )
2348  {
2349  independent = FALSE;
2350  break;
2351  }
2352  }
2353 
2354  /* setting the independent flag */
2355  SCIPbendersSetSubproblemIsIndependent(benders, i, independent);
2356  }
2357  }
2358 
2359  return SCIP_OKAY;
2360 }
2361 
2362 /** informs the Benders' decomposition that the presolving process is being started */
2364  SCIP_BENDERS* benders, /**< Benders' decomposition */
2365  SCIP_SET* set, /**< global SCIP settings */
2366  SCIP_STAT* stat /**< dynamic problem statistics */
2367  )
2368 {
2369  assert(benders != NULL);
2370  assert(set != NULL);
2371  assert(stat != NULL);
2372 
2373  /* if the Benders' decomposition is the original, then the auxiliary variables need to be created. If the Benders'
2374  * decomposition is a copy, then the auxiliary variables already exist. The assignment of the auxiliary variables
2375  * occurs in bendersInit
2376  */
2377  if( !benders->iscopy )
2378  {
2379  /* check the subproblem independence. This check is only performed if the user has not implemented a solve
2380  * subproblem function.
2381  */
2382  if( benders->benderssolvesubconvex == NULL && benders->benderssolvesub == NULL )
2383  SCIP_CALL( checkSubproblemIndependence(set->scip, benders) );
2384 
2385  /* adding the auxiliary variables to the master problem */
2386  SCIP_CALL( addAuxiliaryVariablesToMaster(set->scip, benders) );
2387  }
2388 
2389  /* call presolving initialization method of Benders' decomposition */
2390  if( benders->bendersinitpre != NULL )
2391  {
2392  /* start timing */
2393  SCIPclockStart(benders->setuptime, set);
2394 
2395  SCIP_CALL( benders->bendersinitpre(set->scip, benders) );
2396 
2397  /* stop timing */
2398  SCIPclockStop(benders->setuptime, set);
2399  }
2400 
2401  return SCIP_OKAY;
2402 }
2403 
2404 
2405 /** informs the Benders' decomposition that the presolving process has completed */
2407  SCIP_BENDERS* benders, /**< Benders' decomposition */
2408  SCIP_SET* set, /**< global SCIP settings */
2409  SCIP_STAT* stat /**< dynamic problem statistics */
2410  )
2411 {
2412  assert(benders != NULL);
2413  assert(set != NULL);
2414  assert(stat != NULL);
2415 
2416  /* call presolving deinitialization method of Benders' decomposition */
2417  if( benders->bendersexitpre != NULL )
2418  {
2419  /* start timing */
2420  SCIPclockStart(benders->setuptime, set);
2421 
2422  SCIP_CALL( benders->bendersexitpre(set->scip, benders) );
2423 
2424  /* stop timing */
2425  SCIPclockStop(benders->setuptime, set);
2426  }
2427 
2428  return SCIP_OKAY;
2429 }
2430 
2431 /** informs Benders' decomposition that the branch and bound process is being started */
2433  SCIP_BENDERS* benders, /**< Benders' decomposition */
2434  SCIP_SET* set /**< global SCIP settings */
2435  )
2436 {
2437  int i;
2438 
2439  assert(benders != NULL);
2440  assert(set != NULL);
2441 
2442  /* call solving process initialization method of Benders' decomposition */
2443  if( benders->bendersinitsol != NULL )
2444  {
2445  /* start timing */
2446  SCIPclockStart(benders->setuptime, set);
2447 
2448  SCIP_CALL( benders->bendersinitsol(set->scip, benders) );
2449 
2450  /* stop timing */
2451  SCIPclockStop(benders->setuptime, set);
2452  }
2453 
2454  /* calling the initsol method for the Benders' cuts */
2455  SCIPbendersSortBenderscuts(benders);
2456  for( i = 0; i < benders->nbenderscuts; i++ )
2457  {
2458  SCIP_CALL( SCIPbenderscutInitsol(benders->benderscuts[i], set) );
2459  }
2460 
2461  return SCIP_OKAY;
2462 }
2463 
2464 /** informs Benders' decomposition that the branch and bound process data is being freed */
2466  SCIP_BENDERS* benders, /**< Benders' decomposition */
2467  SCIP_SET* set /**< global SCIP settings */
2468  )
2469 {
2470  int nsubproblems;
2471  int i;
2472 
2473  assert(benders != NULL);
2474  assert(set != NULL);
2475 
2476  nsubproblems = SCIPbendersGetNSubproblems(benders);
2477  /* freeing all subproblems that are independent, this is because they have not bee freed during the subproblem
2478  * solving loop.
2479  */
2480  for( i = 0; i < nsubproblems; i++ )
2481  {
2482  if( SCIPbendersSubproblemIsIndependent(benders, i) )
2483  {
2484  /* disabling the independence of the subproblem so that it can be freed */
2486 
2487  /* freeing the independent subproblem */
2488  SCIP_CALL( SCIPbendersFreeSubproblem(benders, set, i) );
2489  }
2490  }
2491 
2492  /* call solving process deinitialization method of Benders' decomposition */
2493  if( benders->bendersexitsol != NULL )
2494  {
2495  /* start timing */
2496  SCIPclockStart(benders->setuptime, set);
2497 
2498  SCIP_CALL( benders->bendersexitsol(set->scip, benders) );
2499 
2500  /* stop timing */
2501  SCIPclockStop(benders->setuptime, set);
2502  }
2503 
2504  /* sorting the Benders' decomposition cuts in order of priority. Only a single cut is generated for each subproblem
2505  * per solving iteration. This is particularly important in the case of the optimality and feasibility cuts. Since
2506  * these work on two different solutions to the subproblem, it is not necessary to generate both cuts. So, once the
2507  * feasibility cut is generated, then no other cuts will be generated.
2508  */
2509  SCIPbendersSortBenderscuts(benders);
2510 
2511  /* calling the exitsol method for the Benders' cuts */
2512  for( i = 0; i < benders->nbenderscuts; i++ )
2513  {
2514  SCIP_CALL( SCIPbenderscutExitsol(benders->benderscuts[i], set) );
2515  }
2516 
2517  return SCIP_OKAY;
2518 }
2519 
2520 /** activates Benders' decomposition such that it is called in LP solving loop */
2522  SCIP_BENDERS* benders, /**< the Benders' decomposition structure */
2523  SCIP_SET* set, /**< global SCIP settings */
2524  int nsubproblems /**< the number subproblems used in this decomposition */
2525  )
2526 {
2527  SCIP_EVENTHDLR* eventhdlr;
2528  SCIP_EVENTHDLRDATA* eventhdlrdata;
2529  int i;
2530 
2531  assert(benders != NULL);
2532  assert(set != NULL);
2533  assert(set->stage == SCIP_STAGE_INIT || set->stage == SCIP_STAGE_PROBLEM);
2534 
2535  if( !benders->active )
2536  {
2537  benders->active = TRUE;
2538  set->nactivebenders++;
2539  set->benderssorted = FALSE;
2540 
2541  benders->nsubproblems = nsubproblems;
2542  benders->nactivesubprobs = nsubproblems;
2543  benders->prevlowerbound = -SCIPsetInfinity(set);
2544  benders->strengthenround = FALSE;
2545 
2546  /* allocating memory for the subproblems arrays */
2547  SCIP_ALLOC( BMSallocMemoryArray(&benders->subproblems, benders->nsubproblems) );
2548  SCIP_ALLOC( BMSallocMemoryArray(&benders->auxiliaryvars, benders->nsubproblems) );
2549  SCIP_ALLOC( BMSallocMemoryArray(&benders->solvestat, benders->nsubproblems) );
2550  SCIP_ALLOC( BMSallocMemoryArray(&benders->subprobobjval, benders->nsubproblems) );
2553  SCIP_ALLOC( BMSallocMemoryArray(&benders->subprobtype, benders->nsubproblems) );
2556  SCIP_ALLOC( BMSallocMemoryArray(&benders->subprobsetup, benders->nsubproblems) );
2557  SCIP_ALLOC( BMSallocMemoryArray(&benders->indepsubprob, benders->nsubproblems) );
2560 
2561  /* creating the priority queue for the subproblem solving status */
2562  SCIP_CALL( SCIPpqueueCreate(&benders->subprobqueue, benders->nsubproblems, 1.1,
2563  benders->benderssubcomp == NULL ? benderssubcompdefault : benders->benderssubcomp, NULL) );
2564 
2565  for( i = 0; i < benders->nsubproblems; i++ )
2566  {
2567  SCIP_SUBPROBLEMSOLVESTAT* solvestat;
2568 
2569  benders->subproblems[i] = NULL;
2570  benders->auxiliaryvars[i] = NULL;
2571  benders->subprobobjval[i] = SCIPsetInfinity(set);
2572  benders->bestsubprobobjval[i] = SCIPsetInfinity(set);
2573  benders->subproblowerbound[i] = -SCIPsetInfinity(set);
2575  benders->subprobisconvex[i] = FALSE;
2576  benders->subprobisnonlinear[i] = FALSE;
2577  benders->subprobsetup[i] = FALSE;
2578  benders->indepsubprob[i] = FALSE;
2579  benders->subprobenabled[i] = TRUE;
2580  benders->mastervarscont[i] = FALSE;
2581 
2582  /* initialising the subproblem solving status */
2583  SCIP_ALLOC( BMSallocMemory(&solvestat) );
2584  solvestat->idx = i;
2585  solvestat->ncalls = 0;
2586  solvestat->avgiter = 0;
2587  benders->solvestat[i] = solvestat;
2588 
2589  /* inserting the initial elements into the priority queue */
2590  SCIP_CALL( SCIPpqueueInsert(benders->subprobqueue, benders->solvestat[i]) );
2591  }
2592 
2593  /* adding an eventhandler for updating the lower bound when the root node is solved. */
2594  eventhdlrdata = (SCIP_EVENTHDLRDATA*)benders;
2595 
2596  /* include event handler into SCIP */
2598  eventExecBendersNodesolved, eventhdlrdata) );
2599  SCIP_CALL( SCIPsetEventhdlrInitsol(set->scip, eventhdlr, eventInitsolBendersNodesolved) );
2600  assert(eventhdlr != NULL);
2601  }
2602 
2603  return SCIP_OKAY;
2604 }
2605 
2606 /** deactivates Benders' decomposition such that it is no longer called in LP solving loop */
2608  SCIP_BENDERS* benders, /**< the Benders' decomposition structure */
2609  SCIP_SET* set /**< global SCIP settings */
2610  )
2611 {
2612  int i;
2613 
2614  assert(benders != NULL);
2615  assert(set != NULL);
2616  assert(set->stage == SCIP_STAGE_INIT || set->stage == SCIP_STAGE_PROBLEM);
2617 
2618  if( benders->active )
2619  {
2620  int nsubproblems;
2621 
2622  nsubproblems = SCIPbendersGetNSubproblems(benders);
2623 
2624 #ifndef NDEBUG
2625  /* checking whether the auxiliary variables and subproblems are all NULL */
2626  for( i = 0; i < nsubproblems; i++ )
2627  assert(benders->auxiliaryvars[i] == NULL);
2628 #endif
2629 
2630  /* if the subproblems were created by the Benders' decomposition core, then they need to be freed */
2631  if( benders->freesubprobs )
2632  {
2633  for( i = SCIPbendersGetNSubproblems(benders) - 1; i >= 0; i-- )
2634  {
2635  SCIP* subproblem = SCIPbendersSubproblem(benders, i);
2636  SCIP_CALL( SCIPfree(&subproblem) );
2637  }
2638  }
2639 
2640  benders->active = FALSE;
2641  set->nactivebenders--;
2642  set->benderssorted = FALSE;
2643 
2644  /* freeing the priority queue memory */
2645  SCIPpqueueFree(&benders->subprobqueue);
2646 
2647  for( i = nsubproblems - 1; i >= 0; i-- )
2648  BMSfreeMemory(&benders->solvestat[i]);
2649 
2650  /* freeing the memory allocated during the activation of the Benders' decomposition */
2653  BMSfreeMemoryArray(&benders->indepsubprob);
2654  BMSfreeMemoryArray(&benders->subprobsetup);
2657  BMSfreeMemoryArray(&benders->subprobtype);
2660  BMSfreeMemoryArray(&benders->subprobobjval);
2661  BMSfreeMemoryArray(&benders->auxiliaryvars);
2662  BMSfreeMemoryArray(&benders->solvestat);
2663  BMSfreeMemoryArray(&benders->subproblems);
2664  }
2665 
2666  return SCIP_OKAY;
2667 }
2668 
2669 /** returns whether the given Benders' decomposition is in use in the current problem */
2671  SCIP_BENDERS* benders /**< the Benders' decomposition structure */
2672  )
2673 {
2674  assert(benders != NULL);
2675 
2676  return benders->active;
2677 }
2678 
2679 /** updates the lower bound for all auxiliary variables. This is called if the first LP enforced is unbounded. */
2680 static
2682  SCIP_BENDERS* benders, /**< Benders' decomposition */
2683  SCIP_SET* set, /**< global SCIP settings */
2684  SCIP_RESULT* result /**< the result from updating the auxiliary variable lower bound */
2685  )
2686 {
2687  int nsubproblems;
2688  int i;
2689 
2690  assert(benders != NULL);
2691  assert(set != NULL);
2692 
2693  (*result) = SCIP_DIDNOTRUN;
2694 
2695  nsubproblems = SCIPbendersGetNSubproblems(benders);
2696 
2697  for( i = 0; i < nsubproblems; i++ )
2698  {
2699  SCIP_VAR* auxiliaryvar;
2700  SCIP_Real lowerbound;
2701  SCIP_Bool infeasible;
2702 
2703  infeasible = FALSE;
2704 
2705  /* computing the lower bound of the subproblem by solving it without any variable fixings */
2706  SCIP_CALL( SCIPbendersComputeSubproblemLowerbound(benders, set, i, &lowerbound, &infeasible) );
2707 
2708  /* if the subproblem is infeasible, then the original problem is infeasible */
2709  if( infeasible )
2710  {
2711  (*result) = SCIP_INFEASIBLE;
2712  break;
2713  }
2714 
2715  /* retrieving the auxiliary variable */
2716  auxiliaryvar = SCIPbendersGetAuxiliaryVar(benders, i);
2717 
2718  /* only update the lower bound if it is greater than the current lower bound */
2719  if( SCIPsetIsGT(set, lowerbound, SCIPvarGetLbGlobal(auxiliaryvar)) )
2720  {
2721  SCIPsetDebugMsg(set, "Tightened lower bound of <%s> to %g\n", SCIPvarGetName(auxiliaryvar), lowerbound);
2722  /* updating the lower bound of the auxiliary variable */
2723  SCIP_CALL( SCIPchgVarLb(set->scip, auxiliaryvar, lowerbound) );
2724  (*result) = SCIP_REDUCEDDOM;
2725  }
2726 
2727  /* stores the lower bound for the subproblem */
2728  SCIPbendersUpdateSubproblemLowerbound(benders, i, lowerbound);
2729  }
2730 
2731  return SCIP_OKAY;
2732 }
2733 
2734 /** sets the core point used for cut strengthening. If the strenghtenintpoint is set to 'i', then the core point is
2735  * reinitialised each time the incumbent is updated
2736  */
2737 static
2739  SCIP* scip, /**< the SCIP data structure */
2740  SCIP_BENDERS* benders /**< Benders' decomposition */
2741  )
2742 {
2743  SCIP_SOL* bestsol;
2744 
2745  assert(scip != NULL);
2746  assert(benders != NULL);
2747 
2748  /* if the core point is not NULL and the interior point is not reinitialised, then nothing is done */
2749  if( benders->corepoint != NULL && benders->strengthenintpoint != 'i' )
2750  return SCIP_OKAY;
2751 
2752  bestsol = SCIPgetBestSol(scip);
2753 
2754  /* if the core point should be updated, then this only happens if the incumbent solution has been updated */
2755  if( benders->strengthenintpoint == 'i' && benders->initcorepoint == bestsol )
2756  return SCIP_OKAY;
2757 
2758  /* if a corepoint has been used for cut strengthening, then this needs to be freed */
2759  if( benders->corepoint != NULL )
2760  {
2761  SCIP_CALL( SCIPfreeSol(scip, &benders->corepoint) );
2762  }
2763 
2764  switch( benders->strengthenintpoint )
2765  {
2766  SCIP_VAR** vars;
2767  SCIP_Real timelimit;
2768  int nvars;
2769  int i;
2770 
2771  case 'l':
2772  SCIP_CALL( SCIPcreateLPSol(scip, &benders->corepoint, NULL) );
2773  SCIP_CALL( SCIPunlinkSol(scip, benders->corepoint) );
2774  break;
2775  case 'f':
2776  case 'i':
2777  SCIP_CALL( SCIPcreateSolCopy(scip, &benders->corepoint, bestsol) );
2778  SCIP_CALL( SCIPunlinkSol(scip, benders->corepoint) );
2779  benders->initcorepoint = bestsol;
2780  break;
2781  case 'r':
2782  /* prepare time limit */
2783  SCIP_CALL( SCIPgetRealParam(scip, "limits/time", &timelimit) );
2784  if ( ! SCIPisInfinity(scip, timelimit) )
2785  timelimit -= SCIPgetSolvingTime(scip);
2786 
2787  /* if there is time remaining, then compute the relative interior point. Otherwise, return the LP solution */
2788  if ( timelimit > 0.0 )
2789  {
2790  SCIPverbMessage(scip, SCIP_VERBLEVEL_MINIMAL, 0, "Computing relative interior point (time limit: %g, iter limit: %d) ...\n", timelimit, INT_MAX);
2791  SCIP_CALL( SCIPcomputeLPRelIntPoint(scip, TRUE, FALSE, timelimit, INT_MAX, &benders->corepoint) );
2792  }
2793  else
2794  {
2795  SCIP_CALL( SCIPcreateLPSol(scip, &benders->corepoint, NULL) );
2796  SCIP_CALL( SCIPunlinkSol(scip, benders->corepoint) );
2797  }
2798  break;
2799  case 'z':
2800  SCIP_CALL( SCIPcreateSol(scip, &benders->corepoint, NULL) );
2801  break;
2802  case 'o':
2803  SCIP_CALL( SCIPcreateSol(scip, &benders->corepoint, NULL) );
2804 
2805  /* getting the variable data so that the */
2806  SCIP_CALL( SCIPgetVarsData(scip, &vars, &nvars, NULL, NULL, NULL, NULL) );
2807 
2808  /* setting all variable values to 1.0 */
2809  for( i = 0; i < nvars; i++ )
2810  {
2811  SCIP_CALL( SCIPsetSolVal(scip, benders->corepoint, vars[i], 1.0) );
2812  }
2813  break;
2814  default:
2815  SCIP_CALL( SCIPcreateLPSol(scip, &benders->corepoint, NULL) );
2816  SCIP_CALL( SCIPunlinkSol(scip, benders->corepoint) );
2817  }
2818 
2819  return SCIP_OKAY;
2820 }
2821 
2822 /** performs cut strengthening by using an interior solution to generate cuts */
2823 static
2825  SCIP_BENDERS* benders, /**< Benders' decomposition */
2826  SCIP_SET* set, /**< global SCIP settings */
2827  SCIP_SOL* sol, /**< primal CIP solution */
2828  SCIP_BENDERSENFOTYPE type, /**< the type of solution being enforced */
2829  SCIP_Bool checkint, /**< are the subproblems called during a check/enforce of integer sols? */
2830  SCIP_Bool perturbsol, /**< should the solution be perturbed to escape infeasibility? */
2831  SCIP_Bool* auxviol, /**< set to TRUE only if the solution is feasible but the aux vars are violated */
2832  SCIP_Bool* infeasible, /**< is the master problem infeasible with respect to the Benders' cuts? */
2833  SCIP_Bool* skipsolve, /**< should the main solve be skipped as a result of this strengthening? */
2834  SCIP_RESULT* result /**< result of the pricing process */
2835  )
2836 {
2837  SCIP_SOL* sepapoint;
2838  SCIP_VAR** vars;
2839  int prevcutsfound;
2840  int nvars;
2841  int i;
2842 
2843  assert(benders != NULL);
2844  assert(set != NULL);
2845 
2846  (*result) = SCIP_DIDNOTRUN;
2847  (*skipsolve) = FALSE;
2848 
2849  /* the cut stabilisation is only performed when enforcing LP solutions. The solution is not NULL if the stabilisation
2850  * is currently being performed. It is important to avoid recursion
2851  */
2852  if( type != SCIP_BENDERSENFOTYPE_LP || sol != NULL )
2853  return SCIP_OKAY;
2854 
2855  /* checking if a change to the lower bound has occurred */
2856  if( SCIPsetIsGT(set, SCIPgetLowerbound(set->scip), benders->prevlowerbound)
2857  || SCIPgetCurrentNode(set->scip) != benders->prevnode )
2858  {
2859  benders->prevnode = SCIPgetCurrentNode(set->scip);
2860  benders->prevlowerbound = SCIPgetLowerbound(set->scip);
2861  benders->noimprovecount = 0;
2862  }
2863  else
2864  benders->noimprovecount++;
2865 
2866  /* if the number of iterations without improvement exceeds 3*noimprovelimit, then the no stabilisation is performed
2867  */
2868  if( benders->noimprovecount > 3*benders->noimprovelimit )
2869  return SCIP_OKAY;
2870 
2871  /* if there is no incumbent solution, then it is not possible to create the core point and hence the strengthening
2872  * can not be performed
2873  */
2874  if( SCIPgetBestSol(set->scip) == NULL )
2875  return SCIP_OKAY;
2876 
2877  /* if no LP iterations have been performed since the last call of the cut strenghtening, then the strengthening is
2878  * aborted
2879  */
2880  if( benders->prevnlpiter == SCIPgetNLPIterations(set->scip) )
2881  return SCIP_OKAY;
2882 
2883  benders->prevnlpiter = SCIPgetNLPIterations(set->scip);
2884 
2885  /* if the separation point solution is NULL, then we create the solution using the current LP relaxation. */
2886  SCIP_CALL( setAndUpdateCorePoint(set->scip, benders) );
2887 
2888  /* creating the separation point
2889  * TODO: This could be a little to memory heavy, it may be better just to create the separation point once and then
2890  * update it each time.
2891  */
2892  SCIP_CALL( SCIPcreateLPSol(set->scip, &sepapoint, NULL) );
2893  SCIP_CALL( SCIPunlinkSol(set->scip, sepapoint) );
2894 
2895  SCIP_CALL( SCIPgetVarsData(set->scip, &vars, &nvars, NULL, NULL, NULL, NULL) );
2896  assert(vars != NULL);
2897 
2898  /* creating a solution that is a convex combination of the LP solution and the separation point */
2899  for( i = 0; i < nvars; i++ )
2900  {
2901  SCIP_VAR* subvar;
2902  SCIP_Real corepointval;
2903  SCIP_Real lpsolval;
2904  SCIP_Real newsolval;
2905  int j;
2906 
2907  corepointval = SCIPgetSolVal(set->scip, benders->corepoint, vars[i]);
2908  lpsolval = SCIPgetSolVal(set->scip, sol, vars[i]);
2909  newsolval = lpsolval;
2910 
2911  /* checking whether the master variable is mapped to any subproblem variables */
2912  subvar = NULL;
2913  j = 0;
2914  while( subvar == NULL && j < SCIPgetBendersNSubproblems(set->scip, benders) )
2915  {
2916  SCIP_CALL( SCIPgetBendersSubproblemVar(set->scip, benders, vars[i], &subvar, j) );
2917  j++;
2918  }
2919 
2920  /* if the variable is a linking variable and it is not fixed, then a convex combination with the corepoint is
2921  * computed.
2922  */
2923  if( subvar != NULL && SCIPvarGetStatus(vars[i]) != SCIP_VARSTATUS_FIXED )
2924  {
2925  /* if the number of iterations without improvement exceeds noimprovelimit, then no convex combination is
2926  * created
2927  */
2928  if( !perturbsol && benders->noimprovecount <= benders->noimprovelimit )
2929  {
2930  newsolval = lpsolval*benders->convexmult + corepointval*(1 - benders->convexmult);
2931 
2932  /* updating the core point */
2933  SCIP_CALL( SCIPsetSolVal(set->scip, benders->corepoint, vars[i], newsolval) );
2934  }
2935 
2936  /* if the number of iterations without improvement is less than 2*noimprovelimit, then perturbation is
2937  * performed
2938  * TODO: This should be a random vector!!!!
2939  */
2940  if( perturbsol || benders->noimprovecount <= 2*benders->noimprovelimit )
2941  newsolval += benders->perturbeps;
2942  }
2943 
2944  /* updating the separation point */
2945  SCIP_CALL( SCIPsetSolVal(set->scip, sepapoint, vars[i], newsolval) );
2946  }
2947 
2948  /* storing the number of cuts found */
2949  prevcutsfound = SCIPbendersGetNCutsFound(benders);
2950 
2951  SCIPsetDebugMsg(set, "solving Benders' decomposition subproblems with stabilised point.\n");
2952 
2953  /* calling the subproblem solving method to generate cuts from the separation solution */
2954  SCIP_CALL( SCIPsolveBendersSubproblems(set->scip, benders, sepapoint, result, infeasible, auxviol, type, checkint) );
2955 
2956  SCIPsetDebugMsg(set, "solved Benders' decomposition subproblems with stabilised point. noimprovecount %d result %d\n",
2957  benders->noimprovecount, (*result));
2958 
2959  /* if constraints were added, then the main Benders' solving loop is skipped. */
2960  if( !(*infeasible) && ((*result) == SCIP_CONSADDED || (*result) == SCIP_SEPARATED) )
2961  (*skipsolve) = TRUE;
2962 
2963  /* capturing cut strengthening statistics */
2964  benders->nstrengthencalls++;
2965  benders->nstrengthencuts += (SCIPbendersGetNCutsFound(benders) - prevcutsfound);
2966 
2967  /* if no cuts were added, then the strengthening round is marked as failed */
2968  if( SCIPbendersGetNCutsFound(benders) == prevcutsfound )
2969  benders->nstrengthenfails++;
2970 
2971  /* freeing the sepapoint solution */
2972  SCIP_CALL( SCIPfreeSol(set->scip, &sepapoint) );
2973 
2974  return SCIP_OKAY;
2975 }
2976 
2977 
2978 /** Returns whether only the convex relaxations will be checked in this solve loop
2979  * when Benders' is used in the LNS heuristics, only the convex relaxations of the master/subproblems are checked,
2980  * i.e. no integer cuts are generated. In this case, then Benders' decomposition is performed under the assumption
2981  * that all subproblems are convex relaxations.
2982  */
2984  SCIP_BENDERS* benders, /**< Benders' decomposition */
2985  SCIP_Bool subscipsoff /**< flag indicating whether plugins using sub-SCIPs are deactivated */
2986  )
2987 {
2988  return benders->iscopy && benders->lnscheck && subscipsoff;
2989 }
2990 
2991 /** returns the number of subproblems that will be checked in this iteration */
2992 static
2994  SCIP_BENDERS* benders, /**< Benders' decomposition */
2995  SCIP_SET* set, /**< global SCIP settings */
2996  SCIP_BENDERSENFOTYPE type /**< the type of solution being enforced */
2997  )
2998 {
2999  if( benders->ncalls == 0 || type == SCIP_BENDERSENFOTYPE_CHECK
3001  return SCIPbendersGetNSubproblems(benders);
3002  else
3003  return (int) SCIPsetCeil(set, (SCIP_Real) SCIPbendersGetNSubproblems(benders)*benders->subprobfrac);
3004 }
3005 
3006 /** returns whether the solving of the given subproblem needs to be executed */
3007 static
3009  SCIP_BENDERS* benders, /**< Benders' decomposition */
3010  int probnumber /**< the subproblem index */
3011  )
3012 {
3013  return (!SCIPbendersSubproblemIsIndependent(benders, probnumber)
3014  && SCIPbendersSubproblemIsEnabled(benders, probnumber));
3015 }
3016 
3017 /** creates an ordered list of subproblem indices to be solved */
3018 static
3020  SCIP_BENDERS* benders, /**< Benders' decomposition */
3021  SCIP_SET* set, /**< global SCIP settings */
3022  SCIP_BENDERSENFOTYPE type, /**< the type of solution being enforced */
3023  int** solveidx, /**< a list of subproblem indices to the solved in the current iteration */
3024  int* nsolveidx /**< the number of subproblem indices in the list */
3025  )
3026 {
3027  int nsubproblems;
3028  int numtocheck;
3029  int subproblemcount;
3030 
3031  assert(benders != NULL);
3032  assert(set != NULL);
3033  assert((*solveidx) != NULL);
3034  assert(nsolveidx != NULL);
3035  assert(SCIPpqueueNElems(benders->subprobqueue) <= SCIPbendersGetNSubproblems(benders));
3036 
3037  nsubproblems = SCIPbendersGetNSubproblems(benders);
3038 
3039  /* it is possible to only solve a subset of subproblems. This is given by a parameter. */
3040  numtocheck = numSubproblemsToCheck(benders, set, type);
3041 
3042  (*nsolveidx) = 0;
3043 
3044  subproblemcount = 0;
3045  while( subproblemcount < nsubproblems && subproblemcount < numtocheck )
3046  {
3047  SCIP_SUBPROBLEMSOLVESTAT* solvestat;
3048 
3049  solvestat = (SCIP_SUBPROBLEMSOLVESTAT*)SCIPpqueueRemove(benders->subprobqueue);
3050  (*solveidx)[(*nsolveidx)] = solvestat->idx;
3051  (*nsolveidx)++;
3052 
3053  subproblemcount++;
3054  }
3055 }
3056 
3057 /** updates the subproblem solving statistics and inserts the indices into the queue */
3058 static
3060  SCIP_BENDERS* benders, /**< Benders' decomposition */
3061  int* solveidx, /**< the list of indices of subproblems that were solved */
3062  int nsolveidx, /**< the number of subproblem indices */
3063  SCIP_Bool updatestat /**< should the statistics be updated */
3064  )
3065 {
3066  int i;
3067 
3068  assert(benders != NULL);
3069  assert(solveidx != NULL);
3070 
3071  for( i = 0; i < nsolveidx; i++ )
3072  {
3073  SCIP* subproblem;
3074  SCIP_SUBPROBLEMSOLVESTAT* solvestat;
3075 
3076  subproblem = SCIPbendersSubproblem(benders, solveidx[i]);
3077  solvestat = benders->solvestat[solveidx[i]];
3078  assert(solvestat->idx == solveidx[i]);
3079 
3080  /* updating the solving statistics */
3081  if( updatestat )
3082  {
3083  if( subproblem == NULL )
3084  solvestat->avgiter = 1;
3085  else
3086  solvestat->avgiter = (SCIP_Real)(solvestat->avgiter*solvestat->ncalls + SCIPgetNLPIterations(subproblem))
3087  /(SCIP_Real)(solvestat->ncalls + 1);
3088  solvestat->ncalls++;
3089  }
3090 
3091  /* inserting the solving statistics into the priority queue */
3092  SCIP_CALL( SCIPpqueueInsert(benders->subprobqueue, solvestat) );
3093  }
3094 
3095  assert(SCIPpqueueNElems(benders->subprobqueue) == SCIPbendersGetNSubproblems(benders));
3096 
3097  return SCIP_OKAY;
3098 }
3099 
3100 /** Solves each of the Benders' decomposition subproblems for the given solution. All, or a fraction, of subproblems are
3101  * solved before the Benders' decomposition cuts are generated.
3102  * Since a convex relaxation of the subproblem could be solved to generate cuts, a parameter nverified is used to
3103  * identified the number of subproblems that have been solved in their "original" form. For example, if the subproblem
3104  * is a MIP, then if the LP is solved to generate cuts, this does not constitute a verification. The verification is
3105  * only performed when the MIP is solved.
3106  */
3107 static
3109  SCIP_BENDERS* benders, /**< Benders' decomposition */
3110  SCIP_SET* set, /**< global SCIP settings */
3111  SCIP_SOL* sol, /**< primal CIP solution */
3112  SCIP_BENDERSENFOTYPE type, /**< the type of solution being enforced */
3113  SCIP_BENDERSSOLVELOOP solveloop, /**< the current solve loop */
3114  SCIP_Bool checkint, /**< are the subproblems called during a check/enforce of integer sols? */
3115  int* nverified, /**< the number of subproblems verified in the current loop */
3116  int* solveidx, /**< the indices of subproblems to be solved in this loop */
3117  int nsolveidx, /**< the number of subproblems to be solved in this loop */
3118  SCIP_Bool** subprobsolved, /**< an array indicating the subproblems that were solved in this loop. */
3119  SCIP_BENDERSSUBSTATUS** substatus, /**< array to store the status of the subsystem */
3120  SCIP_Bool* infeasible, /**< is the master problem infeasible with respect to the Benders' cuts? */
3121  SCIP_Bool* optimal, /**< is the current solution optimal? */
3122  SCIP_Bool* stopped /**< was the solving process stopped? */
3123  )
3124 {
3125  SCIP_Bool onlyconvexcheck;
3126 #ifdef _OPENMP
3127  int numthreads;
3128  int maxnthreads;
3129 #endif
3130  int i;
3131  int j;
3132 
3133  /* local variables for parallelisation of the solving loop */
3134  int locnverified = *nverified;
3135  SCIP_Bool locinfeasible = *infeasible;
3136  SCIP_Bool locoptimal = *optimal;
3137  SCIP_Bool locstopped = *stopped;
3138 
3139  SCIP_RETCODE retcode = SCIP_OKAY;
3140 
3141  assert(benders != NULL);
3142  assert(set != NULL);
3143 
3144  /* getting the number of threads to use when solving the subproblems. This will be either be
3145  * min(numthreads, maxnthreads).
3146  * NOTE: This may not be correct. The Benders' decomposition parallelisation should not take all minimum threads if
3147  * they are specified. The number of threads should be specified with the Benders' decomposition parameters.
3148  */
3149 #ifdef _OPENMP
3150  SCIP_CALL( SCIPsetGetIntParam(set, "parallel/maxnthreads", &maxnthreads) );
3151  numthreads = MIN(benders->numthreads, maxnthreads);
3152 #endif
3153 
3154  /* in the case of an LNS check, only the convex relaxations of the subproblems will be solved. This is a performance
3155  * feature, since solving the convex relaxation is typically much faster than solving the corresponding CIP. While
3156  * the CIP is not solved during the LNS check, the solutions are still of higher quality than when Benders' is not
3157  * employed.
3158  */
3159  onlyconvexcheck = SCIPbendersOnlyCheckConvexRelax(benders, SCIPsetGetSubscipsOff(set));
3160 
3161  SCIPsetDebugMsg(set, "Performing the subproblem solving process. Number of subproblems to check %d\n", nsolveidx);
3162 
3163  SCIPsetDebugMsg(set, "Benders' decomposition - solve loop %d\n", solveloop);
3164 
3165  if( type == SCIP_BENDERSENFOTYPE_CHECK && sol == NULL )
3166  {
3167  /* TODO: Check whether this is absolutely necessary. I think that this if statment can be removed. */
3168  locinfeasible = TRUE;
3169  }
3170  else
3171  {
3172  /* solving each of the subproblems for Benders' decomposition */
3173  /* TODO: ensure that the each of the subproblems solve and update the parameters with the correct return values
3174  */
3175 #ifndef __INTEL_COMPILER
3176  #pragma omp parallel for num_threads(numthreads) private(i) reduction(&&:locoptimal) reduction(||:locinfeasible) reduction(+:locnverified) reduction(||:locstopped) reduction(min:retcode)
3177 #endif
3178  for( j = 0; j < nsolveidx; j++ )
3179  {
3180  SCIP_Bool subinfeas = FALSE;
3181  SCIP_Bool convexsub;
3182  SCIP_Bool solvesub = TRUE;
3183  SCIP_Bool solved;
3184 
3185  i = solveidx[j];
3187 
3188  /* the subproblem is initially flagged as not solved for this solving loop */
3189  (*subprobsolved)[i] = FALSE;
3190 
3191  /* setting the subsystem status to UNKNOWN at the start of each solve loop */
3192  (*substatus)[i] = SCIP_BENDERSSUBSTATUS_UNKNOWN;
3193 
3194  /* for the second solving loop, if the problem is an LP, it is not solved again. If the problem is a MIP,
3195  * then the subproblem objective function value is set to infinity. However, if the subproblem is proven
3196  * infeasible from the LP, then the IP loop is not performed.
3197  * If the solve loop is SCIP_BENDERSSOLVELOOP_USERCIP, then nothing is done. It is assumed that the user will
3198  * correctly update the objective function within the user-defined solving function.
3199  */
3200  if( solveloop == SCIP_BENDERSSOLVELOOP_CIP )
3201  {
3202  if( convexsub || (*substatus)[i] == SCIP_BENDERSSUBSTATUS_INFEAS )
3203  solvesub = FALSE;
3204  else
3205  {
3206  SCIPbendersSetSubproblemObjval(benders, i, SCIPbendersSubproblem(benders, i) != NULL ?
3207  SCIPinfinity(SCIPbendersSubproblem(benders, i)) : SCIPsetInfinity(set));
3208  }
3209  }
3210 
3211  /* if the subproblem is independent, then it does not need to be solved. In this case, the nverified flag will
3212  * increase by one. When the subproblem is not independent, then it needs to be checked.
3213  */
3214  if( !subproblemIsActive(benders, i) )
3215  {
3216  /* NOTE: There is no need to update the optimal flag. This is because optimal is always TRUE until a
3217  * non-optimal subproblem is found.
3218  */
3219  /* if the auxiliary variable value is infinity, then the subproblem has not been solved yet. Currently the
3220  * subproblem statue is unknown. */
3221  if( SCIPsetIsInfinity(set, SCIPbendersGetAuxiliaryVarVal(benders, set, sol, i))
3222  || SCIPsetIsInfinity(set, -SCIPbendersGetAuxiliaryVarVal(benders, set, sol, i))
3224  {
3225  SCIPbendersSetSubproblemObjval(benders, i, SCIPbendersSubproblem(benders, i) != NULL ?
3226  SCIPinfinity(SCIPbendersSubproblem(benders, i)) : SCIPsetInfinity(set));
3227 
3228  (*substatus)[i] = SCIP_BENDERSSUBSTATUS_UNKNOWN;
3229  locoptimal = FALSE;
3230 
3231  SCIPsetDebugMsg(set, "Benders' decomposition: subproblem %d is not active, but has not been solved."
3232  " setting status to UNKNOWN\n", i);
3233  }
3234  else
3235  {
3237  SCIPbendersGetAuxiliaryVarVal(benders, set, sol, i)) < benders->solutiontol )
3238  {
3239  SCIPbendersSetSubproblemObjval(benders, i, SCIPbendersGetAuxiliaryVarVal(benders, set, sol, i));
3240  (*substatus)[i] = SCIP_BENDERSSUBSTATUS_OPTIMAL;
3241  }
3242  else
3243  {
3245  (*substatus)[i] = SCIP_BENDERSSUBSTATUS_AUXVIOL;
3246  }
3247 
3248  SCIPsetDebugMsg(set, "Benders' decomposition: subproblem %d is not active, setting status to OPTIMAL\n", i);
3249  }
3250 
3251  (*subprobsolved)[i] = TRUE;
3252 
3253  /* the nverified counter is only increased in the convex solving loop */
3254  if( solveloop == SCIP_BENDERSSOLVELOOP_CONVEX || solveloop == SCIP_BENDERSSOLVELOOP_USERCONVEX )
3255  locnverified++;
3256  }
3257  else if( solvesub )
3258  {
3259  retcode = SCIPbendersExecSubproblemSolve(benders, set, sol, i, solveloop, FALSE, &solved, &subinfeas, type);
3260 
3261  /* the solution for the subproblem is only processed if the return code is SCIP_OKAY */
3262  if( retcode == SCIP_OKAY )
3263  {
3264 #ifdef SCIP_DEBUG
3265  if( type == SCIP_BENDERSENFOTYPE_LP )
3266  {
3267  SCIPsetDebugMsg(set, "Enfo LP: Subproblem %d Type %d (%f < %f)\n", i,
3268  SCIPbendersGetSubproblemType(benders, i), SCIPbendersGetAuxiliaryVarVal(benders, set, sol, i),
3269  SCIPbendersGetSubproblemObjval(benders, i));
3270  }
3271 #endif
3272  (*subprobsolved)[i] = solved;
3273 
3274  locinfeasible = locinfeasible || subinfeas;
3275  if( subinfeas )
3276  (*substatus)[i] = SCIP_BENDERSSUBSTATUS_INFEAS;
3277 
3278  /* if the subproblems are solved to check integer feasibility, then the optimality check must be performed.
3279  * This will only be performed if checkint is TRUE and the subproblem was solved. The subproblem may not be
3280  * solved if the user has defined a solving function
3281  */
3282  if( checkint && (*subprobsolved)[i] )
3283  {
3284  /* if the subproblem is feasible, then it is necessary to update the value of the auxiliary variable to the
3285  * objective function value of the subproblem.
3286  */
3287  if( !subinfeas )
3288  {
3289  SCIP_Bool subproboptimal;
3290 
3291  subproboptimal = SCIPbendersSubproblemIsOptimal(benders, set, sol, i);
3292 
3293  if( subproboptimal )
3294  (*substatus)[i] = SCIP_BENDERSSUBSTATUS_OPTIMAL;
3295  else
3296  (*substatus)[i] = SCIP_BENDERSSUBSTATUS_AUXVIOL;
3297 
3298  /* It is only possible to determine the optimality of a solution within a given subproblem in four
3299  * different cases:
3300  * i) solveloop == SCIP_BENDERSSOLVELOOP_CONVEX or USERCONVEX and the subproblem is convex.
3301  * ii) solveloop == SCIP_BENDERSOLVELOOP_CONVEX and only the convex relaxations will be checked.
3302  * iii) solveloop == SCIP_BENDERSSOLVELOOP_USERCIP and the subproblem was solved, since the user has
3303  * defined a solve function, it is expected that the solving is correctly executed.
3304  * iv) solveloop == SCIP_BENDERSSOLVELOOP_CIP and the MIP for the subproblem has been solved.
3305  */
3306  if( convexsub || onlyconvexcheck
3307  || solveloop == SCIP_BENDERSSOLVELOOP_CIP
3308  || solveloop == SCIP_BENDERSSOLVELOOP_USERCIP )
3309  locoptimal = locoptimal && subproboptimal;
3310 
3311 #ifdef SCIP_DEBUG
3312  if( convexsub || solveloop >= SCIP_BENDERSSOLVELOOP_CIP )
3313  {
3314  if( subproboptimal )
3315  {
3316  SCIPsetDebugMsg(set, "Subproblem %d is Optimal (%f >= %f)\n", i,
3317  SCIPbendersGetAuxiliaryVarVal(benders, set, sol, i), SCIPbendersGetSubproblemObjval(benders, i));
3318  }
3319  else
3320  {
3321  SCIPsetDebugMsg(set, "Subproblem %d is NOT Optimal (%f < %f)\n", i,
3322  SCIPbendersGetAuxiliaryVarVal(benders, set, sol, i), SCIPbendersGetSubproblemObjval(benders, i));
3323  }
3324  }
3325 #endif
3326 
3327  /* the nverified variable is only incremented when the original form of the subproblem has been solved.
3328  * What is meant by "original" is that the LP relaxation of CIPs are solved to generate valid cuts. So
3329  * if the subproblem is defined as a CIP, then it is only classified as checked if the CIP is solved.
3330  * There are three cases where the "original" form is solved are:
3331  * i) solveloop == SCIP_BENDERSSOLVELOOP_CONVEX or USERCONVEX and the subproblem is an LP
3332  * - the original form has been solved.
3333  * ii) solveloop == SCIP_BENDERSSOLVELOOP_CIP or USERCIP and the CIP for the subproblem has been
3334  * solved.
3335  * iii) or, only a convex check is performed.
3336  */
3337  if( ((solveloop == SCIP_BENDERSSOLVELOOP_CONVEX || solveloop == SCIP_BENDERSSOLVELOOP_USERCONVEX)
3338  && convexsub)
3339  || ((solveloop == SCIP_BENDERSSOLVELOOP_CIP || solveloop == SCIP_BENDERSSOLVELOOP_USERCIP)
3340  && !convexsub)
3341  || onlyconvexcheck )
3342  locnverified++;
3343  }
3344  }
3345  }
3346  }
3347 
3348  /* checking whether the limits have been exceeded in the master problem */
3349  locstopped = SCIPisStopped(set->scip);
3350  }
3351  }
3352 
3353  /* setting the input parameters to the local variables */
3354  SCIPsetDebugMsg(set, "Local variable values: nverified %d infeasible %u optimal %u stopped %u\n", locnverified,
3355  locinfeasible, locoptimal, locstopped);
3356  *nverified = locnverified;
3357  *infeasible = locinfeasible;
3358  *optimal = locoptimal;
3359  *stopped = locstopped;
3360 
3361  return retcode;
3362 }
3363 
3364 /** Calls the Benders' decompsition cuts for the given solve loop. There are four cases:
3365  * i) solveloop == SCIP_BENDERSSOLVELOOP_CONVEX - only the LP Benders' cuts are called
3366  * ii) solveloop == SCIP_BENDERSSOLVELOOP_CIP - only the CIP Benders' cuts are called
3367  * iii) solveloop == SCIP_BENDERSSOLVELOOP_USERCONVEX - only the LP Benders' cuts are called
3368  * iv) solveloop == SCIP_BENDERSSOLVELOOP_USERCIP - only the CIP Benders' cuts are called
3369  *
3370  * The priority of the results are: SCIP_CONSADDED (SCIP_SEPARATED), SCIP_DIDNOTFIND, SCIP_FEASIBLE, SCIP_DIDNOTRUN. In
3371  * this function, there are four levels of results that need to be assessed. These are:
3372  * i) The result from the individual cut for the subproblem
3373  * ii) The overall result for the subproblem from all cuts
3374  * iii) the overall result for the solve loop from all cuts
3375  * iv) the over all result from all solve loops.
3376  * In each level, the priority of results must be adhered to.
3377  */
3378 static
3380  SCIP_BENDERS* benders, /**< Benders' decomposition */
3381  SCIP_SET* set, /**< global SCIP settings */
3382  SCIP_SOL* sol, /**< primal CIP solution */
3383  SCIP_RESULT* result, /**< result of the pricing process */
3384  SCIP_BENDERSENFOTYPE type, /**< the type of solution being enforced */
3385  SCIP_BENDERSSOLVELOOP solveloop, /**< the current solve loop */
3386  SCIP_Bool checkint, /**< are the subproblems called during a check/enforce of integer sols? */
3387  SCIP_Bool* subprobsolved, /**< an array indicating the subproblems that were solved in this loop. */
3388  SCIP_BENDERSSUBSTATUS* substatus, /**< array to store the status of the subsystem */
3389  int* solveidx, /**< the indices of subproblems to be solved in this loop */
3390  int nsolveidx, /**< the number of subproblems to be solved in this loop */
3391  int** mergecands, /**< the subproblems that are merge candidates */
3392  int* npriomergecands, /**< the number of priority merge candidates. */
3393  int* nmergecands, /**< the number of merge candidates. */
3394  int* nsolveloops /**< the number of solve loops, is updated w.r.t added cuts */
3395  )
3396 {
3397  SCIP_BENDERSCUT** benderscuts;
3398  SCIP_RESULT solveloopresult;
3399  int nbenderscuts;
3400  SCIP_Longint addedcuts = 0;
3401  int i;
3402  int j;
3403  int k;
3404  SCIP_Bool onlyconvexcheck;
3405 
3406  assert(benders != NULL);
3407  assert(set != NULL);
3408 
3409  /* getting the Benders' decomposition cuts */
3410  benderscuts = SCIPbendersGetBenderscuts(benders);
3411  nbenderscuts = SCIPbendersGetNBenderscuts(benders);
3412 
3413  solveloopresult = SCIP_DIDNOTRUN;
3414 
3415  /* in the case of an LNS check, only the convex relaxations of the subproblems will be solved. This is a performance
3416  * feature, since solving the convex relaxation is typically much faster than solving the corresponding CIP. While
3417  * the CIP is not solved during the LNS check, the solutions are still of higher quality than when Benders' is not
3418  * employed.
3419  */
3420  onlyconvexcheck = SCIPbendersOnlyCheckConvexRelax(benders, SCIPsetGetSubscipsOff(set));
3421 
3422  /* It is only possible to add cuts to the problem if it has not already been solved */
3425  && (benders->cutcheck || type != SCIP_BENDERSENFOTYPE_CHECK) )
3426  {
3427  /* This is done in two loops. The first is by subproblem and the second is by cut type. */
3428  for( k = 0; k < nsolveidx; k++ )
3429  {
3430  SCIP_RESULT subprobresult;
3431  SCIP_Bool convexsub;
3432 
3433  i = solveidx[k];
3434 
3436 
3437  /* cuts can only be generated if the subproblem is not independent and if it has been solved. Additionally, the
3438  * status of the subproblem solving must not be INFEASIBLE while in a cut strengthening round.
3439  * The subproblem solved flag is important for the user-defined subproblem solving methods
3440  */
3441  if( subproblemIsActive(benders, i) && subprobsolved[i]
3442  && !(substatus[i] == SCIP_BENDERSSUBSTATUS_INFEAS && benders->strengthenround) )
3443  {
3444  subprobresult = SCIP_DIDNOTRUN;
3445  for( j = 0; j < nbenderscuts; j++ )
3446  {
3447  SCIP_RESULT cutresult;
3448  SCIP_Longint prevaddedcuts;
3449 
3450  assert(benderscuts[j] != NULL);
3451 
3452  prevaddedcuts = SCIPbenderscutGetNFound(benderscuts[j]);
3453  cutresult = SCIP_DIDNOTRUN;
3454 
3455  /* the result is updated only if a Benders' cut is generated or one was not found. However, if a cut has
3456  * been found in a previous iteration, then the result is returned as SCIP_CONSADDED or SCIP_SEPARATED.
3457  * This result is permitted because if a constraint was added, the solution that caused the error in the cut
3458  * generation will be cutoff from the master problem.
3459  */
3460  if( (SCIPbenderscutIsLPCut(benderscuts[j]) && (solveloop == SCIP_BENDERSSOLVELOOP_CONVEX
3461  || solveloop == SCIP_BENDERSSOLVELOOP_USERCONVEX))
3462  || (!SCIPbenderscutIsLPCut(benderscuts[j]) && ((solveloop == SCIP_BENDERSSOLVELOOP_CIP && !convexsub)
3463  || solveloop == SCIP_BENDERSSOLVELOOP_USERCIP)) )
3464  SCIP_CALL( SCIPbenderscutExec(benderscuts[j], set, benders, sol, i, type, &cutresult) );
3465 
3466  addedcuts += (SCIPbenderscutGetNFound(benderscuts[j]) - prevaddedcuts);
3467 
3468  /* the result is updated only if a Benders' cut is generated */
3469  if( cutresult == SCIP_CONSADDED || cutresult == SCIP_SEPARATED )
3470  {
3471  subprobresult = cutresult;
3472 
3473  benders->ncutsfound++;
3474 
3475  /* at most a single cut is generated for each subproblem */
3476  break;
3477  }
3478  else
3479  {
3480  /* checking from lowest priority result */
3481  if( subprobresult == SCIP_DIDNOTRUN )
3482  subprobresult = cutresult;
3483  else if( subprobresult == SCIP_FEASIBLE && cutresult == SCIP_DIDNOTFIND )
3484  subprobresult = cutresult;
3485  /* if the subprobresult is SCIP_DIDNOTFIND, then it can't be updated. */
3486  }
3487  }
3488 
3489  /* the highest priority for the results is CONSADDED and SEPARATED. The solveloopresult will always be
3490  * updated if the subprobresult is either of these.
3491  */
3492  if( subprobresult == SCIP_CONSADDED || subprobresult == SCIP_SEPARATED )
3493  {
3494  solveloopresult = subprobresult;
3495  }
3496  else if( subprobresult == SCIP_FEASIBLE )
3497  {
3498  /* updating the solve loop result based upon the priority */
3499  if( solveloopresult == SCIP_DIDNOTRUN )
3500  solveloopresult = subprobresult;
3501  }
3502  else if( subprobresult == SCIP_DIDNOTFIND )
3503  {
3504  /* updating the solve loop result based upon the priority */
3505  if( solveloopresult == SCIP_DIDNOTRUN || solveloopresult == SCIP_FEASIBLE )
3506  solveloopresult = subprobresult;
3507 
3508  /* since a cut was not found, then merging could be useful to avoid this in subsequent iterations. The
3509  * candidate is labelled as a non-priority merge candidate
3510  */
3511  if( substatus[i] != SCIP_BENDERSSUBSTATUS_OPTIMAL )
3512  {
3513  (*mergecands)[(*nmergecands)] = i;
3514  (*nmergecands)++;
3515  }
3516  }
3517  else if( subprobresult == SCIP_DIDNOTRUN )
3518  {
3519  /* if the subproblem is infeasible and no cut generation methods were run, then the infeasibility will
3520  * never be resolved. As such, the subproblem will be merged into the master problem. If the subproblem
3521  * was not infeasible, then it is added as a possible merge candidate
3522  */
3523  if( substatus[i] == SCIP_BENDERSSUBSTATUS_INFEAS )
3524  {
3525  (*mergecands)[(*nmergecands)] = (*mergecands)[(*npriomergecands)];
3526  (*mergecands)[(*npriomergecands)] = i;
3527  (*npriomergecands)++;
3528  (*nmergecands)++;
3529  }
3530  else if( substatus[i] != SCIP_BENDERSSUBSTATUS_OPTIMAL )
3531  {
3532  (*mergecands)[(*nmergecands)] = i;
3533  (*nmergecands)++;
3534  }
3535  }
3536  }
3537  }
3538  }
3539 
3540  /* updating the overall result based upon the priorities */
3541  if( solveloopresult == SCIP_CONSADDED || solveloopresult == SCIP_SEPARATED )
3542  {
3543  (*result) = solveloopresult;
3544  }
3545  else if( solveloopresult == SCIP_FEASIBLE )
3546  {
3547  /* updating the solve loop result based upon the priority */
3548  if( (*result) == SCIP_DIDNOTRUN )
3549  (*result) = solveloopresult;
3550  }
3551  else if( solveloopresult == SCIP_DIDNOTFIND )
3552  {
3553  /* updating the solve loop result based upon the priority */
3554  if( (*result) == SCIP_DIDNOTRUN || (*result) == SCIP_FEASIBLE )
3555  (*result) = solveloopresult;
3556  }
3557 
3558  /* if no cuts were added, then the number of solve loops is increased */
3559  if( addedcuts == 0 && SCIPbendersGetNConvexSubproblems(benders) < SCIPbendersGetNSubproblems(benders)
3560  && checkint && !onlyconvexcheck )
3561  (*nsolveloops) = 2;
3562 
3563  return SCIP_OKAY;
3564 }
3565 
3566 /** Solves the subproblem using the current master problem solution.
3567  *
3568  * The checkint flag indicates whether integer feasibility can be assumed. If it is not assumed, i.e. checkint ==
3569  * FALSE, then only the convex relaxations of the subproblems are solved. If integer feasibility is assumed, i.e.
3570  * checkint == TRUE, then the convex relaxations and the full CIP are solved to generate Benders' cuts and check
3571  * solution feasibility.
3572  *
3573  * TODO: consider allowing the possibility to pass solution information back from the subproblems instead of the scip
3574  * instance. This would allow the use of different solvers for the subproblems, more importantly allowing the use of an
3575  * LP solver for LP subproblems.
3576  */
3578  SCIP_BENDERS* benders, /**< Benders' decomposition */
3579  SCIP_SET* set, /**< global SCIP settings */
3580  SCIP_SOL* sol, /**< primal CIP solution */
3581  SCIP_RESULT* result, /**< result of the pricing process */
3582  SCIP_Bool* infeasible, /**< is the master problem infeasible with respect to the Benders' cuts? */
3583  SCIP_Bool* auxviol, /**< set to TRUE only if the solution is feasible but the aux vars are violated */
3584  SCIP_BENDERSENFOTYPE type, /**< the type of solution being enforced */
3585  SCIP_Bool checkint /**< should the integer solution be checked by the subproblems */
3586  )
3587 {
3588  int nsubproblems;
3589  int subproblemcount;
3590  int nsolveloops;
3591  int nverified;
3592  int nsolved;
3593  int* mergecands;
3594  int npriomergecands;
3595  int nmergecands;
3596  int* solveidx;
3597  int* executedidx;
3598  int nsolveidx;
3599  int nexecutedidx;
3600  int nfree;
3601  SCIP_Bool* subprobsolved;
3602  SCIP_BENDERSSUBSTATUS* substatus;
3603  SCIP_Bool optimal;
3604  SCIP_Bool allverified;
3605  SCIP_Bool success;
3606  SCIP_Bool stopped;
3607  int i;
3608  int l;
3609 
3610  success = TRUE;
3611  stopped = FALSE;
3612 
3613  SCIPsetDebugMsg(set, "Starting Benders' decomposition subproblem solving. type %d checkint %u\n", type, checkint);
3614 
3615 #ifdef SCIP_MOREDEBUG
3616  SCIP_CALL( SCIPprintSol(set->scip, sol, NULL, FALSE) );
3617 #endif
3618 
3619  /* start timing */
3620  SCIPclockStart(benders->bendersclock, set);
3621 
3622  nsubproblems = SCIPbendersGetNSubproblems(benders);
3623 
3624  (*auxviol) = FALSE;
3625  (*infeasible) = FALSE;
3626 
3627  /* It is assumed that the problem is optimal, until a subproblem is found not to be optimal. However, not all
3628  * subproblems could be checked in each iteration. As such, it is not possible to state that the problem is optimal
3629  * if not all subproblems are checked. Situations where this may occur is when a subproblem is a MIP and only the LP
3630  * is solved. Also, in a distributed computation, then it may be advantageous to only solve some subproblems before
3631  * resolving the master problem. As such, for a problem to be optimal, then (optimal && allverified) == TRUE
3632  */
3633  optimal = TRUE;
3634  nverified = 0;
3635  nsolved = 0;
3636 
3637  assert(benders != NULL);
3638  assert(result != NULL);
3639  assert(infeasible != NULL);
3640  assert(auxviol != NULL);
3641 
3642  /* if the Benders' decomposition is called from a sub-SCIP and the sub-SCIPs have been deactivated, then it is
3643  * assumed that this is an LNS heuristic. As such, the check is not performed and the solution is assumed to be
3644  * feasible
3645  */
3646  if( benders->iscopy && set->subscipsoff
3647  && (!benders->lnscheck
3648  || (benders->lnsmaxdepth > -1 && SCIPgetDepth(benders->sourcescip) >= benders->lnsmaxdepth)
3649  || (benders->lnsmaxcalls > -1 && SCIPbendersGetNCalls(benders) >= benders->lnsmaxcalls)
3650  || (type != SCIP_BENDERSENFOTYPE_CHECK && SCIPgetDepth(set->scip) == 0 && benders->lnsmaxcallsroot > -1
3651  && SCIPbendersGetNCalls(benders) >= benders->lnsmaxcallsroot)) )
3652  {
3653  (*result) = SCIP_DIDNOTRUN;
3654  return SCIP_OKAY;
3655  }
3656 
3657  /* it is not necessary to check all primal solutions by solving the Benders' decomposition subproblems.
3658  * Only the improving solutions are checked to improve efficiency of the algorithm.
3659  * If the solution is non-improving, the result FEASIBLE is returned. While this may be incorrect w.r.t to the
3660  * Benders' subproblems, this solution will never be the optimal solution. A non-improving solution may be used
3661  * within LNS primal heuristics. If this occurs, the improving solution, if found, will be checked by the solving
3662  * the Benders' decomposition subproblems.
3663  * TODO: Add a parameter to control this behaviour.
3664  */
3665  if( checkint && SCIPsetIsLE(set, SCIPgetPrimalbound(set->scip)*(int)SCIPgetObjsense(set->scip),
3666  SCIPgetSolOrigObj(set->scip, sol)*(int)SCIPgetObjsense(set->scip)) )
3667  {
3668  (*result) = SCIP_DIDNOTRUN;
3669  return SCIP_OKAY;
3670  }
3671 
3672  /* if the enforcement type is SCIP_BENDERSENFOTYPE_LP and the LP is currently unbounded. This could mean that there
3673  * is no lower bound on the auxiliary variables. In this case, we try to update the lower bound for the auxiliary
3674  * variables.
3675  */
3677  && benders->updateauxvarbound )
3678  {
3679  SCIP_CALL( updateAuxiliaryVarLowerbound(benders, set, result) );
3680 
3681  /* the auxiliary variable bound will only be updated once. */
3682  benders->updateauxvarbound = FALSE;
3683  }
3684 
3685  /* sets the stored objective function values of the subproblems to infinity */
3686  resetSubproblemObjectiveValue(benders, set);
3687 
3688  *result = SCIP_DIDNOTRUN;
3689 
3690  if( benders->benderspresubsolve != NULL && !benders->strengthenround )
3691  {
3692  SCIP_Bool skipsolve;
3693 
3694  skipsolve = FALSE;
3695  SCIP_CALL( benders->benderspresubsolve(set->scip, benders, sol, type, checkint, infeasible, auxviol, &skipsolve,
3696  result) );
3697 
3698  /* evaluate result */
3699  if( (*result) != SCIP_DIDNOTRUN
3700  && (*result) != SCIP_FEASIBLE
3701  && (*result) != SCIP_INFEASIBLE
3702  && (*result) != SCIP_CONSADDED
3703  && (*result) != SCIP_SEPARATED )
3704  {
3705  SCIPerrorMessage("the user-defined pre subproblem solving method for the Benders' decomposition <%s> returned "
3706  "invalid result <%d>\n", benders->name, *result);
3707  return SCIP_INVALIDRESULT;
3708  }
3709 
3710  /* if the solve must be skipped, then the solving loop is exited and the user defined result is returned */
3711  if( skipsolve )
3712  {
3713  SCIPsetDebugMsg(set, "skipping the subproblem solving for Benders' decomposition <%s>. "
3714  "returning result <%d>\n", benders->name, *result);
3715  return SCIP_OKAY;
3716  }
3717  }
3718 
3719  /* the cut strengthening is performed before the regular subproblem solve is called. To avoid recursion, the flag
3720  * strengthenround is set to TRUE when the cut strengthening is performed. The cut strengthening is not performed as
3721  * part of the large neighbourhood Benders' search.
3722  *
3723  * NOTE: cut strengthening is only applied for fractional solutions and integer solutions if there are no CIP
3724  * subproblems.
3725  */
3726  if( benders->strengthenenabled && !benders->strengthenround && !benders->iscopy
3727  && (!checkint || SCIPbendersGetNConvexSubproblems(benders) == SCIPbendersGetNSubproblems(benders)) )
3728  {
3729  SCIP_Bool skipsolve;
3730 
3731  benders->strengthenround = TRUE;
3732  /* if the user has not requested the solve to be skipped, then the cut strengthening is performed */
3733  SCIP_CALL( performInteriorSolCutStrengthening(benders, set, sol, type, checkint, FALSE, infeasible, auxviol,
3734  &skipsolve, result) );
3735  benders->strengthenround = FALSE;
3736 
3737  /* if the solve must be skipped, then the solving loop is exited and the user defined result is returned */
3738  if( skipsolve )
3739  {
3740  SCIPsetDebugMsg(set, "skipping the subproblem solving because cut strengthening found a cut "
3741  "for Benders' decomposition <%s>. Returning result <%d>\n", benders->name, *result);
3742  return SCIP_OKAY;
3743  }
3744 
3745  /* the result flag need to be reset to DIDNOTRUN for the main subproblem solve */
3746  (*result) = SCIP_DIDNOTRUN;
3747  }
3748 
3749  /* allocating memory for the infeasible subproblem array */
3750  SCIP_CALL( SCIPallocClearBlockMemoryArray(set->scip, &subprobsolved, nsubproblems) );
3751  SCIP_CALL( SCIPallocClearBlockMemoryArray(set->scip, &substatus, nsubproblems) );
3752  SCIP_CALL( SCIPallocClearBlockMemoryArray(set->scip, &mergecands, nsubproblems) );
3753  npriomergecands = 0;
3754  nmergecands = 0;
3755 
3756  /* allocating the memory for the subproblem solving and cut generation indices */
3757  SCIP_CALL( SCIPallocClearBlockMemoryArray(set->scip, &solveidx, nsubproblems) );
3758  SCIP_CALL( SCIPallocClearBlockMemoryArray(set->scip, &executedidx, nsubproblems) );
3759  nsolveidx = 0;
3760  nexecutedidx = 0;
3761 
3762  /* only a subset of the subproblems are initially solved. Both solving loops are executed for the subproblems to
3763  * check whether any cuts are generated. If a cut is generated, then no further subproblems are solved. If a cut is
3764  * not generated, then an additional set of subproblems are solved.
3765  */
3766  while( nsolved < nsubproblems )
3767  {
3768  /* getting the indices for the subproblems that will be solved */
3769  createSolveSubproblemIndexList(benders, set, type, &solveidx, &nsolveidx);
3770 
3771  /* by default the number of solve loops is 1. This is the case if all subproblems are LP or the user has defined a
3772  * benderssolvesub callback. If there is a subproblem that is not an LP, then 2 solve loops are performed. The first
3773  * loop is the LP solving loop, the second solves the subproblem to integer optimality.
3774  */
3775  nsolveloops = 1;
3776 
3777  for( l = 0; l < nsolveloops; l++ )
3778  {
3779  SCIP_BENDERSSOLVELOOP solveloop; /* identifies what problem type is solve in this solve loop */
3780 
3781  /* if either benderssolvesubconvex or benderssolvesub are implemented, then the user callbacks are invoked */
3782  if( benders->benderssolvesubconvex != NULL || benders->benderssolvesub != NULL )
3783  {
3784  if( l == 0 )
3786  else
3787  solveloop = SCIP_BENDERSSOLVELOOP_USERCIP;
3788  }
3789  else
3790  solveloop = (SCIP_BENDERSSOLVELOOP) l;
3791 
3792  /* solving the subproblems for this round of enforcement/checking. */
3793  SCIP_CALL( solveBendersSubproblems(benders, set, sol, type, solveloop, checkint, &nverified,
3794  solveidx, nsolveidx, &subprobsolved, &substatus, infeasible, &optimal, &stopped) );
3795 
3796  /* if the solving has been stopped, then the subproblem solving and cut generation must terminate */
3797  if( stopped )
3798  break;
3799 
3800  /* Generating cuts for the subproblems. Cuts are only generated when the solution is from primal heuristics,
3801  * relaxations or the LP
3802  */
3803  if( type != SCIP_BENDERSENFOTYPE_PSEUDO )
3804  {
3805  SCIP_CALL( generateBendersCuts(benders, set, sol, result, type, solveloop, checkint, subprobsolved,
3806  substatus, solveidx, nsolveidx, &mergecands, &npriomergecands, &nmergecands, &nsolveloops) );
3807  }
3808  else
3809  {
3810  /* The first solving loop solves the convex subproblems and the convex relaxations of the CIP subproblems. The
3811  * second solving loop solves the CIP subproblems. The second solving loop is only called if the integer
3812  * feasibility is being checked and if the convex subproblems and convex relaxations are not infeasible.
3813  */
3814  if( !(*infeasible) && checkint && !SCIPbendersOnlyCheckConvexRelax(benders, SCIPsetGetSubscipsOff(set))
3816  nsolveloops = 2;
3817  }
3818  }
3819 
3820  nsolved += nsolveidx;
3821 
3822  /* storing the indices of the subproblems for which the solving loop was executed */
3823  for( i = 0; i < nsolveidx; i++ )
3824  executedidx[nexecutedidx++] = solveidx[i];
3825 
3826  /* if the result is CONSADDED or SEPARATED, then a cut is generated and no further subproblem processing is
3827  * required
3828  */
3829  if( (*result) == SCIP_CONSADDED || (*result) == SCIP_SEPARATED )
3830  break;
3831  }
3832 
3833  /* inserting the subproblems into the priority queue for the next solve call */
3834  SCIP_CALL( updateSubproblemStatQueue(benders, executedidx, nexecutedidx, TRUE) );
3835 
3836  if( stopped )
3837  goto TERMINATE;
3838 
3839  allverified = (nverified == nsubproblems);
3840 
3841  SCIPsetDebugMsg(set, "End Benders' decomposition subproblem solve. result %d infeasible %u auxviol %u nverified %d\n",
3842  *result, *infeasible, *auxviol, nverified);
3843 
3844 #ifdef SCIP_DEBUG
3845  if( (*result) == SCIP_CONSADDED )
3846  {
3847  SCIPsetDebugMsg(set, "Benders' decomposition: Cut added\n");
3848  }
3849 #endif
3850 
3851  /* if the number of checked pseudo solutions exceeds a set limit, then all subproblems are passed as merge
3852  * candidates. Currently, merging subproblems into the master problem is the only method for resolving numerical
3853  * troubles.
3854  *
3855  * We are only interested in the pseudo solutions that have been checked completely for integrality. This is
3856  * identified by checkint == TRUE. This means that the Benders' decomposition constraint is one of the last
3857  * constraint handlers that must resolve the infeasibility. If the Benders' decomposition framework can't resolve the
3858  * infeasibility, then this will result in an error.
3859  */
3860  if( type == SCIP_BENDERSENFOTYPE_PSEUDO && checkint )
3861  {
3862  benders->npseudosols++;
3863 
3864  if( benders->npseudosols > BENDERS_MAXPSEUDOSOLS )
3865  {
3866  /* if a priority merge candidate already exists, then no other merge candidates need to be added.*/
3867  if( npriomergecands == 0 )
3868  {
3869  /* all subproblems are added to the merge candidate list. The first active subproblem is added as a
3870  * priority merge candidate
3871  */
3872  nmergecands = 0;
3873  npriomergecands = 1;
3874  for( i = 0; i < nsubproblems; i++ )
3875  {
3876  /* only active subproblems are added to the merge candidate list */
3877  if( subproblemIsActive(benders, i) )
3878  {
3879  mergecands[nmergecands] = i;
3880  nmergecands++;
3881  }
3882  }
3883 
3884  SCIPverbMessage(set->scip, SCIP_VERBLEVEL_HIGH, NULL, " The number of checked pseudo solutions exceeds the "
3885  "limit of %d. All active subproblems are merge candidates, with subproblem %d a priority candidate.\n",
3886  BENDERS_MAXPSEUDOSOLS, mergecands[0]);
3887  }
3888  }
3889  }
3890  else
3891  benders->npseudosols = 0;
3892 
3893  /* if the result is SCIP_DIDNOTFIND, then there was a error in generating cuts in all subproblems that are not
3894  * optimal. This result does not cutoff any solution, so the Benders' decomposition algorithm will fail.
3895  *
3896  * It could happen that the cut strengthening approach causes an error the cut generation. In this case, an error
3897  * should not be thrown. So, this check will be skipped when in a strengthening round.
3898  * TODO: Work out a way to ensure Benders' decomposition does not terminate due to a SCIP_DIDNOTFIND result.
3899  */
3900  if( (*result) == SCIP_DIDNOTFIND && !benders->strengthenround )
3901  {
3902  if( type == SCIP_BENDERSENFOTYPE_PSEUDO )
3903  (*result) = SCIP_SOLVELP;
3904  else
3905  (*result) = SCIP_INFEASIBLE;
3906 
3907  SCIPerrorMessage("An error was found when generating cuts for non-optimal subproblems of Benders' "
3908  "decomposition <%s>. Consider merging the infeasible subproblems into the master problem.\n", SCIPbendersGetName(benders));
3909 
3910  /* since no other cuts are generated, then this error will result in a crash. It is possible to avoid the error,
3911  * by merging the affected subproblem into the master problem.
3912  *
3913  * NOTE: If the error occurs while checking solutions, i.e. SCIP_BENDERSENFOTYPE_CHECK, then it is valid to set
3914  * the result to SCIP_INFEASIBLE and the success flag to TRUE
3915  */
3916  if( type != SCIP_BENDERSENFOTYPE_CHECK )
3917  success = FALSE;
3918 
3919  goto POSTSOLVE;
3920  }
3921 
3922  if( type == SCIP_BENDERSENFOTYPE_PSEUDO )
3923  {
3924  if( (*infeasible) || !allverified )
3925  (*result) = SCIP_SOLVELP;
3926  else
3927  {
3928  (*result) = SCIP_FEASIBLE;
3929 
3930  /* if the subproblems are not infeasible, but they are also not optimal. This means that there is a violation
3931  * in the auxiliary variable values. In this case, a feasible result is returned with the auxviol flag set to
3932  * TRUE.
3933  */
3934  (*auxviol) = !optimal;
3935  }
3936  }
3937  else if( checkint && (type == SCIP_BENDERSENFOTYPE_CHECK
3938  || ((*result) != SCIP_CONSADDED && (*result) != SCIP_SEPARATED)) )
3939  {
3940  /* if the subproblems are being solved as part of conscheck, then the results flag must be returned after the solving
3941  * has completed.
3942  */
3943  if( (*infeasible) || !allverified )
3944  (*result) = SCIP_INFEASIBLE;
3945  else
3946  {
3947  (*result) = SCIP_FEASIBLE;
3948 
3949  /* if the subproblems are not infeasible, but they are also not optimal. This means that there is a violation
3950  * in the auxiliary variable values. In this case, a feasible result is returned with the auxviol flag set to
3951  * TRUE.
3952  */
3953  (*auxviol) = !optimal;
3954  }
3955  }
3956 
3957 POSTSOLVE:
3958  /* calling the post-solve call back for the Benders' decomposition algorithm. This allows the user to work directly
3959  * with the solved subproblems and the master problem */
3960  if( benders->benderspostsolve != NULL )
3961  {
3962  SCIP_Bool merged;
3963 
3964  merged = FALSE;
3965 
3966  SCIP_CALL( benders->benderspostsolve(set->scip, benders, sol, type, mergecands, npriomergecands, nmergecands,
3967  checkint, (*infeasible), &merged) );
3968 
3969  if( merged )
3970  {
3971  (*result) = SCIP_CONSADDED;
3972 
3973  /* since subproblems have been merged, then constraints have been added. This could resolve the unresolved
3974  * infeasibility, so the error has been corrected.
3975  */
3976  success = TRUE;
3977  }
3978  else if( !success )
3979  {
3980  SCIPerrorMessage("An error occurred during Benders' decomposition cut generations and no merging had been "
3981  "performed. It is not possible to continue solving the problem by Benders' decomposition\n");
3982  }
3983  }
3984 
3985 TERMINATE:
3986  /* if the solving process has stopped, then all subproblems need to be freed */
3987  if( stopped )
3988  nfree = nsubproblems;
3989  else
3990  nfree = nexecutedidx;
3991 
3992  /* freeing the subproblems after the cuts are generated */
3993  subproblemcount = 0;
3994  while( subproblemcount < nfree )
3995  {
3996  int subidx;
3997 
3998  if( stopped )
3999  subidx = subproblemcount;
4000  else
4001  subidx = executedidx[subproblemcount];
4002 
4003  SCIP_CALL( SCIPbendersFreeSubproblem(benders, set, subidx) );
4004 
4005  subproblemcount++;
4006  }
4007 
4008 #ifndef NDEBUG
4009  for( i = 0; i < nsubproblems; i++ )
4010  assert(SCIPbendersSubproblem(benders, i) == NULL
4012  || !SCIPinProbing(SCIPbendersSubproblem(benders, i))
4013  || !subproblemIsActive(benders, i));
4014 #endif
4015 
4016  /* increment the number of calls to the Benders' decomposition subproblem solve */
4017  benders->ncalls++;
4018 
4019  SCIPsetDebugMsg(set, "End Benders' decomposition execution method. result %d infeasible %u auxviol %u\n", *result,
4020  *infeasible, *auxviol);
4021 
4022  /* end timing */
4023  SCIPclockStop(benders->bendersclock, set);
4024 
4025  /* freeing memory */
4026  SCIPfreeBlockMemoryArray(set->scip, &executedidx, nsubproblems);
4027  SCIPfreeBlockMemoryArray(set->scip, &solveidx, nsubproblems);
4028  SCIPfreeBlockMemoryArray(set->scip, &mergecands, nsubproblems);
4029  SCIPfreeBlockMemoryArray(set->scip, &substatus, nsubproblems);
4030  SCIPfreeBlockMemoryArray(set->scip, &subprobsolved, nsubproblems);
4031 
4032  /* if there was an error in generating cuts and merging was not performed, then the solution is perturbed in an
4033  * attempt to generate a cut and correct the infeasibility
4034  */
4035  if( !success && !stopped )
4036  {
4037  SCIP_Bool skipsolve;
4038  SCIP_RESULT perturbresult;
4039 
4040  skipsolve = FALSE;
4041 
4042  benders->strengthenround = TRUE;
4043  /* if the user has not requested the solve to be skipped, then the cut strengthening is performed */
4044  SCIP_CALL( performInteriorSolCutStrengthening(benders, set, sol, type, checkint, TRUE, infeasible, auxviol,
4045  &skipsolve, &perturbresult) );
4046  benders->strengthenround = FALSE;
4047 
4048  if( perturbresult == SCIP_CONSADDED || perturbresult == SCIP_SEPARATED )
4049  (*result) = perturbresult;
4050 
4051  success = skipsolve;
4052  }
4053 
4054  /* if the Benders' decomposition subproblem check stopped, then we don't have a valid result. In this case, the
4055  * safest thing to do is report INFEASIBLE.
4056  */
4057  if( stopped )
4058  (*result) = SCIP_INFEASIBLE;
4059 
4060  /* if the subproblem verification identifies the solution as feasible, then a check whether slack variables have been
4061  * used is necessary. If any slack variables are non-zero, then the solution is reverified after the objective
4062  * coefficient for the slack variables is increased.
4063  */
4064  if( (*result) == SCIP_FEASIBLE )
4065  {
4066  SCIP_Bool activeslack;
4067 
4068  SCIP_CALL( SCIPbendersSolSlackVarsActive(benders, &activeslack) );
4069  SCIPsetDebugMsg(set, "Type: %d Active slack: %u Feasibility Phase: %u\n", type, activeslack,
4070  benders->feasibilityphase);
4071  if( activeslack )
4072  {
4073  if( type == SCIP_BENDERSENFOTYPE_CHECK )
4074  (*result) = SCIP_INFEASIBLE;
4075  else
4076  {
4077  /* increasing the value of the slack variable by a factor of 10 */
4078  benders->slackvarcoef *= 10.0;
4079 
4080  if( benders->slackvarcoef <= benders->maxslackvarcoef )
4081  {
4082  SCIPmessagePrintVerbInfo(SCIPgetMessagehdlr(set->scip), set->disp_verblevel, SCIP_VERBLEVEL_HIGH, "Increasing the slack variable coefficient to %g.\n", benders->slackvarcoef);
4083  }
4084  else
4085  {
4086  SCIPmessagePrintVerbInfo(SCIPgetMessagehdlr(set->scip), set->disp_verblevel, SCIP_VERBLEVEL_HIGH, "Fixing the slack variables to zero.\n");
4087  }
4088 
4089  /* resolving the subproblems with an increased slack variable */
4090  SCIP_CALL( SCIPsolveBendersSubproblems(set->scip, benders, sol, result, infeasible, auxviol, type, checkint) );
4091  }
4092  }
4093  else if( benders->feasibilityphase )
4094  {
4095  if( type != SCIP_BENDERSENFOTYPE_CHECK )
4096  {
4097  /* disabling the feasibility phase */
4098  benders->feasibilityphase = FALSE;
4099 
4100  /* resolving the subproblems with the slack variables set to zero */
4101  SCIP_CALL( SCIPsolveBendersSubproblems(set->scip, benders, sol, result, infeasible, auxviol, type, checkint) );
4102  }
4103  }
4104  }
4105 
4106  if( !success )
4107  return SCIP_ERROR;
4108  else
4109  return SCIP_OKAY;
4110 }
4111 
4112 /** solves the user-defined subproblem solving function */
4113 static
4115  SCIP_BENDERS* benders, /**< Benders' decomposition */
4116  SCIP_SET* set, /**< global SCIP settings */
4117  SCIP_SOL* sol, /**< primal CIP solution */
4118  int probnumber, /**< the subproblem number */
4119  SCIP_BENDERSSOLVELOOP solveloop, /**< the solve loop iteration. The first iter is for LP, the second for IP */
4120  SCIP_Bool* infeasible, /**< returns whether the current subproblem is infeasible */
4121  SCIP_Real* objective, /**< the objective function value of the subproblem */
4122  SCIP_RESULT* result /**< the result from solving the subproblem */
4123  )
4124 {
4125  assert(benders != NULL);
4126  assert(probnumber >= 0 && probnumber < benders->nsubproblems);
4127  assert(benders->benderssolvesubconvex != NULL || benders->benderssolvesub != NULL);
4128 
4129  assert(solveloop == SCIP_BENDERSSOLVELOOP_USERCONVEX || solveloop == SCIP_BENDERSSOLVELOOP_USERCIP);
4130 
4131  (*objective) = -SCIPsetInfinity(set);
4132 
4133  /* calls the user defined subproblem solving method. Only the convex relaxations are solved during the Large
4134  * Neighbourhood Benders' Search. */
4135  if( solveloop == SCIP_BENDERSSOLVELOOP_USERCONVEX )
4136  {
4137  if( benders->benderssolvesubconvex != NULL )
4138  {
4139  SCIP_CALL( benders->benderssolvesubconvex(set->scip, benders, sol, probnumber,
4140  SCIPbendersOnlyCheckConvexRelax(benders, SCIPsetGetSubscipsOff(set)), objective, result) );
4141  }
4142  else
4143  (*result) = SCIP_DIDNOTRUN;
4144  }
4145  else if( solveloop == SCIP_BENDERSSOLVELOOP_USERCIP )
4146  {
4147  if( benders->benderssolvesub != NULL )
4148  {
4149  SCIP_CALL( benders->benderssolvesub(set->scip, benders, sol, probnumber, objective, result) );
4150  }
4151  else
4152  (*result) = SCIP_DIDNOTRUN;
4153  }
4154 
4155  /* evaluate result */
4156  if( (*result) != SCIP_DIDNOTRUN
4157  && (*result) != SCIP_FEASIBLE
4158  && (*result) != SCIP_INFEASIBLE
4159  && (*result) != SCIP_UNBOUNDED )
4160  {
4161  SCIPerrorMessage("the user-defined solving method for the Benders' decomposition <%s> returned invalid result <%d>\n",
4162  benders->name, *result);
4163  return SCIP_INVALIDRESULT;
4164  }
4165 
4166  if( (*result) == SCIP_INFEASIBLE )
4167  (*infeasible) = TRUE;
4168 
4169  if( (*result) == SCIP_FEASIBLE
4170  && (SCIPsetIsInfinity(set, -(*objective)) || SCIPsetIsInfinity(set, (*objective))) )
4171  {
4172  SCIPerrorMessage("the user-defined solving method for the Benders' decomposition <%s> returned objective value %g\n",
4173  benders->name, (*objective));
4174  return SCIP_ERROR;
4175  }
4176 
4177  /* if the result is SCIP_DIDNOTFIND, then an error is returned and SCIP will terminate. */
4178  if( (*result) == SCIP_DIDNOTFIND )
4179  return SCIP_ERROR;
4180  else
4181  return SCIP_OKAY;
4182 }
4183 
4184 /** executes the subproblem solving process */
4186  SCIP_BENDERS* benders, /**< Benders' decomposition */
4187  SCIP_SET* set, /**< global SCIP settings */
4188  SCIP_SOL* sol, /**< primal CIP solution */
4189  int probnumber, /**< the subproblem number */
4190  SCIP_BENDERSSOLVELOOP solveloop, /**< the solve loop iteration. The first iter is for LP, the second for IP */
4191  SCIP_Bool enhancement, /**< is the solve performed as part of and enhancement? */
4192  SCIP_Bool* solved, /**< flag to indicate whether the subproblem was solved */
4193  SCIP_Bool* infeasible, /**< returns whether the current subproblem is infeasible */
4194  SCIP_BENDERSENFOTYPE type /**< the enforcement type calling this function */
4195  )
4196 { /*lint --e{715}*/
4197  SCIP* subproblem;
4198  SCIP_RESULT result;
4199  SCIP_Real objective;
4200  SCIP_STATUS solvestatus = SCIP_STATUS_UNKNOWN;
4201 
4202  assert(benders != NULL);
4203  assert(probnumber >= 0 && probnumber < benders->nsubproblems);
4204 
4205  SCIPsetDebugMsg(set, "Benders' decomposition: solving subproblem %d\n", probnumber);
4206 
4207  result = SCIP_DIDNOTRUN;
4208  objective = SCIPsetInfinity(set);
4209 
4210  subproblem = SCIPbendersSubproblem(benders, probnumber);
4211 
4212  if( subproblem == NULL && (benders->benderssolvesubconvex == NULL || benders->benderssolvesub == NULL) )
4213  {
4214  SCIPerrorMessage("The subproblem %d is set to NULL, but both bendersSolvesubconvex%s and bendersSolvesub%s "
4215  "are not defined.\n", probnumber, benders->name, benders->name);
4216  SCIPABORT();
4217  return SCIP_ERROR;
4218  }
4219 
4220  /* initially setting the solved flag to FALSE */
4221  (*solved) = FALSE;
4222 
4223  /* if the subproblem solve callback is implemented, then that is used instead of the default setup */
4224  if( solveloop == SCIP_BENDERSSOLVELOOP_USERCONVEX || solveloop == SCIP_BENDERSSOLVELOOP_USERCIP )
4225  {
4226  /* calls the user defined subproblem solving method. Only the convex relaxations are solved during the Large
4227  * Neighbourhood Benders' Search. */
4228  SCIP_CALL( executeUserDefinedSolvesub(benders, set, sol, probnumber, solveloop, infeasible, &objective, &result) );
4229 
4230  /* if the result is DIDNOTRUN, then the subproblem was not solved */
4231  (*solved) = (result != SCIP_DIDNOTRUN);
4232  }
4233  else if( subproblem != NULL )
4234  {
4235  /* setting up the subproblem */
4236  if( solveloop == SCIP_BENDERSSOLVELOOP_CONVEX )
4237  {
4238  SCIP_CALL( SCIPbendersSetupSubproblem(benders, set, sol, probnumber, type) );
4239 
4240  /* if the limits of the master problem were hit during the setup process, then the subproblem will not have
4241  * been setup. In this case, the solving function must be exited.
4242  */
4243  if( !SCIPbendersSubproblemIsSetup(benders, probnumber) )
4244  {
4245  SCIPbendersSetSubproblemObjval(benders, probnumber, SCIPsetInfinity(set));
4246  (*solved) = FALSE;
4247  return SCIP_OKAY;
4248  }
4249  }
4250  else
4251  {
4252  SCIP_CALL( updateEventhdlrUpperbound(benders, probnumber, SCIPbendersGetAuxiliaryVarVal(benders, set, sol, probnumber)) );
4253  }
4254 
4255  /* solving the subproblem
4256  * the LP of the subproblem is solved in the first solveloop.
4257  * In the second solve loop, the MIP problem is solved */
4258  if( solveloop == SCIP_BENDERSSOLVELOOP_CONVEX
4260  {
4261  SCIP_CALL( SCIPbendersSolveSubproblemLP(set->scip, benders, probnumber, &solvestatus, &objective) );
4262 
4263  /* if the (N)LP was solved without error, then the subproblem is labelled as solved */
4264  if( solvestatus == SCIP_STATUS_OPTIMAL || solvestatus == SCIP_STATUS_INFEASIBLE )
4265  (*solved) = TRUE;
4266 
4267  if( solvestatus == SCIP_STATUS_INFEASIBLE )
4268  (*infeasible) = TRUE;
4269  }
4270  else
4271  {
4272  SCIP_SOL* bestsol;
4273 
4274  SCIP_CALL( SCIPbendersSolveSubproblemCIP(set->scip, benders, probnumber, &solvestatus, FALSE) );
4275 
4276  if( solvestatus == SCIP_STATUS_INFEASIBLE )
4277  (*infeasible) = TRUE;
4278 
4279  /* if the generic subproblem solving methods are used, then the CIP subproblems are always solved. */
4280  (*solved) = TRUE;
4281 
4282  bestsol = SCIPgetBestSol(subproblem);
4283  if( bestsol != NULL )
4284  objective = SCIPgetSolOrigObj(subproblem, bestsol)*(int)SCIPgetObjsense(set->scip);
4285  else
4286  objective = SCIPsetInfinity(set);
4287  }
4288  }
4289  else
4290  {
4291  SCIPABORT();
4292  }
4293 
4294  if( !enhancement )
4295  {
4296  /* The following handles the cases when the subproblem is OPTIMAL, INFEASIBLE and UNBOUNDED.
4297  * If a subproblem is unbounded, then the auxiliary variables are set to -infinity and the unbounded flag is
4298  * returned as TRUE. No cut will be generated, but the result will be set to SCIP_FEASIBLE.
4299  */
4300  if( solveloop == SCIP_BENDERSSOLVELOOP_CONVEX || solveloop == SCIP_BENDERSSOLVELOOP_CIP )
4301  {
4302  /* TODO: Consider whether other solutions status should be handled */
4303  if( solvestatus == SCIP_STATUS_OPTIMAL )
4304  SCIPbendersSetSubproblemObjval(benders, probnumber, objective);
4305  else if( solvestatus == SCIP_STATUS_INFEASIBLE )
4306  SCIPbendersSetSubproblemObjval(benders, probnumber, SCIPsetInfinity(set));
4307  else if( solvestatus == SCIP_STATUS_USERINTERRUPT || solvestatus == SCIP_STATUS_BESTSOLLIMIT )
4308  SCIPbendersSetSubproblemObjval(benders, probnumber, objective);
4309  else if( solvestatus == SCIP_STATUS_MEMLIMIT || solvestatus == SCIP_STATUS_TIMELIMIT
4310  || solvestatus == SCIP_STATUS_UNKNOWN )
4311  {
4312  SCIPverbMessage(set->scip, SCIP_VERBLEVEL_FULL, NULL, " Benders' decomposition: Error solving "
4313  "subproblem %d. No cut will be generated for this subproblem.\n", probnumber);
4314  SCIPbendersSetSubproblemObjval(benders, probnumber, SCIPsetInfinity(set));
4315  }
4316  else if( solvestatus == SCIP_STATUS_UNBOUNDED )
4317  {
4318  SCIPerrorMessage("The Benders' decomposition subproblem %d is unbounded. This should not happen.\n",
4319  probnumber);
4320  SCIPABORT();
4321  }
4322  else
4323  {
4324  SCIPerrorMessage("Invalid status returned from solving Benders' decomposition subproblem %d. Solution status: %d\n",
4325  probnumber, solvestatus);
4326  SCIPABORT();
4327  }
4328  }
4329  else
4330  {
4331  assert(solveloop == SCIP_BENDERSSOLVELOOP_USERCONVEX || solveloop == SCIP_BENDERSSOLVELOOP_USERCIP);
4332  if( result == SCIP_FEASIBLE )
4333  SCIPbendersSetSubproblemObjval(benders, probnumber, objective);
4334  else if( result == SCIP_INFEASIBLE )
4335  SCIPbendersSetSubproblemObjval(benders, probnumber, SCIPsetInfinity(set));
4336  else if( result == SCIP_UNBOUNDED )
4337  {
4338  SCIPerrorMessage("The Benders' decomposition subproblem %d is unbounded. This should not happen.\n",
4339  probnumber);
4340  SCIPABORT();
4341  }
4342  else if( result != SCIP_DIDNOTRUN )
4343  {
4344  SCIPerrorMessage("Invalid result <%d> from user-defined subproblem solving method. This should not happen.\n",
4345  result);
4346  }
4347  }
4348  }
4349 
4350  return SCIP_OKAY;
4351 }
4352 
4353 /** sets up the subproblem using the solution to the master problem */
4355  SCIP_BENDERS* benders, /**< Benders' decomposition */
4356  SCIP_SET* set, /**< global SCIP settings */
4357  SCIP_SOL* sol, /**< primal CIP solution */
4358  int probnumber, /**< the subproblem number */
4359  SCIP_BENDERSENFOTYPE type /**< the enforcement type calling this function */
4360  )
4361 {
4362  SCIP* subproblem;
4363  SCIP_VAR** vars;
4364  SCIP_VAR* mastervar;
4365  SCIP_Real solval;
4366  int nvars;
4367  int i;
4368 
4369  assert(benders != NULL);
4370  assert(set != NULL);
4371  assert(probnumber >= 0 && probnumber < SCIPbendersGetNSubproblems(benders));
4372 
4373  subproblem = SCIPbendersSubproblem(benders, probnumber);
4374 
4375  /* the subproblem setup can only be performed if the subproblem is not NULL */
4376  if( subproblem == NULL )
4377  {
4378  SCIPerrorMessage("The subproblem %d is NULL. Thus, the subproblem setup must be performed manually in either "
4379  "bendersSolvesubconvex%s or bendersSolvesub%s.\n", probnumber, benders->name, benders->name);
4380  return SCIP_ERROR;
4381  }
4382  assert(subproblem != NULL);
4383 
4384  /* changing all of the master problem variable to continuous. */
4385  SCIP_CALL( SCIPbendersChgMastervarsToCont(benders, set, probnumber) );
4386 
4387  /* if the Benders' decomposition subproblem is convex and has continuous variables, then probing mode
4388  * must be started.
4389  * If the subproblem contains non-convex constraints or discrete variables, then the problem must be initialised,
4390  * and then put into SCIP_STAGE_SOLVING to be able to change the variable bounds. The probing mode is entered once
4391  * the variable bounds are set.
4392  * In the latter case, the transformed problem is freed after each subproblem solve round. */
4393  if( SCIPbendersGetSubproblemType(benders, probnumber) == SCIP_BENDERSSUBTYPE_CONVEXCONT )
4394  {
4395  SCIP_CALL( SCIPstartProbing(subproblem) );
4396  }
4397  else
4398  {
4399  SCIP_Bool success;
4400 
4401  SCIP_CALL( initialiseSubproblem(benders, set, probnumber, &success) );
4402 
4403  if( !success )
4404  {
4405  /* set the flag to indicate that the subproblems have been set up */
4406  SCIPbendersSetSubproblemIsSetup(benders, probnumber, FALSE);
4407 
4408  return SCIP_OKAY;
4409  }
4410  }
4411 
4412  vars = SCIPgetVars(subproblem);
4413  nvars = SCIPgetNVars(subproblem);
4414 
4415  /* looping over all variables in the subproblem to find those corresponding to the master problem variables. */
4416  /* TODO: It should be possible to store the pointers to the master variables to speed up the subproblem setup */
4417  for( i = 0; i < nvars; i++ )
4418  {
4419  SCIP_CALL( SCIPbendersGetVar(benders, set, vars[i], &mastervar, -1) );
4420 
4421  if( mastervar != NULL )
4422  {
4423  /* It is possible due to numerics that the solution value exceeds the upper or lower bounds. When this
4424  * happens, it causes an error in the LP solver as a result of inconsistent bounds. So the following statements
4425  * are used to ensure that the bounds are not exceeded when applying the fixings for the Benders'
4426  * decomposition subproblems
4427  */
4428  solval = SCIPgetSolVal(set->scip, sol, mastervar);
4429  if( !SCIPisLT(set->scip, solval, SCIPvarGetUbLocal(vars[i])) )
4430  solval = SCIPvarGetUbLocal(vars[i]);
4431  else if( !SCIPisGT(set->scip, solval, SCIPvarGetLbLocal(vars[i])) )
4432  solval = SCIPvarGetLbLocal(vars[i]);
4433 
4434  /* fixing the variable in the subproblem */
4435  if( !SCIPisEQ(subproblem, SCIPvarGetLbLocal(vars[i]), SCIPvarGetUbLocal(vars[i])) )
4436  {
4437  if( SCIPisGT(subproblem, solval, SCIPvarGetLbLocal(vars[i])) )
4438  {
4439  SCIP_CALL( SCIPchgVarLb(subproblem, vars[i], solval) );
4440  }
4441  if( SCIPisLT(subproblem, solval, SCIPvarGetUbLocal(vars[i])) )
4442  {
4443  SCIP_CALL( SCIPchgVarUb(subproblem, vars[i], solval) );
4444  }
4445  }
4446 
4447  assert(SCIPisEQ(subproblem, SCIPvarGetLbLocal(vars[i]), SCIPvarGetUbLocal(vars[i])));
4448  }
4449  else if( strstr(SCIPvarGetName(vars[i]), SLACKVAR_NAME) != NULL )
4450  {
4451  /* if the slack variables have been added to help improve feasibility, then they remain unfixed with a large
4452  * objective coefficient. Once the root node has been solved to optimality, then the slack variables are
4453  * fixed to zero.
4454  */
4455  if( benders->feasibilityphase && SCIPgetDepth(set->scip) == 0 && type != SCIP_BENDERSENFOTYPE_CHECK )
4456  {
4457  /* The coefficient update or variable fixing can only be performed if the subproblem is in probing mode.
4458  * If the slack var coef gets very large, then we fix the slack variable to 0 instead.
4459  */
4460  if( SCIPinProbing(subproblem) )
4461  {
4462  if( benders->slackvarcoef <= benders->maxslackvarcoef )
4463  {
4464  SCIP_CALL( SCIPchgVarObjProbing(subproblem, vars[i], benders->slackvarcoef) );
4465  }
4466  else
4467  {
4468  SCIP_CALL( SCIPchgVarUbProbing(subproblem, vars[i], 0.0) );
4469  }
4470  }
4471  }
4472  else
4473  {
4474  /* if the subproblem is non-linear and convex, then slack variables have been added to the subproblem. These
4475  * need to be fixed to zero when first solving the subproblem. However, if the slack variables have been added
4476  * by setting the execfeasphase runtime parameter, then they must not get fixed to zero
4477  */
4478  assert( !SCIPisEQ(subproblem, SCIPvarGetLbLocal(vars[i]), SCIPvarGetUbLocal(vars[i])) );
4479  assert( SCIPisZero(subproblem, SCIPvarGetLbLocal(vars[i])) );
4480 
4481  if( SCIPisLT(subproblem, 0.0, SCIPvarGetUbLocal(vars[i])) )
4482  {
4483  SCIP_CALL( SCIPchgVarUb(subproblem, vars[i], 0.0) );
4484  }
4485  }
4486  }
4487  }
4488 
4489  /* if the subproblem contain non-convex constraints or discrete variables, then the probing mode is entered after
4490  * setting up the subproblem
4491  */
4492  if( SCIPbendersGetSubproblemType(benders, probnumber) != SCIP_BENDERSSUBTYPE_CONVEXCONT )
4493  {
4494  SCIP_CALL( SCIPstartProbing(subproblem) );
4495  }
4496 
4497  /* set the flag to indicate that the subproblems have been set up */
4498  SCIPbendersSetSubproblemIsSetup(benders, probnumber, TRUE);
4499 
4500  return SCIP_OKAY;
4501 }
4502 
4503 /** Solve a Benders' decomposition subproblems. This will either call the user defined method or the generic solving
4504  * methods. If the generic method is called, then the subproblem must be set up before calling this method. */
4506  SCIP_BENDERS* benders, /**< Benders' decomposition */
4507  SCIP_SET* set, /**< global SCIP settings */
4508  SCIP_SOL* sol, /**< primal CIP solution, can be NULL */
4509  int probnumber, /**< the subproblem number */
4510  SCIP_Bool* infeasible, /**< returns whether the current subproblem is infeasible */
4511  SCIP_Bool solvecip, /**< directly solve the CIP subproblem */
4512  SCIP_Real* objective /**< the objective function value of the subproblem, can be NULL */
4513  )
4514 {
4515  assert(benders != NULL);
4516  assert(set != NULL);
4517  assert(probnumber >= 0 && probnumber < SCIPbendersGetNSubproblems(benders));
4518 
4519  (*infeasible) = FALSE;
4520 
4521  /* the subproblem must be set up before this function is called. */
4522  if( SCIPbendersSubproblem(benders, probnumber) != NULL && !SCIPbendersSubproblemIsSetup(benders, probnumber)
4523  && !SCIPbendersSubproblemIsIndependent(benders, probnumber) )
4524  {
4525  SCIPerrorMessage("Benders' decomposition subproblem %d must be set up before calling SCIPbendersSolveSubproblem(). Call SCIPsetupSubproblem() first.\n", probnumber);
4526  return SCIP_ERROR;
4527  }
4528 
4529  /* if the subproblem solve callback is implemented, then that is used instead of the default setup */
4530  if( benders->benderssolvesubconvex != NULL || benders->benderssolvesub != NULL)
4531  {
4532  SCIP_BENDERSSOLVELOOP solveloop;
4533  SCIP_RESULT result;
4534  SCIP_Real subobj;
4535 
4536  if( solvecip )
4537  solveloop = SCIP_BENDERSSOLVELOOP_USERCIP;
4538  else
4540 
4541  SCIP_CALL( executeUserDefinedSolvesub(benders, set, sol, probnumber, solveloop, infeasible, &subobj, &result) );
4542 
4543  if( objective != NULL )
4544  (*objective) = subobj;
4545  }
4546  else
4547  {
4548  SCIP* subproblem;
4549 
4550  subproblem = SCIPbendersSubproblem(benders, probnumber);
4551  assert(subproblem != NULL);
4552 
4553  /* solving the subproblem */
4554  if( solvecip && SCIPbendersGetSubproblemType(benders, probnumber) != SCIP_BENDERSSUBTYPE_CONVEXCONT )
4555  {
4556  SCIP_STATUS solvestatus;
4557 
4558  SCIP_CALL( SCIPbendersSolveSubproblemCIP(set->scip, benders, probnumber, &solvestatus, solvecip) );
4559 
4560  if( solvestatus == SCIP_STATUS_INFEASIBLE )
4561  (*infeasible) = TRUE;
4562  if( objective != NULL )
4563  (*objective) = SCIPgetSolOrigObj(subproblem, SCIPgetBestSol(subproblem))*(int)SCIPgetObjsense(subproblem);
4564  }
4565  else
4566  {
4567  SCIP_Bool success;
4568 
4569  /* if the subproblem has convex constraints and continuous variables, then it should have been initialised and
4570  * in SCIP_STAGE_SOLVING. In this case, the subproblem only needs to be put into probing mode.
4571  */
4572  if( SCIPbendersGetSubproblemType(benders, probnumber) == SCIP_BENDERSSUBTYPE_CONVEXCONT )
4573  {
4574  /* if the subproblem is not in probing mode, then it must be put into that mode for the LP solve. */
4575  if( !SCIPinProbing(subproblem) )
4576  {
4577  SCIP_CALL( SCIPstartProbing(subproblem) );
4578  }
4579 
4580  success = TRUE;
4581  }
4582  else
4583  {
4584  SCIP_CALL( initialiseSubproblem(benders, set, probnumber, &success) );
4585  }
4586 
4587  /* if setting up the subproblem was successful */
4588  if( success )
4589  {
4590  SCIP_STATUS solvestatus;
4591  SCIP_Real lpobjective;
4592 
4593  SCIP_CALL( SCIPbendersSolveSubproblemLP(set->scip, benders, probnumber, &solvestatus, &lpobjective) );
4594 
4595  if( solvestatus == SCIP_STATUS_INFEASIBLE )
4596  (*infeasible) = TRUE;
4597  else if( objective != NULL )
4598  (*objective) = lpobjective;
4599  }
4600  else
4601  {
4602  if( objective != NULL )
4603  (*objective) = SCIPinfinity(subproblem);
4604  }
4605  }
4606  }
4607 
4608  return SCIP_OKAY;
4609 }
4610 
4611 /** copies the time and memory limit from the master problem to the subproblem */
4612 static
4614  SCIP* scip, /**< the SCIP data structure */
4615  SCIP* subproblem /**< the Benders' decomposition subproblem */
4616  )
4617 {
4618  SCIP_Real mastertimelimit;
4619  SCIP_Real subtimelimit;
4620  SCIP_Real maxsubtimelimit;
4621  SCIP_Real mastermemorylimit;
4622  SCIP_Real submemorylimit;
4623  SCIP_Real maxsubmemorylimit;
4624 
4625  assert(scip != NULL);
4626 
4627  /* setting the time limit for the Benders' decomposition subproblems. It is set to 102% of the remaining time. */
4628  SCIP_CALL( SCIPgetRealParam(scip, "limits/time", &mastertimelimit) );
4629  maxsubtimelimit = SCIPparamGetRealMax(SCIPgetParam(subproblem, "limits/time"));
4630  subtimelimit = (mastertimelimit - SCIPgetSolvingTime(scip)) * 1.02;
4631  subtimelimit = MIN(subtimelimit, maxsubtimelimit);
4632  SCIP_CALL( SCIPsetRealParam(subproblem, "limits/time", MAX(0.0, subtimelimit)) );
4633 
4634  /* setting the memory limit for the Benders' decomposition subproblems. */
4635  SCIP_CALL( SCIPgetRealParam(scip, "limits/memory", &mastermemorylimit) );
4636  maxsubmemorylimit = SCIPparamGetRealMax(SCIPgetParam(subproblem, "limits/memory"));
4637  submemorylimit = mastermemorylimit - (SCIPgetMemUsed(scip) + SCIPgetMemExternEstim(scip))/1048576.0;
4638  submemorylimit = MIN(submemorylimit, maxsubmemorylimit);
4639  SCIP_CALL( SCIPsetRealParam(subproblem, "limits/memory", MAX(0.0, submemorylimit)) );
4640 
4641  return SCIP_OKAY;
4642 }
4643 
4644 /** stores the original parameters from the subproblem */
4645 static
4647  SCIP* subproblem, /**< the SCIP data structure */
4648  SCIP_SUBPROBPARAMS* origparams /**< the original subproblem parameters */
4649  )
4650 {
4651  assert(subproblem != NULL);
4652  assert(origparams != NULL);
4653 
4654  SCIP_CALL( SCIPgetRealParam(subproblem, "limits/memory", &origparams->limits_memory) );
4655  SCIP_CALL( SCIPgetRealParam(subproblem, "limits/time", &origparams->limits_time) );
4656  SCIP_CALL( SCIPgetBoolParam(subproblem, "conflict/enable", &origparams->conflict_enable) );
4657  SCIP_CALL( SCIPgetIntParam(subproblem, "lp/disablecutoff", &origparams->lp_disablecutoff) );
4658  SCIP_CALL( SCIPgetIntParam(subproblem, "lp/scaling", &origparams->lp_scaling) );
4659  SCIP_CALL( SCIPgetCharParam(subproblem, "lp/initalgorithm", &origparams->lp_initalg) );
4660  SCIP_CALL( SCIPgetCharParam(subproblem, "lp/resolvealgorithm", &origparams->lp_resolvealg) );
4661  SCIP_CALL( SCIPgetBoolParam(subproblem, "lp/alwaysgetduals", &origparams->lp_alwaysgetduals) );
4662  SCIP_CALL( SCIPgetBoolParam(subproblem, "misc/scaleobj", &origparams->misc_scaleobj) );
4663  SCIP_CALL( SCIPgetBoolParam(subproblem, "misc/catchctrlc", &origparams->misc_catchctrlc) );
4664  SCIP_CALL( SCIPgetIntParam(subproblem, "propagating/maxrounds", &origparams->prop_maxrounds) );
4665  SCIP_CALL( SCIPgetIntParam(subproblem, "propagating/maxroundsroot", &origparams->prop_maxroundsroot) );
4666  SCIP_CALL( SCIPgetIntParam(subproblem, "constraints/linear/propfreq", &origparams->cons_linear_propfreq) );
4667 
4668  return SCIP_OKAY;
4669 }
4670 
4671 /** sets the parameters for the subproblem */
4672 static
4674  SCIP* scip, /**< the SCIP data structure */
4675  SCIP* subproblem /**< the subproblem SCIP instance */
4676  )
4677 {
4678  assert(scip != NULL);
4679  assert(subproblem != NULL);
4680 
4681  /* copying memory and time limits */
4682  SCIP_CALL( copyMemoryAndTimeLimits(scip, subproblem) );
4683 
4684  /* Do we have to disable presolving? If yes, we have to store all presolving parameters. */
4686 
4687  /* Disabling heuristics so that the problem is not trivially solved */
4689 
4690  /* store parameters that are changed for the generation of the subproblem cuts */
4691  SCIP_CALL( SCIPsetBoolParam(subproblem, "conflict/enable", FALSE) );
4692 
4693  SCIP_CALL( SCIPsetIntParam(subproblem, "lp/disablecutoff", 1) );
4694  SCIP_CALL( SCIPsetIntParam(subproblem, "lp/scaling", 0) );
4695 
4696  SCIP_CALL( SCIPsetCharParam(subproblem, "lp/initalgorithm", 'd') );
4697  SCIP_CALL( SCIPsetCharParam(subproblem, "lp/resolvealgorithm", 'd') );
4698 
4699  SCIP_CALL( SCIPsetBoolParam(subproblem, "lp/alwaysgetduals", TRUE) );
4700  SCIP_CALL( SCIPsetBoolParam(subproblem, "misc/scaleobj", FALSE) );
4701 
4702  /* do not abort subproblem on CTRL-C */
4703  SCIP_CALL( SCIPsetBoolParam(subproblem, "misc/catchctrlc", FALSE) );
4704 
4705 #ifndef SCIP_MOREDEBUG
4706  SCIP_CALL( SCIPsetIntParam(subproblem, "display/verblevel", (int)SCIP_VERBLEVEL_NONE) );
4707 #endif
4708 
4709  SCIP_CALL( SCIPsetIntParam(subproblem, "propagating/maxrounds", 0) );
4710  SCIP_CALL( SCIPsetIntParam(subproblem, "propagating/maxroundsroot", 0) );
4711 
4712  SCIP_CALL( SCIPsetIntParam(subproblem, "constraints/linear/propfreq", -1) );
4713 
4714  SCIP_CALL( SCIPsetIntParam(subproblem, "heuristics/alns/freq", -1) );
4715 
4716  SCIP_CALL( SCIPsetIntParam(subproblem, "separating/aggregation/freq", -1) );
4717  SCIP_CALL( SCIPsetIntParam(subproblem, "separating/gomory/freq", -1) );
4718 
4719  return SCIP_OKAY;
4720 }
4721 
4722 /** resets the original parameters from the subproblem */
4723 static
4725  SCIP* subproblem, /**< the SCIP data structure */
4726  SCIP_SUBPROBPARAMS* origparams /**< the original subproblem parameters */
4727  )
4728 {
4729  assert(subproblem != NULL);
4730  assert(origparams != NULL);
4731 
4732  SCIP_CALL( SCIPsetRealParam(subproblem, "limits/memory", origparams->limits_memory) );
4733  SCIP_CALL( SCIPsetRealParam(subproblem, "limits/time", origparams->limits_time) );
4734  SCIP_CALL( SCIPsetBoolParam(subproblem, "conflict/enable", origparams->conflict_enable) );
4735  SCIP_CALL( SCIPsetIntParam(subproblem, "lp/disablecutoff", origparams->lp_disablecutoff) );
4736  SCIP_CALL( SCIPsetIntParam(subproblem, "lp/scaling", origparams->lp_scaling) );
4737  SCIP_CALL( SCIPsetCharParam(subproblem, "lp/initalgorithm", origparams->lp_initalg) );
4738  SCIP_CALL( SCIPsetCharParam(subproblem, "lp/resolvealgorithm", origparams->lp_resolvealg) );
4739  SCIP_CALL( SCIPsetBoolParam(subproblem, "lp/alwaysgetduals", origparams->lp_alwaysgetduals) );
4740  SCIP_CALL( SCIPsetBoolParam(subproblem, "misc/scaleobj", origparams->misc_scaleobj) );
4741  SCIP_CALL( SCIPsetBoolParam(subproblem, "misc/catchctrlc", origparams->misc_catchctrlc) );
4742  SCIP_CALL( SCIPsetIntParam(subproblem, "propagating/maxrounds", origparams->prop_maxrounds) );
4743  SCIP_CALL( SCIPsetIntParam(subproblem, "propagating/maxroundsroot", origparams->prop_maxroundsroot) );
4744  SCIP_CALL( SCIPsetIntParam(subproblem, "constraints/linear/propfreq", origparams->cons_linear_propfreq) );
4745 
4746  return SCIP_OKAY;
4747 }
4748 
4749 /** solves the LP of the Benders' decomposition subproblem
4750  *
4751  * This requires that the subproblem is in probing mode.
4752  */
4754  SCIP* scip, /**< the SCIP data structure */
4755  SCIP_BENDERS* benders, /**< the Benders' decomposition data structure */
4756  int probnumber, /**< the subproblem number */
4757  SCIP_STATUS* solvestatus, /**< status of subproblem solve */
4758  SCIP_Real* objective /**< optimal value of subproblem, if solved to optimality */
4759  )
4760 {
4761  SCIP* subproblem;
4762  SCIP_SUBPROBPARAMS* origparams;
4763  SCIP_Bool solvenlp;
4764 
4765  assert(benders != NULL);
4766  assert(solvestatus != NULL);
4767  assert(objective != NULL);
4768  assert(SCIPbendersSubproblemIsSetup(benders, probnumber));
4769 
4770  /* TODO: This should be solved just as an LP, so as a MIP. There is too much overhead with the MIP.
4771  * Need to change status check for checking the LP. */
4772  subproblem = SCIPbendersSubproblem(benders, probnumber);
4773  assert(subproblem != NULL);
4774 
4775  /* only solve the NLP relaxation if the NLP has been constructed and there exists an NLPI. If it is not possible to
4776  * solve the NLP relaxation, then the LP relaxation is used to generate Benders' cuts
4777  */
4778  solvenlp = FALSE;
4779  if( SCIPisNLPConstructed(subproblem) && SCIPgetNNlpis(subproblem) > 0
4780  && SCIPbendersGetSubproblemType(benders, probnumber) <= SCIP_BENDERSSUBTYPE_CONVEXDIS )
4781  solvenlp = TRUE;
4782 
4783  *objective = SCIPinfinity(subproblem);
4784 
4785  assert(SCIPisNLPConstructed(subproblem) || SCIPisLPConstructed(subproblem));
4786  assert(SCIPinProbing(subproblem));
4787 
4788  /* allocating memory for the parameter storage */
4789  SCIP_CALL( SCIPallocBlockMemory(subproblem, &origparams) );
4790 
4791  /* store the original parameters of the subproblem */
4792  SCIP_CALL( storeOrigSubproblemParams(subproblem, origparams) );
4793 
4794  /* setting the subproblem parameters */
4795  SCIP_CALL( setSubproblemParams(scip, subproblem) );
4796 
4797  if( solvenlp )
4798  {
4799  SCIP_NLPSOLSTAT nlpsolstat;
4800  SCIP_NLPTERMSTAT nlptermstat;
4801 #ifdef SCIP_MOREDEBUG
4802  SCIP_SOL* nlpsol;
4803 #endif
4804 
4805  SCIP_CALL( SCIPsolveNLP(subproblem) ); /*lint !e666*/
4806 
4807  nlpsolstat = SCIPgetNLPSolstat(subproblem);
4808  nlptermstat = SCIPgetNLPTermstat(subproblem);
4809  SCIPdebugMsg(scip, "NLP solstat %d termstat %d\n", nlpsolstat, nlptermstat);
4810 
4811  if( nlptermstat == SCIP_NLPTERMSTAT_OKAY && (nlpsolstat == SCIP_NLPSOLSTAT_LOCINFEASIBLE || nlpsolstat == SCIP_NLPSOLSTAT_GLOBINFEASIBLE) )
4812  {
4813  /* trust infeasible only if terminated "okay" */
4814  (*solvestatus) = SCIP_STATUS_INFEASIBLE;
4815  }
4816  else if( nlpsolstat == SCIP_NLPSOLSTAT_LOCOPT || nlpsolstat == SCIP_NLPSOLSTAT_GLOBOPT
4817  || nlpsolstat == SCIP_NLPSOLSTAT_FEASIBLE )
4818  {
4819 #ifdef SCIP_MOREDEBUG
4820  SCIP_CALL( SCIPcreateNLPSol(subproblem, &nlpsol, NULL) );
4821  SCIP_CALL( SCIPprintSol(subproblem, nlpsol, NULL, FALSE) );
4822  SCIP_CALL( SCIPfreeSol(subproblem, &nlpsol) );
4823 #endif
4824 
4825  (*solvestatus) = SCIP_STATUS_OPTIMAL;
4826  (*objective) = SCIPretransformObj(subproblem, SCIPgetNLPObjval(subproblem));
4827  }
4828  else if( nlpsolstat == SCIP_NLPSOLSTAT_UNBOUNDED )
4829  {
4830  (*solvestatus) = SCIP_STATUS_UNBOUNDED;
4831  SCIPerrorMessage("The NLP of Benders' decomposition subproblem %d is unbounded. This should not happen.\n",
4832  probnumber);
4833  SCIPABORT();
4834  }
4835  else if( nlptermstat == SCIP_NLPTERMSTAT_TIMELIMIT )
4836  {
4837  (*solvestatus) = SCIP_STATUS_TIMELIMIT;
4838  }
4839  else if( nlptermstat == SCIP_NLPTERMSTAT_INTERRUPT )
4840  {
4841  (*solvestatus) = SCIP_STATUS_USERINTERRUPT;
4842  }
4843  else
4844  {
4845  SCIPerrorMessage("Invalid solution status: %d. Termination status: %d. Solving the NLP relaxation of Benders' decomposition subproblem %d.\n",
4846  nlpsolstat, nlptermstat, probnumber);
4847  SCIPABORT();
4848  }
4849  }
4850  else
4851  {
4852  SCIP_Bool lperror;
4853  SCIP_Bool cutoff;
4854 
4855  SCIP_CALL( SCIPsolveProbingLP(subproblem, -1, &lperror, &cutoff) );
4856 
4857  switch( SCIPgetLPSolstat(subproblem) )
4858  {
4860  {
4861  (*solvestatus) = SCIP_STATUS_INFEASIBLE;
4862  break;
4863  }
4864 
4865  case SCIP_LPSOLSTAT_OPTIMAL :
4866  {
4867  (*solvestatus) = SCIP_STATUS_OPTIMAL;
4868  (*objective) = SCIPgetSolOrigObj(subproblem, NULL)*(int)SCIPgetObjsense(scip);
4869  break;
4870  }
4871 
4873  {
4874  (*solvestatus) = SCIP_STATUS_UNBOUNDED;
4875  SCIPerrorMessage("The LP of Benders' decomposition subproblem %d is unbounded. This should not happen.\n",
4876  probnumber);
4877  SCIPABORT();
4878  break;
4879  }
4880 
4881  case SCIP_LPSOLSTAT_ERROR :
4884  {
4885  if( SCIPgetLPSolstat(subproblem) == SCIP_LPSOLSTAT_TIMELIMIT )
4886  (*solvestatus) = SCIP_STATUS_TIMELIMIT;
4887  else
4888  (*solvestatus) = SCIP_STATUS_UNKNOWN;
4889 
4890  SCIPverbMessage(scip, SCIP_VERBLEVEL_FULL, NULL, " Benders' decomposition: Error solving LP "
4891  "relaxation of subproblem %d. No cut will be generated for this subproblem.\n", probnumber);
4892  break;
4893  }
4894 
4897  default:
4898  {
4899  SCIPerrorMessage("Invalid status: %d. Solving the LP relaxation of Benders' decomposition subproblem %d.\n",
4900  SCIPgetLPSolstat(subproblem), probnumber);
4901  SCIPABORT();
4902  break;
4903  }
4904  }
4905  }
4906 
4907  /* resetting the subproblem parameters */
4908  SCIP_CALL( resetOrigSubproblemParams(subproblem, origparams) );
4909 
4910  /* freeing the parameter storage */
4911  SCIPfreeBlockMemory(subproblem, &origparams);
4912 
4913  return SCIP_OKAY;
4914 }
4915 
4916 /** solves the Benders' decomposition subproblem */
4918  SCIP* scip, /**< the SCIP data structure */
4919  SCIP_BENDERS* benders, /**< the Benders' decomposition data structure */
4920  int probnumber, /**< the subproblem number */
4921  SCIP_STATUS* solvestatus, /**< status of subproblem solve */
4922  SCIP_Bool solvecip /**< directly solve the CIP subproblem */
4923  )
4924 {
4925  SCIP* subproblem;
4926  SCIP_SUBPROBPARAMS* origparams;
4927 
4928  assert(benders != NULL);
4929  assert(solvestatus != NULL);
4930 
4931  subproblem = SCIPbendersSubproblem(benders, probnumber);
4932  assert(subproblem != NULL);
4933 
4934  /* allocating memory for the parameter storage */
4935  SCIP_CALL( SCIPallocBlockMemory(subproblem, &origparams) );
4936 
4937  /* store the original parameters of the subproblem */
4938  SCIP_CALL( storeOrigSubproblemParams(subproblem, origparams) );
4939 
4940  /* If the solve has been stopped for the subproblem, then we need to restart it to complete the solve. The subproblem
4941  * is stopped when it is a MIP so that LP cuts and IP cuts can be generated. */
4942  if( SCIPgetStage(subproblem) == SCIP_STAGE_SOLVING )
4943  {
4944  /* the subproblem should be in probing mode. Otherwise, the event handler did not work correctly */
4945  assert( SCIPinProbing(subproblem) );
4946 
4947  /* the probing mode needs to be stopped so that the MIP can be solved */
4948  SCIP_CALL( SCIPendProbing(subproblem) );
4949 
4950  /* the problem was interrupted in the event handler, so SCIP needs to be informed that the problem is to be restarted */
4951  SCIP_CALL( SCIPrestartSolve(subproblem) );
4952  }
4953  else if( solvecip )
4954  {
4955  /* if the MIP will be solved directly, then the probing mode needs to be skipped.
4956  * This is achieved by setting the solvecip flag in the event handler data to TRUE
4957  */
4958  SCIP_EVENTHDLR* eventhdlr;
4959  SCIP_EVENTHDLRDATA* eventhdlrdata;
4960 
4961  eventhdlr = SCIPfindEventhdlr(subproblem, MIPNODEFOCUS_EVENTHDLR_NAME);
4962  eventhdlrdata = SCIPeventhdlrGetData(eventhdlr);
4963 
4964  eventhdlrdata->solvecip = TRUE;
4965  }
4966  else
4967  {
4968  /* if the problem is not in probing mode, then we need to solve the LP. That requires all methods that will
4969  * modify the structure of the problem need to be deactivated */
4970 
4971  /* setting the subproblem parameters */
4972  SCIP_CALL( setSubproblemParams(scip, subproblem) );
4973 
4974 #ifdef SCIP_EVENMOREDEBUG
4975  SCIP_CALL( SCIPsetBoolParam(subproblem, "display/lpinfo", TRUE) );
4976 #endif
4977  }
4978 
4979 #ifdef SCIP_MOREDEBUG
4980  SCIP_CALL( SCIPsetIntParam(subproblem, "display/verblevel", (int)SCIP_VERBLEVEL_FULL) );
4981  SCIP_CALL( SCIPsetIntParam(subproblem, "display/freq", 1) );
4982 #endif
4983 
4984  SCIP_CALL( SCIPsolve(subproblem) );
4985 
4986  *solvestatus = SCIPgetStatus(subproblem);
4987 
4988  if( *solvestatus != SCIP_STATUS_OPTIMAL && *solvestatus != SCIP_STATUS_UNBOUNDED
4989  && *solvestatus != SCIP_STATUS_INFEASIBLE && *solvestatus != SCIP_STATUS_USERINTERRUPT
4990  && *solvestatus != SCIP_STATUS_BESTSOLLIMIT && *solvestatus != SCIP_STATUS_TIMELIMIT
4991  && *solvestatus != SCIP_STATUS_MEMLIMIT )
4992  {
4993  SCIPerrorMessage("Invalid status: %d. Solving the CIP of Benders' decomposition subproblem %d.\n",
4994  *solvestatus, probnumber);
4995  SCIPABORT();
4996  }
4997 
4998  /* resetting the subproblem parameters */
4999  SCIP_CALL( resetOrigSubproblemParams(subproblem, origparams) );
5000 
5001  /* freeing the parameter storage */
5002  SCIPfreeBlockMemory(subproblem, &origparams);
5003 
5004  return SCIP_OKAY;
5005 }
5006 
5007 /** frees the subproblems */
5009  SCIP_BENDERS* benders, /**< Benders' decomposition */
5010  SCIP_SET* set, /**< global SCIP settings */
5011  int probnumber /**< the subproblem number */
5012  )
5013 {
5014  assert(benders != NULL);
5015  assert(benders->bendersfreesub != NULL
5016  || (benders->bendersfreesub == NULL && benders->benderssolvesubconvex == NULL && benders->benderssolvesub == NULL));
5017  assert(probnumber >= 0 && probnumber < benders->nsubproblems);
5018 
5019  if( benders->bendersfreesub != NULL )
5020  {
5021  SCIP_CALL( benders->bendersfreesub(set->scip, benders, probnumber) );
5022  }
5023  else
5024  {
5025  /* the subproblem is only freed if it is not independent */
5026  if( subproblemIsActive(benders, probnumber) )
5027  {
5028  SCIP* subproblem = SCIPbendersSubproblem(benders, probnumber);
5029 
5030  if( SCIPbendersGetSubproblemType(benders, probnumber) == SCIP_BENDERSSUBTYPE_CONVEXCONT )
5031  {
5032  /* ending probing mode to reset the current node. The probing mode will be restarted at the next solve */
5033  if( SCIPinProbing(subproblem) )
5034  {
5035  SCIP_CALL( SCIPendProbing(subproblem) );
5036  }
5037  }
5038  else
5039  {
5040  /* if the subproblems were solved as part of an enforcement stage, then they will still be in probing mode. The
5041  * probing mode must first be finished and then the problem can be freed */
5042  if( SCIPgetStage(subproblem) >= SCIP_STAGE_TRANSFORMED && SCIPinProbing(subproblem) )
5043  {
5044  SCIP_CALL( SCIPendProbing(subproblem) );
5045  }
5046 
5047  SCIP_CALL( SCIPfreeTransform(subproblem) );
5048  }
5049  }
5050  }
5051 
5052  /* setting the setup flag for the subproblem to FALSE */
5053  SCIPbendersSetSubproblemIsSetup(benders, probnumber, FALSE);
5054  return SCIP_OKAY;
5055 }
5056 
5057 /** compares the subproblem objective value with the auxiliary variable value for optimality */
5059  SCIP_BENDERS* benders, /**< the benders' decomposition structure */
5060  SCIP_SET* set, /**< global SCIP settings */
5061  SCIP_SOL* sol, /**< primal CIP solution */
5062  int probnumber /**< the subproblem number */
5063  )
5064 {
5065  SCIP_Real auxiliaryvarval;
5066  SCIP_Bool optimal;
5067 
5068  assert(benders != NULL);
5069  assert(set != NULL);
5070  assert(probnumber >= 0 && probnumber < benders->nsubproblems);
5071 
5072  optimal = FALSE;
5073 
5074  auxiliaryvarval = SCIPbendersGetAuxiliaryVarVal(benders, set, sol, probnumber);
5075 
5076  SCIPsetDebugMsg(set, "Subproblem %d - Auxiliary Variable: %g Subproblem Objective: %g Reldiff: %g Soltol: %g\n",
5077  probnumber, auxiliaryvarval, SCIPbendersGetSubproblemObjval(benders, probnumber),
5078  SCIPrelDiff(SCIPbendersGetSubproblemObjval(benders, probnumber), auxiliaryvarval), benders->solutiontol);
5079 
5080  if( SCIPrelDiff(SCIPbendersGetSubproblemObjval(benders, probnumber), auxiliaryvarval) < benders->solutiontol )
5081  optimal = TRUE;
5082 
5083  return optimal;
5084 }
5085 
5086 /** returns the value of the auxiliary variable value in a master problem solution */
5088  SCIP_BENDERS* benders, /**< the benders' decomposition structure */
5089  SCIP_SET* set, /**< global SCIP settings */
5090  SCIP_SOL* sol, /**< primal CIP solution */
5091  int probnumber /**< the subproblem number */
5092  )
5093 {
5094  SCIP_VAR* auxiliaryvar;
5095 
5096  assert(benders != NULL);
5097  assert(set != NULL);
5098 
5099  auxiliaryvar = SCIPbendersGetAuxiliaryVar(benders, probnumber);
5100  assert(auxiliaryvar != NULL);
5101 
5102  return SCIPgetSolVal(set->scip, sol, auxiliaryvar);
5103 }
5104 
5105 /** Solves an independent subproblem to identify its lower bound. The lower bound is then used to update the bound on
5106  * the auxiliary variable.
5107  */
5109  SCIP_BENDERS* benders, /**< Benders' decomposition */
5110  SCIP_SET* set, /**< global SCIP settings */
5111  int probnumber, /**< the subproblem to be evaluated */
5112  SCIP_Real* lowerbound, /**< the lowerbound for the subproblem */
5113  SCIP_Bool* infeasible /**< was the subproblem found to be infeasible? */
5114  )
5115 {
5116  SCIP* subproblem;
5117  SCIP_Real dualbound;
5118  SCIP_Real memorylimit;
5119  SCIP_Real timelimit;
5120  SCIP_Longint totalnodes;
5121  int disablecutoff;
5122  int verblevel;
5123  SCIP_Bool lperror;
5124  SCIP_Bool cutoff;
5125 
5126  assert(benders != NULL);
5127  assert(set != NULL);
5128 
5129  if( benders->benderssolvesub != NULL || benders->benderssolvesubconvex != NULL )
5130  {
5131  (*lowerbound) = SCIPvarGetLbGlobal(SCIPbendersGetAuxiliaryVar(benders, probnumber));
5132  (*infeasible) = FALSE;
5133 
5134  SCIPinfoMessage(set->scip, NULL, "Benders' decomposition: a bendersSolvesub or bendersSolvesubconvex has been "
5135  "implemented. SCIPbendersComputeSubproblemLowerbound can not be executed.\n");
5136  SCIPinfoMessage(set->scip, NULL, "Set the auxiliary variable lower bound by calling "
5137  "SCIPbendersUpdateSubproblemLowerbound in bendersCreatesub. The auxiliary variable %d will remain as %g\n",
5138  probnumber, (*lowerbound));
5139 
5140  return SCIP_OKAY;
5141  }
5142  else
5143  {
5144  SCIPverbMessage(set->scip, SCIP_VERBLEVEL_FULL, NULL, "Benders' decomposition: Computing a lower bound for"
5145  " subproblem %d\n", probnumber);
5146  }
5147 
5148  /* getting the subproblem to evaluate */
5149  subproblem = SCIPbendersSubproblem(benders, probnumber);
5150 
5151  (*lowerbound) = -SCIPinfinity(subproblem);
5152  (*infeasible) = FALSE;
5153 
5154  SCIP_CALL( SCIPgetIntParam(subproblem, "display/verblevel", &verblevel) );
5155  SCIP_CALL( SCIPsetIntParam(subproblem, "display/verblevel", (int)SCIP_VERBLEVEL_NONE) );
5156 #ifdef SCIP_MOREDEBUG
5157  SCIP_CALL( SCIPsetIntParam(subproblem, "display/verblevel", (int)SCIP_VERBLEVEL_HIGH) );
5158 #endif
5159 
5160  /* copying memory and time limits */
5161  SCIP_CALL( SCIPgetRealParam(subproblem, "limits/time", &timelimit) );
5162  SCIP_CALL( SCIPgetRealParam(subproblem, "limits/memory", &memorylimit) );
5163  SCIP_CALL( copyMemoryAndTimeLimits(set->scip, subproblem) );
5164 
5165  /* if the subproblem is independent, then the default SCIP settings are used. Otherwise, only the root node is solved
5166  * to compute a lower bound on the subproblem
5167  */
5168  SCIP_CALL( SCIPgetLongintParam(subproblem, "limits/totalnodes", &totalnodes) );
5169  SCIP_CALL( SCIPgetIntParam(subproblem, "lp/disablecutoff", &disablecutoff) );
5170  if( !SCIPbendersSubproblemIsIndependent(benders, probnumber) )
5171  {
5172  SCIP_CALL( SCIPsetLongintParam(subproblem, "limits/totalnodes", 1LL) );
5173  SCIP_CALL( SCIPsetIntParam(subproblem, "lp/disablecutoff", 1) );
5174  }
5175 
5176  /* if the subproblem not independent and is convex, then the probing LP is solved. Otherwise, the MIP is solved */
5177  dualbound = -SCIPinfinity(subproblem);
5178  if( SCIPbendersGetSubproblemType(benders, probnumber) == SCIP_BENDERSSUBTYPE_CONVEXCONT )
5179  {
5180  SCIP_Bool solvenlp = FALSE;
5181 
5182  assert(SCIPisLPConstructed(subproblem) || SCIPisNLPConstructed(subproblem));
5183 
5184  if( SCIPisNLPConstructed(subproblem) && SCIPgetNNlpis(subproblem) > 0
5185  && SCIPbendersGetSubproblemType(benders, probnumber) <= SCIP_BENDERSSUBTYPE_CONVEXDIS )
5186  solvenlp = TRUE;
5187 
5188  SCIP_CALL( SCIPstartProbing(subproblem) );
5189  if( solvenlp )
5190  {
5191  SCIP_NLPSOLSTAT nlpsolstat;
5192  SCIP_NLPTERMSTAT nlptermstat;
5193 
5194  SCIP_CALL( SCIPsolveNLP(subproblem) ); /*lint !e666*/
5195 
5196  nlpsolstat = SCIPgetNLPSolstat(subproblem);
5197  nlptermstat = SCIPgetNLPTermstat(subproblem);
5198  SCIPdebugMsg(set->scip, "NLP solstat %d termstat %d\n", nlpsolstat, nlptermstat);
5199 
5200  if( nlptermstat == SCIP_NLPTERMSTAT_OKAY && (nlpsolstat == SCIP_NLPSOLSTAT_LOCINFEASIBLE || nlpsolstat == SCIP_NLPSOLSTAT_GLOBINFEASIBLE) )
5201  {
5202  /* trust infeasible only if terminated "okay" */
5203  (*infeasible) = TRUE;
5204  }
5205  else if( nlpsolstat == SCIP_NLPSOLSTAT_LOCOPT || nlpsolstat == SCIP_NLPSOLSTAT_GLOBOPT
5206  || nlpsolstat == SCIP_NLPSOLSTAT_FEASIBLE )
5207  {
5208  dualbound = SCIPretransformObj(subproblem, SCIPgetNLPObjval(subproblem));
5209  }
5210  }
5211  else
5212  {
5213  SCIP_CALL( SCIPsolveProbingLP(subproblem, -1, &lperror, &cutoff) );
5214 
5215  if( SCIPgetLPSolstat(subproblem) == SCIP_LPSOLSTAT_INFEASIBLE )
5216  (*infeasible) = TRUE;
5217  else if( SCIPgetLPSolstat(subproblem) == SCIP_LPSOLSTAT_OPTIMAL )
5218  dualbound = SCIPgetSolOrigObj(subproblem, NULL)*(int)SCIPgetObjsense(set->scip);
5219  }
5220  }
5221  else
5222  {
5223  SCIP_EVENTHDLRDATA* eventhdlrdata;
5224 
5225  /* if the subproblem is not convex, then event handlers have been added to interrupt the solve. These must be
5226  * disabled
5227  */
5229  eventhdlrdata->solvecip = TRUE;
5230 
5231  SCIP_CALL( SCIPsolve(subproblem) );
5232 
5233  if( SCIPgetStatus(subproblem) == SCIP_STATUS_INFEASIBLE )
5234  (*infeasible) = TRUE;
5235  else
5236  dualbound = SCIPgetDualbound(subproblem);
5237  }
5238 
5239  /* getting the lower bound value */
5240  (*lowerbound) = dualbound;
5241 
5242  if( !SCIPbendersSubproblemIsIndependent(benders, probnumber) )
5243  {
5244  SCIP_CALL( SCIPsetLongintParam(subproblem, "limits/totalnodes", totalnodes) );
5245  SCIP_CALL( SCIPsetIntParam(subproblem, "lp/disablecutoff", disablecutoff) );
5246  }
5247  SCIP_CALL( SCIPsetIntParam(subproblem, "display/verblevel", verblevel) );
5248  SCIP_CALL( SCIPsetRealParam(subproblem, "limits/memory", memorylimit) );
5249  SCIP_CALL( SCIPsetRealParam(subproblem, "limits/time", timelimit) );
5250 
5251  /* the subproblem must be freed so that it is reset for the subsequent Benders' decomposition solves. If the
5252  * subproblems are independent, they are not freed. SCIPfreeBendersSubproblem must still be called, but in this
5253  * function the independent subproblems are not freed. However, they will still be freed at the end of the
5254  * solving process for the master problem.
5255  */
5256  SCIP_CALL( SCIPbendersFreeSubproblem(benders, set, probnumber) );
5257 
5258  return SCIP_OKAY;
5259 }
5260 
5261 /** Merges a subproblem into the master problem. This process just adds a copy of the subproblem variables and
5262  * constraints to the master problem, but keeps the subproblem stored in the Benders' decomposition data structure. The reason for
5263  * keeping the subproblem available is for when it is queried for solutions after the problem is solved.
5264  *
5265  * Once the subproblem is merged into the master problem, then the subproblem is flagged as disabled. This means that
5266  * it will not be solved in the subsequent subproblem solving loops.
5267  *
5268  * The associated auxiliary variables are kept in the master problem. The objective function of the merged subproblem
5269  * is added as an underestimator constraint.
5270  */
5272  SCIP_BENDERS* benders, /**< Benders' decomposition */
5273  SCIP_SET* set, /**< global SCIP settings */
5274  SCIP_HASHMAP* varmap, /**< a hashmap to store the mapping of subproblem variables corresponding
5275  * to the newly created master variables, or NULL */
5276  SCIP_HASHMAP* consmap, /**< a hashmap to store the mapping of subproblem constraints to the
5277  * corresponding newly created constraints, or NULL */
5278  int probnumber /**< the number of the subproblem that will be merged into the master problem*/
5279  )
5280 {
5281  SCIP* subproblem;
5282  SCIP_HASHMAP* localvarmap;
5283  SCIP_HASHMAP* localconsmap;
5284  SCIP_VAR** vars;
5285  SCIP_VAR* auxiliaryvar;
5286  SCIP_CONS** conss;
5287  SCIP_CONS* objcons;
5288  int nvars;
5289  int nconss;
5290  int i;
5291  SCIP_Bool uselocalvarmap;
5292  SCIP_Bool uselocalconsmap;
5293  char varname[SCIP_MAXSTRLEN];
5294  char consname[SCIP_MAXSTRLEN];
5295  const char* origvarname;
5296 
5297  assert(benders != NULL);
5298  assert(set != NULL);
5299  assert(probnumber >= 0 && probnumber < benders->nsubproblems);
5300 
5301  SCIPverbMessage(set->scip, SCIP_VERBLEVEL_HIGH, NULL, " Benders' decomposition: Infeasibility of subproblem %d can't "
5302  "be resolved. Subproblem %d is being merged into the master problem.\n", probnumber, probnumber);
5303 
5304  /* freeing the subproblem because it will be flagged as independent. Since the subproblem is flagged as independent,
5305  * it will no longer be solved or freed within the solving loop.
5306