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 prop_nlobbt.c
17 * @ingroup DEFPLUGINS_PROP
18 * @brief nlobbt propagator
19 * @author Benjamin Mueller
20 */
21
22 /*---+----1----+----2----+----3----+----4----+----5----+----6----+----7----+----8----+----9----+----0----+----1----+----2*/
23
24 #include "blockmemshell/memory.h"
25 #include "nlpi/nlpi.h"
26 #include "scip/prop_genvbounds.h"
27 #include "scip/prop_nlobbt.h"
28 #include "scip/pub_message.h"
29 #include "scip/pub_misc.h"
30 #include "scip/pub_misc_sort.h"
31 #include "scip/pub_nlp.h"
32 #include "scip/pub_prop.h"
33 #include "scip/pub_tree.h"
34 #include "scip/pub_var.h"
35 #include "scip/scip_general.h"
36 #include "scip/scip_lp.h"
37 #include "scip/scip_mem.h"
38 #include "scip/scip_message.h"
39 #include "scip/scip_nlp.h"
40 #include "scip/scip_nonlinear.h"
41 #include "scip/scip_numerics.h"
42 #include "scip/scip_param.h"
43 #include "scip/scip_prob.h"
44 #include "scip/scip_probing.h"
45 #include "scip/scip_prop.h"
46 #include "scip/scip_randnumgen.h"
47 #include "scip/scip_solvingstats.h"
48 #include "scip/scip_timing.h"
49 #include "scip/scip_tree.h"
50 #include "scip/scip_var.h"
51 #include <string.h>
52
53 #define PROP_NAME "nlobbt"
54 #define PROP_DESC "propagator template"
55 #define PROP_PRIORITY -1100000
56 #define PROP_FREQ -1
57 #define PROP_DELAY TRUE
58 #define PROP_TIMING SCIP_PROPTIMING_AFTERLPLOOP
59
60 #define DEFAULT_MINNONCONVEXFRAC 0.20 /**< default minimum (# convex nlrows)/(# nonconvex nlrows) threshold to apply propagator */
61 #define DEFAULT_MINLINEARFRAC 0.02 /**< default minimum (# convex nlrows)/(# linear nlrows) threshold to apply propagator */
62 #define DEFAULT_FEASTOLFAC 0.01 /**< default factor for NLP feasibility tolerance */
63 #define DEFAULT_RELOBJTOLFAC 0.01 /**< default factor for NLP relative objective tolerance */
64 #define DEFAULT_ADDLPROWS TRUE /**< should (non-initial) LP rows be used? */
65 #define DEFAULT_ITLIMITFACTOR 2.0 /**< multiple of root node LP iterations used as total LP iteration
66 * limit for nlobbt (<= 0: no limit ) */
67 #define DEFAULT_NLPITERLIMIT 500 /**< default iteration limit of NLP solver; 0 for no limit */
68 #define DEFAULT_NLPTIMELIMIT 0.0 /**< default time limit of NLP solver; 0.0 for no limit */
69 #define DEFAULT_NLPVERLEVEL 0 /**< verbosity level of NLP solver */
70 #define DEFAULT_RANDSEED 79 /**< initial random seed */
71
72 /*
73 * Data structures
74 */
75
76 /* status of bound candidates */
77 #define UNSOLVED 1 /**< did not solve LB or UB problem */
78 #define SOLVEDLB 2 /**< solved LB problem */
79 #define SOLVEDUB 4 /**< solved UB problem */
80
81 /** propagator data */
82 struct SCIP_PropData
83 {
84 SCIP_NLPI* nlpi; /**< nlpi used to create the nlpi problem */
85 SCIP_NLPIPROBLEM* nlpiprob; /**< nlpi problem representing the convex NLP relaxation */
86 SCIP_HASHMAP* var2nlpiidx; /**< mapping between variables and nlpi indices */
87 SCIP_VAR** nlpivars; /**< array containing all variables of the nlpi */
88 int nlpinvars; /**< total number of nlpi variables */
89 SCIP_Real* nlscore; /**< score for each nonlinear variable */
90 int* status; /**< array containing a bound status for each candidate */
91 SCIP_PROP* genvboundprop; /**< genvbound propagator */
92 SCIP_RANDNUMGEN* randnumgen; /**< random number generator */
93 SCIP_Bool skipprop; /**< should the propagator be skipped? */
94 SCIP_Longint lastnode; /**< number of last node where obbt was performed */
95 int currpos; /**< current position in the nlpivars array */
96
97 int nlpiterlimit; /**< iteration limit of NLP solver; 0 for no limit */
98 SCIP_Real nlptimelimit; /**< time limit of NLP solver; 0.0 for no limit */
99 int nlpverblevel; /**< verbosity level of NLP solver */
100 SCIP_NLPSTATISTICS* nlpstatistics; /**< statistics from NLP solver */
101
102 SCIP_Real feastolfac; /**< factor for NLP feasibility tolerance */
103 SCIP_Real relobjtolfac; /**< factor for NLP relative objective tolerance */
104 SCIP_Real minnonconvexfrac; /**< minimum (#convex nlrows)/(#nonconvex nlrows) threshold to apply propagator */
105 SCIP_Real minlinearfrac; /**< minimum (#convex nlrows)/(#linear nlrows) threshold to apply propagator */
106 SCIP_Bool addlprows; /**< should (non-initial) LP rows be used? */
107 SCIP_Real itlimitfactor; /**< LP iteration limit for nlobbt will be this factor times total LP
108 * iterations in root node */
109 };
110
111 /*
112 * Local methods
113 */
114
115 /** clears the propagator data */
116 static
propdataClear(SCIP * scip,SCIP_PROPDATA * propdata)117 SCIP_RETCODE propdataClear(
118 SCIP* scip, /**< SCIP data structure */
119 SCIP_PROPDATA* propdata /**< propagator data */
120 )
121 {
122 assert(propdata != NULL);
123
124 if( propdata->nlpiprob != NULL )
125 {
126 assert(propdata->nlpi != NULL);
127
128 SCIPfreeBlockMemoryArray(scip, &propdata->status, propdata->nlpinvars);
129 SCIPfreeBlockMemoryArray(scip, &propdata->nlscore, propdata->nlpinvars);
130 SCIPfreeBlockMemoryArray(scip, &propdata->nlpivars, propdata->nlpinvars);
131 SCIPhashmapFree(&propdata->var2nlpiidx);
132 SCIP_CALL( SCIPnlpiFreeProblem(propdata->nlpi, &propdata->nlpiprob) );
133
134 propdata->nlpinvars = 0;
135 }
136 assert(propdata->nlpinvars == 0);
137
138 propdata->skipprop = FALSE;
139 propdata->currpos = 0;
140 propdata->lastnode = -1;
141
142 return SCIP_OKAY;
143 }
144
145 /** checks whether it is worth to call nonlinear OBBT procedure */
146 static
isNlobbtApplicable(SCIP * scip,SCIP_PROPDATA * propdata)147 SCIP_Bool isNlobbtApplicable(
148 SCIP* scip, /**< SCIP data structure */
149 SCIP_PROPDATA* propdata /**< propagation data */
150 )
151 {
152 SCIP_NLROW** nlrows;
153 int nnonconvexnlrows;
154 int nconvexnlrows;
155 int nlinearnlrows;
156 int nnlrows;
157 int i;
158
159 nlrows = SCIPgetNLPNlRows(scip);
160 nnlrows = SCIPgetNNLPNlRows(scip);
161 nnonconvexnlrows = 0;
162 nconvexnlrows = 0;
163 nlinearnlrows = 0;
164
165 for( i = 0; i < nnlrows; ++i )
166 {
167 if( SCIPnlrowGetNQuadElems(nlrows[i]) == 0 && SCIPnlrowGetExprtree(nlrows[i]) == NULL )
168 ++nlinearnlrows;
169 else if( SCIPnlrowGetCurvature(nlrows[i]) == SCIP_EXPRCURV_CONVEX )
170 {
171 if( !SCIPisInfinity(scip, SCIPnlrowGetRhs(nlrows[i])) )
172 ++nconvexnlrows;
173 if( !SCIPisInfinity(scip, -SCIPnlrowGetLhs(nlrows[i])) )
174 ++nnonconvexnlrows;
175 }
176 else if( SCIPnlrowGetCurvature(nlrows[i]) == SCIP_EXPRCURV_CONCAVE )
177 {
178 if( !SCIPisInfinity(scip, SCIPnlrowGetRhs(nlrows[i])) )
179 ++nnonconvexnlrows;
180 if( !SCIPisInfinity(scip, -SCIPnlrowGetLhs(nlrows[i])) )
181 ++nconvexnlrows;
182 }
183 else
184 {
185 if( !SCIPisInfinity(scip, SCIPnlrowGetRhs(nlrows[i])) )
186 ++nnonconvexnlrows;
187 if( !SCIPisInfinity(scip, -SCIPnlrowGetLhs(nlrows[i])) )
188 ++nnonconvexnlrows;
189 }
190 }
191
192 SCIPdebugMsg(scip, "nconvex=%d nnonconvex=%d nlinear=%d\n", nconvexnlrows, nnonconvexnlrows, nlinearnlrows);
193
194 return nconvexnlrows > 0
195 && nnonconvexnlrows > 0
196 && (SCIPisGE(scip, (SCIP_Real)nconvexnlrows, nnonconvexnlrows * propdata->minnonconvexfrac))
197 && (SCIPisGE(scip, (SCIP_Real)nconvexnlrows, nlinearnlrows * propdata->minlinearfrac));
198 }
199
200 /** filters variables which achieve their lower or dual bound in the current NLP solution */
201 static
filterCands(SCIP * scip,SCIP_PROPDATA * propdata)202 SCIP_RETCODE filterCands(
203 SCIP* scip, /**< SCIP data structure */
204 SCIP_PROPDATA* propdata /**< propagator data */
205 )
206 {
207 SCIP_Real* primal;
208 int i;
209
210 assert(propdata->currpos >= 0 && propdata->currpos < propdata->nlpinvars);
211 assert(SCIPnlpiGetSolstat(propdata->nlpi, propdata->nlpiprob) <= SCIP_NLPSOLSTAT_FEASIBLE);
212
213 SCIP_CALL( SCIPnlpiGetSolution(propdata->nlpi, propdata->nlpiprob, &primal, NULL, NULL, NULL, NULL) );
214 assert(primal != NULL);
215
216 /* we skip all candidates which have been processed already, i.e., starting at propdata->currpos + 1 */
217 for( i = propdata->currpos + 1; i < propdata->nlpinvars; ++i )
218 {
219 SCIP_VAR* var;
220 SCIP_Real val;
221 int varidx;
222
223 /* only uninteresting variables left -> stop filtering */
224 if( SCIPisLE(scip, propdata->nlscore[i], 0.0) )
225 break;
226
227 var = propdata->nlpivars[i];
228 assert(var != NULL && SCIPhashmapExists(propdata->var2nlpiidx, (void*)var));
229
230 varidx = SCIPhashmapGetImageInt(propdata->var2nlpiidx, (void*)var);
231 assert(SCIPgetVars(scip)[varidx] == var);
232 val = primal[varidx];
233
234 if( (propdata->status[i] & SOLVEDLB) == 0 && !SCIPisInfinity(scip, -val) /*lint !e641*/
235 && SCIPisFeasLE(scip, val, SCIPvarGetLbLocal(var)) )
236 {
237 SCIPdebugMsg(scip, "filter LB of %s in [%g,%g] with %g\n", SCIPvarGetName(var), SCIPvarGetLbLocal(var),
238 SCIPvarGetUbLocal(var), val);
239 propdata->status[i] |= SOLVEDLB; /*lint !e641*/
240 assert((propdata->status[i] & SOLVEDLB) != 0); /*lint !e641*/
241 }
242
243 if( (propdata->status[i] & SOLVEDUB) == 0 && !SCIPisInfinity(scip, val) /*lint !e641*/
244 && SCIPisFeasGE(scip, val, SCIPvarGetUbLocal(var)) )
245 {
246 SCIPdebugMsg(scip, "filter UB of %s in [%g,%g] with %g\n", SCIPvarGetName(var), SCIPvarGetLbLocal(var),
247 SCIPvarGetUbLocal(var), val);
248 propdata->status[i] |= SOLVEDUB; /*lint !e641*/
249 assert((propdata->status[i] & SOLVEDUB) != 0); /*lint !e641*/
250 }
251 }
252
253 return SCIP_OKAY;
254 }
255
256 /** tries to add a generalized variable bound by exploiting the dual solution of the last NLP solve (see @ref
257 * prop_nlobbt.h for more information)
258 */
259 static
addGenVBound(SCIP * scip,SCIP_PROPDATA * propdata,SCIP_VAR * var,int varidx,SCIP_BOUNDTYPE boundtype,SCIP_Real cutoffbound)260 SCIP_RETCODE addGenVBound(
261 SCIP* scip, /**< SCIP data structure */
262 SCIP_PROPDATA* propdata, /**< propagator data */
263 SCIP_VAR* var, /**< variable used in last NLP solve */
264 int varidx, /**< variable index in the propdata->nlpivars array */
265 SCIP_BOUNDTYPE boundtype, /**< type of bound provided by the genvbound */
266 SCIP_Real cutoffbound /**< cutoff bound */
267 )
268 {
269 SCIP_VAR** lvbvars;
270 SCIP_Real* lvbcoefs;
271 SCIP_Real* primal;
272 SCIP_Real* dual;
273 SCIP_Real* alpha;
274 SCIP_Real* beta;
275 SCIP_Real constant;
276 SCIP_Real mu;
277 int nlvbvars;
278 int i;
279
280 assert(propdata->genvboundprop != NULL);
281 assert(var != NULL);
282 assert(varidx >= 0 && varidx < propdata->nlpinvars);
283 assert(SCIPnlpiGetSolstat(propdata->nlpi, propdata->nlpiprob) <= SCIP_NLPSOLSTAT_LOCOPT);
284
285 SCIP_CALL( SCIPnlpiGetSolution(propdata->nlpi, propdata->nlpiprob, &primal, &dual, &alpha, &beta, NULL) );
286
287 /* not possible to generate genvbound if the duals for the propagated variable do not disappear */
288 if( !SCIPisFeasZero(scip, alpha[varidx] - beta[varidx]) )
289 return SCIP_OKAY;
290
291 SCIP_CALL( SCIPallocBufferArray(scip, &lvbcoefs, propdata->nlpinvars) );
292 SCIP_CALL( SCIPallocBufferArray(scip, &lvbvars, propdata->nlpinvars) );
293 constant = boundtype == SCIP_BOUNDTYPE_LOWER ? primal[varidx] : -primal[varidx];
294 mu = 0.0;
295 nlvbvars = 0;
296
297 /* collect coefficients of genvbound */
298 for( i = 0; i < propdata->nlpinvars; ++i )
299 {
300 if( !SCIPisZero(scip, beta[i] - alpha[i]) )
301 {
302 lvbvars[nlvbvars] = propdata->nlpivars[i];
303 lvbcoefs[nlvbvars] = beta[i] - alpha[i];
304 ++nlvbvars;
305
306 constant += (alpha[i] - beta[i]) * primal[i];
307 }
308 }
309
310 /* first dual multiplier corresponds to the cutoff row if cutoffbound < SCIPinfinity() */
311 if( !SCIPisInfinity(scip, cutoffbound) && SCIPisGT(scip, dual[0], 0.0) )
312 {
313 mu = dual[0];
314 constant += mu * cutoffbound;
315 }
316
317 /* add genvbound to genvbounds propagator */
318 if( !SCIPisInfinity(scip, REALABS(constant)) && (nlvbvars > 0 || SCIPisFeasGT(scip, mu, 0.0)) )
319 {
320 SCIP_CALL( SCIPgenVBoundAdd(scip, propdata->genvboundprop, lvbvars, var, lvbcoefs, nlvbvars, -mu, constant,
321 boundtype) );
322 SCIPdebugMsg(scip, "add genvbound for %s\n", SCIPvarGetName(var));
323 }
324
325 SCIPfreeBufferArray(scip, &lvbvars);
326 SCIPfreeBufferArray(scip, &lvbcoefs);
327
328 return SCIP_OKAY;
329 }
330
331 /** sets the objective function, solves the NLP, and tightens the given variable; after calling this function, the
332 * objective function is set to zero
333 *
334 * @note function assumes that objective function is zero
335 */
336 static
solveNlp(SCIP * scip,SCIP_PROPDATA * propdata,SCIP_VAR * var,int varidx,SCIP_BOUNDTYPE boundtype,int * nlpiter,SCIP_RESULT * result)337 SCIP_RETCODE solveNlp(
338 SCIP* scip, /**< SCIP data structure */
339 SCIP_PROPDATA* propdata, /**< propagator data */
340 SCIP_VAR* var, /**< variable to propagate */
341 int varidx, /**< variable index in the propdata->nlpivars array */
342 SCIP_BOUNDTYPE boundtype, /**< minimize or maximize var? */
343 int* nlpiter, /**< buffer to store the total number of nlp iterations */
344 SCIP_RESULT* result /**< pointer to store result */
345 )
346 {
347 SCIP_Real timelimit;
348 SCIP_Real* primal;
349 SCIP_Real obj;
350 int iterlimit;
351
352 #ifdef SCIP_DEBUG
353 SCIP_Real oldlb;
354 SCIP_Real oldub;
355
356 oldlb = SCIPvarGetLbLocal(var);
357 oldub = SCIPvarGetUbLocal(var);
358 #endif
359
360 assert(var != NULL);
361 assert(varidx >= 0 && varidx < propdata->nlpinvars);
362 assert(result != NULL && *result != SCIP_CUTOFF);
363
364 *nlpiter = 0;
365
366 /* set time and iteration limit */
367 SCIP_CALL( SCIPgetRealParam(scip, "limits/time", &timelimit) );
368 if( !SCIPisInfinity(scip, timelimit) )
369 {
370 timelimit -= SCIPgetSolvingTime(scip);
371 if( timelimit <= 0.0 )
372 {
373 SCIPdebugMsg(scip, "skip NLP solve; no time left\n");
374 return SCIP_OKAY;
375 }
376 }
377 if( propdata->nlptimelimit > 0.0 )
378 timelimit = MIN(propdata->nlptimelimit, timelimit);
379 iterlimit = propdata->nlpiterlimit > 0 ? propdata->nlpiterlimit : INT_MAX;
380 SCIP_CALL( SCIPnlpiSetRealPar(propdata->nlpi, propdata->nlpiprob, SCIP_NLPPAR_TILIM, timelimit) );
381 SCIP_CALL( SCIPnlpiSetIntPar(propdata->nlpi, propdata->nlpiprob, SCIP_NLPPAR_ITLIM, iterlimit) );
382
383 /* set corresponding objective coefficient and solve NLP */
384 obj = boundtype == SCIP_BOUNDTYPE_LOWER ? 1.0 : -1.0;
385 SCIP_CALL( SCIPnlpiSetObjective(propdata->nlpi, propdata->nlpiprob, 1, &varidx, &obj, 0, NULL, NULL, NULL, 0.0) );
386
387 SCIPdebugMsg(scip, "solve var=%s boundtype=%d nlscore=%g\n", SCIPvarGetName(var), boundtype,
388 propdata->nlscore[propdata->currpos]);
389 SCIP_CALL( SCIPnlpiSolve(propdata->nlpi, propdata->nlpiprob) );
390 SCIPdebugMsg(scip, "NLP solstat = %d\n", SCIPnlpiGetSolstat(propdata->nlpi, propdata->nlpiprob));
391
392 /* collect NLP statistics */
393 assert(propdata->nlpstatistics != NULL);
394 SCIP_CALL( SCIPnlpiGetStatistics(propdata->nlpi, propdata->nlpiprob, propdata->nlpstatistics) );
395 *nlpiter = SCIPnlpStatisticsGetNIterations(propdata->nlpstatistics);
396 SCIPdebugMsg(scip, "iterations %d time %g\n", *nlpiter, SCIPnlpStatisticsGetTotalTime(propdata->nlpstatistics));
397
398 /* filter bound candidates first, otherwise we do not have access to the primal solution values */
399 if( SCIPnlpiGetSolstat(propdata->nlpi, propdata->nlpiprob) <= SCIP_NLPSOLSTAT_FEASIBLE )
400 {
401 SCIP_CALL( filterCands(scip, propdata) );
402 }
403
404 /* try to tighten variable bound */
405 if( SCIPnlpiGetSolstat(propdata->nlpi, propdata->nlpiprob) <= SCIP_NLPSOLSTAT_LOCOPT )
406 {
407 SCIP_Bool tightened;
408 SCIP_Bool infeasible;
409
410 /* try to add a genvbound in the root node */
411 if( propdata->genvboundprop != NULL && SCIPgetDepth(scip) == 0 )
412 {
413 SCIP_CALL( addGenVBound(scip, propdata, var, varidx, boundtype, SCIPgetCutoffbound(scip)) );
414 }
415
416 SCIP_CALL( SCIPnlpiGetSolution(propdata->nlpi, propdata->nlpiprob, &primal, NULL, NULL, NULL, NULL) );
417
418 if( boundtype == SCIP_BOUNDTYPE_LOWER )
419 {
420 SCIP_CALL( SCIPtightenVarLb(scip, var, primal[varidx], FALSE, &infeasible, &tightened) );
421 }
422 else
423 {
424 SCIP_CALL( SCIPtightenVarUb(scip, var, primal[varidx], FALSE, &infeasible, &tightened) );
425 }
426
427 if( infeasible )
428 {
429 SCIPdebugMsg(scip, "detect infeasibility after propagating %s\n", SCIPvarGetName(var));
430 *result = SCIP_CUTOFF;
431 }
432 else if( tightened )
433 {
434 SCIP_Real lb;
435 SCIP_Real ub;
436
437 *result = SCIP_REDUCEDDOM;
438
439 /* update bounds in NLP */
440 lb = SCIPvarGetLbLocal(var);
441 ub = SCIPvarGetUbLocal(var);
442 SCIP_CALL( SCIPnlpiChgVarBounds(propdata->nlpi, propdata->nlpiprob, 1, &varidx, &lb, &ub) );
443
444 #ifdef SCIP_DEBUG
445 SCIPdebugMsg(scip, "tightened bounds of %s from [%g,%g] to [%g,%g]\n", SCIPvarGetName(var), oldlb, oldub, lb, ub);
446 #endif
447 }
448 }
449
450 /* reset objective function */
451 obj = 0.0;
452 SCIP_CALL( SCIPnlpiSetObjective(propdata->nlpi, propdata->nlpiprob, 1, &varidx, &obj, 0, NULL, NULL, NULL, 0.0) );
453
454 return SCIP_OKAY;
455 }
456
457 /** main method of the propagator
458 *
459 * creates a convex NLP relaxation and solves the OBBT-NLPs for each possible candidate;
460 * binary and variables with a small domain will be ignored to reduce the computational cost of the propagator; after
461 * solving each NLP we filter out all variable candidates which are on their lower or upper bound; candidates with a
462 * larger number of occurrences are preferred
463 */
464 static
applyNlobbt(SCIP * scip,SCIP_PROPDATA * propdata,SCIP_RESULT * result)465 SCIP_RETCODE applyNlobbt(
466 SCIP* scip, /**< SCIP data structure */
467 SCIP_PROPDATA* propdata, /**< propagation data */
468 SCIP_RESULT* result /**< pointer to store result */
469 )
470 {
471 int nlpiterleft;
472
473 assert(result != NULL);
474 assert(!propdata->skipprop);
475 assert(SCIPgetNNlpis(scip) > 0);
476
477 *result = SCIP_DIDNOTRUN;
478
479 if( propdata->nlpiprob == NULL && !isNlobbtApplicable(scip, propdata) )
480 {
481 /* do not call the propagator anymore (except after a restart) */
482 SCIPdebugMsg(scip, "nlobbt propagator is not applicable\n");
483 propdata->skipprop = TRUE;
484 return SCIP_OKAY;
485 }
486
487 *result = SCIP_DIDNOTFIND;
488
489 /* compute NLP iteration limit */
490 if( propdata->itlimitfactor > 0.0 )
491 nlpiterleft = (int)(propdata->itlimitfactor * SCIPgetNRootLPIterations(scip));
492 else
493 nlpiterleft = INT_MAX;
494
495 /* recompute NLP relaxation if the variable set changed */
496 if( propdata->nlpiprob != NULL && SCIPgetNVars(scip) != propdata->nlpinvars )
497 {
498 SCIP_CALL( propdataClear(scip, propdata) );
499 assert(propdata->nlpiprob == NULL);
500 }
501
502 /* create or update NLP relaxation */
503 if( propdata->nlpiprob == NULL )
504 {
505 int i;
506
507 propdata->nlpinvars = SCIPgetNVars(scip);
508 propdata->nlpi = SCIPgetNlpis(scip)[0];
509 assert(propdata->nlpi != NULL);
510
511 SCIP_CALL( SCIPnlpiCreateProblem(propdata->nlpi, &propdata->nlpiprob, "nlobbt-nlp") );
512 SCIP_CALL( SCIPhashmapCreate(&propdata->var2nlpiidx, SCIPblkmem(scip), propdata->nlpinvars) );
513 SCIP_CALL( SCIPduplicateBlockMemoryArray(scip, &propdata->nlpivars, SCIPgetVars(scip), propdata->nlpinvars) ); /*lint !e666*/
514 SCIP_CALL( SCIPallocBlockMemoryArray(scip, &propdata->nlscore, propdata->nlpinvars) );
515 SCIP_CALL( SCIPallocBlockMemoryArray(scip, &propdata->status, propdata->nlpinvars) );
516
517 SCIP_CALL( SCIPcreateNlpiProb(scip, propdata->nlpi, SCIPgetNLPNlRows(scip), SCIPgetNNLPNlRows(scip),
518 propdata->nlpiprob, propdata->var2nlpiidx, NULL, propdata->nlscore, SCIPgetCutoffbound(scip), FALSE, TRUE) );
519
520 /* initialize bound status; perturb nlscores by a factor which ensures that zero scores remain zero */
521 assert(propdata->randnumgen != NULL);
522 for( i = 0; i < propdata->nlpinvars; ++i )
523 {
524 propdata->status[i] = UNSOLVED; /*lint !e641*/
525 propdata->nlscore[i] *= 1.0 + SCIPrandomGetReal(propdata->randnumgen, SCIPfeastol(scip), 2.0 * SCIPfeastol(scip));
526 }
527
528 /* add rows of the LP */
529 if( SCIPgetDepth(scip) == 0 )
530 {
531 SCIP_CALL( SCIPaddNlpiProbRows(scip, propdata->nlpi, propdata->nlpiprob, propdata->var2nlpiidx,
532 SCIPgetLPRows(scip), SCIPgetNLPRows(scip)) );
533 }
534 }
535 else
536 {
537 SCIP_CALL( SCIPupdateNlpiProb(scip, propdata->nlpi, propdata->nlpiprob, propdata->var2nlpiidx,
538 propdata->nlpivars, propdata->nlpinvars, SCIPgetCutoffbound(scip)) );
539 }
540
541 assert(propdata->nlpiprob != NULL);
542 assert(propdata->var2nlpiidx != NULL);
543 assert(propdata->nlpivars != NULL);
544 assert(propdata->nlscore != NULL);
545
546 /* sort variables w.r.t. their nlscores if we did not solve any NLP for this node */
547 if( propdata->currpos == 0 )
548 {
549 SCIPsortDownRealIntPtr(propdata->nlscore, propdata->status, (void**)propdata->nlpivars, propdata->nlpinvars);
550 }
551
552 /* set parameters of NLP solver */
553 SCIP_CALL( SCIPnlpiSetRealPar(propdata->nlpi, propdata->nlpiprob, SCIP_NLPPAR_FEASTOL,
554 SCIPfeastol(scip) * propdata->feastolfac) );
555 SCIP_CALL( SCIPnlpiSetRealPar(propdata->nlpi, propdata->nlpiprob, SCIP_NLPPAR_RELOBJTOL,
556 SCIPfeastol(scip) * propdata->relobjtolfac) );
557 SCIP_CALL( SCIPnlpiSetIntPar(propdata->nlpi, propdata->nlpiprob, SCIP_NLPPAR_VERBLEVEL, propdata->nlpverblevel) );
558
559 /* main propagation loop */
560 while( propdata->currpos < propdata->nlpinvars
561 && nlpiterleft > 0
562 && SCIPisGT(scip, propdata->nlscore[propdata->currpos], 0.0)
563 && *result != SCIP_CUTOFF
564 && !SCIPisStopped(scip) )
565 {
566 SCIP_VAR* var;
567 int varidx;
568 int iters;
569
570 var = propdata->nlpivars[propdata->currpos];
571 assert(var != NULL);
572
573 /* skip binary or almost fixed variables */
574 if( SCIPvarGetType(var) == SCIP_VARTYPE_BINARY
575 || SCIPisFeasEQ(scip, SCIPvarGetLbLocal(var), SCIPvarGetUbLocal(var)) )
576 {
577 ++(propdata->currpos);
578 continue;
579 }
580
581 SCIPdebugMsg(scip, "iterations left %d\n", nlpiterleft);
582
583 /* get index of var in the nlpi */
584 assert(SCIPhashmapExists(propdata->var2nlpiidx, (void*)var) );
585 varidx = SCIPhashmapGetImageInt(propdata->var2nlpiidx, (void*)var);
586 assert(var == SCIPgetVars(scip)[varidx]);
587
588 /* case: minimize var */
589 if( (propdata->status[propdata->currpos] & SOLVEDLB) == 0 ) /*lint !e641*/
590 {
591 SCIP_CALL( solveNlp(scip, propdata, var, varidx, SCIP_BOUNDTYPE_LOWER, &iters, result) );
592 nlpiterleft -= iters;
593 }
594
595 /* case: maximize var */
596 if( *result != SCIP_CUTOFF && (propdata->status[propdata->currpos] & SOLVEDUB) == 0 ) /*lint !e641*/
597 {
598 SCIP_CALL( solveNlp(scip, propdata, var, varidx, SCIP_BOUNDTYPE_UPPER, &iters, result) );
599 nlpiterleft -= iters;
600 }
601
602 /* update the current position */
603 ++(propdata->currpos);
604 }
605
606 return SCIP_OKAY;
607 }
608
609 /*
610 * Callback methods of propagator
611 */
612
613 /** destructor of propagator to free user data (called when SCIP is exiting) */
614 static
SCIP_DECL_PROPFREE(propFreeNlobbt)615 SCIP_DECL_PROPFREE(propFreeNlobbt)
616 { /*lint --e{715}*/
617 SCIP_PROPDATA* propdata;
618
619 propdata = SCIPpropGetData(prop);
620 assert(propdata != NULL);
621
622 SCIP_CALL( propdataClear(scip, propdata) );
623 SCIPfreeBlockMemory(scip, &propdata);
624 SCIPpropSetData(prop, NULL);
625
626 return SCIP_OKAY;
627 }
628
629 /** solving process initialization method of propagator (called when branch and bound process is about to begin) */
630 static
SCIP_DECL_PROPINITSOL(propInitsolNlobbt)631 SCIP_DECL_PROPINITSOL(propInitsolNlobbt)
632 { /*lint --e{715}*/
633 SCIP_PROPDATA* propdata;
634
635 assert(scip != NULL);
636 assert(prop != NULL);
637
638 propdata = SCIPpropGetData(prop);
639 assert(propdata != NULL);
640
641 /* if genvbounds propagator is not available, we cannot create genvbounds */
642 propdata->genvboundprop = SCIPfindProp(scip, "genvbounds");
643
644 SCIP_CALL( SCIPcreateRandom(scip, &propdata->randnumgen,
645 DEFAULT_RANDSEED, TRUE) );
646 SCIP_CALL( SCIPnlpStatisticsCreate(SCIPblkmem(scip), &propdata->nlpstatistics) );
647 propdata->lastnode = -1;
648
649 return SCIP_OKAY;
650 }
651
652 /** solving process deinitialization method of propagator (called before branch and bound process data is freed) */
653 static
SCIP_DECL_PROPEXITSOL(propExitsolNlobbt)654 SCIP_DECL_PROPEXITSOL(propExitsolNlobbt)
655 { /*lint --e{715}*/
656 SCIP_PROPDATA* propdata;
657
658 propdata = SCIPpropGetData(prop);
659 assert(propdata != NULL);
660
661 SCIPnlpStatisticsFree(SCIPblkmem(scip), &propdata->nlpstatistics);
662 SCIPfreeRandom(scip, &propdata->randnumgen);
663
664 SCIP_CALL( propdataClear(scip, propdata) );
665
666 return SCIP_OKAY;
667 }
668
669 /** execution method of propagator */
670 static
SCIP_DECL_PROPEXEC(propExecNlobbt)671 SCIP_DECL_PROPEXEC(propExecNlobbt)
672 { /*lint --e{715}*/
673 SCIP_PROPDATA* propdata;
674
675 *result = SCIP_DIDNOTRUN;
676
677 propdata = SCIPpropGetData(prop);
678 assert(propdata != NULL);
679
680 if( propdata->skipprop || SCIPgetStage(scip) != SCIP_STAGE_SOLVING || SCIPinRepropagation(scip)
681 || SCIPinProbing(scip) || SCIPinDive(scip) || !SCIPallowWeakDualReds(scip) || SCIPgetNNlpis(scip) == 0 )
682 {
683 SCIPdebugMsg(scip, "skip nlobbt propagator\n");
684 return SCIP_OKAY;
685 }
686
687 /* only run if LP all columns are in the LP, e.g., do not run if pricers are active */
688 if( !SCIPallColsInLP(scip) )
689 {
690 SCIPdebugMsg(scip, "not all columns in LP, skipping obbt\n");
691 return SCIP_OKAY;
692 }
693
694 /* do not run if SCIP does not have constructed an NLP */
695 if( !SCIPisNLPConstructed(scip) )
696 {
697 SCIPdebugMsg(scip, "NLP not constructed, skipping nlobbt\n");
698 return SCIP_OKAY;
699 }
700
701 /* consider all variables again if we process a new node */
702 if( SCIPnodeGetNumber(SCIPgetCurrentNode(scip)) != propdata->lastnode )
703 {
704 propdata->lastnode = SCIPnodeGetNumber(SCIPgetCurrentNode(scip));
705 propdata->currpos = 0;
706 }
707
708 /* call main procedure of nonlinear OBBT propagator */
709 SCIP_CALL( applyNlobbt(scip, propdata, result) );
710
711 return SCIP_OKAY;
712 }
713
714 /*
715 * propagator specific interface methods
716 */
717
718 /** creates the nlobbt propagator and includes it in SCIP */
SCIPincludePropNlobbt(SCIP * scip)719 SCIP_RETCODE SCIPincludePropNlobbt(
720 SCIP* scip /**< SCIP data structure */
721 )
722 {
723 SCIP_PROPDATA* propdata;
724 SCIP_PROP* prop;
725
726 propdata = NULL;
727 prop = NULL;
728
729 SCIP_CALL( SCIPallocBlockMemory(scip, &propdata) );
730 assert(propdata != NULL);
731 BMSclearMemory(propdata);
732
733 SCIP_CALL( SCIPincludePropBasic(scip, &prop, PROP_NAME, PROP_DESC, PROP_PRIORITY, PROP_FREQ, PROP_DELAY, PROP_TIMING,
734 propExecNlobbt, propdata) );
735 assert(prop != NULL);
736
737 SCIP_CALL( SCIPsetPropFree(scip, prop, propFreeNlobbt) );
738 SCIP_CALL( SCIPsetPropInitsol(scip, prop, propInitsolNlobbt) );
739 SCIP_CALL( SCIPsetPropExitsol(scip, prop, propExitsolNlobbt) );
740
741 SCIP_CALL( SCIPaddRealParam(scip, "propagating/" PROP_NAME "/feastolfac",
742 "factor for NLP feasibility tolerance",
743 &propdata->feastolfac, TRUE, DEFAULT_FEASTOLFAC, 0.0, 1.0, NULL, NULL) );
744
745 SCIP_CALL( SCIPaddRealParam(scip, "propagating/" PROP_NAME "/relobjtolfac",
746 "factor for NLP relative objective tolerance",
747 &propdata->relobjtolfac, TRUE, DEFAULT_RELOBJTOLFAC, 0.0, 1.0, NULL, NULL) );
748
749 SCIP_CALL( SCIPaddRealParam(scip, "propagating/" PROP_NAME "/minnonconvexfrac",
750 "(#convex nlrows)/(#nonconvex nlrows) threshold to apply propagator",
751 &propdata->minnonconvexfrac, TRUE, DEFAULT_MINNONCONVEXFRAC, 0.0, SCIPinfinity(scip), NULL, NULL) );
752
753 SCIP_CALL( SCIPaddRealParam(scip, "propagating/" PROP_NAME "/minlinearfrac",
754 "minimum (#convex nlrows)/(#linear nlrows) threshold to apply propagator",
755 &propdata->minlinearfrac, TRUE, DEFAULT_MINLINEARFRAC, 0.0, SCIPinfinity(scip), NULL, NULL) );
756
757 SCIP_CALL( SCIPaddBoolParam(scip, "propagating/" PROP_NAME "/addlprows",
758 "should non-initial LP rows be used?",
759 &propdata->addlprows, FALSE, DEFAULT_ADDLPROWS, NULL, NULL) );
760
761 SCIP_CALL( SCIPaddIntParam(scip, "propagating/" PROP_NAME "/nlpiterlimit",
762 "iteration limit of NLP solver; 0 for no limit",
763 &propdata->nlpiterlimit, TRUE, DEFAULT_NLPITERLIMIT, 0, INT_MAX, NULL, NULL) );
764
765 SCIP_CALL( SCIPaddRealParam(scip, "propagating/" PROP_NAME "/nlptimelimit",
766 "time limit of NLP solver; 0.0 for no limit",
767 &propdata->nlptimelimit, TRUE, DEFAULT_NLPTIMELIMIT, 0.0, SCIP_REAL_MAX, NULL, NULL) );
768
769 SCIP_CALL( SCIPaddIntParam(scip, "propagating/" PROP_NAME "/nlpverblevel",
770 "verbosity level of NLP solver",
771 &propdata->nlpverblevel, TRUE, DEFAULT_NLPVERLEVEL, 0, 5, NULL, NULL) );
772
773 SCIP_CALL( SCIPaddRealParam(scip, "propagating/" PROP_NAME "/itlimitfactor",
774 "LP iteration limit for nlobbt will be this factor times total LP iterations in root node",
775 &propdata->itlimitfactor, TRUE, DEFAULT_ITLIMITFACTOR, 0.0, SCIP_REAL_MAX, NULL, NULL) );
776
777 return SCIP_OKAY;
778 }
779