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