1 /* * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * */
2 /* */
3 /* This file is part of the program and library */
4 /* SCIP --- Solving Constraint Integer Programs */
5 /* */
6 /* Copyright (c) 2002-2023 Zuse Institute Berlin (ZIB) */
7 /* */
8 /* Licensed under the Apache License, Version 2.0 (the "License"); */
9 /* you may not use this file except in compliance with the License. */
10 /* You may obtain a copy of the License at */
11 /* */
12 /* http://www.apache.org/licenses/LICENSE-2.0 */
13 /* */
14 /* Unless required by applicable law or agreed to in writing, software */
15 /* distributed under the License is distributed on an "AS IS" BASIS, */
16 /* WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. */
17 /* See the License for the specific language governing permissions and */
18 /* limitations under the License. */
19 /* */
20 /* You should have received a copy of the Apache-2.0 license */
21 /* along with SCIP; see the file LICENSE. If not visit scipopt.org. */
22 /* */
23 /* * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * */
24
25 /**@file heur_dps.c
26 * @ingroup DEFPLUGINS_HEUR
27 * @brief dynamic partition search
28 * @author Katrin Halbig
29 *
30 * The dynamic partition search (DPS) is a construction heuristic which additionally needs a
31 * user decomposition with linking constraints only.
32 *
33 * This heuristic splits the problem into several sub-SCIPs according to the given decomposition. Thereby the linking constraints
34 * with their right-hand and left-hand sides are also split. DPS searches for a partition of the sides on the blocks
35 * so that a feasible solution is obtained.
36 * For each block the parts of the original linking constraints are extended by slack variables. Moreover, the objective function
37 * is replaced by the sum of these additional variables weighted by penalty parameters lambda. If all blocks have an optimal solution
38 * of zero, the algorithm terminates with a feasible solution for the main problem. Otherwise, the partition and the penalty parameters
39 * are updated, and the sub-SCIPs are solved again.
40 *
41 * A detailed description can be found in
42 * K. Halbig, A. Göß and D. Weninger (2023). Exploiting user-supplied Decompositions inside Heuristics. https://optimization-online.org/?p=23386
43 */
44
45 /*---+----1----+----2----+----3----+----4----+----5----+----6----+----7----+----8----+----9----+----0----+----1----+----2*/
46
47 #include "scip/heur_dps.h"
48 #include "scip/pub_dcmp.h"
49 #include "scip/pub_heur.h"
50 #include "scip/pub_misc.h"
51 #include "scip/scipdefplugins.h"
52 #include "scip/scip_cons.h"
53 #include "scip/scip_dcmp.h"
54 #include "scip/scip_general.h"
55 #include "scip/scip_heur.h"
56 #include "scip/scip_mem.h"
57 #include "scip/scip_message.h"
58 #include "scip/scip_param.h"
59 #include "scip/scip_prob.h"
60
61
62 #define HEUR_NAME "dps"
63 #define HEUR_DESC "primal heuristic for decomposable MIPs"
64 #define HEUR_DISPCHAR SCIP_HEURDISPCHAR_LNS
65 #define HEUR_PRIORITY 75000
66 #define HEUR_FREQ -1
67 #define HEUR_FREQOFS 0
68 #define HEUR_MAXDEPTH -1
69 #define HEUR_TIMING SCIP_HEURTIMING_BEFORENODE | SCIP_HEURTIMING_AFTERNODE
70 #define HEUR_USESSUBSCIP TRUE /**< does the heuristic use a secondary SCIP instance? */
71
72 #define DEFAULT_MAXIT 50 /**< maximum number of iterations */
73 #define DEFAULT_PENALTY 100.0 /**< multiplier for absolute increase of penalty parameters */
74
75 /* event handler properties */
76 #define EVENTHDLR_NAME "Dps"
77 #define EVENTHDLR_DESC "event handler for " HEUR_NAME " heuristic"
78
79 /*
80 * Data structures
81 */
82
83 /** primal heuristic data */
84 struct SCIP_HeurData
85 {
86 SCIP_CONS** linkingconss; /**< linking constraints */
87 int nlinking; /**< number of linking constraints */
88 int nblocks; /**< number of blocks */
89 int maxit; /**< maximal number of iterations */
90 int timing; /**< should the heuristic run before or after the processing of the node?
91 (0: before, 1: after, 2: both)*/
92 SCIP_Real maxlinkscore; /**< maximal linking score of used decomposition (equivalent to percentage of linking constraints) */
93 SCIP_Real penalty; /**< multiplier for absolute increase of penalty parameters */
94 SCIP_Bool reoptimize; /**< should the problem get reoptimized with the original objective function? */
95 SCIP_Bool reuse; /**< should solutions get reused in subproblems? */
96 SCIP_Bool reoptlimits; /**< should strict limits for reoptimization be set? */
97 };
98
99 /** data related to one block */
100 struct Blockproblem
101 {
102 SCIP* blockscip; /**< SCIP data structure */
103 SCIP_VAR** slackvars; /**< slack variables */
104 SCIP_CONS** linkingconss; /**< linking constraints */
105 int* linkingindices; /**< indices of linking constraints in original problem */
106 int nlinking; /**< number of linking constraints */
107 int nblockvars; /**< number of block variables */
108 int nslackvars; /**< number of slack variables */
109 SCIP_Real* origobj; /**< original objective coefficients */
110 };
111 typedef struct Blockproblem BLOCKPROBLEM;
112
113 /** data related to one linking constraint */
114 struct Linking
115 {
116 SCIP_CONS* linkingcons; /**< corresponding linking constraint of original problem */
117 SCIP_CONS** blockconss; /**< linking constraints of the blocks */
118 SCIP_VAR** slacks; /**< slackvars of the blocks */
119 SCIP_Real* minactivity; /**< minimal activity of constraint for each block */
120 SCIP_Real* maxactivity; /**< maximal activity of constraint for each block */
121 SCIP_Real* currentrhs; /**< current partition of rhs */
122 SCIP_Real* currentlhs; /**< current partition of lhs */
123 int* blocknumbers; /**< number of the blocks */
124 int nblocks; /**< number of blocks in which this linking constraint participates; dimension of arrays */
125 int nslacks; /**< number of slack variables */
126 int nslacksperblock; /**< 2, if ranged constraint; 1, if only rhs or lhs */
127 int lastviolations; /**< number of iterations in which the constraint was violated in succession */
128 SCIP_Bool hasrhs; /**< has linking constraint finite right-hand side? */
129 SCIP_Bool haslhs; /**< has linking constraint finite left-hand side? */
130 };
131 typedef struct Linking LINKING;
132
133 /*
134 * Local methods
135 */
136
137 /** assigns linking variables to last block
138 *
139 * The labels are copied to newdecomp and the linking variables are assigned to the last block (i.e., highest block label).
140 * Constraint labels and statistics are recomputed.
141 */
142 static
143 SCIP_RETCODE assignLinking(
144 SCIP* scip, /**< SCIP data structure */
145 SCIP_DECOMP* newdecomp, /**< decomposition with assigned linking variables */
146 SCIP_VAR** vars, /**< sorted array of variables */
147 SCIP_CONS** conss, /**< sorted array of constraints */
148 int* varlabels, /**< sorted array of variable labels */
149 int* conslabels, /**< sorted array of constraint labels */
150 int nvars, /**< number of variables */
151 int nconss, /**< number of constraints */
152 int nlinkvars /**< number of linking variables */
153 )
154 {
155 int newlabel;
156 int v;
157
158 assert(scip != NULL);
159 assert(newdecomp != NULL);
160 assert(vars != NULL);
161 assert(conss != NULL);
162 assert(varlabels != NULL);
163 assert(conslabels != NULL);
164
165 /* copy the labels */
166 SCIP_CALL( SCIPdecompSetVarsLabels(newdecomp, vars, varlabels, nvars) );
167 SCIP_CALL( SCIPdecompSetConsLabels(newdecomp, conss, conslabels, nconss) );
168
169 /* assign linking variables */
170 newlabel = varlabels[nvars - 1]; /* take always label of last block */
171 assert(newlabel >= 0);
172 for( v = 0; v < nlinkvars; v++ )
173 {
174 SCIP_CALL( SCIPdecompSetVarsLabels(newdecomp, &vars[v], &newlabel, 1) );
175 }
176 SCIPdebugMsg(scip, "assigned %d linking variables\n", nlinkvars);
177
178 /* recompute constraint labels and statistics */
179 SCIP_CALL( SCIPcomputeDecompConsLabels(scip, newdecomp, conss, nconss) );
180 SCIP_CALL( SCIPcomputeDecompStats(scip, newdecomp, TRUE) );
181 nlinkvars = SCIPdecompGetNBorderVars(newdecomp);
182
183 /* get new labels and sort */
184 SCIPdecompGetConsLabels(newdecomp, conss, conslabels, nconss);
185 SCIPdecompGetVarsLabels(newdecomp, vars, varlabels, nvars);
186 SCIPsortIntPtr(conslabels, (void**)conss, nconss);
187 SCIPsortIntPtr(varlabels, (void**)vars, nvars);
188
189 /* After assigning the linking variables, blocks can have zero constraints.
190 * So the remaining variables are labeled as linking in SCIPcomputeDecompStats().
191 * We assign this variables to the same label as above.
192 */
193 if( nlinkvars >= 1 )
194 {
195 assert(varlabels[0] == SCIP_DECOMP_LINKVAR);
196 SCIPdebugMsg(scip, "assign again %d linking variables\n", nlinkvars);
197
198 for( v = 0; v < nlinkvars; v++ )
199 {
200 SCIP_CALL( SCIPdecompSetVarsLabels(newdecomp, &vars[v], &newlabel, 1) );
201 }
202 SCIP_CALL( SCIPcomputeDecompConsLabels(scip, newdecomp, conss, nconss) );
203 SCIP_CALL( SCIPcomputeDecompStats(scip, newdecomp, TRUE) );
204
205 SCIPdecompGetConsLabels(newdecomp, conss, conslabels, nconss);
206 SCIPdecompGetVarsLabels(newdecomp, vars, varlabels, nvars);
207 SCIPsortIntPtr(conslabels, (void**)conss, nconss);
208 SCIPsortIntPtr(varlabels, (void**)vars, nvars);
209 }
210 assert(varlabels[0] != SCIP_DECOMP_LINKVAR);
211
212 return SCIP_OKAY;
213 }
214
215 /** creates a sub-SCIP and sets parameters */
216 static
217 SCIP_RETCODE createSubscip(
218 SCIP* scip, /**< main SCIP data structure */
219 SCIP** subscip /**< pointer to store created sub-SCIP */
220 )
221 {
222 SCIP_Real infvalue;
223
224 assert(scip != NULL);
225 assert(subscip != NULL);
226
227 /* create a new SCIP instance */
228 SCIP_CALL( SCIPcreate(subscip) );
229 SCIP_CALL( SCIPincludeDefaultPlugins(*subscip) );
230
231 /* copy value for infinity */
232 SCIP_CALL( SCIPgetRealParam(scip, "numerics/infinity", &infvalue) );
233 SCIP_CALL( SCIPsetRealParam(*subscip, "numerics/infinity", infvalue) );
234
235 SCIP_CALL( SCIPcopyLimits(scip, *subscip) );
236
237 /* avoid recursive calls */
238 SCIP_CALL( SCIPsetSubscipsOff(*subscip, TRUE) );
239
240 /* disable cutting plane separation */
241 SCIP_CALL( SCIPsetSeparating(*subscip, SCIP_PARAMSETTING_OFF, TRUE) );
242
243 /* disable expensive presolving */
244 SCIP_CALL( SCIPsetPresolving(*subscip, SCIP_PARAMSETTING_FAST, TRUE) );
245
246 /* disable expensive techniques */
247 SCIP_CALL( SCIPsetIntParam(*subscip, "misc/usesymmetry", 0) );
248
249 /* do not abort subproblem on CTRL-C */
250 SCIP_CALL( SCIPsetBoolParam(*subscip, "misc/catchctrlc", FALSE) );
251
252 #ifdef SCIP_DEBUG
253 /* for debugging, enable full output */
254 SCIP_CALL( SCIPsetIntParam(*subscip, "display/verblevel", 5) );
255 SCIP_CALL( SCIPsetIntParam(*subscip, "display/freq", 100000000) );
256 #else
257 /* disable statistic timing inside sub SCIP and output to console */
258 SCIP_CALL( SCIPsetIntParam(*subscip, "display/verblevel", 0) );
259 SCIP_CALL( SCIPsetBoolParam(*subscip, "timing/statistictiming", FALSE) );
260 #endif
261
262 /* speed up sub-SCIP by not checking dual LP feasibility */
263 SCIP_CALL( SCIPsetBoolParam(*subscip, "lp/checkdualfeas", FALSE) );
264
265 return SCIP_OKAY;
266 }
267
268 /** copies the given variables and constraints to the given sub-SCIP */
269 static
270 SCIP_RETCODE copyToSubscip(
271 SCIP* scip, /**< source SCIP */
272 SCIP* subscip, /**< target SCIP */
273 const char* name, /**< name for copied problem */
274 SCIP_VAR** vars, /**< array of variables to copy */
275 SCIP_CONS** conss, /**< array of constraints to copy */
276 SCIP_HASHMAP* varsmap, /**< hashmap for copied variables */
277 SCIP_HASHMAP* conssmap, /**< hashmap for copied constraints */
278 int nvars, /**< number of variables to copy */
279 int nconss, /**< number of constraints to copy */
280 SCIP_Bool* success /**< was copying successful? */
281 )
282 {
283 SCIP_CONS* newcons;
284 SCIP_VAR* newvar;
285 int i;
286
287 assert(scip != NULL);
288 assert(subscip != NULL);
289 assert(vars != NULL);
290 assert(conss != NULL);
291 assert(varsmap != NULL);
292 assert(conssmap != NULL);
293 assert(success != NULL);
294
295 SCIPdebugMsg(scip, "copyToSubscip\n");
296
297 /* create problem in sub-SCIP */
298 SCIP_CALL( SCIPcreateProb(subscip, name, NULL, NULL, NULL, NULL, NULL, NULL, NULL) );
299
300 /* copy variables */
301 for( i = 0; i < nvars; ++i )
302 {
303 SCIP_CALL( SCIPgetVarCopy(scip, subscip, vars[i], &newvar, varsmap, conssmap, FALSE, success) );
304
305 /* abort if variable was not successfully copied */
306 if( !(*success) )
307 {
308 SCIPwarningMessage(scip, "Abort heuristic dps since not all variables were successfully copied.\n");
309 SCIPABORT();
310 return SCIP_OKAY;
311 }
312 }
313 assert(nvars == SCIPgetNOrigVars(subscip));
314
315 /* copy constraints */
316 for( i = 0; i < nconss; ++i )
317 {
318 assert(conss[i] != NULL);
319 assert(!SCIPconsIsModifiable(conss[i]));
320 assert(SCIPconsIsActive(conss[i]));
321 assert(!SCIPconsIsDeleted(conss[i]));
322
323 /* copy the constraint */
324 SCIP_CALL( SCIPgetConsCopy(scip, subscip, conss[i], &newcons, SCIPconsGetHdlr(conss[i]), varsmap, conssmap, NULL,
325 SCIPconsIsInitial(conss[i]), SCIPconsIsSeparated(conss[i]), SCIPconsIsEnforced(conss[i]),
326 SCIPconsIsChecked(conss[i]), SCIPconsIsPropagated(conss[i]), FALSE, FALSE,
327 SCIPconsIsDynamic(conss[i]), SCIPconsIsRemovable(conss[i]), FALSE, FALSE, success) );
328
329 /* abort if constraint was not successfully copied */
330 if( !(*success) )
331 return SCIP_OKAY;
332
333 SCIP_CALL( SCIPaddCons(subscip, newcons) );
334 SCIP_CALL( SCIPreleaseCons(subscip, &newcons) );
335 }
336
337 /* block constraint contains variables which are not part of this block
338 *
339 * todo: maybe they are part of the block, but it is not recognized, because they are, for example, negated or aggregated.
340 */
341 if( nvars != SCIPgetNOrigVars(subscip) )
342 *success = FALSE;
343
344 return SCIP_OKAY;
345 }
346
347 /** creates the subscip for a given block */
348 static
349 SCIP_RETCODE createBlockproblem(
350 SCIP* scip, /**< SCIP data structure */
351 BLOCKPROBLEM* blockproblem, /**< blockproblem that should be created */
352 LINKING** linkings, /**< linkings that will be (partially) initialized */
353 SCIP_CONS** conss, /**< sorted array of constraints of this block */
354 SCIP_VAR** vars, /**< sorted array of variables of this block */
355 int nconss, /**< number of constraints of this block */
356 int nvars, /**< number of variables of this block */
357 SCIP_CONS** linkingconss, /**< linking constraints in the original problem */
358 int nlinking, /**< number of linking constraints in the original problem */
359 int blocknumber, /**< number of block that should be created */
360 SCIP_Bool* success /**< pointer to store whether creation was successful */
361 )
362 {
363 char name[SCIP_MAXSTRLEN];
364 SCIP_HASHMAP* varsmap;
365 SCIP_HASHMAP* conssmap;
366 SCIP_VAR** consvars; /* all vars in original linking cons */
367 SCIP_Real* consvals;
368 int nconsvars;
369 SCIP_VAR** blockvars; /* vars of current linking cons of current block */
370 SCIP_Real* blockvals;
371 int nblockvars;
372 SCIP_VAR** subvars; /* all vars of subscip */
373 int maxnconsvars; /* current size of arrays */
374 int c;
375 int v;
376
377 assert(scip != NULL);
378 assert(blockproblem != NULL);
379 assert(conss != NULL);
380 assert(vars != NULL);
381 assert(blockproblem->blockscip != NULL);
382
383 maxnconsvars = 20; /* start size; increase size if necessary */
384
385 SCIPdebugMsg(scip, "Create blockproblem %d\n", blocknumber);
386
387 /* create the variable/constraint mapping hash map */
388 SCIP_CALL( SCIPhashmapCreate(&varsmap, SCIPblkmem(scip), nvars) );
389 SCIP_CALL( SCIPhashmapCreate(&conssmap, SCIPblkmem(scip), nconss) );
390
391 /* get name of the original problem and add "comp_nr" */
392 (void)SCIPsnprintf(name, SCIP_MAXSTRLEN, "%s_comp_%d", SCIPgetProbName(scip), blocknumber);
393
394 SCIP_CALL( copyToSubscip(scip, blockproblem->blockscip, name, vars, conss, varsmap, conssmap,
395 nvars, nconss, success) );
396 if( !(*success) )
397 {
398 SCIPdebugMsg(scip, "Copy to subscip failed\n");
399 SCIPhashmapFree(&conssmap);
400 SCIPhashmapFree(&varsmap);
401
402 return SCIP_OKAY;
403 }
404
405 /* save number of variables that have a corresponding variable in original problem*/
406 blockproblem->nblockvars = SCIPgetNVars(blockproblem->blockscip);
407
408 /* save original objective and set objective to zero */
409 subvars = SCIPgetVars(blockproblem->blockscip);
410 for( v = 0; v < nvars; v++ )
411 {
412 blockproblem->origobj[v] = SCIPvarGetObj(subvars[v]);
413 SCIP_CALL( SCIPchgVarObj(blockproblem->blockscip, subvars[v], 0.0) );
414 }
415
416 /* allocate memory */
417 SCIP_CALL( SCIPallocBufferArray(blockproblem->blockscip, &blockvars, nvars + 2) ); /* two entries for the slack variables */
418 SCIP_CALL( SCIPallocBufferArray(blockproblem->blockscip, &blockvals, nvars + 2) );
419 SCIP_CALL( SCIPallocBufferArray(blockproblem->blockscip, &consvars, maxnconsvars) );
420 SCIP_CALL( SCIPallocBufferArray(blockproblem->blockscip, &consvals, maxnconsvars) );
421
422 /* find and add parts of linking constraints */
423 SCIPdebugMsg(scip, "add parts of linking constraints\n");
424 for( c = 0; c < nlinking; c++ )
425 {
426 const char* conshdlrname;
427 char consname[SCIP_MAXSTRLEN];
428 SCIP_CONS* newcons;
429 SCIP_Real rhs;
430 SCIP_Real lhs;
431 SCIP_Real minact;
432 SCIP_Real maxact;
433 SCIP_Bool mininfinite;
434 SCIP_Bool maxinfinite;
435
436 assert(linkingconss[c] != NULL);
437
438 newcons = NULL;
439
440 #ifdef SCIP_MORE_DEBUG
441 SCIPdebugMsg(scip, "consider constraint %s\n", SCIPconsGetName(linkingconss[c]));
442 SCIPdebugPrintCons(scip, linkingconss[c], NULL);
443 #endif
444
445 nblockvars = 0;
446
447 /* every constraint with linear representation is allowed */
448 conshdlrname = SCIPconshdlrGetName(SCIPconsGetHdlr(linkingconss[c]));
449 if( !( (strcmp(conshdlrname, "linear") == 0) || (strcmp(conshdlrname, "setppc") == 0)
450 || (strcmp(conshdlrname, "logicor") == 0) || (strcmp(conshdlrname, "knapsack") == 0)
451 || (strcmp(conshdlrname, "varbound") == 0) ) )
452 {
453 SCIPverbMessage(scip, SCIP_VERBLEVEL_FULL, NULL, "Heuristic %s cannot handle linking constraints of type %s\n", HEUR_NAME, conshdlrname);
454 /* TODO which other types can we handle/transform in a linear constraint? */
455
456 *success = FALSE;
457 break; /* releases memory and breaks heuristic */
458 }
459
460 SCIP_CALL( SCIPgetConsNVars(scip, linkingconss[c], &nconsvars, success) );
461
462 /* reallocate memory if we have more variables than maxnconsvars */
463 if( nconsvars > maxnconsvars )
464 {
465 int newsize;
466
467 /* calculate new size */
468 newsize = SCIPcalcMemGrowSize(scip, MAX(2 * maxnconsvars, nconsvars)); /* at least double size */
469
470 SCIP_CALL( SCIPreallocBufferArray(blockproblem->blockscip, &consvars, newsize) );
471 SCIP_CALL( SCIPreallocBufferArray(blockproblem->blockscip, &consvals, newsize) );
472 maxnconsvars = newsize;
473 }
474
475 SCIP_CALL( SCIPgetConsVars(scip, linkingconss[c], consvars, nconsvars, success) );
476 SCIP_CALL( SCIPgetConsVals(scip, linkingconss[c], consvals, nconsvars, success) );
477
478 if( !(*success) )
479 {
480 SCIPdebugMsg(scip, "Create blockproblem failed\n");
481 break; /* releases memory and breaks heuristic */
482 }
483
484 /* check if constraint contains variables of this block */
485 for( v = 0; v < nconsvars; v++ )
486 {
487 if( SCIPhashmapExists(varsmap, (void*)consvars[v]) )
488 {
489 blockvars[nblockvars] = (SCIP_VAR*) SCIPhashmapGetImage(varsmap, (void*)consvars[v]);
490 blockvals[nblockvars] = consvals[v];
491 ++nblockvars;
492 }
493 /* handle negated variables*/
494 else if( SCIPvarGetStatus(consvars[v]) == SCIP_VARSTATUS_NEGATED)
495 {
496 if( SCIPhashmapExists(varsmap, (void*)SCIPvarGetNegationVar(consvars[v])) ) /* negation exists in this block */
497 {
498 /* save negated variable */
499 SCIP_VAR* origblockvar = (SCIP_VAR*) SCIPhashmapGetImage(varsmap, (void*)SCIPvarGetNegationVar(consvars[v]));
500 SCIP_VAR* negblockvar = NULL;
501 SCIP_CALL( SCIPgetNegatedVar(blockproblem->blockscip, origblockvar, &negblockvar) );
502 blockvars[nblockvars] = negblockvar;
503 blockvals[nblockvars] = consvals[v];
504 ++nblockvars;
505 }
506 }
507 }
508
509 /* continue with next linking constraint if it has no part in current block */
510 if( nblockvars == 0 )
511 continue;
512
513 /* get rhs and/or lhs */
514 rhs = SCIPconsGetRhs(scip, linkingconss[c], success);
515 if( !(*success) )
516 {
517 SCIPdebugMsg(scip, "Create blockproblem failed\n");
518 return SCIP_OKAY;
519 }
520 lhs = SCIPconsGetLhs(scip, linkingconss[c], success);
521 if( !(*success) )
522 {
523 SCIPdebugMsg(scip, "Create blockproblem failed\n");
524 return SCIP_OKAY;
525 }
526 assert(!SCIPisInfinity(scip, rhs) || !SCIPisInfinity(scip, -lhs)); /* at least one side bounded */
527 assert(SCIPisLE(scip, lhs, rhs));
528
529 if( !SCIPisInfinity(scip, rhs) )
530 linkings[c]->hasrhs = TRUE;
531 if( !SCIPisInfinity(scip, -lhs) )
532 linkings[c]->haslhs = TRUE;
533 if( !SCIPisInfinity(scip, rhs) && !SCIPisInfinity(scip, -lhs))
534 linkings[c]->nslacksperblock = 2;
535 else
536 linkings[c]->nslacksperblock = 1;
537
538 /* add slack variable for rhs */
539 if( linkings[c]->hasrhs )
540 {
541 /* slack variable z_r >= 0 */
542 char varname[SCIP_MAXSTRLEN];
543 (void)SCIPsnprintf(varname, SCIP_MAXSTRLEN, "z_r_%s", SCIPconsGetName(linkingconss[c]));
544 SCIP_CALL( SCIPcreateVarBasic(blockproblem->blockscip, &blockvars[nblockvars], varname,
545 0.0, SCIPinfinity(scip), 1.0, SCIP_VARTYPE_CONTINUOUS) );
546 blockvals[nblockvars] = -1.0;
547 SCIP_CALL( SCIPaddVar(blockproblem->blockscip, blockvars[nblockvars]) );
548 #ifdef SCIP_MORE_DEBUG
549 SCIPdebugMsg(scip, "Add variable %s\n", SCIPvarGetName(blockvars[nblockvars]));
550 #endif
551 linkings[c]->slacks[linkings[c]->nslacks] = blockvars[nblockvars];
552 blockproblem->slackvars[blockproblem->nslackvars] = blockvars[nblockvars];
553 ++blockproblem->nslackvars;
554 ++linkings[c]->nslacks;
555 ++nblockvars;
556 }
557
558 /* add slack variable for lhs */
559 if( linkings[c]->haslhs )
560 {
561 /* slack variable z_l >= 0 */
562 char varname[SCIP_MAXSTRLEN];
563 (void)SCIPsnprintf(varname, SCIP_MAXSTRLEN, "z_l_%s", SCIPconsGetName(linkingconss[c]));
564 SCIP_CALL( SCIPcreateVarBasic(blockproblem->blockscip, &blockvars[nblockvars], varname,
565 0.0, SCIPinfinity(scip), 1.0, SCIP_VARTYPE_CONTINUOUS) );
566 blockvals[nblockvars] = 1.0;
567 SCIP_CALL( SCIPaddVar(blockproblem->blockscip, blockvars[nblockvars]) );
568 #ifdef SCIP_MORE_DEBUG
569 SCIPdebugMsg(scip, "Add variable %s\n", SCIPvarGetName(blockvars[nblockvars]));
570 #endif
571 linkings[c]->slacks[linkings[c]->nslacks] = blockvars[nblockvars];
572 blockproblem->slackvars[blockproblem->nslackvars] = blockvars[nblockvars];
573 ++blockproblem->nslackvars;
574 ++linkings[c]->nslacks;
575 ++nblockvars;
576 }
577
578 /* add linking constraint with slack variable */
579 (void)SCIPsnprintf(consname, SCIP_MAXSTRLEN, "%s", SCIPconsGetName(linkingconss[c]));
580 SCIP_CALL( SCIPcreateConsBasicLinear(blockproblem->blockscip, &newcons, consname, nblockvars, blockvars, blockvals, lhs, rhs) );
581 SCIP_CALL( SCIPaddCons(blockproblem->blockscip, newcons) );
582 #ifdef SCIP_MORE_DEBUG
583 SCIPdebugMsg(blockproblem->blockscip, "add constraint %s\n", SCIPconsGetName(newcons));
584 SCIPdebugPrintCons(blockproblem->blockscip, newcons, NULL);
585 #endif
586
587 blockproblem->linkingconss[blockproblem->nlinking] = newcons;
588 linkings[c]->blockconss[linkings[c]->nblocks] = newcons;
589 linkings[c]->blocknumbers[linkings[c]->nblocks] = blocknumber;
590 blockproblem->linkingindices[blockproblem->nlinking] = c;
591
592 /* calculate minimal und maximal activity (exclude slack variables) */
593 minact = 0;
594 maxact = 0;
595 mininfinite = FALSE;
596 maxinfinite = FALSE;
597 for( v = 0; v < nblockvars - linkings[c]->nslacksperblock && (!mininfinite || !maxinfinite); v++ )
598 {
599 SCIP_Real lb;
600 SCIP_Real ub;
601 lb = SCIPvarGetLbGlobal(blockvars[v]);
602 ub = SCIPvarGetUbGlobal(blockvars[v]);
603
604 if( blockvals[v] >= 0.0 )
605 {
606 mininfinite = (mininfinite || SCIPisInfinity(scip, -lb));
607 maxinfinite = (maxinfinite || SCIPisInfinity(scip, ub));
608 if( !mininfinite )
609 minact += blockvals[v] * lb;
610 if( !maxinfinite )
611 maxact += blockvals[v] * ub;
612 }
613 else
614 {
615 mininfinite = (mininfinite || SCIPisInfinity(scip, ub));
616 maxinfinite = (maxinfinite || SCIPisInfinity(scip, -lb));
617 if( !mininfinite )
618 minact += blockvals[v] * ub;
619 if( !maxinfinite )
620 maxact += blockvals[v] * lb;
621 }
622 }
623
624 if( mininfinite )
625 linkings[c]->minactivity[linkings[c]->nblocks] = -SCIPinfinity(scip);
626 else
627 linkings[c]->minactivity[linkings[c]->nblocks] = minact;
628 if( maxinfinite )
629 linkings[c]->maxactivity[linkings[c]->nblocks] = SCIPinfinity(scip);
630 else
631 linkings[c]->maxactivity[linkings[c]->nblocks] = maxact;
632 assert(SCIPisLE(scip, linkings[c]->minactivity[linkings[c]->nblocks], linkings[c]->maxactivity[linkings[c]->nblocks]));
633
634 linkings[c]->nblocks++;
635 blockproblem->nlinking++;
636
637 for( v = 1; v <= linkings[c]->nslacksperblock; v++ )
638 {
639 SCIP_CALL( SCIPreleaseVar(blockproblem->blockscip, &blockvars[nblockvars - v]) );
640 }
641
642 SCIP_CALL( SCIPreleaseCons(blockproblem->blockscip, &newcons) );
643 }
644 assert(blockproblem->nlinking <= nlinking);
645
646 /* free memory */
647 SCIPfreeBufferArray(blockproblem->blockscip, &consvals);
648 SCIPfreeBufferArray(blockproblem->blockscip, &consvars);
649 SCIPfreeBufferArray(blockproblem->blockscip, &blockvals);
650 SCIPfreeBufferArray(blockproblem->blockscip, &blockvars);
651
652 SCIPhashmapFree(&conssmap);
653 SCIPhashmapFree(&varsmap);
654
655 return SCIP_OKAY;
656 }
657
658 /** creates data structures and splits problem into blocks */
659 static
660 SCIP_RETCODE createAndSplitProblem(
661 SCIP* scip, /**< SCIP data structure */
662 SCIP_HEURDATA* heurdata, /**< heuristic data */
663 SCIP_DECOMP* decomp, /**< decomposition data structure */
664 BLOCKPROBLEM** blockproblem, /**< array of blockproblem data structures */
665 LINKING** linkings, /**< array of linking data structures */
666 SCIP_VAR** vars, /**< sorted array of variables */
667 SCIP_CONS** conss, /**< sorted array of constraints */
668 SCIP_Bool* success /**< pointer to store whether splitting was successful */
669 )
670 {
671 int* nconssblock;
672 int* nvarsblock;
673 int conssoffset;
674 int varsoffset;
675 int i; /* blocknumber */
676
677 assert(scip != NULL);
678 assert(heurdata != NULL);
679 assert(vars != NULL);
680 assert(conss != NULL);
681
682 SCIP_CALL( SCIPallocBufferArray(scip, &nvarsblock, heurdata->nblocks + 1) );
683 SCIP_CALL( SCIPallocBufferArray(scip, &nconssblock, heurdata->nblocks + 1) );
684 SCIP_CALL( SCIPdecompGetVarsSize(decomp, nvarsblock, heurdata->nblocks + 1) );
685 SCIP_CALL( SCIPdecompGetConssSize(decomp, nconssblock, heurdata->nblocks + 1) );
686 assert(0 == nvarsblock[0]);
687
688 varsoffset = 0;
689 conssoffset = 0;
690
691 for( i = 0; i < heurdata->nblocks; i++)
692 {
693 conssoffset += nconssblock[i];
694 varsoffset += nvarsblock[i];
695
696 SCIP_CALL( createBlockproblem(scip, blockproblem[i], linkings, &conss[conssoffset], &vars[varsoffset], nconssblock[i+1], nvarsblock[i+1],
697 heurdata->linkingconss, heurdata->nlinking, i, success) );
698 if( !(*success) )
699 break;
700 }
701
702 SCIPfreeBufferArray(scip, &nconssblock);
703 SCIPfreeBufferArray(scip, &nvarsblock);
704
705 return SCIP_OKAY;
706 }
707
708 /** rounds partition for one linking constraint to integer value if variables and coefficients are integer
709 *
710 * changes only currentrhs/currentlhs
711 */
712 static
713 SCIP_RETCODE roundPartition(
714 SCIP* scip, /**< SCIP data structure */
715 LINKING* linking, /**< one linking data structure */
716 BLOCKPROBLEM** blockproblem, /**< array of blockproblem data structures */
717 SCIP_Bool roundbyrhs /**< round by right hand side? */
718 )
719 {
720 SCIP_Real* fracPart;
721 int* sorting;
722 int* isinteger;
723 SCIP_Real sumbefor; /* includes value at idx */
724 SCIP_Real sumafter;
725 SCIP_Real diff;
726 int nnonintblocks; /* number of non integer blocks */
727 int idx;
728 int b;
729 int i;
730 int k;
731
732 assert(scip != NULL);
733 assert(linking != NULL);
734 assert(blockproblem != NULL);
735
736 nnonintblocks = 0;
737 idx = 0;
738
739 SCIP_CALL( SCIPallocBufferArray(scip, &fracPart, linking->nblocks) );
740 SCIP_CALL( SCIPallocBufferArray(scip, &sorting, linking->nblocks) );
741 SCIP_CALL( SCIPallocBufferArray(scip, &isinteger, linking->nblocks) );
742
743 /* get integer blocks and fractional parts */
744 for( b = 0; b < linking->nblocks; b++ )
745 {
746 SCIP* subscip;
747 SCIP_CONS* blockcons;
748 SCIP_VAR** blockvars;
749 SCIP_Real* blockvals;
750 int nblockvars;
751 int length; /* number of block variables without slack variables */
752 SCIP_Bool success;
753
754 subscip = blockproblem[linking->blocknumbers[b]]->blockscip;
755 blockcons = linking->blockconss[b];
756 sorting[b] = b; /* store current sorting to sort back */
757
758 SCIP_CALL( SCIPgetConsNVars(subscip, blockcons, &nblockvars, &success) );
759 assert(success);
760 SCIP_CALL( SCIPallocBufferArray(scip, &blockvars, nblockvars) );
761 SCIP_CALL( SCIPallocBufferArray(scip, &blockvals, nblockvars) );
762
763 SCIP_CALL( SCIPgetConsVars(subscip, blockcons, blockvars, nblockvars, &success) );
764 assert(success);
765 SCIP_CALL( SCIPgetConsVals(subscip, blockcons, blockvals, nblockvars, &success) );
766 assert(success);
767
768 /* get number of block variables in this constraint without slack variables */
769 length = nblockvars - linking->nslacksperblock;
770
771 /* get blocks with integer value */
772 isinteger[b] = 1;
773 for( i = 0; i < length; i++ )
774 {
775 if( SCIPvarGetType(blockvars[i]) == SCIP_VARTYPE_CONTINUOUS || !SCIPisIntegral(scip, blockvals[i]) )
776 {
777 isinteger[b] = 0;
778 nnonintblocks++;
779 break;
780 }
781 }
782
783 /* get fractional part of blockconstraints */
784 if( roundbyrhs )
785 fracPart[b] = linking->currentrhs[b] - floor(linking->currentrhs[b]); /* do not use SCIPfrac()! */
786 else
787 fracPart[b] = linking->currentlhs[b] - floor(linking->currentlhs[b]);
788
789 SCIPfreeBufferArray(scip, &blockvals);
790 SCIPfreeBufferArray(scip, &blockvars);
791 }
792
793 /* sort non-integer blocks to the front */
794 SCIPsortIntIntReal(isinteger, sorting, fracPart, linking->nblocks);
795
796 /* sort by fractional part */
797 SCIPsortRealInt(fracPart, sorting, nnonintblocks);
798 SCIPsortRealInt(&fracPart[nnonintblocks], &sorting[nnonintblocks], linking->nblocks - nnonintblocks);
799
800 /* detect blocks for rounding down and rounding up:
801 * integer blocks with small fractional parts are rounded down
802 * integer blocks with big fractional parts are rounded up
803 */
804
805 sumbefor = 0.0;
806 sumafter = 0.0;
807
808 for( i = 0; i < linking->nblocks - nnonintblocks; i++ )
809 sumafter += 1 - fracPart[nnonintblocks + i];
810
811 for( i = 0; i < linking->nblocks - nnonintblocks; i++ )
812 {
813 sumbefor += fracPart[nnonintblocks + i];
814 sumafter -= 1 - fracPart[nnonintblocks + i];
815
816 if( sumbefor >= sumafter )
817 {
818 for( k = 0; k <= i; k++ )
819 fracPart[nnonintblocks + k] = -fracPart[nnonintblocks + k];
820
821 for( k = i + 1; k < linking->nblocks - nnonintblocks; k++ )
822 fracPart[nnonintblocks + k] = 1 - fracPart[nnonintblocks + k];
823
824 idx = i;
825 break;
826 }
827 }
828 diff = sumbefor - sumafter;
829 assert(SCIPisGE(scip, diff, 0.0));
830
831 /* add difference to last non integer block */
832 for( i = nnonintblocks - 1; i >= 0; i-- )
833 {
834 if( SCIPisGT(scip, diff, 0.0) )
835 {
836 fracPart[i] = diff;
837 diff = 0;
838 }
839 else
840 fracPart[i] = 0;
841 }
842
843 /* add difference to last rounded down block if no non integer block exists */
844 if( SCIPisGT(scip, diff, 0.0))
845 {
846 assert(nnonintblocks == 0);
847 fracPart[idx] += diff;
848 }
849
850 /* sort back */
851 SCIPsortIntReal(sorting, fracPart, linking->nblocks);
852
853 /* round partition
854 * if we have a ranged constraint, both sides get rounded in the same way
855 */
856 for( b = 0; b < linking->nblocks; b++ )
857 {
858 if( linking->hasrhs )
859 linking->currentrhs[b] += fracPart[b];
860
861 if( linking->haslhs )
862 linking->currentlhs[b] += fracPart[b];
863 }
864
865 SCIPfreeBufferArray(scip, &isinteger);
866 SCIPfreeBufferArray(scip, &sorting);
867 SCIPfreeBufferArray(scip, &fracPart);
868
869 return SCIP_OKAY;
870 }
871
872 /** calculates initial partition and sets rhs/lhs in blockproblems */
873 static
874 SCIP_RETCODE initCurrent(
875 SCIP* scip, /**< SCIP data structure of main scip */
876 LINKING** linkings, /**< array of linking data structures */
877 BLOCKPROBLEM** blockproblem, /**< array of blockproblem data structures */
878 SCIP_HEURTIMING heurtiming, /**< current point in the node solving process */
879 int nlinking, /**< number of linking constraints */
880 SCIP_Bool* success /**< pointer to store whether initialization was successful */
881 )
882 {
883 LINKING* linking;
884 SCIP_Real rhs;
885 SCIP_Real lhs;
886 SCIP_Real residualrhs;
887 SCIP_Real residuallhs;
888 SCIP_Real goalvalue;
889 int c;
890 int b;
891
892 assert(scip != NULL);
893 assert(linkings != NULL);
894 assert(blockproblem != NULL);
895 assert(nlinking > 0);
896
897 SCIPdebugMsg(scip, "initialize partition\n");
898
899 /* for each linking constraint the rhs/lhs is split between the blocks */
900 for( c = 0; c < nlinking; c++ )
901 {
902 linking = linkings[c];
903 rhs = SCIPconsGetRhs(scip, linking->linkingcons, success);
904 assert(*success);
905 lhs = SCIPconsGetLhs(scip, linking->linkingcons, success);
906 assert(*success);
907 residualrhs = rhs;
908 residuallhs = lhs;
909
910 /* LP value for each block plus part of remaining amount if sum is not equal to rhs/lhs */
911 if( (heurtiming & SCIP_HEURTIMING_AFTERNODE) && (linking->hasrhs || linking->haslhs) )
912 {
913 SCIP_Real sumrhs = 0.0;
914 SCIP_Real sumlhs = 0.0;
915 for( b = 0; b < linking->nblocks; b++ )
916 {
917 SCIP_VAR** consvars;
918 SCIP_Real* consvals;
919 SCIP_CONS* cons;
920 int nconsvars;
921 SCIP_Real lpvalue;
922 int i;
923
924 /* get variables of linking constraint of one block */
925 cons = linking->blockconss[b];
926 SCIP_CALL( SCIPgetConsNVars(blockproblem[linking->blocknumbers[b]]->blockscip, cons, &nconsvars, success) );
927 SCIP_CALL( SCIPallocBufferArray(scip, &consvars, nconsvars) );
928 SCIP_CALL( SCIPallocBufferArray(scip, &consvals, nconsvars) );
929 SCIP_CALL( SCIPgetConsVars(blockproblem[linking->blocknumbers[b]]->blockscip, cons, consvars, nconsvars, success) );
930 SCIP_CALL( SCIPgetConsVals(blockproblem[linking->blocknumbers[b]]->blockscip, cons, consvals, nconsvars, success) );
931 assert(success);
932
933 /* calculate value of partition variable in lp solution */
934 lpvalue = 0.0;
935 for( i = 0; i < nconsvars - linking->nslacksperblock; i++ )
936 {
937 SCIP_VAR* origvar;
938 SCIP_Real varlpvalue;
939
940 origvar = SCIPfindVar(scip, SCIPvarGetName(consvars[i]));
941 if( origvar == NULL ) /* e.g. variable is negated */
942 {
943 *success = FALSE;
944 return SCIP_OKAY;
945 }
946 varlpvalue = SCIPvarGetLPSol(origvar);
947 lpvalue += varlpvalue * consvals[i];
948 }
949 sumrhs += lpvalue;
950 sumlhs += lpvalue;
951
952 /* set currentrhs/lhs at lp value */
953 if( linking->hasrhs )
954 linking->currentrhs[b] = lpvalue;
955 if( linking->haslhs )
956 linking->currentlhs[b] = lpvalue;
957
958 SCIPfreeBufferArray(scip, &consvars);
959 SCIPfreeBufferArray(scip, &consvals);
960 }
961 assert(SCIPisLE(scip, sumrhs, rhs));
962 assert(SCIPisGE(scip, sumlhs, lhs));
963
964 /* distribute remaining amount */
965 if( !SCIPisEQ(scip, rhs, sumrhs) && linking->hasrhs )
966 {
967 SCIP_Real diff;
968 SCIP_Real part;
969 SCIP_Real residual;
970 diff = rhs - sumrhs;
971 residual = 0.0;
972
973 for( b = 0; b < linking->nblocks; b++ )
974 {
975 goalvalue = linking->currentrhs[b] + ( diff / linking->nblocks ) + residual;
976 part = MIN(goalvalue, linking->maxactivity[b]);
977 residual = goalvalue - part;
978 linking->currentrhs[b] = part;
979 }
980 if( !SCIPisZero(scip, residual) )
981 linking->currentrhs[0] += residual;
982 }
983 if( !SCIPisEQ(scip, lhs, sumlhs) && linking->haslhs )
984 {
985 SCIP_Real diff;
986 SCIP_Real part;
987 SCIP_Real residual;
988 diff = sumlhs - lhs; /* always positive*/
989 residual = 0.0;
990
991 for( b = 0; b < linking->nblocks; b++ )
992 {
993 goalvalue = linking->currentlhs[b] - ( diff / linking->nblocks ) + residual;
994 part = MAX(goalvalue, linking->minactivity[b]);
995 residual = goalvalue - part;
996 linking->currentlhs[b] = part;
997 }
998 if( !SCIPisZero(scip, residual) )
999 linking->currentlhs[0] += residual;
1000 }
1001 }
1002 /* equal parts for each block with respect to minimal/maximal activity */
1003 else if( linking->hasrhs || linking->haslhs )
1004 {
1005 if( linking->hasrhs )
1006 {
1007 for( b = 0; b < linking->nblocks; b++ )
1008 {
1009 goalvalue = residualrhs / (linking->nblocks - b);
1010 linking->currentrhs[b] = MIN(MAX(goalvalue, linking->minactivity[b]), linking->maxactivity[b]);
1011 residualrhs -= linking->currentrhs[b];
1012 }
1013 /* add residual partition to first block */
1014 linking->currentrhs[0] += residualrhs;
1015 }
1016 if( linking->haslhs )
1017 {
1018 for( b = 0; b < linking->nblocks; b++ )
1019 {
1020 goalvalue = residuallhs / (linking->nblocks - b);
1021 linking->currentlhs[b] = MIN(MAX(goalvalue, linking->minactivity[b]), linking->maxactivity[b]);
1022 residuallhs -= linking->currentlhs[b];
1023 }
1024 /* add residual partition to first block */
1025 linking->currentlhs[0] += residuallhs;
1026 }
1027 }
1028 else
1029 {
1030 assert(linking->nblocks == 0 && !SCIPconsIsChecked(linking->linkingcons));
1031 }
1032
1033 SCIP_CALL( roundPartition(scip, linking, blockproblem, linking->hasrhs) );
1034
1035 /* set sides in blockproblem at initial partition */
1036 for( b = 0; b < linking->nblocks; b++ )
1037 {
1038 if( linking->hasrhs )
1039 {
1040 SCIP_CALL( SCIPchgRhsLinear(blockproblem[linking->blocknumbers[b]]->blockscip,
1041 linking->blockconss[b], linking->currentrhs[b]) );
1042 #ifdef SCIP_MORE_DEBUG
1043 SCIPdebugMsg(scip, "change rhs of %s in block %d to %f\n",
1044 SCIPconsGetName(linking->linkingcons), linking->blocknumbers[b], linking->currentrhs[b]);
1045 #endif
1046 }
1047 if( linking->haslhs )
1048 {
1049 SCIP_CALL( SCIPchgLhsLinear(blockproblem[linking->blocknumbers[b]]->blockscip,
1050 linking->blockconss[b], linking->currentlhs[b]) );
1051 #ifdef SCIP_MORE_DEBUG
1052 SCIPdebugMsg(scip, "change lhs of %s in block %d to %f\n",
1053 SCIPconsGetName(linking->linkingcons), linking->blocknumbers[b], linking->currentlhs[b]);
1054 #endif
1055 }
1056 }
1057 }
1058
1059 return SCIP_OKAY;
1060 }
1061
1062 /** calculates shift */
1063 static
1064 SCIP_RETCODE calculateShift(
1065 SCIP* scip, /**< SCIP data structure of main scip */
1066 BLOCKPROBLEM** blockproblem, /**< array of blockproblem data structures */
1067 LINKING* linking, /**< one linking data structure */
1068 SCIP_Real** shift, /**< pointer to store direction vector */
1069 int* nviolatedblocksrhs, /**< pointer to store number of blocks which violate rhs */
1070 int* nviolatedblockslhs, /**< pointer to store number of blocks which violate lhs */
1071 SCIP_Bool* update /**< should we update the partition? */
1072 )
1073 {
(1) Event var_decl: |
Declaring variable "nonviolatedblocksrhs" without initializer. |
Also see events: |
[uninit_use] |
1074 int* nonviolatedblocksrhs;
1075 int* nonviolatedblockslhs;
1076 SCIP_Real sumviols = 0.0;
1077 int v;
1078
(2) Event cond_false: |
Condition "linking->hasrhs", taking false branch. |
1079 if( linking->hasrhs )
1080 {
1081 SCIP_CALL( SCIPallocBufferArray(scip, &nonviolatedblocksrhs, linking->nblocks) );
(3) Event if_end: |
End of if statement. |
1082 }
(4) Event cond_true: |
Condition "linking->haslhs", taking true branch. |
1083 if( linking->haslhs )
1084 {
(5) Event cond_false: |
Condition "(nonviolatedblockslhs = BMSallocBufferMemoryArray_call(SCIPbuffer(scip), (size_t)(ptrdiff_t)linking->nblocks, 4UL /* sizeof (*nonviolatedblockslhs) */, "/sapmnt/home1/d029903/my_work/SCIPSoPlex_coverity/scip-soplex-9.0.0_rc1/scip/src/scip/heur_dps.c", 1085)) == NULL", taking false branch. |
(6) Event cond_false: |
Condition "(_restat_ = (((nonviolatedblockslhs = BMSallocBufferMemoryArray_call(SCIPbuffer(scip), (size_t)(ptrdiff_t)linking->nblocks, 4UL /* sizeof (*nonviolatedblockslhs) */, "/sapmnt/home1/d029903/my_work/SCIPSoPlex_coverity/scip-soplex-9.0.0_rc1/scip/src/scip/heur_dps.c", 1085)) == NULL) ? SCIP_NOMEMORY : SCIP_OKAY)) != SCIP_OKAY", taking false branch. |
(7) Event if_end: |
End of if statement. |
1085 SCIP_CALL( SCIPallocBufferArray(scip, &nonviolatedblockslhs, linking->nblocks) );
1086 }
1087
1088 /* get violated blocks and calculate shift */
(8) Event cond_true: |
Condition "v < linking->nblocks", taking true branch. |
(16) Event loop_begin: |
Jumped back to beginning of loop. |
(17) Event cond_false: |
Condition "v < linking->nblocks", taking false branch. |
1089 for( v = 0; v < linking->nblocks; v++ )
1090 {
1091 SCIP* subscip;
1092 SCIP_SOL* subsol;
1093 SCIP_Real slackval;
1094
1095 subscip = blockproblem[linking->blocknumbers[v]]->blockscip;
1096 subsol = SCIPgetBestSol(subscip);
1097
1098 /* if we have ranged constraints, the slack variables of the rhs are in front;
1099 * get slack variable of block; save violated blocks
1100 */
(9) Event cond_false: |
Condition "linking->hasrhs", taking false branch. |
1101 if( linking->hasrhs )
1102 {
1103 slackval = SCIPgetSolVal(subscip, subsol, linking->slacks[v * linking->nslacksperblock]);
1104
1105 /* block is violated */
1106 if( SCIPisPositive(scip, slackval) )
1107 {
1108 (*nviolatedblocksrhs)++;
1109
1110 (*shift)[v] += slackval;
1111 sumviols += slackval;
1112 }
1113 else
1114 {
1115 nonviolatedblocksrhs[v - *nviolatedblocksrhs] = v; /*lint !e644*/
1116 }
(10) Event if_end: |
End of if statement. |
1117 }
(11) Event cond_true: |
Condition "linking->haslhs", taking true branch. |
1118 if( linking->haslhs )
1119 {
1120 slackval = SCIPgetSolVal(subscip, subsol, linking->slacks[(v * linking->nslacksperblock) + linking->nslacksperblock - 1]);
1121
1122 /* block is violated */
(12) Event cond_true: |
Condition "slackval > scip->set->num_epsilon", taking true branch. |
1123 if( SCIPisPositive(scip, slackval) )
1124 {
1125 (*nviolatedblockslhs)++;
1126
1127 (*shift)[v] -= slackval;
1128 sumviols -= slackval;
(13) Event if_fallthrough: |
Falling through to end of if statement. |
1129 }
1130 else
1131 {
1132 nonviolatedblockslhs[v - *nviolatedblockslhs] = v; /*lint !e644*/
(14) Event if_end: |
End of if statement. |
1133 }
1134 }
(15) Event loop: |
Jumping back to the beginning of the loop. |
(18) Event loop_end: |
Reached end of loop. |
1135 }
1136
1137 /* linking constraint is in no block violated or always violated -> do not update partition */
(19) Event cond_false: |
Condition "*nviolatedblocksrhs + *nviolatedblockslhs == 0", taking false branch. |
(20) Event cond_false: |
Condition "linking->nblocks == *nviolatedblocksrhs", taking false branch. |
(21) Event cond_false: |
Condition "linking->nblocks == *nviolatedblockslhs", taking false branch. |
1138 if( *nviolatedblocksrhs + *nviolatedblockslhs == 0 ||
1139 linking->nblocks == *nviolatedblocksrhs || linking->nblocks == *nviolatedblockslhs )
1140 {
1141 *update = FALSE;
1142 /* free memory */
1143 if( linking->haslhs )
1144 SCIPfreeBufferArray(scip, &nonviolatedblockslhs);
1145 if( linking->hasrhs )
1146 SCIPfreeBufferArray(scip, &nonviolatedblocksrhs);
1147 return SCIP_OKAY;
(22) Event if_end: |
End of if statement. |
1148 }
1149
1150 /* set values of non violated blocks */
(23) Event cond_true: |
Condition "sumviols > scip->set->num_epsilon", taking true branch. |
1151 if( SCIPisPositive(scip, sumviols) )
1152 {
1153 /* rhs of original linking constraint is violated */
1154 SCIP_Real residual = sumviols;
1155 SCIP_Real part;
1156 SCIP_Real shift_tmp;
1157
1158 assert(linking->hasrhs);
1159 assert(*nviolatedblocksrhs != 0);
1160
1161 /* substract from each non violated block the same amount with respect to minimal/maximal activity,
1162 * so that the shift is zero in sum
1163 */
(24) Event cond_true: |
Condition "v < linking->nblocks - *nviolatedblocksrhs", taking true branch. |
1164 for( v = 0; v < linking->nblocks - *nviolatedblocksrhs; v++ )
1165 {
(25) Event uninit_use: |
Using uninitialized value "nonviolatedblocksrhs". |
Also see events: |
[var_decl] |
1166 part = linking->currentrhs[nonviolatedblocksrhs[v]] - residual/(linking->nblocks - *nviolatedblocksrhs - v);
1167 part = MIN(MAX(part, linking->minactivity[nonviolatedblocksrhs[v]]), linking->maxactivity[nonviolatedblocksrhs[v]]);
1168 shift_tmp = part - linking->currentrhs[nonviolatedblocksrhs[v]];
1169 residual += shift_tmp;
1170 (*shift)[nonviolatedblocksrhs[v]] += shift_tmp;
1171 }
1172 if( !SCIPisZero(scip, residual) )
1173 {
1174 /* assign residual to ... */
1175 if( linking->nblocks - *nviolatedblocksrhs == 1 )
1176 (*shift)[nonviolatedblocksrhs[0] == 0 ? 1 : 0] -= residual; /* first violated block */
1177 else
1178 (*shift)[nonviolatedblocksrhs[0]] -= residual; /* first nonviolated block */
1179 }
1180 }
1181 if( SCIPisNegative(scip, sumviols) )
1182 {
1183 /* lhs of original linking constraint is violated */
1184 SCIP_Real residual = sumviols;
1185 SCIP_Real part;
1186 SCIP_Real shift_tmp;
1187
1188 assert(linking->haslhs);
1189 assert(*nviolatedblockslhs != 0);
1190
1191 /* add to each non violated block the same amount with respect to minimal/maximal activity,
1192 * so that the shift is zero in sum
1193 */
1194 for( v = 0; v < linking->nblocks - *nviolatedblockslhs; v++ )
1195 {
1196 part = linking->currentlhs[nonviolatedblockslhs[v]] - residual/(linking->nblocks - *nviolatedblockslhs - v);
1197 part = MIN(MAX(part, linking->minactivity[nonviolatedblockslhs[v]]), linking->maxactivity[nonviolatedblockslhs[v]]);
1198 shift_tmp = part - linking->currentlhs[nonviolatedblockslhs[v]];
1199 residual += shift_tmp;
1200 (*shift)[nonviolatedblockslhs[v]] += shift_tmp;
1201 }
1202 if( !SCIPisZero(scip, residual) )
1203 {
1204 /* assign residual to ... */
1205 if( linking->nblocks - *nviolatedblockslhs == 1 )
1206 (*shift)[nonviolatedblockslhs[0] == 0 ? 1 : 0] -= residual; /* first violated block */
1207 else
1208 (*shift)[nonviolatedblockslhs[0]] -= residual; /* first nonviolated block */
1209 }
1210 }
1211
1212 *update = TRUE;
1213
1214 /* free memory */
1215 if( linking->haslhs )
1216 SCIPfreeBufferArray(scip, &nonviolatedblockslhs);
1217 if( linking->hasrhs )
1218 SCIPfreeBufferArray(scip, &nonviolatedblocksrhs);
1219
1220 return SCIP_OKAY;
1221 }
1222
1223 /** update partition */
1224 static
1225 SCIP_RETCODE updatePartition(
1226 SCIP* scip, /**< SCIP data structure of main scip */
1227 LINKING** linkings, /**< linking data structure */
1228 BLOCKPROBLEM** blockproblem, /**< array of blockproblem data structures */
1229 int** nviolatedblocksrhs, /**< pointer to store number of blocks which violate rhs */
1230 int** nviolatedblockslhs, /**< pointer to store number of blocks which violate lhs */
1231 int nlinking, /**< number of linking constraints */
1232 int nblocks, /**< number of blocks */
1233 int iteration, /**< number of iteration */
1234 SCIP_Bool* oneupdate /**< is at least partition of one constraint updated? */
1235 )
1236 {
1237 SCIP_Real* shift; /* direction vector */
1238 int v;
1239 int c;
1240 SCIP_Bool update; /* is partition of current constraint updated? */
1241
1242 assert(scip != NULL);
1243 assert(linkings != NULL);
1244 assert(blockproblem != NULL);
1245 assert(iteration >= 0);
1246 assert(!*oneupdate);
1247
1248 SCIP_CALL( SCIPallocBufferArray(scip, &shift, nblocks) ); /* allocate memory only once */
1249
1250 for( c = 0; c < nlinking; c++ )
1251 {
1252 LINKING* linking; /* one linking data structure */
1253
1254 linking = linkings[c];
1255 (*nviolatedblocksrhs)[c] = 0;
1256 (*nviolatedblockslhs)[c] = 0;
1257 update = TRUE;
1258 BMSclearMemoryArray(shift, linking->nblocks);
1259 assert(nblocks >= linking->nblocks);
1260
1261 SCIP_CALL( calculateShift(scip, blockproblem, linking, &shift, &(*nviolatedblocksrhs)[c], &(*nviolatedblockslhs)[c], &update) );
1262
1263 if( !update )
1264 continue;
1265
1266 *oneupdate = TRUE;
1267
1268 #ifdef SCIP_DEBUG
1269 SCIP_Real sum = 0.0;
1270 /* check sum of shift; must be zero */
1271 for( v = 0; v < linking->nblocks; v++ )
1272 sum += shift[v];
1273 assert(SCIPisFeasZero(scip, sum));
1274 #endif
1275
1276 /* add shift to both sides */
1277 for( v = 0; v < linking->nblocks; v++)
1278 {
1279 if( linking->hasrhs )
1280 linking->currentrhs[v] += shift[v];
1281
1282 if( linking->haslhs )
1283 linking->currentlhs[v] += shift[v];
1284 }
1285
1286 SCIP_CALL( roundPartition(scip, linking, blockproblem, ((linking->hasrhs && ((*nviolatedblocksrhs)[c] != 0)) || !linking->haslhs)) );
1287
1288 /* set sides in blockproblems to new partition */
1289 for( v = 0; v < linking->nblocks; v++ )
1290 {
1291 SCIP* subscip;
1292 subscip = blockproblem[linking->blocknumbers[v]]->blockscip;
1293
1294 if( linking->hasrhs )
1295 {
1296 SCIP_CALL( SCIPchgRhsLinear(subscip, linking->blockconss[v], linking->currentrhs[v]) );
1297 }
1298 if( linking->haslhs )
1299 {
1300 SCIP_CALL( SCIPchgLhsLinear(subscip, linking->blockconss[v], linking->currentlhs[v]) );
1301 }
1302 }
1303 }
1304
1305 /* free memory */
1306 SCIPfreeBufferArray(scip, &shift);
1307
1308 return SCIP_OKAY;
1309 }
1310
1311 /** update penalty parameters lambda
1312 *
1313 * if a linking constraint is violated two times in succession, the corresponding penalty parameter is increased in each block
1314 */
1315 static
1316 SCIP_RETCODE updateLambda(
1317 SCIP* scip, /**< SCIP data structure */
1318 SCIP_HEURDATA* heurdata, /**< heuristic data */
1319 LINKING** linkings, /**< array of linking data structures */
1320 BLOCKPROBLEM** blockproblem, /**< array of blockproblem data structures */
1321 int* nviolatedblocksrhs, /**< number of blocks which violate rhs */
1322 int* nviolatedblockslhs, /**< number of blocks which violate lhs */
1323 int nlinking /**< number of linking constraints */
1324 )
1325 {
1326 SCIP_VAR* slackvar;
1327 SCIP_Real old_obj;
1328 SCIP_Real new_obj;
1329 int c;
1330 int b;
1331
1332 assert(scip != NULL);
1333 assert(linkings != NULL);
1334 assert(blockproblem != NULL);
1335
1336 for( c = 0; c < nlinking; c++ )
1337 {
1338 assert(nviolatedblocksrhs[c] >= 0);
1339 assert(nviolatedblockslhs[c] >= 0);
1340
1341 /* skip constraint if it is not violated */
1342 if( nviolatedblocksrhs[c] + nviolatedblockslhs[c] == 0 )
1343 {
1344 linkings[c]->lastviolations = 0; /* reset flag */
1345 continue;
1346 }
1347
1348 /* add number of violated blocks multiplied with parameter "penalty" to lambda (initial value is 1) */
1349 for( b = 0; b < linkings[c]->nblocks; b++ )
1350 {
1351 SCIP* subscip;
1352 subscip = blockproblem[linkings[c]->blocknumbers[b]]->blockscip;
1353 assert(subscip != NULL);
1354
1355 if( linkings[c]->hasrhs && (nviolatedblocksrhs[c] >= 1) && (linkings[c]->lastviolations >= 1) )
1356 {
1357 slackvar = linkings[c]->slacks[b * linkings[c]->nslacksperblock];
1358 old_obj = SCIPvarGetObj(slackvar);
1359 new_obj = old_obj + heurdata->penalty * nviolatedblocksrhs[c];
1360
1361 SCIP_CALL( SCIPchgVarObj(subscip, slackvar, new_obj) );
1362 }
1363 if( linkings[c]->haslhs && (nviolatedblockslhs[c] >= 1) && (linkings[c]->lastviolations >= 1) )
1364 {
1365 slackvar = linkings[c]->slacks[b * linkings[c]->nslacksperblock + linkings[c]->nslacksperblock - 1];
1366 old_obj = SCIPvarGetObj(slackvar);
1367 new_obj = old_obj + heurdata->penalty * nviolatedblockslhs[c];
1368
1369 SCIP_CALL( SCIPchgVarObj(subscip, slackvar, new_obj) );
1370 }
1371 }
1372
1373 /* update number of violations in the last iterations */
1374 linkings[c]->lastviolations += 1;
1375 }
1376
1377 return SCIP_OKAY;
1378 }
1379
1380 /** computes feasible solution from last stored solution for each block and adds it to the solution storage */
1381 static
1382 SCIP_RETCODE reuseSolution(
1383 LINKING** linkings, /**< array of linking data structures */
1384 BLOCKPROBLEM** blockproblem, /**< array of blockproblem data structures */
1385 int nblocks /**< number of blocks */
1386 )
1387 {
1388 SCIP_SOL** sols;
1389 SCIP_SOL* sol; /* solution of block that will be repaired */
1390 SCIP_SOL* newsol;
1391 SCIP_VAR** blockvars;
1392 SCIP_Real* blockvals;
1393 int nsols;
1394 int nvars;
1395 int b;
1396 int c;
1397 int i;
1398 SCIP_Bool success;
1399
1400 assert(linkings != NULL);
1401 assert(blockproblem != NULL);
1402
1403 for( b = 0; b < nblocks; b++ )
1404 {
1405 SCIP* subscip;
1406
1407 subscip = blockproblem[b]->blockscip;
1408 nsols = SCIPgetNSols(subscip);
1409
1410 /* no solution in solution candidate storage found */
1411 if( nsols == 0 )
1412 return SCIP_OKAY;
1413
1414 /* take last solution */
1415 sols = SCIPgetSols(subscip);
1416 sol = sols[nsols - 1];
1417
1418 /* copy the solution */
1419 nvars = SCIPgetNVars(subscip);
1420 blockvars = SCIPgetVars(subscip);
1421 SCIP_CALL( SCIPallocBufferArray(subscip, &blockvals, nvars) );
1422 SCIP_CALL( SCIPgetSolVals(subscip, sol, nvars, blockvars, blockvals) );
1423 SCIP_CALL( SCIPcreateOrigSol(subscip, &newsol, NULL) );
1424 SCIP_CALL( SCIPsetSolVals(subscip, newsol, nvars, blockvars, blockvals) );
1425
1426 /* correct each coupling constraint:
1427 * lhs <= orig_var - z_r + z_l <= rhs
1428 * adapt slack variables so that constraint is feasible
1429 */
1430 for( c = 0; c < blockproblem[b]->nlinking; c++ )
1431 {
1432 LINKING* linking; /* linking data structure of this constraint */
1433 SCIP_VAR* rvar; /* slack variable z_r */
1434 SCIP_VAR* lvar; /* slack variable z_l */
1435 SCIP_Real rval; /* value of slack variable z_r */
1436 SCIP_Real lval; /* value of slack variable z_l */
1437 SCIP_Real activitycons; /* activity of constraint*/
1438 SCIP_Real activity; /* activity of constraint without slack variables */
1439 SCIP_Real rhs; /* current right hand side */
1440 SCIP_Real lhs; /* current left hand side */
1441
1442 linking = linkings[blockproblem[b]->linkingindices[c]];
1443 rhs = SCIPgetRhsLinear(subscip, blockproblem[b]->linkingconss[c]);
1444 lhs = SCIPgetLhsLinear(subscip, blockproblem[b]->linkingconss[c]);
1445 assert(SCIPisGE(subscip, rhs, lhs));
1446
1447 activitycons = SCIPgetActivityLinear(subscip, blockproblem[b]->linkingconss[c], sol);
1448
1449 /* get slack variables and subtract their value from the activity;
1450 * calculate and set values of slack variables
1451 */
1452 for( i = 0; i < linking->nblocks; i++ )
1453 {
1454 if( linking->blocknumbers[i] == b )
1455 {
1456 if( linking->hasrhs && linking->haslhs )
1457 {
1458 rvar = linking->slacks[2 * i];
1459 lvar = linking->slacks[2 * i + 1];
1460 rval = SCIPgetSolVal(subscip, sol, rvar);
1461 lval = SCIPgetSolVal(subscip, sol, lvar);
1462 activity = activitycons + rval - lval;
1463 SCIP_CALL( SCIPsetSolVal(subscip, newsol, rvar, MAX(0.0, activity - rhs)) );
1464 SCIP_CALL( SCIPsetSolVal(subscip, newsol, lvar, MAX(0.0, lhs - activity)) );
1465 }
1466 else if( linking->hasrhs )
1467 {
1468 rvar = linking->slacks[i];
1469 rval = SCIPgetSolVal(subscip, sol, rvar);
1470 activity = activitycons + rval;
1471 SCIP_CALL( SCIPsetSolVal(subscip, newsol, rvar, MAX(0.0, activity - rhs)) );
1472 }
1473 else /* linking->haslhs */
1474 {
1475 assert(linking->haslhs);
1476 lvar = linking->slacks[i];
1477 lval = SCIPgetSolVal(subscip, sol, lvar);
1478 activity = activitycons - lval;
1479 SCIP_CALL( SCIPsetSolVal(subscip, newsol, lvar, MAX(0.0, lhs - activity)) );
1480 }
1481 break;
1482 }
1483 }
1484 }
1485
1486 SCIPdebugMsg(subscip, "Try adding solution with objective value %.2f\n", SCIPgetSolOrigObj(subscip, newsol));
1487 SCIP_CALL( SCIPaddSolFree(subscip, &newsol, &success) );
1488
1489 if( !success )
1490 SCIPdebugMsg(subscip, "Correcting solution failed\n"); /* maybe not better than old solutions */
1491 else
1492 SCIPdebugMsg(subscip, "Correcting solution successful\n");
1493
1494 SCIPfreeBufferArray(subscip, &blockvals);
1495 }
1496
1497 return SCIP_OKAY;
1498 }
1499
1500 /** reoptimizes the heuristic solution with original objective function */
1501 static
1502 SCIP_RETCODE reoptimize(
1503 SCIP* scip, /**< SCIP data structure */
1504 SCIP_HEUR* heur, /**< pointer to heuristic */
1505 SCIP_SOL* sol, /**< heuristic solution */
1506 BLOCKPROBLEM** blockproblem, /**< array of blockproblem data structures */
1507 int nblocks, /**< number of blockproblems */
1508 SCIP_Bool limits, /**< should strict limits be set? */
1509 SCIP_SOL** newsol, /**< pointer to store improved solution */
1510 SCIP_Bool* success /**< pointer to store whether reoptimization was successful */
1511 )
1512 {
1513 SCIP_Real time;
1514 SCIP_Real timesubscip;
1515 SCIP_Bool check;
1516 int b;
1517 int v;
1518
1519 assert(scip != NULL);
1520 assert(heur != NULL);
1521 assert(sol != NULL);
1522 assert(blockproblem != NULL);
1523
1524 *success = FALSE;
1525 check = FALSE;
1526
1527 /* get used time */
1528 timesubscip = SCIPgetTotalTime(blockproblem[0]->blockscip);
1529
1530 /* for each blockproblem:
1531 * - change back to original objective function
1532 * - fix slack variables to zero
1533 * - set limits and solve problem
1534 */
1535 for( b = 0; b < nblocks; b++ )
1536 {
1537 SCIP* subscip;
1538 SCIP_VAR** blockvars;
1539 SCIP_Real infvalue;
1540 int nvars;
1541 int nsols;
1542
1543 subscip = blockproblem[b]->blockscip;
1544 blockvars = SCIPgetOrigVars(subscip);
1545 nvars = SCIPgetNOrigVars(subscip);
1546 nsols = 0;
1547
1548 /* in order to change objective function */
1549 SCIP_CALL( SCIPfreeTransform(subscip) );
1550
1551 /* change back to original objective function */
1552 for( v = 0; v < blockproblem[b]->nblockvars; v++ )
1553 {
1554 SCIP_CALL( SCIPchgVarObj(subscip, blockvars[v], blockproblem[b]->origobj[v]) );
1555 }
1556
1557 /* fix slack variables to zero */
1558 for( v = blockproblem[b]->nblockvars; v < nvars; v++ )
1559 {
1560 SCIP_CALL( SCIPchgVarUb(subscip, blockvars[v], 0.0) );
1561 SCIP_CALL( SCIPchgVarLb(subscip, blockvars[v], 0.0) );
1562 }
1563
1564 /* do not abort subproblem on CTRL-C */
1565 SCIP_CALL( SCIPsetBoolParam(subscip, "misc/catchctrlc", FALSE) );
1566
1567 /* forbid recursive call of heuristics and separators solving sub-SCIPs */
1568 SCIP_CALL( SCIPsetSubscipsOff(subscip, TRUE) );
1569
1570 #ifdef SCIP_DEBUG
1571 /* for debugging, enable full output */
1572 SCIP_CALL( SCIPsetIntParam(subscip, "display/verblevel", 5) );
1573 SCIP_CALL( SCIPsetIntParam(subscip, "display/freq", 100000000) );
1574 #else
1575 /* disable statistic timing inside sub SCIP and output to console */
1576 SCIP_CALL( SCIPsetIntParam(subscip, "display/verblevel", 0) );
1577 SCIP_CALL( SCIPsetBoolParam(subscip, "timing/statistictiming", FALSE) );
1578 #endif
1579
1580 /* disable cutting plane separation */
1581 SCIP_CALL( SCIPsetSeparating(subscip, SCIP_PARAMSETTING_OFF, TRUE) );
1582
1583 /* disable expensive presolving */
1584 SCIP_CALL( SCIPsetPresolving(subscip, SCIP_PARAMSETTING_FAST, TRUE) );
1585
1586 /* disable expensive techniques */
1587 SCIP_CALL( SCIPsetIntParam(subscip, "misc/usesymmetry", 0) );
1588
1589 /* speed up sub-SCIP by not checking dual LP feasibility */
1590 SCIP_CALL( SCIPsetBoolParam(subscip, "lp/checkdualfeas", FALSE) );
1591
1592 /* copy value for infinity */
1593 SCIP_CALL( SCIPgetRealParam(scip, "numerics/infinity", &infvalue) );
1594 SCIP_CALL( SCIPsetRealParam(subscip, "numerics/infinity", infvalue) );
1595
1596 /* set limits; do not use more time in each subproblem than the heuristic has already used for first solution */
1597 SCIP_CALL( SCIPcopyLimits(scip, subscip) );
1598 if( limits )
1599 {
1600 SCIP_CALL( SCIPsetLongintParam(subscip, "limits/nodes", 1LL) );
1601 SCIP_CALL( SCIPgetRealParam(subscip, "limits/time", &time) );
1602 if( timesubscip < time - 1.0 )
1603 {
1604 SCIP_CALL( SCIPsetRealParam(subscip, "limits/time", timesubscip + 1.0) );
1605 }
1606 SCIP_CALL( SCIPtransformProb(subscip) );
1607 nsols = SCIPgetNSols(subscip);
1608 /* one improving solution */
1609 SCIP_CALL( SCIPsetIntParam(subscip, "limits/bestsol", nsols + 1) );
1610 }
1611 else
1612 {
1613 /* only 50% of remaining time */
1614 SCIP_CALL( SCIPgetRealParam(subscip, "limits/time", &time) );
1615 SCIP_CALL( SCIPsetRealParam(subscip, "limits/time", time * 0.5) );
1616 /* set big node limits */
1617 SCIP_CALL( SCIPsetLongintParam(subscip, "limits/nodes", 1000LL) );
1618 SCIP_CALL( SCIPsetLongintParam(subscip, "limits/stallnodes", 100LL) );
1619 }
1620
1621 /* reoptimize problem */
1622 SCIP_CALL_ABORT( SCIPsolve(subscip) );
1623
1624 if( SCIPgetNSols(subscip) == 0 )
1625 {
1626 /* we found no solution */
1627 return SCIP_OKAY;
1628 }
1629 else if( SCIPgetStatus(subscip) == SCIP_STATUS_BESTSOLLIMIT ||
1630 SCIPgetStatus(subscip) == SCIP_STATUS_OPTIMAL ||
1631 SCIPgetNSols(subscip) > nsols )
1632 {
1633 check = TRUE;
1634 }
1635 }
1636
1637 if( !check )
1638 {
1639 /* we have no better solution */
1640 return SCIP_OKAY;
1641 }
1642
1643 /* create sol of main scip */
1644 SCIP_CALL( SCIPcreateSol(scip, newsol, heur) );
1645
1646 /* copy solution to main scip */
1647 for( b = 0; b < nblocks; b++ )
1648 {
1649 SCIP_SOL* blocksol;
1650 SCIP_VAR** blockvars;
1651 SCIP_Real* blocksolvals;
1652 int nblockvars;
1653
1654 /* get solution of block variables (without slack variables) */
1655 blocksol = SCIPgetBestSol(blockproblem[b]->blockscip);
1656 blockvars = SCIPgetOrigVars(blockproblem[b]->blockscip);
1657 nblockvars = blockproblem[b]->nblockvars;
1658
1659 SCIP_CALL( SCIPallocBufferArray(scip, &blocksolvals, nblockvars) );
1660 SCIP_CALL( SCIPgetSolVals(blockproblem[b]->blockscip, blocksol, nblockvars, blockvars, blocksolvals) );
1661
1662 for( v = 0; v < nblockvars; v++ )
1663 {
1664 SCIP_VAR* origvar;
1665
1666 origvar = SCIPfindVar(scip, SCIPvarGetName(blockvars[v]));
1667 SCIP_CALL( SCIPsetSolVal(scip, *newsol, origvar, blocksolvals[v]) );
1668 }
1669
1670 SCIPfreeBufferArray(scip, &blocksolvals);
1671 }
1672
1673 *success = TRUE;
1674
1675 return SCIP_OKAY;
1676 }
1677
1678
1679 /* ---------------- Callback methods of event handler ---------------- */
1680
1681 /* exec the event handler
1682 *
1683 * interrupt solution process of sub-SCIP if dual bound is greater than zero and a solution is available
1684 */
1685 static
1686 SCIP_DECL_EVENTEXEC(eventExecDps)
1687 {
1688 assert(eventhdlr != NULL);
1689 assert(eventdata != NULL);
1690 assert(strcmp(SCIPeventhdlrGetName(eventhdlr), EVENTHDLR_NAME) == 0);
1691 assert(event != NULL);
1692 assert(SCIPeventGetType(event) & SCIP_EVENTTYPE_LPSOLVED);
1693
1694 SCIPdebugMsg(scip, "dual bound: %0.2f\n", SCIPgetDualbound(scip));
1695
1696 if( SCIPisFeasGT(scip, SCIPgetDualbound(scip), 0.0) && SCIPgetNSols(scip) >= 1 )
1697 {
1698 SCIPdebugMsg(scip, "DPS: interrupt subscip\n");
1699 SCIP_CALL( SCIPinterruptSolve(scip) );
1700 }
1701
1702 return SCIP_OKAY;
1703 }
1704
1705
1706 /*
1707 * Callback methods of primal heuristic
1708 */
1709
1710 /** copy method for primal heuristic plugins (called when SCIP copies plugins) */
1711 static
1712 SCIP_DECL_HEURCOPY(heurCopyDps)
1713 { /*lint --e{715}*/
1714 assert(scip != NULL);
1715 assert(heur != NULL);
1716 assert(strcmp(SCIPheurGetName(heur), HEUR_NAME) == 0);
1717
1718 /* call inclusion method of primal heuristic */
1719 SCIP_CALL( SCIPincludeHeurDps(scip) );
1720
1721 return SCIP_OKAY;
1722 }
1723
1724 /** destructor of primal heuristic to free user data (called when SCIP is exiting) */
1725 static
1726 SCIP_DECL_HEURFREE(heurFreeDps)
1727 { /*lint --e{715}*/
1728 SCIP_HEURDATA* heurdata;
1729
1730 assert(heur != NULL);
1731 assert(scip != NULL);
1732
1733 /* free heuristic data */
1734 heurdata = SCIPheurGetData(heur);
1735 assert(heurdata != NULL);
1736
1737 SCIPfreeBlockMemory(scip, &heurdata);
1738 SCIPheurSetData(heur, NULL);
1739
1740 return SCIP_OKAY;
1741 }
1742
1743 /** execution method of primal heuristic */
1744 static
1745 SCIP_DECL_HEUREXEC(heurExecDps)
1746 { /*lint --e{715}*/
1747 SCIP_DECOMP** alldecomps;
1748 SCIP_DECOMP* decomp;
1749 SCIP_DECOMP* assigneddecomp;
1750 SCIP_VAR** vars;
1751 SCIP_CONS** conss;
1752 SCIP_VAR** sortedvars;
1753 SCIP_CONS** sortedconss;
1754 SCIP_HEURDATA* heurdata;
1755 BLOCKPROBLEM** blockproblem;
1756 LINKING** linkings;
1757 int* sortedvarlabels;
1758 int* sortedconslabels;
1759 SCIP_EVENTHDLR* eventhdlr; /* event handler */
1760 SCIP_Real memory; /* in MB */
1761 SCIP_Real timelimit;
1762 SCIP_Real allslacksval;
1763 SCIP_Real blocksolval;
1764 SCIP_STATUS status;
1765 SCIP_Bool avoidmemout;
1766 SCIP_Bool disablemeasures;
1767 int maxgraphedge;
1768 int ndecomps;
1769 int nvars;
1770 int nconss;
1771 int nblocks;
1772 SCIP_Bool success;
1773 int b;
1774 int c;
1775 int k;
1776
1777 assert( heur != NULL );
1778 assert( scip != NULL );
1779 assert( result != NULL );
1780
1781 heurdata = SCIPheurGetData(heur);
1782 assert(heurdata != NULL);
1783
1784 assigneddecomp = NULL;
1785 blockproblem = NULL;
1786 linkings = NULL;
1787 eventhdlr = NULL;
1788
1789 *result = SCIP_DIDNOTRUN;
1790
1791 if( !((heurtiming & SCIP_HEURTIMING_BEFORENODE) && (heurdata->timing == 0 || heurdata->timing == 2))
1792 && !((heurtiming & SCIP_HEURTIMING_AFTERNODE) && (heurdata->timing == 1 || heurdata->timing == 2)) )
1793 {
1794 return SCIP_OKAY;
1795 }
1796
1797 /* call heuristic only once */
1798 if( SCIPheurGetNCalls(heur) > 0 )
1799 return SCIP_OKAY;
1800
1801 /* -------------------------------------------------------------------- */
1802 SCIPdebugMsg(scip, "initialize dps heuristic\n");
1803
1804 /* take the first transformed decomposition */
1805 SCIPgetDecomps(scip, &alldecomps, &ndecomps, FALSE);
1806 if( ndecomps == 0 )
1807 return SCIP_OKAY;
1808
1809 decomp = alldecomps[0];
1810 assert(decomp != NULL);
1811 SCIPdebugMsg(scip, "First transformed decomposition is selected\n");
1812
1813 nblocks = SCIPdecompGetNBlocks(decomp);
1814 nconss = SCIPgetNConss(scip);
1815 nvars = SCIPgetNVars(scip);
1816
1817 /* if problem has no constraints, no variables or less than two blocks, return */
1818 if( nconss == 0 || nvars == 0 || nblocks <= 1 )
1819 {
1820 SCIPdebugMsg(scip, "problem has no constraints, no variables or less than two blocks\n");
1821 return SCIP_OKAY;
1822 }
1823
1824 /* estimate required memory for all blocks and terminate if not enough memory is available */
1825 SCIP_CALL( SCIPgetRealParam(scip, "limits/memory", &memory) );
1826 SCIP_CALL( SCIPgetBoolParam(scip, "misc/avoidmemout", &avoidmemout) );
1827 if( avoidmemout && (((SCIPgetMemUsed(scip) + SCIPgetMemExternEstim(scip))/1048576.0) * (nblocks/4.0 + 2) >= memory) )
1828 {
1829 SCIPdebugMsg(scip, "The estimated memory usage for %d blocks is too large.\n", nblocks);
1830 return SCIP_OKAY;
1831 }
1832
1833 /* we do not need the block decomposition graph and expensive measures of the decomposition statistics */
1834 SCIP_CALL( SCIPgetIntParam(scip, "decomposition/maxgraphedge", &maxgraphedge) );
1835 if( !SCIPisParamFixed(scip, "decomposition/maxgraphedge") )
1836 {
1837 SCIP_CALL( SCIPsetIntParam(scip, "decomposition/maxgraphedge", 0) );
1838 }
1839 SCIP_CALL( SCIPgetBoolParam(scip, "decomposition/disablemeasures", &disablemeasures) );
1840 if( !SCIPisParamFixed(scip, "decomposition/disablemeasures") )
1841 {
1842 SCIP_CALL( SCIPsetBoolParam(scip, "decomposition/disablemeasures", TRUE) );
1843 }
1844
1845 vars = SCIPgetVars(scip);
1846 conss = SCIPgetConss(scip);
1847 SCIP_CALL( SCIPduplicateBufferArray(scip, &sortedvars, vars, nvars) );
1848 SCIP_CALL( SCIPduplicateBufferArray(scip, &sortedconss, conss, nconss) );
1849 SCIP_CALL( SCIPallocBufferArray(scip, &sortedvarlabels, nvars) );
1850 SCIP_CALL( SCIPallocBufferArray(scip, &sortedconslabels, nconss) );
1851
1852 /* get labels and sort in increasing order */
1853 SCIPdecompGetVarsLabels(decomp, sortedvars, sortedvarlabels, nvars);
1854 SCIPdecompGetConsLabels(decomp, sortedconss, sortedconslabels, nconss);
1855 SCIPsortIntPtr(sortedvarlabels, (void**)sortedvars, nvars);
1856 SCIPsortIntPtr(sortedconslabels, (void**)sortedconss, nconss);
1857
1858 if( sortedvarlabels[0] == SCIP_DECOMP_LINKVAR )
1859 {
1860 /* create new decomposition; don't change the decompositions in the decompstore */
1861 SCIP_CALL( SCIPdecompCreate(&assigneddecomp, SCIPblkmem(scip), nblocks, SCIPdecompIsOriginal(decomp), SCIPdecompUseBendersLabels(decomp)) );
1862
1863 SCIP_CALL( assignLinking(scip, assigneddecomp, sortedvars, sortedconss, sortedvarlabels, sortedconslabels, nvars, nconss, SCIPdecompGetNBorderVars(decomp)) );
1864 assert(SCIPdecompGetNBlocks(decomp) >= SCIPdecompGetNBlocks(assigneddecomp));
1865 decomp = assigneddecomp;
1866
1867 /* number of blocks can get smaller */
1868 nblocks = SCIPdecompGetNBlocks(decomp);
1869 }
1870 else
1871 {
1872 /* The decomposition statistics were computed during transformation of the decomposition store.
1873 * Since propagators can have changed the number of constraints/variables,
1874 * the statistics are no longer up-to-date and have to be recomputed.
1875 */
1876 SCIP_CALL( SCIPcomputeDecompStats(scip, decomp, TRUE) );
1877 nblocks = SCIPdecompGetNBlocks(decomp);
1878 }
1879
1880 /* reset parameters */
1881 SCIP_CALL( SCIPsetIntParam(scip, "decomposition/maxgraphedge", maxgraphedge) );
1882 SCIP_CALL( SCIPsetBoolParam(scip, "decomposition/disablemeasures", disablemeasures) );
1883
1884 #ifdef SCIP_DEBUG
1885 char buffer[SCIP_MAXSTRLEN];
1886 SCIPdebugMsg(scip, "DPS used decomposition:\n%s\n", SCIPdecompPrintStats(decomp, buffer));
1887 #endif
1888
1889 /* check properties of used decomposition */
1890 if( heurdata->maxlinkscore != 1.0 )
1891 {
1892 SCIP_Real linkscore;
1893 int nlinkconss;
1894
1895 nlinkconss = SCIPdecompGetNBorderConss(decomp);
1896
1897 linkscore = (SCIP_Real)nlinkconss / (SCIP_Real)nconss;
1898
1899 if( linkscore > heurdata->maxlinkscore )
1900 {
1901 SCIPdebugMsg(scip, "decomposition has not the required properties\n");
1902 goto TERMINATE;
1903 }
1904 }
1905
1906 if( sortedvarlabels[0] == SCIP_DECOMP_LINKVAR ||
1907 sortedconslabels[0] != SCIP_DECOMP_LINKCONS ||
1908 nblocks <= 1 )
1909 {
1910 SCIPdebugMsg(scip, "Problem has linking variables or no linking constraints or less than two blocks\n");
1911 goto TERMINATE;
1912 }
1913
1914 /* initialize heurdata */
1915 heurdata->linkingconss = sortedconss;
1916 heurdata->nlinking = SCIPdecompGetNBorderConss(decomp);
1917 heurdata->nblocks = nblocks;
1918
1919 /* allocate memory for blockproblems and initialize partially */
1920 SCIP_CALL( SCIPallocBufferArray(scip, &blockproblem, nblocks) );
1921 for( b = 0; b < nblocks; b++ )
1922 {
1923 SCIP_CALL( SCIPallocBlockMemory(scip, &(blockproblem[b])) ); /*lint !e866*/
1924 SCIP_CALL( createSubscip(scip, &blockproblem[b]->blockscip) );
1925
1926 SCIP_CALL( SCIPallocBufferArray(scip, &blockproblem[b]->linkingconss, heurdata->nlinking) );
1927 SCIP_CALL( SCIPallocBufferArray(scip, &blockproblem[b]->linkingindices, heurdata->nlinking) );
1928 SCIP_CALL( SCIPallocBufferArray(scip, &blockproblem[b]->slackvars, heurdata->nlinking * 2) ); /* maximum two slacks per linking constraint */
1929 SCIP_CALL( SCIPallocBufferArray(scip, &blockproblem[b]->origobj, nvars) );
1930 blockproblem[b]->nblockvars = 0;
1931 blockproblem[b]->nlinking = 0;
1932 blockproblem[b]->nslackvars = 0;
1933 }
1934
1935 /* allocate memory for linkings and initialize partially */
1936 SCIP_CALL( SCIPallocBufferArray(scip, &linkings, heurdata->nlinking) );
1937 for( c = 0; c < heurdata->nlinking; c++ )
1938 {
1939 SCIP_CALL( SCIPallocBlockMemory(scip, &(linkings[c])) ); /*lint !e866*/
1940 SCIP_CALL( SCIPallocBufferArray(scip, &(linkings[c])->blockconss, heurdata->nblocks) );
1941 SCIP_CALL( SCIPallocBufferArray(scip, &(linkings[c])->slacks, heurdata->nblocks*2) ); /* maximum two slacks per block */
1942 SCIP_CALL( SCIPallocBufferArray(scip, &(linkings[c])->blocknumbers, heurdata->nblocks) );
1943 SCIP_CALL( SCIPallocBufferArray(scip, &(linkings[c])->minactivity, heurdata->nblocks) );
1944 SCIP_CALL( SCIPallocBufferArray(scip, &(linkings[c])->maxactivity, heurdata->nblocks) );
1945
1946 linkings[c]->linkingcons = heurdata->linkingconss[c];
1947 linkings[c]->currentrhs = NULL;
1948 linkings[c]->currentlhs = NULL;
1949 linkings[c]->nblocks = 0;
1950 linkings[c]->nslacks = 0;
1951 linkings[c]->nslacksperblock = 0;
1952 linkings[c]->lastviolations = 0;
1953 linkings[c]->hasrhs = FALSE;
1954 linkings[c]->haslhs = FALSE;
1955 }
1956
1957 SCIP_CALL( createAndSplitProblem(scip, heurdata, decomp, blockproblem, linkings, sortedvars, sortedconss, &success) );
1958 if( !success )
1959 {
1960 SCIPdebugMsg(scip, "Create and split problem failed\n");
1961 goto TERMINATE;
1962 }
1963
1964 /* allocate memory for current partition*/
1965 for( c = 0; c < heurdata->nlinking; c++ )
1966 {
1967 if( linkings[c]->hasrhs )
1968 {
1969 SCIP_CALL( SCIPallocBufferArray(scip, &(linkings[c])->currentrhs, linkings[c]->nblocks ) );
1970 }
1971
1972 if( linkings[c]->haslhs )
1973 {
1974 SCIP_CALL( SCIPallocBufferArray(scip, &(linkings[c])->currentlhs, linkings[c]->nblocks ) );
1975 }
1976 }
1977
1978 /* initialize partition */
1979 SCIP_CALL( initCurrent(scip, linkings, blockproblem, heurtiming, heurdata->nlinking, &success) );
1980 if( !success )
1981 {
1982 SCIPdebugMsg(scip, "Initialization of partition failed\n");
1983 goto TERMINATE;
1984 }
1985
1986 /* ------------------------------------------------------------------------ */
1987 SCIPdebugMsg(scip, "Start heuristik DPS\n");
1988 *result = SCIP_DIDNOTFIND;
1989
1990 for( k = 0; k < heurdata->maxit; k++ )
1991 {
1992 /* do not exceed the timelimit */
1993 SCIP_CALL( SCIPgetRealParam(scip, "limits/time", &timelimit) );
1994 if( (timelimit - SCIPgetSolvingTime(scip)) <= 0 )
1995 goto TERMINATE;
1996
1997 /* solve the subproblems */
1998 allslacksval = 0.0;
1999 for( b = 0; b < heurdata->nblocks; b++ )
2000 {
2001 SCIP* subscip;
2002 subscip = blockproblem[b]->blockscip;
2003
2004 /* update time and memory limit of subproblem */
2005 SCIP_CALL( SCIPcopyLimits(scip, subscip) );
2006
2007 /* create event handler for LP events */
2008 if( k == 0 )
2009 {
2010 SCIP_CALL( SCIPincludeEventhdlrBasic(subscip, &eventhdlr, EVENTHDLR_NAME, EVENTHDLR_DESC, eventExecDps, NULL) );
2011 if( eventhdlr == NULL )
2012 {
2013 SCIPerrorMessage("event handler for " HEUR_NAME " heuristic not found.\n");
2014 return SCIP_PLUGINNOTFOUND;
2015 }
2016 }
2017
2018 /* catch LP events of sub-SCIP */
2019 SCIP_CALL( SCIPtransformProb(subscip) );
2020 SCIP_CALL( SCIPcatchEvent(subscip, SCIP_EVENTTYPE_LPSOLVED, eventhdlr, (SCIP_EVENTDATA*) heurdata, NULL) );
2021
2022 SCIPdebugMsg(scip, "Solve blockproblem %d\n", b);
2023 SCIP_CALL_ABORT( SCIPsolve(subscip) );
2024
2025 /* drop LP events of sub-SCIP */
2026 SCIP_CALL( SCIPdropEvent(subscip, SCIP_EVENTTYPE_LPSOLVED, eventhdlr, (SCIP_EVENTDATA*) heurdata, -1) );
2027
2028 /* get status and objective value if available */
2029 status = SCIPgetStatus(subscip);
2030 if( status == SCIP_STATUS_INFEASIBLE )
2031 {
2032 SCIPdebugMsg(scip, "Subproblem is infeasible\n");
2033 goto TERMINATE;
2034 }
2035 else if( status == SCIP_STATUS_UNBOUNDED )
2036 {
2037 SCIPdebugMsg(scip, "Subproblem is unbounded\n");
2038 goto TERMINATE;
2039 }
2040 else if( SCIPgetNSols(subscip) >= 1 )
2041 {
2042 blocksolval = SCIPgetPrimalbound(subscip);
2043
2044 if( status == SCIP_STATUS_TIMELIMIT && !SCIPisZero(scip, blocksolval) )
2045 {
2046 SCIPdebugMsg(scip, "Subproblem reached timelimit without optimal solution\n");
2047 goto TERMINATE;
2048 }
2049 SCIPdebugMsg(scip, "Solution value: %f\n", blocksolval);
2050 allslacksval += blocksolval;
2051 }
2052 else
2053 {
2054 SCIPdebugMsg(scip, "No subproblem solution available\n");
2055 goto TERMINATE;
2056 }
2057 }
2058
2059 /* all slack variables are zero -> we found a feasible solution */
2060 if( SCIPisZero(scip, allslacksval) )
2061 {
2062 SCIP_SOL* newsol;
2063
2064 SCIPdebugMsg(scip, "Feasible solution found after %i iterations\n", k);
2065
2066 /* create new solution */
2067 SCIP_CALL( SCIPcreateSol(scip, &newsol, heur) );
2068 for( b = 0; b < heurdata->nblocks; b++ )
2069 {
2070 SCIP_SOL* blocksol;
2071 SCIP_VAR** blockvars;
2072 SCIP_Real* blocksolvals;
2073 int nblockvars;
2074
2075 /* get solution of block variables (without slack variables) */
2076 blocksol = SCIPgetBestSol(blockproblem[b]->blockscip);
2077 blockvars = SCIPgetOrigVars(blockproblem[b]->blockscip);
2078 nblockvars = blockproblem[b]->nblockvars;
2079
2080 SCIP_CALL( SCIPallocBufferArray(scip, &blocksolvals, nblockvars) );
2081 SCIP_CALL( SCIPgetSolVals(blockproblem[b]->blockscip, blocksol, nblockvars, blockvars, blocksolvals) );
2082
2083 for( c = 0; c < nblockvars; c++ )
2084 {
2085 SCIP_VAR* origvar;
2086
2087 origvar = SCIPfindVar(scip, SCIPvarGetName(blockvars[c]));
2088 SCIP_CALL( SCIPsetSolVal(scip, newsol, origvar, blocksolvals[c]) );
2089 }
2090
2091 SCIPfreeBufferArray(scip, &blocksolvals);
2092 }
2093
2094 /* if reoptimization is activated, fix partition and reoptimize with original objective function */
2095 if( heurdata->reoptimize )
2096 {
2097 SCIP_SOL* improvedsol = NULL;
2098 SCIP_CALL( reoptimize(scip, heur, newsol, blockproblem, heurdata->nblocks, heurdata->reoptlimits, &improvedsol, &success) );
2099 assert(improvedsol != NULL || success == FALSE);
2100
2101 if( success )
2102 {
2103 SCIP_CALL( SCIPtrySolFree(scip, &improvedsol, TRUE, FALSE, TRUE, TRUE, TRUE, &success) );
2104 if( success )
2105 {
2106 SCIPdebugMsg(scip, "Reoptimizing solution successful\n");
2107 *result = SCIP_FOUNDSOL;
2108 }
2109 }
2110 }
2111
2112 /* if reoptimization is turned off or reoptimization found no solution, try initial solution */
2113 if( *result != SCIP_FOUNDSOL )
2114 {
2115 SCIPdebugMsg(scip, "Solution has value: %.2f\n", SCIPgetSolOrigObj(scip, newsol));
2116 SCIP_CALL( SCIPtrySolFree(scip, &newsol, TRUE, FALSE, TRUE, TRUE, TRUE, &success) );
2117 if( success )
2118 {
2119 SCIPdebugMsg(scip, "Solution copy successful\n");
2120 *result = SCIP_FOUNDSOL;
2121 }
2122 }
2123 else
2124 {
2125 SCIP_CALL( SCIPfreeSol(scip, &newsol) );
2126 }
2127
2128 goto TERMINATE;
2129 }
2130 /* current partition is not feasible:
2131 * - update partition
2132 * - update penalty parameters lambda
2133 * - reuse last solution
2134 */
2135 else
2136 {
2137 int* nviolatedblocksrhs; /* number of blocks which violate rhs for each linking constraint */
2138 int* nviolatedblockslhs; /* number of blocks which violate lhs for each linking constraint */
2139 SCIP_Bool update = FALSE;
2140
2141 SCIP_CALL( SCIPallocBufferArray(scip, &nviolatedblocksrhs, heurdata->nlinking) );
2142 SCIP_CALL( SCIPallocBufferArray(scip, &nviolatedblockslhs, heurdata->nlinking) );
2143 BMSclearMemoryArray(nviolatedblocksrhs, heurdata->nlinking);
2144 BMSclearMemoryArray(nviolatedblockslhs, heurdata->nlinking);
2145
2146 SCIPdebugMsg(scip, "update partition\n");
2147 SCIP_CALL( updatePartition(scip, linkings, blockproblem, &nviolatedblocksrhs, &nviolatedblockslhs,
2148 heurdata->nlinking, nblocks, k, &update) );
2149
2150 /* terminate if partition is not updated */
2151 if( !update )
2152 {
2153 SCIPfreeBufferArray(scip, &nviolatedblockslhs);
2154 SCIPfreeBufferArray(scip, &nviolatedblocksrhs);
2155 goto TERMINATE;
2156 }
2157
2158 /* in order to change objective function */
2159 for( b = 0; b < heurdata->nblocks; b++ )
2160 {
2161 SCIP_CALL( SCIPfreeTransform(blockproblem[b]->blockscip) );
2162 }
2163
2164 SCIPdebugMsg(scip, "update lambda\n");
2165 SCIP_CALL( updateLambda(scip, heurdata, linkings, blockproblem, nviolatedblocksrhs, nviolatedblockslhs, heurdata->nlinking) );
2166
2167 if( heurdata->reuse )
2168 {
2169 /* reuse old solution in each block if available */
2170 SCIP_CALL( reuseSolution(linkings, blockproblem, nblocks) );
2171 }
2172
2173 SCIPfreeBufferArray(scip, &nviolatedblockslhs);
2174 SCIPfreeBufferArray(scip, &nviolatedblocksrhs);
2175 }
2176 }
2177 SCIPdebugMsg(scip, "maximum number of iterations reached\n");
2178
2179 /* ------------------------------------------------------------------------ */
2180 /* free memory */
2181 TERMINATE:
2182 if( linkings != NULL )
2183 {
2184 for( c = heurdata->nlinking - 1; c >= 0; c-- )
2185 {
2186 if( linkings[c]->currentlhs != NULL )
2187 SCIPfreeBufferArray(scip, &(linkings[c])->currentlhs);
2188
2189 if( linkings[c]->currentrhs != NULL )
2190 SCIPfreeBufferArray(scip, &(linkings[c])->currentrhs);
2191 }
2192
2193 for( c = heurdata->nlinking - 1; c >= 0; c-- )
2194 {
2195 linkings[c]->linkingcons = NULL;
2196 SCIPfreeBufferArray(scip, &(linkings[c])->maxactivity);
2197 SCIPfreeBufferArray(scip, &(linkings[c])->minactivity);
2198 SCIPfreeBufferArray(scip, &(linkings[c])->blocknumbers);
2199 SCIPfreeBufferArray(scip, &(linkings[c])->slacks);
2200 SCIPfreeBufferArray(scip, &(linkings[c])->blockconss);
2201 SCIPfreeBlockMemory(scip, &(linkings[c])); /*lint !e866*/
2202 }
2203 SCIPfreeBufferArray(scip, &linkings);
2204 }
2205
2206 if( blockproblem != NULL )
2207 {
2208 for( b = nblocks - 1; b >= 0; b-- )
2209 {
2210 SCIPfreeBufferArray(scip, &(blockproblem[b])->origobj);
2211 SCIPfreeBufferArray(scip, &(blockproblem[b])->slackvars);
2212 SCIPfreeBufferArray(scip, &(blockproblem[b])->linkingindices);
2213 SCIPfreeBufferArray(scip, &(blockproblem[b])->linkingconss);
2214 SCIP_CALL( SCIPfree(&blockproblem[b]->blockscip) );
2215 SCIPfreeBlockMemory(scip, &(blockproblem[b])); /*lint !e866*/
2216 }
2217 SCIPfreeBufferArray(scip, &blockproblem);
2218 }
2219
2220 if( assigneddecomp != NULL )
2221 SCIPdecompFree(&assigneddecomp, SCIPblkmem(scip));
2222
2223 if( sortedconslabels != NULL )
2224 SCIPfreeBufferArray(scip, &sortedconslabels);
2225
2226 if( sortedvarlabels != NULL )
2227 SCIPfreeBufferArray(scip, &sortedvarlabels);
2228
2229 if( sortedconss != NULL )
2230 SCIPfreeBufferArray(scip, &sortedconss);
2231
2232 if( sortedvars != NULL )
2233 SCIPfreeBufferArray(scip, &sortedvars);
2234
2235 SCIPdebugMsg(scip, "Leave DPS heuristic\n");
2236
2237 return SCIP_OKAY;
2238 }
2239
2240
2241 /*
2242 * primal heuristic specific interface methods
2243 */
2244
2245 /** creates the dps primal heuristic and includes it in SCIP */
2246 SCIP_RETCODE SCIPincludeHeurDps(
2247 SCIP* scip /**< SCIP data structure */
2248 )
2249 {
2250 SCIP_HEURDATA* heurdata;
2251 SCIP_HEUR* heur;
2252
2253 /* create dps primal heuristic data */
2254 SCIP_CALL( SCIPallocBlockMemory(scip, &heurdata) );
2255
2256 heur = NULL;
2257
2258 /* include primal heuristic */
2259 SCIP_CALL( SCIPincludeHeurBasic(scip, &heur,
2260 HEUR_NAME, HEUR_DESC, HEUR_DISPCHAR, HEUR_PRIORITY, HEUR_FREQ, HEUR_FREQOFS,
2261 HEUR_MAXDEPTH, HEUR_TIMING, HEUR_USESSUBSCIP, heurExecDps, heurdata) );
2262
2263 assert(heur != NULL);
2264
2265 /* set non fundamental callbacks via setter functions */
2266 SCIP_CALL( SCIPsetHeurCopy(scip, heur, heurCopyDps) );
2267 SCIP_CALL( SCIPsetHeurFree(scip, heur, heurFreeDps) );
2268
2269 /* add dps primal heuristic parameters */
2270 SCIP_CALL( SCIPaddIntParam(scip, "heuristics/" HEUR_NAME "/maxiterations",
2271 "maximal number of iterations", &heurdata->maxit, FALSE, DEFAULT_MAXIT, 1, INT_MAX, NULL, NULL) );
2272
2273 SCIP_CALL( SCIPaddRealParam(scip, "heuristics/" HEUR_NAME "/maxlinkscore",
2274 "maximal linking score of used decomposition (equivalent to percentage of linking constraints)",
2275 &heurdata->maxlinkscore, FALSE, 1.0, 0.0, 1.0, NULL, NULL) );
2276
2277 SCIP_CALL( SCIPaddRealParam(scip, "heuristics/" HEUR_NAME "/penalty",
2278 "multiplier for absolute increase of penalty parameters (0: no increase)",
2279 &heurdata->penalty, FALSE, DEFAULT_PENALTY, 0.0, SCIP_REAL_MAX, NULL, NULL) );
2280
2281 SCIP_CALL( SCIPaddBoolParam(scip, "heuristics/" HEUR_NAME "/reoptimize",
2282 "should the problem get reoptimized with the original objective function?", &heurdata->reoptimize, FALSE, FALSE, NULL, NULL) );
2283
2284 SCIP_CALL( SCIPaddBoolParam(scip, "heuristics/" HEUR_NAME "/reuse",
2285 "should solutions get reused in subproblems?", &heurdata->reuse, FALSE, FALSE, NULL, NULL) );
2286
2287 SCIP_CALL( SCIPaddBoolParam(scip, "heuristics/" HEUR_NAME "/reoptlimits",
2288 "should strict limits for reoptimization be set?", &heurdata->reoptlimits, FALSE, TRUE, NULL, NULL) );
2289
2290 SCIP_CALL( SCIPaddIntParam(scip, "heuristics/" HEUR_NAME "/timing",
2291 "should the heuristic run before or after the processing of the node? (0: before, 1: after, 2: both)",
2292 &heurdata->timing, FALSE, 0, 0, 2, NULL, NULL) );
2293
2294 return SCIP_OKAY;
2295 }
2296