Scippy

SCIP

Solving Constraint Integer Programs

sepa_closecuts.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-2017 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_closecuts.c
17  * @brief closecuts meta separator
18  * @author Marc Pfetsch
19  *
20  * This separator generates a convex combination of the current LP solution and either the best
21  * primal feasible solution or an interior point of the LP relaxation. If the convex combination is
22  * proper, the new point is closer to the convex hull of the feasible points. The separator then
23  * calls all other separators to separate this point. The idea is that in this way possibly "deeper"
24  * cuts are generated. Note, however, that the new point is not a basic solution, i.e., separators
25  * relying basis information, e.g., Gomory cuts, will not work.
26  *
27  * The other cuts are generated via the sepasol() callbacks in constraints handlers or separators.
28  *
29  * This separator stops after a certain number (parameter @p maxunsuccessful) of unsuccessful
30  * calls. It also inhibits the separation of the ordinary LP solution if it already generated enough
31  * (parameter @p sepathreshold) cuts. The convex combination is determined via the parameter @p
32  * sepacombvalue.
33  *
34  * In general, this separator makes sense if it is expected that there will be many separation
35  * rounds and many cuts will be again deleted, because they are not active after a certain number of
36  * rounds. In particular, branch-and-cut algorithms for combinatorial optimization problems form
37  * good candidates.
38  *
39  * The idea seems to be first proposed in the context of the travelling salesman problem, see@par
40  * The Traveling Salesman Problem: A Computational Study@n
41  * David L. Applegate, Robert E. Bixby, Vasek Chvatal & William J. Cook@n
42  * Princeton University Press 2006@n
43  *
44  * for more details. See also@par
45  * Acceleration of cutting-plane and column generation algorithms: Applications to network design.@n
46  * Walid Ben-Ameur, Jose Neto@n
47  * Networks 49(1): 3-17 (2007).
48  */
49 
50 /*---+----1----+----2----+----3----+----4----+----5----+----6----+----7----+----8----+----9----+----0----+----1----+----2*/
51 
52 #include <assert.h>
53 #include <string.h>
54 
55 #include "scip/sepa_closecuts.h"
56 
57 
58 #define SEPA_NAME "closecuts"
59 #define SEPA_DESC "closecuts meta separator"
60 #define SEPA_PRIORITY 1000000
61 #define SEPA_FREQ -1
62 #define SEPA_MAXBOUNDDIST 1.0
63 #define SEPA_USESSUBSCIP FALSE /**< does the separator use a secondary SCIP instance? */
64 #define SEPA_DELAY FALSE /**< should separation method be delayed, if other separators found cuts? */
65 
66 
67 /* default values for parameters */
68 #define SCIP_DEFAULT_SEPARELINT TRUE /**< generate close cuts w.r.t. relative interior point (best solution otherwise)? */
69 #define SCIP_DEFAULT_SEPACOMBVALUE 0.30 /**< convex combination value for close cuts */
70 #define SCIP_DEFAULT_SEPATHRESHOLD 50 /**< threshold on number of generated cuts below which the ordinary separation is started */
71 #define SCIP_DEFAULT_INCLOBJCUTOFF FALSE /**< include the objective cutoff when computing the relative interior? */
72 #define SCIP_DEFAULT_RECOMPUTERELINT FALSE /**< recompute relative interior in each separation call? */
73 #define SCIP_DEFAULT_MAXUNSUCCESSFUL 0 /**< turn off separation in current node after unsuccessful calls (-1 never turn off) */
74 #define SCIP_DEFAULT_MAXLPITERFACTOR 10.0 /**< factor for maximal LP iterations in relative interior computation compared to node LP iterations */
75 
76 #define SCIP_MIN_LPITERS 100 /**< minimum number of allowed LP iterations in relative interior computation */
77 
78 
79 /** separator data */
80 struct SCIP_SepaData
81 {
82  SCIP_Bool separelint; /**< generate close cuts w.r.t. relative interior point (best solution otherwise)? */
83  SCIP_Bool triedRelint; /**< tried to compute relative interior */
84  SCIP_Real sepacombvalue; /**< convex combination value for close cuts */
85  int sepathreshold; /**< threshold on number of generated cuts below which the ordinary separation is started */
86  SCIP_Bool inclobjcutoff; /**< include the objective cutoff when computing the relative interior? */
87  SCIP_Bool recomputerelint; /**< recompute relative interior in each separation call? */
88  int maxunsuccessful; /**< turn off separation in current node after unsuccessful calls (-1 never turn off) */
89  SCIP_SOL* sepasol; /**< solution that can be used for generating close cuts */
90  SCIP_Longint discardnode; /**< number of node for which separation is discarded */
91  SCIP_Real maxlpiterfactor; /**< factor for maximal LP iterations in relative interior computation compared to node LP iterations */
92  int nunsuccessful; /**< number of consecutive unsuccessful calls */
93 };
94 
95 
96 /** generate point for close cut separation
97  *
98  * The constructed point is the convex combination of the point stored in set->closesol and the
99  * current LP solution. The convexity parameter is set->sepa_closecombvalue. If this parameter is
100  * 0, the point coincides with the LP solution.
101  */
102 static
104  SCIP* scip, /**< SCIP data structure */
105  SCIP_SEPADATA* sepadata, /**< separator data */
106  SCIP_SOL** point /**< point to be generated (or NULL if unsuccessful) */
107  )
108 {
109  SCIP_VAR** vars;
110  SCIP_VAR* var;
111  SCIP_Real val;
112  SCIP_Real alpha;
113  SCIP_Real onealpha;
114  SCIP_Real lb;
115  SCIP_Real ub;
116  int nvars;
117  int i;
118 
119  assert( scip != NULL );
120  assert( point != NULL );
121 
122  *point = NULL;
123  if ( sepadata->sepasol == NULL )
124  return SCIP_OKAY;
125 
126  alpha = sepadata->sepacombvalue;
127  if ( alpha < 0.001 )
128  return SCIP_OKAY;
129  onealpha = 1.0 - alpha;
130 
131  /* create solution */
132  SCIP_CALL( SCIPcreateSol(scip, point, NULL) );
133 
134  /* generate convex combination */
135  vars = SCIPgetVars(scip);
136  nvars = SCIPgetNVars(scip);
137  for (i = 0; i < nvars; ++i)
138  {
139  var = vars[i];
140  val = alpha * SCIPgetSolVal(scip, sepadata->sepasol, var) + onealpha * SCIPvarGetLPSol(var);
141 
142  /* If both the LP relaxation and the base point respect the variable bounds, the computed point will satisfy them
143  * as well. However, variables might be fixed (e.g. by branching) since the time of the computation of the base
144  * point. Thus, we adapt the value to lie inside the bounds in optimized mode. */
145  lb = SCIPvarGetLbLocal(var);
146  ub = SCIPvarGetUbLocal(var);
147  val = MAX(val, lb);
148  val = MIN(val, ub);
149 
150  if ( ! SCIPisZero(scip, val) )
151  {
152  SCIP_CALL( SCIPsetSolVal(scip, *point, var, val) );
153  }
154  }
155 
156  return SCIP_OKAY;
157 }
158 
159 
160 /*
161  * Callback methods of separator
162  */
163 
164 
165 /** copy method for separator plugins (called when SCIP copies plugins) */
166 static
167 SCIP_DECL_SEPACOPY(sepaCopyClosecuts)
168 { /*lint --e{715}*/
169  assert( scip != NULL );
170  assert( sepa != NULL );
171  assert( strcmp(SCIPsepaGetName(sepa), SEPA_NAME) == 0 );
172 
173  /* call inclusion method of constraint handler */
175 
176  return SCIP_OKAY;
177 }
178 
179 /** destructor of separator to free user data (called when SCIP is exiting) */
180 static
181 SCIP_DECL_SEPAFREE(sepaFreeClosecuts)
182 { /*lint --e{715}*/
183  SCIP_SEPADATA* sepadata;
184 
185  assert( sepa != NULL );
186  assert( strcmp(SCIPsepaGetName(sepa), SEPA_NAME) == 0 );
187 
188  /* free separator data */
189  sepadata = SCIPsepaGetData(sepa);
190  assert( sepadata != NULL );
191 
192  SCIPfreeBlockMemory(scip, &sepadata);
193 
194  SCIPsepaSetData(sepa, NULL);
195 
196  return SCIP_OKAY;
197 }
198 
199 
200 /** solving process deinitialization method of separator (called before branch and bound process data is freed) */
201 static
202 SCIP_DECL_SEPAEXITSOL(sepaExitsolClosecuts)
203 { /*lint --e{715}*/
204  SCIP_SEPADATA* sepadata;
205 
206  assert( sepa != NULL );
207  assert( strcmp(SCIPsepaGetName(sepa), SEPA_NAME) == 0 );
208 
209  sepadata = SCIPsepaGetData(sepa);
210  assert( sepadata != NULL );
211 
212  if ( sepadata->separelint && sepadata->sepasol != NULL )
213  {
214  SCIP_CALL( SCIPfreeSol(scip, &sepadata->sepasol) );
215  sepadata->triedRelint = FALSE;
216  }
217 
218  return SCIP_OKAY;
219 }
220 
221 
222 /** LP solution separation method of separator */
223 static
224 SCIP_DECL_SEPAEXECLP(sepaExeclpClosecuts)
225 { /*lint --e{715}*/
226  SCIP_SEPADATA* sepadata;
227  SCIP_Longint currentnodenumber;
228  SCIP_SOL* point = NULL;
229  SCIP_Bool isroot;
230 
231  assert( sepa != NULL );
232  assert( strcmp(SCIPsepaGetName(sepa), SEPA_NAME) == 0 );
233  assert( result != NULL );
234 
235  *result = SCIP_DIDNOTRUN;
236 
237  /* only call separator, if LP has been solved (need LP to compute separation point) */
239  return SCIP_OKAY;
240 
241  /* only call separator, if there are fractional variables */
242  if ( SCIPgetNLPBranchCands(scip) == 0 )
243  return SCIP_OKAY;
244 
245  /* exit if we stopped ... */
246  if ( SCIPisStopped(scip) )
247  return SCIP_OKAY;
248 
249  /* get separation data */
250  sepadata = SCIPsepaGetData(sepa);
251  assert( sepadata != NULL );
252 
253  /* exit if we already decided to discard the current node */
254  currentnodenumber = SCIPnodeGetNumber(SCIPgetCurrentNode(scip));
255  if ( sepadata->discardnode == currentnodenumber )
256  return SCIP_OKAY;
257 
258  SCIPdebugMsg(scip, "Separation method of closecuts separator.\n");
259 
260  /* check whether we have to compute a relative interior point */
261  if ( sepadata->separelint )
262  {
263  if ( sepadata->recomputerelint )
264  {
265  /* check if previous relative interior point should be forgotten, otherwise it is computed only once and the
266  * same point is used for all nodes */
267  if ( sepadata->sepasol != NULL )
268  {
269  SCIP_CALL( SCIPfreeSol(scip, &sepadata->sepasol) );
270  sepadata->triedRelint = FALSE;
271  }
272  }
273  else
274  {
275  /* skip execution, if we unsuccessfully tried to compute a relative interior point */
276  if ( sepadata->sepasol == NULL && sepadata->triedRelint )
277  return SCIP_OKAY;
278  }
279 
280  /* if relative interior point is not available ... */
281  if ( sepadata->sepasol == NULL )
282  {
283  SCIP_Longint nlpiters;
284  SCIP_Real timelimit;
285  int iterlimit;
286 
287  /* prepare time limit */
288  SCIP_CALL( SCIPgetRealParam(scip, "limits/time", &timelimit) );
289  if ( ! SCIPisInfinity(scip, timelimit) )
290  timelimit -= SCIPgetSolvingTime(scip);
291  /* exit if no time left */
292  if ( timelimit <= 0.0 )
293  return SCIP_OKAY;
294 
295  /* determine iteration limit */
296  if ( sepadata->maxlpiterfactor < 0.0 || SCIPisInfinity(scip, sepadata->maxlpiterfactor) )
297  iterlimit = INT_MAX;
298  else
299  {
300  /* determine iteration limit; the number of iterations in the root is only set after its solution, but the
301  * total number of LP iterations is always updated. */
302  if ( SCIPgetDepth(scip) == 0 )
303  nlpiters = SCIPgetNLPIterations(scip);
304  else
305  nlpiters = SCIPgetNRootLPIterations(scip);
306  iterlimit = (int)(sepadata->maxlpiterfactor * nlpiters);
307  iterlimit = MAX(iterlimit, SCIP_MIN_LPITERS);
308  assert(iterlimit > 0);
309  }
310 
311  SCIPverbMessage(scip, SCIP_VERBLEVEL_MINIMAL, 0, "Computing relative interior point (time limit: %g, iter limit: %d) ...\n", timelimit, iterlimit);
312  SCIP_CALL( SCIPcomputeLPRelIntPoint(scip, TRUE, sepadata->inclobjcutoff, timelimit, iterlimit, &sepadata->sepasol) );
313  sepadata->triedRelint = TRUE;
314  }
315  }
316  else
317  {
318  /* get best solution (NULL if not present) */
319  sepadata->sepasol = SCIPgetBestSol(scip);
320  }
321 
322  /* separate close cuts */
323  if ( sepadata->sepasol != NULL )
324  {
325  SCIPdebugMsg(scip, "Generating close cuts ... (combination value: %f)\n", sepadata->sepacombvalue);
326  *result = SCIP_DIDNOTFIND;
327 
328  /* generate point to be separated */
329  SCIP_CALL( generateCloseCutPoint(scip, sepadata, &point) );
330 
331  /* apply a separation round to generated point */
332  if ( point != NULL )
333  {
334  int noldcuts;
335  SCIP_Bool delayed;
336  SCIP_Bool cutoff;
337 
338  noldcuts = SCIPgetNCuts(scip);
339  isroot = (SCIP_Bool) (SCIPgetDepth(scip) == 0);
340 
341  /* separate solution via other separators */
342  SCIP_CALL( SCIPseparateSol(scip, point, isroot, TRUE, FALSE, &delayed, &cutoff) );
343 
344  SCIP_CALL( SCIPfreeSol(scip, &point) );
345  assert( point == NULL );
346 
347  /* the cuts might not violated by the current LP if the computed point is strange */
349 
350  if ( cutoff )
351  *result = SCIP_CUTOFF;
352  else
353  {
354  if ( SCIPgetNCuts(scip) - noldcuts > sepadata->sepathreshold )
355  {
356  sepadata->nunsuccessful = 0;
357  *result = SCIP_NEWROUND;
358  }
359  else
360  {
361  if ( SCIPgetNCuts(scip) > noldcuts )
362  {
363  sepadata->nunsuccessful = 0;
364  *result = SCIP_SEPARATED;
365  }
366  else
367  ++sepadata->nunsuccessful;
368  }
369  }
370 
371  SCIPdebugMsg(scip, "Separated close cuts: %d (enoughcuts: %d, unsuccessful: %d).\n", SCIPgetNCuts(scip) - noldcuts,
372  SCIPgetNCuts(scip) - noldcuts > sepadata->sepathreshold, sepadata->nunsuccessful);
373 
374  if ( sepadata->maxunsuccessful >= 0 && sepadata->nunsuccessful > sepadata->maxunsuccessful )
375  {
376  SCIPdebugMsg(scip, "Turn off close cut separation, because of %d unsuccessful calls.\n", sepadata->nunsuccessful);
377  sepadata->discardnode = currentnodenumber;
378  sepadata->nunsuccessful = 0;
379  }
380  }
381  }
382 
383  return SCIP_OKAY;
384 }
385 
386 
387 /*
388  * separator specific interface methods
389  */
390 
391 /** creates the closecuts separator and includes it in SCIP */
393  SCIP* scip /**< SCIP data structure */
394  )
395 {
396  SCIP_SEPADATA* sepadata;
397  SCIP_SEPA* sepa;
398 
399  /* create closecuts separator data */
400  SCIP_CALL( SCIPallocBlockMemory(scip, &sepadata) );
401  sepadata->sepasol = NULL;
402  sepadata->discardnode = -1;
403  sepadata->nunsuccessful = 0;
404  sepadata->triedRelint = FALSE;
405 
406  /* include separator */
408  sepaExeclpClosecuts, NULL, sepadata) );
409 
410  assert(sepa != NULL);
411 
412  /* set non-NULL pointers to callback methods */
413  SCIP_CALL( SCIPsetSepaCopy(scip, sepa, sepaCopyClosecuts) );
414  SCIP_CALL( SCIPsetSepaFree(scip, sepa, sepaFreeClosecuts) );
415  SCIP_CALL( SCIPsetSepaExitsol(scip, sepa, sepaExitsolClosecuts) );
416 
417  /* add closecuts separator parameters */
419  "separating/closecuts/separelint",
420  "generate close cuts w.r.t. relative interior point (best solution otherwise)?",
421  &sepadata->separelint, TRUE, SCIP_DEFAULT_SEPARELINT, NULL, NULL) );
422 
424  "separating/closecuts/sepacombvalue",
425  "convex combination value for close cuts",
426  &sepadata->sepacombvalue, TRUE, SCIP_DEFAULT_SEPACOMBVALUE, 0.0, 1.0,
427  NULL, NULL) );
428 
430  "separating/closecuts/closethres",
431  "threshold on number of generated cuts below which the ordinary separation is started",
432  &sepadata->sepathreshold, TRUE, SCIP_DEFAULT_SEPATHRESHOLD, -1, INT_MAX, NULL, NULL) );
433 
435  "separating/closecuts/inclobjcutoff",
436  "include an objective cutoff when computing the relative interior?",
437  &sepadata->inclobjcutoff, TRUE, SCIP_DEFAULT_INCLOBJCUTOFF, NULL, NULL) );
438 
440  "separating/closecuts/recomputerelint",
441  "recompute relative interior point in each separation call?",
442  &sepadata->recomputerelint, TRUE, SCIP_DEFAULT_RECOMPUTERELINT, NULL, NULL) );
443 
445  "separating/closecuts/maxunsuccessful",
446  "turn off separation in current node after unsuccessful calls (-1 never turn off)",
447  &sepadata->maxunsuccessful, TRUE, SCIP_DEFAULT_MAXUNSUCCESSFUL, -1, INT_MAX, NULL, NULL) );
448 
450  "separating/closecuts/maxlpiterfactor",
451  "factor for maximal LP iterations in relative interior computation compared to node LP iterations (negative for no limit)",
452  &sepadata->maxlpiterfactor, TRUE, SCIP_DEFAULT_MAXLPITERFACTOR, -1.0, SCIP_REAL_MAX, NULL, NULL) );
453 
454  return SCIP_OKAY;
455 }
456 
457 /** sets point to be used as base point for computing the point to be separated
458  *
459  * The point is only stored if separation of relative interior points is used. The solution is copied.
460  */
462  SCIP* scip, /**< SCIP data structure */
463  SCIP_SOL* sol /**< base point solution */
464  )
465 {
466  SCIP_SEPA* sepa;
467  SCIP_SEPADATA* sepadata;
468 
469  assert( scip != NULL );
470 
471  /* find separator */
472  sepa = SCIPfindSepa(scip, SEPA_NAME);
473  if ( sepa == NULL )
474  {
475  SCIPerrorMessage("Could not find separator <%s>.\n", SEPA_NAME);
476  return SCIP_PLUGINNOTFOUND;
477  }
478  assert( strcmp(SCIPsepaGetName(sepa), SEPA_NAME) == 0 );
479 
480  /* get sepadata */
481  sepadata = SCIPsepaGetData(sepa);
482  assert( sepadata != NULL );
483 
484  /* store point if we have to separate relative interior points */
485  if ( sepadata->separelint )
486  {
487  /* possibly free solution */
488  if ( sepadata->sepasol != NULL )
489  {
490  SCIP_CALL( SCIPfreeSol(scip, &sepadata->sepasol) );
491  }
492 
493  /* copy and store solution */
494  SCIP_CALL( SCIPcreateSolCopy(scip, &sepadata->sepasol, sol) );
495  sepadata->triedRelint = TRUE;
496  }
497 
498  return SCIP_OKAY;
499 }
#define SCIP_DEFAULT_INCLOBJCUTOFF
SCIP_Real SCIPgetSolvingTime(SCIP *scip)
Definition: scip.c:46300
SCIP_Longint SCIPgetNRootLPIterations(SCIP *scip)
Definition: scip.c:42371
SCIP_RETCODE SCIPcomputeLPRelIntPoint(SCIP *scip, SCIP_Bool relaxrows, SCIP_Bool inclobjcutoff, SCIP_Real timelimit, int iterlimit, SCIP_SOL **point)
Definition: scip.c:30134
SCIP_NODE * SCIPgetCurrentNode(SCIP *scip)
Definition: scip.c:41398
SCIP_Longint SCIPgetNLPIterations(SCIP *scip)
Definition: scip.c:42327
SCIP_RETCODE SCIPgetRealParam(SCIP *scip, const char *name, SCIP_Real *value)
Definition: scip.c:4489
SCIP_SEPA * SCIPfindSepa(SCIP *scip, const char *name)
Definition: scip.c:7523
SCIP_Real SCIPvarGetLbLocal(SCIP_VAR *var)
Definition: var.c:17332
#define FALSE
Definition: def.h:64
#define TRUE
Definition: def.h:63
const char * SCIPsepaGetName(SCIP_SEPA *sepa)
Definition: sepa.c:646
enum SCIP_Retcode SCIP_RETCODE
Definition: type_retcode.h:53
#define SEPA_DELAY
static SCIP_DECL_SEPACOPY(sepaCopyClosecuts)
SCIP_RETCODE SCIPsetBasePointClosecuts(SCIP *scip, SCIP_SOL *sol)
#define SCIPfreeBlockMemory(scip, ptr)
Definition: scip.h:22602
#define SCIPallocBlockMemory(scip, ptr)
Definition: scip.h:22585
int SCIPgetNLPBranchCands(SCIP *scip)
Definition: scip.c:37028
SCIP_RETCODE SCIPsetSepaCopy(SCIP *scip, SCIP_SEPA *sepa, SCIP_DECL_SEPACOPY((*sepacopy)))
Definition: scip.c:7427
#define SCIPdebugMsg
Definition: scip.h:455
SCIP_RETCODE SCIPaddIntParam(SCIP *scip, const char *name, const char *desc, int *valueptr, SCIP_Bool isadvanced, int defaultvalue, int minvalue, int maxvalue, SCIP_DECL_PARAMCHGD((*paramchgd)), SCIP_PARAMDATA *paramdata)
Definition: scip.c:4265
SCIP_SEPADATA * SCIPsepaGetData(SCIP_SEPA *sepa)
Definition: sepa.c:557
#define SCIP_DEFAULT_RECOMPUTERELINT
SCIP_Longint SCIPnodeGetNumber(SCIP_NODE *node)
Definition: tree.c:7340
SCIP_RETCODE SCIPcreateSolCopy(SCIP *scip, SCIP_SOL **sol, SCIP_SOL *sourcesol)
Definition: scip.c:38162
#define SCIPerrorMessage
Definition: pub_message.h:45
closecuts meta separator
void SCIPsepaSetData(SCIP_SEPA *sepa, SCIP_SEPADATA *sepadata)
Definition: sepa.c:567
SCIP_RETCODE SCIPseparateSol(SCIP *scip, SCIP_SOL *sol, SCIP_Bool pretendroot, SCIP_Bool allowlocal, SCIP_Bool onlydelayed, SCIP_Bool *delayed, SCIP_Bool *cutoff)
Definition: scip.c:35139
SCIP_Real SCIPvarGetLPSol(SCIP_VAR *var)
Definition: var.c:17650
#define SEPA_NAME
#define SCIP_CALL(x)
Definition: def.h:350
SCIP_RETCODE SCIPremoveInefficaciousCuts(SCIP *scip)
Definition: scip.c:35226
void SCIPverbMessage(SCIP *scip, SCIP_VERBLEVEL msgverblevel, FILE *file, const char *formatstr,...)
Definition: scip.c:1360
static SCIP_RETCODE generateCloseCutPoint(SCIP *scip, SCIP_SEPADATA *sepadata, SCIP_SOL **point)
SCIP_RETCODE SCIPincludeSepaBasic(SCIP *scip, SCIP_SEPA **sepa, const char *name, const char *desc, int priority, int freq, SCIP_Real maxbounddist, SCIP_Bool usessubscip, SCIP_Bool delay, SCIP_DECL_SEPAEXECLP((*sepaexeclp)), SCIP_DECL_SEPAEXECSOL((*sepaexecsol)), SCIP_SEPADATA *sepadata)
Definition: scip.c:7385
static SCIP_DECL_SEPAEXITSOL(sepaExitsolClosecuts)
SCIP_RETCODE SCIPsetSolVal(SCIP *scip, SCIP_SOL *sol, SCIP_VAR *var, SCIP_Real val)
Definition: scip.c:38765
#define SCIP_DEFAULT_MAXLPITERFACTOR
#define SCIP_DEFAULT_SEPATHRESHOLD
SCIP_RETCODE SCIPsetSepaExitsol(SCIP *scip, SCIP_SEPA *sepa, SCIP_DECL_SEPAEXITSOL((*sepaexitsol)))
Definition: scip.c:7507
#define SCIP_Bool
Definition: def.h:61
SCIP_LPSOLSTAT SCIPgetLPSolstat(SCIP *scip)
Definition: scip.c:29287
int SCIPgetDepth(SCIP *scip)
Definition: scip.c:43039
#define MAX(x, y)
Definition: tclique_def.h:75
SCIP_RETCODE SCIPfreeSol(SCIP *scip, SCIP_SOL **sol)
Definition: scip.c:38529
#define SCIP_DEFAULT_SEPARELINT
SCIP_Bool SCIPisInfinity(SCIP *scip, SCIP_Real val)
Definition: scip.c:47033
#define SEPA_PRIORITY
int SCIPgetNVars(SCIP *scip)
Definition: scip.c:11806
SCIP_RETCODE SCIPincludeSepaClosecuts(SCIP *scip)
#define SCIP_REAL_MAX
Definition: def.h:150
SCIP_RETCODE SCIPsetSepaFree(SCIP *scip, SCIP_SEPA *sepa, SCIP_DECL_SEPAFREE((*sepafree)))
Definition: scip.c:7443
SCIP_SOL * SCIPgetBestSol(SCIP *scip)
Definition: scip.c:39876
#define SEPA_FREQ
int SCIPgetNCuts(SCIP *scip)
Definition: scip.c:35190
SCIP_VAR ** SCIPgetVars(SCIP *scip)
Definition: scip.c:11761
#define SCIP_DEFAULT_MAXUNSUCCESSFUL
#define SCIP_Real
Definition: def.h:149
SCIP_Bool SCIPisStopped(SCIP *scip)
Definition: scip.c:1145
#define SEPA_DESC
#define SEPA_USESSUBSCIP
#define SCIP_Longint
Definition: def.h:134
static SCIP_DECL_SEPAFREE(sepaFreeClosecuts)
#define SEPA_MAXBOUNDDIST
SCIP_Bool SCIPisZero(SCIP *scip, SCIP_Real val)
Definition: scip.c:47070
SCIP_Real SCIPvarGetUbLocal(SCIP_VAR *var)
Definition: var.c:17342
#define SCIP_MIN_LPITERS
SCIP_Real SCIPgetSolVal(SCIP *scip, SCIP_SOL *sol, SCIP_VAR *var)
Definition: scip.c:38905
#define SCIP_DEFAULT_SEPACOMBVALUE
SCIP_RETCODE SCIPaddRealParam(SCIP *scip, const char *name, const char *desc, SCIP_Real *valueptr, SCIP_Bool isadvanced, SCIP_Real defaultvalue, SCIP_Real minvalue, SCIP_Real maxvalue, SCIP_DECL_PARAMCHGD((*paramchgd)), SCIP_PARAMDATA *paramdata)
Definition: scip.c:4321
struct SCIP_SepaData SCIP_SEPADATA
Definition: type_sepa.h:38
SCIP_RETCODE SCIPaddBoolParam(SCIP *scip, const char *name, const char *desc, SCIP_Bool *valueptr, SCIP_Bool isadvanced, SCIP_Bool defaultvalue, SCIP_DECL_PARAMCHGD((*paramchgd)), SCIP_PARAMDATA *paramdata)
Definition: scip.c:4239
static SCIP_DECL_SEPAEXECLP(sepaExeclpClosecuts)
SCIP_RETCODE SCIPcreateSol(SCIP *scip, SCIP_SOL **sol, SCIP_HEUR *heur)
Definition: scip.c:37872