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-2024 Zuse Institute Berlin (ZIB) */
7/* */
8/* Licensed under the Apache License, Version 2.0 (the "License"); */
9/* you may not use this file except in compliance with the License. */
10/* You may obtain a copy of the License at */
11/* */
12/* http://www.apache.org/licenses/LICENSE-2.0 */
13/* */
14/* Unless required by applicable law or agreed to in writing, software */
15/* distributed under the License is distributed on an "AS IS" BASIS, */
16/* WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. */
17/* See the License for the specific language governing permissions and */
18/* limitations under the License. */
19/* */
20/* You should have received a copy of the Apache-2.0 license */
21/* along with SCIP; see the file LICENSE. If not visit scipopt.org. */
22/* */
23/* * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * */
24
25/**@file sepa_closecuts.c
26 * @ingroup DEFPLUGINS_SEPA
27 * @brief closecuts meta separator
28 * @author Marc Pfetsch
29 *
30 * This separator generates a convex combination of the current LP solution and either the best
31 * primal feasible solution or an interior point of the LP relaxation. If the convex combination is
32 * proper, the new point is closer to the convex hull of the feasible points. The separator then
33 * calls all other separators to separate this point. The idea is that in this way possibly "deeper"
34 * cuts are generated. Note, however, that the new point is not a basic solution, i.e., separators
35 * relying basis information, e.g., Gomory cuts, will not work.
36 *
37 * The other cuts are generated via the sepasol() callbacks in constraints handlers or separators.
38 *
39 * This separator stops after a certain number (parameter @p maxunsuccessful) of unsuccessful
40 * calls. It also inhibits the separation of the ordinary LP solution if it already generated enough
41 * (parameter @p sepathreshold) cuts. The convex combination is determined via the parameter @p
42 * sepacombvalue.
43 *
44 * In general, this separator makes sense if it is expected that there will be many separation
45 * rounds and many cuts will be again deleted, because they are not active after a certain number of
46 * rounds. In particular, branch-and-cut algorithms for combinatorial optimization problems form
47 * good candidates.
48 *
49 * The idea seems to be first proposed in the context of the travelling salesman problem, see@par
50 * The Traveling Salesman Problem: A Computational Study@n
51 * David L. Applegate, Robert E. Bixby, Vasek Chvatal & William J. Cook@n
52 * Princeton University Press 2006@n
53 *
54 * for more details. See also@par
55 * Acceleration of cutting-plane and column generation algorithms: Applications to network design.@n
56 * Walid Ben-Ameur, Jose Neto@n
57 * Networks 49(1): 3-17 (2007).
58 */
59
60/*---+----1----+----2----+----3----+----4----+----5----+----6----+----7----+----8----+----9----+----0----+----1----+----2*/
61
62#include "scip/pub_message.h"
63#include "scip/pub_sepa.h"
64#include "scip/pub_tree.h"
65#include "scip/pub_var.h"
66#include "scip/scip_branch.h"
67#include "scip/scip_cut.h"
68#include "scip/scip_general.h"
69#include "scip/scip_lp.h"
70#include "scip/scip_mem.h"
71#include "scip/scip_message.h"
72#include "scip/scip_numerics.h"
73#include "scip/scip_param.h"
74#include "scip/scip_prob.h"
75#include "scip/scip_sepa.h"
76#include "scip/scip_sol.h"
78#include "scip/scip_timing.h"
79#include "scip/scip_tree.h"
80#include "scip/sepa_closecuts.h"
81#include <string.h>
82
83
84
85#define SEPA_NAME "closecuts"
86#define SEPA_DESC "closecuts meta separator"
87#define SEPA_PRIORITY 1000000
88#define SEPA_FREQ -1
89#define SEPA_MAXBOUNDDIST 1.0
90#define SEPA_USESSUBSCIP FALSE /**< does the separator use a secondary SCIP instance? */
91#define SEPA_DELAY FALSE /**< should separation method be delayed, if other separators found cuts? */
92
93
94/* default values for parameters */
95#define SCIP_DEFAULT_SEPARELINT TRUE /**< generate close cuts w.r.t. relative interior point (best solution otherwise)? */
96#define SCIP_DEFAULT_SEPACOMBVALUE 0.30 /**< convex combination value for close cuts */
97#define SCIP_DEFAULT_SEPATHRESHOLD 50 /**< threshold on number of generated cuts below which the ordinary separation is started */
98#define SCIP_DEFAULT_INCLOBJCUTOFF FALSE /**< include the objective cutoff when computing the relative interior? */
99#define SCIP_DEFAULT_RECOMPUTERELINT FALSE /**< recompute relative interior in each separation call? */
100#define SCIP_DEFAULT_MAXUNSUCCESSFUL 0 /**< turn off separation in current node after unsuccessful calls (-1 never turn off) */
101#define SCIP_DEFAULT_MAXLPITERFACTOR 10.0 /**< factor for maximal LP iterations in relative interior computation compared to node LP iterations */
102
103#define SCIP_MIN_LPITERS 100 /**< minimum number of allowed LP iterations in relative interior computation */
104
105
106/** separator data */
107struct SCIP_SepaData
108{
109 SCIP_Bool separelint; /**< generate close cuts w.r.t. relative interior point (best solution otherwise)? */
110 SCIP_Bool triedRelint; /**< tried to compute relative interior */
111 SCIP_Real sepacombvalue; /**< convex combination value for close cuts */
112 int sepathreshold; /**< threshold on number of generated cuts below which the ordinary separation is started */
113 SCIP_Bool inclobjcutoff; /**< include the objective cutoff when computing the relative interior? */
114 SCIP_Bool recomputerelint; /**< recompute relative interior in each separation call? */
115 int maxunsuccessful; /**< turn off separation in current node after unsuccessful calls (-1 never turn off) */
116 SCIP_SOL* sepasol; /**< solution that can be used for generating close cuts */
117 SCIP_Longint discardnode; /**< number of node for which separation is discarded */
118 SCIP_Real maxlpiterfactor; /**< factor for maximal LP iterations in relative interior computation compared to node LP iterations */
119 int nunsuccessful; /**< number of consecutive unsuccessful calls */
120};
121
122
123/** generate point for close cut separation
124 *
125 * The constructed point is the convex combination of the point stored in set->closesol and the
126 * current LP solution. The convexity parameter is set->sepa_closecombvalue. If this parameter is
127 * 0, the point coincides with the LP solution.
128 */
129static
131 SCIP* scip, /**< SCIP data structure */
132 SCIP_SEPADATA* sepadata, /**< separator data */
133 SCIP_SOL** point /**< point to be generated (or NULL if unsuccessful) */
134 )
135{
136 SCIP_VAR** vars;
137 SCIP_VAR* var;
138 SCIP_Real val;
139 SCIP_Real alpha;
140 SCIP_Real onealpha;
141 SCIP_Real lb;
142 SCIP_Real ub;
143 int nvars;
144 int i;
145
146 assert( scip != NULL );
147 assert( point != NULL );
148
149 *point = NULL;
150 if ( sepadata->sepasol == NULL )
151 return SCIP_OKAY;
152
153 alpha = sepadata->sepacombvalue;
154 if ( alpha < 0.001 )
155 return SCIP_OKAY;
156 onealpha = 1.0 - alpha;
157
158 /* create solution */
159 SCIP_CALL( SCIPcreateSol(scip, point, NULL) );
160
161 /* generate convex combination */
162 vars = SCIPgetVars(scip);
163 nvars = SCIPgetNVars(scip);
164 for (i = 0; i < nvars; ++i)
165 {
166 var = vars[i];
167 val = alpha * SCIPgetSolVal(scip, sepadata->sepasol, var) + onealpha * SCIPvarGetLPSol(var);
168
169 /* If both the LP relaxation and the base point respect the variable bounds, the computed point will satisfy them
170 * as well. However, variables might be fixed (e.g. by branching) since the time of the computation of the base
171 * point. Thus, we adapt the value to lie inside the bounds in optimized mode. */
172 lb = SCIPvarGetLbLocal(var);
173 ub = SCIPvarGetUbLocal(var);
174 val = MAX(val, lb);
175 val = MIN(val, ub);
176
177 if ( ! SCIPisZero(scip, val) )
178 {
179 SCIP_CALL( SCIPsetSolVal(scip, *point, var, val) );
180 }
181 }
182
183 return SCIP_OKAY;
184}
185
186
187/*
188 * Callback methods of separator
189 */
190
191
192/** copy method for separator plugins (called when SCIP copies plugins) */
193static
194SCIP_DECL_SEPACOPY(sepaCopyClosecuts)
195{ /*lint --e{715}*/
196 assert( scip != NULL );
197 assert( sepa != NULL );
198 assert( strcmp(SCIPsepaGetName(sepa), SEPA_NAME) == 0 );
199
200 /* call inclusion method of constraint handler */
202
203 return SCIP_OKAY;
204}
205
206/** destructor of separator to free user data (called when SCIP is exiting) */
207static
208SCIP_DECL_SEPAFREE(sepaFreeClosecuts)
209{ /*lint --e{715}*/
210 SCIP_SEPADATA* sepadata;
211
212 assert( sepa != NULL );
213 assert( strcmp(SCIPsepaGetName(sepa), SEPA_NAME) == 0 );
214
215 /* free separator data */
216 sepadata = SCIPsepaGetData(sepa);
217 assert( sepadata != NULL );
218
219 SCIPfreeBlockMemory(scip, &sepadata);
220
221 SCIPsepaSetData(sepa, NULL);
222
223 return SCIP_OKAY;
224}
225
226
227/** solving process deinitialization method of separator (called before branch and bound process data is freed) */
228static
229SCIP_DECL_SEPAEXITSOL(sepaExitsolClosecuts)
230{ /*lint --e{715}*/
231 SCIP_SEPADATA* sepadata;
232
233 assert( sepa != NULL );
234 assert( strcmp(SCIPsepaGetName(sepa), SEPA_NAME) == 0 );
235
236 sepadata = SCIPsepaGetData(sepa);
237 assert( sepadata != NULL );
238
239 if ( sepadata->separelint && sepadata->sepasol != NULL )
240 {
241 SCIP_CALL( SCIPfreeSol(scip, &sepadata->sepasol) );
242 sepadata->triedRelint = FALSE;
243 }
244 sepadata->discardnode = -1;
245 sepadata->nunsuccessful = 0;
246
247 return SCIP_OKAY;
248}
249
250
251/** LP solution separation method of separator */
252static
253SCIP_DECL_SEPAEXECLP(sepaExeclpClosecuts)
254{ /*lint --e{715}*/
255 SCIP_SEPADATA* sepadata;
256 SCIP_Longint currentnodenumber;
257 SCIP_SOL* point = NULL;
258 SCIP_Bool isroot;
259
260 assert( sepa != NULL );
261 assert( strcmp(SCIPsepaGetName(sepa), SEPA_NAME) == 0 );
262 assert( result != NULL );
263
264 *result = SCIP_DIDNOTRUN;
265
266 /* only call separator, if LP has been solved (need LP to compute separation point) */
268 return SCIP_OKAY;
269
270 /* only call separator, if there are fractional variables */
271 if ( SCIPgetNLPBranchCands(scip) == 0 )
272 return SCIP_OKAY;
273
274 /* exit if we stopped ... */
275 if ( SCIPisStopped(scip) )
276 return SCIP_OKAY;
277
278 /* get separation data */
279 sepadata = SCIPsepaGetData(sepa);
280 assert( sepadata != NULL );
281
282 /* exit if we already decided to discard the current node */
283 currentnodenumber = SCIPnodeGetNumber(SCIPgetCurrentNode(scip));
284 if ( sepadata->discardnode == currentnodenumber )
285 return SCIP_OKAY;
286
287 SCIPdebugMsg(scip, "Separation method of closecuts separator.\n");
288
289 /* check whether we have to compute a relative interior point */
290 if ( sepadata->separelint )
291 {
292 if ( sepadata->recomputerelint )
293 {
294 /* check if previous relative interior point should be forgotten, otherwise it is computed only once and the
295 * same point is used for all nodes */
296 if ( sepadata->sepasol != NULL )
297 {
298 SCIP_CALL( SCIPfreeSol(scip, &sepadata->sepasol) );
299 sepadata->triedRelint = FALSE;
300 }
301 }
302 else
303 {
304 /* skip execution, if we unsuccessfully tried to compute a relative interior point */
305 if ( sepadata->sepasol == NULL && sepadata->triedRelint )
306 return SCIP_OKAY;
307 }
308
309 /* if relative interior point is not available ... */
310 if ( sepadata->sepasol == NULL )
311 {
312 SCIP_Longint nlpiters;
313 SCIP_Real timelimit;
314 int iterlimit;
315
316 /* prepare time limit */
317 SCIP_CALL( SCIPgetRealParam(scip, "limits/time", &timelimit) );
318 if ( ! SCIPisInfinity(scip, timelimit) )
319 timelimit -= SCIPgetSolvingTime(scip);
320 /* exit if no time left */
321 if ( timelimit <= 0.0 )
322 return SCIP_OKAY;
323
324 /* determine iteration limit */
325 if ( sepadata->maxlpiterfactor < 0.0 || SCIPisInfinity(scip, sepadata->maxlpiterfactor) )
326 iterlimit = INT_MAX;
327 else
328 {
329 /* determine iteration limit; the number of iterations in the root is only set after its solution, but the
330 * total number of LP iterations is always updated.
331 * here we use SCIPgetDepth instead of the depth argument passed to the callback because if we are not in
332 * the root node but depth is 0 (i.e. if we want us to behave as if we are in the root node regarding
333 * limits) then using the total number of iterations so far is a gross overestimation
334 */
335 if ( SCIPgetDepth(scip) == 0 )
336 nlpiters = SCIPgetNLPIterations(scip);
337 else
338 nlpiters = SCIPgetNRootLPIterations(scip);
339 iterlimit = (int)(sepadata->maxlpiterfactor * nlpiters);
340 iterlimit = MAX(iterlimit, SCIP_MIN_LPITERS);
341 assert(iterlimit > 0);
342 }
343
344 SCIPverbMessage(scip, SCIP_VERBLEVEL_MINIMAL, 0, "Computing relative interior point (time limit: %g, iter limit: %d) ...\n", timelimit, iterlimit);
345 SCIP_CALL( SCIPcomputeLPRelIntPoint(scip, TRUE, sepadata->inclobjcutoff, timelimit, iterlimit, &sepadata->sepasol) );
346 sepadata->triedRelint = TRUE;
347 }
348 }
349 else
350 {
351 /* get best solution (NULL if not present) */
352 sepadata->sepasol = SCIPgetBestSol(scip);
353 }
354
355 /* separate close cuts */
356 if ( sepadata->sepasol != NULL )
357 {
358 SCIPdebugMsg(scip, "Generating close cuts ... (combination value: %f)\n", sepadata->sepacombvalue);
359 *result = SCIP_DIDNOTFIND;
360
361 /* generate point to be separated */
362 SCIP_CALL( generateCloseCutPoint(scip, sepadata, &point) );
363
364 /* apply a separation round to generated point */
365 if ( point != NULL )
366 {
367 int noldcuts;
368 SCIP_Bool delayed;
369 SCIP_Bool cutoff;
370
371 noldcuts = SCIPgetNCuts(scip);
372 isroot = (SCIP_Bool) (depth == 0);
373
374 /* separate solution via other separators */
375 SCIP_CALL( SCIPseparateSol(scip, point, isroot, TRUE, FALSE, &delayed, &cutoff) );
376
377 SCIP_CALL( SCIPfreeSol(scip, &point) );
378 assert( point == NULL );
379
380 /* the cuts might not violated by the current LP if the computed point is strange */
382
383 if ( cutoff )
384 *result = SCIP_CUTOFF;
385 else
386 {
387 if ( SCIPgetNCuts(scip) - noldcuts > sepadata->sepathreshold )
388 {
389 sepadata->nunsuccessful = 0;
390 *result = SCIP_NEWROUND;
391 }
392 else
393 {
394 if ( SCIPgetNCuts(scip) > noldcuts )
395 {
396 sepadata->nunsuccessful = 0;
397 *result = SCIP_SEPARATED;
398 }
399 else
400 ++sepadata->nunsuccessful;
401 }
402 }
403
404 SCIPdebugMsg(scip, "Separated close cuts: %d (enoughcuts: %d, unsuccessful: %d).\n", SCIPgetNCuts(scip) - noldcuts,
405 SCIPgetNCuts(scip) - noldcuts > sepadata->sepathreshold, sepadata->nunsuccessful);
406
407 if ( sepadata->maxunsuccessful >= 0 && sepadata->nunsuccessful > sepadata->maxunsuccessful )
408 {
409 SCIPdebugMsg(scip, "Turn off close cut separation, because of %d unsuccessful calls.\n", sepadata->nunsuccessful);
410 sepadata->discardnode = currentnodenumber;
411 sepadata->nunsuccessful = 0;
412 }
413 }
414 }
415
416 return SCIP_OKAY;
417}
418
419
420/*
421 * separator specific interface methods
422 */
423
424/** creates the closecuts separator and includes it in SCIP */
426 SCIP* scip /**< SCIP data structure */
427 )
428{
429 SCIP_SEPADATA* sepadata;
430 SCIP_SEPA* sepa;
431
432 /* create closecuts separator data */
433 SCIP_CALL( SCIPallocBlockMemory(scip, &sepadata) );
434 sepadata->sepasol = NULL;
435 sepadata->discardnode = -1;
436 sepadata->nunsuccessful = 0;
437 sepadata->triedRelint = FALSE;
438
439 /* include separator */
441 sepaExeclpClosecuts, NULL, sepadata) );
442
443 assert(sepa != NULL);
444
445 /* set non-NULL pointers to callback methods */
446 SCIP_CALL( SCIPsetSepaCopy(scip, sepa, sepaCopyClosecuts) );
447 SCIP_CALL( SCIPsetSepaFree(scip, sepa, sepaFreeClosecuts) );
448 SCIP_CALL( SCIPsetSepaExitsol(scip, sepa, sepaExitsolClosecuts) );
449
450 /* add closecuts separator parameters */
452 "separating/closecuts/separelint",
453 "generate close cuts w.r.t. relative interior point (best solution otherwise)?",
454 &sepadata->separelint, TRUE, SCIP_DEFAULT_SEPARELINT, NULL, NULL) );
455
457 "separating/closecuts/sepacombvalue",
458 "convex combination value for close cuts",
459 &sepadata->sepacombvalue, TRUE, SCIP_DEFAULT_SEPACOMBVALUE, 0.0, 1.0,
460 NULL, NULL) );
461
463 "separating/closecuts/closethres",
464 "threshold on number of generated cuts below which the ordinary separation is started",
465 &sepadata->sepathreshold, TRUE, SCIP_DEFAULT_SEPATHRESHOLD, -1, INT_MAX, NULL, NULL) );
466
468 "separating/closecuts/inclobjcutoff",
469 "include an objective cutoff when computing the relative interior?",
470 &sepadata->inclobjcutoff, TRUE, SCIP_DEFAULT_INCLOBJCUTOFF, NULL, NULL) );
471
473 "separating/closecuts/recomputerelint",
474 "recompute relative interior point in each separation call?",
475 &sepadata->recomputerelint, TRUE, SCIP_DEFAULT_RECOMPUTERELINT, NULL, NULL) );
476
478 "separating/closecuts/maxunsuccessful",
479 "turn off separation in current node after unsuccessful calls (-1 never turn off)",
480 &sepadata->maxunsuccessful, TRUE, SCIP_DEFAULT_MAXUNSUCCESSFUL, -1, INT_MAX, NULL, NULL) );
481
483 "separating/closecuts/maxlpiterfactor",
484 "factor for maximal LP iterations in relative interior computation compared to node LP iterations (negative for no limit)",
485 &sepadata->maxlpiterfactor, TRUE, SCIP_DEFAULT_MAXLPITERFACTOR, -1.0, SCIP_REAL_MAX, NULL, NULL) );
486
487 return SCIP_OKAY;
488}
489
490/** sets point to be used as base point for computing the point to be separated
491 *
492 * The point is only stored if separation of relative interior points is used. The solution is copied.
493 */
495 SCIP* scip, /**< SCIP data structure */
496 SCIP_SOL* sol /**< base point solution */
497 )
498{
499 SCIP_SEPA* sepa;
500 SCIP_SEPADATA* sepadata;
501
502 assert( scip != NULL );
503
504 /* find separator */
505 sepa = SCIPfindSepa(scip, SEPA_NAME);
506 if ( sepa == NULL )
507 {
508 SCIPerrorMessage("Could not find separator <%s>.\n", SEPA_NAME);
509 return SCIP_PLUGINNOTFOUND;
510 }
511 assert( strcmp(SCIPsepaGetName(sepa), SEPA_NAME) == 0 );
512
513 /* get sepadata */
514 sepadata = SCIPsepaGetData(sepa);
515 assert( sepadata != NULL );
516
517 /* store point if we have to separate relative interior points */
518 if ( sepadata->separelint )
519 {
520 /* possibly free solution */
521 if ( sepadata->sepasol != NULL )
522 {
523 SCIP_CALL( SCIPfreeSol(scip, &sepadata->sepasol) );
524 }
525
526 /* copy and store solution */
527 SCIP_CALL( SCIPcreateSolCopy(scip, &sepadata->sepasol, sol) );
528 sepadata->triedRelint = TRUE;
529 }
530
531 return SCIP_OKAY;
532}
#define NULL
Definition: def.h:267
#define SCIP_Longint
Definition: def.h:158
#define SCIP_REAL_MAX
Definition: def.h:174
#define SCIP_Bool
Definition: def.h:91
#define MIN(x, y)
Definition: def.h:243
#define SCIP_Real
Definition: def.h:173
#define TRUE
Definition: def.h:93
#define FALSE
Definition: def.h:94
#define MAX(x, y)
Definition: def.h:239
#define SCIP_CALL(x)
Definition: def.h:374
SCIP_Bool SCIPisStopped(SCIP *scip)
Definition: scip_general.c:724
int SCIPgetNVars(SCIP *scip)
Definition: scip_prob.c:1992
SCIP_VAR ** SCIPgetVars(SCIP *scip)
Definition: scip_prob.c:1947
void SCIPverbMessage(SCIP *scip, SCIP_VERBLEVEL msgverblevel, FILE *file, const char *formatstr,...)
Definition: scip_message.c:225
#define SCIPdebugMsg
Definition: scip_message.h:78
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:83
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:139
SCIP_RETCODE SCIPgetRealParam(SCIP *scip, const char *name, SCIP_Real *value)
Definition: scip_param.c:307
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:57
int SCIPgetNLPBranchCands(SCIP *scip)
Definition: scip_branch.c:428
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:735
SCIP_RETCODE SCIPremoveInefficaciousCuts(SCIP *scip)
Definition: scip_cut.c:823
int SCIPgetNCuts(SCIP *scip)
Definition: scip_cut.c:787
SCIP_RETCODE SCIPcomputeLPRelIntPoint(SCIP *scip, SCIP_Bool relaxrows, SCIP_Bool inclobjcutoff, SCIP_Real timelimit, int iterlimit, SCIP_SOL **point)
Definition: scip_lp.c:1094
SCIP_LPSOLSTAT SCIPgetLPSolstat(SCIP *scip)
Definition: scip_lp.c:168
#define SCIPfreeBlockMemory(scip, ptr)
Definition: scip_mem.h:108
#define SCIPallocBlockMemory(scip, ptr)
Definition: scip_mem.h:89
SCIP_Longint SCIPnodeGetNumber(SCIP_NODE *node)
Definition: tree.c:7493
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:109
const char * SCIPsepaGetName(SCIP_SEPA *sepa)
Definition: sepa.c:743
SCIP_RETCODE SCIPsetSepaFree(SCIP *scip, SCIP_SEPA *sepa, SCIP_DECL_SEPAFREE((*sepafree)))
Definition: scip_sepa.c:167
SCIP_SEPA * SCIPfindSepa(SCIP *scip, const char *name)
Definition: scip_sepa.c:247
SCIP_RETCODE SCIPsetSepaExitsol(SCIP *scip, SCIP_SEPA *sepa, SCIP_DECL_SEPAEXITSOL((*sepaexitsol)))
Definition: scip_sepa.c:231
SCIP_SEPADATA * SCIPsepaGetData(SCIP_SEPA *sepa)
Definition: sepa.c:633
void SCIPsepaSetData(SCIP_SEPA *sepa, SCIP_SEPADATA *sepadata)
Definition: sepa.c:643
SCIP_RETCODE SCIPsetSepaCopy(SCIP *scip, SCIP_SEPA *sepa, SCIP_DECL_SEPACOPY((*sepacopy)))
Definition: scip_sepa.c:151
SCIP_SOL * SCIPgetBestSol(SCIP *scip)
Definition: scip_sol.c:2169
SCIP_RETCODE SCIPcreateSol(SCIP *scip, SCIP_SOL **sol, SCIP_HEUR *heur)
Definition: scip_sol.c:184
SCIP_RETCODE SCIPcreateSolCopy(SCIP *scip, SCIP_SOL **sol, SCIP_SOL *sourcesol)
Definition: scip_sol.c:474
SCIP_RETCODE SCIPfreeSol(SCIP *scip, SCIP_SOL **sol)
Definition: scip_sol.c:841
SCIP_RETCODE SCIPsetSolVal(SCIP *scip, SCIP_SOL *sol, SCIP_VAR *var, SCIP_Real val)
Definition: scip_sol.c:1077
SCIP_Real SCIPgetSolVal(SCIP *scip, SCIP_SOL *sol, SCIP_VAR *var)
Definition: scip_sol.c:1217
SCIP_Longint SCIPgetNRootLPIterations(SCIP *scip)
SCIP_Longint SCIPgetNLPIterations(SCIP *scip)
SCIP_Real SCIPgetSolvingTime(SCIP *scip)
Definition: scip_timing.c:378
SCIP_Bool SCIPisInfinity(SCIP *scip, SCIP_Real val)
SCIP_Bool SCIPisZero(SCIP *scip, SCIP_Real val)
int SCIPgetDepth(SCIP *scip)
Definition: scip_tree.c:670
SCIP_NODE * SCIPgetCurrentNode(SCIP *scip)
Definition: scip_tree.c:91
SCIP_Real SCIPvarGetUbLocal(SCIP_VAR *var)
Definition: var.c:18144
SCIP_Real SCIPvarGetLPSol(SCIP_VAR *var)
Definition: var.c:18452
SCIP_Real SCIPvarGetLbLocal(SCIP_VAR *var)
Definition: var.c:18134
SCIP_RETCODE SCIPsetBasePointClosecuts(SCIP *scip, SCIP_SOL *sol)
SCIP_RETCODE SCIPincludeSepaClosecuts(SCIP *scip)
public methods for message output
#define SCIPerrorMessage
Definition: pub_message.h:64
public methods for separators
public methods for branch and bound tree
public methods for problem variables
public methods for branching rule plugins and branching
public methods for cuts and aggregation rows
general public methods
public methods for the LP relaxation, rows and columns
public methods for memory management
public methods for message handling
public methods for numerical tolerances
public methods for SCIP parameter handling
public methods for global and local (sub)problems
public methods for separator plugins
public methods for solutions
public methods for querying solving statistics
public methods for timing
public methods for the branch-and-bound tree
#define SEPA_PRIORITY
#define SEPA_DELAY
#define SEPA_DESC
#define SCIP_DEFAULT_SEPACOMBVALUE
#define SCIP_DEFAULT_INCLOBJCUTOFF
#define SCIP_MIN_LPITERS
#define SEPA_USESSUBSCIP
#define SCIP_DEFAULT_RECOMPUTERELINT
#define SCIP_DEFAULT_SEPARELINT
static SCIP_DECL_SEPACOPY(sepaCopyClosecuts)
#define SEPA_MAXBOUNDDIST
#define SCIP_DEFAULT_MAXUNSUCCESSFUL
#define SEPA_FREQ
#define SEPA_NAME
#define SCIP_DEFAULT_MAXLPITERFACTOR
static SCIP_DECL_SEPAEXECLP(sepaExeclpClosecuts)
static SCIP_DECL_SEPAEXITSOL(sepaExitsolClosecuts)
static SCIP_RETCODE generateCloseCutPoint(SCIP *scip, SCIP_SEPADATA *sepadata, SCIP_SOL **point)
static SCIP_DECL_SEPAFREE(sepaFreeClosecuts)
#define SCIP_DEFAULT_SEPATHRESHOLD
closecuts meta separator
@ SCIP_LPSOLSTAT_OPTIMAL
Definition: type_lp.h:43
@ SCIP_VERBLEVEL_MINIMAL
Definition: type_message.h:54
@ SCIP_DIDNOTRUN
Definition: type_result.h:42
@ SCIP_CUTOFF
Definition: type_result.h:48
@ SCIP_DIDNOTFIND
Definition: type_result.h:44
@ SCIP_SEPARATED
Definition: type_result.h:49
@ SCIP_NEWROUND
Definition: type_result.h:50
@ SCIP_PLUGINNOTFOUND
Definition: type_retcode.h:54
@ SCIP_OKAY
Definition: type_retcode.h:42
enum SCIP_Retcode SCIP_RETCODE
Definition: type_retcode.h:63
struct SCIP_SepaData SCIP_SEPADATA
Definition: type_sepa.h:52