1 /* * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * */
2 /* */
3 /* This file is part of the program and library */
4 /* SCIP --- Solving Constraint Integer Programs */
5 /* */
6 /* Copyright (C) 2002-2022 Konrad-Zuse-Zentrum */
7 /* fuer Informationstechnik Berlin */
8 /* */
9 /* SCIP is distributed under the terms of the ZIB Academic License. */
10 /* */
11 /* You should have received a copy of the ZIB Academic License */
12 /* along with SCIP; see the file COPYING. If not email to scip@zib.de. */
13 /* */
14 /* * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * */
15
16 /**@file nlhdlr_perspective.c
17 * @ingroup DEFPLUGINS_NLHDLR
18 * @brief perspective nonlinear handler
19 * @author Ksenia Bestuzheva
20 */
21
22 #include <string.h>
23
24 #include "scip/nlhdlr_perspective.h"
25 #include "scip/cons_nonlinear.h"
26 #include "scip/scip_sol.h"
27 #include "scip/pub_misc_rowprep.h"
28 #include "scip/nlhdlr.h"
29
30 /* fundamental nonlinear handler properties */
31 #define NLHDLR_NAME "perspective"
32 #define NLHDLR_DESC "perspective handler for expressions"
33 #define NLHDLR_DETECTPRIORITY -20 /**< detect last so that to make use of what other handlers detected */
34 #define NLHDLR_ENFOPRIORITY 125 /**< enforce first because perspective cuts are always stronger */
35
36 #define DEFAULT_MAXPROPROUNDS 1 /**< maximal number of propagation rounds in probing */
37 #define DEFAULT_MINDOMREDUCTION 0.1 /**< minimal relative reduction in a variable's domain for applying probing */
38 #define DEFAULT_MINVIOLPROBING 1e-05 /**< minimal violation w.r.t. auxiliary variables for applying probing */
39 #define DEFAULT_PROBINGONLYINSEPA TRUE /**< whether to do probing only in separation loop */
40 #define DEFAULT_PROBINGFREQ 1 /**< probing frequency (-1 - no probing, 0 - root node only) */
41 #define DEFAULT_CONVEXONLY FALSE /**< whether perspective cuts are added only for convex expressions */
42 #define DEFAULT_TIGHTENBOUNDS TRUE /**< whether variable semicontinuity is used to tighten variable bounds */
43 #define DEFAULT_ADJREFPOINT TRUE /**< whether to adjust the reference point if indicator is not 1 */
44
45 /*
46 * Data structures
47 */
48
49 /** data structure to store information of a semicontinuous variable
50 *
51 * For a variable x (not stored in the struct), this stores the data of nbnds implications
52 * bvars[i] = 0 -> x = vals[i]
53 * bvars[i] = 1 -> lbs[i] <= x <= ubs[i]
54 * where bvars[i] are binary variables.
55 */
56 struct SCVarData
57 {
58 SCIP_Real* vals0; /**< values of the variable when the corresponding bvars[i] = 0 */
59 SCIP_Real* lbs1; /**< global lower bounds of the variable when the corresponding bvars[i] = 1 */
60 SCIP_Real* ubs1; /**< global upper bounds of the variable when the corresponding bvars[i] = 1 */
61 SCIP_VAR** bvars; /**< the binary variables on which the variable domain depends */
62 int nbnds; /**< number of suitable on/off bounds the var has */
63 int bndssize; /**< size of the arrays */
64 };
65 typedef struct SCVarData SCVARDATA;
66
67 /** nonlinear handler expression data
68 *
69 * For an expression expr (not stored in the struct), this stores the data of nindicators implications
70 * indicators[i] = 0 -> expr = exprvals[0]
71 * where indicators[i] is an indicator (binary) variable, corresponding to some bvars entry in SCVarData.
72 *
73 * Also stores the variables the expression depends on.
74 */
75 struct SCIP_NlhdlrExprData
76 {
77 SCIP_Real* exprvals0; /**< 'off' values of the expression for each indicator variable */
78 SCIP_VAR** vars; /**< expression variables (both original and auxiliary) */
79 int nvars; /**< total number of variables in the expression */
80 int varssize; /**< size of the vars array */
81 SCIP_VAR** indicators; /**< all indicator variables for the expression */
82 int nindicators; /**< number of indicator variables */
83 };
84
85 /** nonlinear handler data */
86 struct SCIP_NlhdlrData
87 {
88 SCIP_HASHMAP* scvars; /**< maps semicontinuous variables to their on/off bounds (SCVarData) */
89
90 /* parameters */
91 int maxproprounds; /**< maximal number of propagation rounds in probing */
92 SCIP_Real mindomreduction; /**< minimal relative reduction in a variable's domain for applying probing */
93 SCIP_Real minviolprobing; /**< minimal violation w.r.t. auxiliary variables for applying probing */
94 SCIP_Bool probingonlyinsepa; /**< whether to do probing only in separation loop */
95 int probingfreq; /**< if and when to do probing */
96 SCIP_Bool convexonly; /**< whether perspective cuts are added only for convex expressions */
97 SCIP_Bool tightenbounds; /**< whether variable semicontinuity is used to tighten variable bounds */
98 SCIP_Bool adjrefpoint; /**< whether to adjust the reference point if indicator is not 1 */
99 };
100
101 /*
102 * Local methods
103 */
104
105 /*
106 * Helper methods for working with nlhdlrExprData
107 */
108
109 /** frees nlhdlrexprdata structure */
110 static
111 SCIP_RETCODE freeNlhdlrExprData(
112 SCIP* scip, /**< SCIP data structure */
113 SCIP_NLHDLREXPRDATA* nlhdlrexprdata /**< nlhdlr expression data */
114 )
115 {
116 int v;
117
118 if( nlhdlrexprdata->nindicators != 0 )
119 {
120 assert(nlhdlrexprdata->indicators != NULL);
121 for( v = nlhdlrexprdata->nindicators - 1; v >= 0; --v )
122 {
123 SCIP_CALL( SCIPreleaseVar(scip, &(nlhdlrexprdata->indicators[v])) );
124 }
125 SCIPfreeBlockMemoryArray(scip, &(nlhdlrexprdata->indicators), nlhdlrexprdata->nindicators);
126 SCIPfreeBlockMemoryArrayNull(scip, &(nlhdlrexprdata->exprvals0), nlhdlrexprdata->nindicators);
127 }
128
129 for( v = nlhdlrexprdata->nvars - 1; v >= 0; --v )
130 {
131 SCIP_CALL( SCIPreleaseVar(scip, &(nlhdlrexprdata->vars[v])) );
132 }
133 SCIPfreeBlockMemoryArrayNull(scip, &nlhdlrexprdata->vars, nlhdlrexprdata->varssize);
134
135 return SCIP_OKAY;
136 }
137
138 /* remove an indicator from nlhdlr expression data */
139 static
140 SCIP_RETCODE removeIndicator(
141 SCIP* scip, /**< SCIP data structure */
142 SCIP_NLHDLREXPRDATA* nlexprdata, /**< nlhdlr expression data */
143 int pos /**< position of the indicator */
144 )
145 {
146 int i;
147
148 assert(pos >= 0 && pos < nlexprdata->nindicators);
149
150 SCIP_CALL( SCIPreleaseVar(scip, &nlexprdata->indicators[pos]) );
151 for( i = pos; i < nlexprdata->nindicators - 1; ++i )
152 {
153 nlexprdata->indicators[i] = nlexprdata->indicators[i+1];
154 }
155
156 --nlexprdata->nindicators;
157
158 return SCIP_OKAY;
159 }
160
161 /** adds an auxiliary variable to the vars array in nlhdlrexprdata */
162 static
163 SCIP_RETCODE addAuxVar(
164 SCIP* scip, /**< SCIP data structure */
165 SCIP_NLHDLREXPRDATA* nlhdlrexprdata, /**< nlhdlr expression data */
166 SCIP_HASHMAP* auxvarmap, /**< hashmap linking auxvars to positions in nlhdlrexprdata->vars */
167 SCIP_VAR* auxvar /**< variable to be added */
168 )
169 {
170 int pos;
171 int newsize;
172
173 assert(nlhdlrexprdata != NULL);
174 assert(auxvar != NULL);
175
176 pos = SCIPhashmapGetImageInt(auxvarmap, (void*) auxvar);
177
178 if( pos != INT_MAX )
179 return SCIP_OKAY;
180
181 /* ensure size */
182 if( nlhdlrexprdata->nvars + 1 > nlhdlrexprdata->varssize )
183 {
184 newsize = SCIPcalcMemGrowSize(scip, nlhdlrexprdata->nvars + 1);
185 SCIP_CALL( SCIPreallocBlockMemoryArray(scip, &nlhdlrexprdata->vars, nlhdlrexprdata->varssize, newsize) );
186 nlhdlrexprdata->varssize = newsize;
187 }
188 assert(nlhdlrexprdata->nvars + 1 <= nlhdlrexprdata->varssize);
189
190 nlhdlrexprdata->vars[nlhdlrexprdata->nvars] = auxvar;
191 SCIP_CALL( SCIPcaptureVar(scip, auxvar) );
192 SCIP_CALL( SCIPhashmapSetImageInt(auxvarmap, (void*) auxvar, nlhdlrexprdata->nvars) );
193 ++(nlhdlrexprdata->nvars);
194
195 return SCIP_OKAY;
196 }
197
198 /*
199 * Semicontinuous variable methods
200 */
201
202 /** adds an indicator to the data of a semicontinuous variable */
203 static
204 SCIP_RETCODE addSCVarIndicator(
205 SCIP* scip, /**< SCIP data structure */
206 SCVARDATA* scvdata, /**< semicontinuous variable data */
207 SCIP_VAR* indicator, /**< indicator to be added */
208 SCIP_Real val0, /**< value of the variable when indicator == 0 */
209 SCIP_Real lb1, /**< lower bound of the variable when indicator == 1 */
210 SCIP_Real ub1 /**< upper bound of the variable when indicator == 1 */
211 )
212 {
213 int newsize;
214 int i;
215 SCIP_Bool found;
216 int pos;
217
218 assert(scvdata != NULL);
219 assert(indicator != NULL);
220
221 /* find the position where to insert */
222 if( scvdata->bvars == NULL )
223 {
224 assert(scvdata->nbnds == 0 && scvdata->bndssize == 0);
225 found = FALSE;
226 pos = 0;
227 }
228 else
229 {
230 found = SCIPsortedvecFindPtr((void**)scvdata->bvars, SCIPvarComp, (void*)indicator, scvdata->nbnds, &pos);
231 }
232
233 if( found )
234 return SCIP_OKAY;
235
236 /* ensure sizes */
237 if( scvdata->nbnds + 1 > scvdata->bndssize )
238 {
239 newsize = SCIPcalcMemGrowSize(scip, scvdata->nbnds + 1);
240 SCIP_CALL( SCIPreallocBlockMemoryArray(scip, &scvdata->bvars, scvdata->bndssize, newsize) );
241 SCIP_CALL( SCIPreallocBlockMemoryArray(scip, &scvdata->vals0, scvdata->bndssize, newsize) );
242 SCIP_CALL( SCIPreallocBlockMemoryArray(scip, &scvdata->lbs1, scvdata->bndssize, newsize) );
243 SCIP_CALL( SCIPreallocBlockMemoryArray(scip, &scvdata->ubs1, scvdata->bndssize, newsize) );
244 scvdata->bndssize = newsize;
245 }
246 assert(scvdata->nbnds + 1 <= scvdata->bndssize);
247 assert(scvdata->bvars != NULL);
248
249 /* move entries if needed */
250 for( i = scvdata->nbnds; i > pos; --i )
251 {
252 scvdata->bvars[i] = scvdata->bvars[i-1];
253 scvdata->vals0[i] = scvdata->vals0[i-1];
254 scvdata->lbs1[i] = scvdata->lbs1[i-1];
255 scvdata->ubs1[i] = scvdata->ubs1[i-1];
256 }
257
258 scvdata->bvars[pos] = indicator;
259 scvdata->vals0[pos] = val0;
260 scvdata->lbs1[pos] = lb1;
261 scvdata->ubs1[pos] = ub1;
262 ++scvdata->nbnds;
263
264 return SCIP_OKAY;
265 }
266
267 /** find scvardata of var and position of indicator in it
268 *
269 * If indicator is not there, returns NULL.
270 */
271 static
272 SCVARDATA* getSCVarDataInd(
273 SCIP_HASHMAP* scvars, /**< hashmap linking variables to scvardata */
274 SCIP_VAR* var, /**< variable */
275 SCIP_VAR* indicator, /**< indicator variable */
276 int* pos /**< pointer to store the position of indicator */
277 )
278 {
279 SCIP_Bool exists;
280 SCVARDATA* scvdata;
281
282 assert(var != NULL);
283 assert(scvars != NULL);
284 assert(indicator != NULL);
285
286 scvdata = (SCVARDATA*) SCIPhashmapGetImage(scvars, (void*)var);
287 if( scvdata != NULL )
288 {
289 /* look for the indicator variable */
290 exists = SCIPsortedvecFindPtr((void**)scvdata->bvars, SCIPvarComp, (void*)indicator, scvdata->nbnds, pos);
291 if( !exists )
292 return NULL;
293
294 return scvdata;
295 }
296
297 return NULL;
298 }
299
300 /** checks if a variable is semicontinuous and, if needed, updates the scvars hashmap
301 *
302 * A variable \f$x\f$ is semicontinuous if its bounds depend on at least one binary variable called the indicator,
303 * and indicator = 0 ⇒ \f$x = x^0\f$ for some real constant \f$x^0\f$.
304 */
305 static
306 SCIP_RETCODE varIsSemicontinuous(
307 SCIP* scip, /**< SCIP data structure */
308 SCIP_VAR* var, /**< the variable to check */
309 SCIP_HASHMAP* scvars, /**< semicontinuous variable information */
310 SCIP_Bool* result /**< buffer to store whether var is semicontinuous */
311 )
312 {
313 SCIP_Real lb0;
314 SCIP_Real ub0;
315 SCIP_Real lb1;
316 SCIP_Real ub1;
317 SCIP_Real glb;
318 SCIP_Real gub;
319 SCIP_Bool exists;
320 int c;
321 int pos;
322 SCIP_VAR** vlbvars;
323 SCIP_VAR** vubvars;
324 SCIP_Real* vlbcoefs;
325 SCIP_Real* vubcoefs;
326 SCIP_Real* vlbconstants;
327 SCIP_Real* vubconstants;
328 int nvlbs;
329 int nvubs;
330 SCVARDATA* scvdata;
331 SCIP_VAR* bvar;
332
333 assert(scip != NULL);
334 assert(var != NULL);
335 assert(scvars != NULL);
336 assert(result != NULL);
337
338 scvdata = (SCVARDATA*) SCIPhashmapGetImage(scvars, (void*)var);
339 if( scvdata != NULL )
340 {
341 *result = TRUE;
342 return SCIP_OKAY;
343 }
344
345 vlbvars = SCIPvarGetVlbVars(var);
346 vubvars = SCIPvarGetVubVars(var);
347 vlbcoefs = SCIPvarGetVlbCoefs(var);
348 vubcoefs = SCIPvarGetVubCoefs(var);
349 vlbconstants = SCIPvarGetVlbConstants(var);
350 vubconstants = SCIPvarGetVubConstants(var);
351 nvlbs = SCIPvarGetNVlbs(var);
352 nvubs = SCIPvarGetNVubs(var);
353 glb = SCIPvarGetLbGlobal(var);
354 gub = SCIPvarGetUbGlobal(var);
355
356 *result = FALSE;
357
358 /* Scan through lower bounds; for each binary vlbvar save the corresponding lb0 and lb1.
359 * Then check if there is an upper bound with this vlbvar and save ub0 and ub1.
360 * If the found bounds imply that the var value is fixed to some val0 when vlbvar = 0,
361 * save vlbvar and val0 to scvdata.
362 */
363 for( c = 0; c < nvlbs; ++c )
364 {
365 if( SCIPvarGetType(vlbvars[c]) != SCIP_VARTYPE_BINARY )
366 continue;
367
368 SCIPdebugMsg(scip, "var <%s>[%f, %f] lower bound: %f <%s> %+f", SCIPvarGetName(var), glb, gub, vlbcoefs[c], SCIPvarGetName(vlbvars[c]), vlbconstants[c]);
369
370 bvar = vlbvars[c];
371
372 lb0 = MAX(vlbconstants[c], glb);
373 lb1 = MAX(vlbconstants[c] + vlbcoefs[c], glb);
374
375 /* look for bvar in vubvars */
376 if( vubvars != NULL )
377 exists = SCIPsortedvecFindPtr((void**)vubvars, SCIPvarComp, bvar, nvubs, &pos);
378 else
379 exists = FALSE;
380 if( exists )
381 { /*lint --e{644}*/
382 SCIPdebugMsgPrint(scip, ", upper bound: %f <%s> %+f", vubcoefs[pos], SCIPvarGetName(vubvars[pos]), vubconstants[pos]); /*lint !e613*/
383
384 /* save the upper bounds */
385 ub0 = MIN(vubconstants[pos], gub);
386 ub1 = MIN(vubconstants[pos] + vubcoefs[pos], gub);
387 }
388 else
389 {
390 /* if there is no upper bound with vubvar = bvar, use global var bounds */
391 ub0 = gub;
392 ub1 = gub;
393 }
394
395 /* the 'off' domain of a semicontinuous var should reduce to a single point and be different from the 'on' domain */
396 SCIPdebugMsgPrint(scip, " -> <%s> in [%f, %f] (off), [%f, %f] (on)\n", SCIPvarGetName(var), lb0, ub0, lb1, ub1);
397 if( SCIPisEQ(scip, lb0, ub0) && (!SCIPisEQ(scip, lb0, lb1) || !SCIPisEQ(scip, ub0, ub1)) )
398 {
399 if( scvdata == NULL )
400 {
401 SCIP_CALL( SCIPallocClearBlockMemory(scip, &scvdata) );
402 }
403
404 SCIP_CALL( addSCVarIndicator(scip, scvdata, bvar, lb0, lb1, ub1) );
405 }
406 }
407
408 /* look for vubvars that have not been processed yet */
409 assert(vubvars != NULL || nvubs == 0);
410 for( c = 0; c < nvubs; ++c )
411 {
412 if( SCIPvarGetType(vubvars[c]) != SCIP_VARTYPE_BINARY ) /*lint !e613*/
413 continue;
414
415 bvar = vubvars[c]; /*lint !e613*/
416
417 /* skip vars that are in vlbvars */
418 if( vlbvars != NULL && SCIPsortedvecFindPtr((void**)vlbvars, SCIPvarComp, bvar, nvlbs, &pos) )
419 continue;
420
421 SCIPdebugMsg(scip, "var <%s>[%f, %f] upper bound: %f <%s> %+f",
422 SCIPvarGetName(var), glb, gub, vubcoefs[c], SCIPvarGetName(vubvars[c]), vubconstants[c]); /*lint !e613*/
423
424 lb0 = glb;
425 lb1 = glb;
426 ub0 = MIN(vubconstants[c], gub);
427 ub1 = MIN(vubconstants[c] + vubcoefs[c], gub);
428
429 /* the 'off' domain of a semicontinuous var should reduce to a single point and be different from the 'on' domain */
430 SCIPdebugMsgPrint(scip, " -> <%s> in [%f, %f] (off), [%f, %f] (on)\n", SCIPvarGetName(var), lb0, ub0, lb1, ub1);
431 if( SCIPisEQ(scip, lb0, ub0) && (!SCIPisEQ(scip, lb0, lb1) || !SCIPisEQ(scip, ub0, ub1)) )
432 {
433 if( scvdata == NULL )
434 {
435 SCIP_CALL( SCIPallocClearBlockMemory(scip, &scvdata) );
436 }
437
438 SCIP_CALL( addSCVarIndicator(scip, scvdata, bvar, lb0, lb1, ub1) );
439 }
440 }
441
442 if( scvdata != NULL )
443 {
444 #ifdef SCIP_DEBUG
445 SCIPdebugMsg(scip, "var <%s> has global bounds [%f, %f] and the following on/off bounds:\n", SCIPvarGetName(var), glb, gub);
446 for( c = 0; c < scvdata->nbnds; ++c )
447 {
448 SCIPdebugMsg(scip, " c = %d, bvar <%s>: val0 = %f\n", c, SCIPvarGetName(scvdata->bvars[c]), scvdata->vals0[c]);
449 }
450 #endif
451 SCIP_CALL( SCIPhashmapInsert(scvars, var, scvdata) );
452 *result = TRUE;
453 }
454
455 return SCIP_OKAY;
456 }
457
458 /*
459 * Semicontinuous expression methods
460 */
461
462 /* checks if an expression is semicontinuous
463 *
464 * An expression is semicontinuous if all of its nonlinear variables are semicontinuous
465 * and share at least one common indicator variable
466 */
467 static
468 SCIP_RETCODE exprIsSemicontinuous(
469 SCIP* scip, /**< SCIP data structure */
470 SCIP_NLHDLRDATA* nlhdlrdata, /**< nonlinear handler data */
471 SCIP_NLHDLREXPRDATA* nlhdlrexprdata, /**< nlhdlr expression data */
472 SCIP_EXPR* expr, /**< expression */
473 SCIP_Bool* res /**< buffer to store whether the expression is semicontinuous */
474 )
475 {
476 int v;
477 SCIP_Bool var_is_sc;
478 SCVARDATA* scvdata;
479 SCIP_VAR* var;
480 int nindicators;
481 int nbnds0;
482 int c;
483 SCIP_VAR** indicators;
484 SCIP_Bool* nonlinear;
485
486 *res = FALSE;
487
488 /* constant expression is not semicontinuous; variable expressions are of no interest here */
489 if( nlhdlrexprdata->nvars == 0 )
490 return SCIP_OKAY;
491
492 indicators = NULL;
493 nindicators = 0;
494 nbnds0 = 0;
495
496 if( SCIPisExprSum(scip, expr) )
497 {
498 SCIP_EXPRITER* it;
499 SCIP_EXPR* child;
500 SCIP_EXPR* curexpr;
501 int pos;
502 SCIP_Bool issc;
503
504 /* sums are treated separately because if there are variables that are non-semicontinuous but
505 * appear only linearly, we still want to apply perspective to expr
506 */
507
508 SCIP_CALL( SCIPallocClearBufferArray(scip, &nonlinear, nlhdlrexprdata->nvars) );
509 SCIP_CALL( SCIPcreateExpriter(scip, &it) );
510
511 for( c = 0; c < SCIPexprGetNChildren(expr); ++c )
512 {
513 child = SCIPexprGetChildren(expr)[c];
514
515 if( SCIPisExprVar(scip, child) )
516 {
517 var = SCIPgetVarExprVar(child);
518
519 /* save information on semicontinuity of child */
520 SCIP_CALL( varIsSemicontinuous(scip, var, nlhdlrdata->scvars, &var_is_sc) );
521
522 /* since child is a variable, go on regardless of the value of var_is_sc */
523 continue;
524 }
525
526 issc = TRUE;
527
528 SCIP_CALL( SCIPexpriterInit(it, child, SCIP_EXPRITER_DFS, FALSE) );
529 curexpr = SCIPexpriterGetCurrent(it);
530
531 /* all nonlinear terms of a sum should be semicontinuous in original variables */
532 while( !SCIPexpriterIsEnd(it) )
533 {
534 assert(curexpr != NULL);
535
536 if( SCIPisExprVar(scip, curexpr) )
537 {
538 var = SCIPgetVarExprVar(curexpr);
539
540 if( !SCIPvarIsRelaxationOnly(var) )
541 {
542 SCIP_CALL( varIsSemicontinuous(scip, var, nlhdlrdata->scvars, &var_is_sc) );
543
544 /* mark the variable as nonlinear */
545 (void) SCIPsortedvecFindPtr((void**) nlhdlrexprdata->vars, SCIPvarComp, (void*) var, nlhdlrexprdata->nvars,
546 &pos);
547 assert(0 <= pos && pos < nlhdlrexprdata->nvars);
548 nonlinear[pos] = TRUE;
549
550 if( !var_is_sc )
551 {
552 /* non-semicontinuous child which is (due to a previous check) not a var ->
553 * expr is non-semicontinuous
554 */
555 issc = FALSE;
556 break;
557 }
558 }
559 }
560 curexpr = SCIPexpriterGetNext(it);
561 }
562
563 if( !issc )
564 {
565 SCIPfreeExpriter(&it);
566 goto TERMINATE;
567 }
568 }
569 SCIPfreeExpriter(&it);
570 }
571 else
572 {
573 /* non-sum expression */
574 nonlinear = NULL;
575
576 /* all variables of a non-sum on/off expression should be semicontinuous */
577 for( v = 0; v < nlhdlrexprdata->nvars; ++v )
578 {
579 SCIP_CALL( varIsSemicontinuous(scip, nlhdlrexprdata->vars[v], nlhdlrdata->scvars, &var_is_sc) );
580 if( !var_is_sc )
581 return SCIP_OKAY;
582 }
583 }
584
585 /* look for common binary variables for all variables of the expression */
586
587 SCIPdebugMsg(scip, "Array intersection for var <%s>\n", SCIPvarGetName(nlhdlrexprdata->vars[0]));
588 for( v = 0; v < nlhdlrexprdata->nvars; ++v )
589 {
590 SCIPdebugMsg(scip, "%s; \n", SCIPvarGetName(nlhdlrexprdata->vars[v]));
591
592 if( nonlinear != NULL && !nonlinear[v] )
593 continue;
594
595 scvdata = (SCVARDATA*)SCIPhashmapGetImage(nlhdlrdata->scvars, (void*) nlhdlrexprdata->vars[v]);
596
597 /* we should have exited earlier if there is a nonlinear non-semicontinuous variable */
598 assert(scvdata != NULL);
599
600 if( indicators == NULL )
601 {
602 SCIP_CALL( SCIPduplicateBlockMemoryArray(scip, &indicators, scvdata->bvars, scvdata->nbnds) );
603 nbnds0 = scvdata->nbnds;
604 nindicators = nbnds0;
605 }
606 else
607 {
608 SCIPcomputeArraysIntersectionPtr((void**)indicators, nindicators, (void**)scvdata->bvars, scvdata->nbnds,
609 SCIPvarComp, (void**)indicators, &nindicators);
610 }
611
612 /* if we have found out that the intersection is empty, expr is not semicontinuous */
613 if( indicators != NULL && nindicators == 0 )
614 {
615 SCIPfreeBlockMemoryArray(scip, &indicators, nbnds0);
616 goto TERMINATE;
617 }
618 }
619
620 /* this can happen if all children are linear vars and none are semicontinuous */
621 if( indicators == NULL )
622 {
623 goto TERMINATE;
624 }
625 assert(nindicators > 0 && nindicators <= nbnds0);
626
627 if( nindicators < nbnds0 )
628 {
629 SCIP_CALL( SCIPreallocBlockMemoryArray(scip, &indicators, nbnds0, nindicators) );
630 }
631
632 for( v = 0; v < nindicators; ++v )
633 {
634 SCIP_CALL( SCIPcaptureVar(scip, indicators[v]) );
635 }
636 nlhdlrexprdata->indicators = indicators;
637 nlhdlrexprdata->nindicators = nindicators;
638 *res = TRUE;
639
640 TERMINATE:
641 SCIPfreeBufferArrayNull(scip, &nonlinear);
642
643 return SCIP_OKAY;
644 }
645
646 /** computes the 'off' value of the expression and the 'off' values of
647 * semicontinuous auxiliary variables for each indicator variable
648 */
649 static
650 SCIP_RETCODE computeOffValues(
651 SCIP* scip, /**< SCIP data structure */
652 SCIP_NLHDLRDATA* nlhdlrdata, /**< nonlinear handler data */
653 SCIP_NLHDLREXPRDATA* nlhdlrexprdata, /**< nlhdlr expression data */
654 SCIP_EXPR* expr /**< expression */
655 )
656 {
657 SCIP_EXPRITER* it;
658 SCIP_SOL* sol;
659 int i;
660 int v;
661 int norigvars;
662 SCIP_Real* origvals0;
663 SCIP_VAR** origvars;
664 SCVARDATA* scvdata;
665 SCIP_VAR* auxvar;
666 SCIP_EXPR* curexpr;
667 SCIP_HASHMAP* auxvarmap;
668 SCIP_Bool hasnonsc;
669 int pos;
670
671 assert(expr != NULL);
672
673 SCIP_CALL( SCIPallocBlockMemoryArray(scip, &(nlhdlrexprdata->exprvals0), nlhdlrexprdata->nindicators) );
674 SCIP_CALL( SCIPcreateSol(scip, &sol, NULL) );
675 SCIP_CALL( SCIPallocBufferArray(scip, &origvals0, nlhdlrexprdata->nvars) );
676 SCIP_CALL( SCIPhashmapCreate(&auxvarmap, SCIPblkmem(scip), 10) );
677 SCIP_CALL( SCIPcreateExpriter(scip, &it) );
678 SCIP_CALL( SCIPduplicateBufferArray(scip, &origvars, nlhdlrexprdata->vars, nlhdlrexprdata->nvars) );
679 norigvars = nlhdlrexprdata->nvars;
680
681 for( i = nlhdlrexprdata->nindicators - 1; i >= 0; --i )
682 {
683 hasnonsc = FALSE;
684
685 /* set sol to the off value of all expr vars for this indicator */
686 for( v = 0; v < norigvars; ++v )
687 {
688 /* set vals0[v] = 0 if var is non-sc with respect to indicators[i] - then it will not
689 * contribute to exprvals0[i] since any non-sc var must be linear
690 */
691 scvdata = getSCVarDataInd(nlhdlrdata->scvars, origvars[v], nlhdlrexprdata->indicators[i], &pos);
692 if( scvdata == NULL )
693 {
694 origvals0[v] = 0.0;
695 hasnonsc = TRUE;
696 }
697 else
698 {
699 origvals0[v] = scvdata->vals0[pos];
700 }
701 }
702 SCIP_CALL( SCIPsetSolVals(scip, sol, norigvars, origvars, origvals0) );
703 SCIP_CALL( SCIPevalExpr(scip, expr, sol, 0L) );
704
705 if( SCIPexprGetEvalValue(expr) == SCIP_INVALID ) /*lint !e777*/
706 {
707 SCIPdebugMsg(scip, "expression evaluation failed for %p, removing indicator %s\n",
708 (void*)expr, SCIPvarGetName(nlhdlrexprdata->indicators[i]));
709 /* TODO should we fix the indicator variable to 1? */
710 /* since the loop is backwards, this only modifies the already processed part of nlhdlrexprdata->indicators */
711 SCIP_CALL( removeIndicator(scip, nlhdlrexprdata, i) );
712 continue;
713 }
714
715 nlhdlrexprdata->exprvals0[i] = SCIPexprGetEvalValue(expr);
716
717 /* iterate through the expression and create scvdata for aux vars */
718 SCIP_CALL( SCIPexpriterInit(it, expr, SCIP_EXPRITER_DFS, FALSE) );
719 curexpr = SCIPexpriterGetCurrent(it);
720
721 while( !SCIPexpriterIsEnd(it) )
722 {
723 auxvar = SCIPgetExprAuxVarNonlinear(curexpr);
724
725 if( auxvar != NULL )
726 {
727 SCIP_Bool issc = TRUE;
728 #ifndef NDEBUG
729 SCIP_EXPR** childvarexprs;
730 int nchildvarexprs;
731 SCIP_VAR* var;
732 #endif
733
734 if( hasnonsc )
735 {
736 /* expr is a sum with non-semicontinuous linear terms. Therefore, curexpr might be
737 * non-semicontinuous. In that case the auxvar is also non-semicontinuous, so
738 * we will skip on/off bounds computation.
739 */
740 if( SCIPisExprVar(scip, curexpr) )
741 {
742 /* easy case: curexpr is a variable, can check semicontinuity immediately */
743 scvdata = getSCVarDataInd(nlhdlrdata->scvars, SCIPgetVarExprVar(curexpr),
744 nlhdlrexprdata->indicators[i], &pos);
745 issc = scvdata != NULL;
746 }
747 else if( SCIPisExprSum(scip, curexpr) && curexpr == expr )
748 {
749 /* if expr itself is a sum, this is an exception since a sum with nonlinear terms is
750 * allowed to have both semicontinuous and non-semicontinuous variables; we skip it here
751 * and then analyse it term by term
752 */
753 issc = FALSE;
754 }
755
756 #ifndef NDEBUG
757 if( !SCIPisExprVar(scip, curexpr) && (!SCIPisExprSum(scip, curexpr) || curexpr != expr) )
758 {
759 /* curexpr is a non-variable expression and does not fit the sum special case,
760 * so it belongs to the non-linear part of expr.
761 * Since the non-linear part of expr must be semicontinuous with respect to
762 * nlhdlrexprdata->indicators[i], curexpr must be semicontinuous
763 */
764
765 SCIP_CALL( SCIPallocBufferArray(scip, &childvarexprs, norigvars) );
766 SCIP_CALL( SCIPgetExprVarExprs(scip, curexpr, childvarexprs, &nchildvarexprs) );
767
768 /* all nonlinear variables of a sum on/off term should be semicontinuous */
769 for( v = 0; v < nchildvarexprs; ++v )
770 {
771 var = SCIPgetVarExprVar(childvarexprs[v]);
772 scvdata = getSCVarDataInd(nlhdlrdata->scvars, var, nlhdlrexprdata->indicators[i], &pos);
773 assert(scvdata != NULL);
774
775 SCIP_CALL( SCIPreleaseExpr(scip, &childvarexprs[v]) );
776 }
777
778 SCIPfreeBufferArray(scip, &childvarexprs);
779 }
780 #endif
781 }
782
783 if( issc )
784 {
785 /* we know that all vars are semicontinuous with respect to exprdata->indicators; it remains to:
786 * - get or create the scvardata structure for auxvar
787 * - if had to create scvardata, add it to scvars hashmap
788 * - add the indicator and the off value (= curexpr's off value) to scvardata
789 */
790 scvdata = (SCVARDATA*) SCIPhashmapGetImage(nlhdlrdata->scvars, (void*)auxvar);
791 if( scvdata == NULL )
792 {
793 SCIP_CALL( SCIPallocClearBlockMemory(scip, &scvdata) );
794 SCIP_CALL( SCIPallocBlockMemoryArray(scip, &scvdata->bvars, nlhdlrexprdata->nindicators) );
795 SCIP_CALL( SCIPallocBlockMemoryArray(scip, &scvdata->vals0, nlhdlrexprdata->nindicators) );
796 SCIP_CALL( SCIPallocBlockMemoryArray(scip, &scvdata->lbs1, nlhdlrexprdata->nindicators) );
797 SCIP_CALL( SCIPallocBlockMemoryArray(scip, &scvdata->ubs1, nlhdlrexprdata->nindicators) );
798 scvdata->bndssize = nlhdlrexprdata->nindicators;
799 SCIP_CALL( SCIPhashmapInsert(nlhdlrdata->scvars, auxvar, scvdata) );
800 }
801
802 SCIP_CALL( addSCVarIndicator(scip, scvdata, nlhdlrexprdata->indicators[i],
803 SCIPexprGetEvalValue(curexpr), SCIPvarGetLbGlobal(auxvar), SCIPvarGetUbGlobal(auxvar)) );
804 }
805
806 SCIP_CALL( addAuxVar(scip, nlhdlrexprdata, auxvarmap, auxvar) );
807 }
808
809 curexpr = SCIPexpriterGetNext(it);
810 }
811 }
812
813 SCIPfreeExpriter(&it);
814 SCIPhashmapFree(&auxvarmap);
815 SCIPfreeBufferArray(scip, &origvars);
816 SCIPfreeBufferArray(scip, &origvals0);
817 SCIP_CALL( SCIPfreeSol(scip, &sol) );
818
819 return SCIP_OKAY;
820 }
821
822 /*
823 * Probing and bound tightening methods
824 */
825
826 /** go into probing and set some variable bounds */
827 static
828 SCIP_RETCODE startProbing(
829 SCIP* scip, /**< SCIP data structure */
830 SCIP_NLHDLRDATA* nlhdlrdata, /**< nonlinear handler data */
831 SCIP_NLHDLREXPRDATA* nlhdlrexprdata, /**< nlhdlr expression data */
832 SCIP_VAR* indicator, /**< indicator variable */
833 SCIP_VAR** probingvars, /**< array of vars whose bounds we will change in probing */
834 SCIP_INTERVAL* probingdoms, /**< array of intervals to which bounds of probingvars will be changed in probing */
835 int nprobingvars, /**< number of probing vars */
836 SCIP_SOL* sol, /**< solution to be separated */
837 SCIP_SOL** solcopy, /**< buffer for a copy of sol before going into probing; if *solcopy == sol, then copy is created */
838 SCIP_Bool* cutoff_probing /**< pointer to store whether indicator == 1 is infeasible */
839 )
840 {
841 int v;
842 SCIP_Real newlb;
843 SCIP_Real newub;
844 SCIP_Bool propagate;
845
846 propagate = SCIPgetDepth(scip) == 0;
847
848 /* if a copy of sol has not been created yet, then create one now and copy the relevant var values from sol,
849 * because sol can change after SCIPstartProbing, e.g., when linked to the LP solution
850 */
851 if( *solcopy == sol )
852 {
853 SCIP_CALL( SCIPcreateSol(scip, solcopy, NULL) );
854 for( v = 0; v < nlhdlrexprdata->nvars; ++v )
855 {
856 SCIP_CALL( SCIPsetSolVal(scip, *solcopy, nlhdlrexprdata->vars[v], SCIPgetSolVal(scip, sol, nlhdlrexprdata->vars[v])) );
857 }
858 for( v = 0; v < nlhdlrexprdata->nindicators; ++v )
859 {
860 SCIP_CALL( SCIPsetSolVal(scip, *solcopy, nlhdlrexprdata->indicators[v], SCIPgetSolVal(scip, sol, nlhdlrexprdata->indicators[v])) );
861 }
862 }
863
864 /* go into probing */
865 SCIP_CALL( SCIPstartProbing(scip) );
866
867 /* create a probing node */
868 SCIP_CALL( SCIPnewProbingNode(scip) );
869
870 /* set indicator to 1 */
871 SCIP_CALL( SCIPchgVarLbProbing(scip, indicator, 1.0) );
872
873 /* apply stored bounds */
874 for( v = 0; v < nprobingvars; ++v )
875 {
876 newlb = SCIPintervalGetInf(probingdoms[v]);
877 newub = SCIPintervalGetSup(probingdoms[v]);
878
879 if( SCIPisGT(scip, newlb, SCIPvarGetLbLocal(probingvars[v])) || (newlb >= 0.0 && SCIPvarGetLbLocal(probingvars[v]) < 0.0) )
880 {
881 SCIP_CALL( SCIPchgVarLbProbing(scip, probingvars[v], newlb) );
882 }
883 if( SCIPisLT(scip, newub, SCIPvarGetUbLocal(probingvars[v])) || (newub <= 0.0 && SCIPvarGetUbLocal(probingvars[v]) > 0.0) )
884 {
885 SCIP_CALL( SCIPchgVarUbProbing(scip, probingvars[v], newub) );
886 }
887 }
888
889 if( propagate )
890 {
891 SCIP_Longint ndomreds;
892
893 SCIP_CALL( SCIPpropagateProbing(scip, nlhdlrdata->maxproprounds, cutoff_probing, &ndomreds) );
894 }
895
896 return SCIP_OKAY;
897 }
898
899 /** analyse on/off bounds on a variable
900 *
901 * analyses for
902 * 1. tightening bounds in probing for indicator = 1,
903 * 2. fixing indicator / detecting cutoff if one or both states are infeasible,
904 * 3. tightening local bounds if indicator is fixed.
905 *
906 * `probinglb` and `probingub` are only set if `doprobing` is TRUE.
907 * They are either set to bounds that should be used in probing or to `SCIP_INVALID` if bounds on
908 * `var` shouldn't be changed in probing.
909 */
910 static
911 SCIP_RETCODE analyseVarOnoffBounds(
912 SCIP* scip, /**< SCIP data structure */
913 SCIP_NLHDLRDATA* nlhdlrdata, /**< nonlinear handler data */
914 SCIP_VAR* var, /**< variable */
915 SCIP_VAR* indicator, /**< indicator variable */
916 SCIP_Bool indvalue, /**< indicator value for which the bounds are applied */
917 SCIP_Bool* infeas, /**< pointer to store whether infeasibility has been detected */
918 SCIP_Real* probinglb, /**< pointer to store the lower bound to be applied in probing */
919 SCIP_Real* probingub, /**< pointer to store the upper bound to be applied in probing */
920 SCIP_Bool doprobing, /**< whether we currently consider to go into probing */
921 SCIP_Bool* reduceddom /**< pointer to store whether any variables were fixed */
922 )
923 {
924 SCVARDATA* scvdata;
925 int pos;
926 SCIP_Real sclb;
927 SCIP_Real scub;
928 SCIP_Real loclb;
929 SCIP_Real locub;
930 SCIP_Bool bndchgsuccess;
931
932 assert(var != NULL);
933 assert(indicator != NULL);
934 assert(infeas != NULL);
935 assert(reduceddom != NULL);
936
937 /* shouldn't be called if indicator is fixed to !indvalue */
938 assert((indvalue && SCIPvarGetUbLocal(indicator) > 0.5) || (!indvalue && SCIPvarGetLbLocal(indicator) < 0.5));
939
940 *infeas = FALSE;
941 *reduceddom = FALSE;
942 scvdata = getSCVarDataInd(nlhdlrdata->scvars, var, indicator, &pos);
943 if( doprobing )
944 {
945 assert(probinglb != NULL);
946 assert(probingub != NULL);
947
948 *probinglb = SCIP_INVALID;
949 *probingub = SCIP_INVALID;
950 }
951
952 /* nothing to do for non-semicontinuous variables */
953 if( scvdata == NULL )
954 {
955 return SCIP_OKAY;
956 }
957
958 sclb = indvalue ? scvdata->lbs1[pos] : scvdata->vals0[pos];
959 scub = indvalue ? scvdata->ubs1[pos] : scvdata->vals0[pos];
960 loclb = SCIPvarGetLbLocal(var);
961 locub = SCIPvarGetUbLocal(var);
962
963 /* nothing to do for fixed variables */
964 if( SCIPisEQ(scip, loclb, locub) )
965 return SCIP_OKAY;
966
967 /* use a non-redundant lower bound */
968 if( SCIPisGT(scip, sclb, SCIPvarGetLbLocal(var)) || (sclb >= 0.0 && loclb < 0.0) )
969 {
970 /* first check for infeasibility */
971 if( SCIPisFeasGT(scip, sclb, SCIPvarGetUbLocal(var)) )
972 {
973 SCIP_CALL( SCIPfixVar(scip, indicator, indvalue ? 0.0 : 1.0, infeas, &bndchgsuccess) );
974 *reduceddom += bndchgsuccess;
975 if( *infeas )
976 {
977 return SCIP_OKAY;
978 }
979 }
980 else if( nlhdlrdata->tightenbounds &&
981 (SCIPvarGetUbLocal(indicator) <= 0.5 || SCIPvarGetLbLocal(indicator) >= 0.5) )
982 {
983 /* indicator is fixed; due to a previous check, here it can only be fixed to indvalue;
984 * therefore, sclb is valid for the current node
985 */
986
987 if( indvalue == 0 )
988 {
989 assert(sclb == scub); /*lint !e777*/
990 SCIP_CALL( SCIPfixVar(scip, var, sclb, infeas, &bndchgsuccess) );
991 }
992 else
993 {
994 SCIP_CALL( SCIPtightenVarLb(scip, var, sclb, FALSE, infeas, &bndchgsuccess) );
995 }
996 *reduceddom += bndchgsuccess;
997 if( *infeas )
998 {
999 return SCIP_OKAY;
1000 }
1001 }
1002 }
1003
1004 /* use a non-redundant upper bound */
1005 if( SCIPisLT(scip, scub, SCIPvarGetUbLocal(var)) || (scub <= 0.0 && locub > 0.0) )
1006 {
1007 /* first check for infeasibility */
1008 if( SCIPisFeasLT(scip, scub, SCIPvarGetLbLocal(var)) )
1009 {
1010 SCIP_CALL( SCIPfixVar(scip, indicator, indvalue ? 0.0 : 1.0, infeas, &bndchgsuccess) );
1011 *reduceddom += bndchgsuccess;
1012 if( *infeas )
1013 {
1014 return SCIP_OKAY;
1015 }
1016 }
1017 else if( nlhdlrdata->tightenbounds &&
1018 (SCIPvarGetUbLocal(indicator) <= 0.5 || SCIPvarGetLbLocal(indicator) >= 0.5) )
1019 {
1020 /* indicator is fixed; due to a previous check, here it can only be fixed to indvalue;
1021 * therefore, scub is valid for the current node
1022 */
1023
1024 if( indvalue == 0 )
1025 {
1026 assert(sclb == scub); /*lint !e777*/
1027 SCIP_CALL( SCIPfixVar(scip, var, sclb, infeas, &bndchgsuccess) );
1028 }
1029 else
1030 {
1031 SCIP_CALL( SCIPtightenVarUb(scip, var, scub, FALSE, infeas, &bndchgsuccess) );
1032 }
1033 *reduceddom += bndchgsuccess;
1034 if( *infeas )
1035 {
1036 return SCIP_OKAY;
1037 }
1038 }
1039 }
1040
1041 /* If a bound change has been found and indvalue == TRUE, try to use the new bounds.
1042 * This is only done for indvalue == TRUE since this is where enfo asks other nlhdlrs to estimate,
1043 * and at indicator == FALSE we already only have a single point
1044 */
1045 if( doprobing && indvalue && (((scub - sclb) / (locub - loclb)) <= 1.0 - nlhdlrdata->mindomreduction ||
1046 (sclb >= 0.0 && loclb < 0.0) || (scub <= 0.0 && locub > 0.0)) )
1047 {
1048 *probinglb = sclb;
1049 *probingub = scub;
1050 }
1051
1052 SCIPdebugMsg(scip, "%s in [%g, %g] instead of [%g, %g] (vals0 = %g)\n", SCIPvarGetName(var), sclb, scub,
1053 SCIPvarGetLbLocal(var), SCIPvarGetUbLocal(var), scvdata->vals0[pos]);
1054
1055 return SCIP_OKAY;
1056 }
1057
1058 /** looks for bound tightenings to be applied either in the current node or in probing
1059 *
1060 * Loops through both possible values of indicator and calls analyseVarOnoffBounds().
1061 * Might update the `*doprobing` flag by setting it to `FALSE` if:
1062 * - indicator is fixed or
1063 * - analyseVarOnoffBounds() hasn't found a sufficient improvement at indicator==1.
1064 *
1065 * If `*doprobing==TRUE`, stores bounds suggested by analyseVarOnoffBounds() in order to apply them in probing together
1066 * with the fixing `indicator=1`.
1067 */
1068 static
1069 SCIP_RETCODE analyseOnoffBounds(
1070 SCIP* scip, /**< SCIP data structure */
1071 SCIP_NLHDLRDATA* nlhdlrdata, /**< nonlinear handler data */
1072 SCIP_NLHDLREXPRDATA* nlhdlrexprdata, /**< nlhdlr expression data */
1073 SCIP_VAR* indicator, /**< indicator variable */
1074 SCIP_VAR*** probingvars, /**< array to store variables whose bounds will be changed in probing */
1075 SCIP_INTERVAL** probingdoms, /**< array to store bounds to be applied in probing */
1076 int* nprobingvars, /**< pointer to store number of vars whose bounds will be changed in probing */
1077 SCIP_Bool* doprobing, /**< pointer to the flag telling whether we want to do probing */
1078 SCIP_RESULT* result /**< pointer to store the result */
1079 )
1080 {
1081 int v;
1082 SCIP_VAR* var;
1083 SCIP_Bool infeas;
1084 int b;
1085 SCIP_Real probinglb = SCIP_INVALID;
1086 SCIP_Real probingub = SCIP_INVALID;
1087 SCIP_Bool changed;
1088 SCIP_Bool reduceddom;
1089
1090 assert(indicator != NULL);
1091 assert(nprobingvars != NULL);
1092 assert(doprobing != NULL);
1093 assert(result != NULL);
1094
1095 changed = FALSE;
1096
1097 /* no probing if indicator already fixed */
1098 if( SCIPvarGetUbLocal(indicator) <= 0.5 || SCIPvarGetLbLocal(indicator) >= 0.5 )
1099 {
1100 *doprobing = FALSE;
1101 }
1102
1103 /* consider each possible value of indicator */
1104 for( b = 0; b < 2; ++b )
1105 {
1106 for( v = 0; v < nlhdlrexprdata->nvars; ++v )
1107 {
1108 /* nothing left to do if indicator is already fixed to !indvalue
1109 * (checked in the inner loop since analyseVarOnoff bounds might fix the indicator)
1110 */
1111 if( (b == 1 && SCIPvarGetUbLocal(indicator) <= 0.5) || (b == 0 && SCIPvarGetLbLocal(indicator) >= 0.5) )
1112 {
1113 *doprobing = FALSE;
1114 break;
1115 }
1116
1117 var = nlhdlrexprdata->vars[v];
1118
1119 SCIP_CALL( analyseVarOnoffBounds(scip, nlhdlrdata, var, indicator, b == 1, &infeas, &probinglb,
1120 &probingub, *doprobing, &reduceddom) );
1121
1122 if( infeas )
1123 {
1124 *result = SCIP_CUTOFF;
1125 *doprobing = FALSE;
1126 return SCIP_OKAY;
1127 }
1128 else if( reduceddom )
1129 {
1130 *result = SCIP_REDUCEDDOM;
1131 }
1132
1133 if( !(*doprobing) )
1134 continue;
1135
1136 /* if bounds to be applied in probing have been found, store them */
1137 if( probinglb != SCIP_INVALID ) /*lint !e777*/
1138 {
1139 assert(probingub != SCIP_INVALID); /*lint !e777*/
1140
1141 SCIP_CALL( SCIPreallocBufferArray(scip, probingvars, *nprobingvars + 1) );
1142 SCIP_CALL( SCIPreallocBufferArray(scip, probingdoms, *nprobingvars + 1) );
1143 (*probingvars)[*nprobingvars] = var;
1144 (*probingdoms)[*nprobingvars].inf = probinglb;
1145 (*probingdoms)[*nprobingvars].sup = probingub;
1146 ++*nprobingvars;
1147
1148 changed = TRUE;
1149 }
1150 }
1151 }
1152
1153 if( !changed )
1154 {
1155 *doprobing = FALSE;
1156 }
1157
1158 return SCIP_OKAY;
1159 }
1160
1161 /** saves local bounds on all expression variables, including auxiliary variables, obtained from propagating
1162 * indicator == 1 to the corresponding SCVARDATA (should only be used in the root node)
1163 */
1164 static
1165 SCIP_RETCODE tightenOnBounds(
1166 SCIP_NLHDLREXPRDATA* nlhdlrexprdata, /**< nlhdlr expression data */
1167 SCIP_HASHMAP* scvars, /**< hashmap with semicontinuous variables */
1168 SCIP_VAR* indicator /**< indicator variable */
1169 )
1170 {
1171 int v;
1172 SCIP_VAR* var;
1173 SCVARDATA* scvdata;
1174 int pos;
1175 SCIP_Real lb;
1176 SCIP_Real ub;
1177
1178 for( v = 0; v < nlhdlrexprdata->nvars; ++v )
1179 {
1180 var = nlhdlrexprdata->vars[v];
1181 lb = SCIPvarGetLbLocal(var);
1182 ub = SCIPvarGetUbLocal(var);
1183 scvdata = getSCVarDataInd(scvars, var, indicator, &pos);
1184
1185 if( scvdata != NULL )
1186 {
1187 scvdata->lbs1[pos] = MAX(scvdata->lbs1[pos], lb);
1188 scvdata->ubs1[pos] = MIN(scvdata->ubs1[pos], ub);
1189 }
1190 }
1191
1192 return SCIP_OKAY;
1193 }
1194
1195 /*
1196 * Callback methods of nonlinear handler
1197 */
1198
1199 /** nonlinear handler copy callback */
1200 static
1201 SCIP_DECL_NLHDLRCOPYHDLR(nlhdlrCopyhdlrPerspective)
1202 { /*lint --e{715}*/
1203 assert(targetscip != NULL);
1204 assert(sourcenlhdlr != NULL);
1205 assert(strcmp(SCIPnlhdlrGetName(sourcenlhdlr), NLHDLR_NAME) == 0);
1206
1207 SCIP_CALL( SCIPincludeNlhdlrPerspective(targetscip) );
1208
1209 return SCIP_OKAY;
1210 }
1211
1212
1213 /** callback to free data of handler */
1214 static
1215 SCIP_DECL_NLHDLRFREEHDLRDATA(nlhdlrFreehdlrdataPerspective)
1216 { /*lint --e{715}*/
1217 SCIPfreeBlockMemory(scip, nlhdlrdata);
1218
1219 return SCIP_OKAY;
1220 }
1221
1222
1223 /** callback to free expression specific data */
1224 static
1225 SCIP_DECL_NLHDLRFREEEXPRDATA(nlhdlrFreeExprDataPerspective)
1226 { /*lint --e{715}*/
1227 SCIP_CALL( freeNlhdlrExprData(scip, *nlhdlrexprdata) );
1228 SCIPfreeBlockMemory(scip, nlhdlrexprdata);
1229
1230 return SCIP_OKAY;
1231 }
1232
1233 /** callback to be called in initialization */
1234 #if 0
1235 static
1236 SCIP_DECL_NLHDLRINIT(nlhdlrInitPerspective)
1237 { /*lint --e{715}*/
1238 return SCIP_OKAY;
1239 }
1240 #endif
1241
1242 /** callback to be called in deinitialization */
1243 static
1244 SCIP_DECL_NLHDLREXIT(nlhdlrExitPerspective)
1245 { /*lint --e{715}*/
1246 SCIP_HASHMAPENTRY* entry;
1247 SCVARDATA* data;
1248 int c;
1249 SCIP_NLHDLRDATA* nlhdlrdata;
1250
1251 nlhdlrdata = SCIPnlhdlrGetData(nlhdlr);
1252 assert(nlhdlrdata != NULL);
1253
1254 if( nlhdlrdata->scvars != NULL )
1255 {
1256 for( c = 0; c < SCIPhashmapGetNEntries(nlhdlrdata->scvars); ++c )
1257 {
1258 entry = SCIPhashmapGetEntry(nlhdlrdata->scvars, c);
1259 if( entry != NULL )
1260 {
1261 data = (SCVARDATA*) SCIPhashmapEntryGetImage(entry);
1262 SCIPfreeBlockMemoryArray(scip, &data->ubs1, data->bndssize);
1263 SCIPfreeBlockMemoryArray(scip, &data->lbs1, data->bndssize);
1264 SCIPfreeBlockMemoryArray(scip, &data->vals0, data->bndssize);
1265 SCIPfreeBlockMemoryArray(scip, &data->bvars, data->bndssize);
1266 SCIPfreeBlockMemory(scip, &data);
1267 }
1268 }
1269 SCIPhashmapFree(&nlhdlrdata->scvars);
1270 assert(nlhdlrdata->scvars == NULL);
1271 }
1272
1273 return SCIP_OKAY;
1274 }
1275
1276 /** callback to detect structure in expression tree
1277 *
1278 * We are looking for expressions g(x), where x is a vector of semicontinuous variables that all share at least one
1279 * indicator variable.
1280 */
1281 static
1282 SCIP_DECL_NLHDLRDETECT(nlhdlrDetectPerspective)
1283 { /*lint --e{715}*/
1284 SCIP_NLHDLRDATA* nlhdlrdata;
1285 SCIP_EXPR** varexprs;
1286 SCIP_Bool success = FALSE;
1287 int i;
1288 SCIP_Bool hassepabelow = FALSE;
1289 SCIP_Bool hassepaabove = FALSE;
1290 SCIP_Bool hasnondefault = FALSE;
1291
1292 nlhdlrdata = SCIPnlhdlrGetData(nlhdlr);
1293
1294 assert(scip != NULL);
1295 assert(nlhdlr != NULL);
1296 assert(expr != NULL);
1297 assert(participating != NULL);
1298 assert(enforcing != NULL);
1299 assert(nlhdlrexprdata != NULL);
1300 assert(nlhdlrdata != NULL);
1301
1302 /* do not run if we will have no auxvar to add a cut for */
1303 if( SCIPgetExprNAuxvarUsesNonlinear(expr) == 0 )
1304 return SCIP_OKAY;
1305
1306 if( SCIPgetNBinVars(scip) == 0 )
1307 {
1308 SCIPdebugMsg(scip, "problem has no binary variables, not running perspective detection\n");
1309 return SCIP_OKAY;
1310 }
1311
1312 for( i = 0; i < SCIPgetExprNEnfosNonlinear(expr); ++i )
1313 {
1314 SCIP_NLHDLR* nlhdlr2;
1315 SCIP_NLHDLR_METHOD nlhdlr2participates;
1316 SCIP_Bool sepabelowusesactivity;
1317 SCIP_Bool sepaaboveusesactivity;
1318 SCIPgetExprEnfoDataNonlinear(expr, i, &nlhdlr2, NULL, &nlhdlr2participates, &sepabelowusesactivity, &sepaaboveusesactivity, NULL);
1319
1320 if( (nlhdlr2participates & SCIP_NLHDLR_METHOD_SEPABOTH) == 0 )
1321 continue;
1322
1323 if( !SCIPnlhdlrHasEstimate(nlhdlr2) )
1324 continue;
1325
1326 if( strcmp(SCIPnlhdlrGetName(nlhdlr2), "default") != 0 )
1327 hasnondefault = TRUE;
1328
1329 /* If we are supposed to run only on convex expressions, than check whether there is a nlhdlr
1330 * that participates in separation without using activity for it. Otherwise, check for
1331 * participation regardless of activity usage.
1332 */
1333 if( (nlhdlr2participates & SCIP_NLHDLR_METHOD_SEPABELOW) && (!nlhdlrdata->convexonly || !sepabelowusesactivity) )
1334 hassepabelow = TRUE;
1335
1336 if( (nlhdlr2participates & SCIP_NLHDLR_METHOD_SEPAABOVE) && (!nlhdlrdata->convexonly || !sepaaboveusesactivity) )
1337 hassepaabove = TRUE;
1338 }
1339
1340 /* If a sum expression is handled only by default nlhdlr, then all the children will have auxiliary vars.
1341 * Since the sum will then be linear in auxiliary variables, perspective can't improve anything for it
1342 */
1343 if( SCIPisExprSum(scip, expr) && !hasnondefault )
1344 {
1345 SCIPdebugMsg(scip, "sum expr only has default exprhdlr, not running perspective detection\n");
1346 return SCIP_OKAY;
1347 }
1348
1349 /* If no other nlhdlr separates, neither does perspective (if convexonly, only separation
1350 * without using activity counts)
1351 */
1352 if( !hassepabelow && !hassepaabove )
1353 {
1354 SCIPdebugMsg(scip, "no nlhdlr separates without using activity, not running perspective detection\n");
1355 return SCIP_OKAY;
1356 }
1357
1358 #ifdef SCIP_DEBUG
1359 SCIPdebugMsg(scip, "Called perspective detect, expr = %p: ", (void*)expr);
1360 SCIPprintExpr(scip, expr, NULL);
1361 SCIPdebugMsgPrint(scip, "\n");
1362 #endif
1363
1364 /* allocate memory */
1365 SCIP_CALL( SCIPallocClearBlockMemory(scip, nlhdlrexprdata) );
1366 if( nlhdlrdata->scvars == NULL )
1367 {
1368 SCIP_CALL( SCIPhashmapCreate(&(nlhdlrdata->scvars), SCIPblkmem(scip), SCIPgetNVars(scip)) );
1369 }
1370
1371 /* save varexprs to nlhdlrexprdata */
1372 SCIP_CALL( SCIPgetExprNVars(scip, expr, &(*nlhdlrexprdata)->nvars) );
1373 SCIP_CALL( SCIPallocBlockMemoryArray(scip, &(*nlhdlrexprdata)->vars, (*nlhdlrexprdata)->nvars) );
1374 SCIP_CALL( SCIPallocBufferArray(scip, &varexprs, (*nlhdlrexprdata)->nvars) );
1375 (*nlhdlrexprdata)->varssize = (*nlhdlrexprdata)->nvars;
1376 SCIP_CALL( SCIPgetExprVarExprs(scip, expr, varexprs, &(*nlhdlrexprdata)->nvars) );
1377 for( i = 0; i < (*nlhdlrexprdata)->nvars; ++i )
1378 {
1379 (*nlhdlrexprdata)->vars[i] = SCIPgetVarExprVar(varexprs[i]);
1380 SCIP_CALL( SCIPreleaseExpr(scip, &varexprs[i]) );
1381 SCIP_CALL( SCIPcaptureVar(scip, (*nlhdlrexprdata)->vars[i]) );
1382 }
1383 SCIPsortPtr((void**) (*nlhdlrexprdata)->vars, SCIPvarComp, (*nlhdlrexprdata)->nvars);
1384 SCIPfreeBufferArray(scip, &varexprs);
1385
1386 /* check if expr is semicontinuous and save indicator variables */
1387 SCIP_CALL( exprIsSemicontinuous(scip, nlhdlrdata, *nlhdlrexprdata, expr, &success) );
1388
1389 if( success )
1390 {
1391 assert(*nlhdlrexprdata != NULL);
1392 assert((*nlhdlrexprdata)->nindicators > 0);
1393
1394 if( hassepaabove )
1395 *participating |= SCIP_NLHDLR_METHOD_SEPAABOVE;
1396 if( hassepabelow )
1397 *participating |= SCIP_NLHDLR_METHOD_SEPABELOW;
1398
1399 #ifdef SCIP_DEBUG
1400 SCIPinfoMessage(scip, NULL, "detected an on/off expr: ");
1401 SCIPprintExpr(scip, expr, NULL);
1402 SCIPinfoMessage(scip, NULL, "\n");
1403 #endif
1404 }
(2) Event check_after_deref: |
Null-checking "*nlhdlrexprdata" suggests that it may be null, but it has already been dereferenced on all paths leading to the check. |
Also see events: |
[deref_ptr_in_call] |
1405 else if( *nlhdlrexprdata != NULL )
1406 {
1407 SCIP_CALL( nlhdlrFreeExprDataPerspective(scip, nlhdlr, expr, nlhdlrexprdata) );
1408 }
1409
1410 return SCIP_OKAY;
1411 }
1412
1413
1414 /** auxiliary evaluation callback of nonlinear handler */
1415 static
1416 SCIP_DECL_NLHDLREVALAUX(nlhdlrEvalauxPerspective)
1417 { /*lint --e{715}*/
1418 int e;
1419 SCIP_Real maxdiff;
1420 SCIP_Real auxvarvalue;
1421 SCIP_Real enfoauxval;
1422
1423 assert(scip != NULL);
1424 assert(expr != NULL);
1425 assert(auxvalue != NULL);
1426
1427 auxvarvalue = SCIPgetSolVal(scip, sol, SCIPgetExprAuxVarNonlinear(expr));
1428 maxdiff = 0.0;
1429 *auxvalue = auxvarvalue;
1430
1431 /* use the auxvalue from one of the other nlhdlrs that estimates for this expr: take the one that is farthest
1432 * from the current value of auxvar
1433 */
1434 for( e = 0; e < SCIPgetExprNEnfosNonlinear(expr); ++e )
1435 {
1436 SCIP_NLHDLR* nlhdlr2;
1437 SCIP_NLHDLREXPRDATA* nlhdlr2exprdata;
1438 SCIP_NLHDLR_METHOD nlhdlr2participation;
1439
1440 SCIPgetExprEnfoDataNonlinear(expr, e, &nlhdlr2, &nlhdlr2exprdata, &nlhdlr2participation, NULL, NULL, NULL);
1441
1442 /* skip nlhdlr that do not participate or do not provide estimate */
1443 if( (nlhdlr2participation & SCIP_NLHDLR_METHOD_SEPABOTH) == 0 || !SCIPnlhdlrHasEstimate(nlhdlr2) )
1444 continue;
1445
1446 SCIP_CALL( SCIPnlhdlrEvalaux(scip, nlhdlr2, expr, nlhdlr2exprdata, &enfoauxval, sol) );
1447
1448 SCIPsetExprEnfoAuxValueNonlinear(expr, e, enfoauxval);
1449
1450 if( REALABS(enfoauxval - auxvarvalue) > maxdiff && enfoauxval != SCIP_INVALID ) /*lint !e777*/
1451 {
1452 maxdiff = REALABS(enfoauxval - auxvarvalue);
1453 *auxvalue = enfoauxval;
1454 }
1455 }
1456
1457 return SCIP_OKAY;
1458 }
1459
1460 /** separation initialization method of a nonlinear handler */
1461 static
1462 SCIP_DECL_NLHDLRINITSEPA(nlhdlrInitSepaPerspective)
1463 { /*lint --e{715}*/
1464 int sindicators;
1465
1466 sindicators = nlhdlrexprdata->nindicators;
1467
1468 /* compute 'off' values of expr and subexprs (and thus auxvars too) */
1469 SCIP_CALL( computeOffValues(scip, SCIPnlhdlrGetData(nlhdlr), nlhdlrexprdata, expr) );
1470
1471 /* some indicator variables might have been removed if evaluation failed, check how many remain */
1472 if( nlhdlrexprdata->nindicators == 0 )
1473 {
1474 SCIPfreeBlockMemoryArray(scip, &nlhdlrexprdata->indicators, sindicators);
1475 SCIPfreeBlockMemoryArray(scip, &nlhdlrexprdata->exprvals0, sindicators);
1476 }
1477 else if( nlhdlrexprdata->nindicators < sindicators )
1478 {
1479 SCIP_CALL( SCIPreallocBlockMemoryArray(scip, &nlhdlrexprdata->indicators, sindicators, nlhdlrexprdata->nindicators) );
1480 SCIP_CALL( SCIPreallocBlockMemoryArray(scip, &nlhdlrexprdata->exprvals0, sindicators, nlhdlrexprdata->nindicators) );
1481 }
1482
1483 return SCIP_OKAY;
1484 }
1485
1486
1487 #if 0
1488 /** separation deinitialization method of a nonlinear handler (called during CONSEXITSOL) */
1489 static
1490 SCIP_DECL_NLHDLREXITSEPA(nlhdlrExitSepaPerspective)
1491 { /*lint --e{715}*/
1492 SCIPerrorMessage("method of perspective nonlinear handler not implemented yet\n");
1493 SCIPABORT(); /*lint --e{527}*/
1494
1495 return SCIP_OKAY;
1496 }
1497 #endif
1498
1499 /** nonlinear handler enforcement callback
1500 *
1501 * "Perspectivies" cuts produced by other nonlinear handlers.
1502 *
1503 * Suppose that we want to separate \f$x\f$ from the set \f$\{ x : g(x) \leq 0\}\f$.
1504 * If \f$g(x) = g^0\f$ if indicator \f$z = 0\f$, and a cut is given by \f$\sum_i a_ix_i + c \leq \text{aux}\f$, where \f$x_i = x_i^0\f$ if \f$z = 0\f$ for all \f$i\f$,
1505 * then the "perspectivied" cut is \f[\sum_i a_ix_i + c + (1 - z)\,(g^0 - c - \sum_i a_ix_i^0) \leq \text{aux}.\f]
1506 * This ensures that at \f$z = 1\f$, the new cut is equivalent to the given cut, and at \f$z = 0\f$ it reduces to \f$g^0 \leq \text{aux}\f$.
1507 */
1508 static
1509 SCIP_DECL_NLHDLRENFO(nlhdlrEnfoPerspective)
1510 { /*lint --e{715}*/
1511 SCIP_ROWPREP* rowprep;
1512 SCIP_VAR* auxvar;
1513 int i;
1514 int j;
1515 SCIP_NLHDLRDATA* nlhdlrdata;
1516 SCIP_Real cst0;
1517 SCIP_VAR* indicator;
1518 SCIP_PTRARRAY* rowpreps2;
1519 SCIP_PTRARRAY* rowpreps;
1520 int nrowpreps;
1521 SCIP_SOL* solcopy;
1522 SCIP_Bool doprobing;
1523 SCIP_BOOLARRAY* addedbranchscores2;
1524 SCIP_Bool stop;
1525 int nenfos;
1526 int* enfoposs;
1527 SCIP_SOL* soladj;
1528 int pos;
1529 SCVARDATA* scvdata;
1530
1531 nlhdlrdata = SCIPnlhdlrGetData(nlhdlr);
1532
1533 #ifdef SCIP_DEBUG
1534 SCIPinfoMessage(scip, NULL, "enforcement method of perspective nonlinear handler called for expr %p: ", (void*)expr);
1535 SCIP_CALL( SCIPprintExpr(scip, expr, NULL) );
1536 SCIPinfoMessage(scip, NULL, " at\n");
1537 for( i = 0; i < nlhdlrexprdata->nvars; ++i )
1538 {
1539 SCIPinfoMessage(scip, NULL, "%s = %g\n", SCIPvarGetName(nlhdlrexprdata->vars[i]),
1540 SCIPgetSolVal(scip, sol, nlhdlrexprdata->vars[i]));
1541 }
1542 SCIPinfoMessage(scip, NULL, "%s = %g", SCIPvarGetName(SCIPgetExprAuxVarNonlinear(expr)),
1543 SCIPgetSolVal(scip, sol, SCIPgetExprAuxVarNonlinear(expr)));
1544 #endif
1545
1546 assert(scip != NULL);
1547 assert(expr != NULL);
1548 assert(conshdlr != NULL);
1549 assert(nlhdlrexprdata != NULL);
1550 assert(nlhdlrdata != NULL);
1551
1552 if( nlhdlrexprdata->nindicators == 0 )
1553 {
1554 /* we might have removed all indicators in initsepa */
1555 *result = SCIP_DIDNOTRUN;
1556 return SCIP_OKAY;
1557 }
1558
1559 auxvar = SCIPgetExprAuxVarNonlinear(expr);
1560 assert(auxvar != NULL);
1561
1562 /* detect should have picked only those expressions for which at least one other nlhdlr can enforce */
1563 assert(SCIPgetExprNEnfosNonlinear(expr) > 1);
1564
1565 SCIP_CALL( SCIPallocBufferArray(scip, &enfoposs, SCIPgetExprNEnfosNonlinear(expr) - 1) );
1566
1567 doprobing = FALSE;
1568 nenfos = 0;
1569 soladj = NULL;
1570
1571 /* find suitable nlhdlrs and check if there is enough violation to do probing */
1572 for( j = 0; j < SCIPgetExprNEnfosNonlinear(expr); ++j )
1573 {
1574 SCIP_NLHDLR* nlhdlr2;
1575 SCIP_NLHDLREXPRDATA* nlhdlr2exprdata;
1576 SCIP_NLHDLR_METHOD nlhdlr2participate;
1577 SCIP_Real nlhdlr2auxvalue;
1578 SCIP_Real violation;
1579 SCIP_Bool violbelow;
1580 SCIP_Bool violabove;
1581 SCIP_Bool sepausesactivity = FALSE;
1582
1583 SCIPgetExprEnfoDataNonlinear(expr, j, &nlhdlr2, &nlhdlr2exprdata, &nlhdlr2participate, !overestimate ? &sepausesactivity : NULL, overestimate ? &sepausesactivity: NULL, &nlhdlr2auxvalue); /*lint !e826*/
1584
1585 if( nlhdlr2 == nlhdlr )
1586 continue;
1587
1588 /* if nlhdlr2 cannot estimate, then cannot use it */
1589 if( !SCIPnlhdlrHasEstimate(nlhdlr2) )
1590 continue;
1591
1592 /* if nlhdlr2 does not participate in the separation on the desired side (overestimate), then skip it */
1593 if( (nlhdlr2participate & (overestimate ? SCIP_NLHDLR_METHOD_SEPAABOVE : SCIP_NLHDLR_METHOD_SEPABELOW)) == 0 )
1594 continue;
1595
1596 /* if only working on convex-looking expressions, then skip nlhdlr if it uses activity for estimates */
1597 if( nlhdlrdata->convexonly && sepausesactivity )
1598 continue;
1599
1600 /* evalaux should have called evalaux of nlhdlr2 by now
1601 * check whether handling the violation for nlhdlr2 requires under- or overestimation and this fits to
1602 * overestimate flag
1603 */
1604 SCIP_CALL( SCIPgetExprAbsAuxViolationNonlinear(scip, expr, nlhdlr2auxvalue, sol, &violation, &violbelow,
1605 &violabove) );
1606 assert(violation >= 0.0);
1607
1608 if( (overestimate && !violabove) || (!overestimate && !violbelow) )
1609 continue;
1610
1611 /* if violation is small, cuts would likely be weak - skip perspectification */
1612 if( !allowweakcuts && violation < SCIPfeastol(scip) )
1613 continue;
1614
1615 enfoposs[nenfos] = j;
1616 ++nenfos;
1617
1618 /* enable probing if tightening the domain could be useful for nlhdlr and violation is above threshold */
1619 if( sepausesactivity && violation >= nlhdlrdata->minviolprobing )
1620 doprobing = TRUE;
1621 }
1622
1623 if( nenfos == 0 )
1624 {
1625 *result = SCIP_DIDNOTRUN;
1626 SCIPfreeBufferArray(scip, &enfoposs);
1627 return SCIP_OKAY;
1628 }
1629
1630 /* check probing frequency against depth in b&b tree */
1631 if( nlhdlrdata->probingfreq == -1 || (nlhdlrdata->probingfreq == 0 && SCIPgetDepth(scip) != 0) ||
1632 (nlhdlrdata->probingfreq > 0 && SCIPgetDepth(scip) % nlhdlrdata->probingfreq != 0) )
1633 doprobing = FALSE;
1634
1635 /* if addbranchscores is TRUE, then we can assume to be in enforcement and not in separation */
1636 if( nlhdlrdata->probingonlyinsepa && addbranchscores )
1637 doprobing = FALSE;
1638
1639 /* disable probing if already being in probing or if in a subscip */
1640 if( SCIPinProbing(scip) || SCIPgetSubscipDepth(scip) != 0 )
1641 doprobing = FALSE;
1642
1643 nrowpreps = 0;
1644 *result = SCIP_DIDNOTFIND;
1645 solcopy = sol;
1646 stop = FALSE;
1647
1648 SCIP_CALL( SCIPcreatePtrarray(scip, &rowpreps2) );
1649 SCIP_CALL( SCIPcreatePtrarray(scip, &rowpreps) );
1650 SCIP_CALL( SCIPcreateBoolarray(scip, &addedbranchscores2) );
1651
1652 /* build cuts for every indicator variable */
1653 for( i = 0; i < nlhdlrexprdata->nindicators && !stop; ++i )
1654 {
1655 int v;
1656 int minidx;
1657 int maxidx;
1658 int r;
1659 SCIP_VAR** probingvars;
1660 SCIP_INTERVAL* probingdoms;
1661 int nprobingvars;
1662 SCIP_Bool doprobingind;
1663 SCIP_Real indval;
1664 SCIP_Real solval;
1665 SCIP_Bool adjrefpoint;
1666
1667 indicator = nlhdlrexprdata->indicators[i];
1668 probingvars = NULL;
1669 probingdoms = NULL;
1670 nprobingvars = 0;
1671 doprobingind = doprobing;
1672 solval = SCIPgetSolVal(scip, solcopy, indicator);
1673 adjrefpoint = nlhdlrdata->adjrefpoint && !SCIPisFeasEQ(scip, solval, 1.0);
1674
1675 SCIP_CALL( analyseOnoffBounds(scip, nlhdlrdata, nlhdlrexprdata, indicator, &probingvars, &probingdoms,
1676 &nprobingvars, &doprobingind, result) );
1677
1678 /* don't add perspective cuts for fixed indicators since there is no use for perspectivy */
1679 if( SCIPvarGetLbLocal(indicator) >= 0.5 )
1680 {
1681 assert(!doprobingind);
1682 continue;
1683 }
1684 if( SCIPvarGetUbLocal(indicator) <= 0.5 )
1685 { /* this case is stronger as it implies that everything is fixed;
1686 * therefore we are now happy
1687 */
1688 assert(!doprobingind);
1689 goto TERMINATE;
1690 }
1691
1692 if( doprobingind )
1693 {
1694 SCIP_Bool propagate;
1695 SCIP_Bool cutoff_probing;
1696 SCIP_Bool cutoff;
1697 SCIP_Bool fixed;
1698
1699 #ifndef NDEBUG
1700 SCIP_Real* solvals;
1701 SCIP_CALL( SCIPallocBufferArray(scip, &solvals, nlhdlrexprdata->nvars) );
1702 for( v = 0; v < nlhdlrexprdata->nvars; ++v )
1703 {
1704 solvals[v] = SCIPgetSolVal(scip, sol, nlhdlrexprdata->vars[v]);
1705 }
1706 #endif
1707
1708 propagate = SCIPgetDepth(scip) == 0;
1709
1710 SCIP_CALL( startProbing(scip, nlhdlrdata, nlhdlrexprdata, indicator, probingvars, probingdoms, nprobingvars,
1711 sol, &solcopy, &cutoff_probing) );
1712
1713 #ifndef NDEBUG
1714 for( v = 0; v < nlhdlrexprdata->nvars; ++v )
1715 {
1716 assert(solvals[v] == SCIPgetSolVal(scip, solcopy, nlhdlrexprdata->vars[v])); /*lint !e777*/
1717 }
1718 SCIPfreeBufferArray(scip, &solvals);
1719 #endif
1720
1721 if( propagate )
1722 { /* we are in the root node and startProbing did propagation */
1723 /* probing propagation might have detected infeasibility */
1724 if( cutoff_probing )
1725 {
1726 /* indicator == 1 is infeasible -> set indicator to 0 */
1727 SCIPfreeBufferArrayNull(scip, &probingvars);
1728 SCIPfreeBufferArrayNull(scip, &probingdoms);
1729
1730 SCIP_CALL( SCIPendProbing(scip) );
1731
1732 SCIP_CALL( SCIPfixVar(scip, indicator, 0.0, &cutoff, &fixed) );
1733
1734 if( cutoff )
1735 {
1736 *result = SCIP_CUTOFF;
1737 goto TERMINATE;
1738 }
1739
1740 continue;
1741 }
1742
1743 /* probing propagation in the root node can provide better on/off bounds */
1744 SCIP_CALL( tightenOnBounds(nlhdlrexprdata, nlhdlrdata->scvars, indicator) );
1745 }
1746 }
1747
1748 if( adjrefpoint )
1749 {
1750 /* make sure that when we adjust the point, we don't divide by something too close to 0.0 */
1751 indval = MAX(solval, 0.1);
1752
1753 /* create an adjusted point x^adj = (x* - x0) / z* + x0 */
1754 SCIP_CALL( SCIPcreateSol(scip, &soladj, NULL) );
1755 for( v = 0; v < nlhdlrexprdata->nvars; ++v )
1756 {
1757 if( SCIPvarGetStatus(nlhdlrexprdata->vars[v]) == SCIP_VARSTATUS_FIXED )
1758 continue;
1759
1760 scvdata = getSCVarDataInd(nlhdlrdata->scvars, nlhdlrexprdata->vars[v], indicator, &pos);
1761
1762 /* a non-semicontinuous variable must be linear in expr; skip it */
1763 if( scvdata == NULL )
1764 continue;
1765
1766 SCIP_CALL( SCIPsetSolVal(scip, soladj, nlhdlrexprdata->vars[v],
1767 (SCIPgetSolVal(scip, solcopy, nlhdlrexprdata->vars[v]) - scvdata->vals0[pos]) / indval
1768 + scvdata->vals0[pos]) );
1769 }
1770 for( v = 0; v < nlhdlrexprdata->nindicators; ++v )
1771 {
1772 if( SCIPvarGetStatus(nlhdlrexprdata->indicators[v]) == SCIP_VARSTATUS_FIXED )
1773 continue;
1774
1775 SCIP_CALL( SCIPsetSolVal(scip, soladj, nlhdlrexprdata->indicators[v],
1776 SCIPgetSolVal(scip, solcopy, nlhdlrexprdata->indicators[v])) );
1777 }
1778 if( SCIPvarGetStatus(auxvar) != SCIP_VARSTATUS_FIXED )
1779 {
1780 SCIP_CALL( SCIPsetSolVal(scip, soladj, auxvar, SCIPgetSolVal(scip, solcopy, auxvar)) );
1781 }
1782 }
1783
1784 /* use cuts from every suitable nlhdlr */
1785 for( j = 0; j < nenfos; ++j )
1786 {
1787 SCIP_Bool addedbranchscores2j;
1788 SCIP_NLHDLR* nlhdlr2;
1789 SCIP_NLHDLREXPRDATA* nlhdlr2exprdata;
1790 SCIP_Real nlhdlr2auxvalue;
1791 SCIP_Bool success2;
1792
1793 SCIPgetExprEnfoDataNonlinear(expr, enfoposs[j], &nlhdlr2, &nlhdlr2exprdata, NULL, NULL, NULL, &nlhdlr2auxvalue);
1794 assert(SCIPnlhdlrHasEstimate(nlhdlr2) && nlhdlr2 != nlhdlr);
1795
1796 SCIPdebugMsg(scip, "asking nonlinear handler %s to %sestimate\n", SCIPnlhdlrGetName(nlhdlr2), overestimate ? "over" : "under");
1797
1798 /* ask the nonlinear handler for an estimator */
1799 if( adjrefpoint )
1800 {
1801 SCIP_CALL( SCIPnlhdlrEvalaux(scip, nlhdlr2, expr, nlhdlr2exprdata, &nlhdlr2auxvalue, soladj) );
1802
1803 SCIP_CALL( SCIPnlhdlrEstimate(scip, conshdlr, nlhdlr2, expr,
1804 nlhdlr2exprdata, soladj,
1805 nlhdlr2auxvalue, overestimate, SCIPgetSolVal(scip, solcopy, auxvar),
1806 addbranchscores, rowpreps2, &success2, &addedbranchscores2j) );
1807 }
1808 else
1809 {
1810 SCIP_CALL( SCIPnlhdlrEstimate(scip, conshdlr, nlhdlr2, expr,
1811 nlhdlr2exprdata, solcopy,
1812 nlhdlr2auxvalue, overestimate, SCIPgetSolVal(scip, solcopy, auxvar),
1813 addbranchscores, rowpreps2, &success2, &addedbranchscores2j) );
1814 }
1815
1816 minidx = SCIPgetPtrarrayMinIdx(scip, rowpreps2);
1817 maxidx = SCIPgetPtrarrayMaxIdx(scip, rowpreps2);
1818
1819 assert((success2 && minidx <= maxidx) || (!success2 && minidx > maxidx));
1820
1821 /* perspectivy all cuts from nlhdlr2 and add them to rowpreps */
1822 for( r = minidx; r <= maxidx; ++r )
1823 {
1824 SCIP_Real maxcoef;
1825 SCIP_Real* rowprepcoefs;
1826 SCIP_VAR** rowprepvars;
1827
1828 rowprep = (SCIP_ROWPREP*) SCIPgetPtrarrayVal(scip, rowpreps2, r);
1829 assert(rowprep != NULL);
1830
1831 #ifdef SCIP_DEBUG
1832 SCIPinfoMessage(scip, NULL, "rowprep for expr ");
1833 SCIPprintExpr(scip, expr, NULL);
1834 SCIPinfoMessage(scip, NULL, "rowprep before perspectivy is: \n");
1835 SCIPprintRowprep(scip, rowprep, NULL);
1836 #endif
1837
1838 /* given a rowprep: sum aixi + sum biyi + c, where xi are semicontinuous variables and yi are
1839 * non-semicontinuous variables (which appear in expr linearly, which detect must have ensured),
1840 * perspectivy the semicontinuous part by adding (1-z)(g0 - c - sum aix0i) (the constant is
1841 * treated as belonging to the semicontinuous part)
1842 */
1843
1844 /* we want cst0 = g0 - c - sum aix0i; first add g0 - c */
1845 cst0 = nlhdlrexprdata->exprvals0[i] + SCIProwprepGetSide(rowprep);
1846
1847 maxcoef = 0.0;
1848 rowprepcoefs = SCIProwprepGetCoefs(rowprep);
1849 rowprepvars = SCIProwprepGetVars(rowprep);
1850
1851 for( v = 0; v < SCIProwprepGetNVars(rowprep); ++v )
1852 {
1853 if( REALABS( rowprepcoefs[v]) > maxcoef )
1854 {
1855 maxcoef = REALABS(rowprepcoefs[v]);
1856 }
1857
1858 scvdata = getSCVarDataInd(nlhdlrdata->scvars, rowprepvars[v], indicator, &pos);
1859
1860 /* a non-semicontinuous variable must be linear in expr; skip it */
1861 if( scvdata == NULL )
1862 continue;
1863
1864 cst0 -= rowprepcoefs[v] * scvdata->vals0[pos];
1865 }
1866
1867 /* only perspectivy when the absolute value of cst0 is not too small
1868 * TODO on ex1252a there was cst0=0 - ok to still use the cut?
1869 */
1870 if( cst0 == 0.0 || maxcoef / REALABS(cst0) <= 10.0 / SCIPfeastol(scip) )
1871 {
1872 /* update the rowprep by adding cst0 - cst0*z */
1873 SCIProwprepAddConstant(rowprep, cst0);
1874 SCIP_CALL(SCIPaddRowprepTerm(scip, rowprep, indicator, -cst0));
1875 }
1876 else
1877 {
1878 SCIPfreeRowprep(scip, &rowprep);
1879 continue;
1880 }
1881
1882 SCIP_CALL(SCIPaddRowprepTerm(scip, rowprep, auxvar, -1.0));
1883
1884 SCIPdebugMsg(scip, "rowprep after perspectivy is: \n");
1885 #ifdef SCIP_DEBUG
1886 SCIPprintRowprep(scip, rowprep, NULL);
1887 #endif
1888
1889 SCIP_CALL( SCIPsetPtrarrayVal(scip, rowpreps, nrowpreps, rowprep) );
1890 SCIP_CALL( SCIPsetBoolarrayVal(scip, addedbranchscores2, nrowpreps, addedbranchscores2j) );
1891 ++nrowpreps;
1892 }
1893
1894 SCIP_CALL( SCIPclearPtrarray(scip, rowpreps2) );
1895 }
1896
1897 if( adjrefpoint )
1898 {
1899 SCIP_CALL( SCIPfreeSol(scip, &soladj) );
1900 }
1901
1902 if( doprobingind )
1903 {
1904 SCIP_CALL( SCIPendProbing(scip) );
1905 }
1906
1907 /* add all cuts found for indicator i */
1908 for( r = SCIPgetPtrarrayMinIdx(scip, rowpreps); r <= SCIPgetPtrarrayMaxIdx(scip, rowpreps) && !stop; ++r )
1909 {
1910 SCIP_RESULT resultr;
1911
1912 #ifdef SCIP_DEBUG
1913 SCIPprintRowprep(scip, rowprep, NULL);
1914 #endif
1915 rowprep = (SCIP_ROWPREP*) SCIPgetPtrarrayVal(scip, rowpreps, r);
1916 resultr = SCIP_DIDNOTFIND;
1917
1918 (void) strcat(SCIProwprepGetName(rowprep), "_persp_indicator_");
1919 (void) strcat(SCIProwprepGetName(rowprep), SCIPvarGetName(indicator));
1920
1921 SCIP_CALL( SCIPprocessRowprepNonlinear(scip, nlhdlr, cons, expr, rowprep, overestimate, auxvar, auxvalue,
1922 allowweakcuts, SCIPgetBoolarrayVal(scip, addedbranchscores2, r), addbranchscores, solcopy, &resultr) );
1923
1924 if( resultr == SCIP_SEPARATED )
1925 *result = SCIP_SEPARATED;
1926 else if( resultr == SCIP_CUTOFF )
1927 {
1928 *result = SCIP_CUTOFF;
1929 stop = TRUE;
1930 }
1931 else if( resultr == SCIP_BRANCHED )
1932 {
1933 if( *result != SCIP_SEPARATED && *result != SCIP_REDUCEDDOM )
1934 *result = SCIP_BRANCHED;
1935 }
1936 else if( resultr != SCIP_DIDNOTFIND )
1937 {
1938 SCIPerrorMessage("estimate called by perspective nonlinear handler returned invalid result <%d>\n", resultr);
1939 return SCIP_INVALIDRESULT;
1940 }
1941 }
1942
1943 /* free all rowpreps for indicator i */
1944 for( r = SCIPgetPtrarrayMinIdx(scip, rowpreps); r <= SCIPgetPtrarrayMaxIdx(scip, rowpreps); ++r )
1945 {
1946 rowprep = (SCIP_ROWPREP*) SCIPgetPtrarrayVal(scip, rowpreps, r);
1947 SCIPfreeRowprep(scip, &rowprep);
1948 }
1949
1950 SCIPfreeBufferArrayNull(scip, &probingvars);
1951 SCIPfreeBufferArrayNull(scip, &probingdoms);
1952 SCIP_CALL( SCIPclearPtrarray(scip, rowpreps) );
1953 }
1954
1955 TERMINATE:
1956 SCIP_CALL( SCIPfreeBoolarray(scip, &addedbranchscores2) );
1957 SCIP_CALL( SCIPfreePtrarray(scip, &rowpreps) );
1958 SCIP_CALL( SCIPfreePtrarray(scip, &rowpreps2) );
1959 if( solcopy != sol )
1960 {
1961 SCIP_CALL( SCIPfreeSol(scip, &solcopy) );
1962 }
1963 SCIPfreeBufferArray(scip, &enfoposs);
1964
1965 return SCIP_OKAY;
1966 }
1967
1968
1969 /*
1970 * nonlinear handler specific interface methods
1971 */
1972
1973 /** includes perspective nonlinear handler in nonlinear constraint handler */
1974 SCIP_RETCODE SCIPincludeNlhdlrPerspective(
1975 SCIP* scip /**< SCIP data structure */
1976 )
1977 {
1978 SCIP_NLHDLRDATA* nlhdlrdata;
1979 SCIP_NLHDLR* nlhdlr;
1980
1981 assert(scip != NULL);
1982
1983 /* create nonlinear handler data */
1984 SCIP_CALL( SCIPallocBlockMemory(scip, &nlhdlrdata) );
1985 BMSclearMemory(nlhdlrdata);
1986
1987 SCIP_CALL( SCIPincludeNlhdlrNonlinear(scip, &nlhdlr, NLHDLR_NAME, NLHDLR_DESC, NLHDLR_DETECTPRIORITY,
1988 NLHDLR_ENFOPRIORITY, nlhdlrDetectPerspective, nlhdlrEvalauxPerspective, nlhdlrdata) );
1989 assert(nlhdlr != NULL);
1990
1991 SCIP_CALL( SCIPaddIntParam(scip, "nlhdlr/" NLHDLR_NAME "/maxproprounds",
1992 "maximal number of propagation rounds in probing",
1993 &nlhdlrdata->maxproprounds, FALSE, DEFAULT_MAXPROPROUNDS, -1, INT_MAX, NULL, NULL) );
1994
1995 SCIP_CALL( SCIPaddRealParam(scip, "nlhdlr/" NLHDLR_NAME "/mindomreduction",
1996 "minimal relative reduction in a variable's domain for applying probing",
1997 &nlhdlrdata->mindomreduction, FALSE, DEFAULT_MINDOMREDUCTION, 0.0, 1.0, NULL, NULL) );
1998
1999 SCIP_CALL( SCIPaddRealParam(scip, "nlhdlr/" NLHDLR_NAME "/minviolprobing",
2000 "minimal violation w.r.t. auxiliary variables for applying probing",
2001 &nlhdlrdata->minviolprobing, FALSE, DEFAULT_MINVIOLPROBING, 0.0, SCIP_REAL_MAX, NULL, NULL) );
2002
2003 SCIP_CALL( SCIPaddBoolParam(scip, "nlhdlr/" NLHDLR_NAME "/probingonlyinsepa",
2004 "whether to do probing only in separation",
2005 &nlhdlrdata->probingonlyinsepa, FALSE, DEFAULT_PROBINGONLYINSEPA, NULL, NULL) );
2006
2007 SCIP_CALL( SCIPaddIntParam(scip, "nlhdlr/" NLHDLR_NAME "/probingfreq",
2008 "probing frequency (-1 - no probing, 0 - root node only)",
2009 &nlhdlrdata->probingfreq, FALSE, DEFAULT_PROBINGFREQ, -1, INT_MAX, NULL, NULL) );
2010
2011 SCIP_CALL( SCIPaddBoolParam(scip, "nlhdlr/" NLHDLR_NAME "/convexonly",
2012 "whether perspective cuts are added only for convex expressions",
2013 &nlhdlrdata->convexonly, FALSE, DEFAULT_CONVEXONLY, NULL, NULL) );
2014
2015 SCIP_CALL( SCIPaddBoolParam(scip, "nlhdlr/" NLHDLR_NAME "/tightenbounds",
2016 "whether variable semicontinuity is used to tighten variable bounds",
2017 &nlhdlrdata->tightenbounds, FALSE, DEFAULT_TIGHTENBOUNDS, NULL, NULL) );
2018
2019 SCIP_CALL( SCIPaddBoolParam(scip, "nlhdlr/" NLHDLR_NAME "/adjrefpoint",
2020 "whether to adjust the reference point",
2021 &nlhdlrdata->adjrefpoint, FALSE, DEFAULT_ADJREFPOINT, NULL, NULL) );
2022
2023 SCIPnlhdlrSetCopyHdlr(nlhdlr, nlhdlrCopyhdlrPerspective);
2024 SCIPnlhdlrSetFreeHdlrData(nlhdlr, nlhdlrFreehdlrdataPerspective);
2025 SCIPnlhdlrSetFreeExprData(nlhdlr, nlhdlrFreeExprDataPerspective);
2026 SCIPnlhdlrSetInitExit(nlhdlr, NULL, nlhdlrExitPerspective);
2027 SCIPnlhdlrSetSepa(nlhdlr, nlhdlrInitSepaPerspective, nlhdlrEnfoPerspective, NULL, NULL);
2028
2029 return SCIP_OKAY;
2030 }
2031