1 /* * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * */
2 /* */
3 /* This file is part of the program and library */
4 /* SCIP --- Solving Constraint Integer Programs */
5 /* */
6 /* Copyright (C) 2002-2021 Konrad-Zuse-Zentrum */
7 /* fuer Informationstechnik Berlin */
8 /* */
9 /* SCIP is distributed under the terms of the ZIB Academic License. */
10 /* */
11 /* You should have received a copy of the ZIB Academic License */
12 /* along with SCIP; see the file COPYING. If not visit scipopt.org. */
13 /* */
14 /* * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * */
15
16 /**@file cons_orbitope.c
17 * @ingroup DEFPLUGINS_CONS
18 * @brief constraint handler for (partitioning/packing/full) orbitope constraints w.r.t. the full symmetric group
19 * @author Timo Berthold
20 * @author Marc Pfetsch
21 * @author Christopher Hojny
22 *
23 * The type of constraints of this constraint handler is described in cons_orbitope.h.
24 *
25 * The details of the method implemented here are described in the following papers.
26 *
27 * Packing and Partitioning Orbitopes@n
28 * Volker Kaibel and Marc E. Pfetsch,@n
29 * Math. Program. 114, No. 1, 1-36 (2008)
30 *
31 * Among other things, this paper describes so-called shifted column inequalities of the following
32 * form \f$x(S) \leq x(B)\f$, where \f$S\f$ is a so-called shifted column and \f$B\f$ is a so-called
33 * bar. These inequalities can be used to handle symmetry and they are separated in this constraint
34 * handler. We use the linear time separation algorithm of the paper.@par
35 *
36 * Orbitopal Fixing@n
37 * Volker Kaibel, Matthias Peinhardt, and Marc E. Pfetsch,@n
38 * Discrete Optimization 8, No. 4, 595-610 (2011)
39 * (A preliminary version appears in Proc. IPCO 2007.)
40 *
41 * In this paper a linear time propagation algorithm is described, a variant of which is implemented
42 * here. The implemented variant does not run in linear time, but is very fast in practice.
43 *
44 * <table>
45 * <caption>translation table</caption>
46 * <tr><td>here</td><td>paper</td></tr>
47 * <tr><td></td><td></td></tr>
48 * <tr><td>nspcons </td><td>p </td></tr>
49 * <tr><td>nblocks </td><td>q </td></tr>
50 * <tr><td>vars </td><td>x </td></tr>
51 * <tr><td>vals </td><td>A^\\star</td></tr>
52 * <tr><td>weights </td><td>\\omega </td></tr>
53 * <tr><td>cases </td><td>\\tau </td></tr>
54 * <tr><td>fixtriangle </td><td>-- </td></tr>
55 * <tr><td>resolveprop </td><td>-- </td></tr>
56 * <tr><td>firstnonzeros</td><td>\\mu </td></tr>
57 * <tr><td>lastones </td><td>\\alpha </td></tr>
58 * <tr><td>frontiersteps</td><td>\\Gamma </td></tr>
59 * </table>
60 *
61 * Orbitopal fixing for the full (sub-)orbitope and application to the Unit Commitment Problem@n
62 * Pascale Bendotti, Pierre Fouilhoux, and Cecile Rottner,@n
63 * Optimization Online: http://www.optimization-online.org/DB_HTML/2017/10/6301.html
64 *
65 * Two linear time propagation algorithms for full orbitopes are described in this paper, a static
66 * version and a dynamic one. While the static version uses a fixed variable order, the dynamic
67 * version determines the variable order during the solving process via branching descisions.
68 * We implemented the static version as well as a modified version of the dynamic one. The reason
69 * for the latter is to simplify the compatibility with full orbitope cutting planes.
70 *
71 * Note, however, that the dynamic version may lead to conflicts if orbitopes are copied to subSCIPs.
72 * Since the dynamic version is based on branching decisions, which may be different in main SCIP
73 * and subSCIPs, orbitopes using the dynamic algorithm are not allowed to be copied. However, as a
74 * user might use orbitopes to enforce a certain variable ordering in a solution, we distinguish
75 * whether an orbitope is a model constraint or not. If it is a model constraint, we assume that
76 * a variable order has already been fixed and disable the dynamic algorithm. In this case, orbitope
77 * constraints are copied to subSCIPs. If it is not a model constraint, the orbitope was added to
78 * handle symmetries but not to enforce a solution to have a certain structure. In this case, the
79 * dynamic algorithm can be used and we do not copy orbitope constraints to subSCIPs.
80 *
81 * Polytopes associated with symmetry handling@n
82 * Christopher Hojny and Marc E. Pfetsch,@n
83 * Math. Program. (2018)
84 *
85 * In this paper, a linear time separation algorithm for orbisacks (full orbitopes with two columnes)
86 * is described. We use this algorithm for every pair of adjacent columns within the orbitope as well
87 * as a version that is adapted to the reordering based on the dynamic full orbitope propagation
88 * algorithm to ensure validity of binary points via cutting planes.
89 */
90
91 /*---+----1----+----2----+----3----+----4----+----5----+----6----+----7----+----8----+----9----+----0----+----1----+----2*/
92
93 #include "blockmemshell/memory.h"
94 #include "scip/cons_orbisack.h"
95 #include "scip/cons_orbitope.h"
96 #include "scip/cons_setppc.h"
97 #include "scip/pub_cons.h"
98 #include "scip/pub_message.h"
99 #include "scip/pub_var.h"
100 #include "scip/scip.h"
101 #include "scip/scip_branch.h"
102 #include "scip/scip_conflict.h"
103 #include "scip/scip_cons.h"
104 #include "scip/scip_copy.h"
105 #include "scip/scip_cut.h"
106 #include "scip/scip_general.h"
107 #include "scip/scip_lp.h"
108 #include "scip/scip_mem.h"
109 #include "scip/scip_message.h"
110 #include "scip/scip_numerics.h"
111 #include "scip/scip_param.h"
112 #include "scip/scip_prob.h"
113 #include "scip/scip_probing.h"
114 #include "scip/scip_sol.h"
115 #include "scip/scip_var.h"
116 #include <ctype.h>
117 #include <string.h>
118
119 /* constraint handler properties */
120 #define CONSHDLR_NAME "orbitope"
121 #define CONSHDLR_DESC "symmetry breaking constraint handler relying on (partitioning/packing) orbitopes"
122 #define CONSHDLR_SEPAPRIORITY +40100 /**< priority of the constraint handler for separation */
123 #define CONSHDLR_ENFOPRIORITY -1005200 /**< priority of the constraint handler for constraint enforcing */
124 #define CONSHDLR_CHECKPRIORITY -1005200 /**< priority of the constraint handler for checking feasibility */
125 #define CONSHDLR_SEPAFREQ -1 /**< frequency for separating cuts; zero means to separate only in the root node */
126 #define CONSHDLR_PROPFREQ 1 /**< frequency for propagating domains; zero means only preprocessing propagation */
127 #define CONSHDLR_EAGERFREQ -1 /**< frequency for using all instead of only the useful constraints in separation,
128 * propagation and enforcement, -1 for no eager evaluations, 0 for first only */
129 #define CONSHDLR_MAXPREROUNDS -1 /**< maximal number of presolving rounds the constraint handler participates in (-1: no limit) */
130 #define CONSHDLR_DELAYSEPA FALSE /**< should separation method be delayed, if other separators found cuts? */
131 #define CONSHDLR_DELAYPROP FALSE /**< should propagation method be delayed, if other propagators found reductions? */
132 #define CONSHDLR_NEEDSCONS TRUE /**< should the constraint handler be skipped, if no constraints are available? */
133
134 #define CONSHDLR_PROP_TIMING SCIP_PROPTIMING_BEFORELP /**< propagation timing mask of the constraint handler */
135 #define CONSHDLR_PRESOLTIMING SCIP_PRESOLTIMING_MEDIUM /**< presolving timing of the constraint handler (fast, medium, or exhaustive) */
136
137 #define DEFAULT_PPORBITOPE TRUE /**< whether we check if full orbitopes can be strengthened to packing/partitioning orbitopes */
138 #define DEFAULT_SEPAFULLORBITOPE FALSE /**< whether we separate inequalities for full orbitopes */
139 #define DEFAULT_USEDYNAMICPROP TRUE /**< whether we use a dynamic version of the propagation routine */
140 #define DEFAULT_FORCECONSCOPY FALSE /**< whether orbitope constraints should be forced to be copied to sub SCIPs */
141
142 /*
143 * Data structures
144 */
145
146 /** constraint handler data */
147 struct SCIP_ConshdlrData
148 {
149 SCIP_Bool checkpporbitope; /**< whether we allow upgrading to packing/partitioning orbitopes */
150 SCIP_Bool sepafullorbitope; /**< whether we separate inequalities for full orbitopes orbitopes */
151 SCIP_Bool usedynamicprop; /**< whether we use a dynamic version of the propagation routine */
152 SCIP_Bool forceconscopy; /**< whether orbitope constraints should be forced to be copied to sub SCIPs */
153 };
154
155 /** constraint data for orbitope constraints */
156 struct SCIP_ConsData
157 {
158 SCIP_VAR*** vars; /**< matrix of variables on which the symmetry acts */
159 SCIP_VAR** tmpvars; /**< temporary storage for variables */
160 SCIP_HASHMAP* rowindexmap; /**< map of variables to row index in orbitope matrix */
161 SCIP_Real** vals; /**< LP-solution for those variables */
162 SCIP_Real* tmpvals; /**< temporary storage for values */
163 SCIP_Real** weights; /**< SC weight table */
164 int** cases; /**< indicator of the SC cases */
165 int nspcons; /**< number of set partitioning/packing constraints <=> p */
166 int nblocks; /**< number of symmetric variable blocks <=> q */
167 SCIP_ORBITOPETYPE orbitopetype; /**< type of orbitope constraint */
168 SCIP_Bool resolveprop; /**< should propagation be resolved? */
169 SCIP_Bool istrianglefixed; /**< has the upper right triangle already globally been fixed to zero? */
170 int* roworder; /**< order of orbitope rows if dynamic propagation for full orbitopes
171 * is used. */
172 SCIP_Bool* rowused; /**< whether a row has been considered in roworder */
173 int nrowsused; /**< number of rows that have already been considered in roworder */
174 SCIP_Bool ismodelcons; /**< whether the orbitope is a model constraint */
175 };
176
177
178 /*
179 * Local methods
180 */
181
182 /** frees an orbitope constraint data */
183 static
consdataFree(SCIP * scip,SCIP_CONSDATA ** consdata,SCIP_Bool usedynamicprop)184 SCIP_RETCODE consdataFree(
185 SCIP* scip, /**< SCIP data structure */
186 SCIP_CONSDATA** consdata, /**< pointer to orbitope constraint data */
187 SCIP_Bool usedynamicprop /**< whether we use a dynamic version of the propagation routine */
188 )
189 {
190 int i;
191 int p;
192 int q;
193
194 assert( consdata != NULL );
195 assert( *consdata != NULL );
196
197 if ( usedynamicprop && (*consdata)->rowindexmap != NULL )
198 {
199 SCIPhashmapFree(&((*consdata)->rowindexmap));
200 }
201
202 p = (*consdata)->nspcons;
203 q = (*consdata)->nblocks;
204 for (i = 0; i < p; ++i)
205 {
206 SCIPfreeBlockMemoryArrayNull(scip, &((*consdata)->cases[i]), q); /*lint !e866*/
207 SCIPfreeBlockMemoryArrayNull(scip, &((*consdata)->vars[i]), q); /*lint !e866*/
208 SCIPfreeBlockMemoryArrayNull(scip, &((*consdata)->weights[i]), q); /*lint !e866*/
209 SCIPfreeBlockMemoryArrayNull(scip, &((*consdata)->vals[i]), q); /*lint !e866*/
210 }
211
212 if ( usedynamicprop )
213 {
214 SCIPfreeBlockMemoryArrayNull(scip, &((*consdata)->rowused), p);
215 }
216 SCIPfreeBlockMemoryArrayNull(scip, &((*consdata)->roworder), p);
217 SCIPfreeBlockMemoryArrayNull(scip, &((*consdata)->cases), p);
218 SCIPfreeBlockMemoryArrayNull(scip, &((*consdata)->vars), p);
219 SCIPfreeBlockMemoryArrayNull(scip, &((*consdata)->weights), p);
220 SCIPfreeBlockMemoryArrayNull(scip, &((*consdata)->vals), p);
221
222 SCIPfreeBlockMemoryArrayNull(scip, &((*consdata)->tmpvals), p + q);
223 SCIPfreeBlockMemoryArrayNull(scip, &((*consdata)->tmpvars), p + q);
224
225 SCIPfreeBlockMemory(scip, consdata);
226
227 return SCIP_OKAY;
228 }
229
230
231 /** creates orbitope constraint data */
232 static
consdataCreate(SCIP * scip,SCIP_CONSDATA ** consdata,SCIP_VAR *** vars,int nspcons,int nblocks,SCIP_ORBITOPETYPE orbitopetype,SCIP_Bool resolveprop,SCIP_Bool usedynamicprop,SCIP_Bool ismodelcons)233 SCIP_RETCODE consdataCreate(
234 SCIP* scip, /**< SCIP data structure */
235 SCIP_CONSDATA** consdata, /**< pointer to store constraint data */
236 SCIP_VAR*** vars, /**< variables array, must have size nspcons x nblocks */
237 int nspcons, /**< number of set partitioning (packing) constraints <=> p */
238 int nblocks, /**< number of symmetric variable blocks <=> q */
239 SCIP_ORBITOPETYPE orbitopetype, /**< type of orbitope constraint */
240 SCIP_Bool resolveprop, /**< should propagation be resolved? */
241 SCIP_Bool usedynamicprop, /**< whether we use a dynamic version of the propagation routine */
242 SCIP_Bool ismodelcons /**< whether the orbitope is a model constraint */
243 )
244 {
245 int i;
246 int j;
247
248 assert(consdata != NULL);
249
250 SCIP_CALL( SCIPallocBlockMemory(scip, consdata) );
251
252 SCIP_CALL( SCIPallocBlockMemoryArray(scip, &(*consdata)->vals, nspcons) );
253 SCIP_CALL( SCIPallocBlockMemoryArray(scip, &(*consdata)->weights, nspcons) );
254 SCIP_CALL( SCIPallocBlockMemoryArray(scip, &(*consdata)->vars, nspcons) );
255 SCIP_CALL( SCIPallocBlockMemoryArray(scip, &(*consdata)->cases, nspcons) );
256 SCIP_CALL( SCIPallocBlockMemoryArray(scip, &(*consdata)->roworder, nspcons) );
257
258 if ( usedynamicprop )
259 {
260 SCIP_CALL( SCIPhashmapCreate(&(*consdata)->rowindexmap, SCIPblkmem(scip), nspcons) );
261 SCIP_CALL( SCIPallocBlockMemoryArray(scip, &(*consdata)->rowused, nspcons) );
262 }
263
264 for (i = 0; i < nspcons; ++i)
265 {
266 SCIP_CALL( SCIPallocBlockMemoryArray(scip, &(*consdata)->vals[i], nblocks) ); /*lint !e866*/
267 SCIP_CALL( SCIPallocBlockMemoryArray(scip, &(*consdata)->weights[i], nblocks) ); /*lint !e866*/
268 SCIP_CALL( SCIPduplicateBlockMemoryArray(scip, &(*consdata)->vars[i], vars[i], nblocks) ); /*lint !e866*/
269 SCIP_CALL( SCIPallocBlockMemoryArray(scip, &(*consdata)->cases[i], nblocks) ); /*lint !e866*/
270 (*consdata)->roworder[i] = i;
271
272 if ( usedynamicprop )
273 {
274 (*consdata)->rowused[i] = FALSE;
275 }
276 }
277 (*consdata)->nrowsused = 0;
278
279 (*consdata)->tmpvals = NULL;
280 (*consdata)->tmpvars = NULL;
281 (*consdata)->nspcons = nspcons;
282 (*consdata)->nblocks = nblocks;
283 (*consdata)->orbitopetype = orbitopetype;
284 (*consdata)->resolveprop = resolveprop;
285 (*consdata)->istrianglefixed = FALSE;
286 (*consdata)->ismodelcons = ismodelcons;
287
288 /* get transformed variables, if we are in the transformed problem */
289 if ( SCIPisTransformed(scip) )
290 {
291 SCIP_CALL( SCIPallocBlockMemoryArray(scip, &(*consdata)->tmpvals, nspcons + nblocks) );
292 SCIP_CALL( SCIPallocBlockMemoryArray(scip, &(*consdata)->tmpvars, nspcons + nblocks) );
293
294 for (i = 0; i < nspcons; ++i)
295 {
296 /* make sure that no variable gets multiaggregated (cannot be handled by cons_orbitope, since one cannot easily
297 * eliminate single variables from an orbitope constraint).
298 */
299 for (j = 0; j < nblocks; ++j)
300 {
301 SCIP_CALL( SCIPgetTransformedVar(scip, (*consdata)->vars[i][j], &(*consdata)->vars[i][j]) );
302 SCIP_CALL( SCIPmarkDoNotMultaggrVar(scip, (*consdata)->vars[i][j]) );
303 if ( usedynamicprop )
304 {
305 SCIP_CALL( SCIPhashmapInsert((*consdata)->rowindexmap, (*consdata)->vars[i][j], (void*) (size_t) i) );
306 }
307 }
308 }
309 }
310
311 return SCIP_OKAY;
312 }
313
314
315 /** strengthen full orbitopes to packing/partitioning orbitopes if possible */
316 static
strengthenOrbitopeConstraint(SCIP * scip,SCIP_VAR *** vars,int * nrows,int ncols,SCIP_ORBITOPETYPE * type)317 SCIP_RETCODE strengthenOrbitopeConstraint(
318 SCIP* scip, /**< SCIP data structure */
319 SCIP_VAR*** vars, /**< variable matrix of orbitope constraint */
320 int* nrows, /**< pointer to number of rows of variable matrix */
321 int ncols, /**< number of columns of variable matrix */
322 SCIP_ORBITOPETYPE* type /**< pointer to store type of orbitope constraint after strengthening */
323 )
324 {
325 SCIP_CONSHDLR* setppcconshdlr;
326 SCIP_CONS** setppcconss;
327 int nsetppcconss;
328 int* covered;
329 int nprobvars;
330 int* rowidxvar;
331 int* rowcoveragesetppc;
332 int* rowsinsetppc;
333 int ncovered;
334 int ncoveredpart;
335 int i;
336 int j;
337 int c;
338
339 assert( scip != NULL );
340 assert( vars != NULL );
341 assert( vars != NULL );
342 assert( *nrows > 0 );
343 assert( ncols > 0 );
344 assert( type != NULL );
345
346 *type = SCIP_ORBITOPETYPE_FULL;
347
348 setppcconshdlr = SCIPfindConshdlr(scip, "setppc");
349 if ( setppcconshdlr == NULL )
350 return SCIP_OKAY;
351
352 setppcconss = SCIPconshdlrGetConss(setppcconshdlr);
353 nsetppcconss = SCIPconshdlrGetNConss(setppcconshdlr);
354
355 if ( nsetppcconss == 0 )
356 return SCIP_OKAY;
357 assert( setppcconss != NULL );
358
359 /* whether a row is contained in packing/partitioning constraint */
360 SCIP_CALL( SCIPallocClearBufferArray(scip, &covered, *nrows) );
361 ncovered = 0;
362 ncoveredpart = 0;
363
364 /* array storing index of orbitope row a variable is contained in */
365 nprobvars = SCIPgetNTotalVars(scip);
366
367 SCIP_CALL( SCIPallocBufferArray(scip, &rowidxvar, nprobvars) );
368
369 for (i = 0; i < nprobvars; ++i)
370 rowidxvar[i] = -1;
371
372 for (i = 0; i < *nrows; ++i)
373 {
374 for (j = 0; j < ncols; ++j)
375 {
376 assert( SCIPvarGetIndex(vars[i][j]) >= 0 && SCIPvarGetIndex(vars[i][j]) < nprobvars );
377 rowidxvar[SCIPvarGetIndex(vars[i][j])] = i;
378 }
379 }
380
381 /* storage for number of vars per row that are contained in current setppc cons and
382 * labels of rows intersecting with current setppc cons
383 */
384 SCIP_CALL( SCIPallocClearBufferArray(scip, &rowcoveragesetppc, *nrows) );
385 SCIP_CALL( SCIPallocClearBufferArray(scip, &rowsinsetppc, *nrows) );
386
387 /* iterate over set packing and partitioning constraints and check whether the constraint's
388 * support is a row r of the orbitope (covered[r] = 2) or contains row r (covered[r] = 1)
389 */
390 for (c = 0; c < nsetppcconss && ncoveredpart < ncols; ++c)
391 {
392 int nsetppcvars;
393 SCIP_VAR** setppcvars;
394 SCIP_VAR* var;
395 int nrowintersect = 0;
396 int nvarsinorbitope;
397
398 /* skip covering constraints */
399 if ( SCIPgetTypeSetppc(scip, setppcconss[c]) == SCIP_SETPPCTYPE_COVERING )
400 continue;
401
402 /* get set packing/partitioning variables */
403 nsetppcvars = SCIPgetNVarsSetppc(scip, setppcconss[c]);
404
405 /* constraint does not contain enough variables */
406 if ( nsetppcvars < ncols )
407 continue;
408
409 setppcvars = SCIPgetVarsSetppc(scip, setppcconss[c]);
410 assert( setppcvars != NULL );
411
412 /* upper bound on variables potentially contained in orbitope */
413 nvarsinorbitope = nsetppcvars;
414
415 /* for each setppc var, check whether it appears in a row of the orbitope and store
416 * for each row the number of such variables; can be terminated early, if less than
417 * ncols variables are contained in the orbitope
418 */
419 for (i = 0; i < nsetppcvars && nvarsinorbitope >= ncols; ++i)
420 {
421 int idx;
422 int rowidx;
423
424 var = setppcvars[i];
425 idx = SCIPvarGetIndex(var);
426
427 assert( idx < nprobvars );
428 assert( idx >= 0 );
429
430 rowidx = rowidxvar[idx];
431
432 /* skip variables not contained in the orbitope */
433 if ( rowidx < 0 )
434 {
435 --nvarsinorbitope;
436 continue;
437 }
438
439 /* skip variables corresponding to already treated rows */
440 if ( covered[rowidx] == 2 || (covered[rowidx] == 1 && (nsetppcvars > ncols || nrowintersect > 1)) )
441 {
442 --nvarsinorbitope;
443 continue;
444 }
445
446 /* store information which rows intersect the setppc cons's support */
447 if ( rowcoveragesetppc[rowidx] == 0 )
448 rowsinsetppc[nrowintersect++] = rowidx;
449 rowcoveragesetppc[rowidx] += 1;
450 }
451
452 /* store whether rows coincide with set partitioning cons's support or whether
453 * row is covered by a set packing/partitioning cons's support
454 */
455 if ( SCIPgetTypeSetppc(scip, setppcconss[c]) == SCIP_SETPPCTYPE_PARTITIONING
456 && nrowintersect == 1 && rowcoveragesetppc[rowsinsetppc[0]] == ncols && nsetppcvars == ncols )
457 {
458 if ( covered[rowsinsetppc[0]] == 1 )
459 --ncovered;
460 covered[rowsinsetppc[0]] = 2;
461 ++ncoveredpart;
462 ++ncovered;
463 }
464 else
465 {
466 for (i = 0; i < nrowintersect; ++i)
467 {
468 if ( covered[rowsinsetppc[i]] == 0 && rowcoveragesetppc[rowsinsetppc[i]] >= ncols )
469 {
470 covered[rowsinsetppc[i]] = 1;
471 ++ncovered;
472 }
473 }
474 }
475
476 /* reset data */
477 for (i = 0; i < nrowintersect; ++i)
478 rowcoveragesetppc[rowsinsetppc[i]] = 0;
479 }
480
481 /* check type of orbitope */
482 if ( ncovered == *nrows )
483 {
484 if ( ncoveredpart == *nrows )
485 *type = SCIP_ORBITOPETYPE_PARTITIONING;
486 else
487 *type = SCIP_ORBITOPETYPE_PACKING;
488 }
489 /* If only some rows are contained in set packing/partitioning constraints, it may still be worth it
490 * to exploit the packing/partitioning structure on these rows, because packing/partitioning orbitopes
491 * or more restrictive than full orbitopes. If at least three rows have this property, we discard
492 * all rows not contained in set packing/partitioning constraints and add the smaller sub packing orbitope.
493 */
494 else if ( ncovered >= 3 )
495 {
496 int r = *nrows - 1;
497 while ( r >= 0 )
498 {
499 if ( ! covered[r] )
500 {
501 for (i = r; i < *nrows - 1; ++i)
502 {
503 SCIP_VAR** row = vars[i];
504 vars[i] = vars[i+1];
505 vars[i+1] = row;
506 }
507 *nrows -= 1;
508 }
509 --r;
510 }
511 *type = SCIP_ORBITOPETYPE_PACKING;
512 }
513
514 SCIPfreeBufferArray(scip, &rowsinsetppc);
515 SCIPfreeBufferArray(scip, &rowcoveragesetppc);
516 SCIPfreeBufferArray(scip, &rowidxvar);
517 SCIPfreeBufferArray(scip, &covered);
518
519 return SCIP_OKAY;
520 }
521
522 #ifdef PRINT_MATRIX
523 /** debug method, prints variable matrix */
524 static
printMatrix(SCIP * scip,SCIP_CONSDATA * consdata)525 void printMatrix(
526 SCIP* scip, /**< SCIP data structure */
527 SCIP_CONSDATA* consdata /**< the constraint data */
528 )
529 {
530 int i;
531 int j;
532
533 assert( consdata != NULL );
534 assert( consdata->nspcons > 0 );
535 assert( consdata->nblocks > 0 );
536 assert( consdata->vars != NULL );
537
538 for (j = 0; j < consdata->nblocks; ++j)
539 SCIPinfoMessage(scip, NULL, "-");
540
541 SCIPinfoMessage(scip, NULL, "\n");
542 for (i = 0; i < consdata->nspcons; ++i)
543 {
544 for (j = 0; j < consdata->nblocks; ++j)
545 {
546 if ( SCIPvarGetUbLocal(consdata->vars[i][j]) - SCIPvarGetLbLocal(consdata->vars[i][j]) < 0.5 )
547 SCIPinfoMessage(scip, NULL, "%1.0f", REALABS(SCIPvarGetUbLocal(consdata->vars[i][j])));
548 else
549 SCIPinfoMessage(scip, NULL, " ");
550 }
551 SCIPinfoMessage(scip, NULL, "|\n");
552 }
553 for (j = 0; j < consdata->nblocks; ++j)
554 SCIPinfoMessage(scip, NULL, "-");
555 SCIPinfoMessage(scip, NULL, "\n");
556 }
557 #endif
558
559
560 #ifdef SHOW_SCI
561 /** Print SCI in nice form for debugging */
562 static
printSCI(SCIP * scip,int p,int q,int ** cases,int i,int j)563 SCIP_RETCODE printSCI(
564 SCIP* scip, /**< SCIP pointer */
565 int p, /**< number of rows */
566 int q, /**< number of columns */
567 int** cases, /**< SCI dynamic programming table */
568 int i, /**< row position of bar */
569 int j /**< column position of bar */
570 )
571 {
572 int k;
573 int l;
574 int** M;
575 int p1;
576 int p2;
577
578 SCIP_CALL( SCIPallocBufferArray(scip, &M, p) );
579 for (k = 0; k < p; ++k)
580 {
581 SCIP_CALL( SCIPallocBufferArray(scip, &M[k], q) ); /*lint !e866*/
582 for (l = 0; l < q; ++l)
583 M[k][l] = 0;
584 }
585
586 /* first add bar */
587 for (l = j; l < q; ++l)
588 {
589 assert( M[i][l] == 0 );
590 M[i][l] = 1;
591 }
592
593 /* then add shifted column */
594 p1 = i-1;
595 p2 = j-1;
596 do
597 {
598 assert( cases[p1][p2] != -1 );
599 assert( p1 >= 0 && p1 < i );
600 assert( p2 >= 0 && p2 < j );
601
602 /* if case 1 */
603 if ( cases[p1][p2] == 1 )
604 --p2; /* decrease column */
605 else
606 {
607 /* case 2 or 3: */
608 assert( cases[p1][p2] == 2 || cases[p1][p2] == 3 );
609 assert( M[p1][p2] == 0 );
610 M[p1][p2] = -1;
611 if ( cases[p1][p2] == 3 )
612 break;
613 }
614 --p1; /* decrease row */
615 }
616 while ( p1 >= 0 ); /* should always be true, i.e. the break should end the loop */
617 assert( cases[p1][p2] == 3 );
618
619 /* now output matrix M */
620 for (l = 0; l < q; ++l)
621 SCIPinfoMessage(scip, NULL, "-");
622 SCIPinfoMessage(scip, NULL, "\n");
623
624 for (k = 0; k < p; ++k)
625 {
626 for (l = 0; l < q; ++l)
627 {
628 if ( l > k )
629 SCIPinfoMessage(scip, NULL, "*");
630 else
631 {
632 switch (M[k][l])
633 {
634 case 1:
635 SCIPinfoMessage(scip, NULL, "+");
636 break;
637 case -1:
638 SCIPinfoMessage(scip, NULL, "-");
639 break;
640 case 0:
641 SCIPinfoMessage(scip, NULL, "#");
642 break;
643 default:
644 SCIPerrorMessage("unexpected matrix entry <%d>: should be -1, 0 or +1\n", M[k][l]);
645 SCIPABORT();
646 }
647 }
648 }
649 SCIPinfoMessage(scip, NULL, "\n");
650 }
651
652 for (l = 0; l < q; ++l)
653 SCIPinfoMessage(scip, NULL, "-");
654 SCIPinfoMessage(scip, NULL, "\n");
655
656 for (k = 0; k < p; ++k)
657 SCIPfreeBufferArray(scip, &M[k]);
658 SCIPfreeBufferArray(scip, &M);
659
660 return SCIP_OKAY;
661 }
662 #endif
663
664
665 /** copies the variables values from the solution to the constraint data structure */
666 static
copyValues(SCIP * scip,SCIP_CONSDATA * consdata,SCIP_SOL * sol)667 void copyValues(
668 SCIP* scip, /**< the SCIP data structure */
669 SCIP_CONSDATA* consdata, /**< the constraint data */
670 SCIP_SOL* sol /**< a primal solution or NULL for the current LP optimum */
671 )
672 {
673 int i;
674 int j;
675
676 assert( scip != NULL );
677 assert( consdata != NULL );
678 assert( consdata->nspcons > 0 );
679 assert( consdata->nblocks > 0 );
680 assert( consdata->vars != NULL );
681 assert( consdata->vals != NULL );
682
683 for (i = 0; i < consdata->nspcons; ++i)
684 {
685 for (j = 0; j < consdata->nblocks; ++j)
686 consdata->vals[i][j] = SCIPgetSolVal(scip, sol, consdata->vars[i][j]);
687 }
688 }
689
690
691 /** compute the dynamic programming table for SC
692 *
693 * Build up dynamic programming table in order to find SCs with minimum weight.
694 *
695 * The values of the minimal SCIs are stored in @a weights.
696 * The array @a cases[i][j] stores which of the cases were applied to get @a weights[i][j].
697 * Here, 3 means that we have reached the upper limit.
698 *
699 * We assume that the upper right triangle is fixed to 0. Hence we can perform the computation a
700 * bit more efficient.
701 */
702 static
computeSCTable(SCIP * scip,int nspcons,int nblocks,SCIP_Real ** weights,int ** cases,SCIP_Real ** vals)703 void computeSCTable(
704 SCIP* scip, /**< SCIP pointer */
705 int nspcons, /**< number of set partitioning (packing) constraints <=> p */
706 int nblocks, /**< number of symmetric variable blocks <=> q */
707 SCIP_Real** weights, /**< SC weight table */
708 int** cases, /**< indicator of the SC cases */
709 SCIP_Real** vals /**< current solution */
710 )
711 {
712 SCIP_Real minvalue;
713 int diagsize;
714 int i;
715 int j;
716
717 assert( weights != NULL );
718 assert( cases != NULL );
719 assert( vals != NULL );
720
721 #ifndef NDEBUG
722 /* for debugging */
723 for (i = 0; i < nspcons; ++i)
724 {
725 for (j = 0; j < nblocks; ++j)
726 {
727 if ( i >= j )
728 {
729 weights[i][j] = -1.0;
730 cases[i][j] = -1;
731 }
732 }
733 }
734 #endif
735
736 /* initialize diagonal */
737 minvalue = vals[0][0];
738 weights[0][0] = minvalue;
739 cases[0][0] = 3;
740
741 /* get last row of triangle */
742 diagsize = nblocks;
743 if ( nspcons < nblocks )
744 diagsize = nspcons;
745
746 for (j = 1; j < diagsize; ++j)
747 {
748 /* use LT to move entry as far to the left as possible */
749 if ( SCIPisLT(scip, vals[j][j], minvalue) )
750 {
751 minvalue = vals[j][j];
752 cases[j][j] = 3;
753 }
754 else
755 cases[j][j] = 1;
756 weights[j][j] = minvalue;
757 }
758
759 /* initialize first column */
760 for (i = 1; i < nspcons; ++i)
761 {
762 weights[i][0] = weights[i-1][0] + vals[i][0];
763 cases[i][0] = 2; /* second case */
764 }
765
766 /* build the table */
767 for (i = 2; i < nspcons; ++i)
768 {
769 for (j = 1; j < nblocks && j < i; ++j)
770 {
771 SCIP_Real weightleft;
772 SCIP_Real weightright;
773
774 assert( cases[i-1][j] != -1 );
775 assert( cases[i-1][j-1] != -1 );
776
777 weightleft = weights[i-1][j-1];
778 weightright = vals[i][j] + weights[i-1][j];
779
780 /* For first column: cannot take left possibility */
781 if ( SCIPisLT(scip, weightleft, weightright) )
782 {
783 weights[i][j] = weightleft;
784 cases[i][j] = 1;
785 }
786 else
787 {
788 weights[i][j] = weightright;
789 cases[i][j] = 2;
790 }
791 }
792 }
793 }
794
795
796 /** fix upper right triangle if necessary */
797 static
fixTriangle(SCIP * scip,SCIP_CONS * cons,SCIP_Bool * infeasible,int * nfixedvars)798 SCIP_RETCODE fixTriangle(
799 SCIP* scip, /**< SCIP data structure */
800 SCIP_CONS* cons, /**< constraint to be processed */
801 SCIP_Bool* infeasible, /**< pointer to store TRUE, if the node can be cut off */
802 int* nfixedvars /**< pointer to add up the number of found domain reductions */
803 )
804 {
805 SCIP_CONSDATA* consdata;
806 SCIP_VAR*** vars;
807 SCIP_Bool fixedglobal;
808 SCIP_Bool fixed;
809 int diagsize;
810 int nspcons;
811 int nblocks;
812 int i;
813 int j;
814
815 assert( scip != NULL );
816 assert( cons != NULL );
817 assert( infeasible != NULL );
818 assert( nfixedvars != NULL );
819
820 consdata = SCIPconsGetData(cons);
821 assert( consdata != NULL );
822 assert( consdata->nspcons > 0 );
823 assert( consdata->nblocks > 0 );
824 assert( consdata->vars != NULL );
825
826 *infeasible = FALSE;
827 *nfixedvars = 0;
828
829 if ( consdata->istrianglefixed )
830 return SCIP_OKAY;
831
832 nspcons = consdata->nspcons;
833 nblocks = consdata->nblocks;
834 vars = consdata->vars;
835 fixedglobal = TRUE;
836
837 /* get last row of triangle */
838 diagsize = nblocks;
839 if ( nspcons < nblocks )
840 diagsize = nspcons;
841
842 /* fix variables to 0 */
843 for (i = 0; i < diagsize; ++i)
844 {
845 for (j = i+1; j < nblocks; ++j)
846 {
847 /* fix variable, if not in the root the fixation is local */
848 SCIP_CALL( SCIPfixVar(scip, vars[i][j], 0.0, infeasible, &fixed) );
849
850 if ( *infeasible )
851 {
852 SCIPdebugMsg(scip, "The problem is infeasible: some variable in the upper right triangle is fixed to 1.\n");
853 return SCIP_OKAY;
854 }
855
856 if ( fixed )
857 ++(*nfixedvars);
858
859 if ( SCIPvarGetUbGlobal(vars[i][j]) > 0.5 )
860 fixedglobal = FALSE;
861 }
862 }
863 if ( *nfixedvars > 0 )
864 {
865 SCIPdebugMsg(scip, "<%s>: %s fixed upper right triangle to 0 (fixed vars: %d).\n", SCIPconsGetName(cons), fixedglobal ? "globally" : "locally", *nfixedvars);
866 }
867 else
868 {
869 SCIPdebugMsg(scip, "<%s>: Upper right triangle already fixed to 0.\n", SCIPconsGetName(cons));
870 }
871
872 if ( fixedglobal )
873 consdata->istrianglefixed = TRUE;
874
875 return SCIP_OKAY;
876 }
877
878
879 /** separates shifted column inequalities according to the solution stored in consdata->vals */
880 static
separateSCIs(SCIP * scip,SCIP_CONSHDLR * conshdlr,SCIP_CONS * cons,SCIP_CONSDATA * consdata,SCIP_Bool * infeasible,int * nfixedvars,int * ncuts)881 SCIP_RETCODE separateSCIs(
882 SCIP* scip, /**< the SCIP data structure */
883 SCIP_CONSHDLR* conshdlr, /**< constraint handler */
884 SCIP_CONS* cons, /**< constraint */
885 SCIP_CONSDATA* consdata, /**< the constraint data */
886 SCIP_Bool* infeasible, /**< whether we detected infeasibility */
887 int* nfixedvars, /**< pointer to store the number of variables fixed */
888 int* ncuts /**< pointer to store number of separated SCIs */
889 )
890 {
891 SCIP_Real** vals;
892 SCIP_Real** weights;
893 SCIP_Real* tmpvals;
894 SCIP_VAR*** vars;
895 SCIP_VAR** tmpvars;
896 int** cases;
897 int nspcons;
898 int nblocks;
899 int i;
900 int j;
901 int l;
902
903 assert( scip != NULL );
904 assert( conshdlr != NULL );
905 assert( cons != NULL );
906 assert( infeasible != NULL);
907 assert( nfixedvars != NULL );
908 assert( ncuts != NULL );
909
910 assert( consdata != NULL );
911 assert( consdata->nspcons > 0 );
912 assert( consdata->nblocks > 0 );
913 assert( consdata->vars != NULL );
914 assert( consdata->vals != NULL );
915 assert( consdata->tmpvars != NULL );
916 assert( consdata->tmpvals != NULL );
917 assert( consdata->weights != NULL );
918 assert( consdata->cases != NULL );
919
920 *infeasible = FALSE;
921 *nfixedvars = 0;
922 *ncuts = 0;
923
924 nspcons = consdata->nspcons;
925 nblocks = consdata->nblocks;
926 vars = consdata->vars;
927 vals = consdata->vals;
928 tmpvars = consdata->tmpvars;
929 tmpvals = consdata->tmpvals;
930 weights = consdata->weights;
931 cases = consdata->cases;
932
933 /* check for upper right triangle */
934 if ( ! consdata->istrianglefixed )
935 {
936 SCIP_CALL( fixTriangle(scip, cons, infeasible, nfixedvars) );
937 if ( *infeasible )
938 return SCIP_OKAY;
939 if ( *nfixedvars > 0 )
940 return SCIP_OKAY;
941 }
942
943 /* compute table if necessary (i.e., not computed before) */
944 computeSCTable(scip, nspcons, nblocks, weights, cases, vals);
945
946 /* loop through rows */
947 for (i = 1; i < nspcons && ! (*infeasible); ++i)
948 {
949 SCIP_Real bar; /* value of bar: */
950 int lastcolumn; /* last column considered as part of the bar */
951
952 bar = 0.0;
953 lastcolumn = nblocks - 1;
954 if ( lastcolumn > i )
955 lastcolumn = i;
956
957 /* traverse row from right to left: */
958 /* j >= 1, since for j = 0, i.e., the bar is a complete row, there does not exist an SCI */
959 for (j = lastcolumn; j > 0; --j)
960 {
961 bar += vals[i][j];
962
963 /* check whether weights[i-1][j-1] < bar (<=> bar - weights[i-1][j-1] > 0), i.e. cut is violated) */
964 if ( SCIPisEfficacious(scip, bar - weights[i-1][j-1]) )
965 {
966 SCIP_Real weight;
967 SCIP_ROW* row;
968 #ifdef SCIP_DEBUG
969 char name[SCIP_MAXSTRLEN];
970 #endif
971 int nvars;
972 int p1;
973 int p2;
974
975 nvars = 0;
976 p1 = i-1;
977 p2 = j-1;
978 weight = 0.0;
979
980 /* first add bar */
981 for (l = j; l <= lastcolumn; ++l)
982 {
983 tmpvars[nvars] = vars[i][l];
984 tmpvals[nvars] = 1.0;
985 nvars++;
986 }
987
988 /* then add shifted column */
989 do
990 {
991 assert( cases[p1][p2] != -1 );
992 assert( p1 >= 0 && p1 < i );
993 assert( p2 >= 0 && p2 < j );
994
995 /* if case 1 */
996 if (cases[p1][p2] == 1)
997 p2--; /* decrease column */
998 else
999 {
1000 /* case 2 or 3: */
1001 assert( cases[p1][p2] == 2 || cases[p1][p2] == 3 );
1002 tmpvars[nvars] = vars[p1][p2];
1003 tmpvals[nvars] = -1.0;
1004 nvars++;
1005 weight += vals[p1][p2];
1006 if ( cases[p1][p2] == 3 )
1007 break;
1008 }
1009 p1--; /* decrease row */
1010 }
1011 while ( p1 >= 0 ); /* should always be true, i.e. the break should end the loop */
1012 assert( cases[p1][p2] == 3 );
1013
1014 /* generate cut */
1015 #ifdef SCIP_DEBUG
1016 (void) SCIPsnprintf(name, SCIP_MAXSTRLEN, "sci_%d_%d", i, j);
1017 SCIP_CALL( SCIPcreateEmptyRowConshdlr(scip, &row, conshdlr, name, -SCIPinfinity(scip), 0.0, FALSE, FALSE, TRUE) );
1018 #else
1019 SCIP_CALL( SCIPcreateEmptyRowConshdlr(scip, &row, conshdlr, "", -SCIPinfinity(scip), 0.0, FALSE, FALSE, TRUE) );
1020 #endif
1021 SCIP_CALL( SCIPaddVarsToRow(scip, row, nvars, tmpvars, tmpvals) );
1022 /*SCIP_CALL( SCIPprintRow(scip, row, NULL) ); */
1023 SCIP_CALL( SCIPaddRow(scip, row, FALSE, infeasible) );
1024 SCIP_CALL( SCIPreleaseRow(scip, &row) );
1025 ++(*ncuts);
1026
1027 #ifdef SHOW_SCI
1028 SCIP_CALL( printSCI(scip, nspcons, nblocks, cases, i, j) );
1029 #endif
1030
1031 assert( SCIPisSumEQ(scip, weights[i-1][j-1], weight) );
1032 }
1033 }
1034 }
1035 return SCIP_OKAY;
1036 }
1037
1038
1039 /** propagation method for a single packing or partitioning orbitope constraint */
1040 static
propagatePackingPartitioningCons(SCIP * scip,SCIP_CONS * cons,SCIP_Bool * infeasible,int * nfixedvars)1041 SCIP_RETCODE propagatePackingPartitioningCons(
1042 SCIP* scip, /**< SCIP data structure */
1043 SCIP_CONS* cons, /**< constraint to be processed */
1044 SCIP_Bool* infeasible, /**< pointer to store TRUE, if the node can be cut off */
1045 int* nfixedvars /**< pointer to add up the number of found domain reductions */
1046 )
1047 {
1048 SCIP_CONSDATA* consdata;
1049 SCIP_VAR*** vars;
1050 SCIP_ORBITOPETYPE orbitopetype;
1051 int* firstnonzeros;
1052 int* lastones;
1053 int* frontiersteps;
1054 int lastoneprevrow;
1055 int nspcons;
1056 int nblocks;
1057 int nsteps;
1058 int i;
1059 int j;
1060
1061 assert( scip != NULL );
1062 assert( cons != NULL );
1063 assert( infeasible != NULL );
1064 assert( nfixedvars != NULL );
1065
1066 *nfixedvars = 0;
1067
1068 if( !SCIPallowStrongDualReds(scip) )
1069 return SCIP_OKAY;
1070
1071 consdata = SCIPconsGetData(cons);
1072 assert( consdata != NULL );
1073 assert( consdata->nspcons > 0 );
1074 assert( consdata->nblocks > 0 );
1075 assert( consdata->vars != NULL );
1076
1077 nspcons = consdata->nspcons;
1078 nblocks = consdata->nblocks;
1079 vars = consdata->vars;
1080 orbitopetype = consdata->orbitopetype;
1081
1082 assert( orbitopetype == SCIP_ORBITOPETYPE_PACKING || orbitopetype == SCIP_ORBITOPETYPE_PARTITIONING );
1083
1084 /* fix upper right triangle if still necessary */
1085 if ( ! consdata->istrianglefixed )
1086 {
1087 int nfixed = 0;
1088 SCIP_CALL( fixTriangle(scip, cons, infeasible, &nfixed) );
1089 *nfixedvars += nfixed;
1090 }
1091
1092 /* prepare further propagation */
1093 SCIP_CALL( SCIPallocBufferArray(scip, &firstnonzeros, nspcons) );
1094 SCIP_CALL( SCIPallocBufferArray(scip, &lastones, nspcons) );
1095 SCIP_CALL( SCIPallocBufferArray(scip, &frontiersteps, nblocks) );
1096
1097 #ifdef PRINT_MATRIX
1098 SCIPdebugMsg(scip, "Matrix:\n");
1099 printMatrix(scip, consdata);
1100 #endif
1101
1102 /* propagate */
1103 lastoneprevrow = 0;
1104 lastones[0] = 0;
1105
1106 if ( orbitopetype == SCIP_ORBITOPETYPE_PACKING )
1107 {
1108 /* packing case: if entry (0,0) is fixed to 0 */
1109 if ( SCIPvarGetUbLocal(vars[0][0]) < 0.5 )
1110 {
1111 lastoneprevrow = -1;
1112 lastones[0] = -1;
1113 }
1114 }
1115 nsteps = 0;
1116
1117 for (i = 1; i < nspcons; ++i)
1118 {
1119 int lastcolumn;
1120 int firstnonzeroinrow;
1121 int lastoneinrow;
1122 SCIP_Bool infrontier;
1123
1124 /* last column considered as part of the bar: */
1125 lastcolumn = nblocks - 1;
1126 if ( lastcolumn > i )
1127 lastcolumn = i;
1128
1129 /* find first position not fixed to 0 (partitioning) or fixed to 1 (packing) */
1130 firstnonzeroinrow = -1;
1131 for (j = 0; j <= lastcolumn; ++j)
1132 {
1133 if ( orbitopetype == SCIP_ORBITOPETYPE_PARTITIONING )
1134 {
1135 /* partitioning case: if variable is not fixed to 0 */
1136 if ( SCIPvarGetUbLocal(vars[i][j]) > 0.5 )
1137 {
1138 firstnonzeroinrow = j;
1139 break;
1140 }
1141 }
1142 else
1143 {
1144 /* packing case: if variable is fixed to 1 */
1145 if ( SCIPvarGetLbLocal(vars[i][j]) > 0.5 )
1146 {
1147 firstnonzeroinrow = j;
1148 break;
1149 }
1150 }
1151 }
1152 /* if all variables are fixed to 0 in the partitioning case - should not happen */
1153 if ( firstnonzeroinrow == -1 && orbitopetype == SCIP_ORBITOPETYPE_PARTITIONING )
1154 {
1155 SCIPdebugMsg(scip, " -> Infeasible node: all variables in row %d are fixed to 0.\n", i);
1156 *infeasible = TRUE;
1157 /* conflict should be analyzed by setppc constraint handler */
1158 goto TERMINATE;
1159 }
1160 firstnonzeros[i] = firstnonzeroinrow;
1161 assert( orbitopetype == SCIP_ORBITOPETYPE_PACKING || firstnonzeroinrow >= 0 );
1162 assert( -1 <= firstnonzeroinrow && firstnonzeroinrow <= lastcolumn );
1163
1164 /* compute rightmost possible position for a 1 */
1165 assert( orbitopetype == SCIP_ORBITOPETYPE_PACKING || 0 <= lastoneprevrow );
1166 assert( lastoneprevrow <= lastcolumn );
1167
1168 /* if we are at right border or if entry in column lastoneprevrow+1 is fixed to 0 */
1169 infrontier = FALSE;
1170 assert( lastoneprevrow + 1 >= 0 );
1171 if ( lastoneprevrow == nblocks-1 || SCIPvarGetUbLocal(vars[i][lastoneprevrow+1]) < 0.5 ) /*lint !e679*/
1172 lastoneinrow = lastoneprevrow;
1173 else
1174 {
1175 lastoneinrow = lastoneprevrow + 1;
1176 frontiersteps[nsteps++] = i;
1177 infrontier = TRUE;
1178 }
1179
1180 /* store lastoneinrow */
1181 assert( orbitopetype == SCIP_ORBITOPETYPE_PACKING || 0 <= lastoneinrow );
1182 assert( lastoneinrow <= lastcolumn );
1183 lastones[i] = lastoneinrow;
1184
1185 /* check whether we are infeasible */
1186 if ( firstnonzeroinrow > lastoneinrow )
1187 {
1188 int k;
1189
1190 #ifdef SCIP_DEBUG
1191 if ( orbitopetype == SCIP_ORBITOPETYPE_PARTITIONING )
1192 {
1193 SCIPdebugMsg(scip, " -> Infeasible node: row %d, leftmost nonzero at %d, rightmost 1 at %d\n",
1194 i, firstnonzeroinrow, lastoneinrow);
1195 }
1196 else
1197 {
1198 SCIPdebugMsg(scip, " -> Infeasible node: row %d, 1 at %d, rightmost position for 1 at %d\n",
1199 i, firstnonzeroinrow, lastoneinrow);
1200 }
1201 #endif
1202 /* check if conflict analysis is applicable */
1203 if ( SCIPisConflictAnalysisApplicable(scip) )
1204 {
1205 /* conflict analysis only applicable in SOLVING stage */
1206 assert( SCIPgetStage(scip) == SCIP_STAGE_SOLVING || SCIPinProbing(scip) );
1207
1208 /* perform conflict analysis */
1209 SCIP_CALL( SCIPinitConflictAnalysis(scip, SCIP_CONFTYPE_PROPAGATION, FALSE) );
1210
1211 if ( orbitopetype == SCIP_ORBITOPETYPE_PARTITIONING )
1212 {
1213 /* add bounds (variables fixed to 0) that result in the first nonzero entry */
1214 for (j = 0; j <= lastcolumn; ++j)
1215 {
1216 /* add varaibles in row up to the first variable fixed to 0 */
1217 if ( SCIPvarGetUbLocal(vars[i][j]) > 0.5 )
1218 break;
1219
1220 assert( SCIPvarGetUbLocal(vars[i][j]) < 0.5 );
1221 SCIP_CALL( SCIPaddConflictBinvar(scip, vars[i][j]) );
1222 }
1223 }
1224 else
1225 {
1226 /* add bounds that result in the last one - check top left entry for packing case */
1227 if ( lastones[0] == -1 )
1228 {
1229 assert( SCIPvarGetUbLocal(vars[0][0]) < 0.5 );
1230 SCIP_CALL( SCIPaddConflictBinvar(scip, vars[0][0]) );
1231 }
1232
1233 /* mark variable fixed to 1 */
1234 assert( SCIPvarGetLbLocal(vars[i][firstnonzeroinrow]) > 0.5 );
1235 SCIP_CALL( SCIPaddConflictBinvar(scip, vars[i][firstnonzeroinrow]) );
1236 }
1237
1238 /* add bounds that result in the last one - pass through rows */
1239 for (k = 1; k < i; ++k)
1240 {
1241 int l;
1242 l = lastones[k] + 1;
1243
1244 /* if the frontier has not moved and we are not beyond the matrix boundaries */
1245 if ( l <= nblocks-1 && l <= k && lastones[k-1] == lastones[k] )
1246 {
1247 assert( SCIPvarGetUbLocal(vars[k][l]) < 0.5 );
1248 SCIP_CALL( SCIPaddConflictBinvar(scip, vars[k][l]) );
1249 }
1250 }
1251 SCIP_CALL( SCIPanalyzeConflictCons(scip, cons, NULL) );
1252 }
1253
1254 *infeasible = TRUE;
1255 goto TERMINATE;
1256 }
1257
1258 /* fix entries beyond the last possible position for a 1 in the row to 0 (see Lemma 1 in the paper) */
1259 for (j = lastoneinrow+1; j <= lastcolumn; ++j)
1260 {
1261 /* if the entry is not yet fixed to 0 */
1262 if ( SCIPvarGetUbLocal(vars[i][j]) > 0.5 )
1263 {
1264 SCIP_Bool tightened;
1265 int inferInfo;
1266
1267 SCIPdebugMsg(scip, " -> Fixing entry (%d,%d) to 0.\n", i, j);
1268
1269 tightened = FALSE;
1270
1271 /* fix variable to 0 and store position of (i,lastoneinrow+1) for conflict resolution */
1272 inferInfo = i * nblocks + lastoneinrow + 1;
1273 /* correction according to Lemma 1 in the paper (second part): store (i,lastoneinrow+2) */
1274 if ( !infrontier )
1275 ++inferInfo;
1276 SCIP_CALL( SCIPinferBinvarCons(scip, vars[i][j], FALSE, cons, inferInfo, infeasible, &tightened) );
1277
1278 /* if entry is fixed to one -> infeasible node */
1279 if ( *infeasible )
1280 {
1281 SCIPdebugMsg(scip, " -> Infeasible node: row %d, 1 in column %d beyond rightmost position %d\n", i, j, lastoneinrow);
1282 /* check if conflict analysis is applicable */
1283 if( SCIPisConflictAnalysisApplicable(scip) )
1284 {
1285 int k;
1286
1287 /* conflict analysis only applicable in SOLVING stage */
1288 assert(SCIPgetStage(scip) == SCIP_STAGE_SOLVING || SCIPinProbing(scip));
1289
1290 /* perform conflict analysis */
1291 SCIP_CALL( SCIPinitConflictAnalysis(scip, SCIP_CONFTYPE_PROPAGATION, FALSE) );
1292
1293 /* add current bound */
1294 SCIP_CALL( SCIPaddConflictBinvar(scip, vars[i][j]) );
1295
1296 /* add bounds that result in the last one - check top left entry for packing case */
1297 if ( orbitopetype == SCIP_ORBITOPETYPE_PACKING && lastones[0] == -1 )
1298 {
1299 assert( SCIPvarGetUbLocal(vars[0][0]) < 0.5 );
1300 SCIP_CALL( SCIPaddConflictBinvar(scip, vars[0][0]) );
1301 }
1302
1303 /* add bounds that result in the last one - pass through rows */
1304 for (k = 1; k < i; ++k)
1305 {
1306 int l;
1307 l = lastones[k] + 1;
1308
1309 /* if the frontier has not moved and we are not beyond the matrix boundaries */
1310 if ( l <= nblocks-1 && l <= k && lastones[k-1] == lastones[k] )
1311 {
1312 assert( SCIPvarGetUbLocal(vars[k][l]) < 0.5 );
1313 SCIP_CALL( SCIPaddConflictBinvar(scip, vars[k][l]) );
1314 }
1315 }
1316 SCIP_CALL( SCIPanalyzeConflictCons(scip, cons, NULL) );
1317 }
1318
1319 goto TERMINATE;
1320 }
1321 if ( tightened )
1322 ++(*nfixedvars);
1323 }
1324 }
1325
1326 lastoneprevrow = lastoneinrow;
1327 }
1328
1329 /* check whether fixing any entry to 0 results in a contradiction -> loop through rows in frontiersteps (a.k.a. gamma) */
1330 for (j = 0; j < nsteps; ++j)
1331 {
1332 int s;
1333 int lastoneinrow;
1334
1335 s = frontiersteps[j];
1336 lastoneinrow = lastones[s];
1337 /* note for packing case: if we are in a frontier step then lastoneinrow >= 0 */
1338 assert( 0 <= lastoneinrow && lastoneinrow < nblocks );
1339
1340 /* if entry is not fixed */
1341 if ( SCIPvarGetLbLocal(vars[s][lastoneinrow]) < 0.5 && SCIPvarGetUbLocal(vars[s][lastoneinrow]) > 0.5 )
1342 {
1343 int betaprev;
1344 betaprev = lastoneinrow - 1;
1345
1346 /* loop through rows below s */
1347 for (i = s+1; i < nspcons; ++i)
1348 {
1349 int beta;
1350
1351 assert( betaprev + 1 >= 0 );
1352 if ( betaprev == nblocks-1 || SCIPvarGetUbLocal(vars[i][betaprev+1]) < 0.5 ) /*lint !e679*/
1353 beta = betaprev;
1354 else
1355 beta = betaprev + 1;
1356 assert( -1 <= beta && beta < nblocks );
1357
1358 if ( firstnonzeros[i] > beta )
1359 {
1360 SCIP_Bool tightened;
1361 int inferInfo;
1362
1363 /* can fix (s,lastoneinrow) (a.k.a (s,alpha)) to 1
1364 * (do not need to fix other entries to 0, since they will be
1365 * automatically fixed by SCIPtightenVarLb.)
1366 */
1367 assert( SCIPvarGetLbLocal(vars[s][lastoneinrow]) < 0.5 );
1368 SCIPdebugMsg(scip, " -> Fixing entry (%d,%d) to 1.\n", s, lastoneinrow);
1369
1370 tightened = FALSE;
1371
1372 /* store position (i,firstnonzeros[i]) */
1373 inferInfo = nblocks * nspcons + i * nblocks + firstnonzeros[i];
1374 SCIP_CALL( SCIPinferBinvarCons(scip, vars[s][lastoneinrow], TRUE, cons, inferInfo, infeasible, &tightened) );
1375
1376 assert( !(*infeasible) );
1377 if ( tightened )
1378 ++(*nfixedvars);
1379 break;
1380 }
1381 betaprev = beta;
1382 }
1383 }
1384 }
1385
1386 TERMINATE:
1387 SCIPfreeBufferArray(scip, &frontiersteps);
1388 SCIPfreeBufferArray(scip, &lastones);
1389 SCIPfreeBufferArray(scip, &firstnonzeros);
1390
1391 return SCIP_OKAY;
1392 }
1393
1394
1395 /* Compute dynamic order of rows based on the branching decisions, i.e., the row of the first branching variable
1396 * determines the first row in the new ordering, the row of the second branching variable determines the second
1397 * row in the new ordering if it differs from the row of the first branching variable, and so on.
1398 *
1399 * The roworder array stores this reordering, where acutally only the first maxrowlabel entries encode the
1400 * reordering.
1401 */
1402 static
computeDynamicRowOrder(SCIP * scip,SCIP_HASHMAP * rowindexmap,SCIP_Bool * rowused,int * roworder,int nrows,int ncols,int * maxrowlabel)1403 SCIP_RETCODE computeDynamicRowOrder(
1404 SCIP* scip, /**< SCIP pointer */
1405 SCIP_HASHMAP* rowindexmap, /**< map of variables to indices in orbitope vars matrix */
1406 SCIP_Bool* rowused, /**< bitset marking whether a row has been considered in the new order */
1407 int* roworder, /**< reordering of the rows w.r.t. branching decisions */
1408 int nrows, /**< number of rows in orbitope */
1409 int ncols, /**< number of columns in orbitope */
1410 int* maxrowlabel /**< maximum row label in ordering */
1411 )
1412 {
1413 int i;
1414 SCIP_NODE* node;
1415 int* branchdecisions;
1416 int nbranchdecision;
1417
1418 assert( scip != NULL );
1419 assert( rowindexmap != NULL );
1420 assert( roworder != NULL );
1421 assert( nrows > 0 );
1422 assert( maxrowlabel != NULL );
1423
1424 SCIP_CALL( SCIPallocBufferArray(scip, &branchdecisions, nrows * ncols) );
1425 nbranchdecision = 0;
1426
1427 /* get current node */
1428 node = SCIPgetCurrentNode(scip);
1429
1430 /* follow path to the root (in the root no domains were changed due to branching) */
1431 while ( SCIPnodeGetDepth(node) != 0 )
1432 {
1433 SCIP_BOUNDCHG* boundchg;
1434 SCIP_DOMCHG* domchg;
1435 SCIP_VAR* branchvar;
1436 int nboundchgs;
1437
1438 /* get domain changes of current node */
1439 domchg = SCIPnodeGetDomchg(node);
1440 assert( domchg != NULL );
1441
1442 /* loop through all bound changes */
1443 nboundchgs = SCIPdomchgGetNBoundchgs(domchg);
1444 for (i = 0; i < nboundchgs; ++i)
1445 {
1446 /* get bound change info */
1447 boundchg = SCIPdomchgGetBoundchg(domchg, i);
1448 assert( boundchg != NULL );
1449
1450 /* branching decisions have to be in the beginning of the bound change array */
1451 if ( SCIPboundchgGetBoundchgtype(boundchg) != SCIP_BOUNDCHGTYPE_BRANCHING )
1452 break;
1453
1454 /* get corresponding branching variable */
1455 branchvar = SCIPboundchgGetVar(boundchg);
1456
1457 /* we only consider binary variables */
1458 if ( SCIPvarGetType(branchvar) == SCIP_VARTYPE_BINARY )
1459 {
1460 int rowidx;
1461
1462 /* make sure that branching variable is present in the orbitope */
1463 if ( ! SCIPhashmapExists(rowindexmap, (void*) branchvar) )
1464 continue;
1465
1466 rowidx = (int) (size_t) SCIPhashmapGetImage(rowindexmap, (void*) branchvar);
1467 branchdecisions[nbranchdecision++] = rowidx;
1468 }
1469 }
1470
1471 node = SCIPnodeGetParent(node);
1472 }
1473
1474 /* Insert branching decisions of current path into global row order.
1475 * Iterate in reverse order over branching decisions to get the path
1476 * from the root to the current node.
1477 */
1478 for (i = nbranchdecision - 1; i >= 0; --i)
1479 {
1480 if ( ! rowused[branchdecisions[i]] )
1481 {
1482 roworder[*maxrowlabel] = branchdecisions[i];
1483 rowused[branchdecisions[i]] = TRUE;
1484 *maxrowlabel += 1;
1485 }
1486 }
1487
1488 SCIPfreeBufferArray(scip, &branchdecisions);
1489
1490 return SCIP_OKAY;
1491 }
1492
1493
1494 /* Compute lexicographically minimal face of the hypercube w.r.t. some coordinate fixing */
1495 static
findLexMinFace(SCIP_VAR *** vars,int ** lexminfixes,int * minfixedrowlexmin,SCIP_Bool * infeasible,int m,int n,int nrowsused,SCIP_Bool resprop)1496 SCIP_RETCODE findLexMinFace(
1497 SCIP_VAR*** vars, /**< variable matrix */
1498 int** lexminfixes, /**< fixings characterzing lex-min face */
1499 int* minfixedrowlexmin, /**< index of minimum fixed row for each column or
1500 * NULL (if in prop) */
1501 SCIP_Bool* infeasible, /**< pointer to store whether infeasibility has been
1502 * detected or NULL (if in resprop) */
1503 int m, /**< number of rows in vars */
1504 int n, /**< number of columns in vars */
1505 int nrowsused, /**< number of rows considered in propagation */
1506 SCIP_Bool resprop /**< whether we are in resprop (TRUE) or prop (FALSE) */
1507 )
1508 {
1509 int i;
1510 int j;
1511 *infeasible = FALSE;
1512
1513 assert( vars != NULL );
1514 assert( lexminfixes != NULL );
1515 assert( !resprop || minfixedrowlexmin != NULL );
1516 assert( m > 0 );
1517 assert( n > 0 );
1518 assert( nrowsused > 0 );
1519 assert( nrowsused <= m );
1520 assert( infeasible != NULL );
1521
1522 /* iterate over columns in reverse order and find the lexicographically minimal face
1523 * of the hypercube containing lexminfixes
1524 */
1525 for (j = n - 2; j >= 0; --j)
1526 {
1527 int maxdiscriminating = m;
1528 int minfixed = -1;
1529
1530 /* fix free entries in column j to the corresponding value in column j + 1 and collect some information */
1531 for (i = 0; i < nrowsused; ++i)
1532 {
1533 /* is row i j-discriminating? */
1534 if ( minfixed == -1 && lexminfixes[i][j] != 0 && lexminfixes[i][j + 1] != 1 )
1535 {
1536 assert( lexminfixes[i][j + 1] == 0 );
1537
1538 maxdiscriminating = i;
1539 }
1540
1541 /* is row i j-fixed? */
1542 if ( minfixed == -1 && lexminfixes[i][j] != lexminfixes[i][j + 1] && lexminfixes[i][j] != 2 )
1543 {
1544 assert( lexminfixes[i][j + 1] != 2 );
1545
1546 minfixed = i;
1547
1548 /* detect infeasibility */
1549 if ( maxdiscriminating > minfixed )
1550 {
1551 *infeasible = TRUE;
1552
1553 return SCIP_OKAY;
1554 }
1555 }
1556 }
1557
1558 /* ensure that column j is lexicographically not smaller than column j + 1 */
1559 for (i = 0; i < nrowsused; ++i)
1560 {
1561 if ( lexminfixes[i][j] == 2 )
1562 {
1563 if ( i < maxdiscriminating || minfixed == -1 )
1564 lexminfixes[i][j] = lexminfixes[i][j + 1];
1565 else if ( i == maxdiscriminating )
1566 lexminfixes[i][j] = 1;
1567 else
1568 lexminfixes[i][j] = 0;
1569 }
1570 }
1571
1572 if ( resprop )
1573 {
1574 assert( minfixedrowlexmin != NULL );
1575
1576 /* store minimum fixed row */
1577 if ( minfixed == -1 )
1578 minfixedrowlexmin[j] = nrowsused - 1;
1579 else
1580 minfixedrowlexmin[j] = minfixed;
1581
1582 /* columns 1, ..., n-2 are contained in two columns (take the minimum) and
1583 * the minimum fixed row of column n-1 is determined by column n-2 */
1584 if ( minfixedrowlexmin[j + 1] < minfixedrowlexmin[j] )
1585 minfixedrowlexmin[j + 1] = minfixedrowlexmin[j];
1586 }
1587 }
1588
1589 return SCIP_OKAY;
1590 }
1591
1592
1593 /* Compute lexicographically maximal face of the hypercube w.r.t. some coordinate fixing */
1594 static
findLexMaxFace(SCIP_VAR *** vars,int ** lexmaxfixes,int * minfixedrowlexmax,SCIP_Bool * infeasible,int m,int n,int nrowsused,SCIP_Bool resprop)1595 SCIP_RETCODE findLexMaxFace(
1596 SCIP_VAR*** vars, /**< variable matrix */
1597 int** lexmaxfixes, /**< fixings characterzing lex-max face */
1598 int* minfixedrowlexmax, /**< index of minimum fixed row for each column or
1599 * NULL (if in prop) */
1600 SCIP_Bool* infeasible, /**< pointer to store whether infeasibility has been
1601 * detected or NULL (if in resprop) */
1602 int m, /**< number of rows in vars */
1603 int n, /**< number of columns in vars */
1604 int nrowsused, /**< number of rows considered in propagation */
1605 SCIP_Bool resprop /**< whether we are in resprop (TRUE) or prop (FALSE) */
1606 )
1607 {
1608 int i;
1609 int j;
1610 *infeasible = FALSE;
1611
1612 assert( vars != NULL );
1613 assert( lexmaxfixes != NULL );
1614 assert( !resprop || minfixedrowlexmax != NULL );
1615 assert( m > 0 );
1616 assert( n > 0 );
1617 assert( nrowsused > 0 );
1618 assert( nrowsused <= m );
1619 assert( infeasible != NULL );
1620
1621 for (j = 1; j < n; ++j)
1622 {
1623 int maxdiscriminating = m;
1624 int minfixed = -1;
1625
1626 /* fix free entries in column j to the corresponding value in column j - 1 and collect some information */
1627 for (i = 0; i < nrowsused; ++i)
1628 {
1629 /* is row i j-discriminating? */
1630 if ( minfixed == -1 && lexmaxfixes[i][j - 1] != 0 && lexmaxfixes[i][j] != 1 )
1631 {
1632 assert( lexmaxfixes[i][j - 1] == 1 );
1633
1634 maxdiscriminating = i;
1635 }
1636
1637 /* is row i j-fixed? */
1638 if ( minfixed == -1 && lexmaxfixes[i][j - 1] != lexmaxfixes[i][j] && lexmaxfixes[i][j] != 2 )
1639 {
1640 assert( lexmaxfixes[i][j - 1] != 2 );
1641
1642 minfixed = i;
1643
1644 /* detect infeasibility */
1645 if ( maxdiscriminating > minfixed )
1646 {
1647 *infeasible = TRUE;
1648
1649 return SCIP_OKAY;
1650 }
1651 }
1652 }
1653
1654 /* ensure that column j is lexicographically not greater than column j - 1 */
1655 for (i = 0; i < nrowsused; ++i)
1656 {
1657 if ( lexmaxfixes[i][j] == 2 )
1658 {
1659 if ( i < maxdiscriminating || minfixed == -1 )
1660 lexmaxfixes[i][j] = lexmaxfixes[i][j - 1];
1661 else if ( i == maxdiscriminating )
1662 lexmaxfixes[i][j] = 0;
1663 else
1664 lexmaxfixes[i][j] = 1;
1665 }
1666 }
1667
1668 if ( resprop )
1669 {
1670 assert( minfixedrowlexmax != NULL );
1671
1672 /* store minimum fixed row */
1673 if ( minfixed == -1 )
1674 minfixedrowlexmax[j] = nrowsused - 1;
1675 else
1676 minfixedrowlexmax[j] = minfixed;
1677
1678 /* columns 1, ..., n-2 are contained in two columns (take the minimum) and
1679 * the minimum fixed row of column 0 is determined by column 1 */
1680 if ( minfixedrowlexmax[j - 1] < minfixedrowlexmax[j] )
1681 minfixedrowlexmax[j - 1] = minfixedrowlexmax[j];
1682 }
1683 }
1684
1685 return SCIP_OKAY;
1686 }
1687
1688
1689 /** propagation method for a single packing or partitioning orbitope constraint */
1690 static
propagateFullOrbitopeCons(SCIP * scip,SCIP_CONS * cons,SCIP_Bool * infeasible,int * nfixedvars,SCIP_Bool dynamic)1691 SCIP_RETCODE propagateFullOrbitopeCons(
1692 SCIP* scip, /**< SCIP data structure */
1693 SCIP_CONS* cons, /**< constraint to be processed */
1694 SCIP_Bool* infeasible, /**< pointer to store TRUE, if the node can be cut off */
1695 int* nfixedvars, /**< pointer to add up the number of found domain reductions */
1696 SCIP_Bool dynamic /**< whether we use a dynamic propagation routine */
1697 )
1698 {
1699 SCIP_CONSDATA* consdata;
1700 SCIP_VAR*** vars;
1701 int** lexminfixes;
1702 int** lexmaxfixes;
1703 int* roworder;
1704 int nrowsused;
1705 int i;
1706 int j;
1707 int m;
1708 int n;
1709
1710 assert( scip != NULL );
1711 assert( cons != NULL );
1712 assert( infeasible != NULL );
1713 assert( nfixedvars != NULL );
1714
1715 *nfixedvars = 0;
1716 *infeasible = FALSE;
1717
1718 /* @todo Can the following be removed? */
1719 if ( ! SCIPallowStrongDualReds(scip) )
1720 return SCIP_OKAY;
1721
1722 /* do nothing if we use dynamic propagation and if we are in a probing node */
1723 if ( dynamic && SCIPinProbing(scip) )
1724 return SCIP_OKAY;
1725
1726 consdata = SCIPconsGetData(cons);
1727 assert( consdata != NULL );
1728 assert( consdata->nspcons > 0 );
1729 assert( consdata->nblocks > 0 );
1730 assert( consdata->vars != NULL );
1731 assert( consdata->orbitopetype == SCIP_ORBITOPETYPE_FULL );
1732
1733 m = consdata->nspcons;
1734 n = consdata->nblocks;
1735 vars = consdata->vars;
1736
1737 /* determine order of orbitope rows dynamically by branching decisions */
1738 if ( dynamic )
1739 {
1740 SCIP_CALL( computeDynamicRowOrder(scip, consdata->rowindexmap, consdata->rowused,
1741 consdata->roworder, m, n, &(consdata->nrowsused)) );
1742
1743 /* if no branching variable is contained in the full orbitope */
1744 if ( consdata->nrowsused == 0 )
1745 return SCIP_OKAY;
1746
1747 nrowsused = consdata->nrowsused;
1748 }
1749 else
1750 nrowsused = m;
1751 roworder = consdata->roworder;
1752
1753 /* Initialize lexicographically minimal matrix by fixed entries at the current node.
1754 * Free entries in the last column are set to 0.
1755 */
1756 SCIP_CALL( SCIPallocBufferArray(scip, &lexminfixes, nrowsused) );
1757 for (i = 0; i < nrowsused; ++i)
1758 {
1759 SCIP_CALL( SCIPallocBufferArray(scip, &lexminfixes[i], n) ); /*lint !e866*/
1760 }
1761
1762 for (i = 0; i < nrowsused; ++i)
1763 {
1764 int origrow;
1765
1766 origrow = roworder[i];
1767
1768 for (j = 0; j < n; ++j)
1769 {
1770 if ( SCIPvarGetLbLocal(vars[origrow][j]) > 0.5 )
1771 lexminfixes[i][j] = 1;
1772 else if ( SCIPvarGetUbLocal(vars[origrow][j]) < 0.5 || j == n - 1 )
1773 lexminfixes[i][j] = 0;
1774 else
1775 lexminfixes[i][j] = 2;
1776 }
1777 }
1778
1779 /* find lexicographically minimal face of hypercube containing lexmin fixes */
1780 SCIP_CALL( findLexMinFace(vars, lexminfixes, NULL, infeasible, m, n, nrowsused, FALSE) );
1781
1782 if ( *infeasible == TRUE )
1783 goto FREELEXMIN;
1784
1785 /* Initialize lexicographically maximal matrix by fixed entries at the current node.
1786 * Free entries in the first column are set to 1.
1787 */
1788 SCIP_CALL( SCIPallocBufferArray(scip, &lexmaxfixes, nrowsused) );
1789 for (i = 0; i < nrowsused; ++i)
1790 {
1791 SCIP_CALL( SCIPallocBufferArray(scip, &lexmaxfixes[i], n) ); /*lint !e866*/
1792 }
1793
1794 for (i = 0; i < nrowsused; ++i)
1795 {
1796 int origrow;
1797
1798 origrow = roworder[i];
1799
1800 for (j = 0; j < n; ++j)
1801 {
1802 if ( SCIPvarGetUbLocal(vars[origrow][j]) < 0.5 )
1803 lexmaxfixes[i][j] = 0;
1804 else if ( SCIPvarGetLbLocal(vars[origrow][j]) > 0.5 || j == 0 )
1805 lexmaxfixes[i][j] = 1;
1806 else
1807 lexmaxfixes[i][j] = 2;
1808 }
1809 }
1810
1811 /* find lexicographically maximal face of hypercube containing lexmax fixes */
1812 SCIP_CALL( findLexMaxFace(vars, lexmaxfixes, NULL, infeasible, m, n, nrowsused, FALSE) );
1813
1814 if ( *infeasible )
1815 goto FREELEXMAX;
1816
1817 /* Find for each column j the minimal row in which lexminfixes and lexmaxfixes differ. Fix all entries above this
1818 * row to the corresponding value in lexminfixes (or lexmaxfixes).
1819 */
1820 for (j = 0; j < n; ++j)
1821 {
1822 for (i = 0; i < nrowsused; ++i)
1823 {
1824 int origrow;
1825
1826 origrow = roworder[i];
1827
1828 if ( lexminfixes[i][j] != lexmaxfixes[i][j] )
1829 break;
1830
1831 if ( SCIPvarGetLbLocal(vars[origrow][j]) < 0.5 && SCIPvarGetUbLocal(vars[origrow][j]) > 0.5 )
1832 {
1833 SCIP_Bool success;
1834
1835 SCIP_CALL( SCIPinferBinvarCons(scip, vars[origrow][j], (SCIP_Bool) lexminfixes[i][j],
1836 cons, 0, infeasible, &success) );
1837
1838 if ( success )
1839 *nfixedvars += 1;
1840 }
1841 }
1842 }
1843
1844 FREELEXMAX:
1845 for (i = 0; i < nrowsused; ++i)
1846 SCIPfreeBufferArray(scip, &lexmaxfixes[i]);
1847 SCIPfreeBufferArray(scip, &lexmaxfixes);
1848
1849 FREELEXMIN:
1850 for (i = 0; i < nrowsused; ++i)
1851 SCIPfreeBufferArray(scip, &lexminfixes[i]);
1852 SCIPfreeBufferArray(scip, &lexminfixes);
1853
1854 return SCIP_OKAY;
1855 }
1856
1857
1858 /** propagation method for a single orbitope constraint */
1859 static
propagateCons(SCIP * scip,SCIP_CONS * cons,SCIP_Bool * infeasible,int * nfixedvars,SCIP_Bool usedynamicprop)1860 SCIP_RETCODE propagateCons(
1861 SCIP* scip, /**< SCIP data structure */
1862 SCIP_CONS* cons, /**< constraint to be processed */
1863 SCIP_Bool* infeasible, /**< pointer to store TRUE, if the node can be cut off */
1864 int* nfixedvars, /**< pointer to add up the number of found domain reductions */
1865 SCIP_Bool usedynamicprop /**< whether we use a dynamic version of the propagation routine */
1866 )
1867 {
1868 SCIP_CONSDATA* consdata;
1869 SCIP_ORBITOPETYPE orbitopetype;
1870
1871 assert( scip != NULL );
1872 assert( cons != NULL );
1873 assert( infeasible != NULL );
1874 assert( nfixedvars != NULL );
1875
1876 consdata = SCIPconsGetData(cons);
1877 assert( consdata != NULL );
1878
1879 orbitopetype = consdata->orbitopetype;
1880
1881 if ( orbitopetype == SCIP_ORBITOPETYPE_FULL )
1882 {
1883 SCIP_CALL( propagateFullOrbitopeCons(scip, cons, infeasible, nfixedvars, usedynamicprop && !consdata->ismodelcons) );
1884 }
1885 else
1886 {
1887 assert( orbitopetype == SCIP_ORBITOPETYPE_PACKING || orbitopetype == SCIP_ORBITOPETYPE_PARTITIONING );
1888 SCIP_CALL( propagatePackingPartitioningCons(scip, cons, infeasible, nfixedvars) );
1889 }
1890
1891 return SCIP_OKAY;
1892 }
1893
1894
1895 /** Propagation conflict resolving method of propagator
1896 *
1897 * In this function we use that all variable reductions that can be found by the propagation algorithm
1898 * are only due to the fixed variables that are in or above the minimum fixed row of each pair of adjacent
1899 * columns of the lexmin and lexmax matrices.
1900 *
1901 * Since the storage of an integer is not enough to store the complete information about the fixing,
1902 * we have to use the linear time algorithm for finding the lexmin and lexmax
1903 * matrices and determine from this the minimum fixed rows.
1904 */
1905 static
resolvePropagationFullOrbitope(SCIP * scip,SCIP_CONSHDLR * conshdlr,SCIP_CONS * cons,int inferinfo,SCIP_BDCHGIDX * bdchgidx,SCIP_RESULT * result)1906 SCIP_RETCODE resolvePropagationFullOrbitope(
1907 SCIP* scip, /**< SCIP data structure */
1908 SCIP_CONSHDLR* conshdlr, /**< constraint handler of the corresponding constraint */
1909 SCIP_CONS* cons, /**< constraint that inferred the bound change */
1910 int inferinfo, /**< inference information */
1911 SCIP_BDCHGIDX* bdchgidx, /**< bound change index (time stamp of bound change), or NULL for current time */
1912 SCIP_RESULT* result /**< pointer to store the result of the propagation conflict resolving call */
1913 )
1914 { /*lint --e{715}*/
1915 SCIP_CONSHDLRDATA* conshdlrdata;
1916 SCIP_CONSDATA* consdata;
1917 SCIP_VAR*** vars;
1918 int** lexminfixes;
1919 int** lexmaxfixes;
1920 int* roworder;
1921 int* minfixedrowlexmin;
1922 int* minfixedrowlexmax;
1923 int i;
1924 int j;
1925 int m;
1926 int n;
1927 int nrowsused;
1928 SCIP_Bool dynamic;
1929 SCIP_Bool terminate;
1930
1931 *result = SCIP_DIDNOTFIND;
1932
1933 assert( scip != NULL );
1934 assert( conshdlr != NULL );
1935 assert( cons != NULL );
1936 assert( result != NULL );
1937
1938 consdata = SCIPconsGetData(cons);
1939 assert( consdata != NULL );
1940 assert( consdata->nspcons > 0 );
1941 assert( consdata->nblocks > 0 );
1942 assert( consdata->vars != NULL );
1943 assert( consdata->orbitopetype == SCIP_ORBITOPETYPE_FULL );
1944
1945 conshdlrdata = SCIPconshdlrGetData(conshdlr);
1946 dynamic = conshdlrdata->usedynamicprop && !consdata->ismodelcons;
1947 m = consdata->nspcons;
1948 n = consdata->nblocks;
1949 vars = consdata->vars;
1950
1951 if ( dynamic )
1952 {
1953 assert( consdata->roworder != NULL );
1954 assert( consdata->nrowsused > 0 );
1955
1956 nrowsused = consdata->nrowsused;
1957 }
1958 else
1959 nrowsused = m;
1960 roworder = consdata->roworder;
1961
1962 assert( inferinfo <= consdata->nspcons );
1963
1964 /* Initialize lexicographically minimal matrix by fixed entries at the current node.
1965 * Free entries in the last column are set to 0.
1966 */
1967 SCIP_CALL( SCIPallocBufferArray(scip, &lexminfixes, nrowsused) );
1968 for (i = 0; i < nrowsused; ++i)
1969 {
1970 SCIP_CALL( SCIPallocBufferArray(scip, &lexminfixes[i], n) ); /*lint !e866*/
1971 }
1972
1973 /* store minimum fixed row for each column */
1974 SCIP_CALL( SCIPallocBufferArray(scip, &minfixedrowlexmin, n) );
1975 minfixedrowlexmin[n - 1] = -1;
1976
1977 for (i = 0; i < nrowsused; ++i)
1978 {
1979 int origrow;
1980
1981 origrow = roworder[i];
1982
1983 for (j = 0; j < n; ++j)
1984 {
1985 if ( SCIPvarGetLbAtIndex(vars[origrow][j], bdchgidx, FALSE) > 0.5 )
1986 lexminfixes[i][j] = 1;
1987 else if ( SCIPvarGetUbAtIndex(vars[origrow][j], bdchgidx, FALSE) < 0.5 || j == n - 1 )
1988 lexminfixes[i][j] = 0;
1989 else
1990 lexminfixes[i][j] = 2;
1991 }
1992 }
1993
1994 /* find lexicographically minimal face of hypercube containing lexmin fixes */
1995 SCIP_CALL( findLexMinFace(vars, lexminfixes, minfixedrowlexmin, &terminate, m, n, nrowsused, TRUE) );
1996
1997 if ( terminate )
1998 goto FREELEXMIN;
1999
2000 /* Initialize lexicographically maximal matrix by fixed entries at the current node.
2001 * Free entries in the first column are set to 1.
2002 */
2003 SCIP_CALL( SCIPallocBufferArray(scip, &lexmaxfixes, nrowsused) );
2004 for (i = 0; i < nrowsused; ++i)
2005 {
2006 SCIP_CALL( SCIPallocBufferArray(scip, &lexmaxfixes[i], n) ); /*lint !e866*/
2007 }
2008
2009 /* store minimum fixed row for each column */
2010 SCIP_CALL( SCIPallocBufferArray(scip, &minfixedrowlexmax, n) );
2011 minfixedrowlexmax[0] = -1;
2012
2013 for (i = 0; i < nrowsused; ++i)
2014 {
2015 int origrow;
2016
2017 origrow = roworder[i];
2018
2019 for (j = 0; j < n; ++j)
2020 {
2021 if ( SCIPvarGetUbAtIndex(vars[origrow][j], bdchgidx, FALSE) < 0.5 )
2022 lexmaxfixes[i][j] = 0;
2023 else if ( SCIPvarGetLbAtIndex(vars[origrow][j], bdchgidx, FALSE) > 0.5 || j == 0 )
2024 lexmaxfixes[i][j] = 1;
2025 else
2026 lexmaxfixes[i][j] = 2;
2027 }
2028 }
2029
2030 /* find lexicographically maximal face of hypercube containing lexmax fixes */
2031 SCIP_CALL( findLexMaxFace(vars, lexmaxfixes, minfixedrowlexmax, &terminate, m, n, nrowsused, TRUE) );
2032
2033 if ( terminate )
2034 goto FREELEXMAX;
2035
2036 /* Find for each column j the minimal row in which lexminfixes and lexmaxfixes differ. Fix all entries above this
2037 * row to the corresponding value in lexminfixes (or lexmaxfixes).
2038 */
2039 for (j = 0; j < n; ++j)
2040 {
2041 int ub = MAX(minfixedrowlexmin[j], minfixedrowlexmax[j]);
2042
2043 for (i = 0; i <= ub; ++i)
2044 {
2045 int origrow;
2046
2047 origrow = roworder[i];
2048
2049 if ( SCIPvarGetLbAtIndex(vars[origrow][j], bdchgidx, FALSE) > 0.5 ||
2050 SCIPvarGetUbAtIndex(vars[origrow][j], bdchgidx, FALSE) < 0.5 )
2051 {
2052 SCIP_CALL( SCIPaddConflictBinvar(scip, vars[origrow][j]) );
2053 *result = SCIP_SUCCESS;
2054 }
2055 }
2056 }
2057
2058 FREELEXMAX:
2059 SCIPfreeBufferArray(scip, &minfixedrowlexmax);
2060 for (i = 0; i < nrowsused; ++i)
2061 SCIPfreeBufferArray(scip, &lexmaxfixes[i]);
2062 SCIPfreeBufferArray(scip, &lexmaxfixes);
2063
2064 FREELEXMIN:
2065 SCIPfreeBufferArray(scip, &minfixedrowlexmin);
2066 for (i = 0; i < nrowsused; ++i)
2067 SCIPfreeBufferArray(scip, &lexminfixes[i]);
2068 SCIPfreeBufferArray(scip, &lexminfixes);
2069
2070 return SCIP_OKAY;
2071 }
2072
2073
2074 /** Propagation conflict resolving method of propagator
2075 *
2076 * In this function we use that the propagation method above implicitly propagates SCIs, i.e., every
2077 * fixing can also be gotten via an SCI-fixing.
2078 *
2079 * Since the storage of an integer is not enough to store the complete information about the fixing
2080 * nor a complete shifted column, we have to use the linear time algorithm for SCIs.
2081 *
2082 * The inferinfo integer is set as follows:
2083 *
2084 * - If a shifted column is fixed to 0 and the corresponding bar does not necessarily has value 1
2085 * then we fix these entries to 0 and inferinfo is i * nblocks + j, where (i,j) is the leader of the
2086 * bar. The SCI depends on whether i is in Gamma or not (see Lemma 1 in the paper and the comments
2087 * above).
2088 *
2089 * - If a bar has value 1 and the shifted column has one entry that is not fixed, it can be fixed to
2090 * 1 and inferinfo is (nspcons*nblocks) + i * nblocks + j, where (i,j) is the leader of the bar; see
2091 * Proposition 1 (2c).
2092 */
2093 static
resolvePropagation(SCIP * scip,SCIP_CONS * cons,int inferinfo,SCIP_BDCHGIDX * bdchgidx,SCIP_RESULT * result)2094 SCIP_RETCODE resolvePropagation(
2095 SCIP* scip, /**< SCIP data structure */
2096 SCIP_CONS* cons, /**< constraint that inferred the bound change */
2097 int inferinfo, /**< inference information */
2098 SCIP_BDCHGIDX* bdchgidx, /**< bound change index (time stamp of bound change), or NULL for current time */
2099 SCIP_RESULT* result /**< pointer to store the result of the propagation conflict resolving call */
2100 )
2101 { /*lint --e{715}*/
2102 SCIP_CONSDATA* consdata;
2103 SCIP_Real** vals;
2104 SCIP_Real** weights;
2105 SCIP_VAR*** vars;
2106 SCIP_ORBITOPETYPE orbitopetype;
2107 int** cases;
2108
2109 int i;
2110 int j;
2111 int nspcons;
2112 int nblocks;
2113
2114 assert( scip != NULL );
2115 assert( cons != NULL );
2116 assert( result != NULL );
2117
2118 consdata = SCIPconsGetData(cons);
2119 assert( consdata != NULL );
2120 assert( consdata->nspcons > 0 );
2121 assert( consdata->nblocks > 0 );
2122 assert( consdata->vars != NULL );
2123 assert( consdata->vals != NULL );
2124 assert( consdata->weights != NULL );
2125 assert( consdata->cases != NULL );
2126 assert( consdata->istrianglefixed );
2127
2128 *result = SCIP_DIDNOTFIND;
2129 if ( ! consdata->resolveprop )
2130 return SCIP_OKAY;
2131
2132 nspcons = consdata->nspcons;
2133 nblocks = consdata->nblocks;
2134 vars = consdata->vars;
2135 vals = consdata->vals;
2136 weights = consdata->weights;
2137 orbitopetype = consdata->orbitopetype;
2138 cases = consdata->cases;
2139
2140 SCIPdebugMsg(scip, "Propagation resolution method of orbitope constraint using orbitopal fixing\n");
2141
2142 /* fill table */
2143 for (i = 0; i < nspcons; ++i)
2144 {
2145 int lastcolumn;
2146
2147 /* last column considered as part of the bar: */
2148 lastcolumn = nblocks - 1;
2149 if ( lastcolumn > i )
2150 lastcolumn = i;
2151 for (j = 0; j <= lastcolumn; ++j)
2152 {
2153 /* if the variable was fixed to zero at conflict time */
2154 if ( SCIPgetVarUbAtIndex(scip, vars[i][j], bdchgidx, FALSE) < 0.5 )
2155 vals[i][j] = 0.0;
2156 else
2157 {
2158 /* if the variable was fixed to one at conflict time */
2159 if ( SCIPgetVarLbAtIndex(scip, vars[i][j], bdchgidx, FALSE) > 0.5 )
2160 vals[i][j] = 2.0;
2161 else
2162 vals[i][j] = 1.0;
2163 }
2164 }
2165 }
2166
2167 #ifdef PRINT_MATRIX
2168 SCIPdebugMsg(scip, "Matrix:\n");
2169 printMatrix(scip, consdata);
2170 #endif
2171
2172 /* computation of table: this now minimizes the value of the shifted column */
2173 assert( consdata->istrianglefixed );
2174 computeSCTable(scip, nspcons, nblocks, weights, cases, vals);
2175
2176 /* if we fixed variables in the bar to zero */
2177 assert( inferinfo >= 0 && inferinfo < 2 * nspcons * nblocks );
2178 if ( inferinfo < nspcons * nblocks )
2179 {
2180 int p1;
2181 int p2;
2182 #ifdef SCIP_DEBUG
2183 char str[SCIP_MAXSTRLEN];
2184 char tmpstr[SCIP_MAXSTRLEN];
2185 #endif
2186
2187 i = (int) (inferinfo / nblocks);
2188 j = inferinfo % nblocks;
2189 assert( 0 <= i && i < nspcons );
2190 assert( 0 <= j && j < nblocks );
2191
2192 /* find SCI with value 0 */
2193 assert( weights[i-1][j-1] < 0.5 );
2194
2195 SCIPdebugMsg(scip, " -> reason for x[%d][%d] = ... = x[%d][%d] = 0 was the following SC:\n", i, j, i, MIN(i,nblocks));
2196 #ifdef SCIP_DEBUG
2197 str[0] = '\0';
2198 #endif
2199
2200 p1 = i-1;
2201 p2 = j-1;
2202 do
2203 {
2204 assert( cases[p1][p2] != -1 );
2205 assert( p1 >= 0 && p1 < i );
2206 assert( p2 >= 0 && p2 < j );
2207
2208 /* if case 1 */
2209 if ( cases[p1][p2] == 1 )
2210 --p2; /* decrease column */
2211 else
2212 {
2213 /* case 2 or 3: */
2214 assert( cases[p1][p2] == 2 || cases[p1][p2] == 3 );
2215 assert( SCIPgetVarUbAtIndex(scip, vars[p1][p2], bdchgidx, FALSE) < 0.5 );
2216 SCIP_CALL( SCIPaddConflictUb(scip, vars[p1][p2], bdchgidx) );
2217 *result = SCIP_SUCCESS;
2218
2219 #ifdef SCIP_DEBUG
2220 (void) SCIPsnprintf(tmpstr, SCIP_MAXSTRLEN, " (%d,%d)", p1, p2);
2221 (void) strncat(str, tmpstr, SCIP_MAXSTRLEN);
2222 #endif
2223
2224 if ( cases[p1][p2] == 3 )
2225 break;
2226 }
2227 --p1; /* decrease row */
2228 }
2229 while ( p1 >= 0 ); /* should always be true, i.e. the break should end the loop */
2230 assert( cases[p1][p2] == 3 );
2231
2232 #ifdef SCIP_DEBUG
2233 SCIPdebugMsg(scip, "%s\n", str);
2234 #endif
2235 }
2236 else
2237 {
2238 int k;
2239 int p1;
2240 int p2;
2241 #ifndef NDEBUG
2242 int pos1;
2243 int pos2;
2244 #endif
2245 #ifdef SCIP_DEBUG
2246 char str[SCIP_MAXSTRLEN];
2247 char tmpstr[SCIP_MAXSTRLEN];
2248 #endif
2249
2250 /* if we fixed a variable in the SC to 1 */
2251 inferinfo -= nspcons * nblocks;
2252 i = (int) inferinfo / nblocks;
2253 j = inferinfo % nblocks;
2254 assert( 0 <= i && i < nspcons );
2255 assert( 0 <= j && j < nblocks );
2256
2257 /* In rare cases it might happen that we fixed a variable to 1, but the node later becomes infeasible by globally
2258 * fixing variables to 0. In this case, it might happen that we find a SC with value 0 instead of 1. We then
2259 * cannot use this SC to repropagate (and do not know how to reconstruct the original reasoning). */
2260 if ( weights[i-1][j-1] > 0.5 && weights[i-1][j-1] < 1.5 )
2261 {
2262 SCIPdebugMsg(scip, " -> reason for x[%d][%d] = 1 was the following SC:\n", i, j);
2263 #ifdef SCIP_DEBUG
2264 (void) SCIPsnprintf(str, SCIP_MAXSTRLEN, "SC:");
2265 #endif
2266
2267 p1 = i-1;
2268 p2 = j-1;
2269 #ifndef NDEBUG
2270 pos1 = -1;
2271 pos2 = -1;
2272 #endif
2273 do
2274 {
2275 assert( cases[p1][p2] != -1 );
2276 assert( p1 >= 0 && p1 < i );
2277 assert( p2 >= 0 && p2 < j );
2278
2279 /* if case 1 */
2280 if ( cases[p1][p2] == 1 )
2281 --p2; /* decrease column */
2282 else
2283 {
2284 /* case 2 or 3: reason are formed by variables in SC fixed to 0 */
2285 assert( cases[p1][p2] == 2 || cases[p1][p2] == 3 );
2286 if ( SCIPgetVarUbAtIndex(scip, vars[p1][p2], bdchgidx, FALSE) < 0.5 )
2287 {
2288 SCIP_CALL( SCIPaddConflictUb(scip, vars[p1][p2], bdchgidx) );
2289 *result = SCIP_SUCCESS;
2290
2291 #ifdef SCIP_DEBUG
2292 (void) SCIPsnprintf(tmpstr, SCIP_MAXSTRLEN, " (%d,%d)", p1, p2);
2293 (void) strncat(str, tmpstr, SCIP_MAXSTRLEN);
2294 #endif
2295 }
2296 #ifndef NDEBUG
2297 else
2298 {
2299 assert( SCIPgetVarLbAtIndex(scip, vars[p1][p2], bdchgidx, FALSE) < 0.5 );
2300 assert( pos1 == -1 && pos2 == -1 );
2301 pos1 = p1;
2302 pos2 = p2;
2303 }
2304 #endif
2305 if ( cases[p1][p2] == 3 )
2306 break;
2307 }
2308 --p1; /* decrease row */
2309 }
2310 while ( p1 >= 0 ); /* should always be true, i.e., the break should end the loop */
2311 assert( cases[p1][p2] == 3 );
2312 assert( pos1 >= 0 && pos2 >= 0 );
2313
2314 /* distinguish partitioning/packing */
2315 if ( orbitopetype == SCIP_ORBITOPETYPE_PARTITIONING )
2316 {
2317 /* partitioning case */
2318 #ifdef SCIP_DEBUG
2319 (void) SCIPsnprintf(tmpstr, SCIP_MAXSTRLEN, " before bar: ");
2320 (void) strncat(str, tmpstr, SCIP_MAXSTRLEN);
2321 #endif
2322 /* add variables before the bar in the partitioning case */
2323 for (k = 0; k < j; ++k)
2324 {
2325 assert( SCIPgetVarUbAtIndex(scip, vars[i][k], bdchgidx, FALSE) < 0.5 );
2326 SCIP_CALL( SCIPaddConflictUb(scip, vars[i][k], bdchgidx) );
2327 *result = SCIP_SUCCESS;
2328 #ifdef SCIP_DEBUG
2329 (void) SCIPsnprintf(tmpstr, SCIP_MAXSTRLEN, " (%d,%d)", i, k);
2330 (void) strncat(str, tmpstr, SCIP_MAXSTRLEN);
2331 #endif
2332 }
2333
2334 #ifdef SCIP_DEBUG
2335 SCIPdebugMsg(scip, "%s\n", str);
2336 #endif
2337 }
2338 else
2339 {
2340 /* packing case */
2341 int lastcolumn;
2342
2343 /* last column considered as part of the bar: */
2344 lastcolumn = nblocks - 1;
2345 if ( lastcolumn > i )
2346 lastcolumn = i;
2347
2348 /* search for variable in the bar that is fixed to 1 in the packing case */
2349 for (k = j; k <= lastcolumn; ++k)
2350 {
2351 if ( SCIPgetVarLbAtIndex(scip, vars[i][k], bdchgidx, FALSE) > 0.5 )
2352 {
2353 SCIP_CALL( SCIPaddConflictLb(scip, vars[i][k], bdchgidx) );
2354 *result = SCIP_SUCCESS;
2355 SCIPdebugMsg(scip, " and variable x[%d][%d] fixed to 1.\n", i, k);
2356 break;
2357 }
2358 }
2359 }
2360 }
2361 }
2362
2363 return SCIP_OKAY;
2364 }
2365
2366
2367 /** check packing/partitioning orbitope solution for feasibility */
2368 static
enfopsPackingPartitioningOrbitopeSolution(SCIP * scip,SCIP_CONS * cons,SCIP_RESULT * result)2369 SCIP_RETCODE enfopsPackingPartitioningOrbitopeSolution(
2370 SCIP* scip, /**< SCIP data structure */
2371 SCIP_CONS* cons, /**< pointer to orbitope constraint */
2372 SCIP_RESULT* result /**< pointer to store the result of the enforcing call */
2373 )
2374 {
2375 SCIP_CONSDATA* consdata;
2376 SCIP_Real** weights;
2377 SCIP_Real** vals;
2378 int** cases;
2379 int nspcons;
2380 int nblocks;
2381 int i;
2382 int j;
2383
2384 assert( cons != NULL );
2385
2386 consdata = SCIPconsGetData(cons);
2387
2388 assert( scip != NULL );
2389 assert( consdata != NULL );
2390 assert( consdata->nspcons > 0 );
2391 assert( consdata->nblocks > 0 );
2392 assert( consdata->vals != NULL );
2393 assert( consdata->weights != NULL );
2394 assert( consdata->cases != NULL );
2395
2396 /* check for upper right triangle */
2397 if ( ! consdata->istrianglefixed )
2398 {
2399 SCIP_Bool infeasible = FALSE;
2400 int nfixedvars = 0;
2401
2402 SCIP_CALL( fixTriangle(scip, cons, &infeasible, &nfixedvars) );
2403 if ( infeasible )
2404 {
2405 *result = SCIP_CUTOFF;
2406 return SCIP_OKAY;
2407 }
2408 if ( nfixedvars > 0 )
2409 {
2410 *result = SCIP_REDUCEDDOM;
2411 return SCIP_OKAY;
2412 }
2413 }
2414
2415 nspcons = consdata->nspcons;
2416 nblocks = consdata->nblocks;
2417 vals = consdata->vals;
2418 weights = consdata->weights;
2419 cases = consdata->cases;
2420
2421 /* get solution */
2422 copyValues(scip, consdata, NULL);
2423 SCIPdebugMsg(scip, "Enforcing (pseudo solutions) for orbitope constraint <%s>\n", SCIPconsGetName(cons));
2424
2425 /* compute table */
2426 assert( consdata->istrianglefixed );
2427 computeSCTable(scip, nspcons, nblocks, weights, cases, vals);
2428
2429 /* loop through rows */
2430 for (i = 1; i < nspcons; ++i)
2431 {
2432 SCIP_Real bar = 0.0;
2433 int lastcolumn;
2434
2435 lastcolumn = nblocks - 1;
2436
2437 /* last column considered as part of the bar: */
2438 if ( lastcolumn > i )
2439 lastcolumn = i;
2440
2441 /* traverse row from right to left */
2442 for (j = lastcolumn; j > 0; --j)
2443 {
2444 bar += vals[i][j];
2445 assert( SCIPisIntegral(scip, vals[i][j]) );
2446
2447 /* check whether weights[i-1][j-1] < bar (<=> bar - weights[i-1][j-1] > 0), i.e. cut is violated) */
2448 if ( SCIPisGT(scip, bar - weights[i-1][j-1], 0.0) )
2449 {
2450 SCIPdebugMsg(scip, "Solution is infeasible.\n");
2451 *result = SCIP_INFEASIBLE;
2452 return SCIP_OKAY;
2453 }
2454 }
2455 }
2456
2457 return SCIP_OKAY;
2458 }
2459
2460
2461 /** check packing/partitioning orbitope solution for feasibility */
2462 static
checkPackingPartitioningOrbitopeSolution(SCIP * scip,SCIP_CONS * cons,SCIP_SOL * sol,SCIP_RESULT * result,SCIP_Bool printreason)2463 SCIP_RETCODE checkPackingPartitioningOrbitopeSolution(
2464 SCIP* scip, /**< SCIP data structure */
2465 SCIP_CONS* cons, /**< pointer to orbitope constraint */
2466 SCIP_SOL* sol, /**< solution to be checked */
2467 SCIP_RESULT* result, /**< pointer to store the result of the enforcing call */
2468 SCIP_Bool printreason /**< whether reason for infeasibility should be printed */
2469 )
2470 {
2471 SCIP_CONSDATA* consdata;
2472 SCIP_VAR*** vars;
2473 SCIP_Real** vals;
2474 SCIP_Real** weights;
2475 int** cases;
2476 int nspcons;
2477 int nblocks;
2478 int i;
2479 int j;
2480
2481 /* get data of constraint */
2482 assert( cons != 0 );
2483 consdata = SCIPconsGetData(cons);
2484
2485 assert( consdata != NULL );
2486 assert( consdata->nspcons > 0 );
2487 assert( consdata->nblocks > 0 );
2488 assert( consdata->vars != NULL );
2489 assert( consdata->vals != NULL );
2490 assert( consdata->weights != NULL );
2491 assert( consdata->cases != NULL );
2492
2493 nspcons = consdata->nspcons;
2494 nblocks = consdata->nblocks;
2495 vars = consdata->vars;
2496 vals = consdata->vals;
2497 weights = consdata->weights;
2498 cases = consdata->cases;
2499
2500 /* get solution */
2501 copyValues(scip, consdata, sol);
2502 SCIPdebugMsg(scip, "Checking orbitope constraint <%s> ...\n", SCIPconsGetName(cons));
2503
2504 /* check upper right triangle (if not yet fixed to zero or in debug mode */
2505 #ifdef NDEBUG
2506 if ( ! consdata->istrianglefixed )
2507 #endif
2508 {
2509 int diagsize;
2510
2511 /* get last row of triangle */
2512 diagsize = nblocks;
2513 if ( nspcons < nblocks )
2514 diagsize = nspcons;
2515
2516 /* check variables */
2517 for (i = 0; i < diagsize; ++i)
2518 {
2519 for (j = i+1; j < nblocks; ++j)
2520 {
2521 if ( ! SCIPisFeasZero(scip, vals[i][j]) )
2522 {
2523 if ( printreason )
2524 SCIPinfoMessage(scip, NULL, "variable x[%d][%d] = %f on upper right nonzero.\n", i, j, vals[i][j]);
2525 *result = SCIP_INFEASIBLE;
2526 }
2527 }
2528 }
2529 }
2530
2531 /* compute table */
2532 computeSCTable(scip, nspcons, nblocks, weights, cases, vals);
2533
2534 /* loop through rows */
2535 for (i = 1; i < nspcons; ++i)
2536 {
2537 SCIP_Real bar;
2538 int lastcolumn;
2539
2540 lastcolumn = nblocks - 1;
2541 bar = 0.0;
2542 /* last column considered as part of the bar: */
2543 if ( lastcolumn > i )
2544 lastcolumn = i;
2545
2546 /* traverse row from right to left */
2547 for (j = lastcolumn; j > 0; --j)
2548 {
2549 bar += vals[i][j];
2550 assert( SCIPisFeasIntegral(scip, vals[i][j]) );
2551
2552 /* check whether weights[i-1][j-1] < bar (<=> bar - weights[i-1][j-1] > 0), i.e. cut is violated) */
2553 if ( SCIPisGT(scip, bar - weights[i-1][j-1], 0.0) )
2554 {
2555 SCIPdebugMsg(scip, "Solution is infeasible.\n");
2556 *result = SCIP_INFEASIBLE;
2557
2558 if ( printreason )
2559 {
2560 int l;
2561 int p1;
2562 int p2;
2563
2564 SCIPinfoMessage(scip, NULL, "violated SCI: bar(");
2565
2566 /* first output bar */
2567 for (l = j; l < nblocks; ++l)
2568 SCIPinfoMessage(scip, NULL, "<%s> (%f)", SCIPvarGetName(vars[i][l]), consdata->vals[i][l]);
2569
2570 SCIPinfoMessage(scip, NULL, ") SC(");
2571
2572 /* output shifted column */
2573 p1 = i-1;
2574 p2 = j-1;
2575 do
2576 {
2577 assert( cases[p1][p2] != -1 );
2578 assert( p1 >= 0 && p1 < i );
2579 assert( p2 >= 0 && p2 < j );
2580
2581 /* if case 1 */
2582 if (cases[p1][p2] == 1)
2583 --p2; /* decrease column */
2584 else
2585 {
2586 /* case 2 or 3: */
2587 assert( cases[p1][p2] == 2 || cases[p1][p2] == 3 );
2588 SCIPinfoMessage(scip, NULL, "<%s> (%f)", SCIPvarGetName(vars[p1][p2]), consdata->vals[p1][p2]);
2589 if ( cases[p1][p2] == 3 )
2590 break;
2591 }
2592 --p1; /* decrease row */
2593 }
2594 while ( p1 >= 0 ); /* should always be true, i.e. the break should end the loop */
2595 assert( cases[p1][p2] == 3 );
2596
2597 SCIPinfoMessage(scip, NULL, ")");
2598 }
2599 }
2600 }
2601 }
2602
2603 return SCIP_OKAY;
2604 }
2605
2606
2607 /** check full orbitope solution for feasibility */
2608 static
checkFullOrbitopeSolution(SCIP * scip,SCIP_CONS * cons,SCIP_SOL * sol,SCIP_Bool printreason,SCIP_Bool * feasible)2609 SCIP_RETCODE checkFullOrbitopeSolution(
2610 SCIP* scip, /**< SCIP data structure */
2611 SCIP_CONS* cons, /**< constraint to process */
2612 SCIP_SOL* sol, /**< solution to be checked */
2613 SCIP_Bool printreason, /**< whether reason for infeasibility should be printed */
2614 SCIP_Bool* feasible /**< memory address to store whether solution is feasible */
2615 )
2616 {
2617 SCIP_CONSDATA* consdata;
2618 SCIP_VAR*** vars;
2619 SCIP_VAR** vars1;
2620 SCIP_VAR** vars2;
2621 int nrows;
2622 int ncols;
2623 int j;
2624 int i;
2625
2626 assert( scip != NULL );
2627 assert( cons != NULL );
2628 assert( feasible != NULL );
2629
2630 consdata = SCIPconsGetData(cons);
2631
2632 assert( consdata != NULL );
2633 assert( consdata->vars != NULL );
2634 assert( consdata->nspcons > 0 );
2635 assert( consdata->nblocks > 0 );
2636 assert( ! consdata->ismodelcons ); /* non-model constraints are never checked */
2637
2638 vars = consdata->vars;
2639 nrows = consdata->nspcons;
2640 ncols = consdata->nblocks;
2641
2642 SCIP_CALL( SCIPallocBufferArray(scip, &vars1, nrows) );
2643 SCIP_CALL( SCIPallocBufferArray(scip, &vars2, nrows) );
2644
2645 /* iterate over adjacent columns of orbitope and check whether the first column in this
2646 * column pair is lexicographically not smaller than the second column in the pair */
2647 *feasible = TRUE;
2648 for (j = 1; j < ncols && *feasible; ++j)
2649 {
2650 for (i = 0; i < nrows; ++i)
2651 {
2652 vars1[i] = vars[i][j - 1];
2653 vars2[i] = vars[i][j];
2654 }
2655
2656 SCIP_CALL( SCIPcheckSolutionOrbisack(scip, sol, vars1, vars2, nrows, printreason, feasible) );
2657 }
2658
2659 SCIPfreeBufferArray(scip, &vars2);
2660 SCIPfreeBufferArray(scip, &vars1);
2661
2662 return SCIP_OKAY;
2663 }
2664
2665
2666 /** separate orbisack cover inequalities */
2667 static
separateCoversOrbisack(SCIP * scip,SCIP_CONS * cons,SCIP_SOL * sol,SCIP_Bool dynamic,int * ngen,SCIP_Bool * infeasible)2668 SCIP_RETCODE separateCoversOrbisack(
2669 SCIP* scip, /**< SCIP data structure */
2670 SCIP_CONS* cons, /**< constraint to process */
2671 SCIP_SOL* sol, /**< solution to separate (NULL for the LP solution) */
2672 SCIP_Bool dynamic, /**< whether we use a dynamic row order */
2673 int* ngen, /**< pointer to store number of generated cuts */
2674 SCIP_Bool* infeasible /**< pointer to store whether infeasibility has been detected */
2675 )
2676 {
2677 SCIP_CONSDATA* consdata;
2678 SCIP_VAR*** vars;
2679 int* roworder;
2680 int nrowsused;
2681 int nrows;
2682 int ncols;
2683 int i;
2684 int j;
2685 int origrow;
2686 SCIP_Real rhs;
2687 SCIP_Real lhs;
2688 SCIP_Real* coeffs1;
2689 SCIP_Real* coeffs2;
2690
2691 assert( scip != NULL );
2692 assert( cons != NULL );
2693 assert( ngen != NULL );
2694 assert( infeasible != NULL );
2695
2696 *ngen = 0;
2697 *infeasible = FALSE;
2698
2699 /* get basic data */
2700 consdata = SCIPconsGetData(cons);
2701 assert( consdata != NULL );
2702
2703 vars = consdata->vars;
2704 nrows = consdata->nspcons;
2705 ncols = consdata->nblocks;
2706 nrowsused = dynamic ? consdata->nrowsused : nrows;
2707 roworder = consdata->roworder;
2708
2709 /* allocate memory for cover inequalities */
2710 SCIP_CALL( SCIPallocBufferArray(scip, &coeffs1, nrowsused) );
2711 SCIP_CALL( SCIPallocBufferArray(scip, &coeffs2, nrowsused) );
2712
2713 lhs = 0.0;
2714 rhs = 0.0;
2715
2716 /* separate orbisack cover inequalities for adjacent columns */
2717 for (j = 0; j < ncols - 1 && ! *infeasible; ++j)
2718 {
2719 SCIP_Real rowval;
2720
2721 for (i = 0; i < nrowsused; ++i)
2722 {
2723 origrow = roworder[i];
2724
2725 assert( origrow >= 0 );
2726 assert( origrow < nrows );
2727
2728 rowval = SCIPgetSolVal(scip, sol, vars[origrow][j + 1]) - SCIPgetSolVal(scip, sol, vars[origrow][j]);
2729
2730 /* check whether cover inequality is violated */
2731 if ( SCIPisEfficacious(scip, rowval + lhs - rhs) )
2732 {
2733 SCIP_ROW* row;
2734 int k;
2735
2736 /* set coefficients for current inequality */
2737 coeffs1[i] = -1.0;
2738 coeffs2[i] = 1.0;
2739
2740 /* add violated orbisack cover inequality */
2741 SCIP_CALL( SCIPcreateEmptyRowCons(scip, &row, cons, "orbisackcover", -SCIPinfinity(scip), rhs, FALSE, FALSE, TRUE) );
2742 SCIP_CALL( SCIPcacheRowExtensions(scip, row) );
2743
2744 for (k = 0; k <= i; ++k)
2745 {
2746 int origrow2;
2747
2748 origrow2 = roworder[k];
2749
2750 SCIP_CALL( SCIPaddVarToRow(scip, row, vars[origrow2][j], coeffs1[k]) );
2751 SCIP_CALL( SCIPaddVarToRow(scip, row, vars[origrow2][j + 1], coeffs2[k]) );
2752 }
2753 SCIP_CALL( SCIPflushRowExtensions(scip, row) );
2754
2755 SCIP_CALL( SCIPaddRow(scip, row, FALSE, infeasible) );
2756 #ifdef SCIP_DEBUG
2757 SCIP_CALL( SCIPprintRow(scip, row, NULL) );
2758 #endif
2759 SCIP_CALL( SCIPreleaseRow(scip, &row) );
2760
2761 *ngen += 1;
2762 if ( *infeasible )
2763 break;
2764
2765 /* reset coefficients for next inequality */
2766 coeffs1[i] = 0.0;
2767 coeffs2[i] = 0.0;
2768 }
2769
2770 /* add argmax( 1 - vals[i][0], vals[i][1] ) as coefficient and ensure that both vars1[0] and vars2[0] are
2771 * contained in the LIFTED cover inequality */
2772 rowval = SCIPgetSolVal(scip, sol, vars[origrow][j]) + SCIPgetSolVal(scip, sol, vars[origrow][j + 1]);
2773 if ( SCIPisEfficacious(scip, 1.0 - rowval) )
2774 {
2775 coeffs1[i] = -1.0;
2776 coeffs2[i] = 0.0;
2777 lhs -= SCIPgetSolVal(scip, sol, vars[origrow][j]);
2778
2779 /* apply lifting? */
2780 if ( i == 0 )
2781 {
2782 coeffs2[i] = 1.0;
2783 lhs += SCIPgetSolVal(scip, sol, vars[origrow][j + 1]);
2784 }
2785 }
2786 else
2787 {
2788 coeffs1[i] = 0.0;
2789 coeffs2[i] = 1.0;
2790 lhs += SCIPgetSolVal(scip, sol, vars[origrow][j]);
2791 rhs += 1.0;
2792
2793 /* apply lifting? */
2794 if ( i == 0 )
2795 {
2796 coeffs1[i] = -1.0;
2797 lhs -= SCIPgetSolVal(scip, sol, vars[origrow][j]);
2798 rhs -= 1.0;
2799 }
2800 }
2801 }
2802 }
2803
2804 SCIPfreeBufferArray(scip, &coeffs1);
2805 SCIPfreeBufferArray(scip, &coeffs2);
2806
2807 return SCIP_OKAY;
2808 }
2809
2810
2811 /** separate or enforce constraints */
2812 static
separateConstraints(SCIP * scip,SCIP_CONSHDLR * conshdlr,SCIP_CONS ** conss,int nconss,int nusefulconss,SCIP_SOL * sol,SCIP_RESULT * result,SCIP_Bool enforce)2813 SCIP_RETCODE separateConstraints(
2814 SCIP* scip, /**< SCIP data structure */
2815 SCIP_CONSHDLR* conshdlr, /**< constraint handler */
2816 SCIP_CONS** conss, /**< constraints to process */
2817 int nconss, /**< number of constraints */
2818 int nusefulconss, /**< number of useful (non-obsolete) constraints to process */
2819 SCIP_SOL* sol, /**< solution to separate (NULL for the LP solution) */
2820 SCIP_RESULT* result, /**< pointer to store the result (should be initialized) */
2821 SCIP_Bool enforce /**< whether we enforce orbitope constraints */
2822 )
2823 {
2824 SCIP_Bool infeasible = FALSE;
2825 int nfixedvars = 0;
2826 int ncuts = 0;
2827 int c;
2828
2829 assert( scip != NULL );
2830 assert( conshdlr != NULL );
2831 assert( strcmp(SCIPconshdlrGetName(conshdlr), CONSHDLR_NAME) == 0 );
2832 assert( result != NULL );
2833
2834 /* loop through constraints */
2835 for (c = 0; c < nconss && ! infeasible; c++)
2836 {
2837 SCIP_CONSHDLRDATA* conshdlrdata;
2838 SCIP_CONSDATA* consdata;
2839 int nconsfixedvars = 0;
2840 int nconscuts = 0;
2841 SCIP_ORBITOPETYPE orbitopetype;
2842
2843 assert( conss[c] != NULL );
2844
2845 /* get data of constraint */
2846 consdata = SCIPconsGetData(conss[c]);
2847 assert( consdata != NULL );
2848
2849 /* do not enforce non-model constraints */
2850 if ( enforce && !consdata->ismodelcons )
2851 continue;
2852
2853 /* get solution */
2854 copyValues(scip, consdata, sol);
2855
2856 /* separate */
2857 orbitopetype = consdata->orbitopetype;
2858
2859 conshdlrdata = SCIPconshdlrGetData(conshdlr);
2860 if ( orbitopetype == SCIP_ORBITOPETYPE_PACKING || orbitopetype == SCIP_ORBITOPETYPE_PARTITIONING )
2861 {
2862 SCIP_CALL( separateSCIs(scip, conshdlr, conss[c], consdata, &infeasible, &nconsfixedvars, &nconscuts) );
2863 nfixedvars += nconsfixedvars;
2864 }
2865 else if ( conshdlrdata->sepafullorbitope )
2866 {
2867 SCIP_CALL( separateCoversOrbisack(scip, conss[c], sol, conshdlrdata->usedynamicprop && !consdata->ismodelcons, &nconscuts, &infeasible) );
2868 }
2869 ncuts += nconscuts;
2870
2871 /* stop after the useful constraints if we found cuts of fixed variables */
2872 if ( c >= nusefulconss && (ncuts > 0 || nfixedvars > 0) )
2873 break;
2874 }
2875
2876 if ( infeasible )
2877 {
2878 SCIPdebugMsg(scip, "Infeasible node.\n");
2879 *result = SCIP_CUTOFF;
2880 }
2881 else if ( nfixedvars > 0 )
2882 {
2883 SCIPdebugMsg(scip, "Fixed %d variables.\n", nfixedvars);
2884 *result = SCIP_REDUCEDDOM;
2885 }
2886 else if ( ncuts > 0 )
2887 {
2888 SCIPdebugMsg(scip, "Separated %dinequalities.\n", ncuts);
2889 *result = SCIP_SEPARATED;
2890 }
2891 else
2892 {
2893 SCIPdebugMsg(scip, "No violated inequality found during separation.\n");
2894 }
2895
2896 return SCIP_OKAY;
2897 }
2898
2899
2900 /** check whether all variables in an orbitope constraint are fixed */
2901 static
checkRedundantCons(SCIP * scip,SCIP_CONS * cons,SCIP_Bool * redundant)2902 SCIP_RETCODE checkRedundantCons(
2903 SCIP* scip, /**< SCIP data structure */
2904 SCIP_CONS* cons, /**< constraint to be processed */
2905 SCIP_Bool* redundant /**< pointer to store whether constraint is redundant (contains no active vars) */
2906 )
2907 {
2908 SCIP_CONSDATA* consdata;
2909 SCIP_VAR*** vars;
2910 int i;
2911 int j;
2912 int nrows;
2913 int ncols;
2914
2915 assert( scip != NULL );
2916 assert( cons != NULL );
2917 assert( redundant != NULL );
2918
2919 *redundant = FALSE;
2920
2921 consdata = SCIPconsGetData(cons);
2922 assert( consdata != NULL );
2923 assert( consdata->vars != NULL );
2924 assert( consdata->nspcons > 0 );
2925 assert( consdata->nblocks > 0 );
2926
2927 vars = consdata->vars;
2928 nrows = consdata->nspcons;
2929 ncols = consdata->nblocks;
2930
2931 /* check whether there exists an active variable in the orbitope */
2932 for (i = 0; i < nrows; ++i)
2933 {
2934 for (j = 0; j < ncols; ++j)
2935 {
2936 if ( SCIPvarIsActive(vars[i][j]) )
2937 return SCIP_OKAY;
2938 }
2939 }
2940
2941 *redundant = TRUE;
2942
2943 return SCIP_OKAY;
2944 }
2945
2946
2947 /*
2948 * Callback methods of constraint handler
2949 */
2950
2951 /** copy method for constraint handler plugins (called when SCIP copies plugins) */
2952 static
SCIP_DECL_CONSHDLRCOPY(conshdlrCopyOrbitope)2953 SCIP_DECL_CONSHDLRCOPY(conshdlrCopyOrbitope)
2954 { /*lint --e{715}*/
2955 assert(scip != NULL);
2956 assert(conshdlr != NULL);
2957 assert(strcmp(SCIPconshdlrGetName(conshdlr), CONSHDLR_NAME) == 0);
2958
2959 /* call inclusion method of constraint handler */
2960 SCIP_CALL( SCIPincludeConshdlrOrbitope(scip) );
2961
2962 *valid = TRUE;
2963
2964 return SCIP_OKAY;
2965 }
2966
2967 /** frees constraint handler */
2968 static
SCIP_DECL_CONSFREE(consFreeOrbitope)2969 SCIP_DECL_CONSFREE(consFreeOrbitope)
2970 { /*lint --e{715}*/
2971 SCIP_CONSHDLRDATA* conshdlrdata;
2972
2973 assert( scip != 0 );
2974 assert( conshdlr != 0 );
2975 assert( strcmp(SCIPconshdlrGetName(conshdlr), CONSHDLR_NAME) == 0 );
2976
2977 conshdlrdata = SCIPconshdlrGetData(conshdlr);
2978 assert( conshdlrdata != NULL );
2979
2980 SCIPfreeBlockMemory(scip, &conshdlrdata);
2981
2982 return SCIP_OKAY;
2983 }
2984
2985 /** frees specific constraint data */
2986 static
SCIP_DECL_CONSDELETE(consDeleteOrbitope)2987 SCIP_DECL_CONSDELETE(consDeleteOrbitope)
2988 { /*lint --e{715}*/
2989 SCIP_CONSHDLRDATA* conshdlrdata;
2990
2991 assert(conshdlr != NULL);
2992 assert(strcmp(SCIPconshdlrGetName(conshdlr), CONSHDLR_NAME) == 0);
2993
2994 conshdlrdata = SCIPconshdlrGetData(conshdlr);
2995 assert(conshdlrdata != NULL);
2996
2997 SCIP_CALL( consdataFree(scip, consdata, conshdlrdata->usedynamicprop) );
2998
2999 return SCIP_OKAY;
3000 }
3001
3002 /** transforms constraint data into data belonging to the transformed problem */
3003 static
SCIP_DECL_CONSTRANS(consTransOrbitope)3004 SCIP_DECL_CONSTRANS(consTransOrbitope)
3005 { /*lint --e{715}*/
3006 SCIP_CONSHDLRDATA* conshdlrdata;
3007 SCIP_CONSDATA* sourcedata;
3008 SCIP_CONSDATA* targetdata;
3009
3010 assert(conshdlr != NULL);
3011 assert(strcmp(SCIPconshdlrGetName(conshdlr), CONSHDLR_NAME) == 0);
3012 assert(SCIPgetStage(scip) == SCIP_STAGE_TRANSFORMING);
3013 assert(sourcecons != NULL);
3014 assert(targetcons != NULL);
3015
3016 conshdlrdata = SCIPconshdlrGetData(conshdlr);
3017 assert( conshdlrdata != NULL );
3018
3019 sourcedata = SCIPconsGetData(sourcecons);
3020 assert(sourcedata != NULL);
3021
3022 /* create linear constraint data for target constraint */
3023 SCIP_CALL( consdataCreate(scip, &targetdata, sourcedata->vars, sourcedata->nspcons, sourcedata->nblocks,
3024 sourcedata->orbitopetype, sourcedata->resolveprop, conshdlrdata->usedynamicprop, sourcedata->ismodelcons) );
3025
3026 /* create target constraint */
3027 SCIP_CALL( SCIPcreateCons(scip, targetcons, SCIPconsGetName(sourcecons), conshdlr, targetdata,
3028 SCIPconsIsInitial(sourcecons), SCIPconsIsSeparated(sourcecons), SCIPconsIsEnforced(sourcecons),
3029 SCIPconsIsChecked(sourcecons), SCIPconsIsPropagated(sourcecons),
3030 SCIPconsIsLocal(sourcecons), SCIPconsIsModifiable(sourcecons),
3031 SCIPconsIsDynamic(sourcecons), SCIPconsIsRemovable(sourcecons), SCIPconsIsStickingAtNode(sourcecons)) );
3032
3033 return SCIP_OKAY;
3034 }
3035
3036 /** separation method of constraint handler for LP solutions */
3037 static
SCIP_DECL_CONSSEPALP(consSepalpOrbitope)3038 SCIP_DECL_CONSSEPALP(consSepalpOrbitope)
3039 { /*lint --e{715}*/
3040 assert( scip != NULL );
3041 assert( result != NULL );
3042
3043 SCIPdebugMsg(scip, "Separation of orbitope constraint handler <%s> for LP solution.\n", SCIPconshdlrGetName(conshdlr));
3044
3045 *result = SCIP_DIDNOTRUN;
3046
3047 /* if solution is integer, skip separation */
3048 if ( SCIPgetNLPBranchCands(scip) <= 0 )
3049 return SCIP_OKAY;
3050
3051 *result = SCIP_DIDNOTFIND;
3052
3053 /* separate constraints */
3054 SCIP_CALL( separateConstraints(scip, conshdlr, conss, nconss, nusefulconss, NULL, result, FALSE) );
3055
3056 return SCIP_OKAY;
3057 }
3058
3059 /** separation method of constraint handler for arbitrary primal solutions */
3060 static
SCIP_DECL_CONSSEPASOL(consSepasolOrbitope)3061 SCIP_DECL_CONSSEPASOL(consSepasolOrbitope)
3062 { /*lint --e{715}*/
3063 assert( scip != NULL );
3064 assert( result != NULL );
3065
3066 SCIPdebugMsg(scip, "Separation of orbitope constraint handler <%s> for primal solution.\n", SCIPconshdlrGetName(conshdlr));
3067
3068 *result = SCIP_DIDNOTFIND;
3069
3070 /* separate constraints */
3071 SCIP_CALL( separateConstraints(scip, conshdlr, conss, nconss, nusefulconss, sol, result, FALSE) );
3072
3073 return SCIP_OKAY;
3074 }
3075
3076
3077 /** constraint enforcing method of constraint handler for LP solutions */
3078 static
SCIP_DECL_CONSENFOLP(consEnfolpOrbitope)3079 SCIP_DECL_CONSENFOLP(consEnfolpOrbitope)
3080 { /*lint --e{715}*/
3081 assert( scip != NULL );
3082 assert( result != NULL );
3083
3084 /* we have a negative priority, so we should come after the integrality conshdlr */
3085 assert( SCIPgetNLPBranchCands(scip) == 0 );
3086
3087 SCIPdebugMsg(scip, "Enforcement for orbitope constraint handler <%s> for LP solution.\n", SCIPconshdlrGetName(conshdlr));
3088
3089 *result = SCIP_FEASIBLE;
3090
3091 /* separate constraints */
3092 SCIP_CALL( separateConstraints(scip, conshdlr, conss, nconss, nusefulconss, NULL, result, TRUE) );
3093
3094 return SCIP_OKAY;
3095 }
3096
3097
3098 /** constraint enforcing method of constraint handler for relaxation solutions */
3099 static
SCIP_DECL_CONSENFORELAX(consEnforelaxOrbitope)3100 SCIP_DECL_CONSENFORELAX(consEnforelaxOrbitope)
3101 { /*lint --e{715}*/
3102 assert( result != NULL );
3103 assert( scip != NULL );
3104
3105 SCIPdebugMsg(scip, "Enforcement for orbitope constraint handler <%s> for relaxation solution.\n", SCIPconshdlrGetName(conshdlr));
3106
3107 *result = SCIP_FEASIBLE;
3108
3109 /* separate constraints */
3110 SCIP_CALL( separateConstraints(scip, conshdlr, conss, nconss, nusefulconss, sol, result, TRUE) );
3111
3112 return SCIP_OKAY;
3113 }
3114
3115
3116 /** constraint enforcing method of constraint handler for pseudo solutions */
3117 static
SCIP_DECL_CONSENFOPS(consEnfopsOrbitope)3118 SCIP_DECL_CONSENFOPS(consEnfopsOrbitope)
3119 { /*lint --e{715}*/
3120 int c;
3121
3122 assert( scip != NULL );
3123 assert( conshdlr != NULL );
3124 assert( strcmp(SCIPconshdlrGetName(conshdlr), CONSHDLR_NAME) == 0 );
3125 assert( result != NULL );
3126
3127 *result = SCIP_FEASIBLE;
3128 if ( objinfeasible || solinfeasible )
3129 return SCIP_OKAY;
3130
3131 /* loop through constraints */
3132 for (c = 0; c < nconss; ++c)
3133 {
3134 SCIP_CONS* cons;
3135 SCIP_CONSDATA* consdata;
3136 SCIP_ORBITOPETYPE orbitopetype;
3137 SCIP_Bool feasible;
3138
3139 /* get data of constraint */
3140 cons = conss[c];
3141 assert( cons != 0 );
3142 consdata = SCIPconsGetData(cons);
3143
3144 assert( consdata != NULL );
3145
3146 /* do not enforce non-model constraints */
3147 if ( !consdata->ismodelcons )
3148 continue;
3149
3150 orbitopetype = consdata->orbitopetype;
3151
3152 if ( orbitopetype == SCIP_ORBITOPETYPE_PACKING || orbitopetype == SCIP_ORBITOPETYPE_PARTITIONING )
3153 {
3154 SCIP_CALL( enfopsPackingPartitioningOrbitopeSolution(scip, cons, result) );
3155 }
3156 else
3157 {
3158 SCIP_CALL( checkFullOrbitopeSolution(scip, cons, NULL, FALSE, &feasible) );
3159
3160 if ( ! feasible )
3161 *result = SCIP_INFEASIBLE;
3162 }
3163
3164 if ( *result == SCIP_INFEASIBLE )
3165 break;
3166 }
3167
3168 return SCIP_OKAY;
3169 }
3170
3171
3172 /** feasibility check method of constraint handler for integral solutions */
3173 static
SCIP_DECL_CONSCHECK(consCheckOrbitope)3174 SCIP_DECL_CONSCHECK(consCheckOrbitope)
3175 { /*lint --e{715}*/
3176 int c;
3177 SCIP_CONSDATA* consdata;
3178 SCIP_ORBITOPETYPE orbitopetype;
3179 SCIP_Bool feasible;
3180
3181 assert( scip != NULL );
3182 assert( conshdlr != NULL );
3183 assert( strcmp(SCIPconshdlrGetName(conshdlr), CONSHDLR_NAME) == 0 );
3184 assert( result != NULL );
3185
3186 *result = SCIP_FEASIBLE;
3187
3188 /* loop through constraints */
3189 for( c = 0; c < nconss && (*result == SCIP_FEASIBLE || completely); ++c )
3190 {
3191 assert( conss[c] != 0 );
3192 consdata = SCIPconsGetData(conss[c]);
3193
3194 assert( consdata != NULL );
3195
3196 /* do not check non-model constraints */
3197 if ( !consdata->ismodelcons )
3198 continue;
3199
3200 orbitopetype = consdata->orbitopetype;
3201
3202 if ( orbitopetype == SCIP_ORBITOPETYPE_PACKING || orbitopetype == SCIP_ORBITOPETYPE_PARTITIONING )
3203 {
3204 SCIP_CALL( checkPackingPartitioningOrbitopeSolution(scip, conss[c], sol, result, printreason) );
3205 }
3206 else
3207 {
3208 SCIP_CALL( checkFullOrbitopeSolution(scip, conss[c], sol, printreason, &feasible) );
3209
3210 if ( ! feasible )
3211 *result = SCIP_INFEASIBLE;
3212 }
3213 }
3214 SCIPdebugMsg(scip, "Solution is feasible.\n");
3215
3216 return SCIP_OKAY;
3217 }
3218
3219
3220 /** domain propagation method of constraint handler */
3221 static
SCIP_DECL_CONSPROP(consPropOrbitope)3222 SCIP_DECL_CONSPROP(consPropOrbitope)
3223 { /*lint --e{715}*/
3224 SCIP_CONSHDLRDATA* conshdlrdata;
3225 SCIP_Bool infeasible = FALSE;
3226 int nfixedvars = 0;
3227 int c;
3228
3229 assert( scip != NULL );
3230 assert( conshdlr != NULL );
3231 assert( strcmp(SCIPconshdlrGetName(conshdlr), CONSHDLR_NAME) == 0 );
3232 assert( result != NULL );
3233
3234 *result = SCIP_DIDNOTRUN;
3235
3236 conshdlrdata = SCIPconshdlrGetData(conshdlr);
3237 assert( conshdlrdata != NULL );
3238
3239 /* propagate all useful constraints */
3240 for (c = 0; c < nusefulconss && !infeasible; ++c)
3241 {
3242 assert( conss[c] != 0 );
3243
3244 SCIPdebugMsg(scip, "Propagation of orbitope constraint <%s> ...\n", SCIPconsGetName(conss[c]));
3245
3246 SCIP_CALL( propagateCons(scip, conss[c], &infeasible, &nfixedvars, conshdlrdata->usedynamicprop) );
3247 }
3248
3249 /* return the correct result */
3250 if ( infeasible )
3251 {
3252 *result = SCIP_CUTOFF;
3253 SCIPdebugMsg(scip, "Propagation via orbitopal fixing proved node to be infeasible.\n");
3254 }
3255 else if ( nfixedvars > 0 )
3256 {
3257 *result = SCIP_REDUCEDDOM;
3258 SCIPdebugMsg(scip, "Propagated %d variables via orbitopal fixing.\n", nfixedvars);
3259 }
3260 else if ( nusefulconss > 0 )
3261 {
3262 *result = SCIP_DIDNOTFIND;
3263 SCIPdebugMsg(scip, "Propagation via orbitopal fixing did not find anything.\n");
3264 }
3265
3266 return SCIP_OKAY;
3267 }
3268
3269
3270 /** presolving method of constraint handler */
3271 static
SCIP_DECL_CONSPRESOL(consPresolOrbitope)3272 SCIP_DECL_CONSPRESOL(consPresolOrbitope)
3273 { /*lint --e{715}*/
3274 SCIP_CONSHDLRDATA* conshdlrdata;
3275 SCIP_Bool infeasible = FALSE;
3276 int noldfixedvars;
3277 int c;
3278 SCIP_Bool redundant;
3279
3280 assert( scip != NULL );
3281 assert( conshdlr != NULL );
3282 assert( strcmp(SCIPconshdlrGetName(conshdlr), CONSHDLR_NAME) == 0 );
3283 assert( result != NULL );
3284
3285 *result = SCIP_DIDNOTRUN;
3286 noldfixedvars = *nfixedvars;
3287
3288 conshdlrdata = SCIPconshdlrGetData(conshdlr);
3289 assert( conshdlrdata != NULL );
3290
3291 /* propagate all useful constraints
3292 *
3293 * @todo use an event handler to only propagate if a variable in the orbitope has been fixed
3294 */
3295 for (c = 0; c < nconss && !infeasible; ++c)
3296 {
3297 int nfixed = 0;
3298
3299 assert( conss[c] != 0 );
3300
3301 SCIPdebugMsg(scip, "Presolving of orbitope constraint <%s> ...\n", SCIPconsGetName(conss[c]));
3302
3303 /* first propagate */
3304 SCIP_CALL( propagateCons(scip, conss[c], &infeasible, &nfixed, conshdlrdata->usedynamicprop) );
3305 *nfixedvars += nfixed;
3306
3307 if ( ! infeasible )
3308 {
3309 SCIP_CALL( checkRedundantCons(scip, conss[c], &redundant) );
3310
3311 if ( redundant )
3312 {
3313 SCIPdebugMsg(scip, "Orbitope constraint <%s> is redundant: it does not contain active variables\n",
3314 SCIPconsGetName(conss[c]));
3315 SCIP_CALL( SCIPdelCons(scip, conss[c]) );
3316 assert( ! SCIPconsIsActive(conss[c]) );
3317 (*ndelconss)++;
3318 continue;
3319 }
3320 }
3321 }
3322
3323 if ( infeasible )
3324 {
3325 *result = SCIP_CUTOFF;
3326 SCIPdebugMsg(scip, "Presolving detected infeasibility.\n");
3327 }
3328 else if ( *nfixedvars > noldfixedvars )
3329 {
3330 *result = SCIP_SUCCESS;
3331 }
3332 else if ( nconss > 0 )
3333 {
3334 *result = SCIP_DIDNOTFIND;
3335 SCIPdebugMsg(scip, "Presolving via orbitopal fixing did not find anything.\n");
3336 }
3337
3338 return SCIP_OKAY;
3339 }
3340
3341
3342 /** propagation conflict resolving method of constraint handler */
3343 static
SCIP_DECL_CONSRESPROP(consRespropOrbitope)3344 SCIP_DECL_CONSRESPROP(consRespropOrbitope)
3345 { /*lint --e{715}*/
3346 SCIP_CONSDATA* consdata;
3347 SCIP_ORBITOPETYPE orbitopetype;
3348
3349 assert( scip != NULL );
3350 assert( cons != NULL );
3351 assert( infervar != NULL );
3352 assert( bdchgidx != NULL );
3353 assert( result != NULL );
3354
3355 consdata = SCIPconsGetData(cons);
3356 assert( consdata != NULL );
3357
3358 orbitopetype = consdata->orbitopetype;
3359
3360 /* resolution for full orbitopes not availabe yet */
3361 if ( orbitopetype == SCIP_ORBITOPETYPE_PACKING || orbitopetype == SCIP_ORBITOPETYPE_PARTITIONING )
3362 {
3363 SCIP_CALL( resolvePropagation(scip, cons, inferinfo, bdchgidx, result) );
3364 }
3365 else
3366 {
3367 SCIP_CALL( resolvePropagationFullOrbitope(scip, conshdlr, cons, inferinfo, bdchgidx, result) );
3368 }
3369
3370 return SCIP_OKAY;
3371 }
3372
3373
3374 /** variable rounding lock method of constraint handler */
3375 static
SCIP_DECL_CONSLOCK(consLockOrbitope)3376 SCIP_DECL_CONSLOCK(consLockOrbitope)
3377 { /*lint --e{715}*/
3378 SCIP_CONSDATA* consdata;
3379 SCIP_VAR*** vars;
3380 int i;
3381 int j;
3382 int nspcons;
3383 int nblocks;
3384
3385 assert( scip != NULL );
3386 assert( conshdlr != NULL );
3387 assert( strcmp(SCIPconshdlrGetName(conshdlr), CONSHDLR_NAME) == 0 );
3388 assert( cons != NULL );
3389 assert( locktype == SCIP_LOCKTYPE_MODEL );
3390
3391 consdata = SCIPconsGetData(cons);
3392 assert( consdata != NULL );
3393 assert( consdata->nspcons > 0 );
3394 assert( consdata->nblocks > 0 );
3395 assert( consdata->vars != NULL );
3396
3397 SCIPdebugMsg(scip, "Locking method for orbitope constraint handler\n");
3398
3399 nspcons = consdata->nspcons;
3400 nblocks = consdata->nblocks;
3401 vars = consdata->vars;
3402
3403 /* add up locks and down locks on each variable */
3404 for (i = 0; i < nspcons; ++i)
3405 {
3406 for (j = 0; j < nblocks; ++j)
3407 SCIP_CALL( SCIPaddVarLocksType(scip, vars[i][j], locktype, nlockspos + nlocksneg, nlockspos + nlocksneg) );
3408 }
3409
3410 return SCIP_OKAY;
3411 }
3412
3413
3414 /** constraint display method of constraint handler */
3415 static
SCIP_DECL_CONSPRINT(consPrintOrbitope)3416 SCIP_DECL_CONSPRINT(consPrintOrbitope)
3417 {
3418 SCIP_CONSDATA* consdata;
3419 SCIP_VAR*** vars;
3420 int i;
3421 int j;
3422 int nspcons;
3423 int nblocks;
3424 SCIP_ORBITOPETYPE orbitopetype;
3425
3426 assert( scip != NULL );
3427 assert( conshdlr != NULL );
3428 assert( strcmp(SCIPconshdlrGetName(conshdlr), CONSHDLR_NAME) == 0 );
3429 assert( cons != NULL );
3430
3431 consdata = SCIPconsGetData(cons);
3432 assert( consdata != NULL );
3433 assert( consdata->nspcons > 0 );
3434 assert( consdata->nblocks > 0 );
3435 assert( consdata->vars != NULL );
3436
3437 nspcons = consdata->nspcons;
3438 nblocks = consdata->nblocks;
3439 vars = consdata->vars;
3440 orbitopetype = consdata->orbitopetype;
3441
3442 SCIPdebugMsg(scip, "Printing method for orbitope constraint handler\n");
3443
3444 switch ( orbitopetype )
3445 {
3446 case SCIP_ORBITOPETYPE_PARTITIONING:
3447 SCIPinfoMessage(scip, file, "partOrbitope(");
3448 break;
3449 case SCIP_ORBITOPETYPE_PACKING:
3450 SCIPinfoMessage(scip, file, "packOrbitope(");
3451 break;
3452 case SCIP_ORBITOPETYPE_FULL:
3453 SCIPinfoMessage(scip, file, "fullOrbitope(");
3454 break;
3455 default:
3456 SCIPABORT();
3457 }
3458
3459 for (i = 0; i < nspcons; ++i)
3460 {
3461 for (j = 0; j < nblocks; ++j)
3462 {
3463 if ( j > 0 )
3464 SCIPinfoMessage(scip, file, ",");
3465 SCIP_CALL( SCIPwriteVarName(scip, file, vars[i][j], TRUE) );
3466 }
3467 if ( i < nspcons-1 )
3468 SCIPinfoMessage(scip, file, ".");
3469 }
3470 SCIPinfoMessage(scip, file, ")");
3471
3472 return SCIP_OKAY;
3473 }
3474
3475
3476 /** constraint copying method of constraint handler */
3477 static
SCIP_DECL_CONSCOPY(consCopyOrbitope)3478 SCIP_DECL_CONSCOPY(consCopyOrbitope)
3479 {
3480 SCIP_CONSHDLRDATA* conshdlrdata;
3481 SCIP_CONSDATA* sourcedata;
3482 SCIP_VAR*** sourcevars;
3483 SCIP_VAR*** vars;
3484 int nspcons;
3485 int nblocks;
3486 int i;
3487 int k;
3488 int j;
3489
3490 assert( scip != NULL );
3491 assert( cons != NULL );
3492 assert( sourcescip != NULL );
3493 assert( sourceconshdlr != NULL );
3494 assert( strcmp(SCIPconshdlrGetName(sourceconshdlr), CONSHDLR_NAME) == 0 );
3495 assert( sourcecons != NULL );
3496 assert( varmap != NULL );
3497 assert( valid != NULL );
3498
3499 *valid = TRUE;
3500
3501 SCIPdebugMsg(scip, "Copying method for orbitope constraint handler.\n");
3502
3503 sourcedata = SCIPconsGetData(sourcecons);
3504 assert( sourcedata != NULL );
3505 assert( sourcedata->nspcons > 0 );
3506 assert( sourcedata->nblocks > 0 );
3507 assert( sourcedata->vars != NULL );
3508
3509 conshdlrdata = SCIPconshdlrGetData(sourceconshdlr);
3510 assert( conshdlrdata != NULL );
3511
3512 /* do not copy non-model constraints */
3513 if ( !sourcedata->ismodelcons && !conshdlrdata->forceconscopy )
3514 {
3515 *valid = FALSE;
3516
3517 return SCIP_OKAY;
3518 }
3519
3520 nspcons = sourcedata->nspcons;
3521 nblocks = sourcedata->nblocks;
3522 sourcevars = sourcedata->vars;
3523
3524 SCIP_CALL( SCIPallocBufferArray(scip, &vars, nspcons) );
3525 for (i = 0; i < nspcons && *valid; ++i)
3526 {
3527 SCIP_CALL( SCIPallocBufferArray(scip, &(vars[i]), nblocks) ); /*lint !e866*/
3528
3529 for (j = 0; j < nblocks && *valid; ++j)
3530 {
3531 SCIP_CALL( SCIPgetVarCopy(sourcescip, scip, sourcevars[i][j], &(vars[i][j]), varmap, consmap, global, valid) );
3532 assert( !(*valid) || vars[i][j] != NULL );
3533 }
3534 }
3535
3536 /* only create the target constraint, if all variables could be copied */
3537 if ( *valid )
3538 {
3539 /* create copied constraint */
3540 if ( name == NULL )
3541 name = SCIPconsGetName(sourcecons);
3542
3543 SCIP_CALL( SCIPcreateConsOrbitope(scip, cons, name,
3544 vars, sourcedata->orbitopetype, nspcons, nblocks, sourcedata->resolveprop, sourcedata->ismodelcons,
3545 initial, separate, enforce, check, propagate, local, modifiable, dynamic, removable, stickingatnode) );
3546 }
3547
3548 /* free space; only up to row i if copying failed */
3549 assert( 0 <= i && i <= nspcons );
3550 for (k = i - 1; k >= 0; --k)
3551 SCIPfreeBufferArray(scip, &vars[k]);
3552 SCIPfreeBufferArray(scip, &vars);
3553
3554 return SCIP_OKAY;
3555 }
3556
3557
3558 /** constraint parsing method of constraint handler */
3559 static
SCIP_DECL_CONSPARSE(consParseOrbitope)3560 SCIP_DECL_CONSPARSE(consParseOrbitope)
3561 { /*lint --e{715}*/
3562 const char* s;
3563 char* endptr;
3564 SCIP_ORBITOPETYPE orbitopetype;
3565 SCIP_VAR*** vars;
3566 SCIP_VAR* var;
3567 int nspcons;
3568 int maxnspcons;
3569 int nblocks;
3570 int maxnblocks;
3571 int k;
3572 int j;
3573
3574 assert( success != NULL );
3575
3576 *success = TRUE;
3577 s = str;
3578
3579 /* skip white space */
3580 while ( *s != '\0' && isspace((unsigned char)*s) )
3581 ++s;
3582
3583 if ( strncmp(s, "partOrbitope(", 13) == 0 )
3584 orbitopetype = SCIP_ORBITOPETYPE_PARTITIONING;
3585 else if ( strncmp(s, "packOrbitope(", 13) == 0 )
3586 orbitopetype = SCIP_ORBITOPETYPE_PACKING;
3587 else
3588 {
3589 if ( strncmp(s, "fullOrbitope(", 13) != 0 )
3590 {
3591 SCIPverbMessage(scip, SCIP_VERBLEVEL_MINIMAL, NULL, "Syntax error - expected \"fullOrbitope(\", \"partOrbitope\" or \"packOrbitope\": %s\n", s);
3592 *success = FALSE;
3593 return SCIP_OKAY;
3594 }
3595 orbitopetype = SCIP_ORBITOPETYPE_FULL;
3596 }
3597 s += 13;
3598
3599 /* loop through string */
3600 nspcons = 0;
3601 nblocks = 0;
3602 maxnspcons = 10;
3603 maxnblocks = 10;
3604
3605 SCIP_CALL( SCIPallocBufferArray(scip, &vars, maxnspcons) );
3606 SCIP_CALL( SCIPallocBufferArray(scip, &(vars[0]), maxnblocks) );
3607
3608 j = 0;
3609 do
3610 {
3611 /* skip whitespace */
3612 while ( isspace((int)*s) )
3613 ++s;
3614
3615 /* parse variable name */
3616 SCIP_CALL( SCIPparseVarName(scip, s, &var, &endptr) );
3617 if ( var == NULL )
3618 {
3619 SCIPverbMessage(scip, SCIP_VERBLEVEL_MINIMAL, NULL, "unknown variable name at '%s'\n", str);
3620 *success = FALSE;
3621 return SCIP_OKAY;
3622 }
3623 vars[nspcons][j++] = var;
3624 s = endptr;
3625
3626 if ( j > nblocks )
3627 {
3628 if ( nspcons > 0 )
3629 {
3630 SCIPverbMessage(scip, SCIP_VERBLEVEL_MINIMAL, NULL, "variables per row do not match.\n");
3631 *success = FALSE;
3632 return SCIP_OKAY;
3633 }
3634 nblocks = j;
3635
3636 if ( nblocks > maxnblocks )
3637 {
3638 int newsize;
3639
3640 newsize = SCIPcalcMemGrowSize(scip, nblocks);
3641 SCIP_CALL( SCIPreallocBufferArray(scip, &(vars[nspcons]), newsize) ); /*lint !e866*/
3642 maxnblocks = newsize;
3643 }
3644 }
3645 assert( nblocks <= maxnblocks );
3646
3647 /* skip white space and ',' */
3648 while ( *s != '\0' && ( isspace((unsigned char)*s) || *s == ',' ) )
3649 ++s;
3650
3651 /* begin new row if required */
3652 if ( *s == '.' )
3653 {
3654 ++nspcons;
3655 ++s;
3656
3657 if ( nspcons >= maxnspcons )
3658 {
3659 int newsize;
3660
3661 newsize = SCIPcalcMemGrowSize(scip, nspcons+1);
3662 SCIP_CALL( SCIPreallocBufferArray(scip, &vars, newsize) );
3663 maxnspcons = newsize;
3664 }
3665 assert(nspcons < maxnspcons);
3666
3667 SCIP_CALL( SCIPallocBufferArray(scip, &(vars[nspcons]), nblocks) ); /*lint !e866*/
3668 j = 0;
3669 }
3670 }
3671 while ( *s != ')' );
3672 ++nspcons;
3673
3674 SCIP_CALL( SCIPcreateConsOrbitope(scip, cons, name, vars, orbitopetype, nspcons, nblocks, TRUE, TRUE,
3675 initial, separate, enforce, check, propagate, local, modifiable, dynamic, removable, stickingatnode) );
3676
3677 for (k = nspcons - 1; k >= 0; --k)
3678 SCIPfreeBufferArray(scip, &vars[k]);
3679 SCIPfreeBufferArray(scip, &vars);
3680
3681 return SCIP_OKAY;
3682 }
3683
3684
3685 /** constraint method of constraint handler which returns the variables (if possible) */
3686 static
SCIP_DECL_CONSGETVARS(consGetVarsOrbitope)3687 SCIP_DECL_CONSGETVARS(consGetVarsOrbitope)
3688 { /*lint --e{715}*/
3689 SCIP_CONSDATA* consdata;
3690
3691 assert( cons != NULL );
3692 assert( success != NULL );
3693 assert( vars != NULL );
3694
3695 consdata = SCIPconsGetData(cons);
3696 assert( consdata != NULL );
3697
3698 if ( varssize < consdata->nblocks * consdata->nspcons )
3699 (*success) = FALSE;
3700 else
3701 {
3702 int cnt = 0;
3703 int i;
3704 int j;
3705
3706 for (i = 0; i < consdata->nspcons; ++i)
3707 {
3708 for (j = 0; j < consdata->nblocks; ++j)
3709 vars[cnt++] = consdata->vars[i][j];
3710 }
3711 (*success) = TRUE;
3712 }
3713
3714 return SCIP_OKAY;
3715 }
3716
3717
3718 /** constraint method of constraint handler which returns the number of variables (if possible) */
3719 static
SCIP_DECL_CONSGETNVARS(consGetNVarsOrbitope)3720 SCIP_DECL_CONSGETNVARS(consGetNVarsOrbitope)
3721 { /*lint --e{715}*/
3722 SCIP_CONSDATA* consdata;
3723
3724 assert( cons != NULL );
3725
3726 consdata = SCIPconsGetData(cons);
3727 assert( consdata != NULL );
3728
3729 (*nvars) = consdata->nblocks * consdata->nspcons;
3730 (*success) = TRUE;
3731
3732 return SCIP_OKAY;
3733 }
3734
3735
3736 /*
3737 * constraint specific interface methods
3738 */
3739
3740 /** creates the handler for orbitope constraints and includes it in SCIP */
SCIPincludeConshdlrOrbitope(SCIP * scip)3741 SCIP_RETCODE SCIPincludeConshdlrOrbitope(
3742 SCIP* scip /**< SCIP data structure */
3743 )
3744 {
3745 SCIP_CONSHDLRDATA* conshdlrdata;
3746 SCIP_CONSHDLR* conshdlr;
3747
3748 /* create orbitope constraint handler data */
3749 SCIP_CALL( SCIPallocBlockMemory(scip, &conshdlrdata) );
3750
3751 /* include constraint handler */
3752 SCIP_CALL( SCIPincludeConshdlrBasic(scip, &conshdlr, CONSHDLR_NAME, CONSHDLR_DESC,
3753 CONSHDLR_ENFOPRIORITY, CONSHDLR_CHECKPRIORITY,
3754 CONSHDLR_EAGERFREQ, CONSHDLR_NEEDSCONS,
3755 consEnfolpOrbitope, consEnfopsOrbitope, consCheckOrbitope, consLockOrbitope,
3756 conshdlrdata) );
3757 assert(conshdlr != NULL);
3758
3759 /* set non-fundamental callbacks via specific setter functions */
3760 SCIP_CALL( SCIPsetConshdlrCopy(scip, conshdlr, conshdlrCopyOrbitope, consCopyOrbitope) );
3761 SCIP_CALL( SCIPsetConshdlrFree(scip, conshdlr, consFreeOrbitope) );
3762 SCIP_CALL( SCIPsetConshdlrDelete(scip, conshdlr, consDeleteOrbitope) );
3763 SCIP_CALL( SCIPsetConshdlrGetVars(scip, conshdlr, consGetVarsOrbitope) );
3764 SCIP_CALL( SCIPsetConshdlrGetNVars(scip, conshdlr, consGetNVarsOrbitope) );
3765 SCIP_CALL( SCIPsetConshdlrParse(scip, conshdlr, consParseOrbitope) );
3766 SCIP_CALL( SCIPsetConshdlrPresol(scip, conshdlr, consPresolOrbitope, CONSHDLR_MAXPREROUNDS, CONSHDLR_PRESOLTIMING) );
3767 SCIP_CALL( SCIPsetConshdlrPrint(scip, conshdlr, consPrintOrbitope) );
3768 SCIP_CALL( SCIPsetConshdlrProp(scip, conshdlr, consPropOrbitope, CONSHDLR_PROPFREQ, CONSHDLR_DELAYPROP,
3769 CONSHDLR_PROP_TIMING) );
3770 SCIP_CALL( SCIPsetConshdlrResprop(scip, conshdlr, consRespropOrbitope) );
3771 SCIP_CALL( SCIPsetConshdlrSepa(scip, conshdlr, consSepalpOrbitope, consSepasolOrbitope, CONSHDLR_SEPAFREQ,
3772 CONSHDLR_SEPAPRIORITY, CONSHDLR_DELAYSEPA) );
3773 SCIP_CALL( SCIPsetConshdlrTrans(scip, conshdlr, consTransOrbitope) );
3774 SCIP_CALL( SCIPsetConshdlrEnforelax(scip, conshdlr, consEnforelaxOrbitope) );
3775
3776 SCIP_CALL( SCIPaddBoolParam(scip, "constraints/" CONSHDLR_NAME "/checkpporbitope",
3777 "Strengthen orbitope constraints to packing/partioning orbitopes?",
3778 &conshdlrdata->checkpporbitope, TRUE, DEFAULT_PPORBITOPE, NULL, NULL) );
3779
3780 SCIP_CALL( SCIPaddBoolParam(scip, "constraints/" CONSHDLR_NAME "/sepafullorbitope",
3781 "Whether we separate inequalities for full orbitopes?",
3782 &conshdlrdata->sepafullorbitope, TRUE, DEFAULT_SEPAFULLORBITOPE, NULL, NULL) );
3783
3784 SCIP_CALL( SCIPaddBoolParam(scip, "constraints/" CONSHDLR_NAME "/usedynamicprop",
3785 "Whether we use a dynamic version of the propagation routine.",
3786 &conshdlrdata->usedynamicprop, TRUE, DEFAULT_USEDYNAMICPROP, NULL, NULL) );
3787
3788 SCIP_CALL( SCIPaddBoolParam(scip, "constraints/" CONSHDLR_NAME "/forceconscopy",
3789 "Whether orbitope constraints should be forced to be copied to sub SCIPs.",
3790 &conshdlrdata->forceconscopy, TRUE, DEFAULT_FORCECONSCOPY, NULL, NULL) );
3791
3792 return SCIP_OKAY;
3793 }
3794
3795
3796 /** creates and captures a orbitope constraint
3797 *
3798 * @pre If packing/partitioning orbitopes are used, this constraint handler assumes that constraints which enforce
3799 * the packing/partitioning constraints are contained in the problem. It does not implement, e.g., separation and
3800 * propagation of set packing/partitioning constraints, since this would just copy large parts of the code of the
3801 * setppc constraint handler.
3802 *
3803 * @note the constraint gets captured, hence at one point you have to release it using the method SCIPreleaseCons()
3804 */
SCIPcreateConsOrbitope(SCIP * scip,SCIP_CONS ** cons,const char * name,SCIP_VAR *** vars,SCIP_ORBITOPETYPE orbitopetype,int nspcons,int nblocks,SCIP_Bool resolveprop,SCIP_Bool ismodelcons,SCIP_Bool initial,SCIP_Bool separate,SCIP_Bool enforce,SCIP_Bool check,SCIP_Bool propagate,SCIP_Bool local,SCIP_Bool modifiable,SCIP_Bool dynamic,SCIP_Bool removable,SCIP_Bool stickingatnode)3805 SCIP_RETCODE SCIPcreateConsOrbitope(
3806 SCIP* scip, /**< SCIP data structure */
3807 SCIP_CONS** cons, /**< pointer to hold the created constraint */
3808 const char* name, /**< name of constraint */
3809 SCIP_VAR*** vars, /**< matrix of variables on which the symmetry acts */
3810 SCIP_ORBITOPETYPE orbitopetype, /**< type of orbitope constraint */
3811 int nspcons, /**< number of set partitioning/packing constraints <=> p */
3812 int nblocks, /**< number of symmetric variable blocks <=> q */
3813 SCIP_Bool resolveprop, /**< should propagation be resolved? */
3814 SCIP_Bool ismodelcons, /**< whether the orbitope is a model constraint */
3815 SCIP_Bool initial, /**< should the LP relaxation of constraint be in the initial LP?
3816 * Usually set to TRUE. Set to FALSE for 'lazy constraints'. */
3817 SCIP_Bool separate, /**< should the constraint be separated during LP processing?
3818 * Usually set to TRUE. */
3819 SCIP_Bool enforce, /**< should the constraint be enforced during node processing?
3820 * TRUE for model constraints, FALSE for additional, redundant constraints. */
3821 SCIP_Bool check, /**< should the constraint be checked for feasibility?
3822 * TRUE for model constraints, FALSE for additional, redundant constraints. */
3823 SCIP_Bool propagate, /**< should the constraint be propagated during node processing?
3824 * Usually set to TRUE. */
3825 SCIP_Bool local, /**< is constraint only valid locally?
3826 * Usually set to FALSE. Has to be set to TRUE, e.g., for branching constraints. */
3827 SCIP_Bool modifiable, /**< is constraint modifiable (subject to column generation)?
3828 * Usually set to FALSE. In column generation applications, set to TRUE if pricing
3829 * adds coefficients to this constraint. */
3830 SCIP_Bool dynamic, /**< is constraint subject to aging?
3831 * Usually set to FALSE. Set to TRUE for own cuts which
3832 * are separated as constraints. */
3833 SCIP_Bool removable, /**< should the relaxation be removed from the LP due to aging or cleanup?
3834 * Usually set to FALSE. Set to TRUE for 'lazy constraints' and 'user cuts'. */
3835 SCIP_Bool stickingatnode /**< should the constraint always be kept at the node where it was added, even
3836 * if it may be moved to a more global node?
3837 * Usually set to FALSE. Set to TRUE to for constraints that represent node data. */
3838 )
3839 {
3840 SCIP_CONSHDLRDATA* conshdlrdata;
3841 SCIP_CONSHDLR* conshdlr;
3842 SCIP_CONSDATA* consdata;
3843
3844 /* find the orbitope constraint handler */
3845 conshdlr = SCIPfindConshdlr(scip, CONSHDLR_NAME);
3846 if ( conshdlr == NULL )
3847 {
3848 SCIPerrorMessage("orbitope constraint handler not found\n");
3849 return SCIP_PLUGINNOTFOUND;
3850 }
3851
3852 assert( nspcons > 0 );
3853 assert( nblocks > 0 );
3854
3855 /* run some checks */
3856 #ifndef NDEBUG
3857 {
3858 SCIP_Real obj;
3859 int i;
3860 int j;
3861 for (i = 0; i < nspcons; ++i)
3862 {
3863 /* initialize obj to infinity */
3864 obj = SCIPinfinity(scip);
3865 for (j = 0; j < nblocks; ++j)
3866 {
3867 SCIP_Bool fixedZero;
3868 SCIP_VAR* var;
3869
3870 var = vars[i][j];
3871 assert(var != NULL);
3872
3873 if ( SCIPvarIsNegated(var) )
3874 var = SCIPvarGetNegatedVar(var);
3875
3876 /* all variables need to be binary */
3877 assert( SCIPvarIsBinary(var) );
3878
3879 /* fixed variables have obj = 0; for variables fixed to 0, we assume that there is no
3880 problem (but we cannot always check it, e.g., when in the original problem
3881 variables were fixed and this problem was copied.) */
3882 fixedZero = ( SCIPisZero(scip, SCIPvarGetLbGlobal(var)) && SCIPisZero(scip, SCIPvarGetUbGlobal(var)) );
3883
3884 /* @todo adapt correctness of the following check for sub-scips */
3885 if ( SCIPgetSubscipDepth(scip) == 0 )
3886 {
3887 /* check whether all variables in a row have the same objective */
3888 if ( ! fixedZero && SCIPisInfinity(scip, obj) )
3889 obj = SCIPvarGetObj(var);
3890 else
3891 {
3892 assert( fixedZero || ! SCIPvarIsActive(var) || SCIPisEQ(scip, obj, SCIPvarGetObj(var)) );
3893 }
3894 }
3895 }
3896 }
3897 }
3898 #endif
3899
3900 conshdlrdata = SCIPconshdlrGetData(conshdlr);
3901 if ( conshdlrdata->checkpporbitope && orbitopetype != SCIP_ORBITOPETYPE_PARTITIONING
3902 && orbitopetype != SCIP_ORBITOPETYPE_PACKING )
3903 {
3904 SCIP_CALL( strengthenOrbitopeConstraint(scip, vars, &nspcons, nblocks, &orbitopetype) );
3905 }
3906
3907 /* create constraint data */
3908 SCIP_CALL( consdataCreate(scip, &consdata, vars, nspcons, nblocks, orbitopetype,
3909 resolveprop, conshdlrdata->usedynamicprop, ismodelcons) );
3910
3911 /* create constraint */
3912 SCIP_CALL( SCIPcreateCons(scip, cons, name, conshdlr, consdata, initial, separate, enforce, check, propagate,
3913 local, modifiable, dynamic, removable, stickingatnode) );
3914
3915 return SCIP_OKAY;
3916 }
3917
3918 /** creates and captures an orbitope constraint
3919 * in its most basic variant, i. e., with all constraint flags set to their default values
3920 *
3921 * @note the constraint gets captured, hence at one point you have to release it using the method SCIPreleaseCons()
3922 */
SCIPcreateConsBasicOrbitope(SCIP * scip,SCIP_CONS ** cons,const char * name,SCIP_VAR *** vars,SCIP_ORBITOPETYPE orbitopetype,int nspcons,int nblocks,SCIP_Bool resolveprop,SCIP_Bool ismodelcons)3923 SCIP_RETCODE SCIPcreateConsBasicOrbitope(
3924 SCIP* scip, /**< SCIP data structure */
3925 SCIP_CONS** cons, /**< pointer to hold the created constraint */
3926 const char* name, /**< name of constraint */
3927 SCIP_VAR*** vars, /**< matrix of variables on which the symmetry acts */
3928 SCIP_ORBITOPETYPE orbitopetype, /**< type of orbitope constraint */
3929 int nspcons, /**< number of set partitioning/packing constraints <=> p */
3930 int nblocks, /**< number of symmetric variable blocks <=> q */
3931 SCIP_Bool resolveprop, /**< should propagation be resolved? */
3932 SCIP_Bool ismodelcons /**< whether the orbitope is a model constraint */
3933 )
3934 {
3935 SCIP_CALL( SCIPcreateConsOrbitope(scip, cons, name, vars, orbitopetype, nspcons, nblocks, resolveprop, ismodelcons,
3936 TRUE, TRUE, TRUE, TRUE, TRUE,
3937 FALSE, FALSE, FALSE, FALSE, FALSE) );
3938
3939 return SCIP_OKAY;
3940 }
3941