Scippy

SCIP

Solving Constraint Integer Programs

sepa_oddcycle.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-2014 Konrad-Zuse-Zentrum */
7 /* fuer Informationstechnik Berlin */
8 /* */
9 /* SCIP is distributed under the terms of the ZIB Academic License. */
10 /* */
11 /* You should have received a copy of the ZIB Academic License */
12 /* along with SCIP; see the file COPYING. If not email to scip@zib.de. */
13 /* */
14 /* * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * */
15 
16 /**@file sepa_oddcycle.c
17  * @brief oddcycle separator
18  * @author Robert Waniek
19  * @author Marc Pfetsch
20  *
21  * We separate odd cycle inequalities in the implication graph. Implemented are the classic method
22  * by Groetschel, Lovasz, and Schrijver (GLS) and the levelgraph method by Hoffman and Padberg (HP)
23  *
24  * Odd cycle inequalities are lifted by a heuristic method based on an idea from Alvarez-Valdes,
25  * Parreno, and Tamarit.
26  *
27  * Some part of this code is based on the odd cycle separator of the program colorbitopt by Marc
28  * Pfetsch.
29  */
30 
31 /*---+----1----+----2----+----3----+----4----+----5----+----6----+----7----+----8----+----9----+----0----+----1----+----2*/
32 
33 #include <assert.h>
34 #include <string.h>
35 
36 #include "scip/sepa_oddcycle.h"
37 #include "scip/pub_misc.h"
38 #include "dijkstra/dijkstra.h"
39 
40 
41 #define SEPA_NAME "oddcycle"
42 #define SEPA_DESC "odd cycle separator"
43 #define SEPA_PRIORITY -15000
44 #define SEPA_FREQ -1
45 #define SEPA_MAXBOUNDDIST 1.0
46 #define SEPA_USESSUBSCIP FALSE /**< does the separator use a secondary SCIP instance? */
47 #define SEPA_DELAY FALSE /**< should separation method be delayed, if other separators found cuts? */
48 
49 
50 /* default values for separator settings */
51 #define DEFAULT_SCALEFACTOR 1000 /**< factor for scaling of the arc-weights in the Dijkstra algorithm */
52 #define DEFAULT_USEGLS TRUE /**< use GLS method, otherwise HP method */
53 #define DEFAULT_LIFTODDCYCLES FALSE /**< lift odd cycle cuts */
54 #define DEFAULT_REPAIRCYCLES TRUE /**< try to repair violated cycles in which a variable and its negated appear */
55 #define DEFAULT_ADDSELFARCS TRUE /**< add links between a variable and its negated */
56 #define DEFAULT_INCLUDETRIANGLES TRUE /**< separate triangles (3-cliques) found as 3-cycles or repaired larger cycles */
57 #define DEFAULT_MULTIPLECUTS FALSE /**< still try variable as start, even if it is already covered by a cut */
58 #define DEFAULT_ALLOWMULTIPLECUTS TRUE /**< allow another inequality to use variable, even if it is already covered */
59 #define DEFAULT_LPLIFTCOEF FALSE /**< TRUE: choose lifting candidate with highest value of coefficient*lpvalue
60  * FALSE: choose lifting candidate with highest coefficient */
61 #define DEFAULT_RECALCLIFTCOEF TRUE /**< whether lifting coefficients should be recomputed */
62 #define DEFAULT_MAXSEPACUTS 5000 /**< maximal number of oddcycle cuts separated per separation round */
63 #define DEFAULT_MAXSEPACUTSROOT 5000 /**< maximal number of oddcycle cuts separated per separation round in root node */
64 #define DEFAULT_PERCENTTESTVARS 0 /**< percent of variables to try the chosen method on [0-100] */
65 #define DEFAULT_OFFSETTESTVARS 100 /**< offset of variables to try the chosen method on */
66 #define DEFAULT_MAXCUTSROOT 1 /**< maximal number of oddcycle cuts generated per root of the levelgraph */
67 #define DEFAULT_SORTSWITCH 3 /**< unsorted (0), maxlp (1), minlp (2), maxfrac (3), minfrac (4) */
68 #define DEFAULT_MAXREFERENCE 0 /**< minimal weight on an edge (in level graph or Dijkstra graph) */
69 #define DEFAULT_MAXROUNDS 10 /**< maximal number of rounds pre node */
70 #define DEFAULT_MAXROUNDSROOT 10 /**< maximal number of rounds in the root node */
71 #define DEFAULT_MAXNLEVELS 20 /**< maximal number of levels in level graph */
72 #define DEFAULT_MAXPERNODESLEVEL 100 /**< maximal percentage of nodes allowed in one level of the levelgraph [0-100] */
73 #define DEFAULT_OFFSETNODESLEVEL 10 /**< additional offset of nodes allowed in one level of the levelgraph */
74 #define DEFAULT_SORTROOTNEIGHBORS TRUE /**< sort neighbors of the root in the level graph */
75 #define DEFAULT_MAXCUTSLEVEL 50 /**< maximal number of cuts produced per level */
76 #define DEFAULT_MAXUNSUCESSFULL 3 /**< maximal number of unsuccessful calls at each node */
77 #define DEFAULT_CUTTHRESHOLD -1 /**< maximal number of other cuts s.t. separation is applied (-1 for direct call) */
78 
79 
80 /*
81  * Data structures
82  */
83 
84 /** Graph structure for level graph
85  *
86  * This graph is tailored to the heuristic search for odd holes, @see separateHeur().
87  *
88  * This undirected graph is represented by a directed graph with forward and backward arcs. Arcs are
89  * forward if they lead from a level l to level l+1, i.e., away from the root; backward arcs
90  * lead from a level l+1 to level l. This distinction enables a fast construction and search
91  * process. In the latter only forward or backward arcs have to be searched.
92  *
93  * Target nodes and weights of the arcs incident to each node (adjacency lists) are stored
94  * consecutively in the arrays targetForward, targetBackward, weightForward, and weightBackward.
95  * The end of each list is marked by a -1 in targetForward and targetBackward.
96  */
97 struct levelGraph
98 {
99  unsigned int nnodes; /**< number of nodes */
100  unsigned int narcs; /**< number of arcs */
101  unsigned int maxnodes; /**< maximal number of nodes of the level graph */
102  unsigned int maxarcs; /**< maximal number of arcs of the level graph */
103  unsigned int nlevels; /**< number of levels completely inserted so far */
104  unsigned int* level; /**< level number for each node */
105  unsigned int lastF; /**< last storage element index in targetForward, weightForward - forward direction */
106  unsigned int lastB; /**< last storage element index in targetBackward, weightBackward - backward direction */
107  int* beginForward; /**< forward adjacency list index in targetForward, weightForward for each node */
108  int* beginBackward; /**< backward adjacency list index in targetBackward, weightBackward for each node */
109  int* targetForward; /**< target nodes of forward arcs */
110  int* targetBackward; /**< target nodes of backward arcs */
111  unsigned int* weightForward; /**< weights of forward arcs */
112  unsigned int* weightBackward; /**< weights of backwards arcs */
113  unsigned int sizeForward; /**< size of targetForward and weightForward */
114  unsigned int sizeBackward; /**< size of targetBackward and weightBackward */
115  int* beginAdj; /**< index of list of arcs inside a level (in sourceAdj) for each node
116  * (the index points at the first arc starting from this node) */
117  unsigned int* sourceAdj; /**< source nodes of arcs inside a level */
118  unsigned int* targetAdj; /**< target nodes of arcs inside a level */
119  unsigned int* weightAdj; /**< weights of arcs inside a level */
120  unsigned int* levelAdj; /**< index of the first arc inside a given level */
121  unsigned int sizeAdj; /**< size of sourceAdj, targetAdj and weightAdj */
122 };
123 
124 typedef struct levelGraph LEVELGRAPH;
125 
126 
127 /** sorting type for starting node or root node iteration order
128  *
129  * If the array should be sorted (1-4), the variable array is sorted every round by the chosen
130  * sorttype and the search method tries the nodes in order of the array. If the array is used
131  * unsorted (0), the search methods tries the nodes in order of the array and stores the last
132  * processed start node or root node and continues from this position in the next separation round.
133  */
134 enum sorttype
135 {
136  UNSORTED = 0, /**< variable array is unsorted */
137  MAXIMAL_LPVALUE = 1, /**< variable array is sorted by maximal lp-value */
138  MINIMAL_LPVALUE = 2, /**< variable array is sorted by minimal fractionality */
139  MAXIMAL_FRACTIONALITY = 3, /**< variable array is sorted by maximal lp-value */
140  MINIMAL_FRACTIONALITY = 4 /**< variable array is sorted by minimal fractionality */
141 };
142 typedef enum sorttype SORTTYPE;
143 
144 
145 /** separator data */
146 struct SCIP_SepaData
147 {
148  int scale; /**< factor for scaling of the arc-weights */
149  unsigned int ncuts; /**< number of cuts, added by the separator so far (in current and past calls) */
150  unsigned int oldncuts; /**< number of cuts at the start the current separation round */
151  int nliftedcuts; /**< number of lifted cuts, added by the separator so far (in current and past calls) */
152  SCIP_Bool usegls; /**< use GLS method, otherwise HP method */
153  SCIP_Bool multiplecuts; /**< an odd cycle cut of length L can be generated L times; forbidding multiple cuts
154  * per node might be faster but might miss some cuts in the current round */
155  SCIP_Bool allowmultiplecuts; /**< allow multiple cuts covering one node */
156  SCIP_Bool liftoddcycles; /**< TRUE, iff we try to lift odd cycle inequalities */
157  SCIP_Bool addselfarcs; /**< add arcs between the nodes of a variable and its negated; since not all implications
158  * are in the graph, this often finds more cycles */
159  SCIP_Bool repaircycles; /**< if a variable and its negated appear in a cycle, we can repair the cycle
160  * by removing both and reconnecting the remaining nodes of the cycle */
161  SCIP_Bool includetriangles; /**< handle triangles found as 3-cycles or repaired larger cycles */
162  LEVELGRAPH* levelgraph; /**< level graph when using HP method, NULL otherwise */
163  DIJKSTRA_GRAPH* dijkstragraph; /**< Dijkstra graph if using method by GLS, NULL otherwise */
164  unsigned int* mapping; /**< mapping for getting the index of a variable in the sorted variable array */
165  SCIP_Bool lpliftcoef; /**< TRUE: choose lifting candidate with highest value of coefficient*lpvalue
166  * FALSE: choose lifting candidate with highest coefficient */
167  SCIP_Bool recalcliftcoef; /**< whether lifting coefficients should be recomputed */
168  int maxsepacuts; /**< max. number of oddcycle cuts separated per separation round */
169  int maxsepacutsroot; /**< max. number of oddcycle cuts separated per separation round in the root node */
170  int maxsepacutsround; /**< max. number of oddcycle cuts separated per separation round in the current node */
171  SORTTYPE sortswitch; /**< sorted type: unsorted (0), maxlp (1), minlp (2), maxfrac (3), minfrac (4) */
172  int lastroot; /**< save root of last GLS-method run */
173  SCIP_Bool sortrootneighbors; /**< sort neighbors of the root in the level graph */
174  int percenttestvars; /**< percentage of variables to try the chosen method on [0-100] */
175  int offsettestvars; /**< offset of variables to try the chosen method on */
176  int maxpernodeslevel; /**< percentage of nodes allowed in the same level of the level graph [0-100] */
177  int offsetnodeslevel; /**< additional offset of nodes allowed in one level of the levelgraph */
178  unsigned int maxlevelsize; /**< maximal number of nodes allowed in the same level of the level graph */
179  int maxcutsroot; /**< maximal number of oddcycle cuts generated per root of the levelgraph */
180  int maxcutslevel; /**< maximal number of oddcycle cuts generated per level of the level graph */
181  int maxrounds; /**< maximal number of oddcycle separation rounds per node (-1: unlimited) */
182  int maxroundsroot; /**< maximal number of oddcycle separation rounds in the root node (-1: unlimited) */
183  int maxreference; /**< minimal weight on an edge (in level graph or Dijkstra graph) */
184  int maxnlevels; /**< maximal number of levels in level graph */
185  int maxunsucessfull; /**< maximal number of unsuccessful calls at each node */
186  int nunsucessfull; /**< number of unsuccessful calls at current node */
187  int cutthreshold; /**< maximal number of other cuts s.t. separation is applied (-1 for direct call) */
188  SCIP_Longint lastnode; /**< number of last node */
189 };
190 
191 
192 /*
193  * debugging methods
194  */
195 
196 #ifdef SCIP_OUTPUT
197 
198 /** displays cycle of pred data structure w.r.t. variable names of the original problem (including
199  * status: original or negated node in graph)
200  */
201 static
202 void printCycle(
203  SCIP_VAR** vars, /**< problem variables */
204  unsigned int* pred, /**< cycle stored as predecessor list */
205  unsigned int nbinvars, /**< number of binary problem variables */
206  unsigned int startnode /**< a node of the cycle */
207  )
208 {
209  unsigned int varsindex;
210  unsigned int counter;
211 
212  assert(vars != NULL);
213  assert(pred != NULL);
214  assert(nbinvars > 0);
215  assert(startnode < 4*nbinvars);
216 
217  counter = 0;
218  varsindex = startnode;
219 
220  /* print start/end node */
221  if( varsindex < nbinvars || ( varsindex >= 2*nbinvars && varsindex < 3*nbinvars ) )
222  {
223  SCIPdebugMessage("+ %s\n",SCIPvarGetName(vars[varsindex%nbinvars]));
224  }
225  else
226  {
227  SCIPdebugMessage("- %s\n",SCIPvarGetName(vars[varsindex%nbinvars]));
228  }
229 
230  /* print inner nodes */
231  for( varsindex = pred[startnode]; varsindex != startnode; varsindex = pred[varsindex] )
232  {
233  if( varsindex < nbinvars || ( varsindex >= 2*nbinvars && varsindex < 3*nbinvars ) )
234  {
235  SCIPdebugMessage("+ %s\n",SCIPvarGetName(vars[varsindex%nbinvars]));
236  }
237  else
238  {
239  SCIPdebugMessage("- %s\n",SCIPvarGetName(vars[varsindex%nbinvars]));
240  }
241  ++counter;
242  }
243 
244  /* print start/end node */
245  if( varsindex < nbinvars || ( varsindex >= 2*nbinvars && varsindex < 3*nbinvars ) )
246  {
247  SCIPdebugMessage("+ %s\n",SCIPvarGetName(vars[varsindex%nbinvars]));
248  }
249  else
250  {
251  SCIPdebugMessage("- %s\n",SCIPvarGetName(vars[varsindex%nbinvars]));
252  }
253 
254  ++counter;
255  SCIPdebugMessage("original cycle has %u variables.\n", counter);
256 }
257 #endif
258 
259 
260 /*
261  * lifting methods
262  */
263 
264 /** using the level graph (if possible) or Dijkstra graph data structure (depending on the used
265  * method) we determine whether two nodes are adjacent
266  */
267 static
268 SCIP_Bool isNeighbor(
269  SCIP_VAR** vars, /**< problem variables */
270  unsigned int nbinvars, /**< number of binary problem variables */
271  SCIP_SEPADATA* sepadata, /**< separator data structure */
272  unsigned int a, /**< node index of first variable */
273  unsigned int b /**< node index of second variable */
274  )
275 {
276  unsigned int i;
277 
278  assert(vars != NULL);
279  assert(nbinvars > 2);
280  assert(sepadata != NULL);
281  assert(sepadata->levelgraph != NULL || sepadata->usegls);
282  assert(sepadata->dijkstragraph != NULL || ! sepadata->usegls);
283  assert(a < 2*nbinvars);
284  assert(b < 2*nbinvars);
285  assert(a != b);
286 
287  /* determine adjacency using the Dijkstra graph */
288  if( sepadata->usegls )
289  {
290  if( sepadata->dijkstragraph->outcnt[a] == 0 || sepadata->dijkstragraph->outcnt[b] == 0 )
291  return FALSE;
292 
293  /* @todo later: if helpful: sort head and weight list once */
294  for( i = sepadata->dijkstragraph->outbeg[a]; i < sepadata->dijkstragraph->outbeg[a] + sepadata->dijkstragraph->outcnt[a]; ++i )
295  {
296  if( sepadata->dijkstragraph->head[i] == b + 2*nbinvars )
297  return TRUE;
298  }
299  }
300  else /* determine adjacency using the level graph */
301  {
302  /* if a and b are contained in the level graph (with their arcs), we can check inside the level graph structure */
303  if( (sepadata->levelgraph->beginForward[a] != -1 || sepadata->levelgraph->beginBackward[a] != -1)
304  && (sepadata->levelgraph->beginForward[b] != -1 || sepadata->levelgraph->beginBackward[b] != -1) )
305  {
306  assert(sepadata->levelgraph->level[a] <= sepadata->levelgraph->nlevels);
307  assert(sepadata->levelgraph->level[b] <= sepadata->levelgraph->nlevels);
308 
309  /* if a and b are not in neighbored levels or the same level, they cannot be adjacent */
310  if( sepadata->levelgraph->level[a] > sepadata->levelgraph->level[b] + 1
311  || sepadata->levelgraph->level[b] > sepadata->levelgraph->level[a] + 1 )
312  return FALSE;
313 
314  assert(sepadata->levelgraph->level[a] == sepadata->levelgraph->level[b]
315  || sepadata->levelgraph->level[a]+1 == sepadata->levelgraph->level[b]
316  || sepadata->levelgraph->level[a] == sepadata->levelgraph->level[b]+1);
317 
318  /* first case of adjacent level */
319  if( sepadata->levelgraph->level[a] == sepadata->levelgraph->level[b]+1 )
320  {
321  if( sepadata->levelgraph->beginBackward[a] >= 0 )
322  {
323  i = (unsigned int) sepadata->levelgraph->beginBackward[a];
324  while( sepadata->levelgraph->targetBackward[i] != -1 )
325  {
326  if( sepadata->levelgraph->targetBackward[i] == (int)b )
327  return TRUE;
328  ++i;
329  }
330  }
331  }
332  else if( sepadata->levelgraph->level[a] == sepadata->levelgraph->level[b]-1 ) /* second case of adjacent level */
333  {
334  if( sepadata->levelgraph->beginForward[a] >= 0 )
335  {
336  i = (unsigned int) sepadata->levelgraph->beginForward[a];
337  while( sepadata->levelgraph->targetForward[i] != -1 )
338  {
339  if( sepadata->levelgraph->targetForward[i] == (int)b )
340  return TRUE;
341  ++i;
342  }
343  }
344  }
345  else /* same level (note that an edge between a and b is stored for a if a < b, otherwise it is stored for b) */
346  {
347  assert(sepadata->levelgraph->level[a] == sepadata->levelgraph->level[b]);
348  assert(sepadata->levelgraph->level[a] > 0); /* root has no neighbor in the same level */
349 
350  if( a < b && sepadata->levelgraph->beginAdj[a] >= 0 )
351  {
352  i = (unsigned int) sepadata->levelgraph->beginAdj[a];
353  assert(i >= sepadata->levelgraph->levelAdj[sepadata->levelgraph->level[a]]);
354 
355  while( sepadata->levelgraph->sourceAdj[i] == a && i < sepadata->levelgraph->levelAdj[sepadata->levelgraph->level[a]+1] )
356  {
357  if( sepadata->levelgraph->targetAdj[i] == b )
358  return TRUE;
359 
360  /* if adjacency list ends we are done and a and b are not adjacent */
361  if( sepadata->levelgraph->sourceAdj[i] == 0 && sepadata->levelgraph->targetAdj[i] == 0 )
362  return FALSE;
363 
364  assert(sepadata->levelgraph->sourceAdj[i] < sepadata->levelgraph->targetAdj[i]);
365  ++i;
366  }
367  }
368  if( b < a && sepadata->levelgraph->beginAdj[b] >= 0 )
369  {
370  i = (unsigned int) sepadata->levelgraph->beginAdj[b];
371  assert(i >= sepadata->levelgraph->levelAdj[sepadata->levelgraph->level[b]]);
372 
373  while( sepadata->levelgraph->sourceAdj[i] == b && i < sepadata->levelgraph->levelAdj[sepadata->levelgraph->level[b]+1])
374  {
375  if( sepadata->levelgraph->targetAdj[i] == a )
376  return TRUE;
377 
378  /* if adjacency list ends we are done and a and b are not adjacent */
379  if( sepadata->levelgraph->sourceAdj[i] == 0 && sepadata->levelgraph->targetAdj[i] == 0 )
380  return FALSE;
381 
382  assert(sepadata->levelgraph->sourceAdj[i] < sepadata->levelgraph->targetAdj[i]);
383  ++i;
384  }
385  }
386  }
387  }
388  /* if a or b is not in the levels already completely inserted in the levelgraph, we check
389  * their adjacency by the SCIP data structures
390  */
391  else
392  {
393  SCIP_Bool originala;
394  SCIP_Bool originalb;
395 
396  unsigned int nbinimpls;
397  SCIP_VAR** implvars;
398  SCIP_BOUNDTYPE* impltypes;
399 #ifndef NDEBUG
400  SCIP_Real* implbounds;
401 #endif
402  unsigned int ncliques;
403  SCIP_CLIQUE** cliques;
404  unsigned int ncliquevars;
405  SCIP_VAR** cliquevars;
406  SCIP_Bool* cliquevals;
407 
408  unsigned int j;
409 
410  /* get original variables */
411  originala = TRUE;
412  if( a >= nbinvars )
413  {
414  a = a - nbinvars;
415  originala = FALSE;
416  }
417  assert(a < nbinvars);
418 
419  originalb = TRUE;
420  if( b >= nbinvars )
421  {
422  b = b - nbinvars;
423  originalb = FALSE;
424  }
425  assert(b < nbinvars);
426 
427  /* nodes cannot be connected by trivial observations */
428  if( ( SCIPvarGetNBinImpls(vars[a], originala) + SCIPvarGetNCliques(vars[a], originala) == 0 )
429  || ( SCIPvarGetNBinImpls(vars[b], originalb) + SCIPvarGetNCliques(vars[b], originalb) == 0 ) )
430  return FALSE;
431  if( ( SCIPvarGetNBinImpls(vars[b], originalb) == 0 && SCIPvarGetNCliques(vars[a], originala) == 0 )
432  || ( SCIPvarGetNBinImpls(vars[a], originala) == 0 && SCIPvarGetNCliques(vars[b], originalb) == 0 ) )
433  return FALSE;
434 
435  /* @todo later: possible improvement: do this test for implications and cliques separately if this here is time consuming */
436  /* one of the nodes seems to have more arcs than the other, we swap them (since adjacency is symmetric) */
437  if( SCIPvarGetNBinImpls(vars[a], originala) + 2 * SCIPvarGetNCliques(vars[a], originala) >
438  SCIPvarGetNBinImpls(vars[b], originalb) + 2 * SCIPvarGetNCliques(vars[b], originalb) )
439  {
440  unsigned int temp;
441  SCIP_Bool varfixingtemp;
442 
443  temp = b;
444  varfixingtemp = originalb;
445  b = a;
446  originalb = originala;
447  a = temp;
448  originala = varfixingtemp;
449  }
450 
451  /* check whether there is an implication a = 1 -> b = 0 */
452  nbinimpls = (unsigned int) SCIPvarGetNBinImpls(vars[a], originala);
453  implvars = SCIPvarGetImplVars(vars[a], originala);
454  impltypes = SCIPvarGetImplTypes(vars[a], originala);
455 #ifndef NDEBUG
456  implbounds = SCIPvarGetImplBounds(vars[a], originala);
457 #endif
458 
459  assert(implvars != NULL || nbinimpls == 0);
460  assert(impltypes != NULL || nbinimpls == 0);
461  assert(implbounds != NULL || nbinimpls == 0);
462 
463  for( i = 0; i < nbinimpls; ++i )
464  {
465  assert( implvars != NULL && impltypes != NULL && implbounds != NULL ); /* for lint */
466  if( SCIPvarGetProbindex(vars[b]) == SCIPvarGetProbindex(implvars[i]) )
467  {
468  if( impltypes[i] == SCIP_BOUNDTYPE_UPPER && originalb == TRUE )
469  {
470  assert(implbounds[i] == 0.0);
471  return TRUE;
472  }
473  if( impltypes[i] == SCIP_BOUNDTYPE_LOWER && originalb == FALSE )
474  {
475  assert(implbounds[i] == 1.0);
476  return TRUE;
477  }
478  }
479  }
480 
481  /* check whether a and b are contained in a clique */
482  ncliques = (unsigned int) SCIPvarGetNCliques(vars[a], originala);
483  cliques = SCIPvarGetCliques(vars[a], originala);
484 
485  assert(cliques != NULL || ncliques == 0);
486 
487  for( i = 0; i < ncliques; ++i )
488  {
489  assert( cliques != NULL ); /* for lint */
490  ncliquevars = (unsigned int) SCIPcliqueGetNVars(cliques[i]);
491  cliquevars = SCIPcliqueGetVars(cliques[i]);
492  cliquevals = SCIPcliqueGetValues(cliques[i]);
493 
494  assert(cliquevars != NULL || ncliquevars == 0);
495  assert(cliquevals != NULL || ncliquevars == 0);
496 
497  for( j = 0; j < ncliquevars; ++j )
498  {
499  assert( cliquevals != NULL && cliquevars != NULL ); /* for lint */
500  if( SCIPvarGetProbindex(vars[b]) == SCIPvarGetProbindex(cliquevars[j]) )
501  {
502  if( (cliquevals[j] == FALSE && originalb == TRUE) || ( cliquevals[j] == TRUE && originalb == FALSE ) )
503  return TRUE;
504  }
505  }
506  }
507  }
508  }
509 
510  return FALSE;
511 }
512 
513 /** inside the lifting heuristic we determine the lifting coefficient by counting the length of
514  * chains adjacent to the lifting candidate.
515  *
516  * since we have to exclude all chains adjacent to an already lifted node which is not adjacent to
517  * the current lifting candidate we check all chains of the cycle of length three and block them if
518  * they are adjacent.
519  */
520 static
521 void checkBlocking(
522  unsigned int a, /**< first node of the checked cycle chain of length 3 */
523  unsigned int b, /**< second node of the checked cycle chain of length 3 */
524  unsigned int c, /**< third node of the checked cycle chain of length 3 */
525  unsigned int i, /**< current lifting candidate */
526  unsigned int* cycle, /**< list of cycle nodes in order of the cycle */
527  unsigned int ncyclevars, /**< number of variables in the odd cycle */
528  SCIP_VAR** vars, /**< problem variables */
529  unsigned int nbinvars, /**< number of binary problem variables */
530  unsigned int* lifted, /**< list of lifted nodes */
531  unsigned int* nlifted, /**< number of lifted nodes */
532  SCIP_SEPADATA* sepadata, /**< separator data structure */
533  SCIP_Bool* myi /**< flag array, if cycle node is inner point of a counted chain */
534  )
535 {
536  unsigned int k;
537 
538  assert(a < ncyclevars);
539  assert(b < ncyclevars);
540  assert(c < ncyclevars);
541  assert(cycle != NULL);
542  assert(ncyclevars % 2 == 1);
543  assert(ncyclevars > 2);
544  assert(ncyclevars <= nbinvars);
545  assert(vars != NULL);
546  assert(nbinvars > 2);
547  assert(lifted != NULL);
548  assert(nlifted != NULL);
549  assert(myi != NULL);
550 
551  k = 0;
552  while( (myi[a] || myi[b] || myi[c]) && k < *nlifted )
553  {
554  /* if all three nodes are adjacent to a node which is already lifted and not adjacent with the
555  * current lifting candidate, they cannot be regarded */
556  if( !isNeighbor(vars, nbinvars, sepadata, i, lifted[k])
557  && isNeighbor(vars, nbinvars, sepadata, cycle[a], lifted[k])
558  && isNeighbor(vars, nbinvars, sepadata, cycle[b], lifted[k])
559  && isNeighbor(vars, nbinvars, sepadata, cycle[c], lifted[k]) )
560  {
561  myi[a] = FALSE;
562  myi[b] = FALSE;
563  myi[c] = FALSE;
564  }
565  ++k;
566  }
567 }
568 
569 /** determine the heuristic lifting coefficient by counting the length of the adjacent chains of the
570  * candidate (we have to exclude all chains that are adjacent to an already lifted node which is
571  * not adjacent to the current candidate)
572  */
573 static
574 unsigned int getCoef(
575  SCIP* scip, /**< SCIP data structure */
576  unsigned int i, /**< current lifting candidate */
577  unsigned int* cycle, /**< list of cycle nodes in order of the cycle */
578  unsigned int ncyclevars, /**< number of variables in the odd cycle */
579  SCIP_VAR** vars, /**< problem variables */
580  unsigned int nbinvars, /**< number of binary problem variables */
581  unsigned int* lifted, /**< list of lifted nodes */
582  unsigned int* nlifted, /**< number of lifted nodes */
583  SCIP_SEPADATA* sepadata, /**< separator data structure */
584  SCIP_Bool* myi /**< flag array, if cycle node is inner point of a counted chain */
585  )
586 {
587  int j;
588  unsigned int k;
589  unsigned int coef; /* coefficient of lifting candidate of the current step */
590  unsigned int end;
591 
592  assert(scip != NULL);
593  assert(i < 2*nbinvars);
594  assert(cycle != NULL);
595  assert(ncyclevars % 2 == 1);
596  assert(ncyclevars > 2);
597  assert(ncyclevars <= 2*nbinvars);
598  assert(vars != NULL);
599  assert(nbinvars > 2);
600  assert(nlifted != NULL);
601  assert(lifted != NULL);
602 
603  coef = 0;
604 
605  /* get inner nodes of adjacent chains in cycle */
606  for( j = 1; j < (int)ncyclevars-1; ++j )
607  {
608  myi[j] = isNeighbor(vars, nbinvars, sepadata, i, cycle[j-1]) && isNeighbor(vars, nbinvars, sepadata, i, cycle[j])
609  && isNeighbor(vars, nbinvars, sepadata, i, cycle[j+1]);
610  }
611 
612  /* the first and last node of the cycle are treated separately */
613  myi[0] = isNeighbor(vars, nbinvars, sepadata, i, cycle[ncyclevars-1]) && isNeighbor(vars, nbinvars, sepadata, i, cycle[0])
614  && isNeighbor(vars, nbinvars, sepadata, i, cycle[1]);
615  myi[ncyclevars-1] = isNeighbor(vars, nbinvars, sepadata, i, cycle[ncyclevars-2])
616  && isNeighbor(vars, nbinvars, sepadata, i, cycle[ncyclevars-1])
617  && isNeighbor(vars, nbinvars, sepadata, i, cycle[0]);
618 
619  /* consider already lifted nodes that are not adjacent to current lifting candidate and
620  * remove all inner cycle nodes that are adjacent to them
621  */
622  for( j = 1; j < (int)ncyclevars-1; ++j )
623  {
624  checkBlocking((unsigned int) (j-1), (unsigned int) j, (unsigned int) (j+1), i, cycle, ncyclevars, vars, nbinvars, lifted, nlifted, sepadata, myi);
625  }
626  checkBlocking(ncyclevars-2, ncyclevars-1, 0, i, cycle, ncyclevars, vars, nbinvars, lifted, nlifted, sepadata, myi);
627  checkBlocking(ncyclevars-1, 0, 1, i, cycle, ncyclevars, vars, nbinvars, lifted, nlifted, sepadata, myi);
628 
629  /* calculate lifting coefficient */
630  k = 0;
631 
632  /* first, handle the special case, that the first node of the cycle list is part of a chain */
633  if( myi[0] )
634  {
635  ++k;
636  end = ncyclevars-1;
637  while( myi[end] && end > 0 )
638  {
639  ++k;
640  --end;
641  }
642  assert(k == ncyclevars || end > 0);
643 
644  /* all cycle nodes build a relevant chain (maximal chain s.t. all inner nodes are in myi) */
645  if( end == 0 )
646  {
647  assert(ncyclevars % 2 == 1);
648  coef = (ncyclevars-1)/2;
649  return coef;
650  }
651  assert(!myi[end]);
652 
653  /* current nonempty relevant chain cannot be extended */
654  if( !myi[1] )
655  {
656  coef = (unsigned int) SCIPfloor(scip,(k+1.0)/2.0);
657  assert(coef <= (ncyclevars-1)/2);
658  k = 0;
659  }
660  }
661  else
662  end = ncyclevars;
663 
664  /* find remaining relevant chains */
665  j = 1;
666  while( j < (int)end )
667  {
668  /* skip all nodes that are not inner node */
669  while( !myi[j] && j<(int)end )
670  ++j;
671 
672  /* collect all inner nodes (chain is extended) */
673  while( myi[j] && j<(int)end )
674  {
675  ++k;
676  ++j;
677  }
678 
679  if( k > 0 )
680  {
681  assert(myi[j-1]);
682  coef += (unsigned int) SCIPfloor(scip,(k+1.0)/2.0);
683  assert(coef <= (ncyclevars-1)/2);
684  k = 0;
685  }
686  }
687 
688  return coef;
689 }
690 
691 /** Lifting Heuristic based on an idea by Alvarez-Valdes, Parreno, Tamarit
692  *
693  * This method is based on the observation, that a non-cycle node can be lifted into the inequality
694  * with coefficient \f$1\f$ if the node is adjacent to the nodes of a 3-chain on the cycle.
695  *
696  * The coefficient can be calculated as
697  * \f$\left\lfloor{\frac{|C|-1}{2}}\right\rfloor\f$
698  * where \f$C\f$ is the chain on the cycle.
699  *
700  * If the node is connected to several chains, the coefficients of the chains can be summed up, resulting
701  * in a feasible lifting coefficient.
702  *
703  * Additionally further variables can be lifted by considering chains connected to the additional lifting node
704  * which are not connected to already lifted nodes.
705  *
706  * This method is a feasible heuristic which gives a valid lifted inequality.
707  * (Furthermore the first lifting coefficient is always smaller or equal to the largest possible lifting coefficient.)
708  */
709 static
710 SCIP_RETCODE liftOddCycleCut(
711  SCIP* scip, /**< SCIP data structure */
712  unsigned int* nlifted, /**< number of lifted variables */
713  unsigned int* lifted, /**< indices of the lifted variables */
714  unsigned int* liftcoef, /**< lifting coefficients */
715  SCIP_SEPADATA* sepadata, /**< separator data structure */
716  SCIP_VAR** vars, /**< problem variables */
717  unsigned int nbinvars, /**< number of binary problem variables */
718  unsigned int startnode, /**< a node of the cycle */
719  unsigned int* pred, /**< predecessor of each node (original and negated) in odd cycle */
720  unsigned int ncyclevars, /**< number of variables in the odd cycle */
721  SCIP_Real* vals, /**< values of the variables in the given solution */
722  SCIP_RESULT* result /**< pointer to store the result of the separation call */
723  )
724 {
725  unsigned int* cycle; /* storage for cycle and lifted nodes (and their coefficients) */
726  unsigned int* coef;
727  SCIP_Bool* candList; /* lifting candidate list */
728  unsigned int i;
729  unsigned int j;
730  unsigned int negated;
731  int bestcand;
732  unsigned int liftround;
733  SCIP_Bool* myi;
734 
735  assert(scip != NULL);
736  assert(sepadata != NULL);
737  assert(sepadata->levelgraph != NULL || sepadata->usegls);
738  assert(sepadata->dijkstragraph != NULL || ! sepadata->usegls);
739  assert(vars != NULL);
740  assert(nbinvars > 2);
741  assert(startnode < 2*nbinvars);
742  assert(pred != NULL);
743  assert(ncyclevars % 2 == 1);
744  assert(ncyclevars > 2);
745  assert(ncyclevars <= nbinvars);
746  assert(result != NULL);
747  assert(nlifted != NULL);
748  assert(lifted != NULL);
749  assert(liftcoef != NULL);
750 
751  /* allocate memory for cycle list */
752  SCIP_CALL( SCIPallocBufferArray(scip, &cycle, (int) ncyclevars) );
753 
754  /* transform cycle from predecessor list to array in order of appearance in cycle */
755  cycle[0] = startnode;
756  j = 1;
757  i = pred[startnode];
758  while( i != startnode )
759  {
760  cycle[j] = i;
761  i = pred[i];
762  ++j;
763  }
764  assert(j == ncyclevars);
765 
766  /* allocate memory for coefficients of the lifting candidates (used in every step) */
767  SCIP_CALL( SCIPallocBufferArray(scip, &coef, (int) (2*nbinvars)) );
768 
769  /* allocate memory candidate list and list of lifted nodes */
770  SCIP_CALL( SCIPallocBufferArray(scip, &candList, (int) (2*nbinvars)) );
771 
772  /* allocate memory for counting of chains in getCoef() */
773  SCIP_CALL( SCIPallocBufferArray(scip, &myi, (int) ncyclevars) );
774 
775  if( SCIPisStopped(scip) )
776  goto TERMINATE;
777 
778  /* initialize candidate list */
779  for( i = 0; i < 2*nbinvars; ++i )
780  candList[i] = TRUE;
781 
782  /* remove cycle variables and their negated from candidate list */
783  for( i = 0; i < ncyclevars; ++i )
784  {
785  candList[cycle[i]] = FALSE;
786  if( cycle[i] >= nbinvars )
787  negated = cycle[i] - nbinvars;
788  else
789  negated = cycle[i] + nbinvars;
790  assert(negated < 2*nbinvars);
791  candList[negated] = FALSE;
792  }
793 
794  /* no candidates lifted so far */
795  *nlifted = 0;
796  bestcand = 0;
797  liftround = 0;
798 
799  /* try lifting as long as we have lifting candidates */
800  while( bestcand >= 0 )
801  {
802  /* in case we use a lifting rule, which does not require the first lifting coefficient of all variables: REMOVE this */
803  if( sepadata->recalcliftcoef || liftround == 0 )
804  {
805  for( i = 0; i < 2*nbinvars; ++i )
806  {
807  if( candList[i] )
808  {
809  coef[i] = getCoef(scip, i, cycle, ncyclevars, vars, nbinvars, lifted, nlifted, sepadata, myi);
810  assert(coef[i] <= (ncyclevars-1)/2);
811  if( coef[i] < 1 )
812  candList[i] = FALSE;
813  }
814  }
815  }
816  ++liftround;
817  bestcand = -1;
818  for( i = 0; i < 2*nbinvars; ++i )
819  {
820  if( candList[i] )
821  {
822  /* we want to weight our choice of the lifting node by the value of the current lp solution */
823  if( sepadata->lpliftcoef )
824  {
825  if( bestcand < 0 || coef[i]*vals[i] > coef[bestcand]*vals[bestcand] )
826  bestcand = (int) i;
827  }
828  /* we only regard the coefficient */
829  else
830  {
831  if( bestcand < 0 || coef[i] > coef[bestcand] )
832  bestcand = (int) i;
833  }
834  }
835  }
836 
837  /* there is at least one lifting variable */
838  if( bestcand >= 0 )
839  {
840  if( !(sepadata->recalcliftcoef) )
841  coef[i] = getCoef(scip, (unsigned int) bestcand, cycle, ncyclevars, vars, nbinvars, lifted, nlifted, sepadata, myi);
842  assert(coef[bestcand] <= (ncyclevars-1)/2);
843  candList[bestcand] = FALSE;
844  if( coef[bestcand] > 0 )
845  {
846  if( bestcand >= (int)nbinvars )
847  negated = (unsigned int) bestcand - nbinvars;
848  else
849  negated = (unsigned int) bestcand + nbinvars;
850  assert(negated < 2*nbinvars);
851 
852  candList[negated] = FALSE;
853 
854  assert(*nlifted < nbinvars-ncyclevars);
855  lifted[*nlifted] = (unsigned int) bestcand;
856  liftcoef[*nlifted] = coef[bestcand];
857  ++(*nlifted);
858  }
859  }
860  }
861 
862  TERMINATE:
863  /* free temporary memory */
864  SCIPfreeBufferArray(scip, &myi);
865  SCIPfreeBufferArray(scip, &candList);
866  SCIPfreeBufferArray(scip, &coef);
867  SCIPfreeBufferArray(scip, &cycle);
868 
869  return SCIP_OKAY;
870 }
871 
872 
873 /*
874  * methods for both techniques
875  */
876 
877 /** add the inequality corresponding to the given odd cycle to the LP (if violated)
878  * after lifting it (if requested by user flag)
879  */
880 static
881 SCIP_RETCODE generateOddCycleCut(
882  SCIP* scip, /**< SCIP data structure */
883  SCIP_SEPA* sepa, /**< separator */
884  SCIP_SOL* sol, /**< given primal solution */
885  SCIP_VAR** vars, /**< problem variables */
886  unsigned int nbinvars, /**< number of binary problem variables */
887  unsigned int startnode, /**< a node of the cycle */
888  unsigned int* pred, /**< predecessor of each node (original and negated) in odd cycle */
889  unsigned int ncyclevars, /**< number of variables in the odd cycle */
890  unsigned int* incut, /**< TRUE iff node is covered already by a cut */
891  SCIP_Real* vals, /**< values of the variables in the given solution */
892  SCIP_SEPADATA* sepadata, /**< separator data structure */
893  SCIP_RESULT* result /**< pointer to store the result of the separation call */
894  )
895 {
896  unsigned int i;
897  unsigned int negatedcount;
898  unsigned int negated;
899 
900  /* memory for lifting */
901  unsigned int nlifted; /* number of lifted variables */
902  unsigned int* lifted; /* index of the lifted variables */
903  unsigned int* liftcoef; /* lifting coefficient */
904 
905  /* memory for cut generation */
906  SCIP_ROW* cut;
907  char cutname[SCIP_MAXSTRLEN];
908 
909  assert(scip != NULL);
910  assert(vars != NULL);
911  assert(startnode < 2*nbinvars);
912  assert(pred != NULL);
913  assert(ncyclevars % 2 == 1);
914  assert(ncyclevars <= nbinvars);
915  assert(incut != NULL);
916  assert(sepadata != NULL);
917  assert(sepadata->levelgraph != NULL || sepadata->usegls);
918  assert(sepadata->dijkstragraph != NULL || ! sepadata->usegls);
919  assert(result != NULL);
920 
921 #ifdef SCIP_OUTPUT
922  /* debug method that prints out all found cycles */
923  printCycle(vars, pred, nbinvars, startnode);
924 #endif
925 
926  /* cycle contains only one node */
927  if( ncyclevars < 3 )
928  {
929  SCIPdebugMessage("fixing variable\n");
930  /* strengthening variable bounds due to single-variable-cycle */
931  if( startnode < nbinvars )
932  {
933  SCIP_CALL( SCIPchgVarUb(scip, vars[startnode], 0.0) );
934  }
935  else
936  {
937  negated = startnode - nbinvars;
938  assert(negated < nbinvars);
939  SCIP_CALL( SCIPchgVarLb(scip, vars[negated], 1.0) );
940  }
941  *result = SCIP_REDUCEDDOM;
942  return SCIP_OKAY;
943  }
944 
945  /* cycle is a triangle (can be excluded by user) */
946  if( ncyclevars < 5 && ! sepadata->includetriangles )
947  return SCIP_OKAY;
948 
949  if( SCIPisStopped(scip) )
950  return SCIP_OKAY;
951 
952  /* lift the cycle inequality */
953  nlifted = 0;
954  lifted = NULL;
955  liftcoef = NULL;
956  if( sepadata->liftoddcycles )
957  {
958  SCIP_CALL( SCIPallocBufferArray(scip, &lifted, (int) (nbinvars - ncyclevars)) );
959  SCIP_CALL( SCIPallocBufferArray(scip, &liftcoef, (int) (nbinvars - ncyclevars)) );
960  SCIP_CALL( liftOddCycleCut(scip, &nlifted, lifted, liftcoef, sepadata, vars, nbinvars, startnode, pred, ncyclevars, vals, result) );
961  }
962  /* if we don't try to lift, we generate and add the cut as is */
963 
964  /* create cut */
965  (void) SCIPsnprintf(cutname, SCIP_MAXSTRLEN, "oddcycle_%d", sepadata->ncuts);
966  SCIP_CALL( SCIPcreateEmptyRowSepa(scip, &cut, sepa, cutname, -SCIPinfinity(scip), (ncyclevars-1)/2.0, FALSE, FALSE, TRUE) );
967  SCIP_CALL( SCIPcacheRowExtensions(scip, cut) );
968  negatedcount = 0;
969 
970  /* add variables of odd cycle to cut inequality */
971  i = pred[startnode];
972  while( i != startnode )
973  {
974  assert(i < 2*nbinvars);
975  if( i < nbinvars )
976  {
977  /* inserting original variable */
978  SCIP_CALL( SCIPaddVarToRow(scip, cut, vars[i], 1.0) );
979  incut[i] = TRUE;
980  }
981  else
982  {
983  negated = i - nbinvars;
984  assert(negated < nbinvars);
985 
986  /* inserting negated variable */
987  SCIP_CALL( SCIPaddVarToRow(scip, cut, vars[negated], -1.0) );
988  ++negatedcount;
989  incut[negated] = TRUE;
990  }
991  i = pred[i];
992  }
993 
994  /* insert startnode */
995  if( startnode < nbinvars )
996  {
997  /* inserting original variable */
998  SCIP_CALL( SCIPaddVarToRow(scip, cut, vars[startnode], 1.0) );
999  incut[i] = TRUE;
1000  }
1001  else
1002  {
1003  negated = startnode - nbinvars;
1004  assert(negated < nbinvars);
1005 
1006  /* inserting negated variable */
1007  SCIP_CALL( SCIPaddVarToRow(scip, cut, vars[negated], -1.0) );
1008  ++negatedcount;
1009  incut[negated] = TRUE;
1010  }
1011 
1012  /* add lifted variables to cut inequality (if existing) */
1013  for( i = 0; i < nlifted; ++i)
1014  {
1015  assert(lifted != NULL);
1016  assert(liftcoef != NULL);
1017  if( lifted[i] < nbinvars )
1018  {
1019  assert(vars[lifted[i]] != NULL);
1020  SCIP_CALL( SCIPaddVarToRow(scip, cut, vars[lifted[i]], (SCIP_Real) liftcoef[i]) );
1021  }
1022  else
1023  {
1024  negated = lifted[i] - nbinvars;
1025  assert(negated < nbinvars);
1026  assert(vars[negated] != NULL);
1027  SCIP_CALL( SCIPaddVarToRow(scip, cut, vars[negated], -1.0*liftcoef[i]) );
1028  negatedcount += liftcoef[i];
1029  }
1030  }
1031 
1032  /* modify right hand side corresponding to number of added negated variables */
1033  SCIP_CALL( SCIPchgRowRhs(scip, cut, SCIProwGetRhs(cut) - negatedcount) );
1034  SCIP_CALL( SCIPflushRowExtensions(scip, cut) );
1035 
1036  /* set cut rank: for oddcycle cuts we always set to 1 */
1037  SCIProwChgRank(cut, 1);
1038 
1039  /* not every odd cycle has to be violated due to incompleteness of the implication graph */
1040  if( SCIPisCutEfficacious(scip, sol, cut) )
1041  {
1042  SCIP_Bool infeasible;
1043 
1044  SCIP_CALL( SCIPaddCut(scip, sol, cut, FALSE, &infeasible) );
1045  ++sepadata->ncuts;
1046  if ( nlifted > 0 )
1047  ++sepadata->nliftedcuts;
1048  if ( infeasible )
1049  *result = SCIP_CUTOFF;
1050  else
1051  {
1052  SCIP_CALL( SCIPaddPoolCut(scip, cut) );
1053 
1054  if( *result == SCIP_DIDNOTFIND )
1055  *result = SCIP_SEPARATED;
1056  }
1057 
1058 #ifdef SCIP_OUTPUT
1059  SCIP_CALL( SCIPprintRow(scip, cut, NULL) );
1060 #endif
1061 
1062  assert(*result == SCIP_CUTOFF || *result == SCIP_SEPARATED || *result == SCIP_REDUCEDDOM);
1063  }
1064 
1065  SCIP_CALL( SCIPreleaseRow(scip, &cut) );
1066 
1067  if( sepadata->liftoddcycles )
1068  {
1069  SCIPfreeBufferArray(scip, &liftcoef);
1070  SCIPfreeBufferArray(scip, &lifted);
1071  }
1072  return SCIP_OKAY;
1073 }
1074 
1075 /** check whether the given object is really a cycle without sub-cycles (sub-cycles may be
1076  * calculated by the GLS algorithm in case there is no violated odd cycle inequality) and removes
1077  * pairs of original and negated variables from the cycle
1078  */
1079 static
1080 SCIP_RETCODE cleanCycle(
1081  SCIP* scip, /**< SCIP data structure */
1082  unsigned int* pred, /**< predecessor list of current cycle segment */
1083  SCIP_Bool* incycle, /**< flag array iff node is in cycle segment */
1084  unsigned int* incut, /**< flag array iff node is already covered by a cut */
1085  unsigned int x, /**< index of current variable */
1086  unsigned int startnode, /**< index of first variable of cycle segment */
1087  unsigned int nbinvars, /**< number of binary problem variables */
1088  unsigned int* ncyclevars, /**< number of nodes in current cycle segment */
1089  SCIP_Bool repaircycles, /**< user flag if we should try to repair damaged cycles */
1090  SCIP_Bool allowmultiplecuts, /**< user flag if we allow multiple cuts per node */
1091  SCIP_Bool* success /**< FALSE iff an irreparable cycle appears */
1092  )
1093 {
1094  unsigned int negx;
1095 
1096  assert(scip != NULL);
1097  assert(pred != NULL);
1098  assert(incycle != NULL);
1099  assert(incut != NULL);
1100  assert(ncyclevars != NULL);
1101  assert(*ncyclevars <= nbinvars);
1102  assert(success != NULL);
1103  assert(*success);
1104 
1105  assert(x < 2*nbinvars);
1106 
1107  /* skip variable if it is already covered by a cut and we do not allow multiple cuts per node */
1108  if( incut[x] && !allowmultiplecuts )
1109  {
1110  *success = FALSE;
1111  return SCIP_OKAY;
1112  }
1113 
1114  /* get index of negated variable of current variable */
1115  if( x < nbinvars )
1116  negx = x + nbinvars;
1117  else
1118  negx = x - nbinvars;
1119  assert(negx < 2*nbinvars);
1120 
1121  /* given object is not an odd cycle (contains sub-cycle) or contains original and negated
1122  * variable pair but we should not repair this
1123  */
1124  if( incycle[x] || (incycle[negx] && !repaircycles) )
1125  {
1126  *success = FALSE;
1127  return SCIP_OKAY;
1128  }
1129 
1130  /* cycle does not contain original and negated variable pair */
1131  if( !incycle[negx] )
1132  {
1133  assert(!incycle[x]);
1134  incycle[x] = TRUE;
1135  ++(*ncyclevars);
1136  return SCIP_OKAY;
1137  }
1138 
1139  /* delete original and negated variable and cross-link their neighbors the following way, if possible:
1140  * suppose the cycle contains segments:
1141  * startnode - ... - a - neg(x) - c1 - c2 - ... - cn-1 - cn - x - z=pred(x)
1142  *
1143  * because of the chain a - neg(x) - x - cn it holds that
1144  * a=1 => x=0 => neg(x)=1 => cn=0 and
1145  * cn=1 => x=0 => neg(x)=1 => a=0
1146  * because of the chain z - x - neg(x) - b it holds that
1147  * z=1 => x=0 => neg(x)=1 => c1=0 and
1148  * c1=1 => x=0 => neg(x)=1 => z=0
1149  *
1150  * in addition to that, in our linked list structure we need to relink the chain c1-...-cn in reverse order.
1151  * so we gain the order: a - cn - cn-1 - ... - c2 - c1 - z
1152  */
1153 
1154  /* if negated variable is first node in cycle,
1155  * cross-linking not possible because there is no successor z of neg(x) contained in cycle yet
1156  */
1157  if( negx == startnode )
1158  {
1159  *success = FALSE;
1160  return SCIP_OKAY;
1161  }
1162 
1163  /* if original and negated variable are neighbors, cross linking is not possible,
1164  * but x and neg(x) can simply be removed
1165  * a - neg(x)=pred[a] - x=pred[neg(x)] - z=pred[x] --> a - z=pred[x]=:pred[a]
1166  */
1167  if( pred[negx] == x )
1168  {
1169  unsigned int a;
1170 
1171  /* find a */
1172  a = startnode;
1173  while( pred[a] != negx )
1174  a = pred[a];
1175 
1176  /* link a and z */
1177  pred[a] = pred[x];
1178  }
1179  /* cross linking as mentioned above */
1180  else
1181  {
1182  unsigned int a;
1183  unsigned int z;
1184 
1185  /* memory for chain reverse */
1186  unsigned int* chain;
1187  unsigned int nchain;
1188 
1189  unsigned int i;
1190 
1191  /* allocate temporary memory */
1192  SCIP_CALL( SCIPallocBufferArray(scip, &chain, (int) *ncyclevars) );
1193 
1194  /* find and store a */
1195  a = startnode;
1196  while( pred[a] != negx )
1197  a = pred[a];
1198 
1199  /* store chain */
1200  i = pred[negx];
1201  nchain = 0;
1202  while( i != x )
1203  {
1204  chain[nchain] = i;
1205  ++nchain;
1206  i = pred[i];
1207  }
1208  assert(nchain > 0);
1209 
1210  /* store z */
1211  z = pred[x];
1212 
1213  /* link a and c1 */
1214  pred[a] = chain[nchain-1];
1215 
1216  /* link cn and z */
1217  pred[chain[0]] = z;
1218 
1219  /* reverse the chain */
1220  for( i = nchain-1; i > 0; --i )
1221  pred[chain[i]] = chain[i-1];
1222 
1223  /* free temporary memory */
1224  SCIPfreeBufferArray(scip, &chain);
1225  }
1226 
1227  /* remove negated variable from cycle */
1228  assert(!incycle[x] && incycle[negx]);
1229  incycle[negx] = FALSE;
1230  --(*ncyclevars);
1231 
1232  return SCIP_OKAY;
1233 }
1234 
1235 
1236 /*
1237  * methods for separateHeur()
1238  */
1239 
1240 /** memory reallocation method (the graph is normally very dense, so we dynamically allocate only the memory we need)
1241  *
1242  * Since the array sizes differ the method can be called for each of the three data structure types:
1243  * - Forward: sizeForward, targetForward, weightForward
1244  * - Backward: sizeBackward, targetBackward, weightBackward
1245  * - Adj (inner level edges): sizeAdj, sourceAdj, targetAdj, weightAdj
1246  */
1247 static
1248 SCIP_RETCODE checkArraySizesHeur(
1249  SCIP* scip, /**< SCIP data structure */
1250  LEVELGRAPH* graph, /**< LEVELGRAPH data structure */
1251  unsigned int* size, /**< given size */
1252  int** targetArray, /**< given target array (or NULL if sourceAdjArray and targetAdjArray given) */
1253  unsigned int** weightArray, /**< given weight array */
1254  unsigned int** sourceAdjArray, /**< given sourceAdj array (or NULL if targetArray given) */
1255  unsigned int** targetAdjArray, /**< given targetAdj array (or NULL if targetArray given) */
1256  SCIP_Bool* success /**< FALSE, iff memory reallocation fails */
1257  )
1258 {
1259  SCIP_Real memorylimit;
1260  unsigned int additional;
1261 
1262  assert(scip != NULL);
1263  assert(graph != NULL);
1264  assert(size != NULL);
1265  assert(targetArray != NULL || (sourceAdjArray != NULL && targetAdjArray != NULL));
1266  assert(weightArray != NULL);
1267  assert(success != NULL);
1268 
1269  SCIPdebugMessage("reallocating...\n");
1270 
1271  additional = MIN(graph->maxarcs + graph->maxnodes - *size, *size) * ((int) sizeof(**weightArray));
1272  if( targetArray != NULL )
1273  {
1274  additional += MIN(graph->maxarcs + graph->maxnodes - *size, *size) * ((int) sizeof(**targetArray));
1275  }
1276  else
1277  {
1278  additional += MIN(graph->maxarcs + graph->maxnodes - *size, *size) * ((int) sizeof(**sourceAdjArray));
1279  additional += MIN(graph->maxarcs + graph->maxnodes - *size, *size) * ((int) sizeof(**targetAdjArray));
1280  }
1281 
1282  SCIP_CALL( SCIPgetRealParam(scip, "limits/memory", &memorylimit) );
1283  if( !SCIPisInfinity(scip, memorylimit) )
1284  {
1285  memorylimit -= SCIPgetMemUsed(scip)/1048576.0;
1286  memorylimit -= SCIPgetMemExternEstim(scip)/1048576.0;
1287  }
1288 
1289 
1290  /* if memorylimit would be exceeded or any other limit is reached free all data and exit */
1291  if( memorylimit <= additional/1048576.0 || SCIPisStopped(scip) )
1292  {
1293  *success = FALSE;
1294  SCIPdebugMessage("...memory limit exceeded\n");
1295  return SCIP_OKAY;
1296  }
1297 
1298  *size = 2 * (*size);
1299 
1300  SCIP_CALL( SCIPreallocBufferArray(scip, weightArray, (int) MIN(graph->maxarcs + graph->maxnodes, *size)) );
1301  if( targetArray != NULL )
1302  {
1303  SCIP_CALL( SCIPreallocBufferArray(scip, targetArray, (int) MIN(graph->maxarcs + graph->maxnodes, *size)) );
1304  }
1305  else
1306  {
1307  SCIP_CALL( SCIPreallocBufferArray(scip, sourceAdjArray, (int) MIN(graph->maxarcs, *size)) );
1308  SCIP_CALL( SCIPreallocBufferArray(scip, targetAdjArray, (int) MIN(graph->maxarcs, *size)) );
1309  }
1310 
1311  /* if memorylimit is exceeded free all data and exit */
1312  SCIP_CALL( SCIPgetRealParam(scip, "limits/memory", &memorylimit) );
1313  if( !SCIPisInfinity(scip, memorylimit) )
1314  {
1315  memorylimit -= SCIPgetMemUsed(scip)/1048576.0;
1316  memorylimit -= SCIPgetMemExternEstim(scip)/1048576.0;
1317  }
1318 
1319  if( memorylimit <= 2.0*SCIPgetMemExternEstim(scip)/1048576.0 )
1320  {
1321  *success = FALSE;
1322  SCIPdebugMessage("...memory limit exceeded\n");
1323  return SCIP_OKAY;
1324  }
1325 
1326  SCIPdebugMessage("...with success\n");
1327 
1328  return SCIP_OKAY;
1329 }
1330 
1331 /** Add arc to level graph */
1332 static
1333 SCIP_RETCODE addArc(
1334  SCIP* scip, /**< SCIP data structure */
1335  LEVELGRAPH* graph, /**< LEVELGRAPH data structure */
1336  unsigned int u, /**< source node */
1337  unsigned int v, /**< target node */
1338  unsigned int level, /**< number of current level */
1339  unsigned int weight, /**< weight of the arc */
1340  unsigned int* nAdj, /**< array of numbers of arcs inside levels */
1341  SCIP_Bool* success /**< FALSE, iff memory reallocation fails */
1342  )
1343 {
1344  /* arc is a forward arc */
1345  if( graph->level[v] == level+1 )
1346  {
1347  graph->targetForward[graph->lastF] = (int) v;
1348  graph->weightForward[graph->lastF] = weight;
1349  ++(graph->lastF);
1350  ++(graph->narcs);
1351  if( graph->lastF == graph->sizeForward )
1352  {
1353  SCIP_CALL( checkArraySizesHeur(scip, graph, &(graph->sizeForward), &(graph->targetForward),
1354  &(graph->weightForward), NULL, NULL, success) );
1355  if( !(*success) )
1356  return SCIP_OKAY;
1357  }
1358  }
1359  else
1360  {
1361  assert(graph->level[v] == level || graph->level[v] == level-1);
1362 
1363  /* arc is a backward arc */
1364  if( graph->level[v] == level-1 )
1365  {
1366  graph->targetBackward[graph->lastB] = (int) v;
1367  graph->weightBackward[graph->lastB] = weight;
1368  ++(graph->lastB);
1369  ++(graph->narcs);
1370 
1371  if( graph->lastB == graph->sizeBackward )
1372  {
1373  SCIP_CALL( checkArraySizesHeur(scip, graph, &(graph->sizeBackward), &(graph->targetBackward),
1374  &(graph->weightBackward), NULL, NULL, success) );
1375  if( !(*success) )
1376  return SCIP_OKAY;
1377  }
1378  }
1379  else /* arc is in the same level */
1380  {
1381  assert(graph->level[v] == level);
1382 
1383  /* add arc only once, i.e., if u < v */
1384  if( u < v )
1385  {
1386  graph->sourceAdj[graph->levelAdj[level+1]+*nAdj] = u;
1387  graph->targetAdj[graph->levelAdj[level+1]+*nAdj] = v;
1388  graph->weightAdj[graph->levelAdj[level+1]+*nAdj] = weight;
1389  ++(*nAdj);
1390  ++(graph->narcs);
1391 
1392  if( graph->levelAdj[level+1]+*nAdj == graph->sizeAdj )
1393  {
1394  SCIP_CALL( checkArraySizesHeur(scip, graph, &(graph->sizeAdj), NULL, &(graph->weightAdj),
1395  &(graph->sourceAdj), &(graph->targetAdj), success) );
1396  if( !(*success) )
1397  return SCIP_OKAY;
1398  }
1399  }
1400  }
1401  }
1402  return SCIP_OKAY;
1403 }
1404 
1405 /** add binary implications of the given node u @see createNextLevel() */
1406 static
1407 SCIP_RETCODE addNextLevelBinImpls(
1408  SCIP* scip, /**< SCIP data structure */
1409  SCIP_SEPADATA* sepadata, /**< separator data structure */
1410  SCIP_VAR** vars, /**< problem variables */
1411  SCIP_Real* vals, /**< values of the binary variables in the current LP relaxation */
1412  unsigned int u, /**< current node */
1413  LEVELGRAPH* graph, /**< LEVELGRAPH data structure */
1414  unsigned int level, /**< number of current level */
1415  SCIP_Bool* inlevelgraph, /**< flag array if node is already inserted in level graph */
1416  unsigned int* newlevel, /**< array of nodes of the next level */
1417  unsigned int* nnewlevel, /**< number of nodes of the next level */
1418  unsigned int* nAdj, /**< array of numbers of arcs inside levels */
1419  SCIP_Bool* success /**< FALSE, iff memory reallocation fails */
1420  )
1421 {
1422  SCIP_Bool varfixing;
1423  unsigned int nbinimpls;
1424  unsigned int nbinvars;
1425  unsigned int varsidx;
1426  SCIP_VAR** implvars;
1427  SCIP_BOUNDTYPE* impltypes;
1428  unsigned int j;
1429 
1430  assert(scip != NULL);
1431  assert(vars != NULL);
1432  assert(vals != NULL);
1433  assert(graph != NULL);
1434  assert(graph->targetForward != NULL);
1435  assert(graph->weightForward != NULL);
1436  assert(graph->targetBackward != NULL);
1437  assert(graph->weightBackward != NULL);
1438  assert(graph->sourceAdj != NULL);
1439  assert(graph->targetAdj != NULL);
1440  assert(graph->weightAdj != NULL);
1441  assert(inlevelgraph != NULL);
1442  assert(newlevel != NULL);
1443  assert(nnewlevel != NULL);
1444  assert(nAdj != NULL);
1445  assert(success != NULL);
1446 
1447  assert(u < graph->maxnodes);
1448 
1449  nbinvars = (graph->maxnodes)/2;
1450 
1451  /* current node signifies a problem variable */
1452  if( u < nbinvars )
1453  {
1454  varfixing = TRUE;
1455  varsidx = u;
1456  }
1457  /* current node signifies a negated variable */
1458  else
1459  {
1460  varfixing = FALSE;
1461  varsidx = u - nbinvars;
1462  }
1463  assert(varsidx < nbinvars);
1464  assert(! SCIPisFeasIntegral(scip, vals[varsidx]) );
1465 
1466  /* get binary implications of the current variable */
1467  nbinimpls = (unsigned int) SCIPvarGetNBinImpls(vars[varsidx], varfixing);
1468  implvars = SCIPvarGetImplVars(vars[varsidx], varfixing);
1469  impltypes = SCIPvarGetImplTypes(vars[varsidx], varfixing);
1470 
1471  assert(implvars != NULL || nbinimpls == 0);
1472  assert(impltypes != NULL || nbinimpls == 0);
1473 
1474  for( j = 0; j < nbinimpls; ++j )
1475  {
1476  unsigned int k;
1477  unsigned int v;
1478  unsigned int weight;
1479 
1480  assert( implvars != NULL && impltypes != NULL );
1481  assert( SCIPvarIsBinary(implvars[j]) );
1482 
1483  k = sepadata->mapping[SCIPvarGetProbindex(implvars[j])];
1484  assert(k < nbinvars);
1485 
1486  /* skip integral neighbors */
1487  if( SCIPisFeasIntegral(scip, vals[k]) )
1488  continue;
1489 
1490  /* consider implication to negated variable (x = 1 -> y >= 1 <=> x = 1 -> neg(y) <= 0) */
1491  if( impltypes[j] == SCIP_BOUNDTYPE_LOWER )
1492  v = k + nbinvars;
1493  /* x = 1 -> y <= 0 */
1494  else
1495  {
1496  assert(impltypes[j] == SCIP_BOUNDTYPE_UPPER);
1497  v = k;
1498  }
1499  assert(v < graph->maxnodes);
1500 
1501  /* if variable is a new node, it will be assigned to the next level, but if the level contains
1502  * more nodes than allowed (defined by percent per level plus offset), we skip the rest of the
1503  * nodes
1504  */
1505  if( !inlevelgraph[v] && (*nnewlevel) <= sepadata->maxlevelsize )
1506  {
1507  ++(graph->nnodes);
1508  graph->level[v] = level+1;
1509  inlevelgraph[v] = TRUE;
1510  newlevel[*nnewlevel] = v;
1511  ++(*nnewlevel);
1512  }
1513  assert((*nnewlevel) > sepadata->maxlevelsize || inlevelgraph[v]);
1514 
1515  /* calculate arc weight and add arc, if the neighbor node is on the same or a neighbor level */
1516  if( inlevelgraph[v] && (graph->level[v] == level+1 || graph->level[v] == level || graph->level[v] == level-1))
1517  {
1518  int tmp;
1519 
1520  /* the computation of 1.0 - vals[v] if v is negated is ensured by the fact that v > nbinvars in this case */
1521  /* set weight of arc (x,y) to 1 - x* -y* */
1522  if( varfixing )
1523  {
1524  /* x = 1 -> y <= 0 or y >= 1, i.e., x + y <= 1 or x + (1-y) <= 1 */
1525  tmp = (int) SCIPfeasCeil(scip, sepadata->scale * (1.0 - vals[varsidx] - vals[v]));
1526  weight = (unsigned int) MAX(tmp, sepadata->maxreference);
1527  }
1528  else
1529  {
1530  /* x = 0 <-> neg(x) = 1 -> y <= 0 or y >= 1, i.e., (1-x) + y <= 1 or (1-x) + (1-y) <= 1 */
1531  tmp = (int) SCIPfeasCeil(scip, sepadata->scale * (1.0 - (1.0 - vals[varsidx]) - vals[v]));
1532  weight = (unsigned int) MAX(tmp, sepadata->maxreference);
1533  }
1534 
1535  /* add arc from current to neighbor node */
1536  SCIP_CALL( addArc(scip, graph, u, v, level, weight, nAdj, success) );
1537  if( !(*success) )
1538  return SCIP_OKAY;
1539  }
1540  }
1541 
1542  return SCIP_OKAY;
1543 }
1544 
1545 /** add implications from cliques of the given node u
1546  *
1547  * @see createNextLevel()
1548  */
1549 static
1550 SCIP_RETCODE addNextLevelCliques(
1551  SCIP* scip, /**< SCIP data structure */
1552  SCIP_SEPADATA* sepadata, /**< separator data structure */
1553  SCIP_VAR** vars, /**< problem variables */
1554  SCIP_Real* vals, /**< values of the binary variables in the current LP relaxation */
1555  unsigned int u, /**< current node */
1556  LEVELGRAPH* graph, /**< LEVELGRAPH data structure */
1557  unsigned int level, /**< number of current level */
1558  SCIP_Bool* inlevelgraph, /**< flag array if node is already inserted in level graph */
1559  unsigned int* newlevel, /**< array of nodes of the next level */
1560  unsigned int* nnewlevel, /**< number of nodes of the next level */
1561  unsigned int* nAdj, /**< array of numbers of arcs inside levels */
1562  SCIP_Bool* success /**< FALSE, iff memory reallocation fails */
1563  )
1564 {
1565  SCIP_Bool varfixing;
1566  unsigned int ncliques;
1567  unsigned int nbinvars;
1568  unsigned int varsidx;
1569  SCIP_CLIQUE** cliques;
1570  unsigned int ncliquevars;
1571  SCIP_VAR** cliquevars;
1572  SCIP_Bool* cliquevals;
1573  unsigned int j;
1574  unsigned int k;
1575 
1576  assert(scip != NULL);
1577  assert(vars != NULL);
1578  assert(vals != NULL);
1579  assert(graph != NULL);
1580  assert(graph->targetForward != NULL);
1581  assert(graph->weightForward != NULL);
1582  assert(graph->targetBackward != NULL);
1583  assert(graph->weightBackward != NULL);
1584  assert(graph->sourceAdj != NULL);
1585  assert(graph->targetAdj != NULL);
1586  assert(graph->weightAdj != NULL);
1587  assert(inlevelgraph != NULL);
1588  assert(newlevel != NULL);
1589  assert(nnewlevel != NULL);
1590  assert(nAdj != NULL);
1591  assert(success != NULL);
1592 
1593  assert(u < graph->maxnodes);
1594 
1595  nbinvars = (graph->maxnodes)/2;
1596 
1597  /* current node signifies a problem variable */
1598  if( u < nbinvars )
1599  {
1600  varfixing = TRUE;
1601  varsidx = u;
1602  }
1603  /* current node signifies a negated variable */
1604  else
1605  {
1606  varfixing = FALSE;
1607  varsidx = u - nbinvars;
1608  }
1609  assert(varsidx < nbinvars);
1610  assert(!SCIPisFeasIntegral(scip, vals[varsidx]));
1611 
1612  /* get cliques of the current variable */
1613  ncliques = (unsigned int) SCIPvarGetNCliques(vars[varsidx], varfixing);
1614  if( ncliques == 0 )
1615  return SCIP_OKAY;
1616 
1617  cliques = SCIPvarGetCliques(vars[varsidx], varfixing);
1618  assert(cliques != NULL);
1619 
1620  for( j = 0; j < ncliques; ++j )
1621  {
1622  ncliquevars = (unsigned int) SCIPcliqueGetNVars(cliques[j]);
1623  cliquevars = SCIPcliqueGetVars(cliques[j]);
1624  cliquevals = SCIPcliqueGetValues(cliques[j]);
1625 
1626  assert(cliquevars != NULL || ncliquevars == 0);
1627  assert(cliquevals != NULL || ncliquevars == 0);
1628 
1629  for( k = 0; k < ncliquevars; ++k )
1630  {
1631  unsigned int l;
1632  unsigned int v;
1633  unsigned int weight;
1634 
1635  assert( cliquevars != NULL && cliquevals != NULL ); /* for lint */
1636 
1637  l = sepadata->mapping[SCIPvarGetProbindex(cliquevars[k])];
1638  assert(l < nbinvars);
1639 
1640  /* skip integral neighbors */
1641  if( SCIPisFeasIntegral(scip, vals[l]) )
1642  continue;
1643 
1644  /* consider clique with negated variable (x = 1 -> y >= 1 <=> x = 1 -> neg(y) <= 0) */
1645  if( cliquevals[k] == FALSE )
1646  v = l + nbinvars;
1647  /* x = 1 -> y <= 0 */
1648  else
1649  v = l;
1650  assert(v < graph->maxnodes);
1651 
1652  /* if variable is a new node, it will be assigned to the next level,
1653  * but if the level contains more nodes than allowed
1654  * (defined by percent per level plus offset),
1655  * we skip the rest of the nodes
1656  */
1657  if( !inlevelgraph[v] && (*nnewlevel) <= sepadata->maxlevelsize )
1658  {
1659  ++(graph->nnodes);
1660  graph->level[v] = level+1;
1661  inlevelgraph[v] = TRUE;
1662  newlevel[*nnewlevel] = v;
1663  ++(*nnewlevel);
1664  }
1665  assert(*nnewlevel > sepadata->maxlevelsize || inlevelgraph[v]);
1666 
1667  /* calculate arc weight and add arc, if the neighbor node is on the same or a neighbor level */
1668  if( inlevelgraph[v] && (graph->level[v] == level+1 || graph->level[v] == level || graph->level[v] == level-1))
1669  {
1670  int tmp;
1671 
1672  /* the computation of 1.0 - vals[v] if v is negated is ensured by the fact that v > nbinvars in this case */
1673  /* set weight of arc (x,y) to 1 - x* -y* */
1674  if( varfixing )
1675  {
1676  /* x = 1 -> y <= 0 or y >= 1 */
1677  tmp = (int) SCIPfeasCeil(scip, sepadata->scale * (1 - vals[varsidx] - vals[v]));
1678  weight = (unsigned int) MAX(tmp, sepadata->maxreference);
1679  }
1680  else
1681  {
1682  /* x = 0 <-> neg(x) = 1 -> y <= 0 or y >= 1 */
1683  tmp = (int) SCIPfeasCeil(scip, sepadata->scale * (1.0 - (1.0 - vals[varsidx]) - vals[v]));
1684  weight = (unsigned int) MAX(tmp, sepadata->maxreference);
1685  }
1686 
1687  /* add arc from current to neighbor node */
1688  SCIP_CALL( addArc(scip, graph, u, v, level, weight, nAdj, success) );
1689  if( !(*success) )
1690  return SCIP_OKAY;
1691  }
1692  }
1693  }
1694  return SCIP_OKAY;
1695 }
1696 
1697 
1698 /** sort level of root neighbors
1699  *
1700  * If we limit the size of nodes of a level, we want to add the best neighbors to the next level.
1701  * Since sorting every level is too expensive, we sort the neighbors of the root (if requested).
1702  *
1703  * Create the first level as follows:
1704  * - create flag array for binary variables and their negated and set their values FALSE
1705  * - iterate over the implication and clique neighbors of the root and set their flag array values to TRUE
1706  * - create variable array and insert all variables with flag value TRUE
1707  * - sort variable array by maximal fractionality
1708  * - add variables from sorted array to levelgraph until first level is full (or all variables are inserted)
1709  *
1710  * Even inserting all variables might help for the following creation of further levels since the neighbors
1711  * of nodes with high fractionality often have high fractionalities themselves and would be inserted first
1712  * when further levels would have been sorted (which actually is not the case).
1713  */
1714 static
1715 SCIP_RETCODE insertSortedRootNeighbors(
1716  SCIP* scip, /**< SCIP data structure */
1717  LEVELGRAPH* graph, /**< LEVELGRAPH data structure */
1718  unsigned int nbinvars, /**< number of binary problem variables */
1719  unsigned int ncurlevel, /**< number of nodes of the current level */
1720  unsigned int u, /**< source node */
1721  SCIP_Real* vals, /**< values of the binary variables in the current LP relaxation */
1722  SCIP_VAR** vars, /**< problem variables */
1723  SCIP_SEPADATA* sepadata, /**< separator data structure */
1724  unsigned int* nnewlevel, /**< number of nodes of the next level */
1725  SCIP_Bool* inlevelgraph, /**< nodes in new graph corr. to old graph (-1 if unassigned) */
1726  unsigned int level, /**< number of current level */
1727  unsigned int* newlevel, /**< array of nodes of the next level */
1728  SCIP_Bool* success /**< FALSE, iff memory reallocation fails */
1729  )
1730 {
1731  /* storage for the neighbors of the root */
1732  unsigned int root;
1733  unsigned int nneighbors;
1734  SCIP_Bool* isneighbor;
1735  int* neighbors;
1736  SCIP_Real* sortvals;
1737 
1738  SCIP_Bool varfixing;
1739  unsigned int varsidx;
1740 
1741  /* storage for implications to the neighbors of the root node */
1742  unsigned int nbinimpls;
1743  SCIP_VAR** implvars;
1744  SCIP_BOUNDTYPE* impltypes;
1745 
1746  /* storage for cliques to the neighbors of the root node */
1747  unsigned int ncliques;
1748  SCIP_CLIQUE** cliques;
1749  unsigned int ncliquevars;
1750  SCIP_VAR** cliquevars;
1751  SCIP_Bool* cliquevals;
1752 
1753  unsigned int j;
1754  unsigned int k;
1755  unsigned int v;
1756 
1757  /* allocate flag array for neighbor detection */
1758  SCIP_CALL( SCIPallocBufferArray(scip, &isneighbor, (int) graph->maxnodes) );
1759  BMSclearMemoryArray(isneighbor, graph->maxnodes);
1760 
1761  nbinvars = (graph->maxnodes)/2;
1762 
1763  assert(ncurlevel == 1);
1764  root = u;
1765 
1766  /* current node signifies a problem variable */
1767  if( root < nbinvars )
1768  {
1769  varfixing = TRUE;
1770  varsidx = root;
1771  }
1772  /* current node signifies a negated variable */
1773  else
1774  {
1775  varfixing = FALSE;
1776  varsidx = root - nbinvars;
1777  }
1778  assert(varsidx < nbinvars);
1779  assert(! SCIPisFeasIntegral(scip, vals[varsidx]));
1780  nneighbors = 0;
1781 
1782  /* count implications of the root */
1783  nbinimpls = (unsigned int) SCIPvarGetNBinImpls(vars[varsidx], varfixing);
1784  if( nbinimpls > 0 )
1785  {
1786  unsigned int jidx;
1787 
1788  implvars = SCIPvarGetImplVars(vars[varsidx], varfixing);
1789  impltypes = SCIPvarGetImplTypes(vars[varsidx], varfixing);
1790  assert(implvars != NULL);
1791  assert(impltypes != NULL);
1792 
1793  for( j = 0; j < nbinimpls ; ++j )
1794  {
1795  jidx = sepadata->mapping[SCIPvarGetProbindex(implvars[j])];
1796  assert(jidx < nbinvars);
1797 
1798  if( SCIPisFeasIntegral(scip, vals[jidx]))
1799  continue;
1800 
1801  /* implication x + y <= 1 or (1-x) + y <= 1 */
1802  if( impltypes[j] == SCIP_BOUNDTYPE_UPPER )
1803  {
1804  if ( ! isneighbor[jidx] )
1805  {
1806  ++nneighbors;
1807  isneighbor[jidx] = TRUE;
1808  }
1809  }
1810  /* implication x + (1-y) <= 1 or (1-x) + (1-y) <= 1 */
1811  else
1812  {
1813  assert(impltypes[j] == SCIP_BOUNDTYPE_LOWER);
1814  if ( ! isneighbor[jidx + nbinvars] )
1815  {
1816  ++nneighbors;
1817  isneighbor[jidx + nbinvars] = TRUE;
1818  }
1819  }
1820  }
1821  }
1822 
1823  /* count cliques of the root */
1824  ncliques = (unsigned int) SCIPvarGetNCliques(vars[varsidx], varfixing);
1825  if( ncliques > 0 )
1826  {
1827  cliques = SCIPvarGetCliques(vars[varsidx], varfixing);
1828  assert(cliques != NULL);
1829 
1830  for( j = 0; j < ncliques; ++j )
1831  {
1832  ncliquevars = (unsigned int) SCIPcliqueGetNVars(cliques[j]);
1833  cliquevars = SCIPcliqueGetVars(cliques[j]);
1834  cliquevals = SCIPcliqueGetValues(cliques[j]);
1835 
1836  assert(cliquevars != NULL || ncliquevars == 0);
1837  assert(cliquevals != NULL || ncliquevars == 0);
1838 
1839  for( k = 0; k < ncliquevars; ++k )
1840  {
1841  unsigned int kidx;
1842 
1843  assert( cliquevars != NULL && cliquevals != NULL ); /* for lint */
1844 
1845  kidx = sepadata->mapping[SCIPvarGetProbindex(cliquevars[k])];
1846  assert(kidx < nbinvars);
1847 
1848  /* skip root */
1849  if( kidx == varsidx )
1850  continue;
1851 
1852  /* skip integral neighbors */
1853  if( SCIPisFeasIntegral(scip, vals[kidx]))
1854  continue;
1855 
1856  if( cliquevals[k] == TRUE )
1857  {
1858  if ( ! isneighbor[kidx] )
1859  {
1860  ++nneighbors;
1861  isneighbor[kidx] = TRUE;
1862  }
1863  }
1864  else
1865  {
1866  assert(cliquevals[k] == FALSE);
1867  if ( ! isneighbor[kidx + nbinvars] )
1868  {
1869  ++nneighbors;
1870  isneighbor[kidx+nbinvars] = TRUE;
1871  }
1872  }
1873  }
1874  }
1875  }
1876 
1877  /* root cannot be part of the next level */
1878  assert(! isneighbor[root]);
1879 
1880  /* allocate memory for sorting of root neighbors */
1881  SCIP_CALL( SCIPallocBufferArray(scip, &neighbors, (int) nneighbors) );
1882  SCIP_CALL( SCIPallocBufferArray(scip, &sortvals, (int) nneighbors) );
1883 
1884  k = 0;
1885  for( j = 0; j < graph->maxnodes; ++j )
1886  {
1887  if( isneighbor[j] )
1888  {
1889  assert(j != root);
1890  assert(!SCIPisFeasIntegral(scip, vals[j]));
1891 
1892  neighbors[k] = (int) j;
1893  sortvals[k] = MIN(1.0 - vals[j], vals[j]);
1894  ++k;
1895  }
1896  }
1897  assert(k == nneighbors);
1898 
1899  /* sort neighbors by fractionality */
1900  SCIPsortDownRealInt(sortvals, neighbors, (int) nneighbors);
1901 
1902  /* free temporary memory */
1903  SCIPfreeBufferArray(scip, &sortvals);
1904 
1905  /* insert sorted neighbors until level size limit is reached (or all neighbors are inserted) */
1906  for( j = 0; j < nneighbors && (*nnewlevel) <= sepadata->maxlevelsize; ++j )
1907  {
1908  int tmp;
1909 
1910  v = (unsigned int) neighbors[j];
1911  assert( v < 2 * nbinvars );
1912 
1913  /* only the root is contained in the levelgraph */
1914  assert(! inlevelgraph[v] || v == root+nbinvars || v == root-nbinvars);
1915 
1916  /* insert neighbor into levelgraph */
1917  ++(graph->nnodes);
1918  graph->level[v] = level + 1;
1919  inlevelgraph[v] = TRUE;
1920  newlevel[*nnewlevel] = v;
1921  ++(*nnewlevel);
1922 
1923  assert(! SCIPisFeasIntegral(scip, vals[varsidx]));
1924  assert(! SCIPisFeasIntegral(scip, vals[v]));
1925 
1926  graph->targetForward[graph->lastF] = (int) v;
1927  /* the computation of 1.0 - vals[v] if v is negated is ensured by the fact that v > nbinvars in this case */
1928  if( varfixing )
1929  {
1930  tmp = (int) SCIPfeasCeil(scip, sepadata->scale * (1.0 - vals[varsidx] - vals[v]));
1931  graph->weightForward[graph->lastF] = (unsigned int) MAX(tmp, sepadata->maxreference);
1932  }
1933  else
1934  {
1935  assert( ! varfixing );
1936  tmp = (int) SCIPfeasCeil(scip, sepadata->scale * (1.0 - (1.0 - vals[varsidx]) - vals[v]));
1937  graph->weightForward[graph->lastF] = (unsigned int) MAX(tmp, sepadata->maxreference);
1938  }
1939  ++(graph->lastF);
1940  ++(graph->narcs);
1941  if( graph->lastF == graph->sizeForward )
1942  {
1943  SCIP_CALL( checkArraySizesHeur(scip, graph, &(graph->sizeForward), &(graph->targetForward),
1944  &(graph->weightForward), NULL, NULL, success) );
1945 
1946  if( !(*success) )
1947  break;
1948  }
1949  }
1950 
1951  /* free temporary memory */
1952  SCIPfreeBufferArray(scip, &neighbors);
1953  SCIPfreeBufferArray(scip, &isneighbor);
1954 
1955  return SCIP_OKAY;
1956 }
1957 
1958 /** Find shortest path from start node to root
1959  *
1960  * We perform a BFS to find the shortest path to the root. D stores the distance to the start
1961  * node, P stores the partent nodes in the shortest path tree (-1 if node has not been reached).
1962  */
1963 static
1964 SCIP_RETCODE findShortestPathToRoot(
1965  SCIP* scip, /**< SCIP data structure */
1966  int scale, /**< scaling factor for edge weights */
1967  LEVELGRAPH* graph, /**< LEVELGRAPH data structure */
1968  unsigned int startnode, /**< start node for path search */
1969  unsigned int* distance, /**< distances of searched nodes from root */
1970  unsigned int* queue, /**< node queue for path search */
1971  SCIP_Bool* inQueue, /**< whether node is in the queue */
1972  int* parentTree /**< parent tree (-1 if no parent) */
1973  )
1974 {
1975  unsigned int i;
1976  int startQueue;
1977  int endQueue;
1978  unsigned int u;
1979  int v;
1980  unsigned int d;
1981 
1982  assert(scip != NULL);
1983  assert(graph != NULL);
1984  assert(graph->beginBackward != NULL);
1985  assert(graph->targetBackward != NULL);
1986  assert(graph->weightBackward != NULL);
1987  assert(distance != NULL);
1988  assert(queue != NULL);
1989  assert(inQueue != NULL);
1990  assert(parentTree != NULL);
1991 
1992  /* initialize distances */
1993  for( i = 0; i < graph->maxnodes; ++i )
1994  {
1995  distance[i] = 2*(graph->nnodes)*scale;
1996  parentTree[i] = -1;
1997  inQueue[i] = FALSE;
1998  }
1999  distance[startnode] = 0;
2000 
2001  /* initialize queue */
2002  startQueue = 0;
2003  endQueue = 0;
2004  queue[0] = startnode;
2005  v = 0;
2006 
2007  /* as long as queue is not empty */
2008  while( startQueue <= endQueue )
2009  {
2010  /* pop first node from queue */
2011  u = queue[startQueue];
2012  ++startQueue;
2013 
2014  /* check adjacent nodes */
2015  assert(graph->beginBackward[u] >= 0);
2016  i = (unsigned int) graph->beginBackward[u];
2017  for( v = graph->targetBackward[i]; v >= 0; v = graph->targetBackward[++i] )
2018  {
2019  /* distance to u via current arc: */
2020  d = distance[u] + graph->weightBackward[i];
2021 
2022  /* if we found a shorter connection */
2023  if( d < distance[v] )
2024  {
2025  distance[v] = d;
2026  parentTree[v] = (int) u;
2027 
2028  /* insert in queue if not already present */
2029  if( !inQueue[v] )
2030  {
2031  ++endQueue;
2032  queue[endQueue] = (unsigned int) v;
2033  inQueue[v] = TRUE;
2034  }
2035  }
2036  }
2037  /* it is not necessary to stop if we found the root (in this case there are no arcs left) and we stop anyway */
2038  }
2039  assert(parentTree[u] != -1);
2040 
2041  return SCIP_OKAY;
2042 }
2043 
2044 
2045 /** Block shortest path
2046  *
2047  * We traverse the shortest path found by findShortestPathToRoot() and block all neighbors (in the
2048  * original graph) of nodes in the path, i.e., we set blocked to TRUE. We do not block neighbors of
2049  * the root node, since they have to be used. For the start node we only block nodes on the
2050  * previous layers,
2051  *
2052  * @see findShortestPathToRoot()
2053  */
2054 static
2055 SCIP_RETCODE blockRootPath(
2056  SCIP* scip, /**< SCIP data structure */
2057  LEVELGRAPH* graph, /**< LEVELGRAPH data structure */
2058  unsigned int startnode, /**< start node */
2059  SCIP_Bool* inlevelgraph, /**< nodes in new graph corr. to old graph (-1 if unassigned) */
2060  SCIP_Bool* blocked, /**< whether node is blocked */
2061  int* parentTree, /**< parent tree */
2062  unsigned int root /**< root of the current level graph */
2063  )
2064 {
2065  unsigned int u;
2066  unsigned int i;
2067  int v;
2068 
2069  assert(scip != NULL);
2070  assert(graph != NULL);
2071  assert(graph->level != NULL);
2072  assert(graph->beginForward != NULL);
2073  assert(graph->targetForward != NULL);
2074  assert(graph->beginBackward != NULL);
2075  assert(graph->targetBackward != NULL);
2076  assert(graph->sourceAdj != NULL);
2077  assert(graph->targetAdj != NULL);
2078  assert(inlevelgraph != NULL);
2079  assert(blocked != NULL);
2080  assert(parentTree != NULL);
2081 
2082  assert(parentTree[root] >= 0);
2083 
2084  /* follow the path from the predecessor of root to the start node and block all neighbors */
2085  u = (unsigned int) parentTree[root];
2086  while( u != startnode )
2087  {
2088  /* block neighbors of u in higher level */
2089  i = (unsigned int) graph->beginForward[u];
2090  for( v = graph->targetForward[i]; v >= 0; v = graph->targetForward[++i] )
2091  {
2092  assert(inlevelgraph[v]);
2093  blocked[v] = TRUE;
2094  }
2095 
2096  /* block neighbors of u in lower level */
2097  i = (unsigned int) graph->beginBackward[u];
2098  for( v = graph->targetBackward[i]; v >= 0; v = graph->targetBackward[++i] )
2099  {
2100  assert(inlevelgraph[v]);
2101  blocked[v] = TRUE;
2102  }
2103 
2104  /* block neighbors of u in same level */
2105  assert(graph->level[u] > 0);
2106  for( i = graph->levelAdj[graph->level[u]]; i < graph->levelAdj[graph->level[u]+1]; ++i )
2107  {
2108  assert(graph->sourceAdj[i] < graph->targetAdj[i]);
2109  assert(graph->level[graph->sourceAdj[i]] == graph->level[graph->targetAdj[i]]);
2110 
2111  /* remember that these arcs are only stored for one direction */
2112  if( graph->sourceAdj[i] == u )
2113  {
2114  blocked[graph->targetAdj[i]] = TRUE;
2115  }
2116  if( graph->targetAdj[i] == u )
2117  {
2118  blocked[graph->sourceAdj[i]] = TRUE;
2119  }
2120  }
2121 
2122  /* get next node on the path */
2123  u = (unsigned int) parentTree[u];
2124  }
2125  assert(u == startnode);
2126 
2127  /* block nodes adjacent to start node on previous level */
2128  assert(graph->beginBackward[u] > 0);
2129  i = (unsigned int) graph->beginBackward[u];
2130  for( v = graph->targetBackward[i]; v >= 0; v = graph->targetBackward[++i] )
2131  blocked[v] = TRUE;
2132 
2133  return SCIP_OKAY;
2134 }
2135 
2136 
2137 /** Find shortest path from root to target node
2138  *
2139  * We perform a BFS to find the shortest path from the root. The only difference to
2140  * findShortestPathToRoot() is that we avoid nodes that are blocked.
2141  */
2142 static
2143 SCIP_RETCODE
2144 findUnblockedShortestPathToRoot(
2145  SCIP* scip, /**< SCIP data structure */
2146  int scale, /**< scaling factor for edge weights */
2147  LEVELGRAPH* graph, /**< LEVELGRAPH data structure */
2148  unsigned int startnode, /**< start node for path search */
2149  unsigned int* distance, /**< distances of searched nodes from root */
2150  unsigned int* queue, /**< node queue for path search */
2151  SCIP_Bool* inQueue, /**< whether node is in the queue */
2152  int* parentTreeBackward, /**< parent tree (-1 if no parent) */
2153  unsigned int root, /**< root of the current level graph */
2154  SCIP_Bool* blocked /**< whether nodes can be used */
2155  )
2156 {
2157  unsigned int i;
2158  int startQueue;
2159  int endQueue;
2160  unsigned int u;
2161  int v;
2162  unsigned int d;
2163  int* parentTree;
2164  int* transform;
2165 
2166  assert(scip != NULL);
2167  assert(graph != NULL);
2168  assert(graph->beginBackward != NULL);
2169  assert(graph->targetBackward != NULL);
2170  assert(graph->weightBackward != NULL);
2171  assert(distance != NULL);
2172  assert(queue != NULL);
2173  assert(inQueue != NULL);
2174 
2175  /* allocate temporary memory */
2176  SCIP_CALL( SCIPallocBufferArray(scip, &parentTree, (int) graph->maxnodes) );
2177  SCIP_CALL( SCIPallocBufferArray(scip, &transform, (int) graph->maxnodes) );
2178 
2179  assert(parentTree != NULL);
2180  assert(transform != NULL);
2181 
2182  /* initialize distances */
2183  for( i = 0; i < graph->maxnodes; ++i )
2184  {
2185  distance[i] = 2*(graph->nnodes)*scale;
2186  parentTree[i] = -1;
2187  parentTreeBackward[i] = -1;
2188  transform[i] = -1;
2189  inQueue[i] = FALSE;
2190  }
2191  distance[startnode] = 0;
2192 
2193  /* initialize queue */
2194  startQueue = 0;
2195  endQueue = 0;
2196  queue[0] = startnode;
2197  v = 0;
2198 
2199  /* as long as queue is not empty */
2200  while( startQueue <= endQueue )
2201  {
2202  /* pop first node from queue */
2203  u = queue[startQueue];
2204  ++startQueue;
2205 
2206  /* check adjacent nodes */
2207  assert(graph->beginBackward[u] >= 0);
2208  i = (unsigned int) graph->beginBackward[u];
2209  for( v = graph->targetBackward[i]; v >= 0; v = graph->targetBackward[++i] )
2210  {
2211  if( blocked[v] && v != (int) root)
2212  continue;
2213 
2214  /* distance to u via current arc: */
2215  d = distance[u] + graph->weightBackward[i];
2216 
2217  /* if we found a shorter connection */
2218  if( d < distance[v] )
2219  {
2220  distance[v] = d;
2221  parentTree[v] = (int) u;
2222 
2223  /* insert in queue if not already present */
2224  if( !inQueue[v] )
2225  {
2226  ++endQueue;
2227  queue[endQueue] = (unsigned int) v;
2228  inQueue[v] = TRUE;
2229  }
2230  }
2231  }
2232  /* it is not necessary to stop if we found the root (in this case there are no arcs left) and we stop anyway */
2233  }
2234 
2235  /* reverse order such that it is a path from the root */
2236  v = (int) root;
2237  transform[0] = (int) root;
2238  i = 1;
2239  while(parentTree[v] >= 0)
2240  {
2241  transform[i] = parentTree[v];
2242  ++i;
2243  v = parentTree[v];
2244  }
2245  --i;
2246  while(i > 0)
2247  {
2248  parentTreeBackward[transform[i]] = transform[i-1];
2249  --i;
2250  }
2251 
2252  /* free temporary memory */
2253  SCIPfreeBufferArray(scip, &transform);
2254  SCIPfreeBufferArray(scip, &parentTree);
2255 
2256  return SCIP_OKAY;
2257 }
2258 
2259 /** create next level of level graph for odd cycle separation
2260  *
2261  * @see separateHeur()
2262  */
2263 static
2264 SCIP_RETCODE createNextLevel(
2265  SCIP* scip, /**< SCIP data structure */
2266  SCIP_SEPADATA* sepadata, /**< separator data structure */
2267  SCIP_VAR** vars, /**< problem variables */
2268  SCIP_Real* vals, /**< values of the binary variables in the current LP relaxation */
2269  LEVELGRAPH* graph, /**< LEVELGRAPH data structure */
2270  unsigned int level, /**< number of current level */
2271  SCIP_Bool* inlevelgraph, /**< flag array if node is already inserted in level graph */
2272  unsigned int* curlevel, /**< array of nodes of the current level */
2273  unsigned int ncurlevel, /**< number of nodes of the current level */
2274  unsigned int* newlevel, /**< array of nodes of the next level */
2275  unsigned int* nnewlevel, /**< number of nodes of the next level */
2276  SCIP_Bool* success /**< FALSE, iff memory reallocation fails */
2277  )
2278 {
2279  unsigned int i;
2280  unsigned int nbinvars;
2281  unsigned int nAdj;
2282 
2283  assert(scip != NULL);
2284  assert(vars != NULL);
2285  assert(vals != NULL);
2286  assert(graph != NULL);
2287  assert(graph->level != NULL);
2288  assert(graph->beginForward != NULL);
2289  assert(graph->targetForward != NULL);
2290  assert(graph->weightForward != NULL);
2291  assert(graph->beginBackward != NULL);
2292  assert(graph->targetBackward != NULL);
2293  assert(graph->weightBackward != NULL);
2294  assert(graph->beginAdj != NULL);
2295  assert(graph->levelAdj != NULL);
2296  assert(graph->sourceAdj != NULL);
2297  assert(graph->targetAdj != NULL);
2298  assert(graph->weightAdj != NULL);
2299  assert(inlevelgraph != NULL);
2300  assert(curlevel != NULL);
2301  assert(newlevel != NULL);
2302  assert(success != NULL);
2303 
2304  *nnewlevel = 0;
2305  nAdj = 0;
2306  assert(graph->maxnodes % 2 == 0);
2307  nbinvars = (graph->maxnodes)/2;
2308 
2309  /* for every node in current level add its implications and assign its neighbors to the next
2310  * level, if neighbor is not already existing in the level graph
2311  */
2312  for( i = 0; i < ncurlevel; ++i )
2313  {
2314  unsigned int negated;
2315  unsigned int u;
2316 
2317  /* get node */
2318  u = curlevel[i];
2319  assert(u < graph->maxnodes);
2320  assert(graph->level[u] == level);
2321  assert(graph->beginForward[u] < 0);
2322  assert(graph->beginBackward[u] < 0);
2323  assert(graph->beginAdj[u] < 0);
2324  assert(inlevelgraph[u]);
2325 
2326  /* get negated */
2327  if( u < nbinvars )
2328  negated = u + nbinvars;
2329  else
2330  negated = u - nbinvars;
2331  assert(negated < graph->maxnodes);
2332  assert(negated < nbinvars || u < nbinvars);
2333  assert(negated >= nbinvars || u >= nbinvars);
2334 
2335  /* initialize adjacency lists for node u */
2336  graph->beginForward[u] = (int) graph->lastF;
2337  graph->beginBackward[u] = (int) graph->lastB;
2338  graph->beginAdj[u] = (int) (graph->levelAdj[level+1] + nAdj);
2339 
2340  /* if we want to add arcs between a variable and its negated */
2341  if( sepadata->addselfarcs )
2342  {
2343  /* add negated variable, if not existing in the levelgraph,
2344  * but if the level contains more nodes than allowed
2345  * (defined by percent per level plus offset),
2346  * we skip the rest of the nodes
2347  */
2348  if( !inlevelgraph[negated] && (*nnewlevel) <= sepadata->maxlevelsize )
2349  {
2350  ++(graph->nnodes);
2351  graph->level[negated] = level+1;
2352  inlevelgraph[negated] = TRUE;
2353  newlevel[*nnewlevel] = negated;
2354  ++(*nnewlevel);
2355  }
2356  assert( *nnewlevel > sepadata->maxlevelsize || inlevelgraph[negated] );
2357 
2358  /* add self-arc if negated variable is on a neighbored level */
2359  if( inlevelgraph[negated] && ((graph->level[negated] == level - 1)
2360  || (graph->level[negated] == level) || (graph->level[negated] == level + 1)) )
2361  {
2362  /* add arc from u to its negated variable */
2363  SCIP_CALL( addArc(scip, graph, u, negated, level, 0, &nAdj, success) );
2364  if( !(*success) )
2365  return SCIP_OKAY;
2366  }
2367  }
2368 
2369  /* insert level of sorted root neighbors (if requested) */
2370  if( graph->nlevels == 0 && sepadata->sortrootneighbors )
2371  {
2372  SCIP_CALL( insertSortedRootNeighbors(scip, graph, nbinvars, ncurlevel, u, vals, vars,
2373  sepadata, nnewlevel, inlevelgraph, level, newlevel, success) );
2374  }
2375  else
2376  {
2377  /* add arc from u to all other neighbors of variable implication graph */
2378  SCIP_CALL( addNextLevelBinImpls(scip, sepadata, vars, vals, u, graph, level, inlevelgraph,
2379  newlevel, nnewlevel, &nAdj, success) );
2380 
2381  if( !(*success) )
2382  return SCIP_OKAY;
2383 
2384  SCIP_CALL( addNextLevelCliques(scip, sepadata, vars, vals, u, graph, level, inlevelgraph,
2385  newlevel, nnewlevel, &nAdj, success) );
2386  }
2387  if( !(*success) )
2388  return SCIP_OKAY;
2389 
2390  /* every node has a backward arc */
2391  assert(graph->lastB > (unsigned int) graph->beginBackward[u] || graph->nlevels == 0 );
2392 
2393  /* root has outgoing arcs otherwise we would have skipped it */
2394  assert(graph->lastF > 0);
2395 
2396  /* close adjacency lists */
2397  graph->targetForward[graph->lastF] = -1;
2398  ++(graph->lastF);
2399  if( graph->lastF == graph->sizeForward )
2400  {
2401  SCIP_CALL( checkArraySizesHeur(scip, graph, &(graph->sizeForward), &(graph->targetForward),
2402  &(graph->weightForward), NULL, NULL, success) );
2403 
2404  if( !(*success) )
2405  return SCIP_OKAY;
2406  }
2407  graph->targetBackward[graph->lastB] = -1;
2408  ++(graph->lastB);
2409  if( graph->lastB == graph->sizeBackward )
2410  {
2411  SCIP_CALL( checkArraySizesHeur(scip, graph, &(graph->sizeBackward), &(graph->targetBackward),
2412  &(graph->weightBackward), NULL, NULL, success) );
2413 
2414  if( !(*success) )
2415  return SCIP_OKAY;
2416  }
2417 
2418  /* terminate adjacency list with 0 for current level lifting */
2419  graph->sourceAdj[graph->levelAdj[level+1]+nAdj] = 0;
2420  graph->targetAdj[graph->levelAdj[level+1]+nAdj] = 0;
2421  }
2422  graph->levelAdj[level+2] = graph->levelAdj[level+1]+nAdj;
2423 
2424  return SCIP_OKAY;
2425 }
2426 
2427 /** The heuristic method for finding odd cycles by Hoffman, Padberg uses a level graph
2428  * which is constructed as follows:
2429  * First we choose a node (i.e. a variable of the problem or its negated) as root
2430  * and assign it to level 0 (and no other node is assigned to level 0).
2431  * All neighbors of the root are assigned to level 1 and the arcs between are added.
2432  *
2433  * In general:
2434  * All neighbors of nodes in level i that are assigned to level i+1, if they do not already appear in levels <= i.
2435  * All arcs between nodes in level i and their neighbors are added.
2436  *
2437  * In the construction we only take nodes that are contained in the fractional graph,
2438  * i.e., their current LP values are not integral.
2439  *
2440  * Since SCIP stores implications between original and negated variables,
2441  * the level graph has at most twice the number of fractional binary variables of the problem.
2442  *
2443  * Since the implication graph of SCIP is (normally) incomplete,
2444  * it is possible to use arcs between an original variable and its negated
2445  * to obtain more cycles which are valid but not found due to missing links.
2446  */
2447 static
2448 SCIP_RETCODE separateHeur(
2449  SCIP* scip, /**< SCIP data structure */
2450  SCIP_SEPA* sepa, /**< separator */
2451  SCIP_SEPADATA* sepadata, /**< separator data structure */
2452  SCIP_SOL* sol, /**< given primal solution */
2453  SCIP_RESULT* result /**< pointer to store the result of the separation call */
2454  )
2455 {
2456  /* memory for variable data */
2457  SCIP_VAR** scipvars; /* variables of the current SCIP (unsorted) */
2458  SCIP_VAR** vars; /* variables of the current SCIP (sorted if requested) */
2459  SCIP_Real* vals; /* LP-values of the variables (and negated variables) */
2460  unsigned int nbinvars; /* number of nodecandidates for implicationgraph */
2461  unsigned int* incut; /* flag array for check if a variable is already covered by a cut */
2462 
2463  /* storage for levelgraph */
2464  LEVELGRAPH graph;
2465  unsigned int* curlevel;
2466  unsigned int* newlevel;
2467  unsigned int ncurlevel;
2468  unsigned int nnewlevel;
2469  SCIP_Bool* inlevelgraph;
2470 
2471  /* storage for path finding */
2472  unsigned int* queue;
2473  SCIP_Bool* inQueue;
2474  int* parentTree;
2475  int* parentTreeBackward;
2476  unsigned int* distance;
2477  SCIP_Bool* blocked;
2478 
2479  /* counter and limits */
2480  unsigned int maxroots; /* maximum of level graph roots */
2481  unsigned int rootcounter; /* counter of tried roots */
2482  unsigned int ncutsroot; /* counter of cuts per root */
2483  unsigned int ncutslevel; /* counter of cuts per level */
2484 
2485  unsigned int i;
2486  unsigned int j;
2487  unsigned int k;
2488 
2489  int nscipbinvars;
2490  int nscipintvars;
2491  int nscipimplvars;
2492  int nintegral;
2493  int l;
2494 
2495  assert(scip != NULL);
2496  assert(sepadata != NULL);
2497  assert(result != NULL);
2498 
2499  SCIP_CALL( SCIPgetVarsData(scip, &scipvars, NULL, &nscipbinvars, &nscipintvars, &nscipimplvars, NULL) );
2500  assert(nscipbinvars >= 0);
2501  assert(nscipintvars >= 0);
2502  assert(nscipimplvars >= 0);
2503 
2504  nintegral = nscipbinvars + nscipintvars + nscipimplvars;
2505  assert(scipvars != NULL || nintegral == 0);
2506 
2507  /* collect binary variables, including implicit binary */
2508  SCIP_CALL( SCIPallocBufferArray(scip, &vars, nintegral) );
2509  for (l = 0; l < nscipbinvars; ++l)
2510  vars[l] = scipvars[l]; /*lint !e613*/
2511 
2512  nbinvars = (unsigned int) nscipbinvars;
2513  for (l = nscipbinvars; l < nintegral; ++l)
2514  {
2515  assert( SCIPvarGetType(scipvars[l]) != SCIP_VARTYPE_CONTINUOUS ); /*lint !e613*/
2516  if ( SCIPvarIsBinary(scipvars[l]) ) /*lint !e613*/
2517  vars[nbinvars++] = scipvars[l]; /*lint !e613*/
2518  }
2519 
2520  if( nbinvars == 0 )
2521  {
2522  SCIPfreeBufferArray(scip, &vars);
2523  return SCIP_OKAY;
2524  }
2525 
2526  /* initialize flag array to avoid multiple cuts per variable, if requested by user-flag */
2527  SCIP_CALL( SCIPallocBufferArray(scip, &vals, (int) (2 * nbinvars)) );
2528 
2529  /* prepare values */
2530  assert( vars != NULL );
2531  switch( sepadata->sortswitch )
2532  {
2533  case UNSORTED :
2534  /* if no sorting is requested, we use the normal variable array */
2535  break;
2536 
2537  case MAXIMAL_LPVALUE :
2538  /* store lp-values */
2539  for( i = 0; i < nbinvars; ++i )
2540  vals[i] = SCIPgetSolVal(scip, sol, vars[i]);
2541 
2542  /* sort by lp-value, maximal first */
2543  SCIPsortDownRealPtr(vals, (void**)vars, (int) nbinvars);
2544  break;
2545 
2546  case MINIMAL_LPVALUE :
2547  /* store lp-values */
2548  for( i = 0; i < nbinvars; ++i )
2549  vals[i] = SCIPgetSolVal(scip, sol, vars[i]);
2550 
2551  /* sort by lp-value, minimal first */
2552  SCIPsortRealPtr(vals, (void**)vars, (int) nbinvars);
2553  break;
2554 
2555  case MAXIMAL_FRACTIONALITY :
2556  /* store lp-values and determine fractionality */
2557  for( i = 0; i < nbinvars; ++i )
2558  {
2559  vals[i] = SCIPgetSolVal(scip, sol, vars[i]);
2560  vals[i] = MIN(1.0 - vals[i], vals[i]);
2561  }
2562 
2563  /* sort by fractionality, maximal first */
2564  SCIPsortDownRealPtr(vals, (void**)vars, (int) nbinvars);
2565  break;
2566 
2567  case MINIMAL_FRACTIONALITY :
2568  /* store lp-values and determine fractionality */
2569  for( i = 0; i < nbinvars; ++i )
2570  {
2571  vals[i] = SCIPgetSolVal(scip, sol, vars[i]);
2572  vals[i] = MIN(1.0 - vals[i], vals[i]);
2573  }
2574 
2575  /* sort by fractionality, minimal first */
2576  SCIPsortRealPtr(vals, (void**)vars, (int) nbinvars);
2577  break;
2578 
2579  default :
2580  SCIPerrorMessage("invalid sortswitch value\n");
2581  SCIPABORT();
2582  return SCIP_INVALIDDATA; /*lint !e527*/
2583  }
2584  assert(vars != NULL);
2585 
2586  /* create mapping for getting the index of a variable via its probindex to the index in the sorted variable array */
2587  SCIP_CALL( SCIPallocBufferArray(scip, &(sepadata->mapping), nintegral) );
2588 
2589  /* initialize LP value and cut flag for all variables */
2590  for( i = 0; i < nbinvars; ++i )
2591  {
2592  assert( 0 <= SCIPvarGetProbindex(vars[i]) && SCIPvarGetProbindex(vars[i]) < nintegral); /* since binary, integer, and implicit variables are first */
2593  sepadata->mapping[SCIPvarGetProbindex(vars[i])] = i;
2594  vals[i] = SCIPgetSolVal(scip, sol, vars[i]); /* need to get new values, since they might be corrupted */
2595  }
2596 
2597  for( i = nbinvars; i < 2*nbinvars; ++i )
2598  vals[i] = 1.0 - vals[i - nbinvars];
2599 
2600  /* determine size of level graph */
2601  graph.maxnodes = 2 * nbinvars;
2602 
2603  /* the implication graph is redundant and therefore more implications and clique arcs may occur than should be possible
2604  * @todo later: filtering of edges which were already added
2605  */
2606  /* graph.maxarcs = nbinvars*(2*nbinvars-1); */ /* = 2*nbinvars*(2*nbinvars-1)/2 */
2607  graph.maxarcs = UINT_MAX;
2608 
2609  /* set sizes for graph memory storage */
2610  graph.sizeForward = 100 * graph.maxnodes;
2611  graph.sizeBackward = 100 * graph.maxnodes;
2612  graph.sizeAdj = 100 * graph.maxnodes;
2613 
2614  /* allocate memory for level graph structure */
2615  SCIP_CALL( SCIPallocBufferArray(scip, &graph.level, (int) graph.maxnodes) );
2616  SCIP_CALL( SCIPallocBufferArray(scip, &graph.beginForward, (int) graph.maxnodes) );
2617  SCIP_CALL( SCIPallocBufferArray(scip, &graph.beginBackward, (int) graph.maxnodes) );
2618  SCIP_CALL( SCIPallocBufferArray(scip, &graph.targetForward, (int) MIN(graph.sizeForward, graph.maxarcs)) );
2619  SCIP_CALL( SCIPallocBufferArray(scip, &graph.targetBackward, (int) MIN(graph.sizeBackward, graph.maxarcs)) );
2620  SCIP_CALL( SCIPallocBufferArray(scip, &graph.weightForward, (int) MIN(graph.sizeForward, graph.maxarcs)) );
2621  SCIP_CALL( SCIPallocBufferArray(scip, &graph.weightBackward, (int) MIN(graph.sizeBackward, graph.maxarcs)) );
2622 
2623  SCIP_CALL( SCIPallocBufferArray(scip, &curlevel, (int) graph.maxnodes) );
2624  SCIP_CALL( SCIPallocBufferArray(scip, &newlevel, (int) graph.maxnodes) );
2625  SCIP_CALL( SCIPallocBufferArray(scip, &graph.beginAdj, (int) graph.maxnodes) );
2626  SCIP_CALL( SCIPallocBufferArray(scip, &graph.sourceAdj, (int) MIN(graph.sizeAdj, graph.maxarcs)) );
2627  SCIP_CALL( SCIPallocBufferArray(scip, &graph.targetAdj, (int) MIN(graph.sizeAdj, graph.maxarcs)) );
2628  SCIP_CALL( SCIPallocBufferArray(scip, &graph.weightAdj, (int) MIN(graph.sizeAdj, graph.maxarcs)) );
2629  SCIP_CALL( SCIPallocBufferArray(scip, &graph.levelAdj, (int) graph.maxnodes) );
2630  SCIP_CALL( SCIPallocBufferArray(scip, &inlevelgraph, (int) graph.maxnodes) );
2631 
2632  SCIP_CALL( SCIPallocBufferArray(scip, &queue, (int) graph.maxnodes) );
2633  SCIP_CALL( SCIPallocBufferArray(scip, &inQueue, (int) graph.maxnodes) );
2634  SCIP_CALL( SCIPallocBufferArray(scip, &parentTree, (int) graph.maxnodes) );
2635  SCIP_CALL( SCIPallocBufferArray(scip, &parentTreeBackward, (int) graph.maxnodes) );
2636  SCIP_CALL( SCIPallocBufferArray(scip, &distance, (int) graph.maxnodes) );
2637  SCIP_CALL( SCIPallocBufferArray(scip, &blocked, (int) graph.maxnodes) );
2638 
2639  SCIP_CALL( SCIPallocBufferArray(scip, &incut, (int) (2 * nbinvars)) );
2640 
2641  /* initialize cut flag for all variables */
2642  BMSclearMemoryArray(incut, 2*nbinvars);
2643 
2644  /* determine the number of level graph roots */
2645  maxroots = (unsigned int) SCIPceil(scip, sepadata->offsettestvars + (0.02 * nbinvars * sepadata->percenttestvars));
2646  sepadata->maxlevelsize = (unsigned int) SCIPceil(scip, sepadata->offsetnodeslevel + 0.01 * sepadata->maxpernodeslevel * graph.maxnodes);
2647  rootcounter = 0;
2648 
2649  /* check each node as root */
2650  for( i = (unsigned int) sepadata->lastroot; i < graph.maxnodes && rootcounter < maxroots
2651  && sepadata->ncuts - sepadata->oldncuts < (unsigned int) sepadata->maxsepacutsround
2652  && !SCIPisStopped(scip) ; ++i )
2653  {
2654  /* skip node if it is already covered by a cut and if we do not want to search cycles starting
2655  * with a node already covered by a cut */
2656  if( incut[i] && ! sepadata->multiplecuts )
2657  continue;
2658 
2659  /* skip variable if its LP-value is not fractional */
2660  if( SCIPisFeasIntegral(scip, vals[i]) )
2661  continue;
2662 
2663  /* consider original and negated variable pair and skip variable if there is only one edge leaving the pair */
2664  if( (SCIPvarGetNBinImpls(vars[i % nbinvars], TRUE) + SCIPvarGetNBinImpls(vars[i % nbinvars], FALSE) < 2)
2665  && (SCIPvarGetNCliques(vars[i % nbinvars], TRUE) + SCIPvarGetNCliques(vars[i % nbinvars], FALSE) < 1) )
2666  continue;
2667 
2668  /* skip variable having too less implics and cliques itself */
2669  if( i < nbinvars )
2670  {
2671  if( SCIPvarGetNBinImpls(vars[i % nbinvars], TRUE ) < 1 && SCIPvarGetNCliques(vars[i % nbinvars], TRUE ) < 1 )
2672  continue;
2673 
2674  if( !(sepadata->addselfarcs) && SCIPvarGetNBinImpls(vars[i % nbinvars], TRUE ) < 2
2675  && SCIPvarGetNCliques(vars[i % nbinvars], TRUE ) < 1 )
2676  continue;
2677  }
2678  else
2679  {
2680  if( SCIPvarGetNBinImpls(vars[i % nbinvars], FALSE) < 1 && SCIPvarGetNCliques(vars[i % nbinvars], FALSE) < 1 )
2681  continue;
2682 
2683  if( !(sepadata->addselfarcs) && SCIPvarGetNBinImpls(vars[i % nbinvars], FALSE) < 2
2684  && SCIPvarGetNCliques(vars[i % nbinvars], FALSE) < 1 )
2685  continue;
2686  }
2687 
2688  /* node is actually considered as root node for the level graph */
2689  ++rootcounter;
2690  ncutsroot = 0;
2691 
2692  /* initialize graph */
2693  for( j = 0; j < graph.maxnodes; ++j)
2694  {
2695  graph.beginForward[j] = -1;
2696  graph.beginBackward[j] = -1;
2697  graph.beginAdj[j] = -1;
2698  inlevelgraph[j] = FALSE;
2699  blocked[j] = FALSE;
2700  }
2701  graph.lastF = 0;
2702  graph.lastB = 0;
2703  graph.nlevels = 0;
2704  graph.narcs = 0;
2705 
2706  /* insert root (first level contains root only) */
2707  inlevelgraph[i] = TRUE;
2708  graph.level[i] = 0;
2709  graph.levelAdj[0] = 0;
2710  graph.nnodes = 1;
2711  curlevel[0] = i;
2712  ncurlevel = 1;
2713 
2714  /* there are no arcs inside the root level */
2715  graph.levelAdj[graph.nlevels+1] = 0;
2716 
2717  /* create new levels until there are not more nodes for a new level */
2718  do
2719  {
2720  SCIP_Bool success;
2721  success = TRUE;
2722 
2723  /* all neighbors of nodes in level i that are assigned to level i+1,
2724  if they do not already appear in levels <= i. */
2725  SCIP_CALL( createNextLevel(scip, sepadata, vars, vals, &graph, graph.nlevels, inlevelgraph,
2726  curlevel, ncurlevel, newlevel, &nnewlevel, &success) );
2727 
2728  if( !success )
2729  goto TERMINATE;
2730 
2731  /* search for odd holes */
2732  if( graph.nlevels > 0 && (sepadata->includetriangles || graph.nlevels > 1) )
2733  {
2734  unsigned int maxcutslevel;
2735 
2736  ncutslevel = 0;
2737 
2738  /* calculate maximal cuts in this level due to cut limitations (per level, per root, per separation round) */
2739  maxcutslevel = (unsigned int) sepadata->maxcutslevel;
2740  maxcutslevel = (unsigned int) MIN(maxcutslevel, ncutsroot-sepadata->maxcutsroot);
2741  maxcutslevel = (unsigned int) MIN(maxcutslevel, sepadata->maxsepacutsround + sepadata->oldncuts - sepadata->ncuts);
2742 
2743  /* for each cross edge in this level find both shortest paths to root (as long as no limits are reached) */
2744  for( j = graph.levelAdj[graph.nlevels+1]; j < graph.levelAdj[graph.nlevels+2]
2745  && ncutslevel < maxcutslevel && !SCIPisStopped(scip); ++j )
2746  {
2747  unsigned int ncyclevars;
2748  unsigned int u;
2749 
2750  /* storage for cut generation */
2751  unsigned int* pred; /* predecessor list */
2752  SCIP_Bool* incycle; /* flag array for check of double variables in found cycle */
2753 
2754  assert(graph.sourceAdj[j] < graph.targetAdj[j]);
2755 
2756  /* find shortest path from source to root and update weight of cycle */
2757  SCIP_CALL( findShortestPathToRoot(scip, sepadata->scale, &graph, graph.sourceAdj[j], distance, queue, inQueue, parentTree) );
2758 
2759 #ifndef NDEBUG
2760  /* check that this path ends in the root node */
2761  u = i;
2762  k = 1;
2763  while( u != graph.sourceAdj[j] )
2764  {
2765  assert(parentTree[u] != -1 && k <= graph.maxnodes);
2766  u = (unsigned int) parentTree[u];
2767  ++k;
2768  }
2769 #endif
2770 
2771  /* block all nodes that are adjacent to nodes of the first path */
2772  for( k = 0; k < graph.nnodes; ++k )
2773  blocked[k] = FALSE;
2774  SCIP_CALL( blockRootPath(scip, &graph, graph.sourceAdj[j], inlevelgraph, blocked, parentTree, i) );
2775 
2776  /* if the target is block, no violated odd hole can be found */
2777  if( blocked[graph.targetAdj[j]] )
2778  continue;
2779 
2780  /* find shortest path from root to target node avoiding blocked nodes */
2781  SCIP_CALL( findUnblockedShortestPathToRoot(scip, sepadata->scale, &graph,
2782  graph.targetAdj[j], distance, queue, inQueue, parentTreeBackward, i, blocked) );
2783 
2784  /* no odd cycle cut found */
2785  if( parentTreeBackward[graph.targetAdj[j]] < 0 )
2786  continue;
2787 
2788  /* allocate and initialize predecessor list and flag array representing odd cycle */
2789  SCIP_CALL( SCIPallocBufferArray(scip, &pred, (int) (2 * nbinvars)) );
2790  SCIP_CALL( SCIPallocBufferArray(scip, &incycle, (int) (2 * nbinvars)) );
2791  for( k = 0; k < 2 * nbinvars; ++k )
2792  {
2793  pred[k] = DIJKSTRA_UNUSED;
2794  incycle[k] = FALSE;
2795  }
2796  ncyclevars = 0;
2797  success = TRUE;
2798 
2799  /* check cycle for x-neg(x)-sub-cycles and clean them
2800  * (note that a variable cannot appear twice in a cycle since it is only once in the graph)
2801  * convert parentTreeBackward and parentTree to pred&incycle structure for generateOddCycleCut
2802  */
2803  u = graph.targetAdj[j];
2804 
2805  /* add path to root to cycle */
2806  while( success && u != i )
2807  {
2808  /* insert u in predecessor list */
2809  pred[u] = (unsigned int) parentTreeBackward[u];
2810 
2811  /* remove pairs of original and negated variable from cycle */
2812  SCIP_CALL( cleanCycle(scip, pred, incycle, incut, u, graph.targetAdj[j], nbinvars, &ncyclevars,
2813  sepadata->repaircycles, sepadata->allowmultiplecuts, &success) );
2814 
2815  assert(parentTreeBackward[u] >= 0 || u == i);
2816 
2817  /* select next node on path */
2818  u = (unsigned int) parentTreeBackward[u];
2819  }
2820 
2821  /* add path from root to cycle */
2822  while( success && u != graph.sourceAdj[j] )
2823  {
2824  /* insert u in predecessor list */
2825  pred[u] = (unsigned int) parentTree[u];
2826 
2827  /* remove pairs of original and negated variable from cycle */
2828  SCIP_CALL( cleanCycle(scip, pred, incycle, incut, u, graph.targetAdj[j], nbinvars, &ncyclevars,
2829  sepadata->repaircycles, sepadata->allowmultiplecuts, &success) );
2830 
2831  /* select next node on path */
2832  u = (unsigned int) parentTree[u];
2833  }
2834  assert(!success || u == graph.sourceAdj[j]);
2835 
2836  /* close the cycle */
2837  if( success )
2838  {
2839  pred[u] = graph.targetAdj[j];
2840 
2841  /* remove pairs of original and negated variable from cycle */
2842  SCIP_CALL( cleanCycle(scip, pred, incycle, incut, u, graph.targetAdj[j], nbinvars, &ncyclevars,
2843  sepadata->repaircycles, sepadata->allowmultiplecuts, &success) );
2844  }
2845 
2846  /* generate cut (if cycle is valid) */
2847  if(success)
2848  {
2849  unsigned int oldncuts = sepadata->ncuts;
2850 
2851  sepadata->levelgraph = &graph;
2852  SCIP_CALL( generateOddCycleCut(scip, sepa, sol, vars, nbinvars, graph.targetAdj[j], pred, ncyclevars,
2853  incut, vals, sepadata, result) );
2854 #ifndef NDEBUG
2855  sepadata->levelgraph = NULL;
2856 #endif
2857  if(oldncuts < sepadata->ncuts)
2858  {
2859  ++ncutsroot;
2860  ++ncutslevel;
2861  }
2862  }
2863 
2864  /* free temporary memory */
2865  SCIPfreeBufferArray(scip, &incycle);
2866  SCIPfreeBufferArray(scip, &pred);
2867 
2868  if ( *result == SCIP_CUTOFF )
2869  break;
2870  }
2871  }
2872 
2873  if ( *result == SCIP_CUTOFF )
2874  break;
2875 
2876  /* copy new level to current one */
2877  ++(graph.nlevels);
2878  for( j = 0; j < nnewlevel; ++j )
2879  curlevel[j] = newlevel[j];
2880  ncurlevel = nnewlevel;
2881  }
2882  /* stop level creation loop if new level is empty or any limit is reached */
2883  while( nnewlevel > 0 && !SCIPisStopped(scip)
2884  && graph.nlevels < (unsigned int) sepadata->maxnlevels
2885  && ncutsroot < (unsigned int) sepadata->maxcutsroot
2886  && sepadata->ncuts - sepadata->oldncuts < (unsigned int) sepadata->maxsepacutsround);
2887  }
2888 
2889  /* store the last tried root (when running without sorting the variable array, we don't want
2890  * to always check the same variables and therefore start next time where we stopped last time)
2891  */
2892  if( sepadata->sortswitch == UNSORTED )
2893  {
2894  if( i == graph.maxnodes )
2895  sepadata->lastroot = 0;
2896  else
2897  sepadata->lastroot = (int) i;
2898  }
2899 
2900  TERMINATE:
2901  /* free memory */
2902  SCIPfreeBufferArray(scip, &incut);
2903 
2904  SCIPfreeBufferArray(scip, &blocked);
2905  SCIPfreeBufferArray(scip, &distance);
2906  SCIPfreeBufferArray(scip, &parentTreeBackward);
2907  SCIPfreeBufferArray(scip, &parentTree);
2908  SCIPfreeBufferArray(scip, &inQueue);
2909  SCIPfreeBufferArray(scip, &queue);
2910 
2911  SCIPfreeBufferArray(scip, &inlevelgraph);
2912  SCIPfreeBufferArray(scip, &graph.levelAdj);
2913  SCIPfreeBufferArray(scip, &graph.weightAdj);
2914  SCIPfreeBufferArray(scip, &graph.targetAdj);
2915  SCIPfreeBufferArray(scip, &graph.sourceAdj);
2916  SCIPfreeBufferArray(scip, &graph.beginAdj);
2917  SCIPfreeBufferArray(scip, &newlevel);
2918  SCIPfreeBufferArray(scip, &curlevel);
2919 
2920  SCIPfreeBufferArray(scip, &graph.weightBackward);
2921  SCIPfreeBufferArray(scip, &graph.weightForward);
2922  SCIPfreeBufferArray(scip, &graph.targetBackward);
2923  SCIPfreeBufferArray(scip, &graph.targetForward);
2924  SCIPfreeBufferArray(scip, &graph.beginBackward);
2925  SCIPfreeBufferArray(scip, &graph.beginForward);
2926  SCIPfreeBufferArray(scip, &graph.level);
2927 
2928  SCIPfreeBufferArray(scip, &(sepadata->mapping));
2929  SCIPfreeBufferArray(scip, &vars);
2930  SCIPfreeBufferArray(scip, &vals);
2931 
2932  return SCIP_OKAY;
2933 }
2934 
2935 
2936 /* methods for separateGLS() */
2937 
2938 /** memory reallocation method (the graph is normally very dense, so we dynamically allocate only the memory we need) */
2939 static
2940 SCIP_RETCODE checkArraySizesGLS(
2941  SCIP* scip, /**< SCIP data structure */
2942  unsigned int maxarcs, /**< maximal size of graph->head and graph->weight */
2943  unsigned int* arraysize, /**< current size of graph->head and graph->weight */
2944  DIJKSTRA_GRAPH* graph, /**< Dijkstra Graph data structure */
2945  SCIP_Bool* success /**< FALSE, iff memory reallocation fails */
2946  )
2947 {
2948  SCIP_Real memorylimit;
2949  unsigned int additional;
2950  unsigned int j;
2951  unsigned int oldarraysize;
2952 
2953  assert(scip != NULL);
2954  assert(arraysize != NULL);
2955  assert(graph != NULL);
2956  assert(graph->head != NULL);
2957  assert(graph->weight != NULL);
2958  assert(success != NULL);
2959 
2960  SCIPdebugMessage("reallocating graph->head and graph->weight...\n");
2961 
2962  additional = (MIN(maxarcs, 2 * (*arraysize)) - (*arraysize)) * ((int) sizeof(*(graph->head)));
2963  additional += (MIN(maxarcs, 2 * (*arraysize)) - (*arraysize)) * ((int) sizeof(*(graph->weight)));
2964 
2965  SCIP_CALL( SCIPgetRealParam(scip, "limits/memory", &memorylimit) );
2966  if( !SCIPisInfinity(scip, memorylimit) )
2967  {
2968  memorylimit -= SCIPgetMemUsed(scip)/1048576.0;
2969  memorylimit -= SCIPgetMemExternEstim(scip)/1048576.0;
2970  }
2971 
2972 
2973  /* if memorylimit would be exceeded or any other limit is reached free all data and exit */
2974  if( memorylimit <= additional/1048576.0 || SCIPisStopped(scip) )
2975  {
2976  *success = FALSE;
2977  SCIPdebugMessage("...memory limit exceeded\n");
2978  return SCIP_OKAY;
2979  }
2980 
2981  oldarraysize = *arraysize;
2982  *arraysize = 2*(*arraysize);
2983 
2984  SCIP_CALL( SCIPreallocBufferArray(scip, &(graph->head), (int) MIN(maxarcs, (*arraysize))) );
2985  SCIP_CALL( SCIPreallocBufferArray(scip, &(graph->weight), (int) MIN(maxarcs, (*arraysize))) );
2986 
2987  /* if memorylimit exceeded, leave the separator */
2988  SCIP_CALL( SCIPgetRealParam(scip, "limits/memory", &memorylimit) );
2989 
2990  if( !SCIPisInfinity(scip, memorylimit) )
2991  {
2992  memorylimit -= SCIPgetMemUsed(scip)/1048576.0;
2993  memorylimit -= SCIPgetMemExternEstim(scip)/1048576.0;
2994  }
2995 
2996 
2997  if( memorylimit <= 2.0*SCIPgetMemExternEstim(scip)/1048576.0 )
2998  {
2999  SCIPdebugMessage("...memory limit exceeded - freeing all arrays\n");
3000  *success = FALSE;
3001  return SCIP_OKAY;
3002  }
3003 
3004  /* initialize new segments of graph as empty graph */
3005  for( j = oldarraysize; j < MIN(maxarcs,(*arraysize)); ++j )
3006  {
3007  (graph->head)[j] = DIJKSTRA_UNUSED;
3008  (graph->weight)[j] = DIJKSTRA_UNUSED;
3009  }
3010 
3011  SCIPdebugMessage("...with success\n");
3012 
3013  return SCIP_OKAY;
3014 }
3015 
3016 /** add binary implications of the given node */
3017 static
3018 SCIP_RETCODE addGLSBinImpls(
3019  SCIP* scip, /**< SCIP data structure */
3020  SCIP_SEPADATA* sepadata, /**< separator data structure */
3021  SCIP_VAR** vars, /**< problem variables */
3022  unsigned int varsidx, /**< index of current variable inside the problem variables */
3023  unsigned int dijkindex, /**< index of current variable inside the Dijkstra Graph */
3024  SCIP_Real* vals, /**< value of the variables in the given solution */
3025  unsigned int nbinvars, /**< number of binary problem variables */
3026  unsigned int nbinimpls, /**< number of binary implications of the current node */
3027  DIJKSTRA_GRAPH* graph, /**< Dijkstra Graph data structure */
3028  unsigned int* narcs, /**< current number of arcs inside the Dijkstra Graph */
3029  unsigned int maxarcs, /**< maximal number of arcs inside the Dijkstra Graph */
3030  SCIP_Bool original, /**< TRUE, iff variable is a problem variable */
3031  SCIP_Bool* emptygraph, /**< TRUE, iff there is no arc in the implication graph of the binary variables of SCIP */
3032  unsigned int* arraysize, /**< current size of graph->head and graph->weight */
3033  SCIP_Bool* success /**< FALSE, iff memory reallocation fails */
3034  )
3035 {
3036  SCIP_VAR** implvars; /* implications of the current variable */
3037  SCIP_BOUNDTYPE* impltypes; /* type of the implications of the current variable */
3038  unsigned int m;
3039  SCIP_VAR* neighbor; /* current neighbor of the current variable */
3040  unsigned int neighindex;
3041 #ifndef NDEBUG
3042  SCIP_Real* implbounds;
3043 #endif
3044 
3045  assert(scip != NULL);
3046  assert(sepadata != NULL);
3047  assert(vars != NULL);
3048  assert(graph != NULL);
3049  assert(graph->head != NULL);
3050  assert(graph->weight != NULL);
3051  assert(narcs != NULL);
3052  assert(emptygraph != NULL);
3053  assert(arraysize != NULL);
3054  assert(success != NULL);
3055 
3056  /* get implication data */
3057  implvars = SCIPvarGetImplVars(vars[varsidx], original);
3058  impltypes = SCIPvarGetImplTypes(vars[varsidx], original);
3059 #ifndef NDEBUG
3060  implbounds = SCIPvarGetImplBounds(vars[varsidx], original);
3061 #endif
3062 
3063  assert(implvars != NULL || nbinimpls == 0);
3064  assert(impltypes != NULL || nbinimpls == 0);
3065  assert(implbounds != NULL || nbinimpls == 0);
3066 
3067  /* add all implications to the graph */
3068  for( m = 0; m < nbinimpls; ++m )
3069  {
3070  int tmp;
3071 
3072  assert( implvars != NULL && impltypes != NULL && implbounds != NULL ); /* for lint */
3073  assert( SCIPvarIsBinary(implvars[m]) );
3074 
3075  neighbor = implvars[m];
3076  neighindex = sepadata->mapping[SCIPvarGetProbindex(neighbor)];
3077  assert(neighindex < nbinvars);
3078 
3079  /* we use only variables with fractional LP-solution values */
3080  if( SCIPisFeasIntegral(scip, vals[neighindex]) )
3081  continue;
3082 
3083  /* forward direction (the backward is created at the occurrence of the current variable in the cliquevars of the neighbor) */
3084  /* add implication for x==1 */
3085  if( original )
3086  {
3087  /* implication to y=0 (I->III), i.e., x + y <= 1 */
3088  if( impltypes[m] == SCIP_BOUNDTYPE_UPPER )
3089  {
3090  assert( SCIPisEQ(scip, implbounds[m], 0.0) );
3091 
3092  tmp = (int) SCIPfeasCeil(scip, sepadata->scale * (1.0 - vals[varsidx] - vals[neighindex]));
3093  graph->weight[*narcs] = (unsigned int) MAX(0, tmp);
3094  graph->head[*narcs] = neighindex + 2 * nbinvars;
3095  }
3096  /* implication to y=1 (I->IV), i.e., x + (1-y) <= 1 */
3097  else
3098  {
3099  assert( impltypes[m] == SCIP_BOUNDTYPE_LOWER && SCIPisEQ(scip, implbounds[m], 1.0) );
3100 
3101  tmp = (int) SCIPfeasCeil(scip, sepadata->scale * (1.0 - vals[varsidx] - (1.0 - vals[neighindex])) );
3102  graph->weight[*narcs] = (unsigned int) MAX(0, tmp);
3103  graph->head[*narcs] = neighindex + 3 * nbinvars;
3104  }
3105  }
3106  else /* add implication for x==0 */
3107  {
3108  /* implication to y=0 (II->III), i.e., (1-x) + y <= 1 */
3109  if( impltypes[m] == SCIP_BOUNDTYPE_UPPER )
3110  {
3111  assert(implbounds[m] == 0.0);
3112 
3113  tmp = (int) SCIPfeasCeil(scip, sepadata->scale * (1.0 - (1.0 - vals[varsidx]) - vals[neighindex]) );
3114  graph->weight[*narcs] = (unsigned int) MAX(0, tmp);
3115  graph->head[*narcs] = neighindex + 2 * nbinvars;
3116  }
3117  /* implication to y=1 (II->IV), i.e., (1-x) + (1-y) <= 1 */
3118  else
3119  {
3120  assert(impltypes[m] == SCIP_BOUNDTYPE_LOWER && implbounds[m] == 1.0 );
3121 
3122  tmp = (int) SCIPfeasCeil(scip, sepadata->scale * (1.0 - (1.0 - vals[varsidx]) - (1.0 - vals[neighindex])) );
3123  graph->weight[*narcs] = (unsigned int) MAX(0, tmp);
3124  graph->head[*narcs] = neighindex + 3 * nbinvars;
3125  }
3126  }
3127 
3128  /* update minimum and maximum weight values */
3129  if( graph->weight[*narcs] < graph->minweight )
3130  graph->minweight = graph->weight[*narcs];
3131 
3132  if( graph->weight[*narcs] > graph->maxweight )
3133  graph->maxweight = graph->weight[*narcs];
3134 
3135  assert(graph->head[*narcs] >= 2*nbinvars);
3136  assert(graph->head[*narcs] < 4*nbinvars);
3137  ++(*narcs);
3138 
3139  if( *arraysize == *narcs )
3140  {
3141  SCIP_CALL( checkArraySizesGLS(scip, maxarcs, arraysize, graph, success) );
3142 
3143  if( !(*success) )
3144  return SCIP_OKAY;
3145  }
3146  assert((*narcs) < maxarcs);
3147  ++(graph->outcnt[dijkindex]);
3148 
3149  *emptygraph = FALSE;
3150  }
3151 
3152  return SCIP_OKAY;
3153 }
3154 
3155 /** add implications from cliques of the given node */
3156 static
3157 SCIP_RETCODE addGLSCliques(
3158  SCIP* scip, /**< SCIP data structure */
3159  SCIP_SEPADATA* sepadata, /**< separator data structure */
3160  SCIP_VAR** vars, /**< problem variables */
3161  unsigned int varsidx, /**< index of current variable inside the problem variables */
3162  unsigned int dijkindex, /**< index of current variable inside the Dijkstra Graph */
3163  SCIP_Real* vals, /**< value of the variables in the given solution */
3164  unsigned int nbinvars, /**< number of binary problem variables */
3165  unsigned int ncliques, /**< number of cliques of the current node */
3166  DIJKSTRA_GRAPH* graph, /**< Dijkstra Graph data structure */
3167  unsigned int* narcs, /**< current number of arcs inside the Dijkstra Graph */
3168  unsigned int maxarcs, /**< maximal number of arcs inside the Dijkstra Graph */
3169  SCIP_Bool original, /**< TRUE, iff variable is a problem variable */
3170  SCIP_Bool* emptygraph, /**< TRUE, iff there is no arc in the implication graph of the binary variables of SCIP */
3171  unsigned int* arraysize, /**< current size of graph->head and graph->weight */
3172  SCIP_Bool* success /**< FALSE, iff memory reallocation fails */
3173  )
3174 {
3175  SCIP_VAR* neighbor; /* current neighbor of the current variable */
3176  unsigned int neighindex;
3177  SCIP_CLIQUE** cliques; /* cliques of the current variable (with x==0/1) */
3178  unsigned int ncliquevars; /* number of variables of a clique */
3179  SCIP_VAR** cliquevars; /* variables of a clique */
3180  SCIP_Bool* cliquevals; /* is the cliquevariable fixed to TRUE or to FALSE */
3181  unsigned int k;
3182  unsigned int m;
3183 
3184  assert(scip != NULL);
3185  assert(sepadata != NULL);
3186  assert(vars != NULL);
3187  assert(graph != NULL);
3188  assert(graph->head != NULL);
3189  assert(graph->weight != NULL);
3190  assert(narcs != NULL);
3191  assert(emptygraph != NULL);
3192  assert(arraysize != NULL);
3193  assert(success != NULL);
3194 
3195  /* if current variable has cliques of current clique-type */
3196  cliques = SCIPvarGetCliques(vars[varsidx], original);
3197  assert(cliques != NULL || ncliques == 0);
3198 
3199  for( k = 0; k < ncliques; ++k )
3200  {
3201  assert( cliques != NULL ); /* for lint */
3202 
3203  /* get clique data */
3204  cliquevars = SCIPcliqueGetVars(cliques[k]);
3205  ncliquevars = (unsigned int) SCIPcliqueGetNVars(cliques[k]);
3206  cliquevals = SCIPcliqueGetValues(cliques[k]);
3207 
3208  assert(cliquevars != NULL || ncliquevars == 0);
3209  assert(cliquevals != NULL || ncliquevars == 0);
3210 
3211  /* add arcs for all fractional variables in clique */
3212  for( m = 0; m < ncliquevars; ++m )
3213  {
3214  int tmp;
3215 
3216  assert( cliquevars != NULL && cliquevals != NULL ); /* for lint */
3217  neighbor = cliquevars[m];
3218 
3219  neighindex = sepadata->mapping[SCIPvarGetProbindex(neighbor)];
3220  assert(neighindex < nbinvars);
3221 
3222  /* ignore current variable */
3223  if( neighindex == varsidx )
3224  continue;
3225 
3226  /* we use only variables with fractional LP-solution values */
3227  if( SCIPisFeasIntegral(scip, vals[neighindex]) )
3228  continue;
3229 
3230  /* forward direction (the backward is created at the occurrence of the current variable in the cliquevars of the neighbor) */
3231  /* x==1 */
3232  if( original )
3233  {
3234  /* implication to y=0 (I->III) */
3235  if( cliquevals[m] )
3236  {
3237  tmp = (int) SCIPfeasCeil(scip, sepadata->scale * (1 - vals[varsidx] - vals[neighindex] ));
3238  graph->weight[*narcs] = (unsigned int) MAX(0, tmp);
3239  graph->head[*narcs] = neighindex + 2 * nbinvars;
3240  }
3241  /* implication to y=1 (I->IV) (cliquevals[m] == FALSE) */
3242  else
3243  {
3244  tmp = (int) SCIPfeasCeil(scip, sepadata->scale * (1 - vals[varsidx] - (1 - vals[neighindex]) ));
3245  graph->weight[*narcs] = (unsigned int) MAX(0, tmp);
3246  graph->head[*narcs] = neighindex + 3 * nbinvars;
3247  }
3248  }
3249  /* x==0 */
3250  else
3251  {
3252  /* implication to y=0 (II->III) */
3253  if( cliquevals[m] )
3254  {
3255  tmp = (int) SCIPfeasCeil(scip, sepadata->scale * (1 - (1 - vals[varsidx]) - vals[neighindex] ));
3256  graph->weight[*narcs] = (unsigned int) MAX(0, tmp);
3257  graph->head[*narcs] = neighindex + 2 * nbinvars;
3258  }
3259  /* implication to y=1 (II->IV) (cliquevals[m] == FALSE) */
3260  else
3261  {
3262  tmp = (int) SCIPfeasCeil(scip, sepadata->scale * (1 - (1 - vals[varsidx]) - (1-vals[neighindex]) )) ;
3263  graph->weight[*narcs] = (unsigned int) MAX(0, tmp);
3264  graph->head[*narcs] = neighindex + 3 * nbinvars;
3265  }
3266  }
3267 
3268  /* update minimum and maximum weight values */
3269  if( graph->weight[*narcs] < graph->minweight )
3270  graph->minweight = graph->weight[*narcs];
3271 
3272  if( graph->weight[*narcs] > graph->maxweight )
3273  graph->maxweight = graph->weight[*narcs];
3274 
3275  ++(*narcs);
3276  if( *arraysize == *narcs )
3277  {
3278  SCIP_CALL( checkArraySizesGLS(scip, maxarcs, arraysize, graph, success) );
3279 
3280  if( !(*success) )
3281  return SCIP_OKAY;
3282  }
3283  assert((*narcs) < maxarcs);
3284  ++(graph->outcnt[dijkindex]);
3285 
3286  *emptygraph = FALSE;
3287  }
3288  }
3289 
3290  return SCIP_OKAY;
3291 }
3292 
3293 /** The classical method for finding odd cycles by Groetschel, Lovasz, Schrijver uses a bipartite graph
3294  * which contains in each partition a node for every node in the original graph.
3295  * All arcs uv of the original graph are copied to arcs from u of the first partition to v' of the second partition
3296  * and from u' of the second partition to v of the first partition.
3297  * A Dijkstra algorithm is used to find a path from a node x to its copy x', if existing.
3298  * The nodes in the original graph corresponding to the nodes on the path form an odd cycle.
3299  *
3300  * Since SCIP stores implications between original and negated variables,
3301  * our original graph has at most twice the number of binary variables of the problem.
3302  * By creating the bipartite graph we gain 4 segments of the graph:
3303  *
3304  * I - nodes of the original variables in the first bipartition \n
3305  * II - nodes of the negated variables in the first bipartition \n
3306  * III - nodes of the original variables in the second bipartition \n
3307  * IV - nodes of the negated variables in the second bipartition
3308  *
3309  * The length of every segment if the number of binary variables in the original problem.
3310  *
3311  * Since the implication graph of SCIP is (normally) incomplete,
3312  * it is possible to use arcs between an original variable and its negated
3313  * to obtain more cycles which are valid but not found due to missing links.
3314  */
3315 static
3316 SCIP_RETCODE separateGLS(
3317  SCIP* scip, /**< SCIP data structure */
3318  SCIP_SEPA* sepa, /**< separator */
3319  SCIP_SEPADATA* sepadata, /**< separator data structure */
3320  SCIP_SOL* sol, /**< given primal solution */
3321  SCIP_RESULT* result /**< pointer to store the result of the separation call */
3322  )
3323 {
3324  SCIP_Bool emptygraph; /* flag if graph contains an arc */
3325 
3326  SCIP_Real* vals; /* values of the variables in the given solution */
3327  unsigned int* incut;
3328 
3329  unsigned int i;
3330  unsigned int j;
3331 
3332  SCIP_VAR** scipvars; /* variables of the current SCIP (unsorted) */
3333  SCIP_VAR** vars; /* variables of the current SCIP (sorted if requested) */
3334  unsigned int nbinvars; /* number of binary problem variables */
3335  SCIP_Bool original; /* flag if the current variable is original or negated */
3336 
3337  unsigned int nbinimpls; /* number of binary implications of the current variable */
3338  unsigned int ncliques; /* number of cliques of the current variable */
3339 
3340  DIJKSTRA_GRAPH graph; /* Dijkstra graph data structure */
3341  unsigned int arraysize; /* current size of graph->head and graph->weight */
3342  unsigned int narcs; /* number of arcs in the Dijkstra graph */
3343  unsigned int maxarcs; /* maximum number of arcs in the Dijkstra graph */
3344  unsigned int maxstarts; /* maximum number of start nodes */
3345  unsigned int startcounter; /* counter of tried start nodes */
3346  unsigned long long cutoff; /* cutoff value for Dijkstra algorithm */
3347 
3348  unsigned int startnode; /* start node for Dijkstra algorithm */
3349  unsigned int endnode; /* target node for Dijkstra algorithm */
3350  unsigned long long* dist; /* distance matrix for Dijkstra algorithm */
3351  unsigned int* pred; /* predecessor list for found cycle */
3352  unsigned int* entry; /* storage for Dijkstra algorithm */
3353  unsigned int* order; /* storage for Dijkstra algorithm */
3354  unsigned int dijkindex;
3355  SCIP_Bool success; /* flag for check for several errors */
3356 
3357  SCIP_Bool* incycle; /* flag array if variable is contained in the found cycle */
3358  unsigned int* pred2; /* temporary predecessor list for backprojection of found cycle */
3359 
3360  int nscipbinvars;
3361  int nscipintvars;
3362  int nscipimplvars;
3363  int nintegral;
3364  int k;
3365 
3366  assert(scip != NULL);
3367  assert(sepadata != NULL);
3368  assert(result != NULL);
3369 
3370  success = TRUE;
3371  emptygraph = TRUE;
3372 
3373  SCIP_CALL( SCIPgetVarsData(scip, &scipvars, NULL, &nscipbinvars, &nscipintvars, &nscipimplvars, NULL) );
3374  assert(nscipbinvars >= 0);
3375  assert(nscipintvars >= 0);
3376  assert(nscipimplvars >= 0);
3377 
3378  nintegral = nscipbinvars + nscipintvars + nscipimplvars;
3379  assert(scipvars != NULL || nintegral == 0);
3380 
3381  /* collect binary variables, including implicit binary */
3382  SCIP_CALL( SCIPallocBufferArray(scip, &vars, nintegral) );
3383  for (k = 0; k < nscipbinvars; ++k)
3384  vars[k] = scipvars[k]; /*lint !e613*/
3385 
3386  nbinvars = (unsigned int) nscipbinvars;
3387  for (k = nscipbinvars; k < nintegral; ++k)
3388  {
3389  assert( SCIPvarGetType(scipvars[k]) != SCIP_VARTYPE_CONTINUOUS ); /*lint !e613*/
3390  if ( SCIPvarIsBinary(scipvars[k]) ) /*lint !e613*/
3391  vars[nbinvars++] = scipvars[k]; /*lint !e613*/
3392  }
3393 
3394  if( nbinvars == 0 )
3395  {
3396  SCIPfreeBufferArray(scip, &vars);
3397  return SCIP_OKAY;
3398  }
3399 
3400  /* initialize flag array to avoid multiple cuts per variable, if requested by user-flag */
3401  SCIP_CALL( SCIPallocBufferArray(scip, &vals, (int) (2 * nbinvars)) );
3402 
3403  /* prepare values */
3404  assert( vars != NULL );
3405  switch( sepadata->sortswitch )
3406  {
3407  case UNSORTED :
3408  /* if no sorting is requested, we use the normal variable array */
3409  break;
3410 
3411  case MAXIMAL_LPVALUE :
3412  /* store lp-values */
3413  for( i = 0; i < nbinvars; ++i )
3414  vals[i] = SCIPgetSolVal(scip, sol, vars[i]);
3415 
3416  /* sort by lp-value, maximal first */
3417  SCIPsortDownRealPtr(vals, (void**)vars, (int) nbinvars);
3418  break;
3419 
3420  case MINIMAL_LPVALUE :
3421  /* store lp-values */
3422  for( i = 0; i < nbinvars; ++i )
3423  vals[i] = SCIPgetSolVal(scip, sol, vars[i]);
3424 
3425  /* sort by lp-value, minimal first */
3426  SCIPsortRealPtr(vals, (void**)vars, (int) nbinvars);
3427  break;
3428 
3429  case MAXIMAL_FRACTIONALITY :
3430  /* store lp-values and determine fractionality */
3431  for( i = 0; i < nbinvars; ++i )
3432  {
3433  vals[i] = SCIPgetSolVal(scip, sol, vars[i]);
3434  vals[i] = MIN(1.0 - vals[i], vals[i]);
3435  }
3436 
3437  /* sort by fractionality, maximal first */
3438  SCIPsortDownRealPtr(vals, (void**)vars, (int) nbinvars);
3439  break;
3440 
3441  case MINIMAL_FRACTIONALITY :
3442  /* store lp-values and determine fractionality */
3443  for( i = 0; i < nbinvars; ++i )
3444  {
3445  vals[i] = SCIPgetSolVal(scip, sol, vars[i]);
3446  vals[i] = MIN(1.0 - vals[i], vals[i]);
3447  }
3448 
3449  /* sort by fractionality, minimal first */
3450  SCIPsortRealPtr(vals, (void**)vars, (int) nbinvars);
3451  break;
3452 
3453  default :
3454  SCIPerrorMessage("invalid sortswitch value\n");
3455  SCIPABORT();
3456  return SCIP_INVALIDDATA; /*lint !e527*/
3457  }
3458  assert(vars != NULL);
3459 
3460  /* create mapping for getting the index of a variable via its probindex to the index in the sorted variable array */
3461  SCIP_CALL( SCIPallocBufferArray(scip, &(sepadata->mapping), nintegral) );
3462  SCIP_CALL( SCIPallocBufferArray(scip, &incut, (int) (4 * nbinvars)) );
3463  BMSclearMemoryArray(incut, 4 * nbinvars);
3464 
3465  /* initialize LP value and cut flag for all variables */
3466  for( i = 0; i < nbinvars; ++i )
3467  {
3468  assert( 0 <= SCIPvarGetProbindex(vars[i]) && SCIPvarGetProbindex(vars[i]) < nintegral); /* since binary, integer, and implicit variables are first */
3469  sepadata->mapping[SCIPvarGetProbindex(vars[i])] = i;
3470  vals[i] = SCIPgetSolVal(scip, sol, vars[i]); /* need to get new values, since they might be corrupted */
3471  }
3472 
3473  for( i = nbinvars; i < 2*nbinvars; ++i )
3474  vals[i] = 1 - vals[i - nbinvars];
3475 
3476  /* initialize number of nodes in Dijkstra graph (2*2*n nodes in a mirrored bipartite graph with negated variables) */
3477  graph.nodes = 4 * nbinvars;
3478 
3479  /* Initialize number of arcs in Dijkstra graph, should be (nbinvars+1) * graph.nodes, but might deviate, because
3480  * there might be parallel arcs:
3481  * nbinvars-1 possible arcs per node (it is not possible to be linked to variable and negated)
3482  * + 1 self-arc (arc to negated variable)
3483  * + 1 dummy arc for Dijkstra data structure
3484  * = nbinvars+1 arcs per node
3485  * * graph.nodes
3486  * = (nbinvars+1)*graph.nodes
3487  * + graph.nodes => separating entries for arclist)
3488  *
3489  * Number is corrected below.
3490  */
3491  graph.arcs = 0;
3492 
3493  /* the implication graph is redundant and therefore more implications and clique arcs may occur than should be possible
3494  * @todo later: filtering of edges which were already added, maxarcs should be graph.arcs rather than INT_MAX;
3495  */
3496  maxarcs = UINT_MAX;
3497 
3498  /* allocate memory for Dijkstra graph arrays */
3499  arraysize = 100 * graph.nodes;
3500  SCIP_CALL( SCIPallocBufferArray(scip, &graph.outbeg, (int) graph.nodes) );
3501  SCIP_CALL( SCIPallocBufferArray(scip, &graph.outcnt, (int) graph.nodes) );
3502  SCIP_CALL( SCIPallocBufferArray(scip, &graph.head, (int) MIN(maxarcs, arraysize)) );
3503  SCIP_CALL( SCIPallocBufferArray(scip, &graph.weight, (int) MIN(maxarcs, arraysize)) );
3504  SCIP_CALL( SCIPallocBufferArray(scip, &dist, (int) graph.nodes) );
3505  SCIP_CALL( SCIPallocBufferArray(scip, &pred, (int) graph.nodes) );
3506  SCIP_CALL( SCIPallocBufferArray(scip, &entry, (int) graph.nodes) );
3507  SCIP_CALL( SCIPallocBufferArray(scip, &order, (int) graph.nodes) );
3508 
3509  /* initialize Dijkstra graph as empty graph */
3510  for( i = 0; i < MIN(arraysize, maxarcs); ++i )
3511  {
3512  graph.head[i] = DIJKSTRA_UNUSED;
3513  graph.weight[i] = DIJKSTRA_UNUSED;
3514  }
3515  graph.minweight = DIJKSTRA_FARAWAY;
3516  graph.maxweight = 0;
3517  narcs = 0;
3518 
3519 #ifndef NDEBUG
3520  for( i = 0; i < graph.nodes; ++i )
3521  {
3522  graph.outbeg[i] = 0;
3523  graph.outcnt[i] = 0;
3524  }
3525 #endif
3526 
3527  /* add arcs from first to second partition to Dijkstra graph (based on the original fractional implication graph) */
3528  for( dijkindex = 0; dijkindex < 2 * nbinvars; ++dijkindex )
3529  {
3530  graph.outbeg[dijkindex] = narcs;
3531  graph.outcnt[dijkindex] = 0;
3532 
3533  /* decide if we have original or negated variable */
3534  if( dijkindex < nbinvars )
3535  {
3536  i = dijkindex;
3537  original = TRUE;
3538  }
3539  else
3540  {
3541  i = dijkindex - nbinvars;
3542  original = FALSE;
3543  }
3544  assert(i < nbinvars);
3545 
3546  /* if the variable has a fractional value we add it to the graph */
3547  if( ! SCIPisFeasIntegral(scip, vals[i]) )
3548  {
3549  nbinimpls = (unsigned int) SCIPvarGetNBinImpls(vars[i], original);
3550  ncliques = (unsigned int) SCIPvarGetNCliques(vars[i], original);
3551 
3552  /* insert arcs for binary implications (take var => getImpl(Bin) => getImplVar => add forward-arc) */
3553  /* add implications of implication-type "original" if current variable has them */
3554  if( nbinimpls >= 1 )
3555  {
3556  /* implications from x = 1/0 to y = 0/1 (I/II -> III/IV) */
3557  SCIP_CALL( addGLSBinImpls(scip, sepadata, vars, i, dijkindex, vals, nbinvars, nbinimpls, &graph,
3558  &narcs, maxarcs, original, &emptygraph, &arraysize, &success) );
3559 
3560  if( !success )
3561  goto TERMINATE;
3562  }
3563 
3564  /* insert arcs for cliques (take var => getCliques => take cliquevar => add forward-arc) */
3565  /* add clique arcs of clique-type "original" if current variable has them */
3566  if( ncliques >= 1 )
3567  {
3568  /* x==1/0 -> y==0/1 (I/II -> III/IV) */
3569  SCIP_CALL( addGLSCliques(scip, sepadata, vars, i, dijkindex, vals, nbinvars, ncliques, &graph,
3570  &narcs, maxarcs, original, &emptygraph, &arraysize, &success) );
3571 
3572  if( !success )
3573  goto TERMINATE;
3574  }
3575  }
3576 
3577  /* add link to copy of negated variable (useful if/because the implication graph is incomplete) */
3578  if( sepadata->addselfarcs && graph.outcnt[dijkindex] > 0 )
3579  {
3580  /* I -> IV */
3581  if( original )
3582  {
3583  assert(dijkindex < nbinvars);
3584  graph.head[narcs] = dijkindex + 3*nbinvars;
3585  }
3586  /* II -> III */
3587  else
3588  {
3589  assert(dijkindex >= nbinvars && dijkindex < 2*nbinvars);
3590  graph.head[narcs] = dijkindex + nbinvars;
3591  }
3592  graph.weight[narcs] = 0;
3593 
3594  /* update minimum and maximum weight values */
3595  if( graph.weight[narcs] < graph.minweight )
3596  graph.minweight = graph.weight[narcs];
3597 
3598  if( graph.weight[narcs] > graph.maxweight )
3599  graph.maxweight = graph.weight[narcs];
3600 
3601  ++narcs;
3602  if( arraysize == narcs )
3603  {
3604  SCIP_CALL( checkArraySizesGLS(scip, maxarcs, &arraysize, &graph, &success) );
3605  if( !success )
3606  goto TERMINATE;
3607  }
3608  assert(narcs < maxarcs);
3609  ++(graph.outcnt[dijkindex]);
3610  }
3611 
3612  /* add separating arc */
3613  graph.head[narcs] = DIJKSTRA_UNUSED;
3614  graph.weight[narcs] = DIJKSTRA_UNUSED;
3615  ++narcs;
3616  if( arraysize == narcs )
3617  {
3618  SCIP_CALL( checkArraySizesGLS(scip, maxarcs, &arraysize, &graph, &success) );
3619  if( !success )
3620  goto TERMINATE;
3621  }
3622  assert(narcs < maxarcs);
3623  }
3624 
3625  /* if the graph is empty, there is nothing to do */
3626  if( emptygraph )
3627  goto TERMINATE;
3628 
3629  /* add arcs from second to first partition to Dijkstra graph */
3630  for( i = 0; i < 2*nbinvars; ++i )
3631  {
3632  graph.outbeg[2 * nbinvars + i] = narcs;
3633  graph.outcnt[2 * nbinvars + i] = 0;
3634 
3635  /* copy all arcs to head from the second to the first bipartition */
3636  for( j = graph.outbeg[i]; j < graph.outbeg[i] + graph.outcnt[i]; ++j )
3637  {
3638  /* there are only arcs from first bipartition to the second */
3639  assert(graph.head[j] >= 2*nbinvars && graph.head[j] < 4*nbinvars);
3640 
3641  /* the backward arcs head from III->I or IV->II */
3642  graph.head[narcs] = graph.head[j] - 2 * nbinvars;
3643  graph.weight[narcs] = graph.weight[j];
3644  ++narcs;
3645  if( arraysize == narcs )
3646  {
3647  SCIP_CALL( checkArraySizesGLS(scip, maxarcs, &arraysize, &graph, &success) );
3648 
3649  if( !success )
3650  goto TERMINATE;
3651  }
3652  assert(narcs < maxarcs);
3653  ++(graph.outcnt[2*nbinvars+i]);
3654  }
3655 
3656  /* add separating arc */
3657  graph.head[narcs] = DIJKSTRA_UNUSED;
3658  graph.weight[narcs] = DIJKSTRA_UNUSED;
3659  ++narcs;
3660 
3661  if( arraysize == narcs )
3662  {
3663  SCIP_CALL( checkArraySizesGLS(scip, maxarcs, &arraysize, &graph, &success) );
3664 
3665  if( !success )
3666  goto TERMINATE;
3667  }
3668  assert(narcs < maxarcs);
3669  }
3670 
3671  /* correct number of arcs */
3672  graph.arcs = narcs;
3673 
3674  SCIPdebugMessage("--- graph successfully created (%u nodes, %u arcs) ---\n", graph.nodes, narcs);
3675 
3676  /* graph is now prepared for Dijkstra methods */
3677  assert( dijkstraGraphIsValid(&graph) );
3678 
3679 #ifdef SCIP_ODDCYCLE_WRITEGRAPH
3680  {
3681  char probname [SCIP_MAXSTRLEN];
3682  char filename [SCIP_MAXSTRLEN];
3683  char* name;
3684 
3685  (void) SCIPsnprintf(probname, SCIP_MAXSTRLEN, "%s", SCIPgetProbName(scip));
3686  SCIPsplitFilename(probname, NULL, &name, NULL, NULL);
3687  (void) SCIPsnprintf(filename, SCIP_MAXSTRLEN, "%s_%d.gml", name, SCIPgetNLPs(scip));
3688  SCIP_CALL( SCIPwriteCliqueGraph(scip, filename, TRUE, TRUE) );
3689  SCIPverbMessage(scip, SCIP_VERBLEVEL_HIGH, NULL, "Wrote clique/implication graph to <%s>.\n", filename);
3690  }
3691 #endif
3692 
3693  /* determine the number of start nodes */
3694  maxstarts = (unsigned int) SCIPceil(scip, sepadata->offsettestvars + (0.02 * nbinvars * sepadata->percenttestvars));
3695  startcounter = 0;
3696 
3697  /* allocate and initialize predecessor list and flag array representing odd cycle */
3698  SCIP_CALL( SCIPallocBufferArray(scip, &pred2, (int) (2 * nbinvars)) );
3699  SCIP_CALL( SCIPallocBufferArray(scip, &incycle, (int) (2 * nbinvars)) );
3700 
3701  /* separate odd cycle inequalities by GLS method */
3702  cutoff = (unsigned long long) (0.5 * sepadata->scale);
3703  for( i = (unsigned int) sepadata->lastroot; i < 2 * nbinvars
3704  && startcounter < maxstarts
3705  && sepadata->ncuts - sepadata->oldncuts < (unsigned int) sepadata->maxsepacutsround
3706  && !SCIPisStopped(scip); ++i )
3707  {
3708  unsigned int ncyclevars; /* cycle length */
3709  SCIP_Bool edgedirection; /* partitionindicator for backprojection from bipartite graph to original graph:
3710  * is the current edge a backwards edge, i.e., from second to first partition? */
3711 
3712  /* skip isolated node */
3713  if( graph.head[graph.outbeg[i]] == DIJKSTRA_UNUSED )
3714  continue;
3715 
3716  /* if node has only one arc, there is no odd cycle containing this node
3717  * (but there are invalid odd circuits containing the only neighbor twice)
3718  */
3719  if( graph.head[graph.outbeg[i]+1] == DIJKSTRA_UNUSED )
3720  continue;
3721 
3722  /* search shortest path from node to its counter part in the other partition */
3723  startnode = i;
3724  endnode = i + 2 * nbinvars;
3725 
3726  /* skip node if it is already covered by a cut and
3727  * we do not want to search cycles starting with a node already covered by a cut
3728  */
3729  if( incut[startnode] && ! sepadata->multiplecuts )
3730  continue;
3731 
3732  ++startcounter;
3733 
3734  if ( sepadata->allowmultiplecuts )
3735  (void) dijkstraPairCutoffIgnore(&graph, startnode, endnode, incut, cutoff, dist, pred, entry, order);
3736  else
3737  (void) dijkstraPairCutoff(&graph, startnode, endnode, cutoff, dist, pred, entry, order);
3738 
3739  /* no odd cycle cut found */
3740  if( dist[endnode] == DIJKSTRA_FARAWAY )
3741  continue;
3742 
3743  /* skip check if cutoff has been exceeded */
3744  if ( dist[endnode] >= cutoff )
3745  continue;
3746 
3747  /* detect cycle including:
3748  * project bipartitioned graph to original graph of variables and their negated
3749  * (pred&incycle-structure for generateOddCycleCut)
3750  * check cycles for double variables and try to clean variable-negated-sub-cycles if existing
3751  */
3752  for( j = 0; j < 2 * nbinvars; ++j )
3753  {
3754  pred2[j] = DIJKSTRA_UNUSED;
3755  incycle[j] = FALSE;
3756  }
3757 
3758  ncyclevars = 0;
3759  edgedirection = TRUE;
3760  success = TRUE;
3761 
3762  /* construct odd cycle in implication graph from shortest path on bipartite graph */
3763  for( dijkindex = endnode; dijkindex != startnode && success; dijkindex = pred[dijkindex], edgedirection = !edgedirection )
3764  {
3765  if( edgedirection )
3766  {
3767  /* check that current node is in second partition and next node is in first partition */
3768  assert(dijkindex >= 2 * nbinvars && dijkindex < 4 * nbinvars);
3769  assert(pred[dijkindex] < 2*nbinvars);
3770 
3771  pred2[dijkindex - 2 * nbinvars] = pred[dijkindex];
3772 
3773  /* check whether the object found is really a cycle without sub-cycles
3774  * (sub-cycles may occur in case there is not violated odd cycle inequality)
3775  * and remove pairs of original and negated variable from cycle
3776  */
3777  SCIP_CALL( cleanCycle(scip, pred2, incycle, incut, dijkindex-2*nbinvars, endnode-2*nbinvars, nbinvars,
3778  &ncyclevars, sepadata->repaircycles, sepadata->allowmultiplecuts, &success) );
3779  }
3780  else
3781  {
3782  /* check that current node is in first partition and next node is in second partition */
3783  assert(dijkindex < 2 * nbinvars);
3784  assert(pred[dijkindex] >= 2 * nbinvars && pred[dijkindex] < 4 * nbinvars);
3785 
3786  pred2[dijkindex] = pred[dijkindex] - 2 * nbinvars;
3787 
3788  /* check whether the object found is really a cycle without sub-cycles
3789  * (sub-cycles may occur in case there is not violated odd cycle inequality)
3790  * and remove pairs of original and negated variable from cycle
3791  */
3792  SCIP_CALL( cleanCycle(scip, pred2, incycle, incut, dijkindex, endnode-2*nbinvars, nbinvars, &ncyclevars,
3793  sepadata->repaircycles, sepadata->allowmultiplecuts, &success) );
3794  }
3795  }
3796 
3797  if( success )
3798  {
3799  /* generate cut */
3800  sepadata->dijkstragraph = &graph;
3801  SCIP_CALL( generateOddCycleCut(scip, sepa, sol, vars, nbinvars, startnode, pred2, ncyclevars, incut, vals, sepadata, result) );
3802 #ifndef NDEBUG
3803  sepadata->dijkstragraph = NULL;
3804 #endif
3805  }
3806  }
3807 
3808  /* free temporary memory */
3809  SCIPfreeBufferArray(scip, &incycle);
3810  SCIPfreeBufferArray(scip, &pred2);
3811 
3812  /* store the last tried root (when running without sorting the variable array, we don't want
3813  * to always check the same variables and therefore start next time where we stopped last time)
3814  */
3815  if( sepadata->sortswitch == UNSORTED )
3816  {
3817  if( i == 2 * nbinvars )
3818  sepadata->lastroot = 0;
3819  else
3820  sepadata->lastroot = (int) i;
3821  }
3822 
3823  TERMINATE:
3824  /* free temporary memory */
3825  SCIPfreeBufferArray(scip, &order);
3826  SCIPfreeBufferArray(scip, &entry);
3827  SCIPfreeBufferArray(scip, &pred);
3828  SCIPfreeBufferArray(scip, &dist);
3829  SCIPfreeBufferArray(scip, &graph.weight);
3830  SCIPfreeBufferArray(scip, &graph.head);
3831  SCIPfreeBufferArray(scip, &graph.outcnt);
3832  SCIPfreeBufferArray(scip, &graph.outbeg);
3833  SCIPfreeBufferArray(scip, &incut);
3834  SCIPfreeBufferArray(scip, &(sepadata->mapping));
3835  SCIPfreeBufferArray(scip, &vars);
3836  SCIPfreeBufferArray(scip, &vals);
3837 
3838  return SCIP_OKAY;
3839 }
3840 
3841 
3842 /*
3843  * Callback methods of separator
3844  */
3845 
3846 /** copy method for separator plugins (called when SCIP copies plugins) */
3847 static
3848 SCIP_DECL_SEPACOPY(sepaCopyOddcycle)
3849 { /*lint --e{715}*/
3850  assert(scip != NULL);
3851  assert(sepa != NULL);
3852  assert(strcmp(SCIPsepaGetName(sepa), SEPA_NAME) == 0);
3853 
3854  /* call inclusion method of constraint handler */
3855  SCIP_CALL( SCIPincludeSepaOddcycle(scip) );
3856 
3857  return SCIP_OKAY;
3858 }
3859 
3860 
3861 /** destructor of separator to free user data (called when SCIP is exiting) */
3862 static
3863 SCIP_DECL_SEPAFREE(sepaFreeOddcycle)
3864 {
3865  SCIP_SEPADATA* sepadata;
3866 
3867  sepadata = SCIPsepaGetData(sepa);
3868  assert(sepadata != NULL);
3869 
3870  SCIPfreeMemory(scip, &sepadata);
3871  SCIPsepaSetData(sepa, NULL);
3872 
3873  return SCIP_OKAY;
3874 }
3875 
3876 
3877 /** initialization method of separator (called after problem was transformed) */
3878 static
3879 SCIP_DECL_SEPAINIT(sepaInitOddcycle)
3880 {
3881  SCIP_SEPADATA* sepadata;
3882 
3883  sepadata = SCIPsepaGetData(sepa);
3884  assert(sepadata != NULL);
3885 
3886  sepadata->ncuts = 0;
3887  sepadata->oldncuts = 0;
3888  sepadata->nliftedcuts = 0;
3889  sepadata->lastroot = 0;
3890  sepadata->levelgraph = NULL;
3891  sepadata->dijkstragraph = NULL;
3892 
3893  return SCIP_OKAY;
3894 }
3895 
3896 
3897 /** solving process initialization method of separator (called when branch and bound process is about to begin) */
3898 static
3899 SCIP_DECL_SEPAINITSOL(sepaInitsolOddcycle)
3900 {
3901  SCIP_SEPADATA* sepadata;
3902 
3903  assert(sepa != NULL);
3904 
3905  sepadata = SCIPsepaGetData(sepa);
3906  assert(sepadata != NULL);
3907 
3908  sepadata->nunsucessfull = 0;
3909  sepadata->lastnode = -1;
3910 
3911  return SCIP_OKAY;
3912 }
3913 
3914 
3915 /** LP solution separation method of separator */
3916 static
3917 SCIP_DECL_SEPAEXECLP(sepaExeclpOddcycle)
3918 {
3919  SCIP_SEPADATA* sepadata;
3920  int depth;
3921  int ncalls;
3922  int oldnliftedcuts;
3923 
3924  *result = SCIP_DIDNOTRUN;
3925 
3926  /* get separator data */
3927  sepadata = SCIPsepaGetData(sepa);
3928  assert(sepadata != NULL);
3929 
3930  depth = SCIPgetDepth(scip);
3931  ncalls = SCIPsepaGetNCallsAtNode(sepa);
3932 
3933  /* only call separator a given number of rounds at each b&b node */
3934  if( (depth == 0 && sepadata->maxroundsroot >= 0 && ncalls >= sepadata->maxroundsroot)
3935  || (depth > 0 && sepadata->maxrounds >= 0 && ncalls >= sepadata->maxrounds) )
3936  return SCIP_OKAY;
3937 
3938  /* only call separator if enough binary variables are present */
3939  if( SCIPgetNBinVars(scip) < 3 || (!(sepadata->includetriangles) && SCIPgetNBinVars(scip) < 5))
3940  {
3941  SCIPdebugMessage("skipping separator: not enough binary variables\n");
3942  return SCIP_OKAY;
3943  }
3944 
3945  /* only call separator if enough fractional variables are present */
3946  if( SCIPgetNLPBranchCands(scip) < 3 || (!(sepadata->includetriangles) && SCIPgetNLPBranchCands(scip) < 5))
3947  {
3948  SCIPdebugMessage("skipping separator: not enough fractional variables\n");
3949  return SCIP_OKAY;
3950  }
3951 
3952  /* only call separator if enough implications (but not all implications should come from cliques) are present */
3953  if( SCIPgetNImplications(scip) < 1 || SCIPgetNImplications(scip) + SCIPgetNCliques(scip) < 3 )
3954  {
3955  SCIPdebugMessage("skipping separator: not enough implications present\n");
3956  return SCIP_OKAY;
3957  }
3958 
3959  /* only run if number of cuts already found is small enough */
3960  if ( sepadata->cutthreshold >= 0 && SCIPgetNCutsFoundRound(scip) >= sepadata->cutthreshold )
3961  return SCIP_OKAY;
3962 
3963  /* store node number and reset number of unsuccessful calls */
3964  if ( sepadata->lastnode != SCIPnodeGetNumber(SCIPgetCurrentNode(scip)) )
3965  {
3966  sepadata->nunsucessfull = 0;
3967  sepadata->lastnode = SCIPnodeGetNumber(SCIPgetCurrentNode(scip));
3968  }
3969  else
3970  {
3971  if ( sepadata->nunsucessfull > sepadata->maxunsucessfull )
3972  {
3973  SCIPdebugMessage("skipping separator: number of unsucessfull calls = %d.\n", sepadata->nunsucessfull);
3974  return SCIP_OKAY;
3975  }
3976  }
3977 
3978  *result = SCIP_DIDNOTFIND;
3979  sepadata->oldncuts = sepadata->ncuts;
3980  oldnliftedcuts = sepadata->nliftedcuts;
3981 
3982  if( depth == 0 )
3983  sepadata->maxsepacutsround = sepadata->maxsepacutsroot;
3984  else
3985  sepadata->maxsepacutsround = sepadata->maxsepacuts;
3986 
3987  /* perform the actual separation routines */
3988  if( sepadata->usegls )
3989  {
3990  SCIPdebugMessage("using GLS method for finding odd cycles\n");
3991  SCIP_CALL( separateGLS(scip, sepa, sepadata, NULL, result) );
3992  }
3993  else
3994  {
3995  SCIPdebugMessage("using level graph heuristic for finding odd cycles\n");
3996  SCIP_CALL( separateHeur(scip, sepa, sepadata, NULL, result) );
3997  }
3998 
3999  if( sepadata->ncuts - sepadata->oldncuts > 0 )
4000  {
4001  SCIPdebugMessage("added %u cuts (%d allowed), %d lifted.\n", sepadata->ncuts - sepadata->oldncuts,
4002  sepadata->maxsepacutsround, sepadata->nliftedcuts - oldnliftedcuts);
4003  sepadata->nunsucessfull = 0;
4004  }
4005  else
4006  {
4007  SCIPdebugMessage("no cuts added (%d allowed)\n", sepadata->maxsepacutsround);
4008  ++sepadata->nunsucessfull;
4009  }
4010  SCIPdebugMessage("total sepatime: %.2f - total number of added cuts: %u\n", SCIPsepaGetTime(sepa), sepadata->ncuts);
4011 
4012  return SCIP_OKAY;
4013 }
4014 
4015 
4016 /*
4017  * separator specific interface methods
4018  */
4019 
4020 /** creates the oddcycle separator and includes it in SCIP */
4021 SCIP_RETCODE SCIPincludeSepaOddcycle(
4022  SCIP* scip /**< SCIP data structure */
4023  )
4024 {
4025  SCIP_SEPADATA* sepadata;
4026  SCIP_SEPA* sepa;
4027 
4028  /* create oddcycle separator data */
4029  SCIP_CALL( SCIPallocMemory(scip, &sepadata) );
4030  sepadata->nunsucessfull = 0;
4031  sepadata->lastnode = -1;
4032 
4033  /* include separator */
4034  SCIP_CALL( SCIPincludeSepaBasic(scip, &sepa, SEPA_NAME, SEPA_DESC, SEPA_PRIORITY, SEPA_FREQ, SEPA_MAXBOUNDDIST,
4035  SEPA_USESSUBSCIP, SEPA_DELAY,
4036  sepaExeclpOddcycle, NULL,
4037  sepadata) );
4038 
4039  assert(sepa != NULL);
4040 
4041  /* set non-NULL pointers to callback methods */
4042  SCIP_CALL( SCIPsetSepaCopy(scip, sepa, sepaCopyOddcycle) );
4043  SCIP_CALL( SCIPsetSepaFree(scip, sepa, sepaFreeOddcycle) );
4044  SCIP_CALL( SCIPsetSepaInit(scip, sepa, sepaInitOddcycle) );
4045  SCIP_CALL( SCIPsetSepaInitsol(scip, sepa, sepaInitsolOddcycle) );
4046 
4047  /* add oddcycle separator parameters */
4048  SCIP_CALL( SCIPaddBoolParam(scip, "separating/oddcycle/usegls",
4049  "should the search method by Groetschel, Lovasz, Schrijver be used? Otherwise use levelgraph method by Hoffman, Padberg.",
4050  &sepadata->usegls, FALSE, DEFAULT_USEGLS, NULL, NULL) );
4051  SCIP_CALL( SCIPaddBoolParam(scip, "separating/oddcycle/liftoddcycles",
4052  "should odd cycle cuts be lifted?",
4053  &sepadata->liftoddcycles, FALSE, DEFAULT_LIFTODDCYCLES, NULL, NULL) );
4054  SCIP_CALL( SCIPaddIntParam(scip, "separating/oddcycle/maxsepacuts",
4055  "maximal number of oddcycle cuts separated per separation round",
4056  &sepadata->maxsepacuts, FALSE, DEFAULT_MAXSEPACUTS, 0, INT_MAX, NULL, NULL) );
4057  SCIP_CALL( SCIPaddIntParam(scip, "separating/oddcycle/maxsepacutsroot",
4058  "maximal number of oddcycle cuts separated per separation round in the root node",
4059  &sepadata->maxsepacutsroot, FALSE, DEFAULT_MAXSEPACUTSROOT, 0, INT_MAX, NULL, NULL) );
4060 
4061  /* add advanced parameters */
4062  SCIP_CALL( SCIPaddIntParam(scip, "separating/oddcycle/maxrounds",
4063  "maximal number of oddcycle separation rounds per node (-1: unlimited)",
4064  &sepadata->maxrounds, FALSE, DEFAULT_MAXROUNDS, -1, INT_MAX, NULL, NULL) );
4065  SCIP_CALL( SCIPaddIntParam(scip, "separating/oddcycle/maxroundsroot",
4066  "maximal number of oddcycle separation rounds in the root node (-1: unlimited)",
4067  &sepadata->maxroundsroot, FALSE, DEFAULT_MAXROUNDSROOT, -1, INT_MAX, NULL, NULL) );
4068  SCIP_CALL( SCIPaddIntParam(scip, "separating/oddcycle/scalingfactor",
4069  "factor for scaling of the arc-weights",
4070  &sepadata->scale, TRUE, DEFAULT_SCALEFACTOR, 1, INT_MAX, NULL, NULL) );
4071  SCIP_CALL( SCIPaddBoolParam(scip, "separating/oddcycle/addselfarcs",
4072  "add links between a variable and its negated",
4073  &sepadata->addselfarcs, TRUE, DEFAULT_ADDSELFARCS, NULL, NULL) );
4074  SCIP_CALL( SCIPaddBoolParam(scip, "separating/oddcycle/repaircycles",
4075  "try to repair violated cycles with double appearance of a variable",
4076  &sepadata->repaircycles, TRUE, DEFAULT_REPAIRCYCLES, NULL, NULL) );
4077  SCIP_CALL( SCIPaddBoolParam(scip, "separating/oddcycle/includetriangles",
4078  "separate triangles found as 3-cycles or repaired larger cycles",
4079  &sepadata->includetriangles, TRUE, DEFAULT_INCLUDETRIANGLES, NULL, NULL) );
4080  SCIP_CALL( SCIPaddBoolParam(scip, "separating/oddcycle/multiplecuts",
4081  "even if a variable is already covered by a cut, still try it as start node for a cycle search",
4082  &sepadata->multiplecuts, TRUE, DEFAULT_MULTIPLECUTS, NULL, NULL) );
4083  SCIP_CALL( SCIPaddBoolParam(scip, "separating/oddcycle/allowmultiplecuts",
4084  "even if a variable is already covered by a cut, still allow another cut to cover it too",
4085  &sepadata->allowmultiplecuts, TRUE, DEFAULT_ALLOWMULTIPLECUTS, NULL, NULL) );
4086  SCIP_CALL( SCIPaddBoolParam(scip, "separating/oddcycle/lpliftcoef",
4087  "choose lifting candidate by coef*lpvalue or only by coef",
4088  &sepadata->lpliftcoef, TRUE, DEFAULT_LPLIFTCOEF, NULL, NULL) );
4089  SCIP_CALL( SCIPaddBoolParam(scip, "separating/oddcycle/recalcliftcoef",
4090  "calculate lifting coefficient of every candidate in every step (or only if its chosen)",
4091  &sepadata->recalcliftcoef, TRUE, DEFAULT_RECALCLIFTCOEF, NULL, NULL) );
4092  SCIP_CALL( SCIPaddIntParam(scip, "separating/oddcycle/sortswitch",
4093  "use sorted variable array (unsorted(0),maxlp(1),minlp(2),maxfrac(3),minfrac(4))",
4094  (int*) &sepadata->sortswitch, TRUE, DEFAULT_SORTSWITCH, 0, 4, NULL, NULL) );
4095  SCIP_CALL( SCIPaddBoolParam(scip, "separating/oddcycle/sortrootneighbors",
4096  "sort level of the root neighbors by fractionality (maxfrac)",
4097  &sepadata->sortrootneighbors, TRUE, DEFAULT_SORTROOTNEIGHBORS, NULL, NULL) );
4098  SCIP_CALL( SCIPaddIntParam(scip, "separating/oddcycle/percenttestvars",
4099  "percentage of variables to try the chosen method on [0-100]",
4100  &sepadata->percenttestvars, TRUE, DEFAULT_PERCENTTESTVARS, 0, 100, NULL, NULL) );
4101  SCIP_CALL( SCIPaddIntParam(scip, "separating/oddcycle/offsettestvars",
4102  "offset of variables to try the chosen method on (additional to the percentage of testvars)",
4103  &sepadata->offsettestvars, TRUE, DEFAULT_OFFSETTESTVARS, 0, INT_MAX, NULL, NULL) );
4104  SCIP_CALL( SCIPaddIntParam(scip, "separating/oddcycle/maxpernodeslevel",
4105  "percentage of nodes allowed in the same level of the level graph [0-100]",
4106  &sepadata->maxpernodeslevel, TRUE, DEFAULT_MAXPERNODESLEVEL, 0, 100, NULL, NULL) );
4107  SCIP_CALL( SCIPaddIntParam(scip, "separating/oddcycle/offsetnodeslevel",
4108  "offset of nodes allowed in the same level of the level graph (additional to the percentage of levelnodes)",
4109  &sepadata->offsetnodeslevel, TRUE, DEFAULT_OFFSETNODESLEVEL, 0, INT_MAX, NULL, NULL) );
4110  SCIP_CALL( SCIPaddIntParam(scip, "separating/oddcycle/maxnlevels",
4111  "maximal number of levels in level graph",
4112  &sepadata->maxnlevels, TRUE, DEFAULT_MAXNLEVELS, 0, INT_MAX, NULL, NULL) );
4113  SCIP_CALL( SCIPaddIntParam(scip, "separating/oddcycle/maxcutsroot",
4114  "maximal number of oddcycle cuts generated per chosen variable as root of the level graph",
4115  &sepadata->maxcutsroot, TRUE, DEFAULT_MAXCUTSROOT, 0, INT_MAX, NULL, NULL) );
4116  SCIP_CALL( SCIPaddIntParam(scip, "separating/oddcycle/maxcutslevel",
4117  "maximal number of oddcycle cuts generated in every level of the level graph",
4118  &sepadata->maxcutslevel, TRUE, DEFAULT_MAXCUTSLEVEL, 0, INT_MAX, NULL, NULL) );
4119  SCIP_CALL( SCIPaddIntParam(scip, "separating/oddcycle/maxreference",
4120  "minimal weight on an edge (in level graph or bipartite graph)",
4121  &sepadata->maxreference, TRUE, DEFAULT_MAXREFERENCE, 0, INT_MAX, NULL, NULL) );
4122  SCIP_CALL( SCIPaddIntParam(scip, "separating/oddcycle/maxunsucessfull",
4123  "number of unsuccessful calls at current node",
4124  &sepadata->maxunsucessfull, TRUE, DEFAULT_MAXUNSUCESSFULL, 0, INT_MAX, NULL, NULL) );
4125  SCIP_CALL( SCIPaddIntParam(scip, "separating/oddcycle/cutthreshold",
4126  "maximal number of other cuts s.t. separation is applied (-1 for direct call)",
4127  &sepadata->cutthreshold, TRUE, DEFAULT_CUTTHRESHOLD, -1, INT_MAX, NULL, NULL) );
4128 
4129  return SCIP_OKAY;
4130 }
4131