Scippy

SCIP

Solving Constraint Integer Programs

prop_symmetry.c
Go to the documentation of this file.
1/* * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * */
2/* */
3/* This file is part of the program and library */
4/* SCIP --- Solving Constraint Integer Programs */
5/* */
6/* Copyright (c) 2002-2024 Zuse Institute Berlin (ZIB) */
7/* */
8/* Licensed under the Apache License, Version 2.0 (the "License"); */
9/* you may not use this file except in compliance with the License. */
10/* You may obtain a copy of the License at */
11/* */
12/* http://www.apache.org/licenses/LICENSE-2.0 */
13/* */
14/* Unless required by applicable law or agreed to in writing, software */
15/* distributed under the License is distributed on an "AS IS" BASIS, */
16/* WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. */
17/* See the License for the specific language governing permissions and */
18/* limitations under the License. */
19/* */
20/* You should have received a copy of the Apache-2.0 license */
21/* along with SCIP; see the file LICENSE. If not visit scipopt.org. */
22/* */
23/* * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * */
24
25/**@file prop_symmetry.c
26 * @ingroup DEFPLUGINS_PROP
27 * @brief propagator for handling symmetries
28 * @author Marc Pfetsch
29 * @author Thomas Rehn
30 * @author Christopher Hojny
31 * @author Fabian Wegscheider
32 * @author Jasper van Doornmalen
33 *
34 * This propagator combines the following symmetry handling functionalities:
35 * - It allows to compute symmetries of the problem and to store this information in adequate form. The symmetry
36 * information can be accessed through external functions.
37 * - It implements various methods to handle the symmetries:
38 * - orbital reduction, which generalizes orbital fixing. See symmetry_orbital.c
39 * - (dynamic) orbitopal reduction, which generalizes (dynamic) orbital fixing. See symmetry_orbitopal.c
40 * - static orbitopal fixing (for binary variable domains) for full orbitopes. See cons_orbitope.c
41 * - static orbitopal fixing (for binary variable domains) for packing-partitioning orbitopes. See cons_orbitope.c
42 * - (dynamic) lexicographic reduction. See symmetry_lexred.c
43 * - static lexicographic fixing for binary variable domains (i.e., symresack propagation). See cons_symresack.c
44 * - static lexicographic fixing for binary variable domains on involutions (i.e., orbisacks). See cons_orbisack.c
45 * - Symmetry breaking inequalities based on the Schreier-Sims Table (i.e., SST cuts).
46 * - Strong and weak symmetry breaking inequalities.
47 *
48 *
49 * @section SYMCOMP Symmetry Computation
50 *
51 * The generic functionality of the compute_symmetry.h interface is used.
52 * We do not copy symmetry information, since it is not clear how this information transfers. Moreover, copying
53 * symmetry might inhibit heuristics. But note that solving a sub-SCIP might then happen without symmetry information!
54 *
55 *
56 * @section SYMBREAK Symmetry handling by the (unified) symmetry handling constraints
57 *
58 * Many common methods are captured by a framework that dynamifies symmetry handling constraints. The ideas are
59 * described in@n
60 * J. van Doornmalen, C. Hojny, "A Unified Framework for Symmetry Handling", preprint, 2023,
61 * https://doi.org/10.48550/arXiv.2211.01295.
62 *
63 * This paper shows that various symmetry handling methods are compatible under certain conditions, and provides
64 * generalizations to common symmetry handling constraints from binary variable domains to arbitrary variable domains.
65 * This includes symresack propagation, orbitopal fixing, and orbital fixing, that are generalized to
66 * lexicographic reduction, orbitopal reduction and orbital reduction, respectively. For a description and
67 * implementation, see symmetry_lexred.c, symmetry_orbitopal.c and symmetry_orbital.c, respectively.
68 * The static counterparts on binary variable domains are cons_symresack.c and cons_orbisack.c for lexicographic
69 * reduction (cf. symresack propagation), and cons_orbitope.c and cons_orbisack.c for orbitopal reduction
70 * (cf. orbitopal fixing). We refer to the description of tryAddSymmetryHandlingMethods for the order in which these
71 * methods are applied.
72 *
73 * @section SST Cuts derived from the Schreier Sims table
74 *
75 * SST cuts have been introduced by@n
76 * D. Salvagnin: Symmetry Breaking Inequalities from the Schreier-Sims table. CPAIOR 2018 Proceedings, 521-529, 2018.
77 *
78 * These inequalities are computed as follows. Throughout these procedure a set of so-called leaders is maintained.
79 * Initially the set of leaders is empty. In a first step, select a variable \f$x_i\f$ and compute its orbit w.r.t.
80 * the symmetry group of the mixed-integer program. For each variable \f$x_j\f$ in the orbit of \f$x_i\f$, the
81 * inequality \f$x_i \geq x_j\f$ is a valid symmetry handling inequality, which can be added to the mixed-integer
82 * program. We call \f$x_i\f$ the leader of this inequality. Add the leader \f$x_i\f$ to the set of leaders and
83 * compute the pointwise stabilizer of the leader set. In the next step, select a new variable, compute its orbit
84 * w.r.t. the stabilizer group of the leaders, add the inequalities based on this orbit, and add the new leader
85 * to the set of leaders. This procedure is iterated until the pointwise stabilizer group of the leaders has become
86 * trivial.
87 *
88 * @todo Possibly turn off propagator in subtrees.
89 * @todo Check application of conflict resolution.
90 * @todo Check whether one should switch the role of 0 and 1
91 * @todo Implement stabilizer computation?
92 * @todo Implement isomorphism pruning?
93 * @todo Implement particular preprocessing rules
94 * @todo Separate permuted cuts (first experiments not successful)
95 * @todo Allow the computation of local symmetries
96 * @todo Order rows of orbitopes (in particular packing/partitioning) w.r.t. cliques in conflict graph.
97 * @todo A dynamic variant for packing-partitioning orbitopal structures
98 * @todo A dynamic variant for suborbitopes
99 */
100/* #define SCIP_OUTPUT */
101/* #define SCIP_OUTPUT_COMPONENT */
102/* #define SCIP_DISPLAY_SYM_CHECK */
103
104/*---+----1----+----2----+----3----+----4----+----5----+----6----+----7----+----8----+----9----+----0----+----1----+----2*/
105
106#include "scip/cons_linear.h"
107#include "scip/cons_knapsack.h"
108#include "scip/cons_varbound.h"
109#include "scip/cons_setppc.h"
110#include "scip/cons_and.h"
111#include "scip/cons_logicor.h"
112#include "scip/cons_or.h"
113#include "scip/cons_orbitope.h"
114#include "scip/cons_symresack.h"
115#include "scip/cons_xor.h"
116#include "scip/cons_linking.h"
118#include "scip/cons_indicator.h"
119#include "scip/cons_nonlinear.h"
120#include "scip/cons_sos1.h"
121#include "scip/cons_sos2.h"
122#include "scip/expr_pow.h"
123#include "scip/expr_product.h"
124#include "scip/pub_expr.h"
125#include "scip/misc.h"
127
128#include "scip/prop_symmetry.h"
131#include "scip/symmetry.h"
132#include "scip/symmetry_graph.h"
135#include "scip/symmetry_lexred.h"
136
137#include <math.h>
138#include <string.h>
139
140/* propagator properties */
141#define PROP_NAME "symmetry"
142#define PROP_DESC "propagator for handling symmetry"
143#define PROP_TIMING SCIP_PROPTIMING_BEFORELP /**< propagation timing mask */
144#define PROP_PRIORITY -1000000 /**< propagator priority */
145#define PROP_FREQ 1 /**< propagator frequency */
146#define PROP_DELAY FALSE /**< should propagation method be delayed, if other propagators found reductions? */
147
148#define PROP_PRESOL_PRIORITY -10000000 /**< priority of the presolving method (>= 0: before, < 0: after constraint handlers) */
149#define PROP_PRESOLTIMING SCIP_PRESOLTIMING_EXHAUSTIVE /* timing of the presolving method (fast, medium, or exhaustive) */
150#define PROP_PRESOL_MAXROUNDS -1 /**< maximal number of presolving rounds the presolver participates in (-1: no limit) */
151
152
153/* default parameter values for symmetry computation */
154#define DEFAULT_MAXGENERATORS 1500 /**< limit on the number of generators that should be produced within symmetry detection (0 = no limit) */
155#define DEFAULT_CHECKSYMMETRIES FALSE /**< Should all symmetries be checked after computation? */
156#define DEFAULT_DISPLAYNORBITVARS FALSE /**< Should the number of variables affected by some symmetry be displayed? */
157#define DEFAULT_USECOLUMNSPARSITY FALSE /**< Should the number of conss a variable is contained in be exploited in symmetry detection? */
158#define DEFAULT_DOUBLEEQUATIONS FALSE /**< Double equations to positive/negative version? */
159#define DEFAULT_COMPRESSSYMMETRIES TRUE /**< Should non-affected variables be removed from permutation to save memory? */
160#define DEFAULT_COMPRESSTHRESHOLD 0.5 /**< Compression is used if percentage of moved vars is at most the threshold. */
161#define DEFAULT_SYMFIXNONBINARYVARS FALSE /**< Disabled parameter */
162#define DEFAULT_ENFORCECOMPUTESYMMETRY FALSE /**< always compute symmetries, even if they cannot be handled */
163#define DEFAULT_SYMTYPE (int) SYM_SYMTYPE_PERM /**< type of symmetries to be computed */
164
165/* default parameters for linear symmetry constraints */
166#define DEFAULT_CONSSADDLP TRUE /**< Should the symmetry breaking constraints be added to the LP? */
167#define DEFAULT_ADDSYMRESACKS TRUE /**< Add inequalities for symresacks for each generator? */
168#define DEFAULT_DETECTDOUBLELEX TRUE /**< Should we check whether the components of the symmetry group can be handled by double lex matrices? */
169#define DEFAULT_DETECTORBITOPES TRUE /**< Should we check whether the components of the symmetry group can be handled by orbitopes? */
170#define DEFAULT_DETECTSUBGROUPS TRUE /**< Should we try to detect orbitopes in subgroups of the symmetry group? */
171#define DEFAULT_ADDWEAKSBCS TRUE /**< Should we add weak SBCs for enclosing orbit of symmetric subgroups? */
172#define DEFAULT_ADDSTRONGSBCS FALSE /**< Should we add strong SBCs for enclosing orbit of symmetric subgroups if orbitopes are not used? */
173#define DEFAULT_ADDCONSSTIMING 2 /**< timing of adding constraints (0 = before presolving, 1 = during presolving, 2 = after presolving) */
174#define DEFAULT_MAXNCONSSSUBGROUP 500000 /**< Maximum number of constraints up to which subgroup structures are detected */
175#define DEFAULT_USEDYNAMICPROP TRUE /**< whether dynamic propagation should be used for full orbitopes */
176#define DEFAULT_PREFERLESSROWS TRUE /**< Shall orbitopes with less rows be preferred in detection? */
177
178/* default parameters for symmetry computation */
179#define DEFAULT_SYMCOMPTIMING 2 /**< timing of symmetry computation (0 = before presolving, 1 = during presolving, 2 = at first call) */
180#define DEFAULT_PERFORMPRESOLVING 0 /**< Run orbital fixing during presolving? (disabled parameter) */
181#define DEFAULT_RECOMPUTERESTART 0 /**< Recompute symmetries after a restart has occurred? (0 = never) */
182
183/* default parameters for Schreier Sims constraints */
184#define DEFAULT_SSTTIEBREAKRULE 1 /**< index of tie break rule for selecting orbit for Schreier Sims constraints? */
185#define DEFAULT_SSTLEADERRULE 0 /**< index of rule for selecting leader variables for Schreier Sims constraints? */
186#define DEFAULT_SSTLEADERVARTYPE 14 /**< bitset encoding which variable types can be leaders (1: binary; 2: integer; 4: impl. int; 8: continuous);
187 * if multiple types are allowed, take the one with most affected vars */
188#define DEFAULT_ADDCONFLICTCUTS TRUE /**< Should Schreier Sims constraints be added if we use a conflict based rule? */
189#define DEFAULT_SSTADDCUTS TRUE /**< Should Schreier Sims constraints be added? */
190#define DEFAULT_SSTMIXEDCOMPONENTS TRUE /**< Should Schreier Sims constraints be added if a symmetry component contains variables of different types? */
191
192/* output table properties */
193#define TABLE_NAME_SYMMETRY "symmetry"
194#define TABLE_DESC_SYMMETRY "symmetry handling statistics"
195#define TABLE_POSITION_SYMMETRY 7001 /**< the position of the statistics table */
196#define TABLE_EARLIEST_SYMMETRY SCIP_STAGE_SOLVING /**< output of the statistics table is only printed from this stage onwards */
197
198
199/* other defines */
200#define MAXGENNUMERATOR 64000000 /**< determine maximal number of generators by dividing this number by the number of variables */
201#define COMPRESSNVARSLB 25000 /**< lower bound on the number of variables above which compression could be performed */
202#define DEFAULT_NAUTYMAXNCELLS 100000 /**< terminate symmetry detection using Nauty when number of cells in color refinment is at least this number
203 * (avoids segfaults due to Nauty for large graphs) */
204#define DEFAULT_NAUTYMAXNNODES 10000000 /**< terminate symmetry detection using Nauty when its search tree has at least this number of nodes */
205/*@todo investigate why the Nauty works well for some large instances (miplib2010/mspp16.mps) but not for PB instances (e.g., normalized-celar6-sub0_wcsp.wbo) */
206
207
208/* macros for getting activeness of symmetry handling methods */
209#define ISSYMRETOPESACTIVE(x) (((unsigned) x & SYM_HANDLETYPE_SYMBREAK) != 0)
210#define ISORBITALREDUCTIONACTIVE(x) (((unsigned) x & SYM_HANDLETYPE_ORBITALREDUCTION) != 0)
211#define ISSSTACTIVE(x) (((unsigned) x & SYM_HANDLETYPE_SST) != 0)
212
213#define ISSSTBINACTIVE(x) (((unsigned) x & SCIP_SSTTYPE_BINARY) != 0)
214#define ISSSTINTACTIVE(x) (((unsigned) x & SCIP_SSTTYPE_INTEGER) != 0)
215#define ISSSTIMPLINTACTIVE(x) (((unsigned) x & SCIP_SSTTYPE_IMPLINT) != 0)
216#define ISSSTCONTACTIVE(x) (((unsigned) x & SCIP_SSTTYPE_CONTINUOUS) != 0)
217
218/* enable symmetry statistics */
219#define SYMMETRY_STATISTICS 1
220
221/** propagator data */
222struct SCIP_PropData
223{
224 /* symmetry group information */
225 int npermvars; /**< number of variables for permutations */
226 int nbinpermvars; /**< number of binary variables for permutations */
227 SCIP_VAR** permvars; /**< variables on which permutations act */
228 int nperms; /**< number of permutations */
229 int nmaxperms; /**< maximal number of permutations (needed for freeing storage) */
230 int** perms; /**< pointer to store permutation generators as (nperms x npermvars) matrix */
231 int** permstrans; /**< pointer to store transposed permutation generators as (npermvars x nperms) matrix */
232 SCIP_HASHMAP* permvarmap; /**< map of variables to indices in permvars array */
233 int nmovedpermvars; /**< number of variables moved by any permutation */
234 int nmovedbinpermvars; /**< number of binary variables moved by any permutation */
235 int nmovedintpermvars; /**< number of integer variables moved by any permutation */
236 int nmovedimplintpermvars; /**< number of implicitly integer variables moved by any permutation */
237 int nmovedcontpermvars; /**< number of continuous variables moved by any permutation */
238 SCIP_HASHMAP* customsymopnodetypes; /**< types of operator nodes introduced
239 * by a user for symmetry detection */
240 int nopnodetypes; /**< current number of operator node types used for symmetry detection */
241 SCIP_Real* permvardomaincenter; /**< center of variable domains (needed for signed permutations) */
242 int symtype; /**< type of symmetries to be computed */
243
244 /* components of symmetry group */
245 int ncomponents; /**< number of components of symmetry group */
246 int ncompblocked; /**< number of components that have been blocked */
247 int* components; /**< array containing the indices of permutations sorted by components */
248 int* componentbegins; /**< array containing in i-th position the first position of
249 * component i in components array */
250 int* vartocomponent; /**< array containing for each permvar the index of the component it is
251 * contained in (-1 if not affected) */
252 unsigned* componentblocked; /**< array to store which symmetry methods have been applied to a component using
253 * the same bitset as for misc/usesymmetry */
254 SCIP_Bool* componenthassignedperm; /**< array to indicate whether a component has a signed permutation */
255
256 /* further symmetry information */
257 int nmovedvars; /**< number of variables moved by some permutation */
258 SCIP_Real log10groupsize; /**< log10 of size of symmetry group */
259 SCIP_Bool binvaraffected; /**< whether binary variables are affected by some symmetry */
260
261 /* for symmetry computation */
262 int maxgenerators; /**< limit on the number of generators that should be produced within symmetry detection (0 = no limit) */
263 SCIP_Bool checksymmetries; /**< Should all symmetries be checked after computation? */
264 SCIP_Bool displaynorbitvars; /**< Whether the number of variables in non-trivial orbits shall be computed */
265 SCIP_Bool compresssymmetries; /**< Should non-affected variables be removed from permutation to save memory? */
266 SCIP_Real compressthreshold; /**< Compression is used if percentage of moved vars is at most the threshold. */
267 SCIP_Bool compressed; /**< Whether symmetry data has been compressed */
268 SCIP_Bool computedsymmetry; /**< Have we already tried to compute symmetries? */
269 int usesymmetry; /**< encoding of active symmetry handling methods (for debugging) */
270 SCIP_Bool usecolumnsparsity; /**< Should the number of conss a variable is contained in be exploited in symmetry detection? */
271 SCIP_Bool doubleequations; /**< Double equations to positive/negative version? */
272 SCIP_Bool enforcecomputesymmetry; /**< always compute symmetries, even if they cannot be handled */
273 int symtiming; /**< timing of computing and handling symmetries (0 = before presolving, 1 = during presolving, 2 = after presolving) */
274
275 /* for symmetry constraints */
276 SCIP_Bool triedaddsymmethods; /**< whether we already tried to add symmetry handling methods */
277 SCIP_Bool conssaddlp; /**< Should the symmetry breaking constraints be added to the LP? */
278 SCIP_Bool addsymresacks; /**< Add symresack constraints for each generator? */
279 SCIP_CONS** genorbconss; /**< list of generated orbitope/orbisack/symresack constraints */
280 SCIP_CONS** genlinconss; /**< list of generated linear constraints */
281 int ngenorbconss; /**< number of generated orbitope/orbisack/symresack constraints */
282 int genorbconsssize; /**< size of generated orbitope/orbisack/symresack constraints array */
283 int ngenlinconss; /**< number of generated linear constraints */
284 int genlinconsssize; /**< size of linear constraints array */
285 int nsymresacks; /**< number of symresack constraints */
286 SCIP_Bool detectdoublelex; /**< Should we check whether the components of the symmetry group can be handled by double lex matrices? */
287 SCIP_Bool detectorbitopes; /**< Should we check whether the components of the symmetry group can be handled by orbitopes? */
288 SCIP_Bool detectsubgroups; /**< Should we try to detect orbitopes in subgroups of the symmetry group? */
289 SCIP_Bool addweaksbcs; /**< Should we add weak SBCs for enclosing orbit of symmetric subgroups? */
290 SCIP_Bool addstrongsbcs; /**< Should we add strong SBCs for enclosing orbit of symmetric subgroups if orbitopes are not used? */
291 int norbitopes; /**< number of orbitope constraints */
292 SCIP_Bool* isnonlinvar; /**< array indicating whether variables appear non-linearly */
293 SCIP_CONSHDLR* conshdlr_nonlinear; /**< nonlinear constraint handler */
294 int maxnconsssubgroup; /**< maximum number of constraints up to which subgroup structures are detected */
295 SCIP_Bool usedynamicprop; /**< whether dynamic propagation should be used for full orbitopes */
296 SCIP_Bool preferlessrows; /**< Shall orbitopes with less rows be preferred in detection? */
297
298 /* data necessary for symmetry computation order */
299 int recomputerestart; /**< Recompute symmetries after a restart has occurred? (0 = never, 1 = always, 2 = if symmetry reduction found) */
300 int lastrestart; /**< last restart for which symmetries have been computed */
301 SCIP_Bool symfoundreduction; /**< whether symmetry handling propagation has found a reduction since the last time computing symmetries */
302
303 /* data necessary for Schreier Sims constraints */
304 SCIP_CONS** sstconss; /**< list of generated schreier sims conss */
305 int nsstconss; /**< number of generated schreier sims conss */
306 int maxnsstconss; /**< maximum number of conss in sstconss */
307 int sstleaderrule; /**< rule to select leader */
308 int ssttiebreakrule; /**< tie break rule for leader selection */
309 int sstleadervartype; /**< bitset encoding which variable types can be leaders;
310 * if multiple types are allowed, take the one with most affected vars */
311 int* leaders; /**< index of orbit leaders in permvars */
312 int nleaders; /**< number of orbit leaders in leaders array */
313 int maxnleaders; /**< maximum number of leaders in leaders array */
314 SCIP_Bool addconflictcuts; /**< Should Schreier Sims constraints be added if we use a conflict based rule? */
315 SCIP_Bool sstaddcuts; /**< Should Schreier Sims constraints be added? */
316 SCIP_Bool sstmixedcomponents; /**< Should Schreier Sims constraints be added if a symmetry component contains variables of different types? */
317
318 SCIP_EVENTHDLR* shadowtreeeventhdlr;/**< pointer to event handler for shadow tree */
319 SCIP_ORBITOPALREDDATA* orbitopalreddata; /**< container for the orbitopal reduction data */
320 SCIP_ORBITALREDDATA* orbitalreddata; /**< container for orbital reduction data */
321 SCIP_LEXREDDATA* lexreddata; /**< container for lexicographic reduction propagation */
322};
323
324/** conflict data structure for SST cuts */
325struct SCIP_ConflictData
326{
327 SCIP_VAR* var; /**< variable belonging to node */
328 int orbitidx; /**< orbit of variable w.r.t. current stabilizer subgroup
329 * or -1 if not affected by symmetry */
330 int nconflictinorbit; /**< number of variables the node's var is in conflict with */
331 int orbitsize; /**< size of the variable's orbit */
332 int posinorbit; /**< position of variable in its orbit */
333 SCIP_Bool active; /**< whether variable has not been fixed by Schreier Sims code */
334 SCIP_CLIQUE** cliques; /**< List of setppc constraints. */
335 int ncliques; /**< Number of setppc constraints. */
336};
337typedef struct SCIP_ConflictData SCIP_CONFLICTDATA;
338
339
340/** compare function for sorting an array by the addresses of its members */
341static
342SCIP_DECL_SORTPTRCOMP(sortByPointerValue)
343{
344 /* @todo move to misc.c? */
345 if ( elem1 < elem2 )
346 return -1;
347 else if ( elem1 > elem2 )
348 return +1;
349 return 0;
350}
351
352
353/** checks whether two arrays that are sorted with the same comparator have a common element */
354static
356 void** arr1, /**< first array */
357 int narr1, /**< number of elements in first array */
358 void** arr2, /**< second array */
359 int narr2, /**< number of elements in second array */
360 SCIP_DECL_SORTPTRCOMP((*compfunc)) /**< comparator function that was used to sort arr1 and arr2; must define a total ordering */
361 )
362{
363 /* @todo move to misc.c? */
364 int it1;
365 int it2;
366 int cmp;
367
368 assert( arr1 != NULL || narr1 == 0 );
369 assert( narr1 >= 0 );
370 assert( arr2 != NULL || narr2 == 0 );
371 assert( narr2 >= 0 );
372 assert( compfunc != NULL );
373
374 /* there is no overlap if one of the two arrays is empty */
375 if ( narr1 <= 0 )
376 return FALSE;
377 if ( narr2 <= 0 )
378 return FALSE;
379
380 it1 = 0;
381 it2 = 0;
382
383 while ( TRUE ) /*lint !e716*/
384 {
385 cmp = compfunc(arr1[it1], arr2[it2]);
386 if ( cmp < 0 )
387 {
388 /* comparison function determines arr1[it1] < arr2[it2]
389 * increase iterator for arr1
390 */
391 if ( ++it1 >= narr1 )
392 break;
393 continue;
394 }
395 else if ( cmp > 0 )
396 {
397 /* comparison function determines arr1[it1] > arr2[it2]
398 * increase iterator for arr2
399 */
400 if ( ++it2 >= narr2 )
401 break;
402 continue;
403 }
404 else
405 {
406 /* the entries arr1[it1] and arr2[it2] are the same with respect to the comparison function */
407 assert( cmp == 0 );
408 return TRUE;
409 }
410 }
411
412 /* no overlap detected */
413 assert( it1 >= narr1 || it2 >= narr2 );
414 return FALSE;
415}
416
417
418/*
419 * Display dialog callback methods
420 */
421
422/** displays the cycle of a symmetry */
423static
425 SCIP* scip, /**< SCIP pointer */
426 int* perm, /**< symmetry */
427 SYM_SYMTYPE symtype, /**< type of symmetry */
428 int baseidx, /**< variable index for which cycle is computed */
429 SCIP_Bool* covered, /**< allocated array to store covered variables */
430 int nvars, /**< number of (non-negated) variables in symmetry */
431 SCIP_VAR** vars /**< variables on which symmetry acts */
432 )
433{
434 char* string;
435 int varidx;
436 int j;
437
438 assert( scip != NULL );
439 assert( perm != NULL );
440 assert( 0 <= baseidx );
441 assert( (symtype == SYM_SYMTYPE_PERM && baseidx < nvars) ||
442 (symtype == SYM_SYMTYPE_SIGNPERM && baseidx < 2 * nvars) );
443 assert( covered != NULL );
444
445 /* skip fixed points or elements already covered in previous cycle */
446 if ( perm[baseidx] == baseidx || covered[baseidx] )
447 return SCIP_OKAY;
448
449 varidx = baseidx >= nvars ? baseidx - nvars : baseidx;
450 string = (char*) SCIPvarGetName(vars[varidx]);
451 SCIPinfoMessage(scip, NULL, " (%s<%s>", baseidx >= nvars ? "negated " : "", string);
452 j = perm[baseidx];
453 covered[baseidx] = TRUE;
454 while ( j != baseidx )
455 {
456 covered[j] = TRUE;
457 varidx = j >= nvars ? j - nvars : j;
458 string = (char*) SCIPvarGetName(vars[varidx]);
459 SCIPinfoMessage(scip, NULL, ",%s<%s>", j >= nvars ? "negated " : "", string);
460 j = perm[j];
461 }
462 SCIPinfoMessage(scip, NULL, ")\n");
463
464 return SCIP_OKAY;
465}
466
467/** displays symmetry information without taking components into account */
468static
470 SCIP* scip, /**< SCIP pointer */
471 SCIP_PROPDATA* propdata /**< propagator data */
472 )
473{
474 SCIP_Bool* covered;
475 SYM_SYMTYPE symtype;
476 int* perm;
477 int permlen;
478 int npermvars;
479 int i;
480 int p;
481
482 assert( scip != NULL );
483 assert( propdata != NULL );
484 assert( propdata->nperms > 0 );
485 assert( propdata->permvars != NULL );
486 assert( propdata->npermvars > 0 );
487
488 symtype = (SYM_SYMTYPE) propdata->symtype;
489 npermvars = propdata->npermvars;
490 permlen = symtype == SYM_SYMTYPE_PERM ? npermvars : 2 * npermvars;
491
492 if ( symtype == SYM_SYMTYPE_SIGNPERM )
493 SCIPinfoMessage(scip, NULL, "Display permutations as signed permutations (allowing translations)\n");
494
495 SCIP_CALL( SCIPallocClearBufferArray(scip, &covered, permlen) );
496
497 for (p = 0; p < propdata->nperms; ++p)
498 {
499 SCIPinfoMessage(scip, NULL, "Permutation %d:\n", p);
500 perm = propdata->perms[p];
501
502 for (i = 0; i < permlen; ++i)
503 {
504 SCIP_CALL( displayCycleOfSymmetry(scip, perm, symtype, i, covered, npermvars, propdata->permvars) );
505 }
506
507 for (i = 0; i < permlen; ++i)
508 covered[i] = FALSE;
509 }
510
511 SCIPfreeBufferArray(scip, &covered);
512
513 return SCIP_OKAY;
514}
515
516/** displays symmetry information taking components into account */
517static
519 SCIP* scip, /**< SCIP pointer */
520 SCIP_PROPDATA* propdata /**< propagator data */
521 )
522{
523 SCIP_Bool* covered;
524 SYM_SYMTYPE symtype;
525 int* perm;
526 int comppermlen;
527 int permlen;
528 int npermvars;
529 int i;
530 int p;
531 int c;
532
533 assert( scip != NULL );
534 assert( propdata != NULL );
535 assert( propdata->nperms > 0 );
536 assert( propdata->permvars != NULL );
537 assert( propdata->npermvars > 0 );
538 assert( propdata->ncomponents > 0 );
539
540 symtype = (SYM_SYMTYPE) propdata->symtype;
541 npermvars = propdata->npermvars;
542 permlen = symtype == SYM_SYMTYPE_PERM ? npermvars : 2 * npermvars;
543
544 SCIP_CALL( SCIPallocClearBufferArray(scip, &covered, permlen) );
545
546 for (c = 0; c < propdata->ncomponents; ++c)
547 {
548 int cnt;
549
550 SCIPinfoMessage(scip, NULL, "Display symmetries of component %d.\n", c);
551 if ( propdata->componenthassignedperm[c] )
552 SCIPinfoMessage(scip, NULL, " Symmetries are displayed as signed permutations (allowing translations).\n");
553 else
554 SCIPinfoMessage(scip, NULL, " Symmetries are displayed as permutations.\n");
555
556 comppermlen = propdata->componenthassignedperm[c] ? 2 * npermvars : npermvars;
557
558 for (p = propdata->componentbegins[c], cnt = 0; p < propdata->componentbegins[c + 1]; ++p, ++cnt)
559 {
560 SCIPinfoMessage(scip, NULL, "Permutation %d:\n", cnt);
561 perm = propdata->perms[propdata->components[p]];
562
563 for (i = 0; i < comppermlen; ++i)
564 {
565 SCIP_CALL( displayCycleOfSymmetry(scip, perm, symtype, i, covered, npermvars, propdata->permvars) );
566 }
567
568 for (i = 0; i < comppermlen; ++i)
569 covered[i] = FALSE;
570 }
571 }
572
573 SCIPfreeBufferArray(scip, &covered);
574
575 return SCIP_OKAY;
576}
577
578/** dialog execution method for the display symmetry information command */
579static
580SCIP_DECL_DIALOGEXEC(dialogExecDisplaySymmetry)
581{ /*lint --e{715}*/
582 SCIP_PROPDATA* propdata;
583
584 /* add your dialog to history of dialogs that have been executed */
585 SCIP_CALL( SCIPdialoghdlrAddHistory(dialoghdlr, dialog, NULL, FALSE) );
586
587 propdata = (SCIP_PROPDATA*)SCIPdialogGetData(dialog);
588 assert( propdata != NULL );
589
590 if ( propdata->nperms == -1 )
591 {
592 SCIPinfoMessage(scip, NULL, "Cannot display symmetries. Symmetries have not been computed yet.\n");
593 }
594 else if ( propdata->nperms == 0 )
595 {
596 SCIPinfoMessage(scip, NULL, "Cannot display symmetries. No symmetries detected.\n");
597 }
598 else if ( propdata->ncomponents < 0 )
599 {
601 }
602 else
603 {
605 }
606
607 /* next dialog will be root dialog again */
608 *nextdialog = SCIPdialoghdlrGetRoot(dialoghdlr);
609
610 return SCIP_OKAY;
611}
612
613
614/*
615 * Table callback methods
616 */
617
618/** table data */
619struct SCIP_TableData
620{
621 SCIP_PROPDATA* propdata; /** pass data of propagator for table output function */
622};
623
624
625/** output method of symmetry propagator statistics table to output file stream 'file' */
626static
627SCIP_DECL_TABLEOUTPUT(tableOutputSymmetry)
628{
629 SCIP_TABLEDATA* tabledata;
630 int nred;
631 int ncutoff;
632 SCIP_Real time;
633
634 assert( scip != NULL );
635 assert( table != NULL );
636
637 tabledata = SCIPtableGetData(table);
638 assert( tabledata != NULL );
639 assert( tabledata->propdata != NULL );
640
641 if ( tabledata->propdata->orbitopalreddata || tabledata->propdata->orbitalreddata
642 || tabledata->propdata->lexreddata )
643 {
644 SCIPverbMessage(scip, SCIP_VERBLEVEL_MINIMAL, file, "Symmetry :\n");
645 if ( tabledata->propdata->orbitopalreddata )
646 {
647 SCIP_CALL( SCIPorbitopalReductionGetStatistics(scip, tabledata->propdata->orbitopalreddata, &nred, &ncutoff) );
648 SCIPverbMessage(scip, SCIP_VERBLEVEL_MINIMAL, file, " orbitopal red. : %10d reductions applied,"
649 " %10d cutoffs\n", nred, ncutoff);
650 }
651 if ( tabledata->propdata->orbitalreddata )
652 {
653 SCIP_CALL( SCIPorbitalReductionGetStatistics(scip, tabledata->propdata->orbitalreddata, &nred, &ncutoff) );
654 SCIPverbMessage(scip, SCIP_VERBLEVEL_MINIMAL, file, " orbital reduction: %10d reductions applied,"
655 " %10d cutoffs\n", nred, ncutoff);
656 }
657 if ( tabledata->propdata->lexreddata )
658 {
659 SCIP_CALL( SCIPlexicographicReductionGetStatistics(scip, tabledata->propdata->lexreddata, &nred, &ncutoff) );
660 SCIPverbMessage(scip, SCIP_VERBLEVEL_MINIMAL, file, " lexicographic red: %10d reductions applied,"
661 " %10d cutoffs\n", nred, ncutoff);
662 }
663 if ( tabledata->propdata->shadowtreeeventhdlr )
664 {
665 time = SCIPgetShadowTreeEventHandlerExecutionTime(scip, tabledata->propdata->shadowtreeeventhdlr);
666 SCIPverbMessage(scip, SCIP_VERBLEVEL_MINIMAL, file, " shadow tree time : %10.2f s\n", time);
667 }
668 }
669
670 return SCIP_OKAY;
671}
672
673
674/** destructor of statistics table to free user data (called when SCIP is exiting) */
675static
676SCIP_DECL_TABLEFREE(tableFreeSymmetry)
677{
678 SCIP_TABLEDATA* tabledata;
679 tabledata = SCIPtableGetData(table);
680 assert( tabledata != NULL );
681
682 SCIPfreeBlockMemory(scip, &tabledata);
683
684 return SCIP_OKAY;
685}
686
687
688
689/*
690 * local data structures
691 */
692
693/** data structure to store arrays used for sorting colored component types */
695{
696 int* components; /**< array of components */
697 int* colors; /**< array of colors */
698};
700
701/** sorts variable indices according to their corresponding component in the graph
702 *
703 * Variables are sorted first by the color of their component and then by the component index.
704 *
705 * result:
706 * < 0: ind1 comes before (is better than) ind2
707 * = 0: both indices have the same value
708 * > 0: ind2 comes after (is worse than) ind2
709 */
710static
711SCIP_DECL_SORTINDCOMP(SYMsortGraphCompVars)
712{
714
715 data = (SYM_SORTGRAPHCOMPVARS*) dataptr;
716
717 if ( data->colors[ind1] < data->colors[ind2] )
718 return -1;
719 else if ( data->colors[ind1] > data->colors[ind2] )
720 return 1;
721
722 if ( data->components[ind1] < data->components[ind2] )
723 return -1;
724 if ( data->components[ind1] > data->components[ind2] )
725 return 1;
726
727 return 0;
728}
729
730/** compares two symmetry detection graphs
731 *
732 * Graphs are compared by their number of consnodes, then their constypes, then by their lhs,
733 * then by their rhs, then by their total number of nodes, then by the number of operator nodes,
734 * then by their number of value nodes, and then by their number of edges.
735 *
736 * result:
737 * < 0: ind1 comes before (is better than) ind2
738 * = 0: both indices have the same value
739 * > 0: ind2 comes after (is worse than) ind2
740 */
741static
743 SCIP* scip, /**< SCIP pointer (or NULL) */
744 SYM_GRAPH* G1, /**< first graph in comparison */
745 SYM_GRAPH* G2 /**< second graph in comparison */
746 )
747{
748 int c;
749 int perm1;
750 int perm2;
751
752 /* compare the number of constraint nodes */
753 if ( G1->nconsnodes < G2->nconsnodes )
754 return -1;
755 if ( G1->nconsnodes > G2->nconsnodes )
756 return 1;
757
758 /* compare the constraint nodes of the two graphs */
759 for (c = 0; c < G1->nconsnodes; ++c)
760 {
761 perm1 = G1->consnodeperm[c];
762 perm2 = G2->consnodeperm[c];
763
764 if ( SCIPconsGetHdlr(G1->conss[perm1]) < SCIPconsGetHdlr(G2->conss[perm2]) )
765 return -1;
766 if ( SCIPconsGetHdlr(G1->conss[perm1]) > SCIPconsGetHdlr(G2->conss[perm2]) )
767 return 1;
768
769 /* compare using SCIP functions when SCIP is available */
770 if ( scip != NULL )
771 {
772 if ( SCIPisLT(scip, G1->lhs[perm1], G2->lhs[perm2]) )
773 return -1;
774 if ( SCIPisGT(scip, G1->lhs[perm1], G2->lhs[perm2]) )
775 return 1;
776
777 if ( SCIPisLT(scip, G1->rhs[perm1], G2->rhs[perm2]) )
778 return -1;
779 if ( SCIPisGT(scip, G1->rhs[perm1], G2->rhs[perm2]) )
780 return 1;
781 }
782 else
783 {
784 if ( G1->lhs[perm1] < G2->lhs[perm2] )
785 return -1;
786 if ( G1->lhs[perm1] > G2->lhs[perm2] )
787 return 1;
788
789 if ( G1->rhs[perm1] < G2->rhs[perm2] )
790 return -1;
791 if ( G1->rhs[perm1] > G2->rhs[perm2] )
792 return 1;
793 }
794 }
795
796 /* compare number of remaining node types */
797 if ( G1->nnodes < G2->nnodes )
798 return -1;
799 if ( G1->nnodes > G2->nnodes )
800 return 1;
801
802 if ( G1->nopnodes < G2->nopnodes )
803 return -1;
804 if ( G1->nopnodes > G2->nopnodes )
805 return 1;
806
807 if ( G1->nvalnodes < G2->nvalnodes )
808 return -1;
809 if ( G1->nvalnodes > G2->nvalnodes )
810 return 1;
811
812 if ( G1->nedges < G2->nedges )
813 return -1;
814 if ( G1->nedges > G2->nedges )
815 return 1;
816
817 return 0;
818}
819
820/** sorts symmetry detection graphs
821 *
822 * Graphs are sorted by their number of consnodes, then their constypes, then by their lhs,
823 * then by their rhs, then by their total number of nodes, then by the number of operator nodes,
824 * then by their number of value nodes, and then by their number of edges.
825 *
826 * result:
827 * < 0: ind1 comes before (is better than) ind2
828 * = 0: both indices have the same value
829 * > 0: ind2 comes after (is worse than) ind2
830 */
831static
832SCIP_DECL_SORTINDCOMP(SYMsortSymgraphs)
833{
834 SYM_GRAPH** data;
835 SYM_GRAPH* G1;
836 SYM_GRAPH* G2;
837
838 data = (SYM_GRAPH**) dataptr;
839 G1 = data[ind1];
840 G2 = data[ind2];
841
842 return compareSymgraphs(NULL, G1, G2);
843}
844
845/*
846 * Local methods
847 */
848
849#ifndef NDEBUG
850/** checks that symmetry data is all freed */
851static
853 SCIP_PROPDATA* propdata /**< propagator data */
854 )
855{
856 assert( propdata->permvarmap == NULL );
857 assert( propdata->genorbconss == NULL );
858 assert( propdata->genlinconss == NULL );
859 assert( propdata->ngenlinconss == 0 );
860 assert( propdata->ngenorbconss == 0 );
861 assert( propdata->genorbconsssize == 0 );
862 assert( propdata->genlinconsssize == 0 );
863 assert( propdata->sstconss == NULL );
864 assert( propdata->leaders == NULL );
865
866 assert( propdata->permvardomaincenter == NULL );
867 assert( propdata->permvars == NULL );
868 assert( propdata->perms == NULL );
869 assert( propdata->permstrans == NULL );
870 assert( propdata->npermvars == 0 );
871 assert( propdata->nbinpermvars == 0 );
872 assert( propdata->nperms == -1 || propdata->nperms == 0 );
873 assert( propdata->nmaxperms == 0 );
874 assert( propdata->nmovedpermvars == -1 );
875 assert( propdata->nmovedbinpermvars == 0 );
876 assert( propdata->nmovedintpermvars == 0 );
877 assert( propdata->nmovedimplintpermvars == 0 );
878 assert( propdata->nmovedcontpermvars == 0 );
879 assert( propdata->nmovedvars == -1 );
880 assert( propdata->binvaraffected == FALSE );
881 assert( propdata->isnonlinvar == NULL );
882
883 assert( propdata->componenthassignedperm == NULL );
884 assert( propdata->componentblocked == NULL );
885 assert( propdata->componentbegins == NULL );
886 assert( propdata->components == NULL );
887 assert( propdata->ncomponents == -1 );
888 assert( propdata->ncompblocked == 0 );
889
890 return TRUE;
891}
892#endif
893
894
895/** resets symmetry handling propagators that depend on the branch-and-bound tree structure */
896static
898 SCIP* scip, /**< SCIP pointer */
899 SCIP_PROPDATA* propdata /**< propagator data */
900 )
901{
902 assert( scip != NULL );
903 assert( propdata != NULL );
904
905 /* propagators managed by a different file */
906 if ( propdata->orbitalreddata != NULL )
907 {
908 SCIP_CALL( SCIPorbitalReductionReset(scip, propdata->orbitalreddata) );
909 }
910 if ( propdata->orbitopalreddata != NULL )
911 {
912 SCIP_CALL( SCIPorbitopalReductionReset(scip, propdata->orbitopalreddata) );
913 }
914 if ( propdata->lexreddata != NULL )
915 {
916 SCIP_CALL( SCIPlexicographicReductionReset(scip, propdata->lexreddata) );
917 }
918
919 return SCIP_OKAY;
920}
921
922
923/** frees symmetry data */
924static
926 SCIP* scip, /**< SCIP pointer */
927 SCIP_PROPDATA* propdata /**< propagator data */
928 )
929{
930 int i;
931
932 assert( scip != NULL );
933 assert( propdata != NULL );
934
936
937 if ( propdata->permvarmap != NULL )
938 {
939 SCIPhashmapFree(&propdata->permvarmap);
940 }
941
942 /* release all variables contained in permvars array */
943 for (i = 0; i < propdata->npermvars; ++i)
944 {
945 assert( propdata->permvars[i] != NULL );
946 SCIP_CALL( SCIPreleaseVar(scip, &propdata->permvars[i]) );
947 }
948
949 /* free permstrans matrix*/
950 if ( propdata->permstrans != NULL )
951 {
952 assert( propdata->nperms > 0 );
953 assert( propdata->permvars != NULL );
954 assert( propdata->npermvars > 0 );
955 assert( propdata->nmaxperms > 0 );
956
957 for (i = 0; i < propdata->npermvars; ++i)
958 {
959 SCIPfreeBlockMemoryArray(scip, &propdata->permstrans[i], propdata->nmaxperms);
960 }
961 SCIPfreeBlockMemoryArray(scip, &propdata->permstrans, propdata->npermvars);
962 }
963
964 /* free data of added orbitope/orbisack/symresack constraints */
965 if ( propdata->genorbconss != NULL )
966 {
967 assert( propdata->ngenorbconss > 0 );
968
969 /* release constraints */
970 while ( propdata->ngenorbconss > 0 )
971 {
972 assert( propdata->genorbconss[propdata->ngenorbconss - 1] != NULL );
973 SCIP_CALL( SCIPreleaseCons(scip, &propdata->genorbconss[--propdata->ngenorbconss]) );
974 }
975 assert( propdata->ngenorbconss == 0 );
976
977 /* free pointers to symmetry group and binary variables */
978 SCIPfreeBlockMemoryArray(scip, &propdata->genorbconss, propdata->genorbconsssize);
979 propdata->genorbconsssize = 0;
980 }
981
982 /* free data of added constraints */
983 if ( propdata->genlinconss != NULL )
984 {
985 /* release constraints */
986 for (i = 0; i < propdata->ngenlinconss; ++i)
987 {
988 assert( propdata->genlinconss[i] != NULL );
989 SCIP_CALL( SCIPreleaseCons(scip, &propdata->genlinconss[i]) );
990 }
991
992 /* free pointers to symmetry group and binary variables */
993 SCIPfreeBlockMemoryArray(scip, &propdata->genlinconss, propdata->genlinconsssize);
994 propdata->ngenlinconss = 0;
995 propdata->genlinconsssize = 0;
996 }
997
998 if ( propdata->sstconss != NULL )
999 {
1000 assert( propdata->nsstconss > 0 );
1001
1002 /* release constraints */
1003 for (i = 0; i < propdata->nsstconss; ++i)
1004 {
1005 assert( propdata->sstconss[i] != NULL );
1006 SCIP_CALL( SCIPreleaseCons(scip, &propdata->sstconss[i]) );
1007 }
1008
1009 /* free pointers to symmetry group and binary variables */
1010 SCIPfreeBlockMemoryArray(scip, &propdata->sstconss, propdata->maxnsstconss);
1011 propdata->sstconss = NULL;
1012 propdata->nsstconss = 0;
1013 propdata->maxnsstconss = 0;
1014 }
1015
1016 if ( propdata->leaders != NULL )
1017 {
1018 assert( propdata->maxnleaders > 0 );
1019
1020 SCIPfreeBlockMemoryArray(scip, &propdata->leaders, propdata->maxnleaders);
1021 propdata->maxnleaders = 0;
1022 propdata->leaders = NULL;
1023 propdata->nleaders = 0;
1024 }
1025
1026 /* free components */
1027 if ( propdata->ncomponents > 0 )
1028 {
1029 assert( propdata->componentblocked != NULL );
1030 assert( propdata->vartocomponent != NULL );
1031 assert( propdata->componentbegins != NULL );
1032 assert( propdata->components != NULL );
1033
1034 SCIPfreeBlockMemoryArray(scip, &propdata->componenthassignedperm, propdata->ncomponents);
1035 SCIPfreeBlockMemoryArray(scip, &propdata->componentblocked, propdata->ncomponents);
1036 SCIPfreeBlockMemoryArray(scip, &propdata->vartocomponent, propdata->npermvars);
1037 SCIPfreeBlockMemoryArray(scip, &propdata->componentbegins, propdata->ncomponents + 1);
1038 SCIPfreeBlockMemoryArray(scip, &propdata->components, propdata->nperms);
1039
1040 propdata->ncomponents = -1;
1041 propdata->ncompblocked = 0;
1042 }
1043
1044 /* free main symmetry data */
1045 if ( propdata->nperms > 0 )
1046 {
1047 int permlen;
1048
1049 assert( propdata->permvars != NULL );
1050
1051 if ( (SYM_SYMTYPE) propdata->symtype == SYM_SYMTYPE_SIGNPERM )
1052 permlen = 2 * propdata->npermvars;
1053 else
1054 permlen = propdata->npermvars;
1055
1056 SCIPfreeBlockMemoryArray(scip, &propdata->permvars, propdata->npermvars);
1057 SCIPfreeBlockMemoryArray(scip, &propdata->permvardomaincenter, propdata->npermvars);
1058
1059 if ( propdata->perms != NULL )
1060 {
1061 for (i = 0; i < propdata->nperms; ++i)
1062 {
1063 SCIPfreeBlockMemoryArray(scip, &propdata->perms[i], permlen);
1064 }
1065 SCIPfreeBlockMemoryArray(scip, &propdata->perms, propdata->nmaxperms);
1066 }
1067
1068 SCIPfreeBlockMemoryArrayNull(scip, &propdata->isnonlinvar, propdata->npermvars);
1069
1070 propdata->npermvars = 0;
1071 propdata->nbinpermvars = 0;
1072 propdata->nperms = -1;
1073 propdata->nmaxperms = 0;
1074 propdata->nmovedpermvars = -1;
1075 propdata->nmovedbinpermvars = 0;
1076 propdata->nmovedintpermvars = 0;
1077 propdata->nmovedimplintpermvars = 0;
1078 propdata->nmovedcontpermvars = 0;
1079 propdata->nmovedvars = -1;
1080 propdata->log10groupsize = -1.0;
1081 propdata->binvaraffected = FALSE;
1082 propdata->isnonlinvar = NULL;
1083 }
1084 propdata->nperms = -1;
1085
1086 assert( checkSymmetryDataFree(propdata) );
1087
1088 propdata->computedsymmetry = FALSE;
1089 propdata->compressed = FALSE;
1090
1091 return SCIP_OKAY;
1092}
1093
1094
1095/** makes sure that the constraint array (potentially NULL) of given array size is sufficiently large */
1096static
1098 SCIP* scip, /**< SCIP pointer */
1099 SCIP_CONS*** consarrptr, /**< constraint array pointer */
1100 int* consarrsizeptr, /**< constraint array size pointer */
1101 int consarrsizereq /**< constraint array size required */
1102 )
1103{
1104 int newsize;
1105
1106 assert( scip != NULL );
1107 assert( consarrptr != NULL );
1108 assert( consarrsizeptr != NULL );
1109 assert( consarrsizereq > 0 );
1110 assert( *consarrsizeptr >= 0 );
1111 assert( (*consarrsizeptr == 0) == (*consarrptr == NULL) );
1112
1113 /* array is already sufficiently large */
1114 if ( consarrsizereq <= *consarrsizeptr )
1115 return SCIP_OKAY;
1116
1117 /* compute new size */
1118 newsize = SCIPcalcMemGrowSize(scip, consarrsizereq);
1119 assert( newsize > *consarrsizeptr );
1120
1121 /* allocate or reallocate */
1122 if ( *consarrptr == NULL )
1123 {
1124 SCIP_CALL( SCIPallocBlockMemoryArray(scip, consarrptr, newsize) );
1125 }
1126 else
1127 {
1128 SCIP_CALL( SCIPreallocBlockMemoryArray(scip, consarrptr, *consarrsizeptr, newsize) );
1129 }
1130
1131 *consarrsizeptr = newsize;
1132
1133 return SCIP_OKAY;
1134}
1135
1136/** set symmetry data */
1137static
1139 SCIP* scip, /**< SCIP pointer */
1140 SYM_SYMTYPE symtype, /**< type of symmetries in perms */
1141 SCIP_VAR** vars, /**< vars present at time of symmetry computation */
1142 int nvars, /**< number of vars present at time of symmetry computation */
1143 int nbinvars, /**< number of binary vars present at time of symmetry computation */
1144 SCIP_VAR*** permvars, /**< pointer to permvars array */
1145 int* npermvars, /**< pointer to store number of permvars */
1146 int* nbinpermvars, /**< pointer to store number of binary permvars */
1147 SCIP_Real** permvardomaincenter, /**< pointer to store center points of variable domains */
1148 int** perms, /**< permutations matrix (nperms x nvars) */
1149 int nperms, /**< number of permutations */
1150 int* nmovedvars, /**< pointer to store number of vars affected by symmetry (if usecompression) or NULL */
1151 SCIP_Bool* binvaraffected, /**< pointer to store whether a binary variable is affected by symmetry */
1152 SCIP_Bool usecompression, /**< whether symmetry data shall be compressed */
1153 SCIP_Real compressthreshold, /**< if percentage of moved vars is at most threshold, compression is done */
1154 SCIP_Bool* compressed /**< pointer to store whether compression has been performed */
1155 )
1156{
1157 SCIP_Real ub;
1158 SCIP_Real lb;
1159 int i;
1160 int p;
1161
1162 assert( scip != NULL );
1163 assert( vars != NULL );
1164 assert( nvars > 0 );
1165 assert( permvars != NULL );
1166 assert( npermvars != NULL );
1167 assert( nbinpermvars != NULL );
1168 assert( perms != NULL );
1169 assert( nperms > 0 );
1170 assert( binvaraffected != NULL );
1171 assert( SCIPisGE(scip, compressthreshold, 0.0) );
1172 assert( SCIPisLE(scip, compressthreshold, 1.0) );
1173 assert( compressed != NULL );
1174
1175 /* set default return values */
1176 *permvars = vars;
1177 *npermvars = nvars;
1178 *nbinpermvars = nbinvars;
1179 *binvaraffected = FALSE;
1180 *compressed = FALSE;
1181
1182 /* if we possibly perform compression */
1183 if ( usecompression && SCIPgetNVars(scip) >= COMPRESSNVARSLB )
1184 {
1185 SCIP_Real percentagemovedvars;
1186 int* labelmovedvars;
1187 int* labeltopermvaridx;
1188 int nbinvarsaffected = 0;
1189
1190 assert( nmovedvars != NULL );
1191
1192 *nmovedvars = 0;
1193
1194 /* detect number of moved vars and label moved vars */
1195 SCIP_CALL( SCIPallocBufferArray(scip, &labelmovedvars, nvars) );
1196 SCIP_CALL( SCIPallocBufferArray(scip, &labeltopermvaridx, nvars) );
1197 for (i = 0; i < nvars; ++i)
1198 {
1199 labelmovedvars[i] = -1;
1200
1201 for (p = 0; p < nperms; ++p)
1202 {
1203 if ( perms[p][i] != i )
1204 {
1205 labeltopermvaridx[*nmovedvars] = i;
1206 labelmovedvars[i] = (*nmovedvars)++;
1207
1208 if ( SCIPvarIsBinary(vars[i]) )
1209 ++nbinvarsaffected;
1210 break;
1211 }
1212 }
1213 }
1214
1215 if ( nbinvarsaffected > 0 )
1216 *binvaraffected = TRUE;
1217
1218 /* check whether compression should be performed */
1219 percentagemovedvars = (SCIP_Real) *nmovedvars / (SCIP_Real) nvars;
1220 if ( *nmovedvars > 0 && SCIPisLE(scip, percentagemovedvars, compressthreshold) )
1221 {
1222 /* remove variables from permutations that are not affected by any permutation */
1223 for (p = 0; p < nperms; ++p)
1224 {
1225 /* iterate over labels and adapt permutation (possibly taking signed permutations into account) */
1226 for (i = 0; i < *nmovedvars; ++i)
1227 {
1228 assert( i <= labeltopermvaridx[i] );
1229 if ( perms[p][labeltopermvaridx[i]] >= nvars )
1230 {
1231 int imgvaridx;
1232
1233 assert( symtype == SYM_SYMTYPE_SIGNPERM );
1234
1235 imgvaridx = perms[p][labeltopermvaridx[i]] - nvars;
1236 perms[p][i] = labelmovedvars[imgvaridx] + *nmovedvars;
1237
1238 assert( 0 <= perms[p][i] && perms[p][i] < 2 * (*nmovedvars) );
1239 }
1240 else
1241 perms[p][i] = labelmovedvars[perms[p][labeltopermvaridx[i]]];
1242 }
1243
1244 if ( symtype == SYM_SYMTYPE_SIGNPERM )
1245 {
1246 SCIP_CALL( SCIPreallocBlockMemoryArray(scip, &perms[p], 2 * nvars, 2 * (*nmovedvars)) );
1247 }
1248 else
1249 {
1250 SCIP_CALL( SCIPreallocBlockMemoryArray(scip, &perms[p], nvars, *nmovedvars) );
1251 }
1252 }
1253
1254 /* remove variables from permvars array that are not affected by any symmetry */
1255 SCIP_CALL( SCIPallocBlockMemoryArray(scip, permvars, *nmovedvars) );
1256 for (i = 0; i < *nmovedvars; ++i)
1257 {
1258 (*permvars)[i] = vars[labeltopermvaridx[i]];
1259 }
1260 *npermvars = *nmovedvars;
1261 *nbinpermvars = nbinvarsaffected;
1262 *compressed = TRUE;
1263
1264 SCIPfreeBlockMemoryArray(scip, &vars, nvars);
1265 }
1266 SCIPfreeBufferArray(scip, &labeltopermvaridx);
1267 SCIPfreeBufferArray(scip, &labelmovedvars);
1268 }
1269 else
1270 {
1271 /* detect whether binary variable is affected by symmetry and count number of binary permvars */
1272 for (i = 0; i < nbinvars; ++i)
1273 {
1274 for (p = 0; p < nperms && ! *binvaraffected; ++p)
1275 {
1276 if ( perms[p][i] != i )
1277 {
1278 if ( SCIPvarIsBinary(vars[i]) )
1279 *binvaraffected = TRUE;
1280 break;
1281 }
1282 }
1283 }
1284 }
1285
1286 /* store center points of variable domains */
1287 SCIP_CALL( SCIPallocBlockMemoryArray(scip, permvardomaincenter, *npermvars) );
1288 for (i = 0; i < *npermvars; ++i)
1289 {
1290 ub = SCIPvarGetUbGlobal((*permvars)[i]);
1291 lb = SCIPvarGetLbGlobal((*permvars)[i]);
1292
1293 (*permvardomaincenter)[i] = 0.5 * (ub + lb);
1294 }
1295
1296 return SCIP_OKAY;
1297}
1298
1299/** returns whether a constraint handler can provide required symmetry information */
1300static
1302 SCIP_CONSHDLR* conshdlr, /**< constraint handler */
1303 SYM_SYMTYPE symtype /**< type of symmetries for which information are needed */
1304 )
1305{
1306 assert( conshdlr != NULL );
1307
1308 switch ( symtype )
1309 {
1310 case SYM_SYMTYPE_PERM:
1311 return SCIPconshdlrSupportsPermsymDetection(conshdlr);
1312 default:
1313 assert( symtype == SYM_SYMTYPE_SIGNPERM );
1315 } /*lint !e788*/
1316}
1317
1318/** returns whether all constraint handlers with constraints can provide symmetry information */
1319static
1321 SCIP* scip, /**< SCIP pointer */
1322 SYM_SYMTYPE symtype /**< type of symmetries for which information are needed */
1323 )
1324{
1325 SCIP_CONSHDLR** conshdlrs;
1326 SCIP_CONSHDLR* conshdlr;
1327 int nconshdlrs;
1328 int c;
1329
1330 conshdlrs = SCIPgetConshdlrs(scip);
1331 assert( conshdlrs != NULL );
1332
1333 nconshdlrs = SCIPgetNConshdlrs(scip);
1334 for (c = 0; c < nconshdlrs; ++c)
1335 {
1336 conshdlr = conshdlrs[c];
1337 assert( conshdlr != NULL );
1338
1339 if ( ! conshdlrCanProvideSymInformation(conshdlr, symtype) && SCIPconshdlrGetNConss(conshdlr) > 0 )
1340 {
1341 char name[SCIP_MAXSTRLEN];
1342
1343 if ( symtype == SYM_SYMTYPE_PERM )
1344 (void) SCIPsnprintf(name, SCIP_MAXSTRLEN, "CONSGETPERMSYMGRAPH");
1345 else
1346 (void) SCIPsnprintf(name, SCIP_MAXSTRLEN, "CONSGETSIGNEDPERMSYMGRAPH");
1347
1349 " Symmetry detection interrupted: constraints of type %s do not provide symmetry information.\n"
1350 " If symmetries shall be detected, implement the %s callback.\n",
1351 SCIPconshdlrGetName(conshdlr), name);
1352
1353 return FALSE;
1354 }
1355 }
1356
1357 /* check whether all expressions provide sufficient symmetry information */
1358 conshdlr = SCIPfindConshdlr(scip, "nonlinear");
1359 if ( conshdlr != NULL && SCIPconshdlrGetNConss(conshdlr) > 0 )
1360 {
1361 SCIP_EXPRHDLR* exprhdlr;
1362
1363 for (c = 0; c < SCIPgetNExprhdlrs(scip); ++c)
1364 {
1365 SCIP_Bool found = FALSE;
1366 exprhdlr = SCIPgetExprhdlrs(scip)[c];
1367
1368 if ( SCIPexprhdlrHasGetSymData(exprhdlr) )
1369 continue;
1370
1371 /* check whether exprhdlr is known by SCIP (and handles symmetries correctly) */
1372 if ( strcmp(SCIPexprhdlrGetName(exprhdlr), "var") == 0
1373 || strcmp(SCIPexprhdlrGetName(exprhdlr), "sum") == 0
1374 || strcmp(SCIPexprhdlrGetName(exprhdlr), "product") == 0
1375 || strcmp(SCIPexprhdlrGetName(exprhdlr), "val") == 0
1376 || strcmp(SCIPexprhdlrGetName(exprhdlr), "pow") == 0
1377 || strcmp(SCIPexprhdlrGetName(exprhdlr), "signpow") == 0
1378 || strcmp(SCIPexprhdlrGetName(exprhdlr), "exp") == 0
1379 || strcmp(SCIPexprhdlrGetName(exprhdlr), "log") == 0
1380 || strcmp(SCIPexprhdlrGetName(exprhdlr), "abs") == 0
1381 || strcmp(SCIPexprhdlrGetName(exprhdlr), "sin") == 0
1382 || strcmp(SCIPexprhdlrGetName(exprhdlr), "cos") == 0
1383 || strcmp(SCIPexprhdlrGetName(exprhdlr), "entropy") == 0
1384 || strcmp(SCIPexprhdlrGetName(exprhdlr), "erf") == 0
1385 || strcmp(SCIPexprhdlrGetName(exprhdlr), "varidx") == 0 )
1386 found = TRUE;
1387
1388 /* there exists an unknown expression handler that does not provide symmetry information */
1389 if ( ! found )
1390 {
1391 SCIPwarningMessage(scip, "Expression handler %s does not implement the EXPRGETSYMDATA callback.\n"
1392 "Computed symmetries might be incorrect if the expression uses different constants or assigns\n"
1393 "different coefficients to its children.\n", SCIPexprhdlrGetName(SCIPgetExprhdlrs(scip)[c]));
1394 }
1395 }
1396 }
1397
1398 return TRUE;
1399}
1400
1401/** provides estimates for the number of nodes and edges in a symmetry detection graph */
1402static
1404 SCIP* scip, /**< SCIP pointer */
1405 int* nopnodes, /**< pointer to store estimate for number of operator nodes */
1406 int* nvalnodes, /**< pointer to store estimate for number of value nodes */
1407 int* nconsnodes, /**< pointer to store estimate for number of constraint nodes */
1408 int* nedges /**< pointer to store estimate for number of edges */
1409 )
1410{
1411 SCIP_CONS** conss;
1412 SCIP_Bool success;
1413 int nvars;
1414 int nconss;
1415 int num;
1416 int c;
1417
1418 assert( scip != NULL );
1419 assert( nopnodes != NULL );
1420 assert( nvalnodes != NULL );
1421 assert( nconsnodes != NULL );
1422 assert( nedges != NULL );
1423
1424 nvars = SCIPgetNVars(scip);
1425 nconss = SCIPgetNConss(scip);
1426 conss = SCIPgetConss(scip);
1427 assert( conss != NULL || nconss == 0 );
1428
1429 *nconsnodes = nconss;
1430
1431 /* get estimate from different types of constraints */
1432 *nopnodes = 0;
1433 *nvalnodes = 0;
1434 for (c = 0; c < nconss; ++c)
1435 {
1436 if ( strcmp(SCIPconshdlrGetName(SCIPconsGetHdlr(conss[c])), "bounddisjunction") == 0 )
1437 {
1438 SCIP_CALL( SCIPgetConsNVars(scip, conss[c], &num, &success) );
1439
1440 if ( success )
1441 {
1442 *nopnodes += num;
1443 *nvalnodes += num;
1444 }
1445 }
1446 else if ( strcmp(SCIPconshdlrGetName(SCIPconsGetHdlr(conss[c])), "indicator") == 0 )
1447 {
1448 *nopnodes += 3;
1449 *nvalnodes += 1;
1450 }
1451 else if ( strcmp(SCIPconshdlrGetName(SCIPconsGetHdlr(conss[c])), "nonlinear") == 0 )
1452 {
1453 SCIP_CALL( SCIPgetConsNVars(scip, conss[c], &num, &success) );
1454
1455 /* use binary trees as a proxy for an expression tree */
1456 if ( success )
1457 {
1458 int depth;
1459 int numnodes;
1460 int expval;
1461
1462 depth = (int) log2((double) num);
1463 expval = (int) exp2((double) (depth + 1));
1464 numnodes = MIN(expval, 100);
1465
1466 *nopnodes += numnodes;
1467 *nvalnodes += MAX((int) 0.1 * numnodes, 1);
1468 }
1469 }
1470 else if ( strcmp(SCIPconshdlrGetName(SCIPconsGetHdlr(conss[c])), "SOS1") == 0 )
1471 {
1472 SCIP_CALL( SCIPgetConsNVars(scip, conss[c], &num, &success) );
1473
1474 if ( success )
1475 *nopnodes += num;
1476 }
1477 else if ( strcmp(SCIPconshdlrGetName(SCIPconsGetHdlr(conss[c])), "SOS2") == 0 )
1478 {
1479 SCIP_CALL( SCIPgetConsNVars(scip, conss[c], &num, &success) );
1480
1481 if ( success )
1482 *nopnodes += num - 1;
1483 }
1484 }
1485
1486 /* use a staggered scheme for the number of edges since this can become large
1487 *
1488 * In most cases, edges represent variable coefficients from linear constraints.
1489 * For this reason, use number of variables as proxy.
1490 */
1491 if ( nvars <= 100000 )
1492 *nedges = 100 * nvars;
1493 else if ( nvars <= 1000000 )
1494 *nedges = 32 * nvars;
1495 else if ( nvars <= 16700000 )
1496 *nedges = 16 * nvars;
1497 else
1498 *nedges = INT_MAX / 10;
1499
1500 return SCIP_OKAY;
1501}
1502
1503/** checks whether computed symmetries are indeed symmetries */
1504static
1506 SCIP* scip, /**< SCIP pointer */
1507 SYM_SYMTYPE symtype, /**< type of symmetries to be checked */
1508 int** perms, /**< array of permutations */
1509 int nperms, /**< number of permutations */
1510 int npermvars, /**< number of variables permutations act on */
1511 SYM_SPEC fixedtype /**< variable types that must be fixed by symmetries */
1512 )
1513{
1514 SYM_GRAPH** graphs;
1515 SCIP_CONS** conss;
1516 SCIP_VAR** symvars;
1517 SCIP_Bool success;
1518 int* graphperm;
1519 int* groupbegins;
1520 int ngroups = 1;
1521 int nsymvars;
1522 int nconss;
1523 int p;
1524 int c;
1525 int g;
1526#ifdef SCIP_DISPLAY_SYM_CHECK
1527 int permlen;
1528 SCIP_Bool* covered;
1529#endif
1530
1531 assert( scip != NULL );
1532 assert( perms != NULL );
1533 assert( nperms > 0 );
1534 assert( npermvars > 0 );
1535
1536 /* get symmetry detection graphs for all constraints */
1537 nconss = SCIPgetNConss(scip);
1538 conss = SCIPgetConss(scip);
1539 assert( conss != NULL );
1540
1541 symvars = SCIPgetVars(scip);
1542 nsymvars = SCIPgetNVars(scip);
1543 assert( nsymvars == npermvars );
1544
1545 SCIP_CALL( SCIPallocBufferArray(scip, &graphs, nconss) );
1546
1547 for (c = 0; c < nconss; ++c)
1548 {
1549 SCIP_CALL( SCIPcreateSymgraph(scip, symtype, &graphs[c], symvars, nsymvars, 10, 10, 1, 100) );
1550
1551 success = FALSE;
1552 switch ( symtype )
1553 {
1554 case SYM_SYMTYPE_PERM:
1555 SCIP_CALL( SCIPgetConsPermsymGraph(scip, conss[c], graphs[c], &success) );
1556 break;
1557 default:
1558 assert( symtype == SYM_SYMTYPE_SIGNPERM );
1559 SCIP_CALL( SCIPgetConsSignedPermsymGraph(scip, conss[c], graphs[c], &success) );
1560 } /*lint !e788*/
1561
1562 SCIP_CALL( SCIPcomputeSymgraphColors(scip, graphs[c], fixedtype) );
1563
1564 assert( success );
1565 }
1566
1567 /* sort graphs for quicker comparisons */
1568 SCIP_CALL( SCIPallocBufferArray(scip, &graphperm, nconss) );
1569 SCIP_CALL( SCIPallocBufferArray(scip, &groupbegins, nconss + 1) );
1570 for (c = 0; c < nconss; ++c)
1571 {
1573 }
1574
1575 SCIPsort(graphperm, SYMsortSymgraphs, graphs, nconss);
1576
1577 groupbegins[0] = 0;
1578 for (c = 1; c < nconss; ++c)
1579 {
1580 if ( compareSymgraphs(scip, graphs[graphperm[c]], graphs[graphperm[c-1]]) != 0 )
1581 groupbegins[ngroups++] = c;
1582 }
1583 groupbegins[ngroups] = nconss;
1584
1585 /* remove information from symmetry detection graph that is not needed anymore */
1586 for (c = 0; c < nconss; ++c)
1587 {
1589 }
1590
1591#ifdef SCIP_DISPLAY_SYM_CHECK
1592 permlen = symtype == SYM_SYMTYPE_SIGNPERM ? 2 * npermvars : npermvars;
1593 SCIP_CALL( SCIPallocClearBufferArray(scip, &covered, permlen) );
1594#endif
1595
1596 /* iterate over all permutations and check whether they define symmetries */
1597 for (p = 0; p < nperms; ++p)
1598 {
1599 SYM_GRAPH* graph;
1600 SCIP_Bool found = TRUE;
1601 int d;
1602#ifdef SCIP_DISPLAY_SYM_CHECK
1603 int i;
1604
1605 SCIPinfoMessage(scip, NULL, "Check whether permutation %d is a symmetry:\n", p);
1606 for (i = 0; i < permlen; ++i)
1607 {
1608 SCIP_CALL( displayCycleOfSymmetry(scip, perms[p], symtype, i, covered, npermvars, SCIPgetVars(scip)) );
1609 }
1610
1611 for (i = 0; i < permlen; ++i)
1612 covered[i] = FALSE;
1613 SCIPinfoMessage(scip, NULL, "Check whether every constraint has a symmetric counterpart.\n");
1614#endif
1615
1616 /* for every constraint, create permuted graph by copying nodes and edges */
1617 for (g = 0; g < ngroups; ++g)
1618 {
1619 for (c = groupbegins[g]; c < groupbegins[g+1]; ++c)
1620 {
1621#ifdef SCIP_DISPLAY_SYM_CHECK
1622 SCIPinfoMessage(scip, NULL, "Check whether constraint %d has a symmetric counterpart:\n",
1623 graphperm[c]);
1624 SCIP_CALL( SCIPprintCons(scip, conss[graphperm[c]], NULL) );
1625 SCIPinfoMessage(scip, NULL, "\n");
1626#endif
1627 SCIP_CALL( SCIPcopySymgraph(scip, &graph, graphs[graphperm[c]], perms[p], fixedtype) );
1628
1629 /* if adapted graph is equivalent to original graph, we don't need to check further graphs */
1630 if ( SYMcheckGraphsAreIdentical(scip, symtype, graph, graphs[graphperm[c]]) )
1631 {
1632#ifdef SCIP_DISPLAY_SYM_CHECK
1633 SCIPinfoMessage(scip, NULL, "\tconstraint is symmetric to itself\n");
1634#endif
1635 SCIP_CALL( SCIPfreeSymgraph(scip, &graph) );
1636 continue;
1637 }
1638
1639 /* check whether graph has an isomorphic counterpart */
1640 found = FALSE;
1641 for (d = groupbegins[g]; d < groupbegins[g+1] && ! found; ++d)
1642 {
1643 found = SYMcheckGraphsAreIdentical(scip, symtype, graph, graphs[graphperm[d]]);
1644
1645#ifdef SCIP_DISPLAY_SYM_CHECK
1646 SCIPinfoMessage(scip, NULL, "\tconstraint is %ssymmetric to constraint %d\n\t", !found ? "not " : "", d);
1647 SCIP_CALL( SCIPprintCons(scip, conss[graphperm[d]], NULL) );
1648 SCIPinfoMessage(scip, NULL, "\n");
1649#endif
1650 }
1651
1652 SCIP_CALL( SCIPfreeSymgraph(scip, &graph) );
1653
1654 if ( ! found )
1655 {
1656#ifdef SCIP_DISPLAY_SYM_CHECK
1657 SCIPfreeBufferArray(scip, &covered);
1658#endif
1659 SCIPerrorMessage("permutation %d is not a symmetry\n", p);
1660 return SCIP_ERROR;
1661 }
1662 }
1663 }
1664 }
1665
1666#ifdef SCIP_DISPLAY_SYM_CHECK
1667 SCIPfreeBufferArray(scip, &covered);
1668#endif
1669
1670 SCIPfreeBufferArray(scip, &groupbegins);
1671 SCIPfreeBufferArray(scip, &graphperm);
1672
1673 for (c = nconss - 1; c >= 0; --c)
1674 {
1675 SCIP_CALL( SCIPfreeSymgraph(scip, &graphs[c]) );
1676 }
1677 SCIPfreeBufferArray(scip, &graphs);
1678
1679 return SCIP_OKAY;
1680}
1681
1682/** computes symmetry group of a CIP */
1683static
1685 SCIP* scip, /**< SCIP pointer */
1686 SYM_SYMTYPE symtype, /**< type of symmetries to be computed */
1687 SCIP_Bool compresssymmetries, /**< Should non-affected variables be removed from permutation to save memory? */
1688 SCIP_Real compressthreshold, /**< if percentage of moved vars is at most threshold, compression is done */
1689 int maxgenerators, /**< maximal number of generators constructed (= 0 if unlimited) */
1690 SYM_SPEC fixedtype, /**< variable types that must be fixed by symmetries */
1691 SCIP_Bool checksymmetries, /**< Should all symmetries be checked after computation? */
1692 SCIP_VAR*** permvars, /**< pointer to permvars array */
1693 int* npermvars, /**< pointer to store number of permvars */
1694 int* nbinpermvars, /**< pointer to store number of binary permvars */
1695 SCIP_Real** permvardomaincenter, /**< pointer to store center points of variable domains */
1696 int*** perms, /**< pointer to store permutation matrix (nperms x nvars) */
1697 int* nperms, /**< pointer to store number of permutations */
1698 int* nmaxperms, /**< pointer to store maximal number of permutations
1699 * (needed for freeing storage) */
1700 int* nmovedvars, /**< pointer to store number of vars affected
1701 * by symmetry (if usecompression) or NULL */
1702 SCIP_Bool* binvaraffected, /**< pointer to store whether a binary variable is affected by symmetry */
1703 SCIP_Bool* compressed, /**< pointer to store whether compression has been performed */
1704 SCIP_Real* log10groupsize, /**< pointer to store log10 of size of group */
1705 SCIP_Real* symcodetime, /**< pointer to store the time for symmetry code */
1706 SCIP_Bool* success /**< pointer to store whether symmetry computation was successful */
1707 )
1708{
1709 SCIP_CONS** conss;
1710 SYM_GRAPH* graph;
1711 int nconsnodes = 0;
1712 int nvalnodes = 0;
1713 int nopnodes = 0;
1714 int nedges = 0;
1715 int nconss;
1716 int c;
1717
1718 assert( scip != NULL );
1719 assert( permvars != NULL );
1720 assert( npermvars != NULL );
1721 assert( nbinpermvars != NULL );
1722 assert( perms != NULL );
1723 assert( nperms != NULL );
1724 assert( nmaxperms != NULL );
1725 assert( nmovedvars != NULL );
1726 assert( binvaraffected != NULL );
1727 assert( compressed != NULL );
1728 assert( log10groupsize != NULL );
1729 assert( symcodetime != NULL );
1730 assert( success != NULL );
1731
1732 /* init pointers */
1733 *permvars = NULL;
1734 *npermvars = 0;
1735 *nbinpermvars = 0;
1736 *perms = NULL;
1737 *nperms = 0;
1738 *nmaxperms = 0;
1739 *nmovedvars = -1;
1740 *binvaraffected = FALSE;
1741 *compressed = FALSE;
1742 *log10groupsize = 0;
1743 *success = FALSE;
1744 *symcodetime = 0.0;
1745
1746 /* check whether all constraints can provide symmetry information */
1747 if ( ! conshdlrsCanProvideSymInformation(scip, symtype) )
1748 return SCIP_OKAY;
1749
1750 /* get symmetry detection graphs from constraints */
1751 conss = SCIPgetConss(scip);
1752 nconss = SCIPgetNConss(scip);
1753
1754 assert( conss != NULL || nconss == 0 );
1755
1756 /* exit if no constraints or no variables are available */
1757 if ( nconss == 0 || SCIPgetNVars(scip) == 0 )
1758 {
1759 *success = TRUE;
1760 return SCIP_OKAY;
1761 }
1762
1763 /* get an estimate for the number of nodes and edges */
1764 SCIP_CALL( estimateSymgraphSize(scip, &nopnodes, &nvalnodes, &nconsnodes, &nedges) );
1765
1766 /* create graph */
1768 nopnodes, nvalnodes, nconsnodes, nedges) );
1769
1770 *success = TRUE;
1771 for (c = 0; c < nconss && *success; ++c)
1772 {
1773 if ( symtype == SYM_SYMTYPE_PERM )
1774 {
1775 SCIP_CALL( SCIPgetConsPermsymGraph(scip, conss[c], graph, success) );
1776 }
1777 else
1778 {
1779 assert( symtype == SYM_SYMTYPE_SIGNPERM );
1780 SCIP_CALL( SCIPgetConsSignedPermsymGraph(scip, conss[c], graph, success) );
1781 }
1782
1783 /* terminate early if graph could not be returned */
1784 if ( ! *success )
1785 {
1786 SCIP_CALL( SCIPfreeSymgraph(scip, &graph) );
1787
1788 return SCIP_OKAY;
1789 }
1790 }
1791
1792 SCIP_CALL( SCIPcomputeSymgraphColors(scip, graph, fixedtype) );
1793
1794 /* terminate early in case all variables are different */
1795 if ( (symtype == SYM_SYMTYPE_PERM && SCIPgetSymgraphNVarcolors(graph) == SCIPgetNVars(scip))
1796 || (symtype == SYM_SYMTYPE_SIGNPERM && SCIPgetSymgraphNVarcolors(graph) == 2 * SCIPgetNVars(scip)) )
1797 {
1798 SCIP_CALL( SCIPfreeSymgraph(scip, &graph) );
1799 return SCIP_OKAY;
1800 }
1801
1802 /*
1803 * actually compute symmetries
1804 */
1805 SCIP_CALL( SYMcomputeSymmetryGenerators(scip, maxgenerators, graph, nperms, nmaxperms,
1806 perms, log10groupsize, symcodetime) );
1807
1808 if ( checksymmetries && *nperms > 0 )
1809 {
1810 SCIP_CALL( checkSymmetriesAreSymmetries(scip, symtype, *perms, *nperms, SCIPgetNVars(scip), fixedtype) );
1811 }
1812
1813 /* potentially store symmetries */
1814 if ( *nperms > 0 )
1815 {
1816 SCIP_VAR** vars;
1817 int nvars;
1818
1819 nvars = SCIPgetNVars(scip);
1820 SCIP_CALL( SCIPduplicateBlockMemoryArray(scip, &vars, SCIPgetVars(scip), nvars) ); /*lint !e666*/
1821
1822 SCIP_CALL( setSymmetryData(scip, symtype, vars, nvars, SCIPgetNBinVars(scip), permvars, npermvars, nbinpermvars,
1823 permvardomaincenter, *perms, *nperms, nmovedvars, binvaraffected,
1824 compresssymmetries, compressthreshold, compressed) );
1825 }
1826
1827 /* free symmetry graph */
1828 SCIP_CALL( SCIPfreeSymgraph(scip, &graph) );
1829
1830 return SCIP_OKAY;
1831}
1832
1833/** returns whether a symmetry is a non-standard permutation */
1834static
1836 SCIP* scip, /**< SCIP instance */
1837 int* symmetry, /**< a symmetry encoded as a signed permutation */
1838 SCIP_VAR** vars, /**< array of variables the symmetry acts on */
1839 int nvars /**< number of variables in vars */
1840 )
1841{
1842 int v;
1843
1844 assert( symmetry != NULL );
1845 assert( vars != NULL );
1846 assert( nvars > 0 );
1847
1848 for (v = 0; v < nvars; ++v)
1849 {
1850 /* the symmetry is signed */
1851 if ( symmetry[v] >= nvars )
1852 return TRUE;
1853
1854 /* the domain of symmetric variables is different */
1855 if ( !SCIPisEQ(scip, SCIPvarGetLbLocal(vars[v]), SCIPvarGetLbLocal(vars[symmetry[v]]))
1856 || !SCIPisEQ(scip, SCIPvarGetUbLocal(vars[v]), SCIPvarGetUbLocal(vars[symmetry[v]])) )
1857 {
1858 assert( SCIPisEQ(scip, SCIPvarGetUbLocal(vars[v]) - SCIPvarGetLbLocal(vars[v]),
1859 SCIPvarGetUbLocal(vars[symmetry[v]]) - SCIPvarGetLbLocal(vars[symmetry[v]])) );
1860 return TRUE;
1861 }
1862 }
1863
1864 return FALSE;
1865}
1866
1867/** checks whether component contains non-standard permutations
1868 *
1869 * If all symmetries are standard permutations, stores them as such.
1870 */
1871static
1873 SCIP* scip, /**< SCIP instance */
1874 SCIP_PROPDATA* propdata /**< propagator data */
1875 )
1876{
1877 int* components;
1878 int* componentbegins;
1879 int ncomponents;
1880 int i;
1881 int c;
1882
1883 assert( scip != NULL );
1884 assert( propdata != NULL );
1885 assert( propdata->ncomponents > 0 );
1886 assert( propdata->components != NULL );
1887 assert( propdata->componentbegins != NULL );
1888
1889 components = propdata->components;
1890 componentbegins = propdata->componentbegins;
1891 ncomponents = propdata->ncomponents;
1892
1893 SCIP_CALL( SCIPallocClearBlockMemoryArray(scip, &(propdata->componenthassignedperm), ncomponents) );
1894
1895 /* stop if no non-standard permutations can exist */
1896 if ( (SYM_SYMTYPE) propdata->symtype == SYM_SYMTYPE_PERM )
1897 return SCIP_OKAY;
1898
1899 /* for each component, check whether it has a non-standard permutation */
1900 for (c = 0; c < ncomponents; ++c)
1901 {
1902 for (i = componentbegins[c]; i < componentbegins[c + 1]; ++i)
1903 {
1904 if ( isNonstandardPerm(scip, propdata->perms[components[i]], propdata->permvars, propdata->npermvars) )
1905 {
1906 propdata->componenthassignedperm[c] = TRUE;
1907 break;
1908 }
1909 }
1910 }
1911
1912 return SCIP_OKAY;
1913}
1914
1915/** ensures that the symmetry components are already computed */
1916static
1918 SCIP* scip, /**< SCIP instance */
1919 SCIP_PROPDATA* propdata /**< propagator data */
1920 )
1921{
1922 assert( scip != NULL );
1923 assert( propdata != NULL );
1924
1925 /* symmetries must have been determined */
1926 assert( propdata->nperms >= 0 );
1927
1928 /* stop if already computed */
1929 if ( propdata->ncomponents >= 0 )
1930 return SCIP_OKAY;
1931
1932 /* compute components */
1933 assert( propdata->ncomponents == -1 );
1934 assert( propdata->components == NULL );
1935 assert( propdata->componentbegins == NULL );
1936 assert( propdata->vartocomponent == NULL );
1937
1938#ifdef SCIP_OUTPUT_COMPONENT
1939 SCIPverbMessage(scip, SCIP_VERBLEVEL_HIGH, NULL, " (%.1fs) component computation started\n", SCIPgetSolvingTime(scip));
1940#endif
1941
1942 SCIP_CALL( SCIPcomputeComponentsSym(scip, (SYM_SYMTYPE) propdata->symtype, propdata->perms, propdata->nperms,
1943 propdata->permvars, propdata->npermvars, FALSE, &propdata->components, &propdata->componentbegins,
1944 &propdata->vartocomponent, &propdata->componentblocked, &propdata->ncomponents) );
1945
1946#ifdef SCIP_OUTPUT_COMPONENT
1947 SCIPverbMessage(scip, SCIP_VERBLEVEL_HIGH, NULL, " (%.1fs) component computation finished\n", SCIPgetSolvingTime(scip));
1948#endif
1949
1950 assert( propdata->components != NULL );
1951 assert( propdata->componentbegins != NULL );
1952 assert( propdata->ncomponents > 0 );
1953
1954 /* structure of symmetries can be simplified if they are standard permutations */
1956 assert( propdata->componenthassignedperm != NULL );
1957
1958 return SCIP_OKAY;
1959}
1960
1961
1962/** ensures that permvarmap is initialized */
1963static
1965 SCIP* scip, /**< SCIP instance */
1966 SCIP_PROPDATA* propdata /**< propagator data */
1967 )
1968{
1969 int v;
1970
1971 assert( scip != NULL );
1972 assert( propdata != NULL );
1973
1974 /* symmetries must have been determined */
1975 assert( propdata->nperms >= 0 );
1976
1977 /* stop if already computed */
1978 if ( propdata->permvarmap != NULL )
1979 return SCIP_OKAY;
1980
1981 /* create hashmap for storing the indices of variables */
1982 SCIP_CALL( SCIPhashmapCreate(&propdata->permvarmap, SCIPblkmem(scip), propdata->npermvars) );
1983
1984 /* insert variables into hashmap */
1985 for (v = 0; v < propdata->npermvars; ++v)
1986 {
1987 SCIP_CALL( SCIPhashmapInsertInt(propdata->permvarmap, propdata->permvars[v], v) );
1988 }
1989
1990 return SCIP_OKAY;
1991}
1992
1993
1994/** ensures that permstrans is initialized */
1995static
1997 SCIP* scip, /**< SCIP instance */
1998 SCIP_PROPDATA* propdata /**< propagator data */
1999 )
2000{
2001 int v;
2002 int p;
2003
2004 assert( scip != NULL );
2005 assert( propdata != NULL );
2006
2007 /* symmetries must have been determined */
2008 assert( propdata->nperms >= 0 );
2009
2010 /* stop if already computed */
2011 if ( propdata->permstrans != NULL )
2012 return SCIP_OKAY;
2013
2014 /* transpose symmetries matrix here */
2015 assert( propdata->permstrans == NULL );
2016 SCIP_CALL( SCIPallocBlockMemoryArray(scip, &propdata->permstrans, propdata->npermvars) );
2017 for (v = 0; v < propdata->npermvars; ++v)
2018 {
2019 SCIP_CALL( SCIPallocBlockMemoryArray(scip, &(propdata->permstrans[v]), propdata->nmaxperms) );
2020 for (p = 0; p < propdata->nperms; ++p)
2021 propdata->permstrans[v][p] = propdata->perms[p][v];
2022 }
2023
2024 return SCIP_OKAY;
2025}
2026
2027
2028/** ensures that movedpermvarscounts is initialized */
2029static
2031 SCIP* scip, /**< SCIP instance */
2032 SCIP_PROPDATA* propdata /**< propagator data */
2033 )
2034{
2035 int v;
2036 int p;
2037
2038 assert( scip != NULL );
2039 assert( propdata != NULL );
2040
2041 /* symmetries must have been determined */
2042 assert( propdata->nperms >= 0 );
2043
2044 /* stop if already computed */
2045 if ( propdata->nmovedpermvars >= 0 )
2046 return SCIP_OKAY;
2047 assert( propdata->nmovedpermvars == -1 );
2048
2049 propdata->nmovedpermvars = 0;
2050 propdata->nmovedbinpermvars = 0;
2051 propdata->nmovedintpermvars = 0;
2052 propdata->nmovedimplintpermvars = 0;
2053 propdata->nmovedcontpermvars = 0;
2054
2055 for (p = 0; p < propdata->nperms; ++p)
2056 {
2057 for (v = 0; v < propdata->npermvars; ++v)
2058 {
2059 if ( propdata->perms[p][v] != v )
2060 {
2061 ++propdata->nmovedpermvars;
2062
2063 switch ( SCIPvarGetType(propdata->permvars[v]) )
2064 {
2066 ++propdata->nmovedbinpermvars;
2067 break;
2069 ++propdata->nmovedintpermvars;
2070 break;
2072 ++propdata->nmovedimplintpermvars;
2073 break;
2075 ++propdata->nmovedcontpermvars;
2076 break;
2077 default:
2078 SCIPerrorMessage("Variable provided with unknown vartype\n");
2079 return SCIP_ERROR;
2080 }
2081 }
2082 }
2083 }
2084
2085 return SCIP_OKAY;
2086}
2087
2088
2089/** returns whether any allowed symmetry handling method is effective for the problem instance */
2090static
2092 SCIP* scip, /**< SCIP instance */
2093 SCIP_PROPDATA* propdata /**< propagator data */
2094 )
2095{
2096 /* must always compute symmetry if it is enforced */
2097 if ( propdata->enforcecomputesymmetry )
2098 return TRUE;
2099
2100 /* for dynamic symmetry handling or orbital reduction, branching must be possible */
2101 if ( propdata->usedynamicprop || ISORBITALREDUCTIONACTIVE(propdata->usesymmetry) )
2102 {
2103 /* @todo a proper test whether variables can be branched on or not */
2104 if ( SCIPgetNBinVars(scip) > 0 )
2105 return TRUE;
2106 if ( SCIPgetNIntVars(scip) > 0 )
2107 return TRUE;
2108 /* continuous variables can be branched on if nonlinear constraints exist */
2109 if ( ( SCIPgetNContVars(scip) > 0 || SCIPgetNImplVars(scip) > 0 )
2110 && SCIPconshdlrGetNActiveConss(propdata->conshdlr_nonlinear) > 0 )
2111 return TRUE;
2112 }
2113
2114 /* for SST, matching leadervartypes */
2115 if ( ISSSTACTIVE(propdata->usesymmetry) )
2116 {
2117 if ( ISSSTBINACTIVE(propdata->sstleadervartype) && SCIPgetNBinVars(scip) > 0 ) /*lint !e641*/
2118 return TRUE;
2119 if ( ISSSTINTACTIVE(propdata->sstleadervartype) && SCIPgetNIntVars(scip) > 0 ) /*lint !e641*/
2120 return TRUE;
2121 if ( ISSSTIMPLINTACTIVE(propdata->sstleadervartype) && SCIPgetNImplVars(scip) > 0 ) /*lint !e641*/
2122 return TRUE;
2123 if ( ISSSTCONTACTIVE(propdata->sstleadervartype) && SCIPgetNContVars(scip) > 0 ) /*lint !e641*/
2124 return TRUE;
2125 }
2126
2127 /* for static symmetry handling constraints, binary variables must be present */
2128 if ( ISSYMRETOPESACTIVE(propdata->usesymmetry) )
2129 {
2130 if ( SCIPgetNBinVars(scip) > 0 )
2131 return TRUE;
2132 }
2133
2134 /* if all tests above fail, then the symmetry handling methods cannot achieve anything */
2135 return FALSE;
2136}
2137
2138/** determines symmetry */
2139static
2141 SCIP* scip, /**< SCIP instance */
2142 SCIP_PROPDATA* propdata, /**< propagator data */
2143 SYM_SPEC symspecrequire, /**< symmetry specification for which we need to compute symmetries */
2144 SYM_SPEC symspecrequirefixed /**< symmetry specification of variables which must be fixed by symmetries */
2145 )
2146{ /*lint --e{641}*/
2147 SCIP_Bool successful;
2148 SCIP_Real symcodetime = 0.0;
2149 int maxgenerators;
2150 unsigned int type = 0;
2151 int nvars;
2152 int i;
2153
2154 assert( scip != NULL );
2155 assert( propdata != NULL );
2156 assert( propdata->usesymmetry >= 0 );
2157
2158 /* do not compute symmetry if reoptimization is enabled */
2159 if ( SCIPisReoptEnabled(scip) )
2160 return SCIP_OKAY;
2161
2162 /* do not compute symmetry if Benders decomposition enabled */
2163 if ( SCIPgetNActiveBenders(scip) > 0 )
2164 return SCIP_OKAY;
2165
2166 /* skip symmetry computation if no graph automorphism code was linked */
2167 if ( ! SYMcanComputeSymmetry() )
2168 {
2170 " Deactivated symmetry handling methods, since SCIP was built without symmetry detector (SYM=none).\n");
2171
2172 return SCIP_OKAY;
2173 }
2174
2175 /* do not compute symmetry if there are active pricers */
2176 if ( SCIPgetNActivePricers(scip) > 0 )
2177 return SCIP_OKAY;
2178
2179 /* avoid trivial cases */
2180 nvars = SCIPgetNVars(scip);
2181 if ( nvars <= 0 )
2182 return SCIP_OKAY;
2183
2184 /* do not compute symmetry if we cannot handle it */
2185 if ( !testSymmetryComputationRequired(scip, propdata) )
2186 return SCIP_OKAY;
2187
2188 /* determine symmetry specification */
2189 if ( SCIPgetNBinVars(scip) > 0 )
2190 type |= (int) SYM_SPEC_BINARY;
2191 if ( SCIPgetNIntVars(scip) > 0 )
2192 type |= (int) SYM_SPEC_INTEGER;
2193 /* count implicit integer variables as real variables, since we cannot currently handle integral variables well */
2194 if ( SCIPgetNContVars(scip) > 0 || SCIPgetNImplVars(scip) > 0 )
2195 type |= (int) SYM_SPEC_REAL;
2196
2197 /* skip symmetry computation if required variables are not present */
2198 if ( ! (type & symspecrequire) )
2199 {
2201 " (%.1fs) symmetry computation skipped: type (bin %c, int %c, cont %c) does not match requirements (bin %c, int %c, cont %c).\n",
2203 SCIPgetNBinVars(scip) > 0 ? '+' : '-',
2204 SCIPgetNIntVars(scip) > 0 ? '+' : '-',
2205 SCIPgetNContVars(scip) + SCIPgetNImplVars(scip) > 0 ? '+' : '-',
2206 (symspecrequire & (int) SYM_SPEC_BINARY) != 0 ? '+' : '-',
2207 (symspecrequire & (int) SYM_SPEC_INTEGER) != 0 ? '+' : '-',
2208 (symspecrequire & (int) SYM_SPEC_REAL) != 0 ? '+' : '-');
2209
2210 return SCIP_OKAY;
2211 }
2212
2213 /* skip computation if symmetry has already been computed */
2214 if ( propdata->computedsymmetry )
2215 return SCIP_OKAY;
2216
2217 assert( propdata->npermvars == 0 );
2218 assert( propdata->permvars == NULL );
2219 assert( propdata->nperms < 0 );
2220 assert( propdata->nmaxperms == 0 );
2221 assert( propdata->perms == NULL );
2222
2223 /* output message */
2225 " (%.1fs) symmetry computation started: requiring (bin %c, int %c, cont %c), (fixed: bin %c, int %c, cont %c)\n",
2227 (symspecrequire & (int) SYM_SPEC_BINARY) != 0 ? '+' : '-',
2228 (symspecrequire & (int) SYM_SPEC_INTEGER) != 0 ? '+' : '-',
2229 (symspecrequire & (int) SYM_SPEC_REAL) != 0 ? '+' : '-',
2230 (symspecrequirefixed & (int) SYM_SPEC_BINARY) != 0 ? '+' : '-',
2231 (symspecrequirefixed & (int) SYM_SPEC_INTEGER) != 0 ? '+' : '-',
2232 (symspecrequirefixed & (int) SYM_SPEC_REAL) != 0 ? '+' : '-');
2233
2234 /* output warning if we want to fix certain symmetry parts that we also want to compute */
2235 if ( symspecrequire & symspecrequirefixed )
2236 SCIPwarningMessage(scip, "Warning: some required symmetries must be fixed.\n");
2237
2238 /* determine maximal number of generators depending on the number of variables */
2239 maxgenerators = propdata->maxgenerators;
2240 maxgenerators = MIN(maxgenerators, MAXGENNUMERATOR / nvars);
2241
2242 /* actually compute (global) symmetry */
2243 SCIP_CALL( computeSymmetryGroup(scip, (SYM_SYMTYPE) propdata->symtype,
2244 propdata->compresssymmetries, propdata->compressthreshold,
2245 maxgenerators, symspecrequirefixed, propdata->checksymmetries, &propdata->permvars,
2246 &propdata->npermvars, &propdata->nbinpermvars, &propdata->permvardomaincenter,
2247 &propdata->perms, &propdata->nperms, &propdata->nmaxperms,
2248 &propdata->nmovedvars, &propdata->binvaraffected, &propdata->compressed,
2249 &propdata->log10groupsize, &symcodetime, &successful) );
2250
2251 /* mark that we have computed the symmetry group */
2252 propdata->computedsymmetry = TRUE;
2253
2254 /* store restart level */
2255 propdata->lastrestart = SCIPgetNRuns(scip);
2256
2257 /* return if not successful */
2258 if ( ! successful )
2259 {
2260 assert( checkSymmetryDataFree(propdata) );
2261 SCIPverbMessage(scip, SCIP_VERBLEVEL_HIGH, NULL, " (%.1fs) could not compute symmetry\n", SCIPgetSolvingTime(scip));
2262
2263 return SCIP_OKAY;
2264 }
2265
2266 /* return if no symmetries found */
2267 if ( propdata->nperms == 0 )
2268 {
2269 assert( checkSymmetryDataFree(propdata) );
2270 SCIPverbMessage(scip, SCIP_VERBLEVEL_HIGH, NULL, " (%.1fs) no symmetry present (symcode time: %.2f)\n", SCIPgetSolvingTime(scip), symcodetime);
2271
2272 return SCIP_OKAY;
2273 }
2274 assert( propdata->nperms > 0 );
2275 assert( propdata->npermvars > 0 );
2276
2277 /* display statistics */
2278 SCIPverbMessage(scip, SCIP_VERBLEVEL_HIGH, NULL, " (%.1fs) symmetry computation finished: %d generators found (max: ",
2279 SCIPgetSolvingTime(scip), propdata->nperms);
2280
2281 /* display statistics: maximum number of generators */
2282 if ( maxgenerators == 0 )
2284 else
2285 SCIPverbMessage(scip, SCIP_VERBLEVEL_HIGH, NULL, "%d", maxgenerators);
2286
2287 /* display statistics: log10 group size, number of affected vars*/
2288 SCIPverbMessage(scip, SCIP_VERBLEVEL_HIGH, NULL, ", log10 of symmetry group size: %.1f", propdata->log10groupsize);
2289
2290 if ( propdata->displaynorbitvars )
2291 {
2292 if ( propdata->nmovedvars == -1 )
2293 {
2294 SCIP_CALL( SCIPdetermineNVarsAffectedSym(scip, propdata->perms, propdata->nperms, propdata->permvars,
2295 propdata->npermvars, &(propdata->nmovedvars)) );
2296 }
2297 SCIPverbMessage(scip, SCIP_VERBLEVEL_HIGH, NULL, ", number of affected variables: %d)\n", propdata->nmovedvars);
2298 }
2299 SCIPverbMessage(scip, SCIP_VERBLEVEL_HIGH, NULL, ") (symcode time: %.2f)\n", symcodetime);
2300
2301 /* capture all variables while they are in the permvars array */
2302 for (i = 0; i < propdata->npermvars; ++i)
2303 {
2304 SCIP_CALL( SCIPcaptureVar(scip, propdata->permvars[i]) );
2305 }
2306
2307 return SCIP_OKAY;
2308}
2309
2310
2311/*
2312 * Functions for symmetry constraints
2313 */
2314
2315
2316/** Checks whether given set of 2-cycle permutations forms an orbitope and if so, builds the variable index matrix.
2317 *
2318 * If @p activevars == NULL, then the function assumes all permutations of the component are active and therefore all
2319 * moved vars are considered.
2320 *
2321 * We need to keep track of the number of generated columns, because we might not be able to detect all orbitopes.
2322 * For example (1,2), (2,3), (3,4), (3,5) defines the symmetric group on {1,2,3,4,5}, but the generators we expect
2323 * in our construction need shape (1,2), (2,3), (3,4), (4,5).
2324 *
2325 * @pre @p orbitopevaridx has to be an initialized 2D array of size @p ntwocycles x @p nperms
2326 * @pre @p columnorder has to be an initialized array of size nperms
2327 * @pre @p nusedelems has to be an initialized array of size npermvars
2328 */
2329static
2331 SCIP* scip, /**< SCIP instance */
2332 SCIP_VAR** permvars, /**< array of all permutation variables */
2333 int npermvars, /**< number of permutation variables */
2334 int** perms, /**< array of all permutations of the symmetry group */
2335 int* activeperms, /**< indices of the relevant permutations in perms */
2336 int ntwocycles, /**< number of 2-cycles in the permutations */
2337 int nactiveperms, /**< number of active permutations */
2338 int** orbitopevaridx, /**< pointer to store variable index matrix */
2339 int* columnorder, /**< pointer to store column order */
2340 int* nusedelems, /**< pointer to store how often each element was used */
2341 int* nusedcols, /**< pointer to store number of columns used in orbitope (or NULL) */
2342 SCIP_Shortbool* rowisbinary, /**< pointer to store which rows are binary (or NULL) */
2343 SCIP_Bool* isorbitope, /**< buffer to store result */
2344 SCIP_Shortbool* activevars /**< bitset to store whether a variable is active (or NULL) */
2345 )
2346{ /*lint --e{571}*/
2347 SCIP_Bool* usedperm;
2348 SCIP_Bool foundperm = FALSE;
2349 int nusedperms = 0;
2350 int nfilledcols;
2351 int coltoextend;
2352 int ntestedperms = 0;
2353 int row = 0;
2354 int j;
2355
2356 assert( scip != NULL );
2357 assert( permvars != NULL );
2358 assert( perms != NULL );
2359 assert( activeperms != NULL );
2360 assert( orbitopevaridx != NULL );
2361 assert( columnorder != NULL );
2362 assert( nusedelems != NULL );
2363 assert( isorbitope != NULL );
2364 assert( nactiveperms > 0 );
2365 assert( ntwocycles > 0 );
2366 assert( npermvars > 0 );
2367 assert( activevars == NULL || (0 <= nactiveperms && nactiveperms < npermvars) );
2368
2369 *isorbitope = TRUE;
2370 if ( nusedcols != NULL )
2371 *nusedcols = 0;
2372
2373 /* whether a permutation was considered to contribute to orbitope */
2374 SCIP_CALL( SCIPallocClearBufferArray(scip, &usedperm, nactiveperms) );
2375
2376 /* fill first two columns of orbitopevaridx matrix */
2377
2378 /* look for the first active permutation which moves an active variable */
2379 while ( ! foundperm )
2380 {
2381 int permidx;
2382
2383 assert( ntestedperms < nactiveperms );
2384
2385 permidx = activeperms[ntestedperms];
2386
2387 for (j = 0; j < npermvars; ++j)
2388 {
2389 if ( activevars != NULL && ! activevars[j] )
2390 continue;
2391
2392 assert( activevars == NULL || activevars[perms[permidx][j]] );
2393
2394 /* avoid adding the same 2-cycle twice */
2395 if ( perms[permidx][j] > j )
2396 {
2397 assert( SCIPvarIsBinary(permvars[j]) == SCIPvarIsBinary(permvars[perms[permidx][j]]) );
2398
2399 if ( rowisbinary != NULL && SCIPvarIsBinary(permvars[j]) )
2400 rowisbinary[row] = TRUE;
2401
2402 orbitopevaridx[row][0] = j;
2403 orbitopevaridx[row++][1] = perms[permidx][j];
2404 ++(nusedelems[j]);
2405 ++(nusedelems[perms[permidx][j]]);
2406
2407 foundperm = TRUE;
2408 }
2409
2410 if ( row == ntwocycles )
2411 break;
2412 }
2413
2414 ++ntestedperms;
2415 }
2416
2417 /* in the subgroup case it might happen that a generator has less than ntwocycles many 2-cyles */
2418 if ( row != ntwocycles )
2419 {
2420 *isorbitope = FALSE;
2421 SCIPfreeBufferArray(scip, &usedperm);
2422 return SCIP_OKAY;
2423 }
2424
2425 usedperm[ntestedperms - 1] = TRUE;
2426 ++nusedperms;
2427 columnorder[0] = 0;
2428 columnorder[1] = 1;
2429 nfilledcols = 2;
2430
2431 /* extend orbitopevaridx matrix to the left, i.e., iteratively find new permutations that
2432 * intersect the last added left column in each row in exactly one entry, starting with
2433 * column 0 */
2434 coltoextend = 0;
2435 for (j = ntestedperms; j < nactiveperms; ++j)
2436 { /* lint --e{850} */
2437 SCIP_Bool success = FALSE;
2438 SCIP_Bool infeasible = FALSE;
2439
2440 if ( nusedperms == nactiveperms )
2441 break;
2442
2443 if ( usedperm[j] )
2444 continue;
2445
2446 SCIP_CALL( SCIPextendSubOrbitope(orbitopevaridx, ntwocycles, nfilledcols, coltoextend,
2447 perms[activeperms[j]], TRUE, &nusedelems, permvars, NULL, &success, &infeasible) );
2448
2449 if ( infeasible )
2450 {
2451 *isorbitope = FALSE;
2452 break;
2453 }
2454 else if ( success )
2455 {
2456 usedperm[j] = TRUE;
2457 ++nusedperms;
2458 coltoextend = nfilledcols;
2459 columnorder[nfilledcols++] = -1; /* mark column to be filled from the left */
2460 j = 0; /*lint !e850*/ /* reset j since previous permutations can now intersect with the latest added column */
2461 }
2462 }
2463
2464 if ( ! *isorbitope ) /*lint !e850*/
2465 {
2466 SCIPfreeBufferArray(scip, &usedperm);
2467 return SCIP_OKAY;
2468 }
2469
2470 coltoextend = 1;
2471 for (j = ntestedperms; j < nactiveperms; ++j)
2472 { /*lint --e(850)*/
2473 SCIP_Bool success = FALSE;
2474 SCIP_Bool infeasible = FALSE;
2475
2476 if ( nusedperms == nactiveperms )
2477 break;
2478
2479 if ( usedperm[j] )
2480 continue;
2481
2482 SCIP_CALL( SCIPextendSubOrbitope(orbitopevaridx, ntwocycles, nfilledcols, coltoextend,
2483 perms[activeperms[j]], FALSE, &nusedelems, permvars, NULL, &success, &infeasible) );
2484
2485 if ( infeasible )
2486 {
2487 *isorbitope = FALSE;
2488 break;
2489 }
2490 else if ( success )
2491 {
2492 usedperm[j] = TRUE;
2493 ++nusedperms;
2494 coltoextend = nfilledcols;
2495 columnorder[nfilledcols] = 1; /* mark column to be filled from the right */
2496 ++nfilledcols;
2497 j = 0; /*lint !e850*/ /* reset j since previous permutations can now intersect with the latest added column */
2498 }
2499 }
2500
2501 if ( activevars == NULL && nusedperms < nactiveperms ) /*lint !e850*/
2502 *isorbitope = FALSE;
2503
2504 if ( nusedcols != NULL )
2505 *nusedcols = nfilledcols;
2506 assert( ! *isorbitope || activevars == NULL || nusedperms < nfilledcols );
2507
2508 SCIPfreeBufferArray(scip, &usedperm);
2509
2510 return SCIP_OKAY;
2511}
2512
2513/** choose an order in which the generators should be added for subgroup detection */
2514static
2516 SCIP* scip, /**< SCIP instance */
2517 SCIP_PROPDATA* propdata, /**< pointer to data of symmetry propagator */
2518 int compidx, /**< index of component */
2519 int** genorder, /**< (initialized) buffer to store the resulting order of generator */
2520 int* ntwocycleperms /**< pointer to store the number of 2-cycle permutations in component compidx */
2521 )
2522{
2523 int** perms;
2524 int* components;
2525 int* componentbegins;
2526 int* ntwocycles;
2527 int npermvars;
2528 int npermsincomp;
2529 int i;
2530
2531 assert( scip != NULL );
2532 assert( propdata != NULL );
2533 assert( compidx >= 0 );
2534 assert( compidx < propdata->ncomponents );
2535 assert( genorder != NULL );
2536 assert( *genorder != NULL );
2537 assert( ntwocycleperms != NULL );
2538 assert( propdata->computedsymmetry );
2539 assert( propdata->nperms > 0 );
2540 assert( propdata->perms != NULL );
2541 assert( propdata->npermvars > 0 );
2542 assert( propdata->ncomponents > 0 );
2543 assert( propdata->components != NULL );
2544 assert( propdata->componentbegins != NULL );
2545
2546 perms = propdata->perms;
2547 npermvars = propdata->npermvars;
2548 components = propdata->components;
2549 componentbegins = propdata->componentbegins;
2550 npermsincomp = componentbegins[compidx + 1] - componentbegins[compidx];
2551 *ntwocycleperms = npermsincomp;
2552
2553 SCIP_CALL( SCIPallocBufferArray(scip, &ntwocycles, npermsincomp) );
2554
2555 for (i = 0; i < npermsincomp; ++i)
2556 {
2557 int* perm;
2558 int nbincycles;
2559
2560 perm = perms[components[componentbegins[compidx] + i]];
2561
2562 SCIP_CALL( SCIPisInvolutionPerm(perm, propdata->permvars, npermvars, &(ntwocycles[i]), &nbincycles, FALSE) );
2563
2564 /* we skip permutations which do not purely consist of 2-cycles */
2565 if ( ntwocycles[i] == 0 )
2566 {
2567 /* we change the number of two cycles for this perm so that it will be sorted to the end */
2568 if ( propdata->preferlessrows )
2569 ntwocycles[i] = npermvars;
2570 else
2571 ntwocycles[i] = 0;
2572 --(*ntwocycleperms);
2573 }
2574 else if ( ! propdata->preferlessrows )
2575 ntwocycles[i] = - ntwocycles[i];
2576 }
2577
2578 SCIPsortIntInt(ntwocycles, *genorder, npermsincomp);
2579
2580 SCIPfreeBufferArray(scip, &ntwocycles);
2581
2582 return SCIP_OKAY;
2583}
2584
2585
2586/** builds the graph for symmetric subgroup detection from the given permutation of generators
2587 *
2588 * After execution, @p graphcomponents contains all permvars sorted by their color and component,
2589 * @p graphcompbegins points to the indices where new components in @p graphcomponents start and
2590 * @p compcolorbegins points to the indices where new colors in @p graphcompbegins start.
2591*/
2592static
2594 SCIP* scip, /**< SCIP instance */
2595 SCIP_PROPDATA* propdata, /**< pointer to data of symmetry propagator */
2596 int* genorder, /**< order in which the generators should be considered */
2597 int ntwocycleperms, /**< number of 2-cycle permutations in this component */
2598 int compidx, /**< index of the component */
2599 int** graphcomponents, /**< buffer to store the components of the graph (ordered var indices) */
2600 int** graphcompbegins, /**< buffer to store the indices of each new graph component */
2601 int** compcolorbegins, /**< buffer to store at which indices a new color begins */
2602 int* ngraphcomponents, /**< pointer to store the number of graph components */
2603 int* ncompcolors, /**< pointer to store the number of different colors */
2604 int** usedperms, /**< buffer to store the indices of permutations that were used */
2605 int* nusedperms, /**< pointer to store the number of used permutations in the graph */
2606 int usedpermssize, /**< initial size of usedperms */
2607 SCIP_Shortbool* permused /**< initialized buffer to store which permutations have been used
2608 * (identified by index in component) */
2609 )
2610{
2611 SCIP_DISJOINTSET* vartocomponent;
2612 SCIP_DISJOINTSET* comptocolor;
2613 int** perms;
2614 int* components;
2615 int* componentbegins;
2616 int* componentslastperm;
2617 SYM_SORTGRAPHCOMPVARS graphcompvartype;
2618 int npermvars;
2619 int nextcolor;
2620 int nextcomp;
2621 int j;
2622 int k;
2623
2624 assert( scip != NULL );
2625 assert( propdata != NULL );
2626 assert( graphcomponents != NULL );
2627 assert( graphcompbegins != NULL );
2628 assert( compcolorbegins != NULL );
2629 assert( ngraphcomponents != NULL );
2630 assert( ncompcolors != NULL );
2631 assert( genorder != NULL );
2632 assert( usedperms != NULL );
2633 assert( nusedperms != NULL );
2634 assert( usedpermssize > 0 );
2635 assert( permused != NULL );
2636 assert( ntwocycleperms >= 0 );
2637 assert( compidx >= 0 );
2638 assert( compidx < propdata->ncomponents );
2639 assert( propdata->computedsymmetry );
2640 assert( propdata->nperms > 0 );
2641 assert( propdata->perms != NULL );
2642 assert( propdata->npermvars > 0 );
2643 assert( propdata->ncomponents > 0 );
2644 assert( propdata->components != NULL );
2645 assert( propdata->componentbegins != NULL );
2646 assert( ! propdata->componentblocked[compidx] );
2647
2648 perms = propdata->perms;
2649 npermvars = propdata->npermvars;
2650 components = propdata->components;
2651 componentbegins = propdata->componentbegins;
2652 *nusedperms = 0;
2653
2654 assert( ntwocycleperms <= componentbegins[compidx + 1] - componentbegins[compidx] );
2655
2656 SCIP_CALL( SCIPcreateDisjointset(scip, &vartocomponent, npermvars) );
2657 SCIP_CALL( SCIPcreateDisjointset(scip, &comptocolor, npermvars) );
2658 SCIP_CALL( SCIPallocBufferArray( scip, &componentslastperm, npermvars) );
2659
2660 for (k = 0; k < npermvars; ++k)
2661 componentslastperm[k] = -1;
2662
2663 for (j = 0; j < ntwocycleperms; ++j)
2664 {
2665 int* perm;
2666 int firstcolor = -1;
2667
2668 /* use given order of generators */
2669 perm = perms[components[componentbegins[compidx] + genorder[j]]];
2670 assert( perm != NULL );
2671
2672 /* iteratively handle each swap of perm until an invalid one is found or all edges have been added */
2673 for (k = 0; k < npermvars; ++k)
2674 {
2675 int comp1;
2676 int comp2;
2677 int color1;
2678 int color2;
2679 int img;
2680
2681 img = perm[k];
2682 assert( perm[img] == k );
2683
2684 if ( img <= k )
2685 continue;
2686
2687 comp1 = SCIPdisjointsetFind(vartocomponent, k);
2688 comp2 = SCIPdisjointsetFind(vartocomponent, img);
2689
2690 if ( comp1 == comp2 )
2691 {
2692 /* another permutation has already merged these variables into one component; store its color */
2693 if ( firstcolor < 0 )
2694 {
2695 assert( SCIPdisjointsetFind(comptocolor, comp1) == SCIPdisjointsetFind(comptocolor, comp2) );
2696 firstcolor = SCIPdisjointsetFind(comptocolor, comp1);
2697 }
2698 componentslastperm[comp1] = j;
2699 continue;
2700 }
2701
2702 /* if it is the second time that the component is used for this generator,
2703 * it is not guaranteed that the group acts like the symmetric group, so skip it
2704 */
2705 if ( componentslastperm[comp1] == j || componentslastperm[comp2] == j )
2706 break;
2707
2708 color1 = SCIPdisjointsetFind(comptocolor, comp1);
2709 color2 = SCIPdisjointsetFind(comptocolor, comp2);
2710
2711 /* a generator is not allowed to connect two components of the same color, since they depend on each other */
2712 if ( color1 == color2 )
2713 break;
2714
2715 componentslastperm[comp1] = j;
2716 componentslastperm[comp2] = j;
2717
2718 if ( firstcolor < 0 )
2719 firstcolor = color1;
2720 }
2721
2722 /* if the generator is invalid, delete the newly added edges, go to next generator */
2723 if ( k < npermvars )
2724 continue;
2725
2726 /* if the generator only acts on already existing components, we don't have to store it */
2727 if ( firstcolor == -1 )
2728 continue;
2729
2730 /* check whether we need to resize */
2731 if ( *nusedperms >= usedpermssize )
2732 {
2733 int newsize = SCIPcalcMemGrowSize(scip, (*nusedperms) + 1);
2734 assert( newsize > usedpermssize );
2735
2736 SCIP_CALL( SCIPreallocBufferArray(scip, usedperms, newsize) );
2737
2738 usedpermssize = newsize;
2739 }
2740
2741 (*usedperms)[*nusedperms] = components[componentbegins[compidx] + genorder[j]];
2742 ++(*nusedperms);
2743 permused[genorder[j]] = TRUE;
2744
2745 /* if the generator can be added, update the datastructures for graph components and colors */
2746 for (k = 0; k < npermvars; ++k)
2747 {
2748 int comp1;
2749 int comp2;
2750 int color1;
2751 int color2;
2752 int img;
2753
2754 img = perm[k];
2755 assert( perm[img] == k );
2756
2757 if ( img <= k )
2758 continue;
2759
2760 comp1 = SCIPdisjointsetFind(vartocomponent, k);
2761 comp2 = SCIPdisjointsetFind(vartocomponent, img);
2762
2763 /* components and colors don't have to be updated if the components are the same */
2764 if ( comp1 == comp2 )
2765 continue;
2766
2767 color1 = SCIPdisjointsetFind(comptocolor, comp1);
2768 color2 = SCIPdisjointsetFind(comptocolor, comp2);
2769
2770 if ( color1 != color2 )
2771 {
2772 SCIPdisjointsetUnion(comptocolor, firstcolor, color1, TRUE);
2773 SCIPdisjointsetUnion(comptocolor, firstcolor, color2, TRUE);
2774 }
2775
2776 SCIPdisjointsetUnion(vartocomponent, comp1, comp2, FALSE);
2777
2778 assert( SCIPdisjointsetFind(vartocomponent, k) == SCIPdisjointsetFind(vartocomponent, img) );
2779 assert( SCIPdisjointsetFind(comptocolor, SCIPdisjointsetFind(vartocomponent, k)) == firstcolor );
2780 assert( SCIPdisjointsetFind(comptocolor, SCIPdisjointsetFind(vartocomponent, img)) == firstcolor );
2781 }
2782 }
2783
2784 SCIP_CALL( SCIPallocBlockMemoryArray(scip, graphcomponents, npermvars) );
2785 SCIP_CALL( SCIPallocBufferArray(scip, &(graphcompvartype.components), npermvars) );
2786 SCIP_CALL( SCIPallocBufferArray(scip, &(graphcompvartype.colors), npermvars) );
2787
2788 /*
2789 * At this point, we have built the colored graph. Now we transform the information in the
2790 * disjoint sets to the arrays graphcomponents, graphcompbegins, and compcolorbegins (see above).
2791 */
2792
2793 /* build the struct graphcompvartype which is used to sort the graphcomponents array */
2794 for (j = 0; j < npermvars; ++j)
2795 {
2796 int comp;
2797
2798 comp = SCIPdisjointsetFind(vartocomponent, j);
2799
2800 graphcompvartype.components[j] = comp;
2801 graphcompvartype.colors[j] = SCIPdisjointsetFind(comptocolor, comp);
2802
2803 (*graphcomponents)[j] = j;
2804 }
2805
2806 /* sort graphcomponents first by color, then by component */
2807 SCIPsort(*graphcomponents, SYMsortGraphCompVars, (void*) &graphcompvartype, npermvars);
2808
2809 *ngraphcomponents = SCIPdisjointsetGetComponentCount(vartocomponent);
2810 *ncompcolors = SCIPdisjointsetGetComponentCount(comptocolor);
2811 SCIP_CALL( SCIPallocBlockMemoryArray(scip, graphcompbegins, (*ngraphcomponents) + 1) );
2812 SCIP_CALL( SCIPallocBlockMemoryArray(scip, compcolorbegins, (*ncompcolors) + 1) );
2813
2814 nextcolor = 1;
2815 nextcomp = 1;
2816 (*graphcompbegins)[0] = 0;
2817 (*compcolorbegins)[0] = 0;
2818
2819 /* find the starting indices of new components and new colors */
2820 for (j = 1; j < npermvars; ++j)
2821 {
2822 int idx1;
2823 int idx2;
2824
2825 idx1 = (*graphcomponents)[j];
2826 idx2 = (*graphcomponents)[j-1];
2827
2828 assert( graphcompvartype.colors[idx1] >= graphcompvartype.colors[idx2] );
2829
2830 if ( graphcompvartype.components[idx1] != graphcompvartype.components[idx2] )
2831 {
2832 (*graphcompbegins)[nextcomp] = j;
2833
2834 if ( graphcompvartype.colors[idx1] > graphcompvartype.colors[idx2] )
2835 {
2836 (*compcolorbegins)[nextcolor] = nextcomp;
2837 ++nextcolor;
2838 }
2839
2840 ++nextcomp;
2841 }
2842 }
2843 assert( nextcomp == *ngraphcomponents );
2844 assert( nextcolor == *ncompcolors );
2845
2846 (*compcolorbegins)[nextcolor] = *ngraphcomponents;
2847 (*graphcompbegins)[nextcomp] = npermvars;
2848
2849 SCIPfreeBufferArray(scip, &(graphcompvartype.colors));
2850 SCIPfreeBufferArray(scip, &(graphcompvartype.components));
2851 SCIPfreeBufferArray(scip, &componentslastperm);
2852 SCIPfreeDisjointset(scip, &comptocolor);
2853 SCIPfreeDisjointset(scip, &vartocomponent);
2854
2855 return SCIP_OKAY;
2856}
2857
2858/** adds an orbitope constraint for a suitable color of the subgroup graph */
2859static
2861 SCIP* scip, /**< SCIP instance */
2862 SCIP_PROPDATA* propdata, /**< pointer to data of symmetry propagator */
2863 int* usedperms, /**< array of the permutations that build the orbitope */
2864 int nusedperms, /**< number of permutations in usedperms */
2865 int* compcolorbegins, /**< array indicating where a new graphcolor begins */
2866 int* graphcompbegins, /**< array indicating where a new graphcomponent begins */
2867 int* graphcomponents, /**< array of all variable indices sorted by color and comp */
2868 int graphcoloridx, /**< index of the graph color */
2869 int nrows, /**< number of rows in the orbitope */
2870 int ncols, /**< number of columns in the orbitope */
2871 int* firstvaridx, /**< buffer to store the index of the largest variable (or NULL) */
2872 int* compidxfirstrow, /**< buffer to store the comp index for the first row (or NULL) */
2873 int** lexorder, /**< pointer to array storing lexicographic order defined by sub orbitopes */
2874 int* nvarslexorder, /**< number of variables in lexicographic order */
2875 int* maxnvarslexorder, /**< maximum number of variables in lexicographic order */
2876 SCIP_Bool mayinteract, /**< whether orbitope's symmetries might interact with other symmetries */
2877 SCIP_Bool* success /**< whether the orbitope could be added */
2878 )
2879{ /*lint --e{571}*/
2880 char name[SCIP_MAXSTRLEN];
2881 SCIP_VAR*** orbitopevarmatrix;
2882 SCIP_Shortbool* activevars;
2883 int** orbitopevaridx;
2884 int* columnorder;
2885 int* nusedelems;
2886 SCIP_CONS* cons;
2887 SCIP_Bool isorbitope;
2888 SCIP_Bool infeasible = FALSE;
2889#ifndef NDEBUG
2890 int nactivevars = 0;
2891#endif
2892 int ngencols = 0;
2893 int k;
2894
2895 assert( scip != NULL );
2896 assert( propdata != NULL );
2897 assert( usedperms != NULL );
2898 assert( compcolorbegins != NULL );
2899 assert( graphcompbegins != NULL );
2900 assert( graphcomponents != NULL );
2901 assert( nusedperms > 0 );
2902 assert( nrows > 0 );
2903 assert( ncols > 0 );
2904 assert( lexorder != NULL );
2905 assert( nvarslexorder != NULL );
2906 assert( maxnvarslexorder != NULL );
2907
2908 *success = FALSE;
2909
2910 /* create hashset to mark variables */
2911 SCIP_CALL( SCIPallocClearBufferArray(scip, &activevars, propdata->npermvars) );
2912
2913 /* orbitope matrix for indices of variables in permvars array */
2914 SCIP_CALL( SCIPallocBufferArray(scip, &orbitopevaridx, nrows) );
2915 for (k = 0; k < nrows; ++k)
2916 {
2917 SCIP_CALL( SCIPallocBufferArray(scip, &orbitopevaridx[k], ncols) ); /*lint !e866*/
2918 }
2919
2920 /* order of columns of orbitopevaridx */
2921 SCIP_CALL( SCIPallocBufferArray(scip, &columnorder, ncols) );
2922 for (k = 0; k < ncols; ++k)
2923 columnorder[k] = ncols + 1;
2924
2925 /* count how often an element was used in the potential orbitope */
2926 SCIP_CALL( SCIPallocClearBufferArray(scip, &nusedelems, propdata->npermvars) );
2927
2928 /* mark variables in this subgroup orbitope */
2929 for (k = compcolorbegins[graphcoloridx]; k < compcolorbegins[graphcoloridx+1]; ++k)
2930 {
2931 SCIP_VAR* firstvar;
2932 int compstart;
2933 int l;
2934
2935 compstart = graphcompbegins[k];
2936 firstvar = propdata->permvars[graphcomponents[compstart]];
2937
2938 if ( ! SCIPvarIsBinary(firstvar) )
2939 continue;
2940
2941 for (l = 0; l < ncols; ++l)
2942 {
2943 int varidx;
2944
2945 varidx = graphcomponents[compstart + l];
2946 assert( ! activevars[varidx] );
2947
2948 activevars[varidx] = TRUE;
2949#ifndef NDEBUG
2950 ++nactivevars;
2951#endif
2952 }
2953 }
2954 assert( nactivevars == nrows * ncols );
2955
2956 /* build the variable index matrix for the orbitope
2957 *
2958 * It is possible that we find an orbitope, but not using all possible columns. For example
2959 * (1,2), (2,3), (3,4), (3,5) defines the symmetric group on {1,2,3,4,5}, but the generators
2960 * we expect in our construction need shape (1,2), (2,3), (3,4), (4,5). For this reason,
2961 * we need to store how many columns have been generated.
2962 *
2963 * @todo ensure compatibility with more general generators
2964 */
2965 SCIP_CALL( checkTwoCyclePermsAreOrbitope(scip, propdata->permvars, propdata->npermvars,
2966 propdata->perms, usedperms, nrows, nusedperms, orbitopevaridx, columnorder,
2967 nusedelems, &ngencols, NULL, &isorbitope, activevars) );
2968
2969 /* it might happen that we cannot detect the orbitope if it is generated by permutations with different
2970 * number of 2-cycles.
2971 */
2972 if ( ! isorbitope )
2973 {
2974 SCIPfreeBufferArray(scip, &nusedelems);
2975 SCIPfreeBufferArray(scip, &columnorder);
2976 for (k = nrows - 1; k >= 0; --k)
2977 {
2978 SCIPfreeBufferArray(scip, &orbitopevaridx[k]);
2979 }
2980 SCIPfreeBufferArray(scip, &orbitopevaridx);
2981 SCIPfreeBufferArray(scip, &activevars);
2982
2983 return SCIP_OKAY;
2984 }
2985
2986 /* There are three possibilities for the structure of columnorder:
2987 * 1) [0, 1, -1, -1, ..., -1]
2988 * 2) [0, 1, 1, 1, ..., 1]
2989 * 3) [0, 1, -1, -1, ...., -1, 1, 1, ..., 1]
2990 *
2991 * The '1'-columns will be added to the matrix first and in the last 2
2992 * cases the method starts from the right. So to store the variable index
2993 * that will be in the upper-left corner, we need either the entryin the
2994 * second column (case 1) or the entry in the last column (cases 2 and 3).
2995 */
2996 if ( firstvaridx != NULL )
2997 {
2998 if ( columnorder[ngencols-1] > -1 )
2999 *firstvaridx = orbitopevaridx[0][ngencols-1];
3000 else
3001 *firstvaridx = orbitopevaridx[0][1];
3002 }
3003
3004 /* find corresponding graphcomponent of first variable (needed for weak sbcs) */
3005 if ( compidxfirstrow != NULL && firstvaridx != NULL )
3006 {
3007 *compidxfirstrow = -1;
3008
3009 for (k = compcolorbegins[graphcoloridx]; k < compcolorbegins[graphcoloridx+1] && (*compidxfirstrow) < 0; ++k)
3010 {
3011 SCIP_VAR* firstvar;
3012 int compstart;
3013 int l;
3014
3015 compstart = graphcompbegins[k];
3016 firstvar = propdata->permvars[graphcomponents[compstart]];
3017
3018 if ( ! SCIPvarIsBinary(firstvar) )
3019 continue;
3020
3021 /* iterate over all columns (elements in orbit), because we cannot see from ngencols which columns
3022 * have been left out
3023 */
3024 for (l = 0; l < ncols; ++l)
3025 {
3026 if ( graphcomponents[compstart + l] == *firstvaridx )
3027 {
3028 *compidxfirstrow = k;
3029 break;
3030 }
3031 }
3032 }
3033 assert( *compidxfirstrow > -1 );
3034 }
3035
3036 /* prepare orbitope variable matrix */
3037 SCIP_CALL( SCIPallocBufferArray(scip, &orbitopevarmatrix, nrows) );
3038 for (k = 0; k < nrows; ++k)
3039 {
3040 SCIP_CALL( SCIPallocBufferArray(scip, &orbitopevarmatrix[k], ngencols) );
3041 }
3042
3043 /* build the matrix containing the actual variables of the orbitope */
3044 SCIP_CALL( SCIPgenerateOrbitopeVarsMatrix(scip, &orbitopevarmatrix, nrows, ngencols,
3045 propdata->permvars, propdata->npermvars, orbitopevaridx, columnorder,
3046 nusedelems, NULL, &infeasible, TRUE, lexorder, nvarslexorder, maxnvarslexorder) );
3047
3048 assert( ! infeasible );
3049 assert( firstvaridx == NULL || propdata->permvars[*firstvaridx] == orbitopevarmatrix[0][0] );
3050
3051 (void) SCIPsnprintf(name, SCIP_MAXSTRLEN, "suborbitope_%d_%d", graphcoloridx, propdata->norbitopes);
3052
3053 SCIP_CALL( SCIPcreateConsOrbitope(scip, &cons, name, orbitopevarmatrix,
3054 SCIP_ORBITOPETYPE_FULL, nrows, ngencols, FALSE, mayinteract, FALSE, FALSE, propdata->conssaddlp,
3056
3057 SCIP_CALL( SCIPaddCons(scip, cons) );
3058 *success = TRUE;
3059
3060 /* do not release constraint here - will be done later */
3062 &propdata->genorbconsssize, propdata->ngenorbconss + 1) );
3063 propdata->genorbconss[propdata->ngenorbconss++] = cons;
3064 ++propdata->norbitopes;
3065
3066 for (k = nrows - 1; k >= 0; --k)
3067 SCIPfreeBufferArray(scip, &orbitopevarmatrix[k]);
3068 SCIPfreeBufferArray(scip, &orbitopevarmatrix);
3069 SCIPfreeBufferArray(scip, &nusedelems);
3070 SCIPfreeBufferArray(scip, &columnorder);
3071 for (k = nrows - 1; k >= 0; --k)
3072 SCIPfreeBufferArray(scip, &orbitopevaridx[k]);
3073 SCIPfreeBufferArray(scip, &orbitopevaridx);
3074 SCIPfreeBufferArray(scip, &activevars);
3075
3076 return SCIP_OKAY;
3077}
3078
3079/** adds strong SBCs for a suitable color of the subgroup graph */
3080static
3082 SCIP* scip, /**< SCIP instance */
3083 SCIP_PROPDATA* propdata, /**< pointer to data of symmetry propagator */
3084 int* graphcompbegins, /**< array indicating where a new graphcomponent begins */
3085 int* graphcomponents, /**< array of all variable indices sorted by color and comp */
3086 int graphcompidx, /**< index of the graph component */
3087 SCIP_Bool storelexorder, /**< whether the lexicographic order induced by the orbitope shall be stored */
3088 int** lexorder, /**< pointer to array storing lexicographic order defined by sub orbitopes */
3089 int* nvarsorder, /**< number of variables in lexicographic order */
3090 int* maxnvarsorder /**< maximum number of variables in lexicographic order */
3091 )
3092{
3093 int k;
3094
3095 assert( scip != NULL );
3096 assert( propdata != NULL );
3097 assert( graphcompbegins != NULL );
3098 assert( graphcomponents != NULL );
3099 assert( graphcompidx >= 0 );
3100 assert( ! storelexorder || lexorder != NULL );
3101 assert( ! storelexorder || nvarsorder != NULL );
3102 assert( ! storelexorder || maxnvarsorder != NULL );
3103
3104 /* possibly store lexicographic order defined by strong SBCs */
3105 if ( storelexorder )
3106 {
3107 if ( *maxnvarsorder == 0 )
3108 {
3109 *maxnvarsorder = graphcompbegins[graphcompidx + 1] - graphcompbegins[graphcompidx + 1];
3110 *nvarsorder = 0;
3111
3112 SCIP_CALL( SCIPallocBlockMemoryArray(scip, lexorder, *maxnvarsorder) );
3113 }
3114 else
3115 {
3116 assert( *nvarsorder == *maxnvarsorder );
3117
3118 *maxnvarsorder += graphcompbegins[graphcompidx + 1] - graphcompbegins[graphcompidx + 1];
3119
3120 SCIP_CALL( SCIPreallocBlockMemoryArray(scip, lexorder, *nvarsorder, *maxnvarsorder) );
3121 }
3122
3123 (*lexorder)[*nvarsorder++] = graphcomponents[graphcompbegins[graphcompidx]];
3124 }
3125
3126 /* add strong SBCs (lex-max order) for chosen graph component */
3127 for (k = graphcompbegins[graphcompidx]+1; k < graphcompbegins[graphcompidx+1]; ++k)
3128 {
3129 char name[SCIP_MAXSTRLEN];
3130 SCIP_CONS* cons;
3131 SCIP_VAR* vars[2];
3132 SCIP_Real vals[2] = {1, -1};
3133
3134 vars[0] = propdata->permvars[graphcomponents[k-1]];
3135 vars[1] = propdata->permvars[graphcomponents[k]];
3136
3137 if ( storelexorder )
3138 (*lexorder)[*nvarsorder++] = graphcomponents[k];
3139
3140 (void) SCIPsnprintf(name, SCIP_MAXSTRLEN, "strong_sbcs_%s_%s", SCIPvarGetName(vars[0]), SCIPvarGetName(vars[1]));
3141
3142 SCIP_CALL( SCIPcreateConsLinear(scip, &cons, name, 2, vars, vals, 0.0,
3143 SCIPinfinity(scip), propdata->conssaddlp, propdata->conssaddlp, TRUE,
3144 TRUE, TRUE, FALSE, FALSE, FALSE, FALSE, FALSE) );
3145
3146 SCIP_CALL( SCIPaddCons(scip, cons) );
3147
3148#ifdef SCIP_MORE_DEBUG
3149 SCIP_CALL( SCIPprintCons(scip, cons, NULL) );
3150 SCIPinfoMessage(scip, NULL, "\n");
3151#endif
3152
3153 /* check whether we need to resize */
3155 &propdata->genlinconsssize, propdata->ngenlinconss + 1) );
3156 propdata->genlinconss[propdata->ngenlinconss] = cons;
3157 ++propdata->ngenlinconss;
3158 }
3159
3160 return SCIP_OKAY;
3161}
3162
3163/** adds weak SBCs for a suitable color of the subgroup graph */
3164static
3166 SCIP* scip, /**< SCIP instance */
3167 SCIP_PROPDATA* propdata, /**< pointer to data of symmetry propagator */
3168 int* compcolorbegins, /**< array indicating where a new graphcolor begins */
3169 int* graphcompbegins, /**< array indicating where a new graphcomponent begins */
3170 int* graphcomponents, /**< array of all variable indices sorted by color and comp */
3171 int ncompcolors, /**< number of colors in the graph */
3172 int* chosencomppercolor, /**< array indicating which comp was handled per color */
3173 int* firstvaridxpercolor,/**< array indicating the largest variable per color */
3174 int symgrpcompidx, /**< index of the component of the symmetry group */
3175 int* naddedconss, /**< buffer to store the number of added constraints */
3176 SCIP_Bool storelexorder, /**< whether the lexicographic order induced by the orbitope shall be stored */
3177 int** lexorder, /**< pointer to array storing lexicographic order defined by sub orbitopes */
3178 int* nvarsorder, /**< number of variables in lexicographic order */
3179 int* maxnvarsorder /**< maximum number of variables in lexicographic order */
3180 )
3181{ /*lint --e{571}*/
3182 SCIP_HASHMAP* varsinlexorder;
3183 SCIP_Shortbool* usedvars;
3184 SCIP_VAR* vars[2];
3185 SCIP_Real vals[2] = {1, -1};
3186 SCIP_Shortbool* varfound;
3187 int* orbit[2];
3188 int orbitsize[2] = {1, 1};
3189 int activeorb = 0;
3190 int chosencolor = -1;
3191 int j;
3192 int k;
3193
3194 assert( scip != NULL );
3195 assert( propdata != NULL );
3196 assert( compcolorbegins != NULL );
3197 assert( graphcompbegins != NULL );
3198 assert( graphcomponents != NULL );
3199 assert( firstvaridxpercolor != NULL );
3200 assert( chosencomppercolor != NULL );
3201 assert( naddedconss != NULL );
3202 assert( symgrpcompidx >= 0 );
3203 assert( symgrpcompidx < propdata->ncomponents );
3204 assert( ! storelexorder || lexorder != NULL );
3205 assert( ! storelexorder || nvarsorder != NULL );
3206 assert( ! storelexorder || maxnvarsorder != NULL );
3207
3208 *naddedconss = 0;
3209
3210 SCIP_CALL( SCIPallocCleanBufferArray(scip, &usedvars, propdata->npermvars) );
3211 SCIP_CALL( SCIPallocClearBufferArray(scip, &varfound, propdata->npermvars) );
3212 SCIP_CALL( SCIPallocBufferArray(scip, &orbit[0], propdata->npermvars) );
3213 SCIP_CALL( SCIPallocBufferArray(scip, &orbit[1], propdata->npermvars) );
3214
3215 /* Store the entries in lexorder in a hashmap, for fast lookups. */
3216 if ( lexorder == NULL || *lexorder == NULL )
3217 {
3218 /* Lexorder does not exist, so do not create hashmap. */
3219 varsinlexorder = NULL;
3220 }
3221 else
3222 {
3223 assert( *maxnvarsorder >= 0 );
3224 assert( *nvarsorder >= 0 );
3225
3226 SCIP_CALL( SCIPhashmapCreate(&varsinlexorder, SCIPblkmem(scip), *maxnvarsorder) );
3227
3228 for (k = 0; k < *nvarsorder; ++k)
3229 {
3230 /* add element from lexorder to hashmap.
3231 * Use insert, as duplicate entries in lexorder is not permitted. */
3232 assert((*lexorder)[k] >= 0);
3233 assert( ! SCIPhashmapExists(varsinlexorder, (void*) (size_t) (*lexorder)[k]) ); /* Use int as pointer */
3234 SCIP_CALL( SCIPhashmapInsertInt(varsinlexorder, (void*) (size_t) (*lexorder)[k], k) );
3235 }
3236 }
3237
3238 /* We will store the newest and the largest orbit and activeorb will be used to mark at which entry of the array
3239 * orbit the newly computed one will be stored. */
3240 if ( ncompcolors > 0 )
3241 {
3243 }
3244 for (j = 0; j < ncompcolors; ++j)
3245 {
3246 int graphcomp;
3247 int graphcompsize;
3248 int varidx;
3249
3250 /* skip color for which we did not add anything */
3251 if ( chosencomppercolor[j] < 0 )
3252 continue;
3253
3254 assert( firstvaridxpercolor[j] >= 0 );
3255
3256 graphcomp = chosencomppercolor[j];
3257 graphcompsize = graphcompbegins[graphcomp+1] - graphcompbegins[graphcomp];
3258 varidx = firstvaridxpercolor[j];
3259 assert(varidx >= 0);
3260
3261 /* if the first variable was already contained in another orbit or if there are no variables left anyway, skip the
3262 * component */
3263 if ( varfound[varidx] || graphcompsize == propdata->npermvars )
3264 continue;
3265
3266 /* If varidx is in lexorder, then it must be the first entry of lexorder. */
3267 if ( varsinlexorder != NULL
3268 && SCIPhashmapExists(varsinlexorder, (void*) (size_t) varidx)
3269 && lexorder != NULL && *lexorder != NULL && *maxnvarsorder > 0 && *nvarsorder > 0
3270 && (*lexorder)[0] != varidx )
3271 continue;
3272
3273 /* mark all variables that have been used in strong SBCs */
3274 for (k = graphcompbegins[graphcomp]; k < graphcompbegins[graphcomp+1]; ++k)
3275 {
3276 assert( 0 <= graphcomponents[k] && graphcomponents[k] < propdata->npermvars );
3277
3278 usedvars[graphcomponents[k]] = TRUE;
3279 }
3280
3281 SCIP_CALL( SCIPcomputeOrbitVar(scip, propdata->npermvars, propdata->perms,
3282 propdata->permstrans, propdata->components, propdata->componentbegins,
3283 usedvars, varfound, varidx, symgrpcompidx,
3284 orbit[activeorb], &orbitsize[activeorb]) );
3285
3286 assert( orbit[activeorb][0] == varidx );
3287
3288 if ( orbitsize[activeorb] > orbitsize[1 - activeorb] ) /*lint !e514*/
3289 {
3290 /* if the new orbit is larger then the old largest one, flip activeorb */
3291 activeorb = 1 - activeorb;
3292 chosencolor = j;
3293 }
3294
3295 /* reset array */
3296 for (k = graphcompbegins[graphcomp]; k < graphcompbegins[graphcomp+1]; ++k)
3297 usedvars[graphcomponents[k]] = FALSE;
3298 }
3299
3300 /* check if we have found at least one non-empty orbit */
3301 if ( chosencolor > -1 )
3302 {
3303 /* flip activeorb again to avoid confusion, it is then at the largest orbit */
3304 activeorb = 1 - activeorb;
3305
3306 assert( orbit[activeorb][0] == firstvaridxpercolor[chosencolor] );
3307 vars[0] = propdata->permvars[orbit[activeorb][0]];
3308
3309 assert( chosencolor > -1 );
3310 SCIPdebugMsg(scip, " adding %d weak sbcs for enclosing orbit of color %d.\n", orbitsize[activeorb]-1, chosencolor);
3311
3312 *naddedconss = orbitsize[activeorb] - 1;
3313
3314 /* add weak SBCs for rest of enclosing orbit */
3315 for (j = 1; j < orbitsize[activeorb]; ++j)
3316 {
3317 SCIP_CONS* cons;
3318 char name[SCIP_MAXSTRLEN];
3319
3320 vars[1] = propdata->permvars[orbit[activeorb][j]];
3321
3322 (void) SCIPsnprintf(name, SCIP_MAXSTRLEN, "weak_sbcs_%d_%s_%s", symgrpcompidx, SCIPvarGetName(vars[0]), SCIPvarGetName(vars[1]));
3323
3324 SCIP_CALL( SCIPcreateConsLinear(scip, &cons, name, 2, vars, vals, 0.0,
3325 SCIPinfinity(scip), propdata->conssaddlp, propdata->conssaddlp, TRUE,
3326 TRUE, TRUE, FALSE, FALSE, FALSE, FALSE, FALSE) );
3327
3328 SCIP_CALL( SCIPaddCons(scip, cons) );
3329
3330#ifdef SCIP_MORE_DEBUG
3331 SCIP_CALL( SCIPprintCons(scip, cons, NULL) );
3332 SCIPinfoMessage(scip, NULL, "\n");
3333#endif
3334
3335 /* check whether we need to resize */
3337 &propdata->genlinconsssize, propdata->ngenlinconss + 1) );
3338 propdata->genlinconss[propdata->ngenlinconss] = cons;
3339 ++propdata->ngenlinconss;
3340 }
3341
3342 /* possibly store lexicographic order defined by weak SBCs */
3343 if ( storelexorder )
3344 {
3345 int varidx;
3346
3347 varidx = orbit[activeorb][0];
3348 assert(varidx >= 0);
3349
3350 if ( *maxnvarsorder == 0 )
3351 {
3352 *maxnvarsorder = 1;
3353 *nvarsorder = 0;
3354
3355 SCIP_CALL( SCIPallocBlockMemoryArray(scip, lexorder, *maxnvarsorder) );
3356 (*lexorder)[(*nvarsorder)++] = varidx;
3357 }
3358 else
3359 {
3360 assert( *nvarsorder == *maxnvarsorder );
3361 assert( varsinlexorder != NULL );
3362 assert( lexorder != NULL );
3363 assert( *lexorder != NULL );
3364
3365 /* the leader of the weak inequalities has to be the first element in the lexicographic order */
3366 if ( varidx == (*lexorder)[0] )
3367 {
3368 /* lexorder is already ok!! */
3369 assert( SCIPhashmapExists(varsinlexorder, (void*) (size_t) varidx) );
3370 }
3371 else
3372 {
3373 /* Then varidx must not be in the lexorder,
3374 * We must add it at the front of the array, and maintain the current order. */
3375 assert( ! SCIPhashmapExists(varsinlexorder, (void*) (size_t) varidx) );
3376
3377 ++(*maxnvarsorder);
3378 ++(*nvarsorder);
3379
3380 SCIP_CALL( SCIPreallocBlockMemoryArray(scip, lexorder, *nvarsorder, *maxnvarsorder) );
3381
3382 /* Shift array by one position to the right */
3383 for (k = *maxnvarsorder - 1; k >= 1; --k)
3384 (*lexorder)[k] = (*lexorder)[k - 1];
3385
3386 (*lexorder)[0] = varidx;
3387 }
3388 }
3389 }
3390 }
3391 else
3392 SCIPdebugMsg(scip, " no further weak sbcs are valid\n");
3393
3394 SCIPfreeBufferArray(scip, &orbit[1]);
3395 SCIPfreeBufferArray(scip, &orbit[0]);
3396 if ( varsinlexorder != NULL )
3397 SCIPhashmapFree(&varsinlexorder);
3398 SCIPfreeBufferArray(scip, &varfound);
3399 SCIPfreeCleanBufferArray(scip, &usedvars);
3400
3401 return SCIP_OKAY;
3402}
3403
3404
3405/** temporarily adapt symmetry data to new variable order given by Schreier Sims */
3406static
3408 SCIP* scip, /**< SCIP instance */
3409 int** origperms, /**< permutation matrix w.r.t. original variable ordering */
3410 int** modifiedperms, /**< memory for permutation matrix w.r.t. new variable ordering */
3411 int nperms, /**< number of permutations */
3412 SCIP_VAR** origpermvars, /**< array of permutation vars w.r.t. original variable ordering */
3413 SCIP_VAR** modifiedpermvars, /**< memory for array of permutation vars w.r.t. new variable ordering */
3414 int npermvars, /**< length or modifiedpermvars array */
3415 int* leaders, /**< leaders of Schreier Sims constraints */
3416 int nleaders /**< number of leaders */
3417 )
3418{
3419 int* permvaridx;
3420 int* posinpermvar;
3421 int leader;
3422 int curposleader;
3423 int varidx;
3424 int lidx;
3425 int i;
3426 int l;
3427 int p;
3428
3429 assert( scip != NULL );
3430 assert( origperms != NULL );
3431 assert( modifiedperms != NULL );
3432 assert( nperms > 0 );
3433 assert( origpermvars != NULL );
3434 assert( modifiedpermvars != NULL );
3435 assert( npermvars > 0 );
3436 assert( leaders != NULL );
3437 assert( nleaders > 0 );
3438
3439 /* initialize map from position in lexicographic order to index of original permvar */
3440 SCIP_CALL( SCIPallocBufferArray(scip, &permvaridx, npermvars) );
3441 for (i = 0; i < npermvars; ++i)
3442 permvaridx[i] = i;
3443
3444 /* initialize map from permvaridx to its current position in the reordered permvars array */
3445 SCIP_CALL( SCIPallocBufferArray(scip, &posinpermvar, npermvars) );
3446 for (i = 0; i < npermvars; ++i)
3447 posinpermvar[i] = i;
3448
3449 /* Iterate over leaders and put the l-th leader to the l-th position of the lexicographic order.
3450 * We do this by swapping the l-th leader with the element at position l of the current permvars array. */
3451 for (l = 0; l < nleaders; ++l)
3452 {
3453 leader = leaders[l];
3454 curposleader = posinpermvar[leader];
3455 varidx = permvaridx[curposleader];
3456 lidx = permvaridx[l];
3457
3458 /* swap the permvar at position l with the l-th leader */
3459 permvaridx[curposleader] = lidx;
3460 permvaridx[l] = varidx;
3461
3462 /* update the position map */
3463 posinpermvar[lidx] = curposleader;
3464 posinpermvar[leader] = l;
3465 }
3466
3467 /* update the permvars array to new variable order */
3468 for (i = 0; i < npermvars; ++i)
3469 modifiedpermvars[i] = origpermvars[permvaridx[i]];
3470
3471 /* update the permutation to the new variable order */
3472 for (p = 0; p < nperms; ++p)
3473 {
3474 for (i = 0; i < npermvars; ++i)
3475 modifiedperms[p][i] = posinpermvar[origperms[p][permvaridx[i]]];
3476 }
3477
3478 SCIPfreeBufferArray(scip, &permvaridx);
3479 SCIPfreeBufferArray(scip, &posinpermvar);
3480
3481 return SCIP_OKAY;
3482}
3483
3484
3485/* returns the number of found orbitopes with at least three columns per graph component or 0
3486 * if the found orbitopes do not satisfy certain criteria for being used
3487 */
3488static
3490 SCIP_VAR** permvars, /**< array of variables affected by symmetry */
3491 int* graphcomponents, /**< array of graph components */
3492 int* graphcompbegins, /**< array indicating starting position of graph components */
3493 int* compcolorbegins, /**< array indicating starting positions of potential orbitopes */
3494 int ncompcolors, /**< number of components encoded in compcolorbegins */
3495 int symcompsize /**< size of symmetry component for that we detect suborbitopes */
3496 )
3497{
3498 SCIP_Bool oneorbitopecriterion = FALSE;
3499 SCIP_Bool multorbitopecriterion = FALSE;
3500 int norbitopes = 0;
3501 int j;
3502
3503 assert( graphcompbegins != NULL );
3504 assert( compcolorbegins != NULL );
3505 assert( ncompcolors >= 0 );
3506 assert( symcompsize > 0 );
3507
3508 for (j = 0; j < ncompcolors; ++j)
3509 {
3510 SCIP_VAR* firstvar;
3511 int largestcompsize = 0;
3512 int nbinrows= 0;
3513 int k;
3514
3515 /* skip trivial components */
3516 if ( graphcompbegins[compcolorbegins[j+1]] - graphcompbegins[compcolorbegins[j]] < 2 )
3517 continue;
3518
3519 /* check whether components of this color build an orbitope (with > 2 columns) */
3520 for (k = compcolorbegins[j]; k < compcolorbegins[j+1]; ++k)
3521 {
3522 int compsize;
3523
3524 compsize = graphcompbegins[k+1] - graphcompbegins[k];
3525
3526 /* the first component that we are looking at for this color */
3527 if ( largestcompsize < 1 )
3528 {
3529 if ( compsize < 3 )
3530 break;
3531
3532 largestcompsize = compsize;
3533 }
3534 else if ( compsize != largestcompsize )
3535 break;
3536
3537 firstvar = permvars[graphcomponents[graphcompbegins[k]]];
3538
3539 /* count number of binary orbits (comps) */
3540 if ( SCIPvarIsBinary(firstvar) )
3541 ++nbinrows;
3542 }
3543
3544 /* we have found an orbitope */
3545 if ( k == compcolorbegins[j+1] )
3546 {
3547 SCIP_Real threshold;
3548 int ncols;
3549
3550 ++norbitopes;
3551 ncols = graphcompbegins[compcolorbegins[j] + 1] - graphcompbegins[compcolorbegins[j]];
3552
3553 threshold = 0.7 * (SCIP_Real) symcompsize;
3554
3555 /* check whether criteria for adding orbitopes are satisfied */
3556 if ( nbinrows <= 2 * ncols || (nbinrows <= 8 * ncols && nbinrows < 100) )
3557 multorbitopecriterion = TRUE;
3558 else if ( nbinrows <= 3 * ncols || (SCIP_Real) nbinrows * ncols >= threshold )
3559 oneorbitopecriterion = TRUE;
3560 }
3561 }
3562
3563 if ( (norbitopes == 1 && oneorbitopecriterion) || (norbitopes >= 2 && multorbitopecriterion) )
3564 return norbitopes;
3565
3566 return 0;
3567}
3568
3569
3570/** checks whether subgroups of the components are symmetric groups and adds SBCs for them */
3571static
3573 SCIP* scip, /**< SCIP instance */
3574 SCIP_PROPDATA* propdata, /**< pointer to data of symmetry propagator */
3575 int cidx /**< index of component which shall be handled */
3576 )
3577{
3578 int* genorder;
3579 int p;
3580#ifdef SCIP_DEBUG
3581 int norbitopes = 0;
3582 int nstrongsbcs = 0;
3583 int nweaksbcs = 0;
3584#endif
3585 int** modifiedperms;
3586 SCIP_VAR** modifiedpermvars;
3587 int* nvarsincomponent;
3588
3589 int* graphcomponents;
3590 int* graphcompbegins;
3591 int* compcolorbegins;
3592 int* chosencomppercolor = NULL;
3593 int* firstvaridxpercolor = NULL;
3594 int* usedperms;
3595 int usedpermssize;
3596 int ngraphcomponents;
3597 int ncompcolors;
3598 int ntwocycleperms;
3599 int npermsincomp;
3600 int nusedperms;
3601 int ntrivialcolors = 0;
3602 int j;
3603 int* lexorder = NULL;
3604 int nvarslexorder = 0;
3605 int maxnvarslexorder = 0;
3606 SCIP_Shortbool* permused;
3607 SCIP_Bool allpermsused = FALSE;
3608 SCIP_Bool handlednonbinarysymmetry = FALSE;
3609 int norbitopesincomp;
3610
3611 assert( scip != NULL );
3612 assert( propdata != NULL );
3613 assert( propdata->computedsymmetry );
3614 assert( propdata->nperms >= 0 );
3615 assert( 0 <= cidx && cidx < propdata->ncomponents );
3616 assert( propdata->components != NULL );
3617 assert( propdata->componentbegins != NULL );
3618
3619 /* exit if no symmetry is present or component is blocked */
3620 if ( propdata->nperms == 0 || propdata->componentblocked[cidx] )
3621 return SCIP_OKAY;
3622
3623 /* exit if instance is too large */
3624 if ( SCIPgetNConss(scip) > propdata->maxnconsssubgroup )
3625 return SCIP_OKAY;
3626
3627 assert( propdata->nperms > 0 );
3628 assert( propdata->perms != NULL );
3629 assert( propdata->npermvars > 0 );
3630 assert( propdata->permvars != NULL );
3631
3632 /* create array for permutation order */
3633 SCIP_CALL( SCIPallocBufferArray(scip, &genorder, propdata->nperms) );
3634
3635 /* create arrays for modified permutations in case we adapt the lexicographic order because of suborbitopes */
3636 SCIP_CALL( SCIPallocBufferArray(scip, &modifiedperms, propdata->nperms) );
3637 for (p = 0; p < propdata->nperms; ++p)
3638 {
3639 SCIP_CALL( SCIPallocBufferArray(scip, &modifiedperms[p], propdata->npermvars) );
3640 }
3641 SCIP_CALL( SCIPallocBufferArray(scip, &modifiedpermvars, propdata->npermvars) );
3642
3643 SCIP_CALL( SCIPallocClearBufferArray(scip, &nvarsincomponent, propdata->npermvars) );
3644 for (p = 0; p < propdata->npermvars; ++p)
3645 {
3646 if ( propdata->vartocomponent[p] >= 0 )
3647 ++nvarsincomponent[propdata->vartocomponent[p]];
3648 }
3649
3650 SCIPdebugMsg(scip, "starting subgroup detection routine for component %d\n", cidx);
3651
3652 npermsincomp = propdata->componentbegins[cidx + 1] - propdata->componentbegins[cidx];
3653
3654 /* set the first npermsincomp entries of genorder; the others are not used for this component */
3655 for (j = 0; j < npermsincomp; ++j)
3656 genorder[j] = j;
3657
3658 SCIP_CALL( chooseOrderOfGenerators(scip, propdata, cidx, &genorder, &ntwocycleperms) );
3659
3660 assert( ntwocycleperms >= 0 );
3661 assert( ntwocycleperms <= npermsincomp );
3662
3663 SCIPdebugMsg(scip, "component %d has %d permutations consisting of 2-cycles\n", cidx, ntwocycleperms);
3664
3665#ifdef SCIP_MORE_DEBUG
3666 SCIP_Bool* used;
3667 int perm;
3668 int p;
3669 int k;
3670
3671 SCIP_CALL( SCIPallocBufferArray(scip, &used, propdata->npermvars) );
3672 for (p = propdata->componentbegins[cidx]; p < propdata->componentbegins[cidx+1]; ++p)
3673 {
3674 perm = propdata->components[p];
3675
3676 for (k = 0; k < propdata->npermvars; ++k)
3677 used[k] = FALSE;
3678
3679 SCIPverbMessage(scip, SCIP_VERBLEVEL_HIGH, NULL, "permutation %d\n", perm);
3680
3681 for (k = 0; k < propdata->npermvars; ++k)
3682 {
3683 if ( used[k] )
3684 continue;
3685
3686 j = propdata->perms[perm][k];
3687
3688 if ( k == j )
3689 continue;
3690
3691 SCIPverbMessage(scip, SCIP_VERBLEVEL_HIGH, NULL, "(%s,", SCIPvarGetName(propdata->permvars[k]));
3692 used[k] = TRUE;
3693 while (j != k)
3694 {
3695 SCIPverbMessage(scip, SCIP_VERBLEVEL_HIGH, NULL, "%s,", SCIPvarGetName(propdata->permvars[j]));
3696 used[j] = TRUE;
3697
3698 j = propdata->perms[perm][j];
3699 }
3701 }
3703 }
3704
3705 SCIPfreeBufferArray(scip, &used);
3706#endif
3707
3708 if ( ntwocycleperms < 2 )
3709 {
3710 SCIPdebugMsg(scip, " --> skip\n");
3711 goto FREEBASICMEM;
3712 }
3713
3714 usedpermssize = ntwocycleperms / 2;
3715 SCIP_CALL( SCIPallocBufferArray(scip, &usedperms, usedpermssize) );
3716 SCIP_CALL( SCIPallocClearBufferArray(scip, &permused, npermsincomp) );
3717
3718 SCIP_CALL( buildSubgroupGraph(scip, propdata, genorder, ntwocycleperms, cidx,
3719 &graphcomponents, &graphcompbegins, &compcolorbegins, &ngraphcomponents,
3720 &ncompcolors, &usedperms, &nusedperms, usedpermssize, permused) );
3721
3722 SCIPdebugMsg(scip, " created subgroup detection graph using %d of the permutations\n", nusedperms);
3723
3724 if ( nusedperms == npermsincomp )
3725 allpermsused = TRUE;
3726
3727 assert( graphcomponents != NULL );
3728 assert( graphcompbegins != NULL );
3729 assert( compcolorbegins != NULL );
3730 assert( ngraphcomponents > 0 );
3731 assert( ncompcolors > 0 );
3732 assert( nusedperms <= ntwocycleperms );
3733 assert( ncompcolors < propdata->npermvars );
3734
3735 if ( nusedperms == 0 )
3736 {
3737 SCIPdebugMsg(scip, " -> skipping component, since less no permutation was used\n");
3738
3739 SCIPfreeBufferArray(scip, &permused);
3740 SCIPfreeBufferArray(scip, &usedperms);
3741
3742 goto FREEBASICMEM;
3743 }
3744
3745 SCIPdebugMsg(scip, " number of different colors in the graph: %d\n", ncompcolors);
3746
3747 if ( propdata->addstrongsbcs || propdata->addweaksbcs )
3748 {
3749 SCIP_CALL( SCIPallocBufferArray(scip, &chosencomppercolor, ncompcolors) );
3750 SCIP_CALL( SCIPallocBufferArray(scip, &firstvaridxpercolor, ncompcolors) );
3751
3752 /* Initialize the arrays with -1 to encode that we have not added orbitopes/strong SBCs
3753 * yet. In case we do not modify this entry, no weak inequalities are added based on
3754 * this component.
3755 */
3756 for (j = 0; j < ncompcolors; ++j)
3757 {
3758 chosencomppercolor[j] = -1;
3759 firstvaridxpercolor[j] = -1;
3760 }
3761 }
3762
3763 norbitopesincomp = getNOrbitopesInComp(propdata->permvars, graphcomponents, graphcompbegins, compcolorbegins,
3764 ncompcolors, nvarsincomponent[cidx]);
3765
3766 /* if there is just one orbitope satisfying the requirements, handle the full component by symresacks */
3767 if ( norbitopesincomp == 1 )
3768 {
3769 int k;
3770
3771 for (k = 0; k < npermsincomp; ++k)
3772 {
3773 SCIP_CONS* cons;
3774 char name[SCIP_MAXSTRLEN];
3775
3776 (void) SCIPsnprintf(name, SCIP_MAXSTRLEN, "symresack_comp%d_perm%d", cidx, k);
3777
3779 propdata->perms[propdata->components[propdata->componentbegins[cidx] + k]],
3780 propdata->permvars, propdata->npermvars, FALSE,
3781 propdata->conssaddlp, TRUE, FALSE, TRUE, TRUE, FALSE, FALSE, FALSE, FALSE, FALSE) );
3782 SCIP_CALL( SCIPaddCons(scip, cons));
3783
3784 /* do not release constraint here - will be done later */
3786 &propdata->genorbconsssize, propdata->ngenorbconss + 1) );
3787 propdata->genorbconss[propdata->ngenorbconss++] = cons;
3788 ++propdata->nsymresacks;
3789
3790 if ( ! propdata->componentblocked[cidx] )
3791 {
3792 propdata->componentblocked[cidx] |= SYM_HANDLETYPE_SYMBREAK;
3793 ++propdata->ncompblocked;
3794 }
3795
3796 SCIPdebugMsg(scip, " add symresack for permutation %d of component %d\n", k, cidx);
3797 }
3798
3799 goto FREEMEMORY;
3800 }
3801
3802 for (j = 0; j < ncompcolors; ++j)
3803 {
3804 int nbinarycomps = 0;
3805 int largestcolorcomp = -1;
3806 int largestcompsize = 0;
3807 int k;
3808 SCIP_Bool isorbitope = TRUE;
3809 SCIP_Bool orbitopeadded = FALSE;
3810 SCIP_Bool useorbitope;
3811#ifdef SCIP_DEBUG
3812 SCIP_Bool binaffected = FALSE;
3813 SCIP_Bool intaffected = FALSE;
3814 SCIP_Bool contaffected = FALSE;
3815#endif
3816
3817 /* skip trivial components */
3818 if ( graphcompbegins[compcolorbegins[j+1]] - graphcompbegins[compcolorbegins[j]] < 2 )
3819 {
3820 if( chosencomppercolor != NULL )
3821 chosencomppercolor[j] = -1;
3822
3823 ++ntrivialcolors;
3824 continue;
3825 }
3826
3827 SCIPdebugMsg(scip, " color %d has %d components with overall %d variables\n", j, compcolorbegins[j+1] - compcolorbegins[j],
3828 graphcompbegins[compcolorbegins[j+1]] - graphcompbegins[compcolorbegins[j]]);
3829
3830 /* check whether components of this color might build an orbitope (with > 2 columns) */
3831 for (k = compcolorbegins[j]; k < compcolorbegins[j+1]; ++k)
3832 {
3833 SCIP_VAR* firstvar;
3834 int compsize;
3835
3836 compsize = graphcompbegins[k+1] - graphcompbegins[k];
3837
3838 /* the first component that we are looking at for this color */
3839 if ( largestcompsize < 1 )
3840 {
3841 if ( compsize < 3 )
3842 {
3843 isorbitope = FALSE;
3844 break;
3845 }
3846
3847 largestcompsize = compsize;
3848 largestcolorcomp = k;
3849 }
3850 else if ( compsize != largestcompsize )
3851 {
3852 /* variable orbits (compsize) have not the same size, cannot define orbitope */
3853 isorbitope = FALSE;
3854 break;
3855 }
3856
3857 firstvar = propdata->permvars[graphcomponents[graphcompbegins[k]]];
3858
3859 /* count number of binary orbits (comps) */
3860 if ( SCIPvarIsBinary(firstvar) )
3861 ++nbinarycomps;
3862 }
3863
3864#ifdef SCIP_DEBUG
3865 for (k = compcolorbegins[j]; k < compcolorbegins[j+1]; ++k)
3866 {
3867 SCIP_VAR* firstvar;
3868
3869 firstvar = propdata->permvars[graphcomponents[graphcompbegins[k]]];
3870
3871 if ( SCIPvarIsBinary(firstvar) )
3872 binaffected = TRUE;
3873 else if (SCIPvarIsIntegral(firstvar) )
3874 intaffected = TRUE;
3875 else
3876 contaffected = TRUE;
3877 }
3878
3879 SCIPdebugMsg(scip, " affected types (bin,int,cont): (%d,%d,%d)\n", binaffected, intaffected, contaffected);
3880#endif
3881
3882 /* only use the orbitope if there are binary rows */
3883 useorbitope = FALSE;
3884 if ( norbitopesincomp > 0 && nbinarycomps > 0 )
3885 useorbitope = TRUE;
3886
3887 if ( isorbitope && useorbitope )
3888 {
3889 int firstvaridx;
3890 int chosencomp;
3891
3892 SCIPdebugMsg(scip, " detected an orbitope with %d rows and %d columns\n", nbinarycomps, largestcompsize);
3893
3894 assert( nbinarycomps > 0 );
3895 assert( largestcompsize > 2 );
3896
3897 /* add the orbitope constraint for this color
3898 *
3899 * It might happen that we cannot generate the orbitope matrix if the orbitope is not generated by permutations
3900 * all having the same number of 2-cycles, e.g., the orbitope generated by (1,2)(4,5), (2,3), (5,6).
3901 */
3902 SCIP_CALL( addOrbitopeSubgroup(scip, propdata, usedperms, nusedperms, compcolorbegins,
3903 graphcompbegins, graphcomponents, j, nbinarycomps, largestcompsize, &firstvaridx, &chosencomp,
3904 &lexorder, &nvarslexorder, &maxnvarslexorder, allpermsused, &orbitopeadded) );
3905
3906 if ( orbitopeadded )
3907 {
3908 if ( propdata->addstrongsbcs || propdata->addweaksbcs )
3909 {
3910 assert( chosencomppercolor != NULL );
3911 assert( firstvaridxpercolor != NULL );
3912
3913 /* adapt the first variable per color to be compatible with the created orbiope (upper left variable) */
3914 assert( compcolorbegins[j] <= chosencomp && chosencomp < compcolorbegins[j+1] );
3915 assert( 0 <= firstvaridx && firstvaridx < propdata->npermvars );
3916
3917 chosencomppercolor[j] = chosencomp;
3918 firstvaridxpercolor[j] = firstvaridx;
3919 }
3920
3921 if ( ! propdata->componentblocked[cidx] )
3922 {
3923 propdata->componentblocked[cidx] |= SYM_HANDLETYPE_SYMBREAK;
3924 ++propdata->ncompblocked;
3925 }
3926
3927#ifdef SCIP_DEBUG
3928 ++norbitopes;
3929#endif
3930 }
3931 }
3932
3933 /* if no (useable) orbitope was found, possibly add strong SBCs */
3934 if ( propdata->addstrongsbcs && ! orbitopeadded )
3935 {
3936 assert( largestcolorcomp >= 0 );
3937 assert( largestcolorcomp < ngraphcomponents );
3938 assert( largestcompsize > 0 );
3939
3940 if( propdata->addweaksbcs )
3941 {
3942 assert( chosencomppercolor != NULL );
3943 assert( firstvaridxpercolor != NULL );
3944
3945 chosencomppercolor[j] = largestcolorcomp;
3946 firstvaridxpercolor[j] = graphcomponents[graphcompbegins[largestcolorcomp]];
3947 }
3948
3949 SCIPdebugMsg(scip, " choosing component %d with %d variables and adding strong SBCs\n",
3950 largestcolorcomp, graphcompbegins[largestcolorcomp+1] - graphcompbegins[largestcolorcomp]);
3951
3952 /* add the strong SBCs for the corresponding component */
3953 SCIP_CALL( addStrongSBCsSubgroup(scip, propdata, graphcompbegins, graphcomponents, largestcolorcomp,
3954 propdata->addsymresacks, &lexorder, &nvarslexorder, &maxnvarslexorder) );
3955
3956 /* store whether symmetries on non-binary symmetries have been handled */
3957 if ( ! SCIPvarIsBinary(propdata->permvars[graphcomponents[graphcompbegins[largestcolorcomp]]]) )
3958 handlednonbinarysymmetry = TRUE;
3959
3960 if ( ! propdata->componentblocked[cidx] )
3961 {
3962 propdata->componentblocked[cidx] |= SYM_HANDLETYPE_SYMBREAK;
3963 ++propdata->ncompblocked;
3964 }
3965
3966#ifdef SCIP_DEBUG
3967 nstrongsbcs += graphcompbegins[largestcolorcomp+1] - graphcompbegins[largestcolorcomp] - 1;
3968#endif
3969 }
3970 else if ( ! orbitopeadded )
3971 {
3972 SCIPdebugMsg(scip, " no useable orbitope found and no SBCs added\n");
3973
3974 /* mark the color as not handled */
3975 if ( propdata->addweaksbcs )
3976 {
3977 assert( chosencomppercolor != NULL );
3978 chosencomppercolor[j] = -1; /*lint !e613*/
3979 }
3980 }
3981 }
3982
3983 SCIPdebugMsg(scip, " skipped %d trivial colors\n", ntrivialcolors);
3984
3985 /* possibly add weak SBCs for enclosing orbit of first component */
3986 if ( propdata->addweaksbcs && propdata->componentblocked[cidx] && nusedperms < npermsincomp )
3987 {
3988 int naddedconss;
3989
3990 assert( firstvaridxpercolor != NULL );
3991 assert( chosencomppercolor != NULL );
3992
3993 SCIP_CALL( addWeakSBCsSubgroup(scip, propdata, compcolorbegins, graphcompbegins,
3994 graphcomponents, ncompcolors, chosencomppercolor, firstvaridxpercolor,
3995 cidx, &naddedconss, propdata->addsymresacks, &lexorder, &nvarslexorder, &maxnvarslexorder) );
3996
3997 assert( naddedconss < propdata->npermvars );
3998
3999#ifdef SCIP_DEBUG
4000 nweaksbcs += naddedconss;
4001#endif
4002 }
4003 else
4004 SCIPdebugMsg(scip, " don't add weak sbcs because all generators were used or the settings forbid it\n");
4005
4006 /* if suborbitopes or strong group actions have been found, potentially add symresacks adapted to
4007 * variable order given by lexorder if no symmetries on non-binary variables have been handled
4008 */
4009 if ( nvarslexorder > 0 && propdata->addsymresacks && ! handlednonbinarysymmetry )
4010 {
4011 int k;
4012
4013 SCIP_CALL( adaptSymmetryDataSST(scip, propdata->perms, modifiedperms, propdata->nperms,
4014 propdata->permvars, modifiedpermvars, propdata->npermvars, lexorder, nvarslexorder) );
4015
4016 for (k = 0; k < npermsincomp; ++k)
4017 {
4018 SCIP_CONS* cons;
4019 char name[SCIP_MAXSTRLEN];
4020 int* symresackperm;
4021 SCIP_Bool actsonbinary = FALSE;
4022
4023 /* skip permutations that have been used to build an orbitope */
4024 if ( permused[k] )
4025 continue;
4026
4027 /* skip permutations that do not act on binary variables */
4028 symresackperm = modifiedperms[propdata->components[propdata->componentbegins[cidx] + k]];
4029 for (j = 0; j < propdata->nperms && !actsonbinary; ++j)
4030 {
4031 if ( symresackperm[j] != j && SCIPvarIsBinary(modifiedpermvars[j]) )
4032 actsonbinary = TRUE;
4033 }
4034
4035 if ( ! actsonbinary )
4036 continue;
4037
4038 (void) SCIPsnprintf(name, SCIP_MAXSTRLEN, "symresack_comp%d_perm%d", cidx, k);
4039
4041 symresackperm, modifiedpermvars, propdata->npermvars, FALSE,
4042 propdata->conssaddlp, TRUE, FALSE, TRUE, TRUE, FALSE, FALSE, FALSE, FALSE, FALSE) );
4043 SCIP_CALL( SCIPaddCons(scip, cons));
4044
4045 /* do not release constraint here - will be done later */
4047 &propdata->genorbconsssize, propdata->ngenorbconss + 1) );
4048 propdata->genorbconss[propdata->ngenorbconss++] = cons;
4049 ++propdata->nsymresacks;
4050
4051 if ( ! propdata->componentblocked[cidx] )
4052 {
4053 propdata->componentblocked[cidx] |= SYM_HANDLETYPE_SYMBREAK;
4054 ++propdata->ncompblocked;
4055 }
4056
4057 SCIPdebugMsg(scip, " add symresack for permutation %d of component %d adapted to suborbitope lexorder\n", k, cidx);
4058 }
4059 }
4060
4061 FREEMEMORY:
4062 SCIPfreeBlockMemoryArrayNull(scip, &lexorder, maxnvarslexorder);
4063
4064 SCIPfreeBufferArrayNull(scip, &firstvaridxpercolor);
4065 SCIPfreeBufferArrayNull(scip, &chosencomppercolor);
4066 SCIPfreeBlockMemoryArrayNull(scip, &compcolorbegins, ncompcolors + 1);
4067 SCIPfreeBlockMemoryArrayNull(scip, &graphcompbegins, ngraphcomponents + 1);
4068 SCIPfreeBlockMemoryArrayNull(scip, &graphcomponents, propdata->npermvars);
4069 SCIPfreeBufferArrayNull(scip, &permused);
4070 SCIPfreeBufferArrayNull(scip, &usedperms);
4071
4072#ifdef SCIP_DEBUG
4073 SCIPdebugMsg(scip, "total number of added (sub-)orbitopes: %d\n", norbitopes);
4074 SCIPdebugMsg(scip, "total number of added strong sbcs: %d\n", nstrongsbcs);
4075 SCIPdebugMsg(scip, "total number of added weak sbcs: %d\n", nweaksbcs);
4076#endif
4077
4078 FREEBASICMEM:
4079 SCIPfreeBufferArray(scip, &nvarsincomponent);
4080
4081 SCIPfreeBufferArray(scip, &modifiedpermvars);
4082 for (p = propdata->nperms - 1; p >= 0; --p)
4083 {
4084 SCIPfreeBufferArray(scip, &modifiedperms[p]);
4085 }
4086 SCIPfreeBufferArray(scip, &modifiedperms);
4087 SCIPfreeBufferArray(scip, &genorder);
4088
4089 return SCIP_OKAY;
4090}
4091
4092
4093/*
4094 * Functions for symmetry constraints
4095 */
4096
4097
4098/** update symmetry information of conflict graph */
4099static
4101 SCIP* scip, /**< SCIP instance */
4102 SCIP_CONFLICTDATA* varconflicts, /**< conflict structure */
4103 SCIP_VAR** conflictvars, /**< variables encoded in conflict structure */
4104 int nconflictvars, /**< number of nodes/vars in conflict structure */
4105 int* orbits, /**< array of non-trivial orbits */
4106 int* orbitbegins, /**< array containing begin positions of new orbits in orbits array */
4107 int norbits /**< number of non-trivial orbits */
4108 )
4109{
4110 int i;
4111 int j;
4112 int ii;
4113 int jj;
4114 int r; /* r from orbit, the orbit index. */
4115
4116 assert( scip != NULL );
4117 assert( varconflicts != NULL );
4118 assert( conflictvars != NULL );
4119 assert( nconflictvars > 0 );
4120 assert( orbits != NULL );
4121 assert( orbitbegins != NULL );
4122 assert( norbits >= 0 );
4123
4124 /* initialize/reset variable information of nodes in conflict graph */
4125 for (i = 0; i < nconflictvars; ++i)
4126 {
4127 /* (re-)set node data */
4128 varconflicts[i].orbitidx = -1;
4129 varconflicts[i].nconflictinorbit = 0;
4130 varconflicts[i].orbitsize = -1;
4131 varconflicts[i].posinorbit = -1;
4132 }
4133
4134 /* add orbit information to nodes of conflict graph */
4135 for (r = 0; r < norbits; ++r)
4136 {
4137 int posinorbit = 0;
4138 int orbitsize;
4139
4140 orbitsize = orbitbegins[r + 1] - orbitbegins[r];
4141 assert( orbitsize >= 0 );
4142
4143 for (i = orbitbegins[r]; i < orbitbegins[r + 1]; ++i)
4144 {
4145 int pos;
4146
4147 /* get variable and position in conflict graph */
4148 pos = orbits[i];
4149 assert( pos < nconflictvars );
4150 assert( varconflicts[pos].var == conflictvars[pos] );
4151
4152 varconflicts[pos].orbitidx = r;
4153 varconflicts[pos].nconflictinorbit = 0;
4154 varconflicts[pos].orbitsize = orbitsize;
4155 varconflicts[pos].posinorbit = posinorbit++;
4156 }
4157
4158 /* determine nconflictsinorbit
4159 *
4160 * For each pair of active variables in this orbit, check if it is part of a conflict clique.
4161 * Use that we store the cliques of this type in varconflicts[pos].cliques.
4162 * These lists are sorted (by the address of the constraint), so we only need to check for each i, j in the orbit
4163 * whether they are contained in the same clique.
4164 */
4165 for (i = orbitbegins[r]; i < orbitbegins[r + 1]; ++i)
4166 {
4167 ii = orbits[i];
4168 assert( varconflicts[ii].orbitidx == r );
4169
4170 /* skip inactive variables */
4171 if ( ! varconflicts[ii].active )
4172 continue;
4173
4174 for (j = i + 1; j < orbitbegins[r + 1]; ++j)
4175 {
4176 jj = orbits[j];
4177 assert( varconflicts[jj].orbitidx == r );
4178
4179 /* skip inactive variables */
4180 if ( ! varconflicts[jj].active )
4181 continue;
4182
4183 /* Check if i and j are overlapping in some clique, where only one of the two could have value 1.
4184 * Use that cliques are sorted by the constraint address.
4185 *
4186 * @todo A better sorted order would be: First constraints with large variables (higher hitting probability)
4187 * and then by a unique constraint identifier (address, or conspos).
4188 */
4189 if ( checkSortedArraysHaveOverlappingEntry((void**)varconflicts[ii].cliques,
4190 varconflicts[ii].ncliques, (void**)varconflicts[jj].cliques, varconflicts[jj].ncliques,
4191 sortByPointerValue) )
4192 {
4193 /* there is overlap! */
4194 ++varconflicts[ii].nconflictinorbit;
4195 ++varconflicts[jj].nconflictinorbit;
4196 }
4197 }
4198 }
4199 }
4200
4201 return SCIP_OKAY;
4202}
4203
4204
4205/** create conflict graph either for symmetric or for all variables
4206 *
4207 * This routine just creates the graph, but does not add (symmetry) information to its nodes.
4208 * This has to be done separately by the routine updateSymInfoConflictGraphSST().
4209 *
4210 * The function returns with varconflicts as NULL when we do not create it.
4211 */
4212static
4214 SCIP* scip, /**< SCIP instance */
4215 SCIP_CONFLICTDATA** varconflicts, /**< pointer to store the variable conflict data */
4216 SCIP_VAR** conflictvars, /**< array of variables to encode in conflict graph */
4217 int nconflictvars, /**< number of vars to encode in conflict graph */
4218 SCIP_HASHMAP* conflictvarmap /**< map of variables to indices in conflictvars array */
4219 )
4220{
4221 SCIP_CLIQUE** cliques;
4222 SCIP_VAR** cliquevars;
4223 SCIP_CLIQUE* clique;
4224 int* tmpncliques;
4225 int ncliques;
4226 int ncliquevars;
4227 int node;
4228 int c;
4229 int i;
4230
4231#ifdef SCIP_DEBUG
4232 int varncliques = 0;
4233#endif
4234
4235 assert( scip != NULL );
4236 assert( varconflicts != NULL );
4237 assert( conflictvars != NULL );
4238 assert( nconflictvars > 0 );
4239
4240 /* we set the pointer of varconflicts to NULL to illustrate that we didn't generate it */
4241 *varconflicts = NULL;
4242
4243 /* get cliques for creating conflict structure */
4244
4245 cliques = SCIPgetCliques(scip);
4246 ncliques = SCIPgetNCliques(scip);
4247 if ( ncliques == 0 )
4248 {
4249 SCIPdebugMsg(scip, "No cliques present --> construction of conflict structure aborted.\n");
4250 return SCIP_OKAY;
4251 }
4252
4253 /* construct variable conflicts */
4254 SCIPdebugMsg(scip, "Construction of conflict structure:\n");
4255 SCIP_CALL( SCIPallocBlockMemoryArray(scip, varconflicts, nconflictvars) );
4256 for (i = 0; i < nconflictvars; ++i)
4257 {
4258 (*varconflicts)[i].ncliques = 0;
4259 (*varconflicts)[i].active = TRUE;
4260 (*varconflicts)[i].var = conflictvars[i];
4261 /* set remaining variable conflictdata at neutral entries */
4262 (*varconflicts)[i].cliques = NULL;
4263 (*varconflicts)[i].orbitidx = -1;
4264 (*varconflicts)[i].nconflictinorbit = 0;
4265 (*varconflicts)[i].orbitsize = -1;
4266 (*varconflicts)[i].posinorbit = -1;
4267 }
4268
4269 /* Store, for each variable, the conflict cliques it is contained in.
4270 * In three steps:
4271 * (1.) Count the number of cliques it's contained in, per var, then
4272 * (2.) Create the array of this size, and
4273 * (3.) Fill the array with the cliques.
4274 * Starting with (1.):
4275 */
4276 for (c = 0; c < ncliques; ++c)
4277 {
4278 clique = cliques[c];
4279 assert( clique != NULL );
4280
4281 cliquevars = SCIPcliqueGetVars(clique);
4282 ncliquevars = SCIPcliqueGetNVars(clique);
4283 assert( cliquevars != NULL );
4284 assert( ncliquevars > 0 );
4285
4286 SCIPdebugMsg(scip, "\tIdentify edges for clique ID: %u; Index: %d).\n", SCIPcliqueGetId(clique),
4287 SCIPcliqueGetIndex(clique));
4288
4289 /* for all variables, list which cliques it is part of */
4290 for (i = 0; i < ncliquevars; ++i)
4291 {
4292 node = SCIPhashmapGetImageInt(conflictvarmap, cliquevars[i]);
4293
4294 /* skip variables not in the conflictvars array (so not in hashmap, too) */
4295 if ( node == INT_MAX )
4296 continue;
4297 assert( node >= 0 );
4298 assert( node < nconflictvars );
4299
4300 assert( (*varconflicts)[node].var == cliquevars[i] );
4301 (*varconflicts)[node].active = TRUE;
4302 (*varconflicts)[node].ncliques++;
4303 }
4304 }
4305
4306 /* (2.) allocate the arrays */
4307 for (i = 0; i < nconflictvars; ++i)
4308 {
4309 assert( (*varconflicts)[i].ncliques >= 0 );
4310 assert( (*varconflicts)[i].cliques == NULL );
4311 if ( (*varconflicts)[i].ncliques > 0 )
4312 {
4313 SCIP_CALL( SCIPallocBlockMemoryArray(scip, &(*varconflicts)[i].cliques, (*varconflicts)[i].ncliques) );
4314 }
4315 }
4316
4317 /* (3.) fill the clique constraints */
4318 SCIP_CALL( SCIPallocClearBufferArray(scip, &tmpncliques, nconflictvars) );
4319 for (c = 0; c < ncliques; ++c)
4320 {
4321 clique = cliques[c];
4322 assert( clique != NULL );
4323
4324 cliquevars = SCIPcliqueGetVars(clique);
4325 ncliquevars = SCIPcliqueGetNVars(clique);
4326 assert( cliquevars != NULL );
4327 assert( ncliquevars > 0 );
4328
4329 SCIPdebugMsg(scip, "\tAdd edges for clique ID: %u; Index: %d).\n", SCIPcliqueGetId(clique),
4330 SCIPcliqueGetIndex(clique));
4331
4332 /* for all variables, list which cliques it is part of */
4333 for (i = 0; i < ncliquevars; ++i)
4334 {
4335 node = SCIPhashmapGetImageInt(conflictvarmap, cliquevars[i]);
4336
4337 /* skip variables not in the conflictvars array (so not in hashmap, too) */
4338 if ( node == INT_MAX )
4339 continue;
4340
4341 assert( node >= 0 );
4342 assert( node < nconflictvars );
4343 assert( (*varconflicts)[node].var == cliquevars[i] );
4344
4345 /* add clique to the cliques */
4346 assert( tmpncliques[node] < (*varconflicts)[node].ncliques );
4347 assert( (*varconflicts)[node].cliques != NULL );
4348 (*varconflicts)[node].cliques[tmpncliques[node]++] = clique;
4349
4350#ifdef SCIP_DEBUG
4351 varncliques++;
4352#endif
4353 }
4354 }
4355
4356 /* sort the variable cliques by the address, so checkSortedArraysHaveOverlappingEntry can detect intersections */
4357 for (i = 0; i < nconflictvars; ++i)
4358 {
4359 SCIPsortPtr((void**)(*varconflicts)[i].cliques, sortByPointerValue, (*varconflicts)[i].ncliques);
4360 }
4361
4362#ifndef NDEBUG
4363 for (i = 0; i < nconflictvars; ++i)
4364 {
4365 assert( tmpncliques[i] == (*varconflicts)[i].ncliques );
4366 }
4367#endif
4368
4369 SCIPfreeBufferArray(scip, &tmpncliques);
4370
4371#ifdef SCIP_DEBUG
4372 SCIPdebugMsg(scip, "Construction of conflict graph terminated; %d variable-clique combinations detected.\n",
4373 varncliques);
4374#endif
4375
4376 return SCIP_OKAY;
4377}
4378
4379/** frees conflict graph */
4380static
4382 SCIP* scip, /**< SCIP instance */
4383 SCIP_CONFLICTDATA** varconflicts, /**< conflict graph */
4384 int nvars /**< number of nodes in conflict graph */
4385 )
4386{
4387 int i;
4388 int n;
4389
4390 assert( scip != NULL );
4391 assert( varconflicts != NULL );
4392 assert( *varconflicts != NULL );
4393 assert( nvars >= 0 );
4394
4395 for (i = nvars - 1; i >= 0; --i)
4396 {
4397 n = (*varconflicts)[i].ncliques;
4398 SCIPfreeBlockMemoryArray(scip, &(*varconflicts)[i].cliques, n);
4399 }
4400 SCIPfreeBlockMemoryArray(scip, varconflicts, nvars);
4401
4402 return SCIP_OKAY;
4403}
4404
4405
4406/** adds symresack constraints */
4407static
4409 SCIP* scip, /**< SCIP instance */
4410 SCIP_PROPDATA* propdata, /**< data of symmetry propagator */
4411 int cidx /**< index of component to be handled */
4412 )
4413{ /*lint --e{641}*/
4414 int* components;
4415 int* componentbegins;
4416 SCIP_VAR** permvars;
4417 SCIP_Bool conssaddlp;
4418 int** modifiedperms = NULL;
4419 SCIP_VAR** modifiedpermvars = NULL;
4420 int** perms;
4421 int nsymresackcons = 0;
4422 int npermvars;
4423 int nperms;
4424 int i;
4425 int p;
4426
4427 assert( scip != NULL );
4428 assert( propdata != NULL );
4429 assert( propdata->npermvars >= 0 );
4430 assert( propdata->nbinpermvars >= 0 );
4431
4432 /* if no symmetries on binary variables are present */
4433 if ( propdata->nbinpermvars == 0 )
4434 {
4435 assert( propdata->binvaraffected == 0 );
4436 return SCIP_OKAY;
4437 }
4438
4439 perms = propdata->perms;
4440 nperms = propdata->nperms;
4441 permvars = propdata->permvars;
4442 npermvars = propdata->npermvars;
4443 conssaddlp = propdata->conssaddlp;
4444 components = propdata->components;
4445 componentbegins = propdata->componentbegins;
4446
4447 assert( nperms <= 0 || perms != NULL );
4448 assert( permvars != NULL );
4449 assert( npermvars > 0 );
4450 assert( components != NULL );
4451 assert( componentbegins != NULL );
4452 assert( 0 <= cidx && cidx < propdata->ncomponents );
4453
4454 /* exit if component is already blocked by incompatible methods */
4455 if ( propdata->componentblocked[cidx] & (~SYM_HANDLETYPE_SST) )
4456 return SCIP_OKAY;
4457 if ( (propdata->componentblocked[cidx] & SYM_HANDLETYPE_SST) )
4458 {
4459 /* the leader must be binary for compatability */
4460 if ( (ISSSTINTACTIVE(propdata->sstleadervartype)
4461 || ISSSTIMPLINTACTIVE(propdata->sstleadervartype)
4462 || ISSSTCONTACTIVE(propdata->sstleadervartype)) )
4463 return SCIP_OKAY;
4464 }
4465
4466 /* skip component if it has signed permutations */
4467 if ( propdata->componenthassignedperm[cidx] )
4468 return SCIP_OKAY;
4469
4470 /* adapt natural variable order to a variable order that is compatible with Schreier Sims constraints */
4471 if ( propdata->nleaders > 0 && ISSSTBINACTIVE(propdata->sstleadervartype) )
4472 {
4473 SCIP_CALL( SCIPallocBufferArray(scip, &modifiedperms, nperms) );
4474 for (p = 0; p < nperms; ++p)
4475 {
4476 SCIP_CALL( SCIPallocBufferArray(scip, &modifiedperms[p], npermvars) );
4477 }
4478 SCIP_CALL( SCIPallocBufferArray(scip, &modifiedpermvars, npermvars) );
4479
4480 for (i = 0; i < npermvars; ++i)
4481 modifiedpermvars[i] = permvars[i];
4482
4483 SCIP_CALL( adaptSymmetryDataSST(scip, perms, modifiedperms, nperms, permvars, modifiedpermvars, npermvars,
4484 propdata->leaders, propdata->nleaders) );
4485 }
4486
4487 /* loop through perms in component cidx and add symresack constraints */
4488 for (p = componentbegins[cidx]; p < componentbegins[cidx + 1]; ++p)
4489 {
4490 SCIP_CONS* cons;
4491 int permidx;
4492 char name[SCIP_MAXSTRLEN];
4493
4494 permidx = components[p];
4495
4496 (void) SCIPsnprintf(name, SCIP_MAXSTRLEN, "symbreakcons_component%d_perm%d", cidx, permidx);
4497
4498 /* adapt permutation to leader */
4499 if ( propdata->nleaders > 0 && ISSSTBINACTIVE(propdata->sstleadervartype) )
4500 {
4501 assert( (propdata->componentblocked[cidx] & SYM_HANDLETYPE_SST) != 0 );
4502 assert( modifiedperms != NULL );
4503 assert( modifiedpermvars != NULL );
4504
4505 SCIP_CALL( SCIPcreateSymbreakCons(scip, &cons, name, modifiedperms[permidx], modifiedpermvars, npermvars, FALSE,
4506 conssaddlp, TRUE, FALSE, TRUE, TRUE, FALSE, FALSE, FALSE, FALSE, FALSE) );
4507 }
4508 else
4509 {
4510 SCIP_CALL( SCIPcreateSymbreakCons(scip, &cons, name, perms[permidx], permvars, npermvars, FALSE,
4511 conssaddlp, TRUE, FALSE, TRUE, TRUE, FALSE, FALSE, FALSE, FALSE, FALSE) );
4512 }
4513 propdata->componentblocked[cidx] |= SYM_HANDLETYPE_SYMBREAK;
4514 SCIP_CALL( SCIPaddCons(scip, cons) );
4515
4516 /* do not release constraint here - will be done later */
4518 &propdata->genorbconsssize, propdata->ngenorbconss + 1) );
4519 propdata->genorbconss[propdata->ngenorbconss++] = cons;
4520 ++propdata->nsymresacks;
4521 ++nsymresackcons;
4522 }
4523
4524 if ( propdata->nleaders > 0 && ISSSTBINACTIVE(propdata->sstleadervartype) )
4525 {
4526 assert( modifiedperms != NULL );
4527 assert( modifiedpermvars != NULL );
4528
4529 SCIPfreeBufferArray(scip, &modifiedpermvars);
4530 for (p = nperms - 1; p >= 0; --p)
4531 {
4532 SCIPfreeBufferArray(scip, &modifiedperms[p]);
4533 }
4534 SCIPfreeBufferArray(scip, &modifiedperms);
4535 }
4536
4537 SCIPdebugMsg(scip, "Added %d symresack constraints.\n", nsymresackcons);
4538
4539 return SCIP_OKAY;
4540}
4541
4542
4543/** add Schreier Sims constraints for a specific orbit and update Schreier Sims table */
4544static
4546 SCIP* scip, /**< SCIP instance */
4547 SCIP_CONFLICTDATA* varconflicts, /**< conflict graph or NULL if useconflictgraph == FALSE */
4548 SCIP_PROPDATA* propdata, /**< data of symmetry propagator */
4549 SCIP_VAR** permvars, /**< permvars array */
4550 int* orbits, /**< symmetry orbits */
4551 int* orbitbegins, /**< array storing begin position for each orbit */
4552 int orbitidx, /**< index of orbit for Schreier Sims constraints */
4553 int orbitleaderidx, /**< index of leader variable for Schreier Sims constraints */
4554 SCIP_Shortbool* orbitvarinconflict, /**< indicator whether orbitvar is in conflict with orbit leader */
4555 int norbitvarinconflict,/**< number of variables in conflict with orbit leader */
4556 int* nchgbds /**< pointer to store number of bound changes (or NULL) */
4557 )
4558{ /*lint --e{613,641}*/
4559 SCIP_CONS* cons;
4560 char name[SCIP_MAXSTRLEN];
4561 SCIP_VAR* vars[2];
4562 SCIP_Real vals[2];
4563 int orbitsize;
4564 int posleader;
4565 int poscur;
4566 int ncuts = 0;
4567 SCIP_Bool addcuts = FALSE;
4568 int i;
4569#ifndef NDEBUG
4570 int j;
4571#endif
4572
4573 assert( scip != NULL );
4574 assert( propdata != NULL );
4575 assert( permvars != NULL );
4576 assert( orbits != NULL );
4577 assert( orbitbegins != NULL );
4578 assert( orbitidx >= 0 );
4579 assert( orbitleaderidx >= 0 );
4580 assert( orbitvarinconflict != NULL || varconflicts == NULL );
4581 assert( norbitvarinconflict >= 0 );
4582 assert( nchgbds != NULL );
4583
4584 orbitsize = orbitbegins[orbitidx + 1] - orbitbegins[orbitidx];
4585
4586 /* variables in conflict with leader are fixed and not treated by a cut; trailing -1 to not count the leader */
4587 if ( propdata->sstaddcuts )
4588 addcuts = TRUE;
4589 else if ( propdata->sstleaderrule == SCIP_LEADERRULE_MAXCONFLICTSINORBIT
4590 || propdata->ssttiebreakrule == SCIP_LEADERTIEBREAKRULE_MAXCONFLICTSINORBIT )
4591 addcuts = propdata->addconflictcuts;
4592
4593 if ( addcuts )
4594 ncuts = orbitsize - norbitvarinconflict - 1;
4595
4596 /* (re-)allocate memory for Schreier Sims constraints and leaders */
4597 if ( ncuts > 0 )
4598 {
4599 if ( propdata->nsstconss == 0 )
4600 {
4601 assert( propdata->sstconss == NULL );
4602 assert( propdata->maxnsstconss == 0 );
4603 propdata->maxnsstconss = 2 * ncuts;
4604 SCIP_CALL( SCIPallocBlockMemoryArray(scip, &(propdata->sstconss), propdata->maxnsstconss) );
4605 }
4606 else if ( propdata->nsstconss + ncuts > propdata->maxnsstconss )
4607 {
4608 int newsize;
4609
4610 newsize = SCIPcalcMemGrowSize(scip, propdata->maxnsstconss + 2 * ncuts);
4611 SCIP_CALL( SCIPreallocBlockMemoryArray(scip, &(propdata->sstconss),
4612 propdata->maxnsstconss, newsize) );
4613 propdata->maxnsstconss = newsize;
4614 }
4615 }
4616
4617 if ( propdata->nleaders == 0 )
4618 {
4619 propdata->maxnleaders = MIN(propdata->nperms, propdata->npermvars);
4620 SCIP_CALL( SCIPallocBlockMemoryArray(scip, &(propdata->leaders), propdata->maxnleaders) );
4621 }
4622 assert( propdata->nleaders < propdata->maxnleaders );
4623
4624 /* add Schreier Sims constraints vars[0] >= vars[1], where vars[0] is always the leader */
4625 posleader = orbitbegins[orbitidx] + orbitleaderidx;
4626 vars[0] = permvars[orbits[posleader]];
4627 vals[0] = -1.0;
4628 vals[1] = 1.0;
4629 propdata->leaders[propdata->nleaders++] = orbits[posleader];
4630 *nchgbds = 0;
4631 for (i = 0, poscur = orbitbegins[orbitidx]; i < orbitsize; ++i, ++poscur)
4632 {
4633 if ( i == orbitleaderidx )
4634 {
4635 assert( orbitvarinconflict == NULL || ! orbitvarinconflict[i] );
4636 continue;
4637 }
4638
4639 vars[1] = permvars[orbits[poscur]];
4640#ifndef NDEBUG
4641 for (j = 0; j < propdata->nleaders - 1; ++j)
4642 {
4643 assert( propdata->leaders[j] != orbits[poscur] );
4644 }
4645#endif
4646
4647 /* if the i-th variable in the orbit is in a conflict with the leader, fix it to 0 */
4648 if ( varconflicts != NULL )
4649 {
4650 if ( orbitvarinconflict[i] )
4651 {
4652 assert( SCIPvarIsBinary(vars[1]) );
4653 assert( SCIPvarGetLbLocal(vars[1]) < 0.5 );
4654 assert( varconflicts != NULL );
4655
4656 /* if variable is fixed */
4657 if ( SCIPvarGetUbLocal(vars[1]) > 0.5 )
4658 {
4659 SCIP_CALL( SCIPchgVarUb(scip, vars[1], 0.0) );
4660 ++(*nchgbds);
4661
4662 /* deactivate the fixed variable (cannot contribute to a conflict anymore) */
4663 assert( varconflicts[orbits[poscur]].active );
4664 varconflicts[orbits[poscur]].active = FALSE;
4665 }
4666
4667 /* reset value */
4668 orbitvarinconflict[i] = FALSE;
4669 }
4670 else if ( addcuts )
4671 {
4672 (void) SCIPsnprintf(name, SCIP_MAXSTRLEN, "SSTcut_%d_%d", orbits[posleader], orbits[poscur]);
4673 SCIP_CALL( SCIPcreateConsLinear(scip, &cons, name, 2, vars, vals, - SCIPinfinity(scip), 0.0,
4675
4676 SCIP_CALL( SCIPaddCons(scip, cons) );
4677 propdata->sstconss[propdata->nsstconss++] = cons;
4678 }
4679 }
4680 else if ( addcuts )
4681 {
4682 (void) SCIPsnprintf(name, SCIP_MAXSTRLEN, "SSTcut_%d_%d", orbits[posleader], orbits[poscur]);
4683 SCIP_CALL( SCIPcreateConsLinear(scip, &cons, name, 2, vars, vals, - SCIPinfinity(scip), 0.0,
4685
4686 SCIP_CALL( SCIPaddCons(scip, cons) );
4687 propdata->sstconss[propdata->nsstconss++] = cons;
4688 }
4689 }
4690
4691 return SCIP_OKAY;
4692}
4693
4694
4695/** selection rule of next orbit/leader in orbit for Schreier Sims constraints */
4696static
4698 SCIP* scip, /**< SCIP instance */
4699 SCIP_CONFLICTDATA* varconflicts, /**< variable conflicts structure, or NULL if we do not use it */
4700 SCIP_VAR** conflictvars, /**< variables encoded in conflict graph */
4701 int nconflictvars, /**< number of variables encoded in conflict graph */
4702 int* orbits, /**< orbits of stabilizer subgroup, expressed in terms of conflictvars */
4703 int* orbitbegins, /**< array storing the begin position of each orbit in orbits */
4704 int norbits, /**< number of orbits */
4705 int leaderrule, /**< rule to select leader */
4706 int tiebreakrule, /**< tie break rule to select leader */
4707 SCIP_VARTYPE leadervartype, /**< variable type of leader */
4708 int* orbitidx, /**< pointer to index of selected orbit */
4709 int* leaderidx, /**< pointer to leader in orbit */
4710 SCIP_Shortbool* orbitvarinconflict, /**< array to store whether a var in the orbit is conflicting with leader */
4711 int* norbitvarinconflict,/**< pointer to store number of vars in the orbit in conflict with leader */
4712 SCIP_Bool* success /**< pointer to store whether orbit cut be selected successfully */
4713 )
4714{
4715 int varidx;
4716 int orbitcriterion;
4717 int curcriterion = INT_MIN;
4718 int orbitsize;
4719 int i;
4720 int leader = -1;
4721
4722 assert( scip != NULL );
4723 assert( conflictvars != NULL );
4724 assert( nconflictvars > 0 );
4725 assert( orbits != NULL );
4726 assert( orbitbegins != NULL );
4727 assert( norbits > 0 );
4728 assert( orbitidx != NULL );
4729 assert( leaderidx != NULL );
4730 assert( orbitvarinconflict != NULL || varconflicts == NULL );
4731 assert( norbitvarinconflict != NULL );
4732 assert( success != NULL );
4733
4734 *orbitidx = 0;
4735 *leaderidx = 0;
4736 *norbitvarinconflict = 0;
4737 *success = FALSE;
4738
4739 /* terminate if leader or tiebreak rule cannot be checked */
4740 if ( varconflicts == NULL && (leaderrule == (int) SCIP_LEADERRULE_MAXCONFLICTSINORBIT
4741 || tiebreakrule == (int) SCIP_LEADERTIEBREAKRULE_MAXCONFLICTSINORBIT) )
4742 return SCIP_OKAY;
4743
4744 /* select the leader and its orbit */
4745 if ( leaderrule == (int) SCIP_LEADERRULE_FIRSTINORBIT || leaderrule == (int) SCIP_LEADERRULE_LASTINORBIT )
4746 {
4747 orbitcriterion = INT_MIN;
4748
4749 /* iterate over orbits and select the first one that meets the tiebreak rule */
4750 for (i = 0; i < norbits; ++i)
4751 {
4752 /* skip orbits containing vars different to the leader's vartype */
4753 /* Conflictvars is permvars! */
4754 if ( SCIPvarGetType(conflictvars[orbits[orbitbegins[i]]]) != leadervartype )
4755 continue;
4756
4757 if ( tiebreakrule == (int) SCIP_LEADERTIEBREAKRULE_MINORBIT )
4758 curcriterion = orbitbegins[i] - orbitbegins[i + 1];
4759 else if ( tiebreakrule == (int) SCIP_LEADERTIEBREAKRULE_MAXORBIT )
4760 curcriterion = orbitbegins[i + 1] - orbitbegins[i];
4761 else
4762 {
4763 assert( tiebreakrule == (int) SCIP_LEADERTIEBREAKRULE_MAXCONFLICTSINORBIT );
4764
4765 /* get first or last active variable in orbit */
4766 if ( leaderrule == (int) SCIP_LEADERRULE_FIRSTINORBIT )
4767 {
4768 int cnt;
4769
4770 cnt = orbitbegins[i];
4771
4772 do
4773 {
4774 varidx = orbits[cnt++];
4775 }
4776 while ( SCIPvarGetProbindex(conflictvars[varidx]) == -1 && cnt < orbitbegins[i + 1]);
4777 }
4778 else
4779 {
4780 int cnt;
4781
4782 cnt = orbitbegins[i + 1] - 1;
4783
4784 do
4785 {
4786 varidx = orbits[cnt--];
4787 }
4788 while ( SCIPvarGetProbindex(conflictvars[varidx]) == -1 && cnt >= orbitbegins[i]);
4789 }
4790
4791 /* skip inactive variables */
4792 if ( SCIPvarGetProbindex(conflictvars[varidx]) == -1 )
4793 continue;
4794
4795 assert( varconflicts[varidx].orbitidx == i );
4796 /* coverity[var_deref_op] */
4797 curcriterion = varconflicts[varidx].nconflictinorbit;
4798 }
4799
4800 /* update selected orbit */
4801 if ( curcriterion > orbitcriterion )
4802 {
4803 orbitcriterion = curcriterion;
4804 *orbitidx = i;
4805 *success = TRUE;
4806
4807 if ( leaderrule == (int) SCIP_LEADERRULE_FIRSTINORBIT )
4808 *leaderidx = 0;
4809 else
4810 *leaderidx = orbitbegins[i + 1] - orbitbegins[i] - 1;
4811 }
4812 }
4813
4814 /* store variables in conflict with leader */
4815 if ( *success && varconflicts != NULL )
4816 {
4817 leader = orbits[orbitbegins[*orbitidx] + *leaderidx];
4818 assert( leader < nconflictvars );
4819
4820 if ( tiebreakrule == (int) SCIP_LEADERTIEBREAKRULE_MAXCONFLICTSINORBIT
4821 && varconflicts[leader].ncliques > 0 )
4822 {
4823 /* count how many active variables in the orbit conflict with "leader"
4824 * This is only needed if there are possible conflicts.
4825 */
4826 int varmapid;
4827
4828 orbitsize = orbitbegins[*orbitidx + 1] - orbitbegins[*orbitidx];
4829 assert( varconflicts != NULL );
4830 assert( leader >= 0 && leader < nconflictvars );
4831
4832 assert( orbitvarinconflict != NULL );
4833
4834 for (i = 0; i < orbitsize; ++i)
4835 {
4836 /* skip the leader */
4837 if ( i == *leaderidx )
4838 continue;
4839
4840 /* get variable index in conflict graph */
4841 varmapid = orbits[orbitbegins[*orbitidx] + i];
4842
4843 /* only active variables */
4844 if ( ! varconflicts[varmapid].active )
4845 continue;
4846
4847 /* check if leader and var have overlap */
4848 if ( checkSortedArraysHaveOverlappingEntry((void**)varconflicts[leader].cliques,
4849 varconflicts[leader].ncliques, (void**)varconflicts[varmapid].cliques,
4850 varconflicts[varmapid].ncliques, sortByPointerValue) )
4851 {
4852 /* there is overlap! */
4853 orbitvarinconflict[i] = TRUE;
4854 ++(*norbitvarinconflict);
4855 }
4856 }
4857 }
4858 }
4859 }
4860 else
4861 {
4862 /* only three possible values for leaderrules, so it must be MAXCONFLICTSINORBIT
4863 * In this case, the code must have computed the conflict graph.
4864 */
4865 assert( leaderrule == (int) SCIP_LEADERRULE_MAXCONFLICTSINORBIT );
4866 assert( varconflicts != NULL );
4867
4868 orbitcriterion = 0;
4869
4870 /* iterate over variables and select the first one that meets the tiebreak rule */
4871 for (i = 0; i < nconflictvars; ++i)
4872 {
4873 /* skip vars different to the leader's vartype */
4874 if ( SCIPvarGetType(conflictvars[i]) != leadervartype )
4875 continue;
4876
4877 /* skip variables not affected by symmetry */
4878 /* coverity[var_deref_op] */
4879 if ( varconflicts[i].orbitidx == -1 )
4880 continue;
4881
4882 curcriterion = varconflicts[i].nconflictinorbit;
4883
4884 if ( curcriterion > orbitcriterion )
4885 {
4886 orbitcriterion = curcriterion;
4887 *orbitidx = varconflicts[i].orbitidx;
4888 *leaderidx = varconflicts[i].posinorbit;
4889 *success = TRUE;
4890 }
4891 }
4892
4893 /* store variables in conflict with leader */
4894 leader = orbits[orbitbegins[*orbitidx] + *leaderidx];
4895 assert( leader < nconflictvars );
4896 assert( norbitvarinconflict != NULL );
4897
4898 if ( *success && varconflicts[leader].ncliques > 0 )
4899 {
4900 /* count how many active variables in the orbit conflict with leader */
4901 int varmapid;
4902
4903 orbitsize = orbitbegins[*orbitidx + 1] - orbitbegins[*orbitidx];
4904 assert( varconflicts != NULL );
4905 assert( leader >= 0 && leader < nconflictvars );
4906
4907 assert( orbitvarinconflict != NULL );
4908
4909 for (i = 0; i < orbitsize; ++i)
4910 {
4911 /* skip the leader */
4912 if ( i == *leaderidx )
4913 continue;
4914
4915 /* get variable index in conflict graph */
4916 varmapid = orbits[orbitbegins[*orbitidx] + i];
4917 /* only active variables */
4918 if ( ! varconflicts[varmapid].active )
4919 continue;
4920
4921 /* check if leader and var have overlap */
4922 if ( checkSortedArraysHaveOverlappingEntry((void**)varconflicts[leader].cliques,
4923 varconflicts[leader].ncliques, (void**)varconflicts[varmapid].cliques,
4924 varconflicts[varmapid].ncliques, sortByPointerValue) )
4925 {
4926 /* there is overlap! */
4927 orbitvarinconflict[i] = TRUE;
4928 ++(*norbitvarinconflict);
4929 }
4930 }
4931 }
4932 }
4933
4934 return SCIP_OKAY;
4935}
4936
4937
4938/** add Schreier Sims constraints to the problem */
4939static
4941 SCIP* scip, /**< SCIP instance */
4942 SCIP_PROPDATA* propdata, /**< data of symmetry propagator */
4943 SCIP_Bool onlywithcontvars, /**< only handle components that contain continuous variables with SST */
4944 int* nchgbds, /**< pointer to store number of bound changes (or NULL) */
4945 int cidx /**< index of component which shall be handled */
4946 )
4947{ /*lint --e{641}*/
4948 SCIP_CONFLICTDATA* varconflicts = NULL;
4949 SCIP_HASHMAP* permvarmap;
4950 SCIP_VAR** permvars;
4951 int** permstrans;
4952 int npermvars;
4953 int nmovedpermvars;
4954 int nmovedbinpermvars;
4955 int nmovedintpermvars;
4956 int nmovedimplintpermvars;
4957 int nmovedcontpermvars;
4958 int nperms;
4959
4960 int* orbits;
4961 int* orbitbegins;
4962 int norbits;
4963 int* components;
4964 int* componentbegins;
4965 int* vartocomponent;
4966 int ncomponents;
4967 unsigned* componentblocked;
4968
4969 int orbitidx;
4970 int orbitleaderidx;
4971 SCIP_Shortbool* orbitvarinconflict = NULL;
4972 int norbitvarinconflict;
4973 SCIP_Shortbool* inactiveperms;
4974 int ninactiveperms;
4975 int posleader;
4976 int leaderrule;
4977 int tiebreakrule;
4978 int leadervartype;
4980 int nvarsselectedtype;
4981 SCIP_Bool conflictgraphcreated = FALSE;
4982 SCIP_Bool mixedcomponents;
4983 int norbitleadercomponent;
4984 int* perm;
4985 SCIP_VARTYPE vartype;
4986
4987 int i;
4988 int c;
4989 int p;
4990 SCIP_Bool success = TRUE;
4991
4992 assert( scip != NULL );
4993 assert( propdata != NULL );
4994 assert( propdata->computedsymmetry );
4995
4996 permvars = propdata->permvars;
4997 npermvars = propdata->npermvars;
4998 nperms = propdata->nperms;
4999 assert( permvars != NULL );
5000 assert( npermvars > 0 );
5001 assert( nperms > 0 );
5002
5004 permvarmap = propdata->permvarmap;
5005 assert( permvarmap != NULL );
5006
5008 permstrans = propdata->permstrans;
5009 assert( permstrans != NULL );
5010
5011 components = propdata->components;
5012 componentbegins = propdata->componentbegins;
5013 componentblocked = propdata->componentblocked;
5014 vartocomponent = propdata->vartocomponent;
5015 ncomponents = propdata->ncomponents;
5016
5017 assert( components != NULL );
5018 assert( componentbegins != NULL );
5019 assert( vartocomponent != NULL );
5020 assert( componentblocked != NULL );
5021 assert( ncomponents > 0 );
5022 assert( 0 <= cidx && cidx < ncomponents );
5023
5024 /* exit if component is blocked */
5025 if ( componentblocked[cidx] )
5026 return SCIP_OKAY;
5027
5028 /* skip component if it has signed permutations */
5029 if ( propdata->componenthassignedperm[cidx] )
5030 return SCIP_OKAY;
5031
5032 leaderrule = propdata->sstleaderrule;
5033 tiebreakrule = propdata->ssttiebreakrule;
5034 leadervartype = propdata->sstleadervartype;
5035 mixedcomponents = propdata->sstmixedcomponents;
5036
5037 /* if not already computed, get number of affected vars */
5039 nmovedpermvars = propdata->nmovedpermvars;
5040 nmovedbinpermvars = propdata->nmovedbinpermvars;
5041 nmovedintpermvars = propdata->nmovedintpermvars;
5042 nmovedimplintpermvars = propdata->nmovedimplintpermvars;
5043 nmovedcontpermvars = propdata->nmovedcontpermvars;
5044 assert( nmovedpermvars > 0 ); /* nperms > 0 implies this */
5045
5046 /* determine the leader's vartype */
5047 nvarsselectedtype = 0;
5048 if ( ISSSTBINACTIVE(leadervartype) && nmovedbinpermvars > nvarsselectedtype )
5049 {
5050 selectedtype = SCIP_VARTYPE_BINARY;
5051 nvarsselectedtype = nmovedbinpermvars;
5052 }
5053
5054 if ( ISSSTINTACTIVE(leadervartype) && nmovedintpermvars > nvarsselectedtype )
5055 {
5056 selectedtype = SCIP_VARTYPE_INTEGER;
5057 nvarsselectedtype = nmovedintpermvars;
5058 }
5059
5060 if ( ISSSTIMPLINTACTIVE(leadervartype) && nmovedimplintpermvars > nvarsselectedtype )
5061 {
5062 selectedtype = SCIP_VARTYPE_IMPLINT;
5063 nvarsselectedtype = nmovedimplintpermvars;
5064 }
5065
5066 if ( ISSSTCONTACTIVE(leadervartype) && nmovedcontpermvars > nvarsselectedtype )
5067 {
5068 selectedtype = SCIP_VARTYPE_CONTINUOUS;
5069 nvarsselectedtype = nmovedcontpermvars;
5070 }
5071
5072 /* terminate if no variables of a possible leader type is affected */
5073 if ( nvarsselectedtype == 0 )
5074 return SCIP_OKAY;
5075
5076 /* ignore this component if no continuous variables are contained */
5077 if ( onlywithcontvars )
5078 {
5079 for (p = componentbegins[cidx]; p < componentbegins[cidx + 1]; ++p)
5080 {
5081 perm = propdata->perms[p];
5082 for (i = 0; i < propdata->npermvars; ++i)
5083 {
5084 if ( perm[i] == i )
5085 continue;
5086 vartype = SCIPvarGetType(propdata->permvars[i]);
5087 if ( vartype == SCIP_VARTYPE_CONTINUOUS || vartype == SCIP_VARTYPE_IMPLINT )
5088 goto COMPONENTOK;
5089 }
5090 }
5091 /* loop terminated naturally, so component does not have continuous or implicitly integer variables. */
5092 return SCIP_OKAY;
5093
5094 COMPONENTOK:
5095 ;
5096 }
5097
5098 /* @todo online create the conflict graph for the variable in the current component */
5099 /* possibly create conflict graph; graph is not created if no cliques are present */
5100 if ( selectedtype == SCIP_VARTYPE_BINARY && (leaderrule == SCIP_LEADERRULE_MAXCONFLICTSINORBIT
5102 {
5103 SCIP_CALL( createConflictGraphSST(scip, &varconflicts, permvars, npermvars, permvarmap) );
5104 conflictgraphcreated = varconflicts != NULL;
5105 }
5106
5107 /* allocate data structures necessary for orbit computations and conflict graph */
5108 SCIP_CALL( SCIPallocBufferArray(scip, &inactiveperms, nperms) );
5109 SCIP_CALL( SCIPallocBufferArray(scip, &orbits, npermvars) );
5110 SCIP_CALL( SCIPallocBufferArray(scip, &orbitbegins, npermvars) );
5111
5112 if ( conflictgraphcreated )
5113 {
5114 SCIP_CALL( SCIPallocClearBufferArray(scip, &orbitvarinconflict, npermvars) );
5115 }
5116
5117 SCIPdebugMsg(scip, "Start selection of orbits and leaders for Schreier Sims constraints.\n");
5118 SCIPdebugMsg(scip, "orbitidx\tleaderidx\torbitsize\n");
5119
5120 if ( nchgbds != NULL )
5121 *nchgbds = 0;
5122
5123 /* initialize array indicating whether permutations shall not be considered for orbit permutations */
5124 for (c = 0; c < ncomponents; ++c)
5125 {
5126 for (p = componentbegins[c]; p < componentbegins[c + 1]; ++p)
5127 {
5128 if ( c == cidx )
5129 inactiveperms[components[p]] = FALSE;
5130 else
5131 inactiveperms[components[p]] = TRUE;
5132 }
5133 }
5134 ninactiveperms = nperms - componentbegins[cidx + 1] + componentbegins[cidx];
5135
5136 /* as long as the stabilizer is non-trivial, add Schreier Sims constraints */
5137 norbitleadercomponent = 0;
5138 while ( ninactiveperms < nperms )
5139 {
5140 int nchanges = 0;
5141
5142 /* compute orbits w.r.t. active perms */
5143 SCIP_CALL( SCIPcomputeOrbitsFilterSym(scip, npermvars, permstrans, nperms, inactiveperms,
5144 orbits, orbitbegins, &norbits, components, componentbegins, vartocomponent,
5145 componentblocked, ncomponents, nmovedpermvars) );
5146
5147 /* stop if we require pure components and a component contains variables of different types */
5148 if ( ! mixedcomponents )
5149 {
5150 for (p = 0; p < norbits; ++p)
5151 {
5152 /* stop if the first element of an orbits has the wrong vartype */
5153 if ( SCIPvarGetType(permvars[orbits[orbitbegins[p]]]) != selectedtype )
5154 {
5155 success = FALSE;
5156 break;
5157 }
5158 }
5159 }
5160
5161 if ( ! success )
5162 break;
5163
5164 /* update symmetry information of conflict graph */
5165 if ( conflictgraphcreated )
5166 {
5167 SCIP_CALL( updateSymInfoConflictGraphSST(scip, varconflicts, permvars, npermvars, orbits, orbitbegins,
5168 norbits) );
5169 }
5170
5171 /* possibly adapt the leader and tie-break rule */
5172 if ( leaderrule == SCIP_LEADERRULE_MAXCONFLICTSINORBIT && ! conflictgraphcreated )
5173 leaderrule = SCIP_LEADERRULE_FIRSTINORBIT;
5174 if ( leaderrule == SCIP_LEADERRULE_MAXCONFLICTSINORBIT && selectedtype != SCIP_VARTYPE_BINARY )
5175 leaderrule = SCIP_LEADERRULE_FIRSTINORBIT;
5176 if ( tiebreakrule == SCIP_LEADERTIEBREAKRULE_MAXCONFLICTSINORBIT && ! conflictgraphcreated )
5177 tiebreakrule = SCIP_LEADERTIEBREAKRULE_MAXORBIT;
5178 if ( tiebreakrule == SCIP_LEADERTIEBREAKRULE_MAXCONFLICTSINORBIT && selectedtype != SCIP_VARTYPE_BINARY )
5179 tiebreakrule = SCIP_LEADERTIEBREAKRULE_MAXORBIT;
5180
5181 /* select orbit and leader */
5182 SCIP_CALL( selectOrbitLeaderSSTConss(scip, varconflicts, permvars, npermvars, orbits, orbitbegins,
5183 norbits, propdata->sstleaderrule, propdata->ssttiebreakrule, selectedtype, &orbitidx, &orbitleaderidx,
5184 orbitvarinconflict, &norbitvarinconflict, &success) );
5185
5186 if ( ! success )
5187 break;
5188
5189 assert( 0 <= orbitidx && orbitidx < norbits );
5190 assert( 0 <= orbitleaderidx && orbitleaderidx < orbitbegins[orbitidx + 1] - orbitbegins[orbitidx] );
5191 SCIPdebugMsg(scip, "%d\t\t%d\t\t%d\n", orbitidx, orbitleaderidx, orbitbegins[orbitidx + 1] - orbitbegins[orbitidx]);
5192
5193 /* add Schreier Sims constraints for the selected orbit and update Schreier Sims table */
5194 SCIP_CALL( addSSTConssOrbitAndUpdateSST(scip, varconflicts, propdata, permvars,
5195 orbits, orbitbegins, orbitidx, orbitleaderidx, orbitvarinconflict, norbitvarinconflict, &nchanges) );
5196
5197 ++norbitleadercomponent;
5198
5199 if ( nchgbds != NULL )
5200 *nchgbds += nchanges;
5201
5202 /* deactivate permutations that move the orbit leader */
5203 posleader = orbits[orbitbegins[orbitidx] + orbitleaderidx];
5204 for (p = 0; p < nperms; ++p)
5205 {
5206 if ( inactiveperms[p] )
5207 continue;
5208
5209 if ( permstrans[posleader][p] != posleader )
5210 {
5211 inactiveperms[p] = TRUE;
5212 ++ninactiveperms;
5213 }
5214 }
5215 }
5216
5217 /* if Schreier Sims constraints have been added, store that Schreier Sims has been used for this component */
5218 if ( norbitleadercomponent > 0 )
5219 componentblocked[cidx] |= SYM_HANDLETYPE_SST;
5220
5221 if ( conflictgraphcreated )
5222 {
5223 SCIPfreeBufferArray(scip, &orbitvarinconflict);
5224 }
5225 SCIPfreeBufferArray(scip, &orbitbegins);
5226 SCIPfreeBufferArray(scip, &orbits);
5227 if ( varconflicts != NULL )
5228 {
5229 /* nconflictvars at construction is npermvars */
5230 SCIP_CALL( freeConflictGraphSST(scip, &varconflicts, npermvars) );
5231 }
5232 SCIPfreeBufferArray(scip, &inactiveperms);
5233
5234 return SCIP_OKAY;
5235}
5236
5237
5238/** orbitopal reduction */
5239static
5241 SCIP* scip, /**< SCIP instance */
5242 SCIP_PROPDATA* propdata, /**< propdata */
5243 int id, /**< ID for orbitope constraint (needed for name) */
5244 int** varidxmatrix, /**< matrix containing variable indices in orbitope matrix */
5245 int nrows, /**< number of rows of orbitope */
5246 int ncols, /**< number of columns of orbitope */
5247 SCIP_Bool* success /**< pointer to store whether orbitope could be added successfully */
5248 )
5249{
5250 char name[SCIP_MAXSTRLEN];
5251 int i;
5252 int j;
5253
5254 SCIP_Bool ispporbitope;
5255 SCIP_VAR*** varmatrix;
5256 SCIP_Bool* pprows;
5257 int npprows;
5258 SCIP_ORBITOPETYPE type;
5259
5260 assert( scip != NULL );
5261 assert( propdata != NULL );
5262 assert( propdata->usedynamicprop );
5263 assert( varidxmatrix != NULL );
5264 assert( nrows > 0 );
5265 assert( ncols > 0 );
5266 assert( success != NULL );
5267
5268 *success = FALSE;
5269
5270 /* add linear constraints x_1 >= x_2 >= ... >= x_ncols for single-row orbitopes */
5271 if ( nrows == 1 )
5272 {
5273 /* restrict to the packing and partitioning rows */
5274 SCIP_CONS* cons;
5275 SCIP_VAR* consvars[2];
5276 SCIP_Real conscoefs[2] = { -1.0, 1.0 };
5277
5278 /* for all adjacent column pairs, add linear constraint */
5280 &propdata->genlinconsssize, propdata->ngenlinconss + ncols - 1) );
5281 for (i = 0; i < ncols - 1; ++i)
5282 {
5283 (void) SCIPsnprintf(name, SCIP_MAXSTRLEN, "orbitope_1row_comp_%d_col%d", id, i);
5284
5285 consvars[0] = propdata->permvars[varidxmatrix[0][i]];
5286 consvars[1] = propdata->permvars[varidxmatrix[0][i + 1]];
5287
5288 SCIP_CALL( SCIPcreateConsLinear(scip, &cons, name, 2, consvars, conscoefs, -SCIPinfinity(scip), 0.0,
5289 propdata->conssaddlp, propdata->conssaddlp, TRUE, TRUE, TRUE, FALSE, FALSE, FALSE, FALSE, FALSE ) );
5290
5291 SCIP_CALL( SCIPaddCons(scip, cons) );
5292 propdata->genlinconss[propdata->ngenlinconss++] = cons;
5293 }
5294
5295 *success = TRUE;
5296 return SCIP_OKAY;
5297 }
5298
5299 /* for only 2 columns, the the component can be completely handled by lexicographic reduction */
5300 if ( ncols == 2 && propdata->lexreddata != NULL )
5301 {
5302 int* orbisackperm;
5303
5304 /* If the component is an orbitope with 2 columns, then there is 1 generator of order 2. */
5305 orbisackperm = propdata->perms[propdata->components[propdata->componentbegins[id]]];
5306
5308 propdata->permvars, propdata->npermvars, orbisackperm, (SYM_SYMTYPE) propdata->symtype,
5309 propdata->permvardomaincenter, TRUE, success) );
5310 if ( *success )
5311 return SCIP_OKAY;
5312 }
5313
5314 /* create orbitope variable matrix */
5315 SCIP_CALL( SCIPallocBufferArray(scip, &varmatrix, nrows) );
5316 for (i = 0; i < nrows; ++i)
5317 {
5318 SCIP_CALL( SCIPallocBufferArray(scip, &varmatrix[i], ncols) );
5319 for (j = 0; j < ncols; ++j)
5320 varmatrix[i][j] = propdata->permvars[varidxmatrix[i][j]];
5321 }
5322
5323 pprows = NULL;
5324 SCIP_CALL( SCIPisPackingPartitioningOrbitope(scip, varmatrix, nrows, ncols, &pprows, &npprows, &type) );
5325
5326 /* does it have at least 3 packing-partitioning rows? */
5327 ispporbitope = npprows >= 3; /* (use same magic number as cons_orbitope.c) */
5328
5329 if ( ispporbitope ) /* @todo if it's a pporbitope, we do it statically right now. */
5330 {
5331 /* restrict to the packing and partitioning rows */
5332 SCIP_CONS* cons;
5333 SCIP_VAR*** ppvarsarrayonlypprows;
5334 int r;
5335
5336 assert( pprows != NULL );
5337
5338 SCIP_CALL( SCIPallocBufferArray(scip, &ppvarsarrayonlypprows, npprows) );
5339
5340 r = 0;
5341 for (i = 0; i < nrows; ++i)
5342 {
5343 if ( pprows[i] )
5344 {
5345 assert( r < npprows );
5346 ppvarsarrayonlypprows[r++] = varmatrix[i];
5347 }
5348 }
5349 assert( r == npprows );
5350
5351 (void) SCIPsnprintf(name, SCIP_MAXSTRLEN, "orbitope_pp_comp_%d", id);
5352 SCIP_CALL( SCIPcreateConsOrbitope(scip, &cons, name, ppvarsarrayonlypprows, SCIP_ORBITOPETYPE_PACKING,
5353 npprows, ncols, FALSE, FALSE, FALSE, FALSE, propdata->conssaddlp,
5355
5356 SCIP_CALL( SCIPaddCons(scip, cons) );
5357
5358 /* check whether we need to resize */
5360 &propdata->genlinconsssize, propdata->ngenlinconss + 1) );
5361 /* @todo we add orbitopes to the dynamically sized array `genlinconss` instead of `genorbconss` to ensure
5362 * compatability with the static orbitope function, which allocates this array statically
5363 */
5364 propdata->genlinconss[propdata->ngenlinconss++] = cons;
5365 *success = TRUE;
5366
5367 SCIPfreeBufferArray(scip, &ppvarsarrayonlypprows);
5368 }
5369 else
5370 {
5371 /* use orbitopal reduction for component */
5372 SCIP_COLUMNORDERING columnordering;
5373 SCIP_VAR** orbitopevarmatrix;
5374 int nelem;
5375 int pos = 0;
5376
5377 /* variable array */
5378 nelem = nrows * ncols;
5379 SCIP_CALL( SCIPallocBufferArray(scip, &orbitopevarmatrix, nelem) );
5380 for (i = 0; i < nrows; ++i)
5381 {
5382 for (j = 0; j < ncols; ++j)
5383 orbitopevarmatrix[pos++] = varmatrix[i][j];
5384 }
5385
5386 /* get column ordering */
5387 columnordering = SCIPorbitopalReductionGetDefaultColumnOrdering(propdata->orbitopalreddata);
5388
5389 (void) SCIPsnprintf(name, SCIP_MAXSTRLEN, "orbitope_full_comp_%d", id);
5390 SCIP_CALL( SCIPorbitopalReductionAddOrbitope(scip, propdata->orbitopalreddata,
5391 SCIP_ROWORDERING_BRANCHING, columnordering,
5392 orbitopevarmatrix, nrows, ncols, success) );
5393 *success = TRUE;
5394
5395 SCIPfreeBufferArray(scip, &orbitopevarmatrix);
5396 }
5397
5398 SCIPfreeBlockMemoryArrayNull(scip, &pprows, nrows);
5399
5400 for (i = nrows - 1; i >= 0; --i)
5401 {
5402 SCIPfreeBufferArray(scip, &varmatrix[i]);
5403 }
5404 SCIPfreeBufferArray(scip, &varmatrix);
5405
5406 return SCIP_OKAY;
5407}
5408
5409
5410/** applies pp-orbitope upgrade if at least 50% of the permutations in a component correspond to pp-orbisacks */
5411static
5413 SCIP* scip, /**< SCIP instance */
5414 SCIP_PROPDATA* propdata, /**< propdata */
5415 int** componentperms, /**< permutations in the component */
5416 int componentsize, /**< number of permutations in the component */
5417 SCIP_Bool hassignedperm, /**< whether the component has a signed permutation */
5418 SCIP_Bool* success /**< whether the packing partitioning upgrade succeeded */
5419 )
5420{
5421 int c;
5422 int i;
5423 int j;
5424 int p;
5425 int* perm;
5426 SCIP_CONSHDLR* setppcconshdlr;
5427 SCIP_CONS** setppcconss;
5428 SCIP_CONS* cons;
5429 SCIP_CONS** setppconsssort;
5430 int nsetppconss;
5431 int nsetppcvars;
5432 SCIP_VAR** setppcvars;
5433 int nsetppcconss;
5434 int** pporbisackperms;
5435 int npporbisackperms;
5436 SCIP_VAR* var;
5437 int varid;
5438 SCIP_CONS*** permvarssetppcconss;
5439 int* npermvarssetppcconss;
5440 int* maxnpermvarssetppcconss;
5441 int maxntwocycles;
5442 int ntwocycles;
5443
5444 assert( scip != NULL );
5445 assert( propdata != NULL );
5446 assert( componentperms != NULL );
5447 assert( componentsize > 0 );
5448 assert( success != NULL );
5449
5450 /* we did not upgrade yet */
5451 *success = FALSE;
5452
5453 /* currently, we cannot handle signed permutations */
5454 if ( hassignedperm )
5455 return SCIP_OKAY;
5456
5457 setppcconshdlr = SCIPfindConshdlr(scip, "setppc");
5458 if ( setppcconshdlr == NULL )
5459 return SCIP_OKAY;
5460
5461 nsetppcconss = SCIPconshdlrGetNConss(setppcconshdlr);
5462 if ( nsetppcconss == 0 )
5463 return SCIP_OKAY;
5464
5465 setppcconss = SCIPconshdlrGetConss(setppcconshdlr);
5466 assert( setppcconss != NULL );
5467
5469
5470 /* collect non-covering constraints and sort by pointer for easy intersection finding */
5471 SCIP_CALL( SCIPallocBufferArray(scip, &setppconsssort, nsetppcconss) );
5472 nsetppconss = 0;
5473 for (c = 0; c < nsetppcconss; ++c)
5474 {
5475 cons = setppcconss[c];
5476
5477 /* only packing or partitioning constraints, no covering types */
5479 continue;
5480
5481 setppconsssort[nsetppconss++] = cons;
5482 }
5483 SCIPsortPtr((void**) setppconsssort, sortByPointerValue, nsetppconss);
5484
5485 /* For each permvar, introduce an array of setppc constraints (initially NULL) for each variable,
5486 * and populate it with the setppc constraints that it contains. This array follows the ordering by cons ptr address.
5487 */
5488 SCIP_CALL( SCIPallocCleanBufferArray(scip, &permvarssetppcconss, propdata->npermvars) );
5489 SCIP_CALL( SCIPallocCleanBufferArray(scip, &npermvarssetppcconss, propdata->npermvars) );
5490 SCIP_CALL( SCIPallocCleanBufferArray(scip, &maxnpermvarssetppcconss, propdata->npermvars) );
5491 for (c = 0; c < nsetppconss; ++c)
5492 {
5493 assert( c >= 0 );
5494 assert( c < nsetppconss );
5495 cons = setppconsssort[c];
5496 assert( cons != NULL );
5497
5498 setppcvars = SCIPgetVarsSetppc(scip, cons);
5499 nsetppcvars = SCIPgetNVarsSetppc(scip, cons);
5500
5501 for (i = 0; i < nsetppcvars; ++i)
5502 {
5503 var = setppcvars[i];
5504 assert( var != NULL );
5505 varid = SCIPhashmapGetImageInt(propdata->permvarmap, (void*) var);
5506 assert( varid == INT_MAX || varid < propdata->npermvars );
5507 assert( varid >= 0 );
5508 if ( varid < propdata->npermvars )
5509 {
5511 &(permvarssetppcconss[varid]), &maxnpermvarssetppcconss[varid], npermvarssetppcconss[varid] + 1) );
5512 assert( npermvarssetppcconss[varid] < maxnpermvarssetppcconss[varid] );
5513 permvarssetppcconss[varid][npermvarssetppcconss[varid]++] = cons;
5514 }
5515 }
5516 }
5517
5518 /* for all permutations, test involutions on binary variables and test if they are captured by setppc conss */
5519 SCIP_CALL( SCIPallocBufferArray(scip, &pporbisackperms, componentsize) );
5520 maxntwocycles = 0;
5521 npporbisackperms = 0;
5522 for (p = 0; p < componentsize; ++p)
5523 {
5524 perm = componentperms[p];
5525 ntwocycles = 0;
5526
5527 /* check if the binary orbits are involutions */
5528 for (i = 0; i < propdata->npermvars; ++i)
5529 {
5530 j = perm[i];
5531
5532 /* ignore fixed points in permutation */
5533 if ( i == j )
5534 continue;
5535 /* only check for situations where i and j are binary variables */
5536 assert( SCIPvarGetType(propdata->permvars[i]) == SCIPvarGetType(propdata->permvars[j]) );
5537 if ( SCIPvarGetType(propdata->permvars[i]) != SCIP_VARTYPE_BINARY )
5538 continue;
5539 /* the permutation must be an involution on binary variables */
5540 if ( perm[j] != i )
5541 goto NEXTPERMITER;
5542 /* i and j are a two-cycle, so we find this once for i and once for j. Only handle this once for i < j. */
5543 if ( i > j )
5544 continue;
5545 /* disqualify permutation if i and j are not in a common set packing constraint */
5546 if ( !checkSortedArraysHaveOverlappingEntry((void**) permvarssetppcconss[i], npermvarssetppcconss[i],
5547 (void**) permvarssetppcconss[j], npermvarssetppcconss[j], sortByPointerValue) )
5548 goto NEXTPERMITER;
5549 ++ntwocycles;
5550 }
5551
5552 /* The permutation qualif