1 /* * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * */
2 /* */
3 /* This file is part of the program and library */
4 /* SCIP --- Solving Constraint Integer Programs */
5 /* */
6 /* Copyright (C) 2002-2021 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 scipopt.org. */
13 /* */
14 /* * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * */
15
16 /**@file sepa_gauge.c
17 * @ingroup DEFPLUGINS_SEPA
18 * @brief gauge separator
19 * @author Felipe Serrano
20 *
21 * @todo should separator only be run when SCIPallColsInLP is true?
22 * @todo add SCIPisStopped(scip) to the condition of time consuming loops
23 * @todo check if it makes sense to implement the copy callback
24 */
25
26 /*---+----1----+----2----+----3----+----4----+----5----+----6----+----7----+----8----+----9----+----0----+----1----+----2*/
27
28 #include <assert.h>
29 #include <string.h>
30
31 #include "blockmemshell/memory.h"
32 #include "nlpi/exprinterpret.h"
33 #include "nlpi/nlpi.h"
34 #include "nlpi/nlpi_ipopt.h"
35 #include "nlpi/nlpioracle.h"
36 #include "nlpi/pub_expr.h"
37 #include "scip/pub_lp.h"
38 #include "scip/pub_message.h"
39 #include "scip/pub_misc.h"
40 #include "scip/pub_nlp.h"
41 #include "scip/pub_sepa.h"
42 #include "scip/scip_cut.h"
43 #include "scip/scip_lp.h"
44 #include "scip/scip_mem.h"
45 #include "scip/scip_message.h"
46 #include "scip/scip_nlp.h"
47 #include "scip/scip_nonlinear.h"
48 #include "scip/scip_numerics.h"
49 #include "scip/scip_param.h"
50 #include "scip/scip_prob.h"
51 #include "scip/scip_sepa.h"
52 #include "scip/scip_sol.h"
53 #include "scip/scip_solvingstats.h"
54 #include "scip/scip_timing.h"
55 #include "scip/sepa_gauge.h"
56 #include <string.h>
57
58
59 #define SEPA_NAME "gauge"
60 #define SEPA_DESC "gauge separator"
61 #define SEPA_PRIORITY 0
62 #define SEPA_FREQ -1
63 #define SEPA_MAXBOUNDDIST 1.0
64 #define SEPA_USESSUBSCIP FALSE /**< does the separator use a secondary SCIP instance? */
65 #define SEPA_DELAY FALSE /**< should separation method be delayed, if other separators found cuts? */
66
67 #define VIOLATIONFAC 100 /* constraints regarded as violated when violation > VIOLATIONFAC*SCIPfeastol */
68 #define MAX_ITER 75 /* maximum number of iterations for the line search */
69
70 #define DEFAULT_NLPTIMELIMIT 0.0 /**< default time limit of NLP solver; 0.0 for no limit */
71 #define DEFAULT_NLPITERLIM 1000 /**< default NLP iteration limit */
72
73 #define NLPFEASFAC 1e-1/**< NLP feasibility tolerance = NLPFEASFAC * SCIP's feasibility tolerance */
74 #define NLPVERBOSITY 0 /**< NLP solver verbosity */
75
76 #define INTERIOROBJVARLB -100 /**< lower bound of the objective variable when computing interior point */
77 /*
78 * Data structures
79 */
80
81 /** side that makes a nlrow convex */
82 enum ConvexSide
83 {
84 LHS = 0, /**< left hand side */
85 RHS = 1 /**< right hand side */
86 };
87 typedef enum ConvexSide CONVEXSIDE;
88
89 /** position of a point */
90 enum Position
91 {
92 INTERIOR = 0, /**< point is in the interior of the region */
93 BOUNDARY = 1, /**< point is in the boundary of the region */
94 EXTERIOR = 2 /**< point is in the exterior of the region */
95 };
96 typedef enum Position POSITION;
97
98 /** separator data */
99 struct SCIP_SepaData
100 {
101 SCIP_NLROW** nlrows; /**< stores convex nlrows */
102 CONVEXSIDE* convexsides; /**< which sides make the nlrows convex */
103 int* nlrowsidx; /**< indices of nlrows that violate the current lp solution */
104 int nnlrowsidx; /**< total number of convex nonlinear nlrows that violate the current lp solution */
105 int nnlrows; /**< total number of convex nonlinear nlrows */
106 int nlrowssize; /**< memory allocated for nlrows, convexsides and nlrowsidx */
107
108 SCIP_Bool isintsolavailable; /**< do we have an interior point available? */
109 SCIP_Bool skipsepa; /**< whether separator should be skipped */
110 SCIP_SOL* intsol; /**< stores interior point */
111
112 SCIP_EXPRINT* exprinterpreter; /**< expression interpreter to compute gradients */
113
114 int ncuts; /**< number of cuts generated */
115
116 /* parameters */
117 SCIP_Real nlptimelimit; /**< time limit of NLP solver; 0.0 for no limit */
118 int nlpiterlimit; /**< iteration limit of NLP solver; 0 for no limit */
119 };
120
121 /*
122 * Local methods
123 */
124
125 /** stores, from the constraints represented by nlrows, the nonlinear convex ones in sepadata */
126 static
storeNonlinearConvexNlrows(SCIP * scip,SCIP_SEPADATA * sepadata,SCIP_NLROW ** nlrows,int nnlrows)127 SCIP_RETCODE storeNonlinearConvexNlrows(
128 SCIP* scip, /**< SCIP data structure */
129 SCIP_SEPADATA* sepadata, /**< separator data */
130 SCIP_NLROW** nlrows, /**< nlrows from which to store convex ones */
131 int nnlrows /**< number of nlrows */
132 )
133 {
134 int i;
135
136 assert(scip != NULL);
137 assert(sepadata != NULL);
138 assert(nlrows != NULL);
139 assert(nnlrows > 0);
140
141 SCIPdebugMsg(scip, "storing convex nlrows\n");
142
143 SCIP_CALL( SCIPallocBlockMemoryArray(scip, &(sepadata->nlrows), nnlrows) );
144 SCIP_CALL( SCIPallocBlockMemoryArray(scip, &(sepadata->convexsides), nnlrows) );
145 SCIP_CALL( SCIPallocBlockMemoryArray(scip, &(sepadata->nlrowsidx), nnlrows) );
146 sepadata->nlrowssize = nnlrows;
147
148 sepadata->nnlrows = 0;
149 for( i = 0; i < nnlrows; ++i )
150 {
151 SCIP_NLROW* nlrow;
152
153 nlrow = nlrows[i];
154 assert(nlrow != NULL);
155
156 /* linear case */
157 if( SCIPnlrowGetCurvature(nlrow) == SCIP_EXPRCURV_LINEAR ||
158 (SCIPnlrowGetNQuadElems(nlrow) == 0 && SCIPnlrowGetExprtree(nlrow) == NULL) )
159 continue;
160
161 /* nonlinear case */
162 if( !SCIPisInfinity(scip, SCIPnlrowGetRhs(nlrow)) && SCIPnlrowGetCurvature(nlrow) == SCIP_EXPRCURV_CONVEX )
163 {
164 sepadata->convexsides[sepadata->nnlrows] = RHS;
165 sepadata->nlrows[sepadata->nnlrows] = nlrow;
166 ++(sepadata->nnlrows);
167 }
168 else if( !SCIPisInfinity(scip, -SCIPnlrowGetLhs(nlrow)) && SCIPnlrowGetCurvature(nlrow) == SCIP_EXPRCURV_CONCAVE )
169 {
170 sepadata->convexsides[sepadata->nnlrows] = LHS;
171 sepadata->nlrows[sepadata->nnlrows] = nlrow;
172 ++(sepadata->nnlrows);
173 }
174 }
175
176 return SCIP_OKAY;
177 }
178
179 /** computes an interior point of a convex NLP relaxation; builds the convex relaxation, modifies it to find an interior
180 * point, solves it and frees it; more details in @ref sepa_gauge.h
181 *
182 * @note the method also counts the number of nonlinear convex constraints and if there are < 2, then the convex
183 * relaxation is not interesting and the separator will not run again
184 */
185 static
computeInteriorPoint(SCIP * scip,SCIP_SEPADATA * sepadata)186 SCIP_RETCODE computeInteriorPoint(
187 SCIP* scip, /**< SCIP data structure */
188 SCIP_SEPADATA* sepadata /**< separator data */
189 )
190 {
191 SCIP_NLPIORACLE* nlpioracle;
192 SCIP_NLPIPROBLEM* nlpiprob;
193 SCIP_NLPI* nlpi;
194 SCIP_HASHMAP* var2nlpiidx;
195 SCIP_Real timelimit;
196 SCIP_Real objvarlb;
197 SCIP_Real minusone;
198 SCIP_Real one;
199 int nconvexnlrows;
200 int iterlimit;
201 int objvaridx;
202 int nconss;
203 int nvars;
204 int i;
205
206 assert(scip != NULL);
207 assert(sepadata != NULL);
208 assert(!sepadata->skipsepa);
209
210 SCIPdebugMsg(scip, "Computing interior point\n");
211
212 /* create convex relaxation NLP */
213 assert(SCIPgetNNlpis(scip) > 0);
214
215 nlpi = SCIPgetNlpis(scip)[0];
216 assert(nlpi != NULL);
217
218 nvars = SCIPgetNVars(scip);
219 SCIP_CALL( SCIPnlpiCreateProblem(nlpi, &nlpiprob, "gauge-interiorpoint-nlp") );
220 SCIP_CALL( SCIPhashmapCreate(&var2nlpiidx, SCIPblkmem(scip), nvars) );
221 SCIP_CALL( SCIPcreateNlpiProb(scip, nlpi, SCIPgetNLPNlRows(scip), SCIPgetNNLPNlRows(scip), nlpiprob, var2nlpiidx,
222 NULL, NULL, SCIPgetCutoffbound(scip), FALSE, TRUE) );
223
224 /* add objective variable; the problem is \min t, s.t. g(x) <= t, l(x) <= 0, where g are nonlinear and l linear */
225 objvaridx = nvars;
226 objvarlb = INTERIOROBJVARLB;
227 one = 1.0;
228 SCIP_CALL( SCIPnlpiAddVars(nlpi, nlpiprob, 1, &objvarlb, NULL, NULL) );
229 SCIP_CALL( SCIPnlpiSetObjective(nlpi, nlpiprob, 1, &objvaridx, &one, 0, NULL, NULL, NULL, 0.0) );
230
231 /* add objective variables to constraints; for this we need to get nlpi oracle to have access to number of
232 * constraints and which constraints are nonlinear
233 */
234 /* @todo: this code is only valid when using IPOPT and needs to be changed when new NLP solvers get interfaced */
235 assert(strcmp(SCIPnlpiGetName(nlpi), "ipopt") == 0);
236 nlpioracle = (SCIP_NLPIORACLE *)SCIPgetNlpiOracleIpopt(nlpiprob);
237 assert(nlpioracle != NULL);
238 assert(SCIPnlpiOracleGetNVars(nlpioracle) == objvaridx + 1);
239
240 minusone = -1.0;
241 nconvexnlrows = 0;
242 nconss = SCIPnlpiOracleGetNConstraints(nlpioracle);
243 for( i = 0; i < nconss; i++ )
244 {
245 if( SCIPnlpiOracleGetConstraintDegree(nlpioracle, i) > 1 )
246 {
247 SCIP_CALL( SCIPnlpiChgLinearCoefs(nlpi, nlpiprob, i, 1, &objvaridx, &minusone) );
248 ++nconvexnlrows;
249 }
250 }
251 SCIPdebug( SCIP_CALL( SCIPnlpiOraclePrintProblem(nlpioracle, SCIPgetMessagehdlr(scip), NULL) ) );
252
253 /* check if convex relaxation is interesting */
254 if( nconvexnlrows < 2 )
255 {
256 SCIPdebugMsg(scip, "convex relaxation is not interesting, only %d nonlinear convex rows; abort\n", nconvexnlrows);
257 sepadata->skipsepa = TRUE;
258 goto CLEANUP;
259 }
260
261 /* add linear rows */
262 SCIP_CALL( SCIPaddNlpiProbRows(scip, nlpi, nlpiprob, var2nlpiidx, SCIPgetLPRows(scip), SCIPgetNLPRows(scip)) );
263
264 /* set parameters in nlpi; time and iterations limit, tolerance, verbosity; for time limit, get time limit of scip;
265 * if scip doesn't have much time left, don't run separator. otherwise, timelimit is the minimum between whats left
266 * for scip and the timelimit setting
267 */
268 SCIP_CALL( SCIPgetRealParam(scip, "limits/time", &timelimit) );
269 if( !SCIPisInfinity(scip, timelimit) )
270 {
271 timelimit -= SCIPgetSolvingTime(scip);
272 if( timelimit <= 1.0 )
273 {
274 SCIPdebugMsg(scip, "skip NLP solve; no time left\n");
275 sepadata->skipsepa = TRUE;
276 goto CLEANUP;
277 }
278 }
279 if( sepadata->nlptimelimit > 0.0 )
280 timelimit = MIN(sepadata->nlptimelimit, timelimit);
281 SCIP_CALL( SCIPnlpiSetRealPar(nlpi, nlpiprob, SCIP_NLPPAR_TILIM, timelimit) );
282
283 iterlimit = sepadata->nlpiterlimit > 0 ? sepadata->nlpiterlimit : INT_MAX;
284 SCIP_CALL( SCIPnlpiSetIntPar(nlpi, nlpiprob, SCIP_NLPPAR_ITLIM, iterlimit) );
285 SCIP_CALL( SCIPnlpiSetRealPar(nlpi, nlpiprob, SCIP_NLPPAR_FEASTOL, NLPFEASFAC * SCIPfeastol(scip)) );
286 SCIP_CALL( SCIPnlpiSetRealPar(nlpi, nlpiprob, SCIP_NLPPAR_RELOBJTOL, MAX(SCIPfeastol(scip), SCIPdualfeastol(scip))) ); /*lint !e666*/
287 SCIP_CALL( SCIPnlpiSetIntPar(nlpi, nlpiprob, SCIP_NLPPAR_VERBLEVEL, NLPVERBOSITY) );
288
289 /* compute interior point */
290 SCIPdebugMsg(scip, "starting interior point computation\n");
291 SCIP_CALL( SCIPnlpiSolve(nlpi, nlpiprob) );
292 SCIPdebugMsg(scip, "finish interior point computation\n");
293
294 #ifdef SCIP_DEBUG
295 {
296 SCIP_NLPSTATISTICS* nlpstatistics;
297
298 /* get statistics */
299 SCIP_CALL( SCIPnlpStatisticsCreate(SCIPblkmem(scip), &nlpstatistics) );
300 SCIP_CALL( SCIPnlpiGetStatistics(nlpi, nlpiprob, nlpstatistics) );
301
302 SCIPdebugMsg(scip, "nlpi took iters %d, time %g searching for an find interior point: solstat %d\n",
303 SCIPnlpStatisticsGetNIterations(nlpstatistics), SCIPnlpStatisticsGetTotalTime(nlpstatistics),
304 SCIPnlpiGetSolstat(nlpi, nlpiprob));
305
306 SCIPnlpStatisticsFree(SCIPblkmem(scip), &nlpstatistics);
307 }
308 #endif
309
310 if( SCIPnlpiGetSolstat(nlpi, nlpiprob) <= SCIP_NLPSOLSTAT_FEASIBLE )
311 {
312 SCIP_Real* nlpisol;
313
314 SCIP_CALL( SCIPnlpiGetSolution(nlpi, nlpiprob, &nlpisol, NULL, NULL, NULL, NULL) );
315
316 assert(nlpisol != NULL);
317 SCIPdebugMsg(scip, "NLP solved: sol found has objvalue = %g\n", nlpisol[objvaridx]);
318
319 /* if we found an interior point store it */
320 if( SCIPisFeasNegative(scip, nlpisol[objvaridx]) )
321 {
322 SCIPdebugMsg(scip, "Interior point found!, storing it\n");
323 SCIP_CALL( SCIPcreateSol(scip, &sepadata->intsol, NULL) );
324 for( i = 0; i < nvars; i ++ )
325 {
326 SCIP_VAR* var;
327
328 var = SCIPgetVars(scip)[i];
329 assert(SCIPhashmapExists(var2nlpiidx, (void*)var) );
330
331 /* @todo: filter zero? */
332 SCIP_CALL( SCIPsetSolVal(scip, sepadata->intsol, var,
333 nlpisol[SCIPhashmapGetImageInt(var2nlpiidx, (void *)var)]) );
334 }
335
336 sepadata->isintsolavailable = TRUE;
337 }
338 else
339 {
340 SCIPdebugMsg(scip, "We got a feasible point but not interior (objval: %g)\n", nlpisol[objvaridx]);
341 sepadata->skipsepa = TRUE;
342 }
343 }
344 else
345 {
346 SCIPdebugMsg(scip, "We couldn't get an interior point (stat: %d)\n", SCIPnlpiGetSolstat(nlpi, nlpiprob));
347 sepadata->skipsepa = TRUE;
348 }
349
350 CLEANUP:
351 /* free memory */
352 SCIPhashmapFree(&var2nlpiidx);
353 SCIP_CALL( SCIPnlpiFreeProblem(nlpi, &nlpiprob) );
354
355 return SCIP_OKAY;
356 }
357
358
359 /** find whether point is in the interior, at the boundary or in the exterior of the region described by the
360 * intersection of nlrows[i] <= rhs if convexsides[i] = RHS or lhs <= nlrows[i] if convexsides[i] = LHS
361 * @note: point corresponds to a convex combination between the lp solution and the interior point
362 */
363 static
findPointPosition(SCIP * scip,SCIP_NLROW ** nlrows,int * nlrowsidx,int nnlrowsidx,CONVEXSIDE * convexsides,SCIP_SOL * point,POSITION * position)364 SCIP_RETCODE findPointPosition(
365 SCIP* scip, /**< SCIP data structure */
366 SCIP_NLROW** nlrows, /**< nlrows defining the region */
367 int* nlrowsidx, /**< indices of nlrows defining the region */
368 int nnlrowsidx, /**< number of nlrows indices */
369 CONVEXSIDE* convexsides, /**< sides of the nlrows involved in the region */
370 SCIP_SOL* point, /**< point for which we want to know its position */
371 POSITION* position /**< buffer to store position of sol */
372 )
373 {
374 int i;
375
376 assert(scip != NULL);
377 assert(nlrows != NULL);
378 assert(convexsides != NULL);
379 assert(nnlrowsidx > 0);
380 assert(point != NULL);
381 assert(position != NULL);
382
383 *position = INTERIOR;
384 for( i = 0; i < nnlrowsidx; i++ )
385 {
386 SCIP_NLROW* nlrow;
387 SCIP_Real activity;
388 CONVEXSIDE convexside;
389
390 nlrow = nlrows[nlrowsidx[i]];
391 convexside = convexsides[nlrowsidx[i]];
392
393 /* compute activity of nlrow at point */
394 SCIP_CALL( SCIPgetNlRowSolActivity(scip, nlrow, point, &activity) );
395
396 if( convexside == RHS )
397 {
398 assert(!SCIPisInfinity(scip, SCIPnlrowGetRhs(nlrow)));
399
400 /* if nlrow <= rhs is violated, then we are in the exterior */
401 if( SCIPisFeasGT(scip, activity, SCIPnlrowGetRhs(nlrow)) )
402 {
403 *position = EXTERIOR;
404 SCIPdebugMsg(scip, "exterior because cons <%s> has activity %g. rhs: %g\n", SCIPnlrowGetName(nlrow),
405 activity, SCIPnlrowGetRhs(nlrow));
406 SCIPdebug( SCIPprintNlRow(scip, nlrow, NULL) );
407
408 return SCIP_OKAY;
409 }
410
411 /* if nlrow(point) == rhs, then we are currently at the boundary */
412 if( SCIPisFeasEQ(scip, activity, SCIPnlrowGetRhs(nlrow)) )
413 *position = BOUNDARY;
414 }
415 else
416 {
417 assert(!SCIPisInfinity(scip, -SCIPnlrowGetLhs(nlrow)));
418 assert(convexside == LHS);
419
420 /* if lhs <= nlrow is violated, then we are in the exterior */
421 if( SCIPisFeasLT(scip, activity, SCIPnlrowGetLhs(nlrow)) )
422 {
423 *position = EXTERIOR;
424 return SCIP_OKAY;
425 }
426
427 /* if lhs == nlrow(point), then we are currently at the boundary */
428 if( SCIPisFeasEQ(scip, activity, SCIPnlrowGetLhs(nlrow)) )
429 *position = BOUNDARY;
430 }
431 }
432
433 return SCIP_OKAY;
434 }
435
436
437 /** returns, in convexcomb, the convex combination
438 * \f$ \lambda \f$ endpoint + (1 - \f$ lambda \f$) startpoint = startpoint + \f$ \lambda \f$ (endpoint - tosepasol) */
439 static
buildConvexCombination(SCIP * scip,SCIP_Real lambda,SCIP_SOL * startpoint,SCIP_SOL * endpoint,SCIP_SOL * convexcomb)440 SCIP_RETCODE buildConvexCombination(
441 SCIP* scip, /**< SCIP data structure */
442 SCIP_Real lambda, /**< convex combination multiplier */
443 SCIP_SOL* startpoint, /**< point corresponding to \f$ \lambda = 0 \f$ */
444 SCIP_SOL* endpoint, /**< point corresponding to \f$ \lambda = 1 \f$ */
445 SCIP_SOL* convexcomb /**< solution to store convex combination of intsol and tosepasol */
446 )
447 {
448 SCIP_VAR** vars;
449 int nvars;
450 int i;
451
452 assert(scip != NULL);
453 assert(startpoint != NULL);
454 assert(endpoint != NULL);
455 assert(convexcomb != NULL);
456
457 vars = SCIPgetVars(scip);
458 nvars = SCIPgetNVars(scip);
459
460 for( i = 0; i < nvars; i++ )
461 {
462 SCIP_Real val;
463 SCIP_VAR* var;
464
465 var = vars[i];
466 val = lambda * SCIPgetSolVal(scip, endpoint, var) + (1.0 - lambda) * SCIPgetSolVal(scip, startpoint, var);
467
468 if( !SCIPisZero(scip, val) )
469 {
470 SCIP_CALL( SCIPsetSolVal(scip, convexcomb, var, val) );
471 }
472 else
473 {
474 SCIP_CALL( SCIPsetSolVal(scip, convexcomb, var, 0.0) );
475 }
476 }
477
478 return SCIP_OKAY;
479 }
480
481
482 /** performs binary search to find the point belonging to the segment [intsol, tosepasol] that intersects the boundary
483 * of the region described by the intersection of nlrows[i] <= rhs if convexsides[i] = RHS or lhs <= nlrows[i] if not,
484 * for i in nlrowsidx
485 */
486 static
findBoundaryPoint(SCIP * scip,SCIP_NLROW ** nlrows,int * nlrowsidx,int nnlrowsidx,CONVEXSIDE * convexsides,SCIP_SOL * intsol,SCIP_SOL * tosepasol,SCIP_SOL * sol,POSITION * position)487 SCIP_RETCODE findBoundaryPoint(
488 SCIP* scip, /**< SCIP data structure */
489 SCIP_NLROW** nlrows, /**< nlrows defining the region */
490 int* nlrowsidx, /**< indices of nlrows defining the region */
491 int nnlrowsidx, /**< number of nlrows indices */
492 CONVEXSIDE* convexsides, /**< sides of the nlrows involved in the region */
493 SCIP_SOL* intsol, /**< point acting as 'interior point' */
494 SCIP_SOL* tosepasol, /**< solution that should be separated */
495 SCIP_SOL* sol, /**< convex combination of intsol and lpsol */
496 POSITION* position /**< buffer to store position of sol */
497 )
498 {
499 SCIP_Real lb;
500 SCIP_Real ub;
501 int i;
502
503 assert(scip != NULL);
504 assert(nlrows != NULL);
505 assert(nlrowsidx != NULL);
506 assert(convexsides != NULL);
507 assert(intsol != NULL);
508 assert(tosepasol != NULL);
509 assert(sol != NULL);
510 assert(position != NULL);
511
512 SCIPdebugMsg(scip, "starting binary search\n");
513 lb = 0.0; /* corresponds to intsol */
514 ub = 1.0; /* corresponds to tosepasol */
515 for( i = 0; i < MAX_ITER; i++ )
516 {
517 /* sol = (ub+lb)/2 * lpsol + (1 - (ub+lb)/2) * intsol */
518 SCIP_CALL( buildConvexCombination(scip, (ub + lb)/2.0, intsol, tosepasol, sol) );
519
520 /* find poisition of point: boundary, interior, exterior */
521 SCIP_CALL( findPointPosition(scip, nlrows, nlrowsidx, nnlrowsidx, convexsides, sol, position) );
522 SCIPdebugMsg(scip, "Position: %d, lambda: %g\n", *position, (ub + lb)/2.0);
523
524 switch( *position )
525 {
526 case BOUNDARY:
527 SCIPdebugMsg(scip, "Done\n");
528 return SCIP_OKAY;
529
530 case INTERIOR:
531 /* want to be closer to tosepasol */
532 lb = (ub + lb)/2.0;
533 break;
534
535 case EXTERIOR:
536 /* want to be closer to intsol */
537 ub = (ub + lb)/2.0;
538 break;
539 }
540 }
541 SCIPdebugMsg(scip, "Done\n");
542 return SCIP_OKAY;
543 }
544
545
546 /** computes gradient of exprtree at sol */
547 static
computeGradient(SCIP * scip,SCIP_EXPRINT * exprint,SCIP_SOL * sol,SCIP_EXPRTREE * exprtree,SCIP_Real * grad)548 SCIP_RETCODE computeGradient(
549 SCIP* scip, /**< SCIP data structure */
550 SCIP_EXPRINT* exprint, /**< expressions interpreter */
551 SCIP_SOL* sol, /**< point where we compute gradient */
552 SCIP_EXPRTREE* exprtree, /**< exprtree for which we compute the gradient */
553 SCIP_Real* grad /**< buffer to store the gradient */
554 )
555 {
556 SCIP_Real* x;
557 SCIP_Real val;
558 int nvars;
559 int i;
560
561 assert(scip != NULL);
562 assert(exprint != NULL);
563 assert(sol != NULL);
564 assert(exprtree != NULL);
565 assert(grad != NULL);
566
567 nvars = SCIPexprtreeGetNVars(exprtree);
568 assert(nvars > 0);
569
570 SCIP_CALL( SCIPallocBufferArray(scip, &x, nvars) );
571
572 /* compile expression exprtree, if not done before */
573 if( SCIPexprtreeGetInterpreterData(exprtree) == NULL )
574 {
575 SCIP_CALL( SCIPexprintCompile(exprint, exprtree) );
576 }
577
578 for( i = 0; i < nvars; ++i )
579 {
580 x[i] = SCIPgetSolVal(scip, sol, SCIPexprtreeGetVars(exprtree)[i]);
581 }
582
583 SCIP_CALL( SCIPexprintGrad(exprint, exprtree, x, TRUE, &val, grad) );
584
585 /*SCIPdebug( for( i = 0; i < nvars; ++i ) printf("%e [%s]\n", grad[i], SCIPvarGetName(SCIPexprtreeGetVars(exprtree)[i])) );*/
586
587 SCIPfreeBufferArray(scip, &x);
588
589 return SCIP_OKAY;
590 }
591
592
593 /** computes gradient cut (linearization) of nlrow at sol */
594 static
generateCut(SCIP * scip,SCIP_SOL * sol,SCIP_EXPRINT * exprint,SCIP_NLROW * nlrow,CONVEXSIDE convexside,SCIP_ROW * row,SCIP_Bool * success)595 SCIP_RETCODE generateCut(
596 SCIP* scip, /**< SCIP data structure */
597 SCIP_SOL* sol, /**< point used to construct gradient cut (x_0) */
598 SCIP_EXPRINT* exprint, /**< expression interpreter */
599 SCIP_NLROW* nlrow, /**< constraint */
600 CONVEXSIDE convexside, /**< whether we use rhs or lhs of nlrow */
601 SCIP_ROW* row, /**< storage for cut */
602 SCIP_Bool* success /**< buffer to store whether the gradient was finite */
603 )
604 {
605 SCIP_Real activity;
606 SCIP_Real gradx0; /* <grad f(x_0), x_0> */
607 int i;
608
609 assert(scip != NULL);
610 assert(exprint != NULL);
611 assert(nlrow != NULL);
612 assert(row != NULL);
613
614 gradx0 = 0.0;
615 *success = TRUE;
616
617 /* an nlrow has a linear part, quadratic part and expression tree; ideally one would just build the gradient but we
618 * do not know if the different parts share variables or not, so we can't just build the gradient; for this reason
619 * we create the row right away and compute the gradients of each part independently and add them to the row; the
620 * row takes care to add coeffs corresponding to the same variable when they appear in different parts of the nlrow
621 */
622
623 SCIP_CALL( SCIPcacheRowExtensions(scip, row) );
624
625 /* linear part */
626 for( i = 0; i < SCIPnlrowGetNLinearVars(nlrow); i++ )
627 {
628 gradx0 += SCIPgetSolVal(scip, sol, SCIPnlrowGetLinearVars(nlrow)[i]) * SCIPnlrowGetLinearCoefs(nlrow)[i];
629 SCIP_CALL( SCIPaddVarToRow(scip, row, SCIPnlrowGetLinearVars(nlrow)[i], SCIPnlrowGetLinearCoefs(nlrow)[i]) );
630 }
631
632 /* quadratic part */
633 for( i = 0; i < SCIPnlrowGetNQuadElems(nlrow); i++ )
634 {
635 SCIP_VAR* var1;
636 SCIP_VAR* var2;
637 SCIP_Real grad1;
638 SCIP_Real grad2;
639 SCIP_Real solval1;
640 SCIP_Real solval2;
641
642 var1 = SCIPnlrowGetQuadVars(nlrow)[SCIPnlrowGetQuadElems(nlrow)[i].idx1];
643 var2 = SCIPnlrowGetQuadVars(nlrow)[SCIPnlrowGetQuadElems(nlrow)[i].idx2];
644 solval1 = SCIPgetSolVal(scip, sol, var1);
645 solval2 = SCIPgetSolVal(scip, sol, var2);
646 grad1 = SCIPnlrowGetQuadElems(nlrow)[i].coef * solval2; /* note that solval2 is correct */
647 grad2 = SCIPnlrowGetQuadElems(nlrow)[i].coef * solval1;
648
649 SCIP_CALL( SCIPaddVarToRow(scip, row, var1, grad1) );
650 SCIP_CALL( SCIPaddVarToRow(scip, row, var2, grad2) );
651
652 gradx0 += grad1 * solval1 + grad2 * solval2;
653 }
654
655 /* expression tree part */
656 {
657 SCIP_Real* grad;
658 SCIP_EXPRTREE* tree;
659
660 tree = SCIPnlrowGetExprtree(nlrow);
661
662 if( tree != NULL && SCIPexprtreeGetNVars(tree) > 0 )
663 {
664 SCIP_CALL( SCIPallocBufferArray(scip, &grad, SCIPexprtreeGetNVars(tree)) );
665
666 SCIP_CALL( computeGradient(scip, exprint, sol, tree, grad) );
667
668 for( i = 0; i < SCIPexprtreeGetNVars(tree); i++ )
669 {
670 /* check gradient entries: function might not be differentiable */
671 if( !SCIPisFinite(grad[i]) || SCIPisInfinity(scip, grad[i]) || SCIPisInfinity(scip, -grad[i]) )
672 {
673 *success = FALSE;
674 break;
675 }
676
677 gradx0 += grad[i] * SCIPgetSolVal(scip, sol, SCIPexprtreeGetVars(tree)[i]);
678 SCIP_CALL( SCIPaddVarToRow(scip, row, SCIPexprtreeGetVars(tree)[i], grad[i]) );
679 }
680
681 SCIPfreeBufferArray(scip, &grad);
682 }
683 }
684
685 SCIP_CALL( SCIPflushRowExtensions(scip, row) );
686
687 /* if there was a problem computing the cut -> return */
688 if( ! *success )
689 return SCIP_OKAY;
690
691 #ifdef CUT_DEBUG
692 SCIPdebugMsg(scip, "gradient: ");
693 SCIPdebug( SCIP_CALL( SCIPprintRow(scip, row, NULL) ) );
694 SCIPdebugMsg(scip, "gradient dot x_0: %g\n", gradx0);
695 #endif
696
697 /* gradient cut is f(x_0) - <grad f(x_0), x_0> + <grad f(x_0), x> <= rhs or >= lhs */
698 SCIP_CALL( SCIPgetNlRowSolActivity(scip, nlrow, sol, &activity) );
699 if( convexside == RHS )
700 {
701 assert(!SCIPisInfinity(scip, SCIPnlrowGetRhs(nlrow)));
702 SCIP_CALL( SCIPchgRowRhs(scip, row, SCIPnlrowGetRhs(nlrow) - activity + gradx0) );
703 }
704 else
705 {
706 assert(convexside == LHS);
707 assert(!SCIPisInfinity(scip, -SCIPnlrowGetLhs(nlrow)));
708 SCIP_CALL( SCIPchgRowLhs(scip, row, SCIPnlrowGetLhs(nlrow) - activity + gradx0) );
709 }
710
711 #ifdef CUT_DEBUG
712 SCIPdebugMsg(scip, "gradient cut: ");
713 SCIPdebug( SCIP_CALL( SCIPprintRow(scip, row, NULL) ) );
714 #endif
715
716 return SCIP_OKAY;
717 }
718
719 /** tries to generate gradient cuts at the point on the segment [intsol, tosepasol] that intersecs the boundary of the
720 * convex relaxation
721 * -# checks that the relative interior of the segment actually intersects the boundary
722 * (this check is needed since intsol is not necessarily an interior point)
723 * -# finds point on the boundary
724 * -# generates gradient cut at point on the boundary
725 */
726 static
separateCuts(SCIP * scip,SCIP_SEPA * sepa,SCIP_SOL * tosepasol,SCIP_RESULT * result)727 SCIP_RETCODE separateCuts(
728 SCIP* scip, /**< SCIP data structure */
729 SCIP_SEPA* sepa, /**< the cut separator itself */
730 SCIP_SOL* tosepasol, /**< solution that should be separated */
731 SCIP_RESULT* result /**< pointer to store the result of the separation call */
732 )
733 {
734 SCIP_SEPADATA* sepadata;
735 SCIP_NLROW** nlrows;
736 CONVEXSIDE* convexsides;
737 SCIP_SOL* sol;
738 SCIP_SOL* intsol;
739 POSITION position;
740 int* nlrowsidx;
741 int nnlrowsidx;
742 int i;
743
744 assert(sepa != NULL);
745
746 sepadata = SCIPsepaGetData(sepa);
747 assert(sepadata != NULL);
748
749 intsol = sepadata->intsol;
750 nlrows = sepadata->nlrows;
751 nlrowsidx = sepadata->nlrowsidx;
752 nnlrowsidx = sepadata->nnlrowsidx;
753 convexsides = sepadata->convexsides;
754
755 assert(intsol != NULL);
756 assert(nlrows != NULL);
757 assert(nlrowsidx != NULL);
758 assert(nnlrowsidx > 0);
759 assert(convexsides != NULL);
760
761 /* to evaluate the nlrow one needs a solution */
762 SCIP_CALL( SCIPcreateSol(scip, &sol, NULL) );
763
764 /* don't separate if, under SCIP tolerances, only a slight perturbation of the interior point in the direction of
765 * tosepasol gives a point that is in the exterior */
766 SCIP_CALL( buildConvexCombination(scip, VIOLATIONFAC * SCIPfeastol(scip), intsol, tosepasol, sol) );
767 SCIP_CALL( findPointPosition(scip, nlrows, nlrowsidx, nnlrowsidx, convexsides, sol, &position) );
768
769 if( position == EXTERIOR )
770 {
771 #ifdef SCIP_DEBUG
772 SCIPdebugMsg(scip, "segment joining intsol and tosepasol seems to be contained in the exterior of the region, can't separate\n");
773 /* move from intsol in the direction of -tosepasol to check if we are really tangent to the region */
774 SCIP_CALL( buildConvexCombination(scip, -1e-3, intsol, tosepasol, sol) );
775 SCIP_CALL( findPointPosition(scip, nlrows, nlrowsidx, nnlrowsidx, convexsides, sol, &position) );
776 if( position == EXTERIOR )
777 {
778 SCIPdebugMsg(scip, "line through intsol and tosepasol is tangent to region; can't separate\n");
779 }
780 SCIP_CALL( findPointPosition(scip, nlrows, nlrowsidx, nnlrowsidx, convexsides, intsol, &position) );
781 printf("Position of intsol is %s\n",
782 position == EXTERIOR ? "exterior" : position == INTERIOR ? "interior": "boundary");
783 SCIP_CALL( findPointPosition(scip, nlrows, nlrowsidx, nnlrowsidx, convexsides, tosepasol, &position) );
784 printf("Position of tosepasol is %s\n",
785 position == EXTERIOR ? "exterior" : position == INTERIOR ? "interior": "boundary");
786
787 /* slightly move from intsol in the direction of +-tosepasol */
788 SCIP_CALL( buildConvexCombination(scip, 1e-5, intsol, tosepasol, sol) );
789 SCIP_CALL( findPointPosition(scip, nlrows, nlrowsidx, nnlrowsidx, convexsides, sol, &position) );
790 printf("Position of intsol + 0.00001(tosepasol - inisol) is %s\n",
791 position == EXTERIOR ? "exterior" : position == INTERIOR ? "interior": "boundary");
792 SCIPdebug( SCIPprintSol(scip, sol, NULL, FALSE) );
793
794 SCIP_CALL( buildConvexCombination(scip, -1e-5, intsol, tosepasol, sol) );
795 SCIP_CALL( findPointPosition(scip, nlrows, nlrowsidx, nnlrowsidx, convexsides, sol, &position) );
796 printf("Position of intsol - 0.00001(tosepasol - inisol) is %s\n",
797 position == EXTERIOR ? "exterior" : position == INTERIOR ? "interior": "boundary");
798 SCIPdebug( SCIPprintSol(scip, sol, NULL, FALSE) );
799 #endif
800 *result = SCIP_DIDNOTFIND;
801 goto CLEANUP;
802 }
803
804 /* find point on boundary */
805 if( position != BOUNDARY )
806 {
807 SCIP_CALL( findBoundaryPoint(scip, nlrows, nlrowsidx, nnlrowsidx, convexsides, intsol, tosepasol, sol,
808 &position) );
809
810 /* if MAX_ITER weren't enough to find a point in the boundary we don't separate */
811 if( position != BOUNDARY )
812 {
813 SCIPdebugMsg(scip, "couldn't find boundary point, don't separate\n");
814 goto CLEANUP;
815 }
816 }
817
818 /* generate cuts at sol */
819 for( i = 0; i < nnlrowsidx; i++ )
820 {
821 SCIP_NLROW* nlrow;
822 SCIP_ROW* row;
823 SCIP_Real activity;
824 CONVEXSIDE convexside;
825 SCIP_Bool success;
826 char rowname[SCIP_MAXSTRLEN];
827
828 nlrow = nlrows[nlrowsidx[i]];
829 convexside = convexsides[nlrowsidx[i]];
830
831 (void) SCIPsnprintf(rowname, SCIP_MAXSTRLEN, "%s_%u", SCIPnlrowGetName(nlrow), ++(sepadata->ncuts));
832
833 /* only separate nlrows that are tight at the boundary point */
834 SCIP_CALL( SCIPgetNlRowSolActivity(scip, nlrow, sol, &activity) );
835 SCIPdebugMsg(scip, "cons <%s> at boundary point has activity: %g\n", SCIPnlrowGetName(nlrow), activity);
836
837 if( (convexside == RHS && !SCIPisFeasEQ(scip, activity, SCIPnlrowGetRhs(nlrow)))
838 || (convexside == LHS && !SCIPisFeasEQ(scip, activity, SCIPnlrowGetLhs(nlrow))) )
839 continue;
840
841 /* cut is globally valid, since we work on nlrows from the NLP built at the root node, which are globally valid */
842 /* @todo: when local nlrows get supported in SCIP, one can think of recomputing the interior point */
843 SCIP_CALL( SCIPcreateEmptyRowSepa(scip, &row, sepa, rowname, -SCIPinfinity(scip), SCIPinfinity(scip),
844 FALSE, FALSE , TRUE) );
845 SCIP_CALL( generateCut(scip, sol, sepadata->exprinterpreter, nlrow, convexside, row, &success) );
846
847 /* add cut */
848 SCIPdebugMsg(scip, "cut <%s> has efficacy %g\n", SCIProwGetName(row), SCIPgetCutEfficacy(scip, NULL, row));
849 if( success && SCIPisCutEfficacious(scip, NULL, row) )
850 {
851 SCIP_Bool infeasible;
852
853 SCIPdebugMsg(scip, "adding cut\n");
854 SCIP_CALL( SCIPaddRow(scip, row, FALSE, &infeasible) );
855
856 if( infeasible )
857 {
858 *result = SCIP_CUTOFF;
859 SCIP_CALL( SCIPreleaseRow(scip, &row) );
860 break;
861 }
862 else
863 {
864 *result = SCIP_SEPARATED;
865 }
866 }
867
868 /* release the row */
869 SCIP_CALL( SCIPreleaseRow(scip, &row) );
870 }
871
872 CLEANUP:
873 SCIP_CALL( SCIPfreeSol(scip, &sol) );
874
875 return SCIP_OKAY;
876 }
877
878 /*
879 * Callback methods of separator
880 */
881
882
883 /** destructor of separator to free user data (called when SCIP is exiting) */
884 static
SCIP_DECL_SEPAFREE(sepaFreeGauge)885 SCIP_DECL_SEPAFREE(sepaFreeGauge)
886 { /*lint --e{715}*/
887 SCIP_SEPADATA* sepadata;
888
889 assert(strcmp(SCIPsepaGetName(sepa), SEPA_NAME) == 0);
890
891 /* free separator data */
892 sepadata = SCIPsepaGetData(sepa);
893 assert(sepadata != NULL);
894
895 SCIPfreeBlockMemory(scip, &sepadata);
896
897 SCIPsepaSetData(sepa, NULL);
898
899 return SCIP_OKAY;
900 }
901
902
903 /** solving process deinitialization method of separator (called before branch and bound process data is freed) */
904 static
SCIP_DECL_SEPAEXITSOL(sepaExitsolGauge)905 SCIP_DECL_SEPAEXITSOL(sepaExitsolGauge)
906 { /*lint --e{715}*/
907 SCIP_SEPADATA* sepadata;
908
909 assert(sepa != NULL);
910
911 sepadata = SCIPsepaGetData(sepa);
912
913 assert(sepadata != NULL);
914
915 /* free memory and reset data */
916 if( sepadata->isintsolavailable )
917 {
918 SCIPfreeBlockMemoryArray(scip, &sepadata->nlrowsidx, sepadata->nlrowssize);
919 SCIPfreeBlockMemoryArray(scip, &sepadata->convexsides, sepadata->nlrowssize);
920 SCIPfreeBlockMemoryArray(scip, &sepadata->nlrows, sepadata->nlrowssize);
921 SCIP_CALL( SCIPfreeSol(scip, &sepadata->intsol) );
922 SCIP_CALL( SCIPexprintFree(&sepadata->exprinterpreter) );
923
924 sepadata->nnlrows = 0;
925 sepadata->nnlrowsidx = 0;
926 sepadata->nlrowssize = 0;
927 sepadata->isintsolavailable = FALSE;
928 }
929 assert(sepadata->nnlrows == 0);
930 assert(sepadata->nnlrowsidx == 0);
931 assert(sepadata->nlrowssize == 0);
932 assert(sepadata->isintsolavailable == FALSE);
933
934 sepadata->skipsepa = FALSE;
935
936 return SCIP_OKAY;
937 }
938
939
940 /** LP solution separation method of separator */
941 static
SCIP_DECL_SEPAEXECLP(sepaExeclpGauge)942 SCIP_DECL_SEPAEXECLP(sepaExeclpGauge)
943 { /*lint --e{715}*/
944 SCIP_SEPADATA* sepadata;
945 SCIP_SOL* lpsol;
946 int i;
947
948 assert(scip != NULL);
949 assert(sepa != NULL);
950
951 sepadata = SCIPsepaGetData(sepa);
952
953 assert(sepadata != NULL);
954
955 *result = SCIP_DIDNOTRUN;
956
957 /* do not run if there is no interesting convex relaxation (with at least two nonlinear convex constraint) */
958 if( sepadata->skipsepa )
959 {
960 SCIPdebugMsg(scip, "not running because convex relaxation is uninteresting\n");
961 return SCIP_OKAY;
962 }
963
964 /* do not run if SCIP has not constructed an NLP */
965 if( !SCIPisNLPConstructed(scip) )
966 {
967 SCIPdebugMsg(scip, "NLP not constructed, skipping gauge separator\n");
968 return SCIP_OKAY;
969 }
970
971 /* do not run if SCIP has no way of solving nonlinear problems */
972 if( SCIPgetNNlpis(scip) == 0 )
973 {
974 SCIPdebugMsg(scip, "Skip gauge separator: no nlpi and SCIP can't solve nonlinear problems without a nlpi\n");
975 return SCIP_OKAY;
976 }
977
978 /* if we don't have an interior point compute one; if we fail to compute one, then separator will not be run again;
979 * otherwise, we also store the convex nlrows in sepadata
980 */
981 if( !sepadata->isintsolavailable )
982 {
983 /* @todo: one could store the convex nonlinear rows inside computeInteriorPoint */
984 SCIP_CALL( computeInteriorPoint(scip, sepadata) );
985 assert(sepadata->skipsepa || sepadata->isintsolavailable);
986
987 if( sepadata->skipsepa )
988 return SCIP_OKAY;
989
990 SCIP_CALL( storeNonlinearConvexNlrows(scip, sepadata, SCIPgetNLPNlRows(scip), SCIPgetNNLPNlRows(scip)) );
991 SCIP_CALL( SCIPexprintCreate(SCIPblkmem(scip), &sepadata->exprinterpreter) );
992 }
993
994 #ifdef SCIP_DISABLED_CODE
995 /* get interior point: try to compute an interior point, otherwise use primal solution, otherwise use NLP solution */
996 /* @todo: - decide order:
997 * - we can also use convex combination of solutions; there is a function SCIPvarGetAvgSol!
998 * - can add an event handler to only update when a new solution has been found
999 */
1000 if( !sepadata->isintsolavailable )
1001 {
1002 if( SCIPgetNSols(scip) > 0 )
1003 {
1004 SCIPdebugMsg(scip, "Using current primal solution as interior point!\n");
1005 SCIP_CALL( SCIPcreateSolCopy(scip, &sepadata->intsol, SCIPgetBestSol(scip)) );
1006 sepadata->isintsolavailable = TRUE;
1007 }
1008 else if( SCIPhasNLPSolution(scip) )
1009 {
1010 SCIPdebugMsg(scip, "Using NLP solution as interior point!\n");
1011 SCIP_CALL( SCIPcreateNLPSol(scip, &sepadata->intsol, NULL) );
1012 sepadata->isintsolavailable = TRUE;
1013 }
1014 else
1015 {
1016 SCIPdebugMsg(scip, "We couldn't find an interior point, don't have a feasible nor an NLP solution; skip separator\n");
1017 return SCIP_OKAY;
1018 }
1019 }
1020 #endif
1021
1022 /* store lp sol (or pseudo sol when lp is not solved) to be able to use it to compute nlrows' activities */
1023 SCIP_CALL( SCIPcreateCurrentSol(scip, &lpsol, NULL) );
1024
1025 /* store indices of relevant constraints, ie, the ones that violate the lp sol */
1026 sepadata->nnlrowsidx = 0;
1027 for( i = 0; i < sepadata->nnlrows; i++ )
1028 {
1029 SCIP_NLROW* nlrow;
1030 SCIP_Real activity;
1031
1032 nlrow = sepadata->nlrows[i];
1033
1034 SCIP_CALL( SCIPgetNlRowSolActivity(scip, nlrow, lpsol, &activity) );
1035
1036 if( sepadata->convexsides[i] == RHS )
1037 {
1038 assert(!SCIPisInfinity(scip, SCIPnlrowGetRhs(nlrow)));
1039
1040 if( activity - SCIPnlrowGetRhs(nlrow) < VIOLATIONFAC * SCIPfeastol(scip) )
1041 continue;
1042 }
1043 else
1044 {
1045 assert(sepadata->convexsides[i] == LHS);
1046 assert(!SCIPisInfinity(scip, -SCIPnlrowGetLhs(nlrow)));
1047
1048 if( SCIPnlrowGetLhs(nlrow) - activity < VIOLATIONFAC * SCIPfeastol(scip) )
1049 continue;
1050 }
1051
1052 sepadata->nlrowsidx[sepadata->nnlrowsidx] = i;
1053 ++(sepadata->nnlrowsidx);
1054 }
1055
1056 /* separate only if there are violated nlrows */
1057 SCIPdebugMsg(scip, "there are %d violated nlrows\n", sepadata->nnlrowsidx);
1058 if( sepadata->nnlrowsidx > 0 )
1059 {
1060 SCIP_CALL( separateCuts(scip, sepa, lpsol, result) );
1061 }
1062
1063 /* free lpsol */
1064 SCIP_CALL( SCIPfreeSol(scip, &lpsol) );
1065
1066 return SCIP_OKAY;
1067 }
1068
1069
1070 /*
1071 * separator specific interface methods
1072 */
1073
1074 /** creates the gauge separator and includes it in SCIP */
SCIPincludeSepaGauge(SCIP * scip)1075 SCIP_RETCODE SCIPincludeSepaGauge(
1076 SCIP* scip /**< SCIP data structure */
1077 )
1078 {
1079 SCIP_SEPADATA* sepadata;
1080 SCIP_SEPA* sepa;
1081
1082 /* create gauge separator data */
1083 SCIP_CALL( SCIPallocBlockMemory(scip, &sepadata) );
1084
1085 /* this sets all data in sepadata to 0 */
1086 BMSclearMemory(sepadata);
1087
1088 /* include separator */
1089 SCIP_CALL( SCIPincludeSepaBasic(scip, &sepa, SEPA_NAME, SEPA_DESC, SEPA_PRIORITY, SEPA_FREQ, SEPA_MAXBOUNDDIST,
1090 SEPA_USESSUBSCIP, SEPA_DELAY,
1091 sepaExeclpGauge, NULL,
1092 sepadata) );
1093
1094 assert(sepa != NULL);
1095
1096 /* set non fundamental callbacks via setter functions */
1097 SCIP_CALL( SCIPsetSepaFree(scip, sepa, sepaFreeGauge) );
1098 SCIP_CALL( SCIPsetSepaExitsol(scip, sepa, sepaExitsolGauge) );
1099
1100 /* add gauge separator parameters */
1101 SCIP_CALL( SCIPaddIntParam(scip, "separating/" SEPA_NAME "/nlpiterlimit",
1102 "iteration limit of NLP solver; 0 for no limit",
1103 &sepadata->nlpiterlimit, TRUE, DEFAULT_NLPITERLIM, 0, INT_MAX, NULL, NULL) );
1104
1105 SCIP_CALL( SCIPaddRealParam(scip, "separating/" SEPA_NAME "/nlptimelimit",
1106 "time limit of NLP solver; 0.0 for no limit",
1107 &sepadata->nlptimelimit, TRUE, DEFAULT_NLPTIMELIMIT, 0.0, SCIP_REAL_MAX, NULL, NULL) );
1108
1109 return SCIP_OKAY;
1110 }
1111