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-2018 Konrad-Zuse-Zentrum */
7 /* fuer Informationstechnik Berlin */
8 /* */
9 /* SCIP is distributed under the terms of the ZIB Academic License. */
10 /* */
11 /* You should have received a copy of the ZIB Academic License */
12 /* along with SCIP; see the file COPYING. If not visit 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 "scip/pub_message.h"
53 #include "scip/pub_sepa.h"
54 #include "scip/pub_tree.h"
55 #include "scip/pub_var.h"
56 #include "scip/scip_branch.h"
57 #include "scip/scip_cut.h"
58 #include "scip/scip_general.h"
59 #include "scip/scip_lp.h"
60 #include "scip/scip_mem.h"
61 #include "scip/scip_message.h"
62 #include "scip/scip_numerics.h"
63 #include "scip/scip_param.h"
64 #include "scip/scip_prob.h"
65 #include "scip/scip_sepa.h"
66 #include "scip/scip_sol.h"
67 #include "scip/scip_solvingstats.h"
68 #include "scip/scip_timing.h"
69 #include "scip/scip_tree.h"
70 #include "scip/sepa_closecuts.h"
71 #include <string.h>
72 
73 
74 
75 #define SEPA_NAME "closecuts"
76 #define SEPA_DESC "closecuts meta separator"
77 #define SEPA_PRIORITY 1000000
78 #define SEPA_FREQ -1
79 #define SEPA_MAXBOUNDDIST 1.0
80 #define SEPA_USESSUBSCIP FALSE /**< does the separator use a secondary SCIP instance? */
81 #define SEPA_DELAY FALSE /**< should separation method be delayed, if other separators found cuts? */
82 
83 
84 /* default values for parameters */
85 #define SCIP_DEFAULT_SEPARELINT TRUE /**< generate close cuts w.r.t. relative interior point (best solution otherwise)? */
86 #define SCIP_DEFAULT_SEPACOMBVALUE 0.30 /**< convex combination value for close cuts */
87 #define SCIP_DEFAULT_SEPATHRESHOLD 50 /**< threshold on number of generated cuts below which the ordinary separation is started */
88 #define SCIP_DEFAULT_INCLOBJCUTOFF FALSE /**< include the objective cutoff when computing the relative interior? */
89 #define SCIP_DEFAULT_RECOMPUTERELINT FALSE /**< recompute relative interior in each separation call? */
90 #define SCIP_DEFAULT_MAXUNSUCCESSFUL 0 /**< turn off separation in current node after unsuccessful calls (-1 never turn off) */
91 #define SCIP_DEFAULT_MAXLPITERFACTOR 10.0 /**< factor for maximal LP iterations in relative interior computation compared to node LP iterations */
92 
93 #define SCIP_MIN_LPITERS 100 /**< minimum number of allowed LP iterations in relative interior computation */
94 
95 
96 /** separator data */
97 struct SCIP_SepaData
98 {
99  SCIP_Bool separelint; /**< generate close cuts w.r.t. relative interior point (best solution otherwise)? */
100  SCIP_Bool triedRelint; /**< tried to compute relative interior */
101  SCIP_Real sepacombvalue; /**< convex combination value for close cuts */
102  int sepathreshold; /**< threshold on number of generated cuts below which the ordinary separation is started */
103  SCIP_Bool inclobjcutoff; /**< include the objective cutoff when computing the relative interior? */
104  SCIP_Bool recomputerelint; /**< recompute relative interior in each separation call? */
105  int maxunsuccessful; /**< turn off separation in current node after unsuccessful calls (-1 never turn off) */
106  SCIP_SOL* sepasol; /**< solution that can be used for generating close cuts */
107  SCIP_Longint discardnode; /**< number of node for which separation is discarded */
108  SCIP_Real maxlpiterfactor; /**< factor for maximal LP iterations in relative interior computation compared to node LP iterations */
109  int nunsuccessful; /**< number of consecutive unsuccessful calls */
110 };
111 
112 
113 /** generate point for close cut separation
114  *
115  * The constructed point is the convex combination of the point stored in set->closesol and the
116  * current LP solution. The convexity parameter is set->sepa_closecombvalue. If this parameter is
117  * 0, the point coincides with the LP solution.
118  */
119 static
121  SCIP* scip, /**< SCIP data structure */
122  SCIP_SEPADATA* sepadata, /**< separator data */
123  SCIP_SOL** point /**< point to be generated (or NULL if unsuccessful) */
124  )
125 {
126  SCIP_VAR** vars;
127  SCIP_VAR* var;
128  SCIP_Real val;
129  SCIP_Real alpha;
130  SCIP_Real onealpha;
131  SCIP_Real lb;
132  SCIP_Real ub;
133  int nvars;
134  int i;
135 
136  assert( scip != NULL );
137  assert( point != NULL );
138 
139  *point = NULL;
140  if ( sepadata->sepasol == NULL )
141  return SCIP_OKAY;
142 
143  alpha = sepadata->sepacombvalue;
144  if ( alpha < 0.001 )
145  return SCIP_OKAY;
146  onealpha = 1.0 - alpha;
147 
148  /* create solution */
149  SCIP_CALL( SCIPcreateSol(scip, point, NULL) );
150 
151  /* generate convex combination */
152  vars = SCIPgetVars(scip);
153  nvars = SCIPgetNVars(scip);
154  for (i = 0; i < nvars; ++i)
155  {
156  var = vars[i];
157  val = alpha * SCIPgetSolVal(scip, sepadata->sepasol, var) + onealpha * SCIPvarGetLPSol(var);
158 
159  /* If both the LP relaxation and the base point respect the variable bounds, the computed point will satisfy them
160  * as well. However, variables might be fixed (e.g. by branching) since the time of the computation of the base
161  * point. Thus, we adapt the value to lie inside the bounds in optimized mode. */
162  lb = SCIPvarGetLbLocal(var);
163  ub = SCIPvarGetUbLocal(var);
164  val = MAX(val, lb);
165  val = MIN(val, ub);
166 
167  if ( ! SCIPisZero(scip, val) )
168  {
169  SCIP_CALL( SCIPsetSolVal(scip, *point, var, val) );
170  }
171  }
172 
173  return SCIP_OKAY;
174 }
175 
176 
177 /*
178  * Callback methods of separator
179  */
180 
181 
182 /** copy method for separator plugins (called when SCIP copies plugins) */
183 static
184 SCIP_DECL_SEPACOPY(sepaCopyClosecuts)
185 { /*lint --e{715}*/
186  assert( scip != NULL );
187  assert( sepa != NULL );
188  assert( strcmp(SCIPsepaGetName(sepa), SEPA_NAME) == 0 );
189 
190  /* call inclusion method of constraint handler */
192 
193  return SCIP_OKAY;
194 }
195 
196 /** destructor of separator to free user data (called when SCIP is exiting) */
197 static
198 SCIP_DECL_SEPAFREE(sepaFreeClosecuts)
199 { /*lint --e{715}*/
200  SCIP_SEPADATA* sepadata;
201 
202  assert( sepa != NULL );
203  assert( strcmp(SCIPsepaGetName(sepa), SEPA_NAME) == 0 );
204 
205  /* free separator data */
206  sepadata = SCIPsepaGetData(sepa);
207  assert( sepadata != NULL );
208 
209  SCIPfreeBlockMemory(scip, &sepadata);
210 
211  SCIPsepaSetData(sepa, NULL);
212 
213  return SCIP_OKAY;
214 }
215 
216 
217 /** solving process deinitialization method of separator (called before branch and bound process data is freed) */
218 static
219 SCIP_DECL_SEPAEXITSOL(sepaExitsolClosecuts)
220 { /*lint --e{715}*/
221  SCIP_SEPADATA* sepadata;
222 
223  assert( sepa != NULL );
224  assert( strcmp(SCIPsepaGetName(sepa), SEPA_NAME) == 0 );
225 
226  sepadata = SCIPsepaGetData(sepa);
227  assert( sepadata != NULL );
228 
229  if ( sepadata->separelint && sepadata->sepasol != NULL )
230  {
231  SCIP_CALL( SCIPfreeSol(scip, &sepadata->sepasol) );
232  sepadata->triedRelint = FALSE;
233  }
234 
235  return SCIP_OKAY;
236 }
237 
238 
239 /** LP solution separation method of separator */
240 static
241 SCIP_DECL_SEPAEXECLP(sepaExeclpClosecuts)
242 { /*lint --e{715}*/
243  SCIP_SEPADATA* sepadata;
244  SCIP_Longint currentnodenumber;
245  SCIP_SOL* point = NULL;
246  SCIP_Bool isroot;
247 
248  assert( sepa != NULL );
249  assert( strcmp(SCIPsepaGetName(sepa), SEPA_NAME) == 0 );
250  assert( result != NULL );
251 
252  *result = SCIP_DIDNOTRUN;
253 
254  /* only call separator, if LP has been solved (need LP to compute separation point) */
256  return SCIP_OKAY;
257 
258  /* only call separator, if there are fractional variables */
259  if ( SCIPgetNLPBranchCands(scip) == 0 )
260  return SCIP_OKAY;
261 
262  /* exit if we stopped ... */
263  if ( SCIPisStopped(scip) )
264  return SCIP_OKAY;
265 
266  /* get separation data */
267  sepadata = SCIPsepaGetData(sepa);
268  assert( sepadata != NULL );
269 
270  /* exit if we already decided to discard the current node */
271  currentnodenumber = SCIPnodeGetNumber(SCIPgetCurrentNode(scip));
272  if ( sepadata->discardnode == currentnodenumber )
273  return SCIP_OKAY;
274 
275  SCIPdebugMsg(scip, "Separation method of closecuts separator.\n");
276 
277  /* check whether we have to compute a relative interior point */
278  if ( sepadata->separelint )
279  {
280  if ( sepadata->recomputerelint )
281  {
282  /* check if previous relative interior point should be forgotten, otherwise it is computed only once and the
283  * same point is used for all nodes */
284  if ( sepadata->sepasol != NULL )
285  {
286  SCIP_CALL( SCIPfreeSol(scip, &sepadata->sepasol) );
287  sepadata->triedRelint = FALSE;
288  }
289  }
290  else
291  {
292  /* skip execution, if we unsuccessfully tried to compute a relative interior point */
293  if ( sepadata->sepasol == NULL && sepadata->triedRelint )
294  return SCIP_OKAY;
295  }
296 
297  /* if relative interior point is not available ... */
298  if ( sepadata->sepasol == NULL )
299  {
300  SCIP_Longint nlpiters;
301  SCIP_Real timelimit;
302  int iterlimit;
303 
304  /* prepare time limit */
305  SCIP_CALL( SCIPgetRealParam(scip, "limits/time", &timelimit) );
306  if ( ! SCIPisInfinity(scip, timelimit) )
307  timelimit -= SCIPgetSolvingTime(scip);
308  /* exit if no time left */
309  if ( timelimit <= 0.0 )
310  return SCIP_OKAY;
311 
312  /* determine iteration limit */
313  if ( sepadata->maxlpiterfactor < 0.0 || SCIPisInfinity(scip, sepadata->maxlpiterfactor) )
314  iterlimit = INT_MAX;
315  else
316  {
317  /* determine iteration limit; the number of iterations in the root is only set after its solution, but the
318  * total number of LP iterations is always updated. */
319  if ( SCIPgetDepth(scip) == 0 )
320  nlpiters = SCIPgetNLPIterations(scip);
321  else
322  nlpiters = SCIPgetNRootLPIterations(scip);
323  iterlimit = (int)(sepadata->maxlpiterfactor * nlpiters);
324  iterlimit = MAX(iterlimit, SCIP_MIN_LPITERS);
325  assert(iterlimit > 0);
326  }
327 
328  SCIPverbMessage(scip, SCIP_VERBLEVEL_MINIMAL, 0, "Computing relative interior point (time limit: %g, iter limit: %d) ...\n", timelimit, iterlimit);
329  SCIP_CALL( SCIPcomputeLPRelIntPoint(scip, TRUE, sepadata->inclobjcutoff, timelimit, iterlimit, &sepadata->sepasol) );
330  sepadata->triedRelint = TRUE;
331  }
332  }
333  else
334  {
335  /* get best solution (NULL if not present) */
336  sepadata->sepasol = SCIPgetBestSol(scip);
337  }
338 
339  /* separate close cuts */
340  if ( sepadata->sepasol != NULL )
341  {
342  SCIPdebugMsg(scip, "Generating close cuts ... (combination value: %f)\n", sepadata->sepacombvalue);
343  *result = SCIP_DIDNOTFIND;
344 
345  /* generate point to be separated */
346  SCIP_CALL( generateCloseCutPoint(scip, sepadata, &point) );
347 
348  /* apply a separation round to generated point */
349  if ( point != NULL )
350  {
351  int noldcuts;
352  SCIP_Bool delayed;
353  SCIP_Bool cutoff;
354 
355  noldcuts = SCIPgetNCuts(scip);
356  isroot = (SCIP_Bool) (SCIPgetDepth(scip) == 0);
357 
358  /* separate solution via other separators */
359  SCIP_CALL( SCIPseparateSol(scip, point, isroot, TRUE, FALSE, &delayed, &cutoff) );
360 
361  SCIP_CALL( SCIPfreeSol(scip, &point) );
362  assert( point == NULL );
363 
364  /* the cuts might not violated by the current LP if the computed point is strange */
366 
367  if ( cutoff )
368  *result = SCIP_CUTOFF;
369  else
370  {
371  if ( SCIPgetNCuts(scip) - noldcuts > sepadata->sepathreshold )
372  {
373  sepadata->nunsuccessful = 0;
374  *result = SCIP_NEWROUND;
375  }
376  else
377  {
378  if ( SCIPgetNCuts(scip) > noldcuts )
379  {
380  sepadata->nunsuccessful = 0;
381  *result = SCIP_SEPARATED;
382  }
383  else
384  ++sepadata->nunsuccessful;
385  }
386  }
387 
388  SCIPdebugMsg(scip, "Separated close cuts: %d (enoughcuts: %d, unsuccessful: %d).\n", SCIPgetNCuts(scip) - noldcuts,
389  SCIPgetNCuts(scip) - noldcuts > sepadata->sepathreshold, sepadata->nunsuccessful);
390 
391  if ( sepadata->maxunsuccessful >= 0 && sepadata->nunsuccessful > sepadata->maxunsuccessful )
392  {
393  SCIPdebugMsg(scip, "Turn off close cut separation, because of %d unsuccessful calls.\n", sepadata->nunsuccessful);
394  sepadata->discardnode = currentnodenumber;
395  sepadata->nunsuccessful = 0;
396  }
397  }
398  }
399 
400  return SCIP_OKAY;
401 }
402 
403 
404 /*
405  * separator specific interface methods
406  */
407 
408 /** creates the closecuts separator and includes it in SCIP */
410  SCIP* scip /**< SCIP data structure */
411  )
412 {
413  SCIP_SEPADATA* sepadata;
414  SCIP_SEPA* sepa;
415 
416  /* create closecuts separator data */
417  SCIP_CALL( SCIPallocBlockMemory(scip, &sepadata) );
418  sepadata->sepasol = NULL;
419  sepadata->discardnode = -1;
420  sepadata->nunsuccessful = 0;
421  sepadata->triedRelint = FALSE;
422 
423  /* include separator */
425  sepaExeclpClosecuts, NULL, sepadata) );
426 
427  assert(sepa != NULL);
428 
429  /* set non-NULL pointers to callback methods */
430  SCIP_CALL( SCIPsetSepaCopy(scip, sepa, sepaCopyClosecuts) );
431  SCIP_CALL( SCIPsetSepaFree(scip, sepa, sepaFreeClosecuts) );
432  SCIP_CALL( SCIPsetSepaExitsol(scip, sepa, sepaExitsolClosecuts) );
433 
434  /* add closecuts separator parameters */
436  "separating/closecuts/separelint",
437  "generate close cuts w.r.t. relative interior point (best solution otherwise)?",
438  &sepadata->separelint, TRUE, SCIP_DEFAULT_SEPARELINT, NULL, NULL) );
439 
441  "separating/closecuts/sepacombvalue",
442  "convex combination value for close cuts",
443  &sepadata->sepacombvalue, TRUE, SCIP_DEFAULT_SEPACOMBVALUE, 0.0, 1.0,
444  NULL, NULL) );
445 
447  "separating/closecuts/closethres",
448  "threshold on number of generated cuts below which the ordinary separation is started",
449  &sepadata->sepathreshold, TRUE, SCIP_DEFAULT_SEPATHRESHOLD, -1, INT_MAX, NULL, NULL) );
450 
452  "separating/closecuts/inclobjcutoff",
453  "include an objective cutoff when computing the relative interior?",
454  &sepadata->inclobjcutoff, TRUE, SCIP_DEFAULT_INCLOBJCUTOFF, NULL, NULL) );
455 
457  "separating/closecuts/recomputerelint",
458  "recompute relative interior point in each separation call?",
459  &sepadata->recomputerelint, TRUE, SCIP_DEFAULT_RECOMPUTERELINT, NULL, NULL) );
460 
462  "separating/closecuts/maxunsuccessful",
463  "turn off separation in current node after unsuccessful calls (-1 never turn off)",
464  &sepadata->maxunsuccessful, TRUE, SCIP_DEFAULT_MAXUNSUCCESSFUL, -1, INT_MAX, NULL, NULL) );
465 
467  "separating/closecuts/maxlpiterfactor",
468  "factor for maximal LP iterations in relative interior computation compared to node LP iterations (negative for no limit)",
469  &sepadata->maxlpiterfactor, TRUE, SCIP_DEFAULT_MAXLPITERFACTOR, -1.0, SCIP_REAL_MAX, NULL, NULL) );
470 
471  return SCIP_OKAY;
472 }
473 
474 /** sets point to be used as base point for computing the point to be separated
475  *
476  * The point is only stored if separation of relative interior points is used. The solution is copied.
477  */
479  SCIP* scip, /**< SCIP data structure */
480  SCIP_SOL* sol /**< base point solution */
481  )
482 {
483  SCIP_SEPA* sepa;
484  SCIP_SEPADATA* sepadata;
485 
486  assert( scip != NULL );
487 
488  /* find separator */
489  sepa = SCIPfindSepa(scip, SEPA_NAME);
490  if ( sepa == NULL )
491  {
492  SCIPerrorMessage("Could not find separator <%s>.\n", SEPA_NAME);
493  return SCIP_PLUGINNOTFOUND;
494  }
495  assert( strcmp(SCIPsepaGetName(sepa), SEPA_NAME) == 0 );
496 
497  /* get sepadata */
498  sepadata = SCIPsepaGetData(sepa);
499  assert( sepadata != NULL );
500 
501  /* store point if we have to separate relative interior points */
502  if ( sepadata->separelint )
503  {
504  /* possibly free solution */
505  if ( sepadata->sepasol != NULL )
506  {
507  SCIP_CALL( SCIPfreeSol(scip, &sepadata->sepasol) );
508  }
509 
510  /* copy and store solution */
511  SCIP_CALL( SCIPcreateSolCopy(scip, &sepadata->sepasol, sol) );
512  sepadata->triedRelint = TRUE;
513  }
514 
515  return SCIP_OKAY;
516 }
#define SCIP_DEFAULT_INCLOBJCUTOFF
SCIP_Real SCIPgetSolvingTime(SCIP *scip)
Definition: scip_timing.c:436
SCIP_Longint SCIPgetNRootLPIterations(SCIP *scip)
#define NULL
Definition: def.h:239
SCIP_RETCODE SCIPcomputeLPRelIntPoint(SCIP *scip, SCIP_Bool relaxrows, SCIP_Bool inclobjcutoff, SCIP_Real timelimit, int iterlimit, SCIP_SOL **point)
Definition: scip_lp.c:1073
public methods for SCIP parameter handling
SCIP_NODE * SCIPgetCurrentNode(SCIP *scip)
Definition: scip_tree.c:158
public methods for branch and bound tree
SCIP_Longint SCIPgetNLPIterations(SCIP *scip)
public methods for memory management
SCIP_RETCODE SCIPgetRealParam(SCIP *scip, const char *name, SCIP_Real *value)
Definition: scip_param.c:379
SCIP_SEPA * SCIPfindSepa(SCIP *scip, const char *name)
Definition: scip_sepa.c:316
SCIP_Real SCIPvarGetLbLocal(SCIP_VAR *var)
Definition: var.c:17399
public methods for timing
#define FALSE
Definition: def.h:65
#define TRUE
Definition: def.h:64
const char * SCIPsepaGetName(SCIP_SEPA *sepa)
Definition: sepa.c:689
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)
public methods for problem variables
#define SCIPfreeBlockMemory(scip, ptr)
Definition: scip_mem.h:114
#define SCIPallocBlockMemory(scip, ptr)
Definition: scip_mem.h:97
int SCIPgetNLPBranchCands(SCIP *scip)
Definition: scip_branch.c:417
SCIP_RETCODE SCIPsetSepaCopy(SCIP *scip, SCIP_SEPA *sepa, SCIP_DECL_SEPACOPY((*sepacopy)))
Definition: scip_sepa.c:220
#define SCIPdebugMsg
Definition: scip_message.h:88
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_param.c:155
public methods for separator plugins
public methods for numerical tolerances
SCIP_SEPADATA * SCIPsepaGetData(SCIP_SEPA *sepa)
Definition: sepa.c:600
public methods for querying solving statistics
#define SCIP_DEFAULT_RECOMPUTERELINT
public methods for the branch-and-bound tree
SCIP_Longint SCIPnodeGetNumber(SCIP_NODE *node)
Definition: tree.c:7344
SCIP_RETCODE SCIPcreateSolCopy(SCIP *scip, SCIP_SOL **sol, SCIP_SOL *sourcesol)
Definition: scip_sol.c:667
#define SCIPerrorMessage
Definition: pub_message.h:45
closecuts meta separator
void SCIPsepaSetData(SCIP_SEPA *sepa, SCIP_SEPADATA *sepadata)
Definition: sepa.c:610
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_cut.c:778
SCIP_Real SCIPvarGetLPSol(SCIP_VAR *var)
Definition: var.c:17717
#define SEPA_NAME
#define SCIP_CALL(x)
Definition: def.h:351
SCIP_RETCODE SCIPremoveInefficaciousCuts(SCIP *scip)
Definition: scip_cut.c:866
void SCIPverbMessage(SCIP *scip, SCIP_VERBLEVEL msgverblevel, FILE *file, const char *formatstr,...)
Definition: scip_message.c:296
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_sepa.c:178
static SCIP_DECL_SEPAEXITSOL(sepaExitsolClosecuts)
SCIP_RETCODE SCIPsetSolVal(SCIP *scip, SCIP_SOL *sol, SCIP_VAR *var, SCIP_Real val)
Definition: scip_sol.c:1270
#define SCIP_DEFAULT_MAXLPITERFACTOR
#define SCIP_DEFAULT_SEPATHRESHOLD
SCIP_RETCODE SCIPsetSepaExitsol(SCIP *scip, SCIP_SEPA *sepa, SCIP_DECL_SEPAEXITSOL((*sepaexitsol)))
Definition: scip_sepa.c:300
#define SCIP_Bool
Definition: def.h:62
SCIP_LPSOLSTAT SCIPgetLPSolstat(SCIP *scip)
Definition: scip_lp.c:226
int SCIPgetDepth(SCIP *scip)
Definition: scip_tree.c:715
#define MIN(x, y)
Definition: def.h:209
public methods for cuts and aggregation rows
SCIP_RETCODE SCIPfreeSol(SCIP *scip, SCIP_SOL **sol)
Definition: scip_sol.c:1034
#define SCIP_DEFAULT_SEPARELINT
SCIP_Bool SCIPisInfinity(SCIP *scip, SCIP_Real val)
#define SEPA_PRIORITY
public methods for the LP relaxation, rows and columns
int SCIPgetNVars(SCIP *scip)
Definition: scip_prob.c:2044
SCIP_RETCODE SCIPincludeSepaClosecuts(SCIP *scip)
#define SCIP_REAL_MAX
Definition: def.h:151
public methods for branching rule plugins and branching
general public methods
#define MAX(x, y)
Definition: def.h:208
SCIP_RETCODE SCIPsetSepaFree(SCIP *scip, SCIP_SEPA *sepa, SCIP_DECL_SEPAFREE((*sepafree)))
Definition: scip_sepa.c:236
SCIP_SOL * SCIPgetBestSol(SCIP *scip)
Definition: scip_sol.c:2379
public methods for solutions
#define SEPA_FREQ
int SCIPgetNCuts(SCIP *scip)
Definition: scip_cut.c:830
public methods for message output
SCIP_VAR ** SCIPgetVars(SCIP *scip)
Definition: scip_prob.c:1999
#define SCIP_DEFAULT_MAXUNSUCCESSFUL
#define SCIP_Real
Definition: def.h:150
SCIP_Bool SCIPisStopped(SCIP *scip)
Definition: scip_general.c:739
#define SEPA_DESC
public methods for message handling
#define SEPA_USESSUBSCIP
#define SCIP_Longint
Definition: def.h:135
static SCIP_DECL_SEPAFREE(sepaFreeClosecuts)
#define SEPA_MAXBOUNDDIST
SCIP_Bool SCIPisZero(SCIP *scip, SCIP_Real val)
SCIP_Real SCIPvarGetUbLocal(SCIP_VAR *var)
Definition: var.c:17409
public methods for separators
#define SCIP_MIN_LPITERS
public methods for global and local (sub)problems
SCIP_Real SCIPgetSolVal(SCIP *scip, SCIP_SOL *sol, SCIP_VAR *var)
Definition: scip_sol.c:1410
#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_param.c:211
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_param.c:129
static SCIP_DECL_SEPAEXECLP(sepaExeclpClosecuts)
SCIP_RETCODE SCIPcreateSol(SCIP *scip, SCIP_SOL **sol, SCIP_HEUR *heur)
Definition: scip_sol.c:377