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 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 */
|
(1) Event cond_true: |
Condition "scvdata->bvars == NULL", taking true branch. |
|
(2) Event var_compare_op: |
Comparing "scvdata->bvars" to null implies that "scvdata->bvars" might be null. |
| Also see events: |
[var_deref_op] |
222 if( scvdata->bvars == NULL )
223 {
224 assert(scvdata->nbnds == 0 && scvdata->bndssize == 0);
225 found = FALSE;
226 pos = 0;
|
(3) Event if_fallthrough: |
Falling through to end of if statement. |
227 }
228 else
229 {
230 found = SCIPsortedvecFindPtr((void**)scvdata->bvars, SCIPvarComp, (void*)indicator, scvdata->nbnds, &pos);
|
(4) Event if_end: |
End of if statement. |
231 }
232
|
(5) Event cond_false: |
Condition "found", taking false branch. |
233 if( found )
|
(6) Event if_end: |
End of if statement. |
234 return SCIP_OKAY;
235
236 /* ensure sizes */
|
(7) Event cond_false: |
Condition "scvdata->nbnds + 1 > scvdata->bndssize", taking false branch. |
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;
|
(8) Event if_end: |
End of if statement. |
245 }
246 assert(scvdata->nbnds + 1 <= scvdata->bndssize);
247 assert(scvdata->bvars != NULL);
248
249 /* move entries if needed */
|
(9) Event cond_true: |
Condition "i > pos", taking true branch. |
250 for( i = scvdata->nbnds; i > pos; --i )
251 {
|
(10) Event var_deref_op: |
Dereferencing null pointer "scvdata->bvars". |
| Also see events: |
[var_compare_op] |
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
748 {
749 /* curexpr is a non-variable expression, so it belongs to the non-linear part of expr
750 * since the non-linear part of expr must be semicontinuous with respect to
751 * nlhdlrexprdata->indicators[i], curexpr must be semicontinuous
752 */
753 issc = TRUE;
754
755 #ifndef NDEBUG
756 SCIP_CALL( SCIPallocBufferArray(scip, &childvarexprs, norigvars) );
757 SCIP_CALL( SCIPgetExprVarExprs(scip, curexpr, childvarexprs, &nchildvarexprs) );
758
759 /* all nonlinear variables of a sum on/off term should be semicontinuous */
760 for( v = 0; v < nchildvarexprs; ++v )
761 {
762 var = SCIPgetVarExprVar(childvarexprs[v]);
763 scvdata = getSCVarDataInd(nlhdlrdata->scvars, var, nlhdlrexprdata->indicators[i], &pos);
764 assert(scvdata != NULL);
765
766 SCIP_CALL( SCIPreleaseExpr(scip, &childvarexprs[v]) );
767 }
768
769 SCIPfreeBufferArray(scip, &childvarexprs);
770 #endif
771 }
772 }
773
774 if( issc )
775 {
776 /* we know that all vars are semicontinuous with respect to exprdata->indicators; it remains to:
777 * - get or create the scvardata structure for auxvar
778 * - if had to create scvardata, add it to scvars hashmap
779 * - add the indicator and the off value (= curexpr's off value) to scvardata
780 */
781 scvdata = (SCVARDATA*) SCIPhashmapGetImage(nlhdlrdata->scvars, (void*)auxvar);
782 if( scvdata == NULL )
783 {
784 SCIP_CALL( SCIPallocClearBlockMemory(scip, &scvdata) );
785 SCIP_CALL( SCIPallocBlockMemoryArray(scip, &scvdata->bvars, nlhdlrexprdata->nindicators) );
786 SCIP_CALL( SCIPallocBlockMemoryArray(scip, &scvdata->vals0, nlhdlrexprdata->nindicators) );
787 SCIP_CALL( SCIPallocBlockMemoryArray(scip, &scvdata->lbs1, nlhdlrexprdata->nindicators) );
788 SCIP_CALL( SCIPallocBlockMemoryArray(scip, &scvdata->ubs1, nlhdlrexprdata->nindicators) );
789 scvdata->bndssize = nlhdlrexprdata->nindicators;
790 SCIP_CALL( SCIPhashmapInsert(nlhdlrdata->scvars, auxvar, scvdata) );
791 }
792
793 SCIP_CALL( addSCVarIndicator(scip, scvdata, nlhdlrexprdata->indicators[i],
794 SCIPexprGetEvalValue(curexpr), SCIPvarGetLbGlobal(auxvar), SCIPvarGetUbGlobal(auxvar)) );
795 }
796
797 SCIP_CALL( addAuxVar(scip, nlhdlrexprdata, auxvarmap, auxvar) );
798 }
799
800 curexpr = SCIPexpriterGetNext(it);
801 }
802 }
803
804 SCIPfreeExpriter(&it);
805 SCIPhashmapFree(&auxvarmap);
806 SCIPfreeBufferArray(scip, &origvars);
807 SCIPfreeBufferArray(scip, &origvals0);
808 SCIP_CALL( SCIPfreeSol(scip, &sol) );
809
810 return SCIP_OKAY;
811 }
812
813 /*
814 * Probing and bound tightening methods
815 */
816
817 /** go into probing and set some variable bounds */
818 static
819 SCIP_RETCODE startProbing(
820 SCIP* scip, /**< SCIP data structure */
821 SCIP_NLHDLRDATA* nlhdlrdata, /**< nonlinear handler data */
822 SCIP_NLHDLREXPRDATA* nlhdlrexprdata, /**< nlhdlr expression data */
823 SCIP_VAR* indicator, /**< indicator variable */
824 SCIP_VAR** probingvars, /**< array of vars whose bounds we will change in probing */
825 SCIP_INTERVAL* probingdoms, /**< array of intervals to which bounds of probingvars will be changed in probing */
826 int nprobingvars, /**< number of probing vars */
827 SCIP_SOL* sol, /**< solution to be separated */
828 SCIP_SOL** solcopy, /**< buffer for a copy of sol before going into probing; if *solcopy == sol, then copy is created */
829 SCIP_Bool* cutoff_probing /**< pointer to store whether indicator == 1 is infeasible */
830 )
831 {
832 int v;
833 SCIP_Real newlb;
834 SCIP_Real newub;
835 SCIP_Bool propagate;
836
837 propagate = SCIPgetDepth(scip) == 0;
838
839 /* if a copy of sol has not been created yet, then create one now and copy the relevant var values from sol,
840 * because sol can change after SCIPstartProbing, e.g., when linked to the LP solution
841 */
842 if( *solcopy == sol )
843 {
844 SCIP_CALL( SCIPcreateSol(scip, solcopy, NULL) );
845 for( v = 0; v < nlhdlrexprdata->nvars; ++v )
846 {
847 SCIP_CALL( SCIPsetSolVal(scip, *solcopy, nlhdlrexprdata->vars[v], SCIPgetSolVal(scip, sol, nlhdlrexprdata->vars[v])) );
848 }
849 for( v = 0; v < nlhdlrexprdata->nindicators; ++v )
850 {
851 SCIP_CALL( SCIPsetSolVal(scip, *solcopy, nlhdlrexprdata->indicators[v], SCIPgetSolVal(scip, sol, nlhdlrexprdata->indicators[v])) );
852 }
853 }
854
855 /* go into probing */
856 SCIP_CALL( SCIPstartProbing(scip) );
857
858 /* create a probing node */
859 SCIP_CALL( SCIPnewProbingNode(scip) );
860
861 /* set indicator to 1 */
862 SCIP_CALL( SCIPchgVarLbProbing(scip, indicator, 1.0) );
863
864 /* apply stored bounds */
865 for( v = 0; v < nprobingvars; ++v )
866 {
867 newlb = SCIPintervalGetInf(probingdoms[v]);
868 newub = SCIPintervalGetSup(probingdoms[v]);
869
870 if( SCIPisGT(scip, newlb, SCIPvarGetLbLocal(probingvars[v])) || (newlb >= 0.0 && SCIPvarGetLbLocal(probingvars[v]) < 0.0) )
871 {
872 SCIP_CALL( SCIPchgVarLbProbing(scip, probingvars[v], newlb) );
873 }
874 if( SCIPisLT(scip, newub, SCIPvarGetUbLocal(probingvars[v])) || (newub <= 0.0 && SCIPvarGetUbLocal(probingvars[v]) > 0.0) )
875 {
876 SCIP_CALL( SCIPchgVarUbProbing(scip, probingvars[v], newub) );
877 }
878 }
879
880 if( propagate )
881 {
882 SCIP_Longint ndomreds;
883
884 SCIP_CALL( SCIPpropagateProbing(scip, nlhdlrdata->maxproprounds, cutoff_probing, &ndomreds) );
885 }
886
887 return SCIP_OKAY;
888 }
889
890 /** analyse on/off bounds on a variable
891 *
892 * analyses for
893 * 1. tightening bounds in probing for indicator = 1,
894 * 2. fixing indicator / detecting cutoff if one or both states are infeasible,
895 * 3. tightening local bounds if indicator is fixed.
896 *
897 * `probinglb` and `probingub` are only set if `doprobing` is TRUE.
898 * They are either set to bounds that should be used in probing or to `SCIP_INVALID` if bounds on
899 * `var` shouldn't be changed in probing.
900 */
901 static
902 SCIP_RETCODE analyseVarOnoffBounds(
903 SCIP* scip, /**< SCIP data structure */
904 SCIP_NLHDLRDATA* nlhdlrdata, /**< nonlinear handler data */
905 SCIP_VAR* var, /**< variable */
906 SCIP_VAR* indicator, /**< indicator variable */
907 SCIP_Bool indvalue, /**< indicator value for which the bounds are applied */
908 SCIP_Bool* infeas, /**< pointer to store whether infeasibility has been detected */
909 SCIP_Real* probinglb, /**< pointer to store the lower bound to be applied in probing */
910 SCIP_Real* probingub, /**< pointer to store the upper bound to be applied in probing */
911 SCIP_Bool doprobing, /**< whether we currently consider to go into probing */
912 SCIP_Bool* reduceddom /**< pointer to store whether any variables were fixed */
913 )
914 {
915 SCVARDATA* scvdata;
916 int pos;
917 SCIP_Real sclb;
918 SCIP_Real scub;
919 SCIP_Real loclb;
920 SCIP_Real locub;
921 SCIP_Bool bndchgsuccess;
922
923 assert(var != NULL);
924 assert(indicator != NULL);
925 assert(infeas != NULL);
926 assert(reduceddom != NULL);
927
928 /* shouldn't be called if indicator is fixed to !indvalue */
929 assert((indvalue && SCIPvarGetUbLocal(indicator) > 0.5) || (!indvalue && SCIPvarGetLbLocal(indicator) < 0.5));
930
931 *infeas = FALSE;
932 *reduceddom = FALSE;
933 scvdata = getSCVarDataInd(nlhdlrdata->scvars, var, indicator, &pos);
934 if( doprobing )
935 {
936 assert(probinglb != NULL);
937 assert(probingub != NULL);
938
939 *probinglb = SCIP_INVALID;
940 *probingub = SCIP_INVALID;
941 }
942
943 /* nothing to do for non-semicontinuous variables */
944 if( scvdata == NULL )
945 {
946 return SCIP_OKAY;
947 }
948
949 sclb = indvalue ? scvdata->lbs1[pos] : scvdata->vals0[pos];
950 scub = indvalue ? scvdata->ubs1[pos] : scvdata->vals0[pos];
951 loclb = SCIPvarGetLbLocal(var);
952 locub = SCIPvarGetUbLocal(var);
953
954 /* nothing to do for fixed variables */
955 if( SCIPisEQ(scip, loclb, locub) )
956 return SCIP_OKAY;
957
958 /* use a non-redundant lower bound */
959 if( SCIPisGT(scip, sclb, SCIPvarGetLbLocal(var)) || (sclb >= 0.0 && loclb < 0.0) )
960 {
961 /* first check for infeasibility */
962 if( SCIPisFeasGT(scip, sclb, SCIPvarGetUbLocal(var)) )
963 {
964 SCIP_CALL( SCIPfixVar(scip, indicator, indvalue ? 0.0 : 1.0, infeas, &bndchgsuccess) );
965 *reduceddom += bndchgsuccess;
966 if( *infeas )
967 {
968 return SCIP_OKAY;
969 }
970 }
971 else if( nlhdlrdata->tightenbounds &&
972 (SCIPvarGetUbLocal(indicator) <= 0.5 || SCIPvarGetLbLocal(indicator) >= 0.5) )
973 {
974 /* indicator is fixed; due to a previous check, here it can only be fixed to indvalue;
975 * therefore, sclb is valid for the current node
976 */
977
978 if( indvalue == 0 )
979 {
980 assert(sclb == scub); /*lint !e777*/
981 SCIP_CALL( SCIPfixVar(scip, var, sclb, infeas, &bndchgsuccess) );
982 }
983 else
984 {
985 SCIP_CALL( SCIPtightenVarLb(scip, var, sclb, FALSE, infeas, &bndchgsuccess) );
986 }
987 *reduceddom += bndchgsuccess;
988 if( *infeas )
989 {
990 return SCIP_OKAY;
991 }
992 }
993 }
994
995 /* use a non-redundant upper bound */
996 if( SCIPisLT(scip, scub, SCIPvarGetUbLocal(var)) || (scub <= 0.0 && locub > 0.0) )
997 {
998 /* first check for infeasibility */
999 if( SCIPisFeasLT(scip, scub, SCIPvarGetLbLocal(var)) )
1000 {
1001 SCIP_CALL( SCIPfixVar(scip, indicator, indvalue ? 0.0 : 1.0, infeas, &bndchgsuccess) );
1002 *reduceddom += bndchgsuccess;
1003 if( *infeas )
1004 {
1005 return SCIP_OKAY;
1006 }
1007 }
1008 else if( nlhdlrdata->tightenbounds &&
1009 (SCIPvarGetUbLocal(indicator) <= 0.5 || SCIPvarGetLbLocal(indicator) >= 0.5) )
1010 {
1011 /* indicator is fixed; due to a previous check, here it can only be fixed to indvalue;
1012 * therefore, scub is valid for the current node
1013 */
1014
1015 if( indvalue == 0 )
1016 {
1017 assert(sclb == scub); /*lint !e777*/
1018 SCIP_CALL( SCIPfixVar(scip, var, sclb, infeas, &bndchgsuccess) );
1019 }
1020 else
1021 {
1022 SCIP_CALL( SCIPtightenVarUb(scip, var, scub, FALSE, infeas, &bndchgsuccess) );
1023 }
1024 *reduceddom += bndchgsuccess;
1025 if( *infeas )
1026 {
1027 return SCIP_OKAY;
1028 }
1029 }
1030 }
1031
1032 /* If a bound change has been found and indvalue == TRUE, try to use the new bounds.
1033 * This is only done for indvalue == TRUE since this is where enfo asks other nlhdlrs to estimate,
1034 * and at indicator == FALSE we already only have a single point
1035 */
1036 if( doprobing && indvalue && (((scub - sclb) / (locub - loclb)) <= 1.0 - nlhdlrdata->mindomreduction ||
1037 (sclb >= 0.0 && loclb < 0.0) || (scub <= 0.0 && locub > 0.0)) )
1038 {
1039 *probinglb = sclb;
1040 *probingub = scub;
1041 }
1042
1043 SCIPdebugMsg(scip, "%s in [%g, %g] instead of [%g, %g] (vals0 = %g)\n", SCIPvarGetName(var), sclb, scub,
1044 SCIPvarGetLbLocal(var), SCIPvarGetUbLocal(var), scvdata->vals0[pos]);
1045
1046 return SCIP_OKAY;
1047 }
1048
1049 /** looks for bound tightenings to be applied either in the current node or in probing
1050 *
1051 * Loops through both possible values of indicator and calls analyseVarOnoffBounds().
1052 * Might update the `*doprobing` flag by setting it to `FALSE` if:
1053 * - indicator is fixed or
1054 * - analyseVarOnoffBounds() hasn't found a sufficient improvement at indicator==1.
1055 *
1056 * If `*doprobing==TRUE`, stores bounds suggested by analyseVarOnoffBounds() in order to apply them in probing together
1057 * with the fixing `indicator=1`.
1058 */
1059 static
1060 SCIP_RETCODE analyseOnoffBounds(
1061 SCIP* scip, /**< SCIP data structure */
1062 SCIP_NLHDLRDATA* nlhdlrdata, /**< nonlinear handler data */
1063 SCIP_NLHDLREXPRDATA* nlhdlrexprdata, /**< nlhdlr expression data */
1064 SCIP_VAR* indicator, /**< indicator variable */
1065 SCIP_VAR*** probingvars, /**< array to store variables whose bounds will be changed in probing */
1066 SCIP_INTERVAL** probingdoms, /**< array to store bounds to be applied in probing */
1067 int* nprobingvars, /**< pointer to store number of vars whose bounds will be changed in probing */
1068 SCIP_Bool* doprobing, /**< pointer to the flag telling whether we want to do probing */
1069 SCIP_RESULT* result /**< pointer to store the result */
1070 )
1071 {
1072 int v;
1073 SCIP_VAR* var;
1074 SCIP_Bool infeas;
1075 int b;
1076 SCIP_Real probinglb = SCIP_INVALID;
1077 SCIP_Real probingub = SCIP_INVALID;
1078 SCIP_Bool changed;
1079 SCIP_Bool reduceddom;
1080
1081 assert(indicator != NULL);
1082 assert(nprobingvars != NULL);
1083 assert(doprobing != NULL);
1084 assert(result != NULL);
1085
1086 changed = FALSE;
1087
1088 /* no probing if indicator already fixed */
1089 if( SCIPvarGetUbLocal(indicator) <= 0.5 || SCIPvarGetLbLocal(indicator) >= 0.5 )
1090 {
1091 *doprobing = FALSE;
1092 }
1093
1094 /* consider each possible value of indicator */
1095 for( b = 0; b < 2; ++b )
1096 {
1097 for( v = 0; v < nlhdlrexprdata->nvars; ++v )
1098 {
1099 /* nothing left to do if indicator is already fixed to !indvalue
1100 * (checked in the inner loop since analyseVarOnoff bounds might fix the indicator)
1101 */
1102 if( (b == 1 && SCIPvarGetUbLocal(indicator) <= 0.5) || (b == 0 && SCIPvarGetLbLocal(indicator) >= 0.5) )
1103 {
1104 *doprobing = FALSE;
1105 break;
1106 }
1107
1108 var = nlhdlrexprdata->vars[v];
1109
1110 SCIP_CALL( analyseVarOnoffBounds(scip, nlhdlrdata, var, indicator, b == 1, &infeas, &probinglb,
1111 &probingub, *doprobing, &reduceddom) );
1112
1113 if( infeas )
1114 {
1115 *result = SCIP_CUTOFF;
1116 *doprobing = FALSE;
1117 return SCIP_OKAY;
1118 }
1119 else if( reduceddom )
1120 {
1121 *result = SCIP_REDUCEDDOM;
1122 }
1123
1124 if( !(*doprobing) )
1125 continue;
1126
1127 /* if bounds to be applied in probing have been found, store them */
1128 if( probinglb != SCIP_INVALID ) /*lint !e777*/
1129 {
1130 assert(probingub != SCIP_INVALID); /*lint !e777*/
1131
1132 SCIP_CALL( SCIPreallocBufferArray(scip, probingvars, *nprobingvars + 1) );
1133 SCIP_CALL( SCIPreallocBufferArray(scip, probingdoms, *nprobingvars + 1) );
1134 (*probingvars)[*nprobingvars] = var;
1135 (*probingdoms)[*nprobingvars].inf = probinglb;
1136 (*probingdoms)[*nprobingvars].sup = probingub;
1137 ++*nprobingvars;
1138
1139 changed = TRUE;
1140 }
1141 }
1142 }
1143
1144 if( !changed )
1145 {
1146 *doprobing = FALSE;
1147 }
1148
1149 return SCIP_OKAY;
1150 }
1151
1152 /** saves local bounds on all expression variables, including auxiliary variables, obtained from propagating
1153 * indicator == 1 to the corresponding SCVARDATA (should only be used in the root node)
1154 */
1155 static
1156 SCIP_RETCODE tightenOnBounds(
1157 SCIP_NLHDLREXPRDATA* nlhdlrexprdata, /**< nlhdlr expression data */
1158 SCIP_HASHMAP* scvars, /**< hashmap with semicontinuous variables */
1159 SCIP_VAR* indicator /**< indicator variable */
1160 )
1161 {
1162 int v;
1163 SCIP_VAR* var;
1164 SCVARDATA* scvdata;
1165 int pos;
1166 SCIP_Real lb;
1167 SCIP_Real ub;
1168
1169 for( v = 0; v < nlhdlrexprdata->nvars; ++v )
1170 {
1171 var = nlhdlrexprdata->vars[v];
1172 lb = SCIPvarGetLbLocal(var);
1173 ub = SCIPvarGetUbLocal(var);
1174 scvdata = getSCVarDataInd(scvars, var, indicator, &pos);
1175
1176 if( scvdata != NULL )
1177 {
1178 scvdata->lbs1[pos] = MAX(scvdata->lbs1[pos], lb);
1179 scvdata->ubs1[pos] = MIN(scvdata->ubs1[pos], ub);
1180 }
1181 }
1182
1183 return SCIP_OKAY;
1184 }
1185
1186 /*
1187 * Callback methods of nonlinear handler
1188 */
1189
1190 /** nonlinear handler copy callback */
1191 static
1192 SCIP_DECL_NLHDLRCOPYHDLR(nlhdlrCopyhdlrPerspective)
1193 { /*lint --e{715}*/
1194 assert(targetscip != NULL);
1195 assert(sourcenlhdlr != NULL);
1196 assert(strcmp(SCIPnlhdlrGetName(sourcenlhdlr), NLHDLR_NAME) == 0);
1197
1198 SCIP_CALL( SCIPincludeNlhdlrPerspective(targetscip) );
1199
1200 return SCIP_OKAY;
1201 }
1202
1203
1204 /** callback to free data of handler */
1205 static
1206 SCIP_DECL_NLHDLRFREEHDLRDATA(nlhdlrFreehdlrdataPerspective)
1207 { /*lint --e{715}*/
1208 SCIPfreeBlockMemory(scip, nlhdlrdata);
1209
1210 return SCIP_OKAY;
1211 }
1212
1213
1214 /** callback to free expression specific data */
1215 static
1216 SCIP_DECL_NLHDLRFREEEXPRDATA(nlhdlrFreeExprDataPerspective)
1217 { /*lint --e{715}*/
1218 SCIP_CALL( freeNlhdlrExprData(scip, *nlhdlrexprdata) );
1219 SCIPfreeBlockMemory(scip, nlhdlrexprdata);
1220
1221 return SCIP_OKAY;
1222 }
1223
1224 /** callback to be called in initialization */
1225 #if 0
1226 static
1227 SCIP_DECL_NLHDLRINIT(nlhdlrInitPerspective)
1228 { /*lint --e{715}*/
1229 return SCIP_OKAY;
1230 }
1231 #endif
1232
1233 /** callback to be called in deinitialization */
1234 static
1235 SCIP_DECL_NLHDLREXIT(nlhdlrExitPerspective)
1236 { /*lint --e{715}*/
1237 SCIP_HASHMAPENTRY* entry;
1238 SCVARDATA* data;
1239 int c;
1240 SCIP_NLHDLRDATA* nlhdlrdata;
1241
1242 nlhdlrdata = SCIPnlhdlrGetData(nlhdlr);
1243 assert(nlhdlrdata != NULL);
1244
1245 if( nlhdlrdata->scvars != NULL )
1246 {
1247 for( c = 0; c < SCIPhashmapGetNEntries(nlhdlrdata->scvars); ++c )
1248 {
1249 entry = SCIPhashmapGetEntry(nlhdlrdata->scvars, c);
1250 if( entry != NULL )
1251 {
1252 data = (SCVARDATA*) SCIPhashmapEntryGetImage(entry);
1253 SCIPfreeBlockMemoryArray(scip, &data->ubs1, data->bndssize);
1254 SCIPfreeBlockMemoryArray(scip, &data->lbs1, data->bndssize);
1255 SCIPfreeBlockMemoryArray(scip, &data->vals0, data->bndssize);
1256 SCIPfreeBlockMemoryArray(scip, &data->bvars, data->bndssize);
1257 SCIPfreeBlockMemory(scip, &data);
1258 }
1259 }
1260 SCIPhashmapFree(&nlhdlrdata->scvars);
1261 assert(nlhdlrdata->scvars == NULL);
1262 }
1263
1264 return SCIP_OKAY;
1265 }
1266
1267 /** callback to detect structure in expression tree
1268 *
1269 * We are looking for expressions g(x), where x is a vector of semicontinuous variables that all share at least one
1270 * indicator variable.
1271 */
1272 static
1273 SCIP_DECL_NLHDLRDETECT(nlhdlrDetectPerspective)
1274 { /*lint --e{715}*/
1275 SCIP_NLHDLRDATA* nlhdlrdata;
1276 SCIP_EXPR** varexprs;
1277 SCIP_Bool success = FALSE;
1278 int i;
1279 SCIP_Bool hassepabelow = FALSE;
1280 SCIP_Bool hassepaabove = FALSE;
1281 SCIP_Bool hasnondefault = FALSE;
1282
1283 nlhdlrdata = SCIPnlhdlrGetData(nlhdlr);
1284
1285 assert(scip != NULL);
1286 assert(nlhdlr != NULL);
1287 assert(expr != NULL);
1288 assert(participating != NULL);
1289 assert(enforcing != NULL);
1290 assert(nlhdlrexprdata != NULL);
1291 assert(nlhdlrdata != NULL);
1292
1293 /* do not run if we will have no auxvar to add a cut for */
1294 if( SCIPgetExprNAuxvarUsesNonlinear(expr) == 0 )
1295 return SCIP_OKAY;
1296
1297 if( SCIPgetNBinVars(scip) == 0 )
1298 {
1299 SCIPdebugMsg(scip, "problem has no binary variables, not running perspective detection\n");
1300 return SCIP_OKAY;
1301 }
1302
1303 for( i = 0; i < SCIPgetExprNEnfosNonlinear(expr); ++i )
1304 {
1305 SCIP_NLHDLR* nlhdlr2;
1306 SCIP_NLHDLR_METHOD nlhdlr2participates;
1307 SCIP_Bool sepabelowusesactivity;
1308 SCIP_Bool sepaaboveusesactivity;
1309 SCIPgetExprEnfoDataNonlinear(expr, i, &nlhdlr2, NULL, &nlhdlr2participates, &sepabelowusesactivity, &sepaaboveusesactivity, NULL);
1310
1311 if( (nlhdlr2participates & SCIP_NLHDLR_METHOD_SEPABOTH) == 0 )
1312 continue;
1313
1314 if( !SCIPnlhdlrHasEstimate(nlhdlr2) )
1315 continue;
1316
1317 if( strcmp(SCIPnlhdlrGetName(nlhdlr2), "default") != 0 )
1318 hasnondefault = TRUE;
1319
1320 /* If we are supposed to run only on convex expressions, than check whether there is a nlhdlr
1321 * that participates in separation without using activity for it. Otherwise, check for
1322 * participation regardless of activity usage.
1323 */
1324 if( (nlhdlr2participates & SCIP_NLHDLR_METHOD_SEPABELOW) && (!nlhdlrdata->convexonly || !sepabelowusesactivity) )
1325 hassepabelow = TRUE;
1326
1327 if( (nlhdlr2participates & SCIP_NLHDLR_METHOD_SEPAABOVE) && (!nlhdlrdata->convexonly || !sepaaboveusesactivity) )
1328 hassepaabove = TRUE;
1329 }
1330
1331 /* If a sum expression is handled only by default nlhdlr, then all the children will have auxiliary vars.
1332 * Since the sum will then be linear in auxiliary variables, perspective can't improve anything for it
1333 */
1334 if( SCIPisExprSum(scip, expr) && !hasnondefault )
1335 {
1336 SCIPdebugMsg(scip, "sum expr only has default exprhdlr, not running perspective detection\n");
1337 return SCIP_OKAY;
1338 }
1339
1340 /* If no other nlhdlr separates, neither does perspective (if convexonly, only separation
1341 * without using activity counts)
1342 */
1343 if( !hassepabelow && !hassepaabove )
1344 {
1345 SCIPdebugMsg(scip, "no nlhdlr separates without using activity, not running perspective detection\n");
1346 return SCIP_OKAY;
1347 }
1348
1349 #ifdef SCIP_DEBUG
1350 SCIPdebugMsg(scip, "Called perspective detect, expr = %p: ", (void*)expr);
1351 SCIPprintExpr(scip, expr, NULL);
1352 SCIPdebugMsgPrint(scip, "\n");
1353 #endif
1354
1355 /* allocate memory */
1356 SCIP_CALL( SCIPallocClearBlockMemory(scip, nlhdlrexprdata) );
1357 if( nlhdlrdata->scvars == NULL )
1358 {
1359 SCIP_CALL( SCIPhashmapCreate(&(nlhdlrdata->scvars), SCIPblkmem(scip), SCIPgetNVars(scip)) );
1360 }
1361
1362 /* save varexprs to nlhdlrexprdata */
1363 SCIP_CALL( SCIPgetExprNVars(scip, expr, &(*nlhdlrexprdata)->nvars) );
1364 SCIP_CALL( SCIPallocBlockMemoryArray(scip, &(*nlhdlrexprdata)->vars, (*nlhdlrexprdata)->nvars) );
1365 SCIP_CALL( SCIPallocBufferArray(scip, &varexprs, (*nlhdlrexprdata)->nvars) );
1366 (*nlhdlrexprdata)->varssize = (*nlhdlrexprdata)->nvars;
1367 SCIP_CALL( SCIPgetExprVarExprs(scip, expr, varexprs, &(*nlhdlrexprdata)->nvars) );
1368 for( i = 0; i < (*nlhdlrexprdata)->nvars; ++i )
1369 {
1370 (*nlhdlrexprdata)->vars[i] = SCIPgetVarExprVar(varexprs[i]);
1371 SCIP_CALL( SCIPreleaseExpr(scip, &varexprs[i]) );
1372 SCIP_CALL( SCIPcaptureVar(scip, (*nlhdlrexprdata)->vars[i]) );
1373 }
1374 SCIPsortPtr((void**) (*nlhdlrexprdata)->vars, SCIPvarComp, (*nlhdlrexprdata)->nvars);
1375 SCIPfreeBufferArray(scip, &varexprs);
1376
1377 /* check if expr is semicontinuous and save indicator variables */
1378 SCIP_CALL( exprIsSemicontinuous(scip, nlhdlrdata, *nlhdlrexprdata, expr, &success) );
1379
1380 if( success )
1381 {
1382 assert(*nlhdlrexprdata != NULL);
1383 assert((*nlhdlrexprdata)->nindicators > 0);
1384
1385 if( hassepaabove )
1386 *participating |= SCIP_NLHDLR_METHOD_SEPAABOVE;
1387 if( hassepabelow )
1388 *participating |= SCIP_NLHDLR_METHOD_SEPABELOW;
1389
1390 #ifdef SCIP_DEBUG
1391 SCIPinfoMessage(scip, NULL, "detected an on/off expr: ");
1392 SCIPprintExpr(scip, expr, NULL);
1393 SCIPinfoMessage(scip, NULL, "\n");
1394 #endif
1395 }
1396 else if( *nlhdlrexprdata != NULL )
1397 {
1398 SCIP_CALL( nlhdlrFreeExprDataPerspective(scip, nlhdlr, expr, nlhdlrexprdata) );
1399 }
1400
1401 return SCIP_OKAY;
1402 }
1403
1404
1405 /** auxiliary evaluation callback of nonlinear handler */
1406 static
1407 SCIP_DECL_NLHDLREVALAUX(nlhdlrEvalauxPerspective)
1408 { /*lint --e{715}*/
1409 int e;
1410 SCIP_Real maxdiff;
1411 SCIP_Real auxvarvalue;
1412 SCIP_Real enfoauxval;
1413
1414 assert(scip != NULL);
1415 assert(expr != NULL);
1416 assert(auxvalue != NULL);
1417
1418 auxvarvalue = SCIPgetSolVal(scip, sol, SCIPgetExprAuxVarNonlinear(expr));
1419 maxdiff = 0.0;
1420 *auxvalue = auxvarvalue;
1421
1422 /* use the auxvalue from one of the other nlhdlrs that estimates for this expr: take the one that is farthest
1423 * from the current value of auxvar
1424 */
1425 for( e = 0; e < SCIPgetExprNEnfosNonlinear(expr); ++e )
1426 {
1427 SCIP_NLHDLR* nlhdlr2;
1428 SCIP_NLHDLREXPRDATA* nlhdlr2exprdata;
1429 SCIP_NLHDLR_METHOD nlhdlr2participation;
1430
1431 SCIPgetExprEnfoDataNonlinear(expr, e, &nlhdlr2, &nlhdlr2exprdata, &nlhdlr2participation, NULL, NULL, NULL);
1432
1433 /* skip nlhdlr that do not participate or do not provide estimate */
1434 if( (nlhdlr2participation & SCIP_NLHDLR_METHOD_SEPABOTH) == 0 || !SCIPnlhdlrHasEstimate(nlhdlr2) )
1435 continue;
1436
1437 SCIP_CALL( SCIPnlhdlrEvalaux(scip, nlhdlr2, expr, nlhdlr2exprdata, &enfoauxval, sol) );
1438
1439 SCIPsetExprEnfoAuxValueNonlinear(expr, e, enfoauxval);
1440
1441 if( REALABS(enfoauxval - auxvarvalue) > maxdiff && enfoauxval != SCIP_INVALID ) /*lint !e777*/
1442 {
1443 maxdiff = REALABS(enfoauxval - auxvarvalue);
1444 *auxvalue = enfoauxval;
1445 }
1446 }
1447
1448 return SCIP_OKAY;
1449 }
1450
1451 /** separation initialization method of a nonlinear handler */
1452 static
1453 SCIP_DECL_NLHDLRINITSEPA(nlhdlrInitSepaPerspective)
1454 { /*lint --e{715}*/
1455 int sindicators;
1456
1457 sindicators = nlhdlrexprdata->nindicators;
1458
1459 /* compute 'off' values of expr and subexprs (and thus auxvars too) */
1460 SCIP_CALL( computeOffValues(scip, SCIPnlhdlrGetData(nlhdlr), nlhdlrexprdata, expr) );
1461
1462 /* some indicator variables might have been removed if evaluation failed, check how many remain */
1463 if( nlhdlrexprdata->nindicators == 0 )
1464 {
1465 SCIPfreeBlockMemoryArray(scip, &nlhdlrexprdata->indicators, sindicators);
1466 SCIPfreeBlockMemoryArray(scip, &nlhdlrexprdata->exprvals0, sindicators);
1467 }
1468 else if( nlhdlrexprdata->nindicators < sindicators )
1469 {
1470 SCIP_CALL( SCIPreallocBlockMemoryArray(scip, &nlhdlrexprdata->indicators, sindicators, nlhdlrexprdata->nindicators) );
1471 SCIP_CALL( SCIPreallocBlockMemoryArray(scip, &nlhdlrexprdata->exprvals0, sindicators, nlhdlrexprdata->nindicators) );
1472 }
1473
1474 return SCIP_OKAY;
1475 }
1476
1477
1478 #if 0
1479 /** separation deinitialization method of a nonlinear handler (called during CONSEXITSOL) */
1480 static
1481 SCIP_DECL_NLHDLREXITSEPA(nlhdlrExitSepaPerspective)
1482 { /*lint --e{715}*/
1483 SCIPerrorMessage("method of perspective nonlinear handler not implemented yet\n");
1484 SCIPABORT(); /*lint --e{527}*/
1485
1486 return SCIP_OKAY;
1487 }
1488 #endif
1489
1490 /** nonlinear handler enforcement callback
1491 *
1492 * "Perspectivies" cuts produced by other nonlinear handlers.
1493 *
1494 * Suppose that we want to separate \f$x\f$ from the set \f$\{ x : g(x) \leq 0\}\f$.
1495 * 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$,
1496 * 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]
1497 * 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$.
1498 */
1499 static
1500 SCIP_DECL_NLHDLRENFO(nlhdlrEnfoPerspective)
1501 { /*lint --e{715}*/
1502 SCIP_ROWPREP* rowprep;
1503 SCIP_VAR* auxvar;
1504 int i;
1505 int j;
1506 SCIP_NLHDLRDATA* nlhdlrdata;
1507 SCIP_Real cst0;
1508 SCIP_VAR* indicator;
1509 SCIP_PTRARRAY* rowpreps2;
1510 SCIP_PTRARRAY* rowpreps;
1511 int nrowpreps;
1512 SCIP_SOL* solcopy;
1513 SCIP_Bool doprobing;
1514 SCIP_BOOLARRAY* addedbranchscores2;
1515 SCIP_Bool stop;
1516 int nenfos;
1517 int* enfoposs;
1518 SCIP_SOL* soladj;
1519 int pos;
1520 SCVARDATA* scvdata;
1521
1522 nlhdlrdata = SCIPnlhdlrGetData(nlhdlr);
1523
1524 #ifdef SCIP_DEBUG
1525 SCIPinfoMessage(scip, NULL, "enforcement method of perspective nonlinear handler called for expr %p: ", (void*)expr);
1526 SCIP_CALL( SCIPprintExpr(scip, expr, NULL) );
1527 SCIPinfoMessage(scip, NULL, " at\n");
1528 for( i = 0; i < nlhdlrexprdata->nvars; ++i )
1529 {
1530 SCIPinfoMessage(scip, NULL, "%s = %g\n", SCIPvarGetName(nlhdlrexprdata->vars[i]),
1531 SCIPgetSolVal(scip, sol, nlhdlrexprdata->vars[i]));
1532 }
1533 SCIPinfoMessage(scip, NULL, "%s = %g", SCIPvarGetName(SCIPgetExprAuxVarNonlinear(expr)),
1534 SCIPgetSolVal(scip, sol, SCIPgetExprAuxVarNonlinear(expr)));
1535 #endif
1536
1537 assert(scip != NULL);
1538 assert(expr != NULL);
1539 assert(conshdlr != NULL);
1540 assert(nlhdlrexprdata != NULL);
1541 assert(nlhdlrdata != NULL);
1542
1543 if( nlhdlrexprdata->nindicators == 0 )
1544 {
1545 /* we might have removed all indicators in initsepa */
1546 *result = SCIP_DIDNOTRUN;
1547 return SCIP_OKAY;
1548 }
1549
1550 auxvar = SCIPgetExprAuxVarNonlinear(expr);
1551 assert(auxvar != NULL);
1552
1553 /* detect should have picked only those expressions for which at least one other nlhdlr can enforce */
1554 assert(SCIPgetExprNEnfosNonlinear(expr) > 1);
1555
1556 SCIP_CALL( SCIPallocBufferArray(scip, &enfoposs, SCIPgetExprNEnfosNonlinear(expr) - 1) );
1557
1558 doprobing = FALSE;
1559 nenfos = 0;
1560 soladj = NULL;
1561
1562 /* find suitable nlhdlrs and check if there is enough violation to do probing */
1563 for( j = 0; j < SCIPgetExprNEnfosNonlinear(expr); ++j )
1564 {
1565 SCIP_NLHDLR* nlhdlr2;
1566 SCIP_NLHDLREXPRDATA* nlhdlr2exprdata;
1567 SCIP_NLHDLR_METHOD nlhdlr2participate;
1568 SCIP_Real nlhdlr2auxvalue;
1569 SCIP_Real violation;
1570 SCIP_Bool violbelow;
1571 SCIP_Bool violabove;
1572 SCIP_Bool sepausesactivity = FALSE;
1573
1574 SCIPgetExprEnfoDataNonlinear(expr, j, &nlhdlr2, &nlhdlr2exprdata, &nlhdlr2participate, !overestimate ? &sepausesactivity : NULL, overestimate ? &sepausesactivity: NULL, &nlhdlr2auxvalue); /*lint !e826*/
1575
1576 if( nlhdlr2 == nlhdlr )
1577 continue;
1578
1579 /* if nlhdlr2 cannot estimate, then cannot use it */
1580 if( !SCIPnlhdlrHasEstimate(nlhdlr2) )
1581 continue;
1582
1583 /* if nlhdlr2 does not participate in the separation on the desired side (overestimate), then skip it */
1584 if( (nlhdlr2participate & (overestimate ? SCIP_NLHDLR_METHOD_SEPAABOVE : SCIP_NLHDLR_METHOD_SEPABELOW)) == 0 )
1585 continue;
1586
1587 /* if only working on convex-looking expressions, then skip nlhdlr if it uses activity for estimates */
1588 if( nlhdlrdata->convexonly && sepausesactivity )
1589 continue;
1590
1591 /* evalaux should have called evalaux of nlhdlr2 by now
1592 * check whether handling the violation for nlhdlr2 requires under- or overestimation and this fits to
1593 * overestimate flag
1594 */
1595 SCIP_CALL( SCIPgetExprAbsAuxViolationNonlinear(scip, expr, nlhdlr2auxvalue, sol, &violation, &violbelow,
1596 &violabove) );
1597 assert(violation >= 0.0);
1598
1599 if( (overestimate && !violabove) || (!overestimate && !violbelow) )
1600 continue;
1601
1602 /* if violation is small, cuts would likely be weak - skip perspectification */
1603 if( !allowweakcuts && violation < SCIPfeastol(scip) )
1604 continue;
1605
1606 enfoposs[nenfos] = j;
1607 ++nenfos;
1608
1609 /* enable probing if tightening the domain could be useful for nlhdlr and violation is above threshold */
1610 if( sepausesactivity && violation >= nlhdlrdata->minviolprobing )
1611 doprobing = TRUE;
1612 }
1613
1614 if( nenfos == 0 )
1615 {
1616 *result = SCIP_DIDNOTRUN;
1617 SCIPfreeBufferArray(scip, &enfoposs);
1618 return SCIP_OKAY;
1619 }
1620
1621 /* check probing frequency against depth in b&b tree */
1622 if( nlhdlrdata->probingfreq == -1 || (nlhdlrdata->probingfreq == 0 && SCIPgetDepth(scip) != 0) ||
1623 (nlhdlrdata->probingfreq > 0 && SCIPgetDepth(scip) % nlhdlrdata->probingfreq != 0) )
1624 doprobing = FALSE;
1625
1626 /* if addbranchscores is TRUE, then we can assume to be in enforcement and not in separation */
1627 if( nlhdlrdata->probingonlyinsepa && addbranchscores )
1628 doprobing = FALSE;
1629
1630 /* disable probing if already being in probing or if in a subscip */
1631 if( SCIPinProbing(scip) || SCIPgetSubscipDepth(scip) != 0 )
1632 doprobing = FALSE;
1633
1634 nrowpreps = 0;
1635 *result = SCIP_DIDNOTFIND;
1636 solcopy = sol;
1637 stop = FALSE;
1638
1639 SCIP_CALL( SCIPcreatePtrarray(scip, &rowpreps2) );
1640 SCIP_CALL( SCIPcreatePtrarray(scip, &rowpreps) );
1641 SCIP_CALL( SCIPcreateBoolarray(scip, &addedbranchscores2) );
1642
1643 /* build cuts for every indicator variable */
1644 for( i = 0; i < nlhdlrexprdata->nindicators && !stop; ++i )
1645 {
1646 int v;
1647 int minidx;
1648 int maxidx;
1649 int r;
1650 SCIP_VAR** probingvars;
1651 SCIP_INTERVAL* probingdoms;
1652 int nprobingvars;
1653 SCIP_Bool doprobingind;
1654 SCIP_Real indval;
1655 SCIP_Real solval;
1656 SCIP_Bool adjrefpoint;
1657
1658 indicator = nlhdlrexprdata->indicators[i];
1659 probingvars = NULL;
1660 probingdoms = NULL;
1661 nprobingvars = 0;
1662 doprobingind = doprobing;
1663 solval = SCIPgetSolVal(scip, solcopy, indicator);
1664 adjrefpoint = nlhdlrdata->adjrefpoint && !SCIPisFeasEQ(scip, solval, 1.0);
1665
1666 SCIP_CALL( analyseOnoffBounds(scip, nlhdlrdata, nlhdlrexprdata, indicator, &probingvars, &probingdoms,
1667 &nprobingvars, &doprobingind, result) );
1668
1669 /* don't add perspective cuts for fixed indicators since there is no use for perspectivy */
1670 if( SCIPvarGetLbLocal(indicator) >= 0.5 )
1671 {
1672 assert(!doprobingind);
1673 continue;
1674 }
1675 if( SCIPvarGetUbLocal(indicator) <= 0.5 )
1676 { /* this case is stronger as it implies that everything is fixed;
1677 * therefore we are now happy
1678 */
1679 assert(!doprobingind);
1680 goto TERMINATE;
1681 }
1682
1683 if( doprobingind )
1684 {
1685 SCIP_Bool propagate;
1686 SCIP_Bool cutoff_probing;
1687 SCIP_Bool cutoff;
1688 SCIP_Bool fixed;
1689
1690 #ifndef NDEBUG
1691 SCIP_Real* solvals;
1692 SCIP_CALL( SCIPallocBufferArray(scip, &solvals, nlhdlrexprdata->nvars) );
1693 for( v = 0; v < nlhdlrexprdata->nvars; ++v )
1694 {
1695 solvals[v] = SCIPgetSolVal(scip, sol, nlhdlrexprdata->vars[v]);
1696 }
1697 #endif
1698
1699 propagate = SCIPgetDepth(scip) == 0;
1700
1701 SCIP_CALL( startProbing(scip, nlhdlrdata, nlhdlrexprdata, indicator, probingvars, probingdoms, nprobingvars,
1702 sol, &solcopy, &cutoff_probing) );
1703
1704 #ifndef NDEBUG
1705 for( v = 0; v < nlhdlrexprdata->nvars; ++v )
1706 {
1707 assert(solvals[v] == SCIPgetSolVal(scip, solcopy, nlhdlrexprdata->vars[v])); /*lint !e777*/
1708 }
1709 SCIPfreeBufferArray(scip, &solvals);
1710 #endif
1711
1712 if( propagate )
1713 { /* we are in the root node and startProbing did propagation */
1714 /* probing propagation might have detected infeasibility */
1715 if( cutoff_probing )
1716 {
1717 /* indicator == 1 is infeasible -> set indicator to 0 */
1718 SCIPfreeBufferArrayNull(scip, &probingvars);
1719 SCIPfreeBufferArrayNull(scip, &probingdoms);
1720
1721 SCIP_CALL( SCIPendProbing(scip) );
1722
1723 SCIP_CALL( SCIPfixVar(scip, indicator, 0.0, &cutoff, &fixed) );
1724
1725 if( cutoff )
1726 {
1727 *result = SCIP_CUTOFF;
1728 goto TERMINATE;
1729 }
1730
1731 continue;
1732 }
1733
1734 /* probing propagation in the root node can provide better on/off bounds */
1735 SCIP_CALL( tightenOnBounds(nlhdlrexprdata, nlhdlrdata->scvars, indicator) );
1736 }
1737 }
1738
1739 if( adjrefpoint )
1740 {
1741 /* make sure that when we adjust the point, we don't divide by something too close to 0.0 */
1742 indval = MAX(solval, 0.1);
1743
1744 /* create an adjusted point x^adj = (x* - x0) / z* + x0 */
1745 SCIP_CALL( SCIPcreateSol(scip, &soladj, NULL) );
1746 for( v = 0; v < nlhdlrexprdata->nvars; ++v )
1747 {
1748 if( SCIPvarGetStatus(nlhdlrexprdata->vars[v]) == SCIP_VARSTATUS_FIXED )
1749 continue;
1750
1751 scvdata = getSCVarDataInd(nlhdlrdata->scvars, nlhdlrexprdata->vars[v], indicator, &pos);
1752
1753 /* a non-semicontinuous variable must be linear in expr; skip it */
1754 if( scvdata == NULL )
1755 continue;
1756
1757 SCIP_CALL( SCIPsetSolVal(scip, soladj, nlhdlrexprdata->vars[v],
1758 (SCIPgetSolVal(scip, solcopy, nlhdlrexprdata->vars[v]) - scvdata->vals0[pos]) / indval
1759 + scvdata->vals0[pos]) );
1760 }
1761 for( v = 0; v < nlhdlrexprdata->nindicators; ++v )
1762 {
1763 if( SCIPvarGetStatus(nlhdlrexprdata->indicators[v]) == SCIP_VARSTATUS_FIXED )
1764 continue;
1765
1766 SCIP_CALL( SCIPsetSolVal(scip, soladj, nlhdlrexprdata->indicators[v],
1767 SCIPgetSolVal(scip, solcopy, nlhdlrexprdata->indicators[v])) );
1768 }
1769 if( SCIPvarGetStatus(auxvar) != SCIP_VARSTATUS_FIXED )
1770 {
1771 SCIP_CALL( SCIPsetSolVal(scip, soladj, auxvar, SCIPgetSolVal(scip, solcopy, auxvar)) );
1772 }
1773 }
1774
1775 /* use cuts from every suitable nlhdlr */
1776 for( j = 0; j < nenfos; ++j )
1777 {
1778 SCIP_Bool addedbranchscores2j;
1779 SCIP_NLHDLR* nlhdlr2;
1780 SCIP_NLHDLREXPRDATA* nlhdlr2exprdata;
1781 SCIP_Real nlhdlr2auxvalue;
1782 SCIP_Bool success2;
1783
1784 SCIPgetExprEnfoDataNonlinear(expr, enfoposs[j], &nlhdlr2, &nlhdlr2exprdata, NULL, NULL, NULL, &nlhdlr2auxvalue);
1785 assert(SCIPnlhdlrHasEstimate(nlhdlr2) && nlhdlr2 != nlhdlr);
1786
1787 SCIPdebugMsg(scip, "asking nonlinear handler %s to %sestimate\n", SCIPnlhdlrGetName(nlhdlr2), overestimate ? "over" : "under");
1788
1789 /* ask the nonlinear handler for an estimator */
1790 if( adjrefpoint )
1791 {
1792 SCIP_CALL( SCIPnlhdlrEvalaux(scip, nlhdlr2, expr, nlhdlr2exprdata, &nlhdlr2auxvalue, soladj) );
1793
1794 SCIP_CALL( SCIPnlhdlrEstimate(scip, conshdlr, nlhdlr2, expr,
1795 nlhdlr2exprdata, soladj,
1796 nlhdlr2auxvalue, overestimate, SCIPgetSolVal(scip, solcopy, auxvar),
1797 addbranchscores, rowpreps2, &success2, &addedbranchscores2j) );
1798 }
1799 else
1800 {
1801 SCIP_CALL( SCIPnlhdlrEstimate(scip, conshdlr, nlhdlr2, expr,
1802 nlhdlr2exprdata, solcopy,
1803 nlhdlr2auxvalue, overestimate, SCIPgetSolVal(scip, solcopy, auxvar),
1804 addbranchscores, rowpreps2, &success2, &addedbranchscores2j) );
1805 }
1806
1807 minidx = SCIPgetPtrarrayMinIdx(scip, rowpreps2);
1808 maxidx = SCIPgetPtrarrayMaxIdx(scip, rowpreps2);
1809
1810 assert((success2 && minidx <= maxidx) || (!success2 && minidx > maxidx));
1811
1812 /* perspectivy all cuts from nlhdlr2 and add them to rowpreps */
1813 for( r = minidx; r <= maxidx; ++r )
1814 {
1815 SCIP_Real maxcoef;
1816 SCIP_Real* rowprepcoefs;
1817 SCIP_VAR** rowprepvars;
1818
1819 rowprep = (SCIP_ROWPREP*) SCIPgetPtrarrayVal(scip, rowpreps2, r);
1820 assert(rowprep != NULL);
1821
1822 #ifdef SCIP_DEBUG
1823 SCIPinfoMessage(scip, NULL, "rowprep for expr ");
1824 SCIPprintExpr(scip, expr, NULL);
1825 SCIPinfoMessage(scip, NULL, "rowprep before perspectivy is: \n");
1826 SCIPprintRowprep(scip, rowprep, NULL);
1827 #endif
1828
1829 /* given a rowprep: sum aixi + sum biyi + c, where xi are semicontinuous variables and yi are
1830 * non-semicontinuous variables (which appear in expr linearly, which detect must have ensured),
1831 * perspectivy the semicontinuous part by adding (1-z)(g0 - c - sum aix0i) (the constant is
1832 * treated as belonging to the semicontinuous part)
1833 */
1834
1835 /* we want cst0 = g0 - c - sum aix0i; first add g0 - c */
1836 cst0 = nlhdlrexprdata->exprvals0[i] + SCIProwprepGetSide(rowprep);
1837
1838 maxcoef = 0.0;
1839 rowprepcoefs = SCIProwprepGetCoefs(rowprep);
1840 rowprepvars = SCIProwprepGetVars(rowprep);
1841
1842 for( v = 0; v < SCIProwprepGetNVars(rowprep); ++v )
1843 {
1844 if( REALABS( rowprepcoefs[v]) > maxcoef )
1845 {
1846 maxcoef = REALABS(rowprepcoefs[v]);
1847 }
1848
1849 scvdata = getSCVarDataInd(nlhdlrdata->scvars, rowprepvars[v], indicator, &pos);
1850
1851 /* a non-semicontinuous variable must be linear in expr; skip it */
1852 if( scvdata == NULL )
1853 continue;
1854
1855 cst0 -= rowprepcoefs[v] * scvdata->vals0[pos];
1856 }
1857
1858 /* only perspectivy when the absolute value of cst0 is not too small
1859 * TODO on ex1252a there was cst0=0 - ok to still use the cut?
1860 */
1861 if( cst0 == 0.0 || maxcoef / REALABS(cst0) <= 10.0 / SCIPfeastol(scip) )
1862 {
1863 /* update the rowprep by adding cst0 - cst0*z */
1864 SCIProwprepAddConstant(rowprep, cst0);
1865 SCIP_CALL(SCIPaddRowprepTerm(scip, rowprep, indicator, -cst0));
1866 }
1867 else
1868 {
1869 SCIPfreeRowprep(scip, &rowprep);
1870 continue;
1871 }
1872
1873 SCIP_CALL(SCIPaddRowprepTerm(scip, rowprep, auxvar, -1.0));
1874
1875 SCIPdebugMsg(scip, "rowprep after perspectivy is: \n");
1876 #ifdef SCIP_DEBUG
1877 SCIPprintRowprep(scip, rowprep, NULL);
1878 #endif
1879
1880 SCIP_CALL( SCIPsetPtrarrayVal(scip, rowpreps, nrowpreps, rowprep) );
1881 SCIP_CALL( SCIPsetBoolarrayVal(scip, addedbranchscores2, nrowpreps, addedbranchscores2j) );
1882 ++nrowpreps;
1883 }
1884
1885 SCIP_CALL( SCIPclearPtrarray(scip, rowpreps2) );
1886 }
1887
1888 if( adjrefpoint )
1889 {
1890 SCIP_CALL( SCIPfreeSol(scip, &soladj) );
1891 }
1892
1893 if( doprobingind )
1894 {
1895 SCIP_CALL( SCIPendProbing(scip) );
1896 }
1897
1898 /* add all cuts found for indicator i */
1899 for( r = SCIPgetPtrarrayMinIdx(scip, rowpreps); r <= SCIPgetPtrarrayMaxIdx(scip, rowpreps) && !stop; ++r )
1900 {
1901 SCIP_RESULT resultr;
1902
1903 #ifdef SCIP_DEBUG
1904 SCIPprintRowprep(scip, rowprep, NULL);
1905 #endif
1906 rowprep = (SCIP_ROWPREP*) SCIPgetPtrarrayVal(scip, rowpreps, r);
1907 resultr = SCIP_DIDNOTFIND;
1908
1909 (void) strcat(SCIProwprepGetName(rowprep), "_persp_indicator_");
1910 (void) strcat(SCIProwprepGetName(rowprep), SCIPvarGetName(indicator));
1911
1912 SCIP_CALL( SCIPprocessRowprepNonlinear(scip, nlhdlr, cons, expr, rowprep, overestimate, auxvar, auxvalue,
1913 allowweakcuts, SCIPgetBoolarrayVal(scip, addedbranchscores2, r), addbranchscores, solcopy, &resultr) );
1914
1915 if( resultr == SCIP_SEPARATED )
1916 *result = SCIP_SEPARATED;
1917 else if( resultr == SCIP_CUTOFF )
1918 {
1919 *result = SCIP_CUTOFF;
1920 stop = TRUE;
1921 }
1922 else if( resultr == SCIP_BRANCHED )
1923 {
1924 if( *result != SCIP_SEPARATED && *result != SCIP_REDUCEDDOM )
1925 *result = SCIP_BRANCHED;
1926 }
1927 else if( resultr != SCIP_DIDNOTFIND )
1928 {
1929 SCIPerrorMessage("estimate called by perspective nonlinear handler returned invalid result <%d>\n", resultr);
1930 return SCIP_INVALIDRESULT;
1931 }
1932 }
1933
1934 /* free all rowpreps for indicator i */
1935 for( r = SCIPgetPtrarrayMinIdx(scip, rowpreps); r <= SCIPgetPtrarrayMaxIdx(scip, rowpreps); ++r )
1936 {
1937 rowprep = (SCIP_ROWPREP*) SCIPgetPtrarrayVal(scip, rowpreps, r);
1938 SCIPfreeRowprep(scip, &rowprep);
1939 }
1940
1941 SCIPfreeBufferArrayNull(scip, &probingvars);
1942 SCIPfreeBufferArrayNull(scip, &probingdoms);
1943 SCIP_CALL( SCIPclearPtrarray(scip, rowpreps) );
1944 }
1945
1946 TERMINATE:
1947 SCIP_CALL( SCIPfreeBoolarray(scip, &addedbranchscores2) );
1948 SCIP_CALL( SCIPfreePtrarray(scip, &rowpreps) );
1949 SCIP_CALL( SCIPfreePtrarray(scip, &rowpreps2) );
1950 if( solcopy != sol )
1951 {
1952 SCIP_CALL( SCIPfreeSol(scip, &solcopy) );
1953 }
1954 SCIPfreeBufferArray(scip, &enfoposs);
1955
1956 return SCIP_OKAY;
1957 }
1958
1959
1960 /*
1961 * nonlinear handler specific interface methods
1962 */
1963
1964 /** includes perspective nonlinear handler in nonlinear constraint handler */
1965 SCIP_RETCODE SCIPincludeNlhdlrPerspective(
1966 SCIP* scip /**< SCIP data structure */
1967 )
1968 {
1969 SCIP_NLHDLRDATA* nlhdlrdata;
1970 SCIP_NLHDLR* nlhdlr;
1971
1972 assert(scip != NULL);
1973
1974 /* create nonlinear handler data */
1975 SCIP_CALL( SCIPallocBlockMemory(scip, &nlhdlrdata) );
1976 BMSclearMemory(nlhdlrdata);
1977
1978 SCIP_CALL( SCIPincludeNlhdlrNonlinear(scip, &nlhdlr, NLHDLR_NAME, NLHDLR_DESC, NLHDLR_DETECTPRIORITY,
1979 NLHDLR_ENFOPRIORITY, nlhdlrDetectPerspective, nlhdlrEvalauxPerspective, nlhdlrdata) );
1980 assert(nlhdlr != NULL);
1981
1982 SCIP_CALL( SCIPaddIntParam(scip, "nlhdlr/" NLHDLR_NAME "/maxproprounds",
1983 "maximal number of propagation rounds in probing",
1984 &nlhdlrdata->maxproprounds, FALSE, DEFAULT_MAXPROPROUNDS, -1, INT_MAX, NULL, NULL) );
1985
1986 SCIP_CALL( SCIPaddRealParam(scip, "nlhdlr/" NLHDLR_NAME "/mindomreduction",
1987 "minimal relative reduction in a variable's domain for applying probing",
1988 &nlhdlrdata->mindomreduction, FALSE, DEFAULT_MINDOMREDUCTION, 0.0, 1.0, NULL, NULL) );
1989
1990 SCIP_CALL( SCIPaddRealParam(scip, "nlhdlr/" NLHDLR_NAME "/minviolprobing",
1991 "minimal violation w.r.t. auxiliary variables for applying probing",
1992 &nlhdlrdata->minviolprobing, FALSE, DEFAULT_MINVIOLPROBING, 0.0, SCIP_REAL_MAX, NULL, NULL) );
1993
1994 SCIP_CALL( SCIPaddBoolParam(scip, "nlhdlr/" NLHDLR_NAME "/probingonlyinsepa",
1995 "whether to do probing only in separation",
1996 &nlhdlrdata->probingonlyinsepa, FALSE, DEFAULT_PROBINGONLYINSEPA, NULL, NULL) );
1997
1998 SCIP_CALL( SCIPaddIntParam(scip, "nlhdlr/" NLHDLR_NAME "/probingfreq",
1999 "probing frequency (-1 - no probing, 0 - root node only)",
2000 &nlhdlrdata->probingfreq, FALSE, DEFAULT_PROBINGFREQ, -1, INT_MAX, NULL, NULL) );
2001
2002 SCIP_CALL( SCIPaddBoolParam(scip, "nlhdlr/" NLHDLR_NAME "/convexonly",
2003 "whether perspective cuts are added only for convex expressions",
2004 &nlhdlrdata->convexonly, FALSE, DEFAULT_CONVEXONLY, NULL, NULL) );
2005
2006 SCIP_CALL( SCIPaddBoolParam(scip, "nlhdlr/" NLHDLR_NAME "/tightenbounds",
2007 "whether variable semicontinuity is used to tighten variable bounds",
2008 &nlhdlrdata->tightenbounds, FALSE, DEFAULT_TIGHTENBOUNDS, NULL, NULL) );
2009
2010 SCIP_CALL( SCIPaddBoolParam(scip, "nlhdlr/" NLHDLR_NAME "/adjrefpoint",
2011 "whether to adjust the reference point",
2012 &nlhdlrdata->adjrefpoint, FALSE, DEFAULT_ADJREFPOINT, NULL, NULL) );
2013
2014 SCIPnlhdlrSetCopyHdlr(nlhdlr, nlhdlrCopyhdlrPerspective);
2015 SCIPnlhdlrSetFreeHdlrData(nlhdlr, nlhdlrFreehdlrdataPerspective);
2016 SCIPnlhdlrSetFreeExprData(nlhdlr, nlhdlrFreeExprDataPerspective);
2017 SCIPnlhdlrSetInitExit(nlhdlr, NULL, nlhdlrExitPerspective);
2018 SCIPnlhdlrSetSepa(nlhdlr, nlhdlrInitSepaPerspective, nlhdlrEnfoPerspective, NULL, NULL);
2019
2020 return SCIP_OKAY;
2021 }
2022