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 visit scip.zib.de. */
13 /* */
14 /* * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * */
15
16 /**@file nlhdlr_convex.c
17 * @ingroup DEFPLUGINS_NLHDLR
18 * @brief nonlinear handlers for convex and concave expressions
19 * @author Benjamin Mueller
20 * @author Stefan Vigerske
21 *
22 * TODO convex: perturb reference point if separation fails due to too large numbers
23 */
24
25 #include <string.h>
26
27 #include "scip/nlhdlr_convex.h"
28 #include "scip/pub_nlhdlr.h"
29 #include "scip/scip_expr.h"
30 #include "scip/cons_nonlinear.h"
31 #include "scip/expr_var.h"
32 #include "scip/pub_misc_rowprep.h"
33 #include "scip/dbldblarith.h"
34
35 /* fundamental nonlinear handler properties */
36 #define CONVEX_NLHDLR_NAME "convex"
37 #define CONVEX_NLHDLR_DESC "handler that identifies and estimates convex expressions"
38 #define CONVEX_NLHDLR_DETECTPRIORITY 50
39 #define CONVEX_NLHDLR_ENFOPRIORITY 50
40
41 #define CONCAVE_NLHDLR_NAME "concave"
42 #define CONCAVE_NLHDLR_DESC "handler that identifies and estimates concave expressions"
43 #define CONCAVE_NLHDLR_DETECTPRIORITY 40
44 #define CONCAVE_NLHDLR_ENFOPRIORITY 40
45
46 #define DEFAULT_DETECTSUM FALSE
47 #define DEFAULT_EXTENDEDFORM TRUE
48 #define DEFAULT_CVXQUADRATIC_CONVEX TRUE
49 #define DEFAULT_CVXQUADRATIC_CONCAVE FALSE
50 #define DEFAULT_CVXSIGNOMIAL TRUE
51 #define DEFAULT_CVXPRODCOMP TRUE
52 #define DEFAULT_HANDLETRIVIAL FALSE
53
54 #define INITLPMAXVARVAL 1000.0 /**< maximal absolute value of variable for still generating a linearization cut at that point in initlp */
55
56 /*lint -e440*/
57 /*lint -e441*/
58 /*lint -e666*/
59 /*lint -e777*/
60
61 /*
62 * Data structures
63 */
64
65 /** nonlinear handler expression data */
66 struct SCIP_NlhdlrExprData
67 {
68 SCIP_EXPR* nlexpr; /**< expression (copy) for which this nlhdlr estimates */
69 SCIP_HASHMAP* nlexpr2origexpr; /**< mapping of our copied expression to original expression */
70
71 int nleafs; /**< number of distinct leafs of nlexpr, i.e., number of distinct (auxiliary) variables handled */
72 SCIP_EXPR** leafexprs; /**< distinct leaf expressions (excluding value-expressions), thus variables */
73 };
74
75 /** nonlinear handler data */
76 struct SCIP_NlhdlrData
77 {
78 SCIP_Bool isnlhdlrconvex; /**< whether this data is used for the convex nlhdlr (TRUE) or the concave one (FALSE) */
79 SCIP_SOL* evalsol; /**< solution used for evaluating expression in a different point,
80 e.g., for facet computation of vertex-polyhedral function */
81
82 /* parameters */
83 SCIP_Bool detectsum; /**< whether to run detection when the root of an expression is a non-quadratic sum */
84 SCIP_Bool extendedform; /**< whether to create extended formulations instead of looking for maximal possible subexpression */
85
86 /* advanced parameters (maybe remove some day) */
87 SCIP_Bool cvxquadratic; /**< whether to use convexity check on quadratics */
88 SCIP_Bool cvxsignomial; /**< whether to use convexity check on signomials */
89 SCIP_Bool cvxprodcomp; /**< whether to use convexity check on product composition f(h)*h */
90 SCIP_Bool handletrivial; /**< whether to handle trivial expressions, i.e., those where all children are variables */
91 };
92
93 /** data struct to be be passed on to vertexpoly-evalfunction (see SCIPcomputeFacetVertexPolyhedralNonlinear) */
94 typedef struct
95 {
96 SCIP_NLHDLREXPRDATA* nlhdlrexprdata;
97 SCIP_SOL* evalsol;
98 SCIP* scip;
99 } VERTEXPOLYFUN_EVALDATA;
100
101 /** stack used in constructExpr to store expressions that need to be investigated ("to do list") */
102 typedef struct
103 {
104 SCIP_EXPR** stack; /**< stack elements */
105 int stacksize; /**< allocated space (in number of pointers) */
106 int stackpos; /**< position of top element of stack */
107 } EXPRSTACK;
108
109 #define DECL_CURVCHECK(x) SCIP_RETCODE x( \
110 SCIP* scip, /**< SCIP data structure */ \
111 SCIP_EXPR* nlexpr, /**< nlhdlr-expr to check */ \
112 SCIP_Bool isrootexpr, /**< whether nlexpr is the root from where detection has been started */ \
113 EXPRSTACK* stack, /**< stack where to add generated leafs */ \
114 SCIP_HASHMAP* nlexpr2origexpr, /**< mapping from our expression copy to original expression */ \
115 SCIP_NLHDLRDATA* nlhdlrdata, /**< data of nlhdlr */ \
116 SCIP_HASHMAP* assumevarfixed, /**< hashmap containing variables that should be assumed to be fixed, or NULL */ \
117 SCIP_Bool* success /**< whether we found something */ \
118 )
119
120 /*
121 * static methods
122 */
123
124 /** create nlhdlr-expression
125 *
126 * does not create children, i.e., assumes that this will be a leaf
127 */
128 static
129 SCIP_RETCODE nlhdlrExprCreate(
130 SCIP* scip, /**< SCIP data structure */
131 SCIP_HASHMAP* nlexpr2origexpr, /**< mapping from copied to original expression */
132 SCIP_EXPR** nlhdlrexpr, /**< buffer to store created expr */
133 SCIP_EXPR* origexpr, /**< original expression to be copied */
134 SCIP_EXPRCURV curv /**< curvature to achieve */
135 )
136 {
137 assert(scip != NULL);
138 assert(nlexpr2origexpr != NULL);
139 assert(nlhdlrexpr != NULL);
140 assert(origexpr != NULL);
141
142 if( SCIPexprGetNChildren(origexpr) == 0 )
143 {
144 /* for leaves, do not copy */
145 *nlhdlrexpr = origexpr;
146 SCIPcaptureExpr(*nlhdlrexpr);
147 if( !SCIPhashmapExists(nlexpr2origexpr, (void*)*nlhdlrexpr) )
148 {
149 SCIP_CALL( SCIPhashmapInsert(nlexpr2origexpr, (void*)*nlhdlrexpr, (void*)origexpr) );
150 }
151 return SCIP_OKAY;
152 }
153
154 /* create copy of expression, but without children */
155 SCIP_CALL( SCIPduplicateExprShallow(scip, origexpr, nlhdlrexpr, NULL, NULL) );
156 assert(*nlhdlrexpr != NULL); /* copies within the same SCIP must always work */
157
158 /* store the curvature we want to get in the curvature flag of the copied expression
159 * it's a bit of a misuse, but once we are done with everything, this is actually correct
160 */
161 SCIPexprSetCurvature(*nlhdlrexpr, curv);
162
163 /* remember which the original expression was */
164 SCIP_CALL( SCIPhashmapInsert(nlexpr2origexpr, (void*)*nlhdlrexpr, (void*)origexpr) );
165
166 return SCIP_OKAY;
167 }
168
169 /** expand nlhdlr-expression by adding children according to original expression */
170 static
171 SCIP_RETCODE nlhdlrExprGrowChildren(
172 SCIP* scip, /**< SCIP data structure */
173 SCIP_HASHMAP* nlexpr2origexpr, /**< mapping from copied to original expression */
174 SCIP_EXPR* nlhdlrexpr, /**< expression for which to create children */
175 SCIP_EXPRCURV* childrencurv /**< curvature required for children, or NULL if to set to UNKNOWN */
176 )
177 {
178 SCIP_EXPR* origexpr;
179 SCIP_EXPR* child;
180 int nchildren;
181 int i;
182
183 assert(scip != NULL);
184 assert(nlhdlrexpr != NULL);
185 assert(SCIPexprGetNChildren(nlhdlrexpr) == 0);
186
187 origexpr = (SCIP_EXPR*)SCIPhashmapGetImage(nlexpr2origexpr, (void*)nlhdlrexpr);
188
189 nchildren = SCIPexprGetNChildren(origexpr);
190 if( nchildren == 0 )
191 return SCIP_OKAY;
192
193 for( i = 0; i < nchildren; ++i )
194 {
195 SCIP_CALL( nlhdlrExprCreate(scip, nlexpr2origexpr, &child, SCIPexprGetChildren(origexpr)[i],
196 childrencurv != NULL ? childrencurv[i] : SCIP_EXPRCURV_UNKNOWN) );
197 SCIP_CALL( SCIPappendExprChild(scip, nlhdlrexpr, child) );
198 /* append captures child, so we can release the capture from nlhdlrExprCreate */
199 SCIP_CALL( SCIPreleaseExpr(scip, &child) );
200 }
201
202 assert(SCIPexprGetNChildren(nlhdlrexpr) == SCIPexprGetNChildren(origexpr));
203
204 return SCIP_OKAY;
205 }
206
207 /** evaluate expression at solution w.r.t. auxiliary variables */
208 static
209 SCIP_DECL_VERTEXPOLYFUN(nlhdlrExprEvalConcave)
210 {
211 VERTEXPOLYFUN_EVALDATA* evaldata = (VERTEXPOLYFUN_EVALDATA*)funcdata;
212 int i;
213
214 assert(args != NULL);
215 assert(nargs == evaldata->nlhdlrexprdata->nleafs);
216 assert(evaldata != NULL);
217
218 #ifdef SCIP_MORE_DEBUG
219 SCIPdebugMsg(evaldata->scip, "eval vertexpolyfun at\n");
220 #endif
221 for( i = 0; i < nargs; ++i )
222 {
223 #ifdef SCIP_MORE_DEBUG
224 SCIPdebugMsg(evaldata->scip, " <%s> = %g\n",
225 SCIPvarGetName(SCIPgetVarExprVar(evaldata->nlhdlrexprdata->leafexprs[i])), args[i]);
226 #endif
227 SCIP_CALL_ABORT( SCIPsetSolVal(evaldata->scip, evaldata->evalsol,
228 SCIPgetVarExprVar(evaldata->nlhdlrexprdata->leafexprs[i]), args[i]) );
229 }
230
231 SCIP_CALL_ABORT( SCIPevalExpr(evaldata->scip, evaldata->nlhdlrexprdata->nlexpr, evaldata->evalsol, 0L) );
232
233 return SCIPexprGetEvalValue(evaldata->nlhdlrexprdata->nlexpr);
234 }
235
236 /** initialize expression stack */
237 static
238 SCIP_RETCODE exprstackInit(
239 SCIP* scip, /**< SCIP data structure */
240 EXPRSTACK* exprstack, /**< stack to initialize */
241 int initsize /**< initial size */
242 )
243 {
244 assert(scip != NULL);
245 assert(exprstack != NULL);
246 assert(initsize > 0);
247
248 SCIP_CALL( SCIPallocBufferArray(scip, &exprstack->stack, initsize) );
249 exprstack->stacksize = initsize;
250 exprstack->stackpos = -1;
251
252 return SCIP_OKAY;
253 }
254
255 /** free expression stack */
256 static
257 void exprstackFree(
258 SCIP* scip, /**< SCIP data structure */
259 EXPRSTACK* exprstack /**< free expression stack */
260 )
261 {
262 assert(scip != NULL);
263 assert(exprstack != NULL);
264
265 SCIPfreeBufferArray(scip, &exprstack->stack);
266 }
267
268 /** add expressions to expression stack */
269 static
270 SCIP_RETCODE exprstackPush(
271 SCIP* scip, /**< SCIP data structure */
272 EXPRSTACK* exprstack, /**< expression stack */
273 int nexprs, /**< number of expressions to push */
274 SCIP_EXPR** exprs /**< expressions to push */
275 )
276 {
277 assert(scip != NULL);
278 assert(exprstack != NULL);
279
280 if( nexprs == 0 )
281 return SCIP_OKAY;
282
283 assert(exprs != NULL);
284
285 if( exprstack->stackpos+1 + nexprs > exprstack->stacksize ) /*lint !e644*/
286 {
287 exprstack->stacksize = SCIPcalcMemGrowSize(scip, exprstack->stackpos+1 + nexprs); /*lint !e644*/
288 SCIP_CALL( SCIPreallocBufferArray(scip, &exprstack->stack, exprstack->stacksize) );
289 }
290
291 memcpy(exprstack->stack + (exprstack->stackpos+1), exprs, nexprs * sizeof(SCIP_EXPR*)); /*lint !e679*/ /*lint !e737*/
292 exprstack->stackpos += nexprs;
293
294 return SCIP_OKAY;
295 }
296
297 /** gives expression from top of expression stack and removes it from stack */
298 static
299 SCIP_EXPR* exprstackPop(
300 EXPRSTACK* exprstack /**< expression stack */
301 )
302 {
303 assert(exprstack != NULL);
304 assert(exprstack->stackpos >= 0);
305
306 return exprstack->stack[exprstack->stackpos--];
307 }
308
309 /** indicate whether expression stack is empty */
310 static
311 SCIP_Bool exprstackIsEmpty(
312 EXPRSTACK* exprstack /**< expression stack */
313 )
314 {
315 assert(exprstack != NULL);
316
317 return exprstack->stackpos < 0;
318 }
319
320 /** looks whether given expression is (proper) quadratic and has a given curvature
321 *
322 * If having a given curvature, currently require all arguments of quadratic to be linear.
323 * Hence, not using this for a simple square term, as curvCheckExprhdlr may provide a better condition on argument curvature then.
324 * Also we wouldn't do anything useful for a single bilinear term.
325 * Thus, run on sum's only.
326 */
327 static
328 DECL_CURVCHECK(curvCheckQuadratic)
329 { /*lint --e{715}*/
330 SCIP_EXPR* expr;
331 SCIP_EXPRCURV presentcurv;
332 SCIP_EXPRCURV wantedcurv;
333 SCIP_HASHSET* lonelysquares = NULL;
334 SCIP_Bool isquadratic;
335 int nbilinexprs;
336 int nquadexprs;
337 int i;
338
339 assert(nlexpr != NULL);
340 assert(stack != NULL);
341 assert(nlexpr2origexpr != NULL);
342 assert(success != NULL);
343
344 *success = FALSE;
345
346 if( !nlhdlrdata->cvxquadratic )
347 return SCIP_OKAY;
348
349 if( !SCIPisExprSum(scip, nlexpr) )
350 return SCIP_OKAY;
351
352 wantedcurv = SCIPexprGetCurvature(nlexpr);
353 if( wantedcurv == SCIP_EXPRCURV_LINEAR )
354 return SCIP_OKAY;
355 assert(wantedcurv == SCIP_EXPRCURV_CONVEX || wantedcurv == SCIP_EXPRCURV_CONCAVE);
356
357 expr = (SCIP_EXPR*)SCIPhashmapGetImage(nlexpr2origexpr, (void*)nlexpr);
358 assert(expr != NULL);
359
360 /* check whether quadratic */
361 SCIP_CALL( SCIPcheckExprQuadratic(scip, expr, &isquadratic) );
362
363 /* if not quadratic, then give up here */
364 if( !isquadratic )
365 return SCIP_OKAY;
366
367 SCIPexprGetQuadraticData(expr, NULL, NULL, NULL, NULL, &nquadexprs, &nbilinexprs, NULL, NULL);
368
369 /* if only single square term (+linear), then give up here (let curvCheckExprhdlr handle this) */
370 if( nquadexprs <= 1 )
371 return SCIP_OKAY;
372
373 /* if root expression is only sum of squares (+linear) and detectsum is disabled, then give up here, too */
374 if( isrootexpr && !nlhdlrdata->detectsum && nbilinexprs == 0 )
375 return SCIP_OKAY;
376
377 /* get curvature of quadratic
378 * TODO as we know what curvature we want, we could first do some simple checks like computing xQx for a random x
379 */
380 SCIP_CALL( SCIPcomputeExprQuadraticCurvature(scip, expr, &presentcurv, assumevarfixed, FALSE) );
381
382 /* if not having desired curvature, return */
383 if( presentcurv != wantedcurv )
384 return SCIP_OKAY;
385
386 *success = TRUE;
387
388 if( !nlhdlrdata->detectsum )
389 {
390 /* first step towards block-decomposition of quadratic term:
391 * collect all square-expressions (in original expr) which have no adjacent bilinear term
392 * we will treat these x^2 as linear, i.e., add an auxvar for them, so x^2 maybe linearized
393 * more efficiently (in particular if x is discrete)
394 */
395 SCIP_CALL( SCIPhashsetCreate(&lonelysquares, SCIPblkmem(scip), nquadexprs) );
396 for( i = 0; i < nquadexprs; ++i )
397 {
398 int nadjbilin;
399 SCIP_EXPR* sqrexpr;
400
401 SCIPexprGetQuadraticQuadTerm(expr, i, NULL, NULL, NULL, &nadjbilin, NULL, &sqrexpr);
402 if( nadjbilin == 0 )
403 {
404 assert(sqrexpr != NULL);
405 SCIP_CALL( SCIPhashsetInsert(lonelysquares, SCIPblkmem(scip), (void*)sqrexpr) );
406 }
407 }
408 }
409
410 /* add immediate children to nlexpr */
411 SCIP_CALL( nlhdlrExprGrowChildren(scip, nlexpr2origexpr, nlexpr, NULL) );
412 assert(SCIPexprGetNChildren(nlexpr) == SCIPexprGetNChildren(expr));
413
414 /* put children that are not square or product on stack
415 * grow child for children that are square or product and put this child on stack
416 * require all children to be linear
417 */
418 for( i = 0; i < SCIPexprGetNChildren(nlexpr); ++i )
419 {
420 SCIP_EXPR* child;
421 SCIP_EXPRCURV curvlinear[2] = { SCIP_EXPRCURV_LINEAR, SCIP_EXPRCURV_LINEAR };
422
423 child = SCIPexprGetChildren(nlexpr)[i];
424 assert(child != NULL);
425
426 assert(SCIPhashmapGetImage(nlexpr2origexpr, (void*)child) == SCIPexprGetChildren(expr)[i]);
427
428 if( SCIPisExprPower(scip, child) && SCIPgetExponentExprPow(child) == 2.0 &&
429 (lonelysquares == NULL || !SCIPhashsetExists(lonelysquares, SCIPexprGetChildren(expr)[i])) )
430 {
431 /* square term that isn't lonely, i.e., orig-version of child is a square-expr and nadjbilin>0 */
432 SCIP_CALL( nlhdlrExprGrowChildren(scip, nlexpr2origexpr, child, curvlinear) );
433 assert(SCIPexprGetNChildren(child) == 1);
434 SCIP_CALL( exprstackPush(scip, stack, 1, SCIPexprGetChildren(child)) );
435 }
436 else if( SCIPisExprProduct(scip, child) && SCIPexprGetNChildren(SCIPexprGetChildren(expr)[i]) == 2 )
437 /* using original version of child here as NChildren(child)==0 atm */
438 {
439 /* bilinear term */
440 SCIP_CALL( nlhdlrExprGrowChildren(scip, nlexpr2origexpr, child, curvlinear) );
441 assert(SCIPexprGetNChildren(child) == 2);
442 SCIP_CALL( exprstackPush(scip, stack, 2, SCIPexprGetChildren(child)) );
443 }
444 else
445 {
446 /* linear term (or term to be considered as linear) or lonely square term
447 * if we want extended formulations, then require linearity, so an auxvar will be introduced if it is nonlinear
448 * if we do not want extended formulations, then the term needs to have curvature "wantedcurv"
449 * thus, if the coef is negative, then the child needs to have the curvature opposite to "wantedcurv"
450 */
451 if( nlhdlrdata->extendedform )
452 SCIPexprSetCurvature(child, SCIP_EXPRCURV_LINEAR);
453 else
454 SCIPexprSetCurvature(child, SCIPexprcurvMultiply(SCIPgetCoefsExprSum(nlexpr)[i], wantedcurv));
455 SCIP_CALL( exprstackPush(scip, stack, 1, &child) );
456 }
457 }
458
459 if( lonelysquares != NULL )
460 SCIPhashsetFree(&lonelysquares, SCIPblkmem(scip));
461
462 return SCIP_OKAY;
463 }
464
465 /** looks whether top of given expression looks like a signomial that can have a given curvature
466 *
467 * e.g., sqrt(x)*sqrt(y) is convex if x,y >= 0 and x and y are convex
468 *
469 * unfortunately, doesn't work for tls, because i) it's originally sqrt(x*y), and ii) it is expanded into some sqrt(z*y+y);
470 * but works for cvxnonsep_nsig
471 */
472 static
473 DECL_CURVCHECK(curvCheckSignomial)
474 { /*lint --e{715}*/
475 SCIP_EXPR* expr;
476 SCIP_EXPR* child;
477 SCIP_Real* exponents;
478 SCIP_INTERVAL* bounds;
479 SCIP_EXPRCURV* curv;
480 int nfactors;
481 int i;
482
483 assert(nlexpr != NULL);
484 assert(stack != NULL);
485 assert(nlexpr2origexpr != NULL);
486 assert(success != NULL);
487
488 *success = FALSE;
489
490 if( !nlhdlrdata->cvxsignomial )
491 return SCIP_OKAY;
492
493 if( !SCIPisExprProduct(scip, nlexpr) )
494 return SCIP_OKAY;
495
496 expr = (SCIP_EXPR*)SCIPhashmapGetImage(nlexpr2origexpr, (void*)nlexpr);
497 assert(expr != NULL);
498
499 nfactors = SCIPexprGetNChildren(expr);
500 if( nfactors <= 1 ) /* boooring */
501 return SCIP_OKAY;
502
503 SCIP_CALL( SCIPallocBufferArray(scip, &exponents, nfactors) );
504 SCIP_CALL( SCIPallocBufferArray(scip, &bounds, nfactors) );
505 SCIP_CALL( SCIPallocBufferArray(scip, &curv, nfactors) );
506
507 for( i = 0; i < nfactors; ++i )
508 {
509 child = SCIPexprGetChildren(expr)[i];
510 assert(child != NULL);
511
512 if( !SCIPisExprPower(scip, child) )
513 {
514 exponents[i] = 1.0;
515 SCIP_CALL( SCIPevalExprActivity(scip, child) );
516 bounds[i] = SCIPexprGetActivity(child);
517 }
518 else
519 {
520 exponents[i] = SCIPgetExponentExprPow(child);
521 SCIP_CALL( SCIPevalExprActivity(scip, SCIPexprGetChildren(child)[0]) );
522 bounds[i] = SCIPexprGetActivity(SCIPexprGetChildren(child)[0]);
523 }
524 }
525
526 if( !SCIPexprcurvMonomialInv(SCIPexprcurvMultiply(SCIPgetCoefExprProduct(expr), SCIPexprGetCurvature(nlexpr)),
527 nfactors, exponents, bounds, curv) )
528 goto TERMINATE;
529
530 /* add immediate children to nlexpr
531 * some entries in curv actually apply to arguments of pow's, will correct this next
532 */
533 SCIP_CALL( nlhdlrExprGrowChildren(scip, nlexpr2origexpr, nlexpr, curv) );
534 assert(SCIPexprGetNChildren(nlexpr) == nfactors);
535
536 /* put children that are not power on stack
537 * grow child for children that are power and put this child on stack
538 * if extendedform, then require children to be linear
539 * unless they are linear, an auxvar will be introduced for them and thus they will be handled as var here
540 */
541 for( i = 0; i < nfactors; ++i )
542 {
543 child = SCIPexprGetChildren(nlexpr)[i];
544 assert(child != NULL);
545
546 if( SCIPisExprPower(scip, child) )
547 {
548 SCIP_CALL( nlhdlrExprGrowChildren(scip, nlexpr2origexpr, child, &curv[i]) );
549 assert(SCIPexprGetNChildren(child) == 1);
550 child = SCIPexprGetChildren(child)[0];
551 }
552 assert(SCIPexprGetNChildren(child) == 0);
553
554 if( nlhdlrdata->extendedform )
555 {
556 SCIPexprSetCurvature(child, SCIP_EXPRCURV_LINEAR);
557 #ifdef SCIP_DEBUG
558 SCIPinfoMessage(scip, NULL, "Extendedform: Require linearity for ");
559 SCIPprintExpr(scip, child, NULL);
560 SCIPinfoMessage(scip, NULL, "\n");
561 #endif
562 }
563
564 SCIP_CALL( exprstackPush(scip, stack, 1, &child) );
565 }
566
567 *success = TRUE;
568
569 TERMINATE:
570 SCIPfreeBufferArray(scip, &curv);
571 SCIPfreeBufferArray(scip, &bounds);
572 SCIPfreeBufferArray(scip, &exponents);
573
574 return SCIP_OKAY;
575 }
576
577 /** looks for \f$f(c h(x)+d) h(x) \cdot \text{constant}\f$ and tries to conclude conditions on curvature
578 *
579 * Assume \f$h\f$ is univariate:
580 * - First derivative is \f$f'(c h + d) c h' h + f(c h + d) h'\f$.
581 * - Second derivative is \f{align}{&f''(c h + d) c h' c h' h + f'(c h + d) (c h'' h + c h' h') + f'(c h + d) c h' h' + f(c h + d) h'' \\
582 * =& f''(c h + d) c^2 h'^2 h + f'(c h + d) c h'' h + 2 f'(c h + d) c h'^2 + f(c h + d) h''.\f}
583 * Remove always positive factors leaves \f[f''(c h + d) h,\quad f'(c h + d) c h'' h,\quad f'(c h + d) c,\quad f(c h + d) h''.\f]
584 * For convexity we want all these terms to be nonnegative. For concavity we want all of them to be nonpositive.
585 * Note, that in each term either both \f$f'(c h + d)\f$ and \f$c\f$ occur, or none of them.
586 * - Thus, \f$f(c h(x) + d)h(x)\f$ is convex if \f$cf\f$ is monotonically increasing \f$(c f' \geq 0)\f$ and either
587 * - \f$f\f$ is convex \f$(f'' \geq 0)\f$ and \f$h\f$ is nonnegative \f$(h \geq 0)\f$ and \f$h\f$ is convex \f$(h'' \geq 0)\f$ and [\f$f\f$ is nonnegative \f$(f \geq 0)\f$ or \f$h\f$ is linear \f$(h''=0)\f$], or
588 * - \f$f\f$ is concave \f$(f'' \leq 0)\f$ and \f$h\f$ is nonpositive \f$(h \leq 0)\f$ and \f$h\f$ is concave \f$(h'' \leq 0)\f$ and [\f$f\f$ is nonpositive \f$(f \leq 0)\f$ or \f$h\f$ is linear \f$(h''=0)\f$].
589 * - Further, \f$f(c h(x) + d)h(x)\f$ is concave if \f$cf\f$ is monotonically decreasing \f$(c f' \leq 0)\f$ and either
590 * - f is convex \f$(f'' \geq 0)\f$ and \f$h\f$ is nonpositive \f$(h \leq 0)\f$ and \f$h\f$ is concave \f$(h'' \leq 0)\f$ and [\f$f\f$ is nonnegative \f$(f \geq 0)\f$ or \f$h\f$ is linear \f$(h''=0)\f$], or
591 * - f is concave \f$(f'' \leq 0)\f$ and \f$h\f$ is nonnegative \f$(h >= 0)\f$ and \f$h\f$ is convex \f$(h'' \geq 0)\f$ and [\f$f\f$ is nonpositive \f$(f \leq 0)\f$ or \f$h\f$ is linear \f$(h''=0)\f$].
592 *
593 * This should hold also for multivariate and linear \f$h\f$, as things are invariant under linear transformations.
594 * Similar to signomial, I'll assume that this will also hold for other multivariate \f$h\f$ (someone has a formal proof?).
595 */
596 static
597 DECL_CURVCHECK(curvCheckProductComposite)
598 { /*lint --e{715}*/
599 SCIP_EXPR* expr;
600 SCIP_EXPR* f;
601 SCIP_EXPR* h = NULL;
602 SCIP_Real c = 0.0;
603 SCIP_EXPR* ch = NULL; /* c * h */
604 SCIP_INTERVAL fbounds;
605 SCIP_INTERVAL hbounds;
606 SCIP_MONOTONE fmonotonicity;
607 SCIP_EXPRCURV desiredcurv;
608 SCIP_EXPRCURV hcurv;
609 SCIP_EXPRCURV dummy;
610 int fidx;
611
612 assert(nlexpr != NULL);
613 assert(stack != NULL);
614 assert(nlexpr2origexpr != NULL);
615 assert(success != NULL);
616
617 *success = FALSE;
618
619 if( !nlhdlrdata->cvxprodcomp )
620 return SCIP_OKAY;
621
622 if( !SCIPisExprProduct(scip, nlexpr) )
623 return SCIP_OKAY;
624
625 expr = (SCIP_EXPR*)SCIPhashmapGetImage(nlexpr2origexpr, (void*)nlexpr);
626 assert(expr != NULL);
627
628 if( SCIPexprGetNChildren(expr) != 2 )
629 return SCIP_OKAY;
630
631 /* check whether we have f(c * h(x)) * h(x) or h(x) * f(c * h(x)) */
632 for( fidx = 0; fidx <= 1; ++fidx )
633 {
634 f = SCIPexprGetChildren(expr)[fidx];
635
636 if( SCIPexprGetNChildren(f) != 1 )
637 continue;
638
639 ch = SCIPexprGetChildren(f)[0];
640 c = 1.0;
641 h = ch;
642
643 /* check whether ch is of the form c*h(x), then switch h to child ch */
644 if( SCIPisExprSum(scip, ch) && SCIPexprGetNChildren(ch) == 1 )
645 {
646 c = SCIPgetCoefsExprSum(ch)[0];
647 h = SCIPexprGetChildren(ch)[0];
648 assert(c != 1.0 || SCIPgetConstantExprSum(ch) != 0.0); /* we could handle this, but it should have been simplified away */
649 }
650
651 #ifndef NLHDLR_CONVEX_UNITTEST
652 /* can assume that duplicate subexpressions have been identified and comparing pointer is sufficient */
653 if( SCIPexprGetChildren(expr)[1-fidx] == h )
654 #else
655 /* called from unittest -> duplicate subexpressions were not identified -> compare more expensively */
656 if( SCIPcompareExpr(scip, SCIPexprGetChildren(expr)[1-fidx], h) == 0 )
657 #endif
658 break;
659 }
660 if( fidx == 2 )
661 return SCIP_OKAY;
662
663 #ifdef SCIP_MORE_DEBUG
664 SCIPinfoMessage(scip, NULL, "f(c*h+d)*h with f = %s, c = %g, d = %g, h = ", SCIPexprhdlrGetName(SCIPexprGetHdlr(f)),
665 c, h != ch ? SCIPgetConstantExprSum(ch) : 0.0);
666 SCIPprintExpr(scip, h, NULL);
667 SCIPinfoMessage(scip, NULL, "\n");
668 #endif
669
670 assert(c != 0.0);
671
672 SCIP_CALL( SCIPevalExprActivity(scip, f) );
673 SCIP_CALL( SCIPevalExprActivity(scip, h) );
674 fbounds = SCIPexprGetActivity(f);
675 hbounds = SCIPexprGetActivity(h);
676
677 /* if h has mixed sign, then cannot conclude anything */
678 if( hbounds.inf < 0.0 && hbounds.sup > 0.0 )
679 return SCIP_OKAY;
680
681 SCIP_CALL( SCIPcallExprMonotonicity(scip, f, 0, &fmonotonicity) );
682
683 /* if f is not monotone, then cannot conclude anything */
684 if( fmonotonicity == SCIP_MONOTONE_UNKNOWN )
685 return SCIP_OKAY;
686
687 /* curvature we want to achieve (negate if product has negative coef) */
688 desiredcurv = SCIPexprcurvMultiply(SCIPgetCoefExprProduct(nlexpr), SCIPexprGetCurvature(nlexpr));
689
690 /* now check the conditions as stated above */
691 if( desiredcurv == SCIP_EXPRCURV_CONVEX )
692 {
693 /* f(c h(x)+d)h(x) is convex if c*f is monotonically increasing (c f' >= 0) and either
694 * - f is convex (f'' >= 0) and h is nonnegative (h >= 0) and h is convex (h'' >= 0) and [f is nonnegative (f >= 0) or h is linear (h''=0)], or
695 * - f is concave (f'' <= 0) and h is nonpositive (h <= 0) and h is concave (h'' <= 0) and [f is nonpositive (f <= 0) or h is linear (h''=0)]
696 * as the curvature requirements on f are on f only and not the composition f(h), we can ignore the requirements returned by SCIPcallExprCurvature (last arg)
697 */
698 if( (c > 0.0 && fmonotonicity != SCIP_MONOTONE_INC) || (c < 0.0 && fmonotonicity != SCIP_MONOTONE_DEC) )
699 return SCIP_OKAY;
700
701 /* check whether f can be convex (h>=0) or concave (h<=0), resp., and derive requirements for h */
702 if( hbounds.inf >= 0 )
703 {
704 SCIP_CALL( SCIPcallExprCurvature(scip, f, SCIP_EXPRCURV_CONVEX, success, &dummy) );
705
706 /* now h also needs to be convex; and if f < 0, then h actually needs to be linear */
707 if( fbounds.inf < 0.0 )
708 hcurv = SCIP_EXPRCURV_LINEAR;
709 else
710 hcurv = SCIP_EXPRCURV_CONVEX;
711 }
712 else
713 {
714 SCIP_CALL( SCIPcallExprCurvature(scip, f, SCIP_EXPRCURV_CONCAVE, success, &dummy) );
715
716 /* now h also needs to be concave; and if f > 0, then h actually needs to be linear */
717 if( fbounds.sup > 0.0 )
718 hcurv = SCIP_EXPRCURV_LINEAR;
719 else
720 hcurv = SCIP_EXPRCURV_CONCAVE;
721 }
722 }
723 else
724 {
725 /* f(c h(x)+d)*h(x) is concave if c*f is monotonically decreasing (c f' <= 0) and either
726 * - f is convex (f'' >= 0) and h is nonpositive (h <= 0) and h is concave (h'' <= 0) and [f is nonnegative (f >= 0) or h is linear (h''=0)], or
727 * - f is concave (f'' <= 0) and h is nonnegative (h >= 0) and h is convex (h'' >= 0) and [f is nonpositive (f <= 0) or h is linear (h''=0)]
728 * as the curvature requirements on f are on f only and not the composition f(h), we can ignore the requirements returned by SCIPcallExprCurvature (last arg)
729 */
730 if( (c > 0.0 && fmonotonicity != SCIP_MONOTONE_DEC) || (c < 0.0 && fmonotonicity != SCIP_MONOTONE_INC) )
731 return SCIP_OKAY;
732
733 /* check whether f can be convex (h<=0) or concave (h>=0), resp., and derive requirements for h */
734 if( hbounds.sup <= 0 )
735 {
736 SCIP_CALL( SCIPcallExprCurvature(scip, f, SCIP_EXPRCURV_CONVEX, success, &dummy) );
737
738 /* now h also needs to be concave; and if f < 0, then h actually needs to be linear */
739 if( fbounds.inf < 0.0 )
740 hcurv = SCIP_EXPRCURV_LINEAR;
741 else
742 hcurv = SCIP_EXPRCURV_CONCAVE;
743 }
744 else
745 {
746 SCIP_CALL( SCIPcallExprCurvature(scip, f, SCIP_EXPRCURV_CONCAVE, success, &dummy) );
747
748 /* now h also needs to be convex; and if f > 0, then h actually needs to be linear */
749 if( fbounds.sup > 0.0 )
750 hcurv = SCIP_EXPRCURV_LINEAR;
751 else
752 hcurv = SCIP_EXPRCURV_CONVEX;
753 }
754 }
755
756 if( !*success )
757 return SCIP_OKAY;
758
759 /* add immediate children (f and ch) to nlexpr; we set required curvature for h further below */
760 SCIP_CALL( nlhdlrExprGrowChildren(scip, nlexpr2origexpr, nlexpr, NULL) );
761 assert(SCIPexprGetNChildren(nlexpr) == 2);
762
763 /* copy of f (and h) should have same child position in nlexpr as f (and h) has on expr (resp) */
764 assert(SCIPhashmapGetImage(nlexpr2origexpr, (void*)SCIPexprGetChildren(nlexpr)[fidx]) == (void*)f);
765 #ifndef NLHDLR_CONVEX_UNITTEST
766 assert(SCIPhashmapGetImage(nlexpr2origexpr, (void*)SCIPexprGetChildren(nlexpr)[1-fidx]) == (void*)h);
767 #endif
768 /* push this h onto stack for further checking */
769 SCIP_CALL( exprstackPush(scip, stack, 1, &(SCIPexprGetChildren(nlexpr)[1-fidx])) );
770
771 /* if we prefer extended formulations, then we always want h() to be linear */
772 if( nlhdlrdata->extendedform )
773 hcurv = SCIP_EXPRCURV_LINEAR;
774
775 /* h-child of product should have curvature hcurv */
776 SCIPexprSetCurvature(SCIPexprGetChildren(nlexpr)[1-fidx], hcurv);
777
778 if( h != ch )
779 {
780 /* add copy of ch as child to copy of f */
781 SCIP_CALL( nlhdlrExprGrowChildren(scip, nlexpr2origexpr, SCIPexprGetChildren(nlexpr)[fidx], NULL) );
782 assert(SCIPexprGetNChildren(SCIPexprGetChildren(nlexpr)[fidx]) == 1);
783 assert(SCIPhashmapGetImage(nlexpr2origexpr, (void*)SCIPexprGetChildren(SCIPexprGetChildren(nlexpr)[fidx])[0]) == (void*)ch);
784
785 /* add copy of h (created above as child of product) as child in copy of ch */
786 SCIP_CALL( SCIPappendExprChild(scip,
787 SCIPexprGetChildren(SCIPexprGetChildren(nlexpr)[fidx])[0] /* copy of ch */,
788 SCIPexprGetChildren(nlexpr)[1-fidx] /* copy of h */) );
789 }
790 else
791 {
792 /* add copy of h (created above as child of product) as child in copy of f */
793 SCIP_CALL( SCIPappendExprChild(scip,
794 SCIPexprGetChildren(nlexpr)[fidx] /* copy of f */,
795 SCIPexprGetChildren(nlexpr)[1-fidx] /* copy of h */) );
796 }
797
798 return SCIP_OKAY;
799 }
800
801 /** use expression handlers curvature callback to check whether given curvature can be achieved */
802 static
803 DECL_CURVCHECK(curvCheckExprhdlr)
804 { /*lint --e{715}*/
805 SCIP_EXPR* origexpr;
806 int nchildren;
807 SCIP_EXPRCURV* childcurv;
808
809 assert(nlexpr != NULL);
810 assert(stack != NULL);
811 assert(nlexpr2origexpr != NULL);
812 assert(success != NULL);
813
814 origexpr = (SCIP_EXPR*)SCIPhashmapGetImage(nlexpr2origexpr, nlexpr);
815 assert(origexpr != NULL);
816 nchildren = SCIPexprGetNChildren(origexpr);
817
818 if( nchildren == 0 )
819 {
820 /* if originally no children, then should be var or value, which should have every curvature,
821 * so should always be success
822 */
823 SCIP_CALL( SCIPcallExprCurvature(scip, origexpr, SCIPexprGetCurvature(nlexpr), success, NULL) );
824 assert(*success);
825
826 return SCIP_OKAY;
827 }
828
829 /* ignore sums if > 1 children
830 * NOTE: this means that for something like 1+f(x), even if f is a trivial convex expression, we would handle 1+f(x)
831 * with this nlhdlr, instead of formulating this as 1+z and handling z=f(x) with the default nlhdlr, i.e., the exprhdlr
832 * today, I prefer handling this here, as it avoids introducing an extra auxiliary variable
833 */
834 if( isrootexpr && !nlhdlrdata->detectsum && SCIPisExprSum(scip, nlexpr) && nchildren > 1 )
835 return SCIP_OKAY;
836
837 SCIP_CALL( SCIPallocBufferArray(scip, &childcurv, nchildren) );
838
839 /* check whether and under which conditions origexpr can have desired curvature */
840 SCIP_CALL( SCIPcallExprCurvature(scip, origexpr, SCIPexprGetCurvature(nlexpr), success, childcurv) );
841 #ifdef SCIP_MORE_DEBUG
842 SCIPprintExpr(scip, origexpr, NULL);
843 SCIPinfoMessage(scip, NULL, " is %s? %d\n", SCIPexprcurvGetName(SCIPexprGetCurvature(nlexpr)), *success);
844 #endif
845 if( !*success )
846 goto TERMINATE;
847
848 /* if origexpr can have curvature curv, then don't treat it as leaf, but include its children */
849 SCIP_CALL( nlhdlrExprGrowChildren(scip, nlexpr2origexpr, nlexpr, childcurv) );
850 assert(SCIPexprGetChildren(nlexpr) != NULL);
851 assert(SCIPexprGetNChildren(nlexpr) == nchildren);
852
853 /* If we prefer extended formulations, then require all children to be linear.
854 * Unless they are, auxvars will be introduced and they will be handles as variables, which can be an
855 * advantage in the context of extended formulations.
856 */
857 if( nlhdlrdata->extendedform )
858 {
859 int i;
860 for( i = 0; i < nchildren; ++i )
861 SCIPexprSetCurvature(SCIPexprGetChildren(nlexpr)[i], SCIP_EXPRCURV_LINEAR);
862 #ifdef SCIP_DEBUG
863 SCIPinfoMessage(scip, NULL, "require linearity for children of ");
864 SCIPprintExpr(scip, origexpr, NULL);
865 SCIPinfoMessage(scip, NULL, "\n");
866 #endif
867 }
868
869 /* add children expressions to to-do list (stack) */
870 SCIP_CALL( exprstackPush(scip, stack, nchildren, SCIPexprGetChildren(nlexpr)) );
871
872 TERMINATE:
873 SCIPfreeBufferArray(scip, &childcurv);
874
875 return SCIP_OKAY;
876 }
877
878 /** curvature check and expression-growing methods
879 *
880 * some day this could be plugins added by users at runtime, but for now we have a fixed list here
881 * @note curvCheckExprhdlr should be last
882 */
883 static DECL_CURVCHECK((*CURVCHECKS[])) = { curvCheckProductComposite, curvCheckSignomial, curvCheckQuadratic, curvCheckExprhdlr };
884 /** number of curvcheck methods */
885 static const int NCURVCHECKS = sizeof(CURVCHECKS) / sizeof(void*);
886
887 /** checks whether expression is a sum with more than one child and each child being a variable or going to be a variable if `expr` is a nlhdlr-specific copy
888 *
889 * Within constructExpr(), we can have an expression of any type which is a copy of an original expression,
890 * but without children. At the end of constructExpr() (after the loop with the stack), these expressions
891 * will remain as leafs and will eventually be turned into variables in collectLeafs(). Thus, we treat
892 * every child that has no children as if it were a variable. Theoretically, there is still the possibility
893 * that it could be a constant (value-expression), but simplify should have removed these.
894 */
895 static
896 SCIP_Bool exprIsMultivarLinear(
897 SCIP* scip, /**< SCIP data structure */
898 SCIP_EXPR* expr /**< expression to check */
899 )
900 {
901 int nchildren;
902 int c;
903
904 assert(expr != NULL);
905
906 if( !SCIPisExprSum(scip, expr) )
907 return FALSE;
908
909 nchildren = SCIPexprGetNChildren(expr);
910 if( nchildren <= 1 )
911 return FALSE;
912
913 for( c = 0; c < nchildren; ++c )
914 /*if( !SCIPisExprVar(scip, SCIPexprGetChildren(expr)[c]) ) */
915 if( SCIPexprGetNChildren(SCIPexprGetChildren(expr)[c]) > 0 )
916 return FALSE;
917
918 return TRUE;
919 }
920
921 /** constructs a subexpression (as nlhdlr-expression) of maximal size that has a given curvature
922 *
923 * If the curvature cannot be achieved for an expression in the original expression graph,
924 * then this expression becomes a leaf in the nlhdlr-expression.
925 *
926 * Sets `*rootnlexpr` to NULL if failed.
927 */
928 static
929 SCIP_RETCODE constructExpr(
930 SCIP* scip, /**< SCIP data structure */
931 SCIP_NLHDLRDATA* nlhdlrdata, /**< nonlinear handler data */
932 SCIP_EXPR** rootnlexpr, /**< buffer to store created expression */
933 SCIP_HASHMAP* nlexpr2origexpr, /**< mapping from our expression copy to original expression */
934 int* nleafs, /**< number of leafs in constructed expression */
935 SCIP_EXPR* rootexpr, /**< expression */
936 SCIP_EXPRCURV curv, /**< curvature to achieve */
937 SCIP_HASHMAP* assumevarfixed, /**< hashmap containing variables that should be assumed to be fixed, or NULL */
938 SCIP_Bool assumecurvature, /**< whether to assume that desired curvature is given (skips curvature checks) */
939 SCIP_Bool* curvsuccess /**< pointer to store whether the curvature could be achieved
940 w.r.t. the original variables (might be NULL) */
941 )
942 {
943 SCIP_EXPR* nlexpr;
944 EXPRSTACK stack; /* to do list: expressions where to check whether they can have the desired curvature when taking their children into account */
945 int oldstackpos;
946 SCIP_Bool isrootexpr = TRUE;
947
948 assert(scip != NULL);
949 assert(nlhdlrdata != NULL);
950 assert(rootnlexpr != NULL);
951 assert(nlexpr2origexpr != NULL);
952 assert(nleafs != NULL);
953 assert(rootexpr != NULL);
954 assert(curv == SCIP_EXPRCURV_CONVEX || curv == SCIP_EXPRCURV_CONCAVE);
955
956 /* create root expression */
957 SCIP_CALL( nlhdlrExprCreate(scip, nlexpr2origexpr, rootnlexpr, rootexpr, curv) );
958
959 *nleafs = 0;
960 if( curvsuccess != NULL )
961 *curvsuccess = TRUE;
962
963 SCIP_CALL( exprstackInit(scip, &stack, 20) );
964 SCIP_CALL( exprstackPush(scip, &stack, 1, rootnlexpr) );
965 while( !exprstackIsEmpty(&stack) )
966 {
967 /* take expression from stack */
968 nlexpr = exprstackPop(&stack);
969 assert(nlexpr != NULL);
970 assert(SCIPexprGetNChildren(nlexpr) == 0);
971
972 oldstackpos = stack.stackpos;
973 if( nlhdlrdata->isnlhdlrconvex && !SCIPexprhdlrHasBwdiff(SCIPexprGetHdlr(nlexpr)) )
974 {
975 /* if bwdiff is not implemented, then we could not generate cuts in the convex nlhdlr, so "stop" (treat nlexpr as variable) */
976 }
977 else if( !nlhdlrdata->isnlhdlrconvex && exprIsMultivarLinear(scip, (SCIP_EXPR*)SCIPhashmapGetImage(nlexpr2origexpr, (void*)nlexpr)) )
978 {
979 /* if we are in the concave handler, we would like to treat linear multivariate subexpressions by a new auxvar always,
980 * e.g., handle log(x+y) as log(z), z=x+y, because the estimation problem will be smaller then without making the estimator worse
981 * (cons_nonlinear does this, too)
982 * this check takes care of this when x and y are original variables
983 * however, it isn't unlikely that we will have sums that become linear after we add auxvars for some children
984 * this will be handled in a postprocessing below
985 * for now, the check is performed on the original expression since there is not enough information in nlexpr yet
986 */
987 #ifdef SCIP_MORE_DEBUG
988 SCIPprintExpr(scip, SCIPhashmapGetImage(nlexpr2origexpr, (void*)nlexpr), NULL);
989 SCIPinfoMessage(scip, NULL, "... is a multivariate linear sum that we'll treat as auxvar\n");
990 #endif
991 }
992 else if( SCIPexprGetCurvature(nlexpr) != SCIP_EXPRCURV_UNKNOWN && !assumecurvature )
993 {
994 /* if we are here, either convexity or concavity is required; try to check for this curvature */
995 SCIP_Bool success;
996 int method;
997
998 /* try through curvature check methods until one succeeds */
999 for( method = 0; method < NCURVCHECKS; ++method )
1000 {
1001 SCIP_CALL( CURVCHECKS[method](scip, nlexpr, isrootexpr, &stack, nlexpr2origexpr, nlhdlrdata, assumevarfixed, &success) );
1002 if( success )
1003 break;
1004 }
1005 }
1006 else
1007 {
1008 /* if we don't care about curvature in this subtree anymore (very unlikely),
1009 * or we are told to assume that the desired curvature is present (assumecurvature==TRUE),
1010 * then only continue iterating this subtree to assemble leaf expressions
1011 */
1012 SCIP_CALL( nlhdlrExprGrowChildren(scip, nlexpr2origexpr, nlexpr, NULL) );
1013
1014 /* add children expressions, if any, to to-do list (stack) */
1015 SCIP_CALL( exprstackPush(scip, &stack, SCIPexprGetNChildren(nlexpr), SCIPexprGetChildren(nlexpr)) );
1016 }
1017 assert(stack.stackpos >= oldstackpos); /* none of the methods above should have removed something from the stack */
1018
1019 isrootexpr = FALSE;
1020
1021 /* if nothing was added, then none of the successors of nlexpr were added to the stack
1022 * this is either because nlexpr was already a variable or value expressions, thus a leaf,
1023 * or because the desired curvature could not be achieved, so it will be handled as variables, thus a leaf
1024 */
1025 if( stack.stackpos == oldstackpos )
1026 {
1027 ++*nleafs;
1028
1029 /* check whether the new leaf is not an original variable (or constant) */
1030 if( curvsuccess != NULL && !SCIPisExprVar(scip, nlexpr) && !SCIPisExprValue(scip, nlexpr) )
1031 *curvsuccess = FALSE;
1032 }
1033 }
1034
1035 exprstackFree(scip, &stack);
1036
1037 if( !nlhdlrdata->isnlhdlrconvex && *rootnlexpr != NULL )
1038 {
1039 /* remove multivariate linear subexpressions, that is, change some f(z1+z2) into f(z3) (z3=z1+z2 will be done by nlhdlr_default)
1040 * this handles the case that was not covered by the above check, which could recognize f(x+y) for x, y original variables
1041 */
1042 SCIP_EXPRITER* it;
1043
1044 SCIP_CALL( SCIPcreateExpriter(scip, &it) );
1045 SCIP_CALL( SCIPexpriterInit(it, *rootnlexpr, SCIP_EXPRITER_DFS, FALSE) );
1046 SCIPexpriterSetStagesDFS(it, SCIP_EXPRITER_VISITINGCHILD);
1047
1048 while( !SCIPexpriterIsEnd(it) )
1049 {
1050 SCIP_EXPR* child;
1051
1052 child = SCIPexpriterGetChildExprDFS(it);
1053 assert(child != NULL);
1054
1055 /* We want to change some f(x+y+z) into just f(), where f is the expression the iterator points to
1056 * and x+y+z is child. A child of a child, e.g., z, may not be a variable yet (these are added in collectLeafs later),
1057 * but an expression of some nonlinear type without children.
1058 */
1059 if( exprIsMultivarLinear(scip, child) )
1060 {
1061 /* turn child (x+y+z) into a sum without children
1062 * collectLeafs() should then replace this by an auxvar
1063 */
1064 #ifdef SCIP_MORE_DEBUG
1065 SCIPprintExpr(scip, child, NULL);
1066 SCIPinfoMessage(scip, NULL, "... is a multivariate linear sum that we'll treat as auxvar instead (postprocess)\n");
1067 #endif
1068
1069 /* TODO remove children from nlexpr2origexpr ?
1070 * should also do this if they are not used somewhere else; we could check nuses for this
1071 * however, it shouldn't matter to have some stray entries in the hashmap either
1072 */
1073 SCIP_CALL( SCIPremoveExprChildren(scip, child) );
1074 assert(SCIPexprGetNChildren(child) == 0);
1075
1076 (void) SCIPexpriterSkipDFS(it);
1077 }
1078 else
1079 {
1080 (void) SCIPexpriterGetNext(it);
1081 }
1082 }
1083
1084 SCIPfreeExpriter(&it);
1085 }
1086
1087 if( *rootnlexpr != NULL )
1088 {
1089 SCIP_Bool istrivial = TRUE;
1090
1091 /* if handletrivial is enabled, then only require that rootnlexpr itself has required curvature (so has children; see below) and
1092 * that we are not a trivial sum (because the previous implementation of this nlhdlr didn't allow this, either)
1093 */
1094 if( !nlhdlrdata->handletrivial || SCIPisExprSum(scip, *rootnlexpr) )
1095 {
1096 /* if all children do not have children, i.e., are variables, or will be replaced by auxvars, then free
1097 * also if rootnlexpr has no children, then free
1098 */
1099 int i;
1100 for( i = 0; i < SCIPexprGetNChildren(*rootnlexpr); ++i )
1101 {
1102 if( SCIPexprGetNChildren(SCIPexprGetChildren(*rootnlexpr)[i]) > 0 )
1103 {
1104 istrivial = FALSE;
1105 break;
1106 }
1107 }
1108 }
1109 else if( SCIPexprGetNChildren(*rootnlexpr) > 0 ) /* if handletrivial, then just require children */
1110 istrivial = FALSE;
1111
1112 if( istrivial )
1113 {
1114 SCIP_CALL( SCIPreleaseExpr(scip, rootnlexpr) );
1115 }
1116 }
1117
1118 return SCIP_OKAY;
1119 }
1120
1121 /** collects (non-value) leaf expressions and ensure that they correspond to a variable (original or auxiliary)
1122 *
1123 * For children where we could not achieve the desired curvature, get the auxvar and replace the child by a
1124 * var-expression that points to this auxvar.
1125 * Collect all leaf expressions (if not a value-expression) and index them.
1126 */
1127 static
1128 SCIP_RETCODE collectLeafs(
1129 SCIP* scip, /**< SCIP data structure */
1130 SCIP_NLHDLREXPRDATA* nlhdlrexprdata /**< nlhdlr expression data */
1131 )
1132 {
1133 SCIP_EXPRITER* it;
1134 SCIP_EXPR* nlexpr;
1135 SCIP_HASHMAP* leaf2index;
1136 int i;
1137
1138 assert(nlhdlrexprdata != NULL);
1139 assert(nlhdlrexprdata->nlexpr != NULL);
1140 assert(nlhdlrexprdata->nlexpr2origexpr != NULL);
1141 /* nleafs should be the upper bound on the number of variables given by constructExpr
1142 * leafexprs should be NULL, as this is what we want to setup here
1143 */
1144 assert(nlhdlrexprdata->nleafs > 0);
1145 assert(nlhdlrexprdata->leafexprs == NULL);
1146
1147 /* collect all auxvars and collect all variables */
1148 SCIP_CALL( SCIPhashmapCreate(&leaf2index, SCIPblkmem(scip), nlhdlrexprdata->nleafs) );
1149 nlhdlrexprdata->nleafs = 0; /* we start a new count, this time skipping value-expressions */
1150
1151 SCIP_CALL( SCIPcreateExpriter(scip, &it) );
1152 SCIP_CALL( SCIPexpriterInit(it, nlhdlrexprdata->nlexpr, SCIP_EXPRITER_DFS, FALSE) );
1153 SCIPexpriterSetStagesDFS(it, SCIP_EXPRITER_VISITINGCHILD);
1154
1155 for( nlexpr = SCIPexpriterGetCurrent(it); !SCIPexpriterIsEnd(it); nlexpr = SCIPexpriterGetNext(it) )
1156 {
1157 SCIP_EXPR* child;
1158 SCIP_EXPR* origexpr;
1159
1160 assert(nlexpr != NULL);
1161
1162 child = SCIPexpriterGetChildExprDFS(it);
1163
1164 /* if the to-be-visited child has children, then it doesn't need to be replaced by a new expression (representing the auxvar) */
1165 if( SCIPexprGetNChildren(child) > 0 )
1166 continue;
1167
1168 origexpr = (SCIP_EXPR*)SCIPhashmapGetImage(nlhdlrexprdata->nlexpr2origexpr, (void*)child);
1169 assert(origexpr != NULL);
1170
1171 if( SCIPexprGetNChildren(origexpr) > 0 )
1172 {
1173 SCIP_EXPR* newchild;
1174 int childidx;
1175 SCIP_VAR* var;
1176
1177 /* having a child that had children in original but not in copy means that we could not achieve the desired curvature
1178 * thus, replace by a new child that points to the auxvar of the original expression
1179 * we registered in createNlhdlrExprData that we need an auxvar, so it should exist now
1180 */
1181 var = SCIPgetExprAuxVarNonlinear(origexpr);
1182 assert(var != NULL);
1183
1184 SCIP_CALL( SCIPcreateExprVar(scip, &newchild, var, NULL, NULL) ); /* this captures newchild once */
1185
1186 childidx = SCIPexpriterGetChildIdxDFS(it);
1187 SCIP_CALL( SCIPreplaceExprChild(scip, nlexpr, childidx, newchild) ); /* this captures newchild again */
1188
1189 /* do not remove child->origexpr from hashmap, as child may appear again due to common subexprs
1190 * (created by curvCheckProductComposite, for example)
1191 * if it doesn't reappear, though, but the memory address is reused, we need to make sure it
1192 * points to the right origexpr
1193 */
1194 /* SCIP_CALL( SCIPhashmapRemove(nlexpr2origexpr, (void*)child) ); */
1195 SCIP_CALL( SCIPhashmapSetImage(nlhdlrexprdata->nlexpr2origexpr, (void*)newchild, (void*)origexpr) );
1196
1197 if( !SCIPhashmapExists(leaf2index, (void*)newchild) )
1198 {
1199 /* new leaf -> new index and remember in hashmap */
1200 SCIP_CALL( SCIPhashmapInsertInt(leaf2index, (void*)newchild, nlhdlrexprdata->nleafs++) );
1201 }
1202
1203 child = newchild;
1204 SCIP_CALL( SCIPreleaseExpr(scip, &newchild) ); /* because it was captured by both create and replace */
1205 }
1206 else if( SCIPisExprVar(scip, child) )
1207 {
1208 /* if variable, then add to hashmap, if not already there */
1209 if( !SCIPhashmapExists(leaf2index, (void*)child) )
1210 {
1211 SCIP_CALL( SCIPhashmapInsertInt(leaf2index, (void*)child, nlhdlrexprdata->nleafs++) );
1212 }
1213 }
1214 /* else: it's probably a value-expression, nothing to do */
1215
1216 /* update integrality flag for future leaf expressions: convex nlhdlr may use this information */
1217 SCIP_CALL( SCIPcomputeExprIntegrality(scip, child) );
1218 }
1219 assert(nlhdlrexprdata->nleafs > 0);
1220
1221 SCIPfreeExpriter(&it);
1222
1223 /* assemble auxvars array */
1224 SCIP_CALL( SCIPallocBlockMemoryArray(scip, &(nlhdlrexprdata->leafexprs), nlhdlrexprdata->nleafs) );
1225 for( i = 0; i < SCIPhashmapGetNEntries(leaf2index); ++i )
1226 {
1227 SCIP_HASHMAPENTRY* entry;
1228 SCIP_EXPR* leaf;
1229 int idx;
1230
1231 entry = SCIPhashmapGetEntry(leaf2index, i);
1232 if( entry == NULL )
1233 continue;
1234
1235 leaf = (SCIP_EXPR*) SCIPhashmapEntryGetOrigin(entry);
1236 assert(leaf != NULL);
1237 assert(SCIPisExprVar(scip, leaf));
1238
1239 idx = SCIPhashmapEntryGetImageInt(entry);
1240 assert(idx >= 0);
1241 assert(idx < nlhdlrexprdata->nleafs);
1242
1243 nlhdlrexprdata->leafexprs[idx] = leaf;
1244
1245 SCIPdebugMsg(scip, "leaf %d: <%s>\n", idx, SCIPvarGetName(SCIPgetVarExprVar(leaf)));
1246 }
1247
1248 SCIPhashmapFree(&leaf2index);
1249
1250 return SCIP_OKAY;
1251 }
1252
1253 /** creates nonlinear handler expression data structure and registers expr usage */
1254 static
1255 SCIP_RETCODE createNlhdlrExprData(
1256 SCIP* scip, /**< SCIP data structure */
1257 SCIP_NLHDLRDATA* nlhdlrdata, /**< nlhdlr data */
1258 SCIP_NLHDLREXPRDATA** nlhdlrexprdata, /**< pointer to store nlhdlr expression data */
1259 SCIP_EXPR* expr, /**< original expression */
1260 SCIP_EXPR* nlexpr, /**< our copy of expression */
1261 SCIP_HASHMAP* nlexpr2origexpr, /**< mapping of expression copy to original */
1262 int nleafs, /**< number of leafs as counted by constructExpr */
1263 SCIP_NLHDLR_METHOD participating /**< the enfo methods in which we plan to participate */
1264 )
1265 {
1266 SCIP_EXPRITER* it;
1267 SCIP_Bool usingaux;
1268
1269 assert(scip != NULL);
1270 assert(expr != NULL);
1271 assert(nlhdlrexprdata != NULL);
1272 assert(*nlhdlrexprdata == NULL);
1273 assert(nlexpr != NULL);
1274 assert(nlexpr2origexpr != NULL);
1275
1276 assert(SCIPexprGetNChildren(nlexpr) > 0);
1277 assert(SCIPexprGetChildren(nlexpr) != NULL);
1278
1279 SCIP_CALL( SCIPallocClearBlockMemory(scip, nlhdlrexprdata) );
1280 (*nlhdlrexprdata)->nlexpr = nlexpr;
1281 (*nlhdlrexprdata)->nlexpr2origexpr = nlexpr2origexpr;
1282 (*nlhdlrexprdata)->nleafs = nleafs;
1283
1284 usingaux = FALSE;
1285
1286 SCIP_CALL( SCIPcreateExpriter(scip, &it) );
1287 SCIP_CALL( SCIPexpriterInit(it, nlexpr, SCIP_EXPRITER_DFS, FALSE) );
1288 SCIPexpriterSetStagesDFS(it, SCIP_EXPRITER_VISITINGCHILD);
1289
1290 for( ; !SCIPexpriterIsEnd(it); (void) SCIPexpriterGetNext(it) )
1291 {
1292 SCIP_EXPR* child;
1293 SCIP_EXPR* origexpr;
1294
1295 /* check whether to-be-visited child needs to be replaced by a new expression (representing the auxvar)
1296 * if child has children, then that is not the case
1297 * if child has no children, but also corresponding origexpr has no chilren, then this is also not the case
1298 */
1299 child = SCIPexpriterGetChildExprDFS(it);
1300 if( SCIPexprGetNChildren(child) > 0 )
1301 continue;
1302
1303 origexpr = (SCIP_EXPR*)SCIPhashmapGetImage(nlexpr2origexpr, (void*)child);
1304 assert(origexpr != NULL);
1305
1306 /* if child had children in original but not in copy means that we could not achieve the desired curvature
1307 * thus, we will later replace by a new child that points to the auxvar of the original expression
1308 * as we do not have the auxvar now, we will only register that we will need the auxvar later (if origexpr isn't a variable or constant)
1309 * if we are working for the concave nlhdlr, then we also indicate interest on the exprs activity for estimate (distinguish below or above)
1310 */
1311 SCIP_CALL( SCIPregisterExprUsageNonlinear(scip, origexpr,
1312 SCIPexprGetNChildren(origexpr) > 0, FALSE,
1313 !nlhdlrdata->isnlhdlrconvex && (participating & SCIP_NLHDLR_METHOD_SEPABELOW),
1314 !nlhdlrdata->isnlhdlrconvex && (participating & SCIP_NLHDLR_METHOD_SEPAABOVE)) );
1315
1316 /* remember that we use an auxvar */
1317 if( SCIPexprGetNChildren(origexpr) > 0 )
1318 usingaux = TRUE;
1319 }
1320
1321 SCIPfreeExpriter(&it);
1322
1323 #ifdef SCIP_DEBUG
1324 SCIPprintExpr(scip, nlexpr, NULL);
1325 SCIPinfoMessage(scip, NULL, " (%p) is handled as %s\n", SCIPhashmapGetImage(nlexpr2origexpr, (void*)nlexpr),
1326 SCIPexprcurvGetName(SCIPexprGetCurvature(nlexpr)));
1327 #endif
1328
1329 /* If we don't work on the extended formulation, then set curvature also in original expression
1330 * (in case someone wants to pick this up; this might be removed again).
1331 * This doesn't ensure that every convex or concave original expression is actually marked here.
1332 * Not only because our tests are incomprehensive, but also because we may not detect on sums,
1333 * prefer extended formulations (in nlhdlr_convex), or introduce auxvars for linear subexpressions
1334 * on purpose (in nlhdlr_concave).
1335 */
1336 if( !usingaux )
1337 SCIPexprSetCurvature(expr, SCIPexprGetCurvature(nlexpr));
1338
1339 return SCIP_OKAY;
1340 }
1341
1342 /** adds an estimator for a vertex-polyhedral (e.g., concave) function to a given rowprep
1343 *
1344 * Calls \ref SCIPcomputeFacetVertexPolyhedralNonlinear() for given function and
1345 * box set to local bounds of auxiliary variables.
1346 */
1347 static
1348 SCIP_RETCODE estimateVertexPolyhedral(
1349 SCIP* scip, /**< SCIP data structure */
1350 SCIP_CONSHDLR* conshdlr, /**< nonlinear constraint handler */
1351 SCIP_NLHDLR* nlhdlr, /**< nonlinear handler */
1352 SCIP_NLHDLREXPRDATA* nlhdlrexprdata, /**< nonlinear handler expression data */
1353 SCIP_SOL* sol, /**< solution to use, unless usemidpoint is TRUE */
1354 SCIP_Bool usemidpoint, /**< whether to use the midpoint of the domain instead of sol */
1355 SCIP_Bool overestimate, /**< whether over- or underestimating */
1356 SCIP_Real targetvalue, /**< a target value to achieve; if not reachable, then can give up early */
1357 SCIP_ROWPREP* rowprep, /**< rowprep where to store estimator */
1358 SCIP_Bool* success /**< buffer to store whether successful */
1359 )
1360 {
1361 SCIP_NLHDLRDATA* nlhdlrdata;
1362 VERTEXPOLYFUN_EVALDATA evaldata;
1363 SCIP_Real* xstar;
1364 SCIP_Real* box;
1365 SCIP_Real facetconstant;
1366 SCIP_VAR* var;
1367 int i;
1368 SCIP_Bool allfixed;
1369
1370 assert(scip != NULL);
1371 assert(nlhdlr != NULL);
1372 assert(nlhdlrexprdata != NULL);
1373 assert(rowprep != NULL);
1374 assert(success != NULL);
1375
1376 *success = FALSE;
1377
1378 /* caller is responsible to have checked whether we can estimate, i.e., expression curvature and overestimate flag match */
1379 assert( overestimate || SCIPexprGetCurvature(nlhdlrexprdata->nlexpr) == SCIP_EXPRCURV_CONCAVE); /* if underestimate, then must be concave */
1380 assert(!overestimate || SCIPexprGetCurvature(nlhdlrexprdata->nlexpr) == SCIP_EXPRCURV_CONVEX); /* if overestimate, then must be convex */
1381
1382 #ifdef SCIP_DEBUG
1383 SCIPinfoMessage(scip, NULL, "%sestimate expression ", overestimate ? "over" : "under");
1384 SCIPprintExpr(scip, nlhdlrexprdata->nlexpr, NULL);
1385 SCIPinfoMessage(scip, NULL, " at point\n");
1386 for( i = 0; i < nlhdlrexprdata->nleafs; ++i )
1387 {
1388 var = SCIPgetVarExprVar(nlhdlrexprdata->leafexprs[i]);
1389 assert(var != NULL);
1390
1391 SCIPinfoMessage(scip, NULL, " <%s> = %g [%g,%g]\n", SCIPvarGetName(var),
1392 usemidpoint ? 0.5 * (SCIPvarGetLbLocal(var) + SCIPvarGetUbLocal(var)) : SCIPgetSolVal(scip, sol, var),
1393 SCIPvarGetLbLocal(var), SCIPvarGetUbLocal(var));
1394 }
1395 #endif
1396
1397 nlhdlrdata = SCIPnlhdlrGetData(nlhdlr);
1398 assert(nlhdlrdata != NULL);
1399
1400 if( nlhdlrdata->evalsol == NULL )
1401 {
1402 SCIP_CALL( SCIPcreateSol(scip, &nlhdlrdata->evalsol, NULL) );
1403 }
1404
1405 evaldata.nlhdlrexprdata = nlhdlrexprdata;
1406 evaldata.evalsol = nlhdlrdata->evalsol;
1407 evaldata.scip = scip;
1408
1409 SCIP_CALL( SCIPallocBufferArray(scip, &xstar, nlhdlrexprdata->nleafs) );
1410 SCIP_CALL( SCIPallocBufferArray(scip, &box, 2*nlhdlrexprdata->nleafs) );
1411
1412 allfixed = TRUE;
1413 for( i = 0; i < nlhdlrexprdata->nleafs; ++i )
1414 {
1415 var = SCIPgetVarExprVar(nlhdlrexprdata->leafexprs[i]);
1416 assert(var != NULL);
1417
1418 box[2*i] = SCIPvarGetLbLocal(var);
1419 if( SCIPisInfinity(scip, -box[2*i]) )
1420 {
1421 SCIPdebugMsg(scip, "lower bound at -infinity, no estimate possible\n");
1422 goto TERMINATE;
1423 }
1424
1425 box[2*i+1] = SCIPvarGetUbLocal(var);
1426 if( SCIPisInfinity(scip, box[2*i+1]) )
1427 {
1428 SCIPdebugMsg(scip, "upper bound at +infinity, no estimate possible\n");
1429 goto TERMINATE;
1430 }
1431
1432 if( !SCIPisRelEQ(scip, box[2*i], box[2*i+1]) )
1433 allfixed = FALSE;
1434
1435 if( usemidpoint )
1436 xstar[i] = 0.5 * (box[2*i] + box[2*i+1]);
1437 else
1438 xstar[i] = SCIPgetSolVal(scip, sol, var);
1439 assert(xstar[i] != SCIP_INVALID);
1440 }
1441
1442 if( allfixed )
1443 {
1444 /* SCIPcomputeFacetVertexPolyhedralNonlinear prints a warning and does not succeed if all is fixed */
1445 SCIPdebugMsg(scip, "all variables fixed, skip estimate\n");
1446 goto TERMINATE;
1447 }
1448
1449 SCIP_CALL( SCIPensureRowprepSize(scip, rowprep, nlhdlrexprdata->nleafs + 1) );
1450
1451 SCIP_CALL( SCIPcomputeFacetVertexPolyhedralNonlinear(scip, conshdlr, overestimate, nlhdlrExprEvalConcave, (void*)&evaldata,
1452 xstar, box, nlhdlrexprdata->nleafs, targetvalue, success, SCIProwprepGetCoefs(rowprep), &facetconstant) );
1453
1454 if( !*success )
1455 {
1456 SCIPdebugMsg(scip, "failed to compute facet of convex hull\n");
1457 goto TERMINATE;
1458 }
1459
1460 SCIProwprepSetLocal(rowprep, TRUE);
1461 SCIProwprepAddConstant(rowprep, facetconstant);
1462 for( i = 0; i < nlhdlrexprdata->nleafs; ++i )
1463 {
1464 SCIP_CALL( SCIPaddRowprepTerm(scip, rowprep, SCIPgetVarExprVar(nlhdlrexprdata->leafexprs[i]), SCIProwprepGetCoefs(rowprep)[i]) );
1465 }
1466
1467 #ifdef SCIP_DEBUG
1468 SCIPinfoMessage(scip, NULL, "computed estimator: ");
1469 SCIPprintRowprep(scip, rowprep, NULL);
1470 #endif
1471
1472 TERMINATE:
1473 SCIPfreeBufferArray(scip, &box);
1474 SCIPfreeBufferArray(scip, &xstar);
1475
1476 return SCIP_OKAY;
1477 }
1478
1479 /** adds an estimator computed via a gradient to a given rowprep */
1480 static
1481 SCIP_RETCODE estimateGradient(
1482 SCIP* scip, /**< SCIP data structure */
1483 SCIP_NLHDLREXPRDATA* nlhdlrexprdata, /**< nonlinear handler expression data */
1484 SCIP_SOL* sol, /**< solution to use */
1485 SCIP_Real auxvalue, /**< value of nlexpr in sol - we may not be able to take this value
1486 from nlexpr if it was evaluated at a different sol recently */
1487 SCIP_ROWPREP* rowprep, /**< rowprep where to store estimator */
1488 SCIP_Bool* success /**< buffer to store whether successful */
1489 )
1490 {
1491 SCIP_EXPR* nlexpr;
1492 SCIP_Real QUAD(constant);
1493 int i;
1494
1495 assert(nlhdlrexprdata != NULL);
1496 assert(rowprep != NULL);
1497 assert(success != NULL);
1498
1499 nlexpr = nlhdlrexprdata->nlexpr;
1500 assert(nlexpr != NULL);
1501
1502 #ifdef SCIP_DEBUG
1503 SCIPinfoMessage(scip, NULL, "estimate expression ");
1504 SCIPprintExpr(scip, nlexpr, NULL);
1505 SCIPinfoMessage(scip, NULL, " by gradient\n");
1506 #endif
1507
1508 *success = FALSE;
1509
1510 /* evaluation error -> skip */
1511 if( auxvalue == SCIP_INVALID )
1512 {
1513 SCIPdebugMsg(scip, "evaluation error / too large value (%g) for %p\n", auxvalue, (void*)nlexpr);
1514 return SCIP_OKAY;
1515 }
1516
1517 /* compute gradient (TODO: this also re-evaluates (soltag=0), which shouldn't be necessary unless we tried ConvexSecant before) */
1518 SCIP_CALL( SCIPevalExprGradient(scip, nlexpr, sol, 0L) );
1519
1520 /* gradient evaluation error -> skip */
1521 if( SCIPexprGetDerivative(nlexpr) == SCIP_INVALID )
1522 {
1523 SCIPdebugMsg(scip, "gradient evaluation error for %p\n", (void*)nlexpr);
1524 return SCIP_OKAY;
1525 }
1526
1527 /* add gradient underestimator to rowprep: f(sol) + (x - sol) \nabla f(sol)
1528 * constant will store f(sol) - sol * \nabla f(sol)
1529 * to avoid some cancellation errors when linear variables take huge values (like 1e20),
1530 * we use double-double arithemtic here
1531 */
1532 QUAD_ASSIGN(constant, SCIPexprGetEvalValue(nlexpr)); /* f(sol) */
1533 for( i = 0; i < nlhdlrexprdata->nleafs; ++i )
1534 {
1535 SCIP_VAR* var;
1536 SCIP_Real deriv;
1537 SCIP_Real varval;
1538
1539 assert(SCIPexprGetDiffTag(nlhdlrexprdata->leafexprs[i]) == SCIPexprGetDiffTag(nlexpr));
1540 deriv = SCIPexprGetDerivative(nlhdlrexprdata->leafexprs[i]);
1541 if( deriv == SCIP_INVALID )
1542 {
1543 SCIPdebugMsg(scip, "gradient evaluation error for component %d of %p\n", i, (void*)nlexpr);
1544 return SCIP_OKAY;
1545 }
1546
1547 var = SCIPgetVarExprVar(nlhdlrexprdata->leafexprs[i]);
1548 assert(var != NULL);
1549
1550 varval = SCIPgetSolVal(scip, sol, var);
1551
1552 SCIPdebugMsg(scip, "add %g * (<%s> - %g) to rowprep\n", deriv, SCIPvarGetName(var), varval);
1553
1554 /* add deriv * var to rowprep and deriv * (-varval) to constant */
1555 SCIP_CALL( SCIPaddRowprepTerm(scip, rowprep, var, deriv) );
1556 SCIPquadprecSumQD(constant, constant, -deriv * varval);
1557 }
1558
1559 SCIProwprepAddConstant(rowprep, QUAD_TO_DBL(constant));
1560 SCIProwprepSetLocal(rowprep, FALSE);
1561
1562 *success = TRUE;
1563
1564 return SCIP_OKAY;
1565 }
1566
1567 /** adds an estimator generated by putting a secant through the coordinates given by the two closest integer points */
1568 static
1569 SCIP_RETCODE estimateConvexSecant(
1570 SCIP* scip, /**< SCIP data structure */
1571 SCIP_NLHDLR* nlhdlr, /**< nonlinear handler */
1572 SCIP_NLHDLREXPRDATA* nlhdlrexprdata, /**< nonlinear handler expression data */
1573 SCIP_SOL* sol, /**< solution to use, unless usemidpoint is TRUE */
1574 SCIP_ROWPREP* rowprep, /**< rowprep where to store estimator */
1575 SCIP_Bool* success /**< buffer to store whether successful */
1576 )
1577 {
1578 SCIP_NLHDLRDATA* nlhdlrdata;
1579 SCIP_EXPR* nlexpr;
1580 SCIP_VAR* var;
1581 SCIP_Real x;
1582 SCIP_Real left, right;
1583 SCIP_Real fleft, fright;
1584
1585 assert(nlhdlrexprdata != NULL);
1586 assert(nlhdlrexprdata->nleafs == 1);
1587 assert(rowprep != NULL);
1588 assert(success != NULL);
1589
1590 nlexpr = nlhdlrexprdata->nlexpr;
1591 assert(nlexpr != NULL);
1592
1593 *success = FALSE;
1594
1595 nlhdlrdata = SCIPnlhdlrGetData(nlhdlr);
1596 assert(nlhdlrdata != NULL);
1597
1598 var = SCIPgetVarExprVar(nlhdlrexprdata->leafexprs[0]);
1599 assert(var != NULL);
1600
1601 x = SCIPgetSolVal(scip, sol, var);
1602
1603 #ifdef SCIP_DEBUG
1604 SCIPinfoMessage(scip, NULL, "estimate expression ");
1605 SCIPprintExpr(scip, nlexpr, NULL);
1606 SCIPinfoMessage(scip, NULL, " by secant\n");
1607 SCIPinfoMessage(scip, NULL, "integral variable <%s> = %g [%g,%g]\n", SCIPvarGetName(var),
1608 x, SCIPvarGetLbGlobal(var), SCIPvarGetUbGlobal(var));
1609 #endif
1610
1611 /* find out coordinates of var left and right to sol */
1612 if( SCIPisIntegral(scip, x) )
1613 {
1614 x = SCIPround(scip, x);
1615 if( SCIPisEQ(scip, x, SCIPvarGetLbGlobal(var)) )
1616 {
1617 left = x;
1618 right = left + 1.0;
1619 }
1620 else
1621 {
1622 right = x;
1623 left = right - 1.0;
1624 }
1625 }
1626 else
1627 {
1628 left = SCIPfloor(scip, x);
1629 right = SCIPceil(scip, x);
1630 }
1631 assert(left != right);
1632
1633 /* now evaluate at left and right */
1634 if( nlhdlrdata->evalsol == NULL )
1635 {
1636 SCIP_CALL( SCIPcreateSol(scip, &nlhdlrdata->evalsol, NULL) );
1637 }
1638
1639 SCIP_CALL( SCIPsetSolVal(scip, nlhdlrdata->evalsol, var, left) );
1640 SCIP_CALL( SCIPevalExpr(scip, nlexpr, nlhdlrdata->evalsol, 0L) );
1641
1642 /* evaluation error or a too large constant -> skip */
1643 fleft = SCIPexprGetEvalValue(nlexpr);
1644 if( SCIPisInfinity(scip, REALABS(fleft)) )
1645 {
1646 SCIPdebugMsg(scip, "evaluation error / too large value (%g) for %p\n", SCIPexprGetEvalValue(nlexpr), (void*)nlexpr);
1647 return SCIP_OKAY;
1648 }
1649
1650 SCIP_CALL( SCIPsetSolVal(scip, nlhdlrdata->evalsol, var, right) );
1651 SCIP_CALL( SCIPevalExpr(scip, nlexpr, nlhdlrdata->evalsol, 0L) );
1652
1653 /* evaluation error or a too large constant -> skip */
1654 fright = SCIPexprGetEvalValue(nlexpr);
1655 if( SCIPisInfinity(scip, REALABS(fright)) )
1656 {
1657 SCIPdebugMsg(scip, "evaluation error / too large value (%g) for %p\n", SCIPexprGetEvalValue(nlexpr), (void*)nlexpr);
1658 return SCIP_OKAY;
1659 }
1660
1661 SCIPdebugMsg(scip, "f(%g)=%g, f(%g)=%g\n", left, fleft, right, fright);
1662
1663 /* skip if too steep
1664 * for clay0204h, this resulted in a wrong cut from f(0)=1e12 f(1)=0.99998,
1665 * since due to limited precision, this was handled as if f(1)=1
1666 */
1667 if( (!SCIPisZero(scip, fleft) && REALABS(fright/fleft)*SCIPepsilon(scip) > 1.0) ||
1668 (!SCIPisZero(scip, fright) && REALABS(fleft/fright)*SCIPepsilon(scip) > 1.0) )
1669 {
1670 SCIPdebugMsg(scip, "function is too steep, abandoning\n");
1671 return SCIP_OKAY;
1672 }
1673
1674 /* now add f(left) + (f(right) - f(left)) * (x - left) as estimator to rowprep */
1675 SCIP_CALL( SCIPaddRowprepTerm(scip, rowprep, var, fright - fleft) );
1676 SCIProwprepAddConstant(rowprep, fleft - (fright - fleft) * left);
1677 SCIProwprepSetLocal(rowprep, FALSE);
1678
1679 *success = TRUE;
1680
1681 return SCIP_OKAY;
1682 }
1683
1684 /*
1685 * Callback methods of convex nonlinear handler
1686 */
1687
1688 /** free handler data of convex or concave nlhdlr */
1689 static
1690 SCIP_DECL_NLHDLRFREEHDLRDATA(nlhdlrfreeHdlrDataConvexConcave)
1691 { /*lint --e{715}*/
1692 assert(scip != NULL);
1693 assert(nlhdlrdata != NULL);
1694 assert(*nlhdlrdata != NULL);
1695 assert((*nlhdlrdata)->evalsol == NULL);
1696
1697 SCIPfreeBlockMemory(scip, nlhdlrdata);
1698
1699 return SCIP_OKAY;
1700 }
1701
1702 /** callback to free expression specific data */
1703 static
1704 SCIP_DECL_NLHDLRFREEEXPRDATA(nlhdlrfreeExprDataConvexConcave)
1705 { /*lint --e{715}*/
1706 assert(scip != NULL);
1707 assert(nlhdlrexprdata != NULL);
1708 assert(*nlhdlrexprdata != NULL);
1709
1710 SCIPfreeBlockMemoryArrayNull(scip, &(*nlhdlrexprdata)->leafexprs, (*nlhdlrexprdata)->nleafs);
1711 SCIP_CALL( SCIPreleaseExpr(scip, &(*nlhdlrexprdata)->nlexpr) );
1712 SCIPhashmapFree(&(*nlhdlrexprdata)->nlexpr2origexpr);
1713
1714 SCIPfreeBlockMemory(scip, nlhdlrexprdata);
1715
1716 return SCIP_OKAY;
1717 }
1718
1719 /** deinitialization of problem-specific data */
1720 static
1721 SCIP_DECL_NLHDLREXIT(nlhdlrExitConvex)
1722 {
1723 SCIP_NLHDLRDATA* nlhdlrdata;
1724
1725 nlhdlrdata = SCIPnlhdlrGetData(nlhdlr);
1726 assert(nlhdlrdata != NULL);
1727
1728 if( nlhdlrdata->evalsol != NULL )
1729 {
1730 SCIP_CALL( SCIPfreeSol(scip, &nlhdlrdata->evalsol) );
1731 }
1732
1733 return SCIP_OKAY;
1734 }
1735
1736 /** checks whether expression (or -expression) is convex, possibly after introducing auxiliary variables */
1737 static
1738 SCIP_DECL_NLHDLRDETECT(nlhdlrDetectConvex)
1739 { /*lint --e{715}*/
1740 SCIP_NLHDLRDATA* nlhdlrdata;
1741 SCIP_EXPR* nlexpr = NULL;
1742 SCIP_HASHMAP* nlexpr2origexpr;
1743 int nleafs = 0;
1744
1745 assert(scip != NULL);
1746 assert(nlhdlr != NULL);
1747 assert(expr != NULL);
1748 assert(enforcing != NULL);
1749 assert(participating != NULL);
1750 assert(nlhdlrexprdata != NULL);
1751
1752 /* we currently do not participate if only activity computation is required */
1753 if( (*enforcing & SCIP_NLHDLR_METHOD_SEPABOTH) == SCIP_NLHDLR_METHOD_SEPABOTH )
1754 return SCIP_OKAY;
1755
1756 /* ignore pure constants and variables */
1757 if( SCIPexprGetNChildren(expr) == 0 )
1758 return SCIP_OKAY;
1759
1760 nlhdlrdata = SCIPnlhdlrGetData(nlhdlr);
1761 assert(nlhdlrdata != NULL);
1762 assert(nlhdlrdata->isnlhdlrconvex);
1763
1764 SCIPdebugMsg(scip, "nlhdlr_convex detect for expr %p\n", (void*)expr);
1765
1766 /* initialize mapping from copied expression to original one
1767 * 20 is not a bad estimate for the size of convex subexpressions that we can usually discover
1768 * when expressions will be allowed to store "user"data, we could get rid of this hashmap (TODO)
1769 */
1770 SCIP_CALL( SCIPhashmapCreate(&nlexpr2origexpr, SCIPblkmem(scip), 20) );
1771
1772 if( (*enforcing & SCIP_NLHDLR_METHOD_SEPABELOW) == 0 ) /* if no separation below yet */
1773 {
1774 SCIP_CALL( constructExpr(scip, nlhdlrdata, &nlexpr, nlexpr2origexpr, &nleafs, expr,
1775 SCIP_EXPRCURV_CONVEX, NULL, SCIPassumeConvexNonlinear(conshdlr), NULL) );
1776 if( nlexpr != NULL )
1777 {
1778 assert(SCIPexprGetNChildren(nlexpr) > 0); /* should not be trivial */
1779
1780 *participating |= SCIP_NLHDLR_METHOD_SEPABELOW;
1781
1782 SCIPdebugMsg(scip, "detected expr %p to be convex -> can enforce expr <= auxvar\n", (void*)expr);
1783 }
1784 else
1785 {
1786 SCIP_CALL( SCIPhashmapRemoveAll(nlexpr2origexpr) );
1787 }
1788 }
1789
1790 if( (*enforcing & SCIP_NLHDLR_METHOD_SEPAABOVE) == 0 && nlexpr == NULL ) /* if no separation above and not convex */
1791 {
1792 SCIP_CALL( constructExpr(scip, nlhdlrdata, &nlexpr, nlexpr2origexpr, &nleafs, expr,
1793 SCIP_EXPRCURV_CONCAVE, NULL, SCIPassumeConvexNonlinear(conshdlr), NULL) );
1794 if( nlexpr != NULL )
1795 {
1796 assert(SCIPexprGetNChildren(nlexpr) > 0); /* should not be trivial */
1797
1798 *participating |= SCIP_NLHDLR_METHOD_SEPAABOVE;
1799
1800 SCIPdebugMsg(scip, "detected expr %p to be concave -> can enforce expr >= auxvar\n", (void*)expr);
1801 }
1802 }
1803
1804 /* everything we participate in we also enforce */
1805 *enforcing |= *participating;
1806
1807 assert(*participating || nlexpr == NULL);
1808 if( !*participating )
1809 {
1810 SCIPhashmapFree(&nlexpr2origexpr);
1811 return SCIP_OKAY;
1812 }
1813
1814 /* create the expression data of the nonlinear handler
1815 * notify conshdlr about expr for which we will require auxiliary variables
1816 */
1817 SCIP_CALL( createNlhdlrExprData(scip, nlhdlrdata, nlhdlrexprdata, expr, nlexpr, nlexpr2origexpr, nleafs, *participating) );
1818
1819 return SCIP_OKAY;
1820 }
1821
1822 /** auxiliary evaluation callback */
1823 static
1824 SCIP_DECL_NLHDLREVALAUX(nlhdlrEvalAuxConvexConcave)
1825 { /*lint --e{715}*/
1826 assert(nlhdlrexprdata != NULL);
1827 assert(nlhdlrexprdata->nlexpr != NULL);
1828 assert(auxvalue != NULL);
1829
1830 SCIP_CALL( SCIPevalExpr(scip, nlhdlrexprdata->nlexpr, sol, 0L) );
1831 *auxvalue = SCIPexprGetEvalValue(nlhdlrexprdata->nlexpr);
1832
1833 return SCIP_OKAY;
1834 }
1835
1836 /** init sepa callback that initializes LP */
1837 static
1838 SCIP_DECL_NLHDLRINITSEPA(nlhdlrInitSepaConvex)
1839 { /*lint --e{715}*/
1840 SCIP_EXPR* nlexpr;
1841 SCIP_EXPRCURV curvature;
1842 SCIP_Bool success;
1843 SCIP_ROWPREP* rowprep = NULL;
1844 SCIP_ROW* row;
1845 SCIP_Real lb;
1846 SCIP_Real ub;
1847 SCIP_Real lambda;
1848 SCIP_SOL* sol;
1849 int k;
1850
1851 assert(scip != NULL);
1852 assert(expr != NULL);
1853 assert(nlhdlrexprdata != NULL);
1854
1855 /* setup nlhdlrexprdata->leafexprs */
1856 SCIP_CALL( collectLeafs(scip, nlhdlrexprdata) );
1857
1858 nlexpr = nlhdlrexprdata->nlexpr;
1859 assert(nlexpr != NULL);
1860 assert(SCIPhashmapGetImage(nlhdlrexprdata->nlexpr2origexpr, (void*)nlexpr) == expr);
1861
1862 curvature = SCIPexprGetCurvature(nlexpr);
1863 assert(curvature == SCIP_EXPRCURV_CONVEX || curvature == SCIP_EXPRCURV_CONCAVE);
1864
1865 /* we can only be estimating on the convex side */
1866 if( curvature == SCIP_EXPRCURV_CONVEX )
1867 overestimate = FALSE;
1868 else if( curvature == SCIP_EXPRCURV_CONCAVE )
1869 underestimate = FALSE;
1870 if( !overestimate && !underestimate )
1871 return SCIP_OKAY;
1872
1873 /* linearizes at 5 different points obtained as convex combination of the lower and upper bound of the variables
1874 * present in the convex expression; whether more weight is given to the lower or upper bound of a variable depends
1875 * on whether the fixing of the variable to that value is better for the objective function
1876 */
1877 SCIP_CALL( SCIPcreateSol(scip, &sol, NULL) );
1878
1879 *infeasible = FALSE;
1880
1881 for( k = 0; k < 5; ++k )
1882 {
1883 int i;
1884 lambda = 0.1 * (k+1); /* lambda = 0.1, 0.2, 0.3, 0.4, 0.5 */
1885
1886 for( i = 0; i < nlhdlrexprdata->nleafs; ++i )
1887 {
1888 SCIP_VAR* var;
1889
1890 var = SCIPgetVarExprVar(nlhdlrexprdata->leafexprs[i]);
1891
1892 lb = SCIPvarGetLbGlobal(var);
1893 ub = SCIPvarGetUbGlobal(var);
1894
1895 /* make sure the absolute values of bounds are not too large */
1896 if( ub > -INITLPMAXVARVAL )
1897 lb = MAX(lb, -INITLPMAXVARVAL);
1898 if( lb < INITLPMAXVARVAL )
1899 ub = MIN(ub, INITLPMAXVARVAL);
1900
1901 /* in the case when ub < -maxabsbnd or lb > maxabsbnd, we still want to at least make bounds finite */
1902 if( SCIPisInfinity(scip, -lb) )
1903 lb = MIN(-10.0, ub - 0.1*REALABS(ub));
1904 if( SCIPisInfinity(scip, ub) )
1905 ub = MAX( 10.0, lb + 0.1*REALABS(lb));
1906
1907 if( SCIPvarGetBestBoundType(var) == SCIP_BOUNDTYPE_LOWER )
1908 SCIP_CALL( SCIPsetSolVal(scip, sol, var, lambda * ub + (1.0 - lambda) * lb) );
1909 else
1910 SCIP_CALL( SCIPsetSolVal(scip, sol, var, lambda * lb + (1.0 - lambda) * ub) );
1911 }
1912
1913 SCIP_CALL( SCIPcreateRowprep(scip, &rowprep, overestimate ? SCIP_SIDETYPE_LEFT : SCIP_SIDETYPE_RIGHT, TRUE) );
1914 SCIP_CALL( estimateGradient(scip, nlhdlrexprdata, sol, 0.0, rowprep, &success) );
1915 if( !success )
1916 {
1917 SCIPdebugMsg(scip, "failed to linearize for k = %d\n", k);
1918 if( rowprep != NULL )
1919 SCIPfreeRowprep(scip, &rowprep);
1920 continue;
1921 }
1922
1923 /* add auxiliary variable */
1924 SCIP_CALL( SCIPaddRowprepTerm(scip, rowprep, SCIPgetExprAuxVarNonlinear(expr), -1.0) );
1925
1926 /* straighten out numerics */
1927 SCIP_CALL( SCIPcleanupRowprep2(scip, rowprep, NULL, SCIPgetHugeValue(scip), &success) );
1928 if( !success )
1929 {
1930 SCIPdebugMsg(scip, "failed to cleanup rowprep numerics for k = %d\n", k);
1931 if( rowprep != NULL )
1932 SCIPfreeRowprep(scip, &rowprep);
1933 continue;
1934 }
1935
1936 (void) SCIPsnprintf(SCIProwprepGetName(rowprep), SCIP_MAXSTRLEN, "%sestimate_gradient%p_initsepa_%d",
1937 overestimate ? "over" : "under", (void*)expr, k);
1938 SCIP_CALL( SCIPgetRowprepRowCons(scip, &row, rowprep, cons) );
1939 SCIPfreeRowprep(scip, &rowprep);
1940
1941 #ifdef SCIP_DEBUG
1942 SCIPinfoMessage(scip, NULL, "initsepa computed row: ");
1943 SCIPprintRow(scip, row, NULL);
1944 #endif
1945
1946 SCIP_CALL( SCIPaddRow(scip, row, FALSE, infeasible) );
1947 SCIP_CALL( SCIPreleaseRow(scip, &row) );
1948
1949 if( *infeasible )
1950 break;
1951 }
1952
1953 SCIP_CALL( SCIPfreeSol(scip, &sol) );
1954
1955 return SCIP_OKAY;
1956 }
1957
1958 /** estimator callback */
1959 static
1960 SCIP_DECL_NLHDLRESTIMATE(nlhdlrEstimateConvex)
1961 { /*lint --e{715}*/
1962 SCIP_ROWPREP* rowprep;
1963
1964 assert(scip != NULL);
1965 assert(expr != NULL);
1966 assert(nlhdlrexprdata != NULL);
1967 assert(nlhdlrexprdata->nlexpr != NULL);
1968 assert(rowpreps != NULL);
1969 assert(success != NULL);
1970
1971 assert(SCIPhashmapGetImage(nlhdlrexprdata->nlexpr2origexpr, (void*)nlhdlrexprdata->nlexpr) == expr);
1972
1973 /* we must be called only for the side that we indicated to participate in during DETECT */
1974 assert(SCIPexprGetCurvature(nlhdlrexprdata->nlexpr) == SCIP_EXPRCURV_CONVEX
1975 || SCIPexprGetCurvature(nlhdlrexprdata->nlexpr) == SCIP_EXPRCURV_CONCAVE);
1976 assert(!overestimate || SCIPexprGetCurvature(nlhdlrexprdata->nlexpr) == SCIP_EXPRCURV_CONCAVE);
1977 assert( overestimate || SCIPexprGetCurvature(nlhdlrexprdata->nlexpr) == SCIP_EXPRCURV_CONVEX);
1978
1979 *success = FALSE;
1980 *addedbranchscores = FALSE;
1981
1982 /* we can skip eval as nlhdlrEvalAux should have been called for same solution before */
1983 /* SCIP_CALL( nlhdlrExprEval(scip, nlexpr, sol) ); */
1984 assert(auxvalue == SCIPexprGetEvalValue(nlhdlrexprdata->nlexpr)); /* given value (originally from
1985 nlhdlrEvalAuxConvexConcave) should coincide with the one stored in nlexpr */
1986
1987 SCIP_CALL( SCIPcreateRowprep(scip, &rowprep, overestimate ? SCIP_SIDETYPE_LEFT : SCIP_SIDETYPE_RIGHT, TRUE) );
1988
1989 if( nlhdlrexprdata->nleafs == 1 && SCIPexprIsIntegral(nlhdlrexprdata->leafexprs[0]) )
1990 {
1991 SCIP_CALL( estimateConvexSecant(scip, nlhdlr, nlhdlrexprdata, sol, rowprep, success) );
1992
1993 (void) SCIPsnprintf(SCIProwprepGetName(rowprep), SCIP_MAXSTRLEN, "%sestimate_convexsecant%p_%s%d",
1994 overestimate ? "over" : "under",
1995 (void*)expr,
1996 sol != NULL ? "sol" : "lp",
1997 sol != NULL ? SCIPsolGetIndex(sol) : SCIPgetNLPs(scip));
1998 }
1999
2000 /* if secant method was not used or failed, then try with gradient */
2001 if( !*success )
2002 {
2003 SCIP_CALL( estimateGradient(scip, nlhdlrexprdata, sol, auxvalue, rowprep, success) );
2004
2005 (void) SCIPsnprintf(SCIProwprepGetName(rowprep), SCIP_MAXSTRLEN, "%sestimate_convexgradient%p_%s%d",
2006 overestimate ? "over" : "under",
2007 (void*)expr,
2008 sol != NULL ? "sol" : "lp",
(1) Event invalid_type: |
Argument "(sol != NULL) ? sol->index : SCIPgetNLPs(scip)" to format specifier "%d" was expected to have type "int" but has type "long long". [details] |
2009 sol != NULL ? SCIPsolGetIndex(sol) : SCIPgetNLPs(scip));
2010 }
2011
2012 if( *success )
2013 {
2014 SCIP_CALL( SCIPsetPtrarrayVal(scip, rowpreps, 0, rowprep) );
2015 }
2016 else
2017 {
2018 SCIPfreeRowprep(scip, &rowprep);
2019 }
2020
2021 return SCIP_OKAY;
2022 }
2023
2024 /** include nlhdlr in another scip instance */
2025 static
2026 SCIP_DECL_NLHDLRCOPYHDLR(nlhdlrCopyhdlrConvex)
2027 { /*lint --e{715}*/
2028 assert(targetscip != NULL);
2029 assert(sourcenlhdlr != NULL);
2030 assert(strcmp(SCIPnlhdlrGetName(sourcenlhdlr), CONVEX_NLHDLR_NAME) == 0);
2031
2032 SCIP_CALL( SCIPincludeNlhdlrConvex(targetscip) );
2033
2034 return SCIP_OKAY;
2035 }
2036
2037 /** includes convex nonlinear handler in nonlinear constraint handler */
2038 SCIP_RETCODE SCIPincludeNlhdlrConvex(
2039 SCIP* scip /**< SCIP data structure */
2040 )
2041 {
2042 SCIP_NLHDLR* nlhdlr;
2043 SCIP_NLHDLRDATA* nlhdlrdata;
2044
2045 assert(scip != NULL);
2046
2047 SCIP_CALL( SCIPallocBlockMemory(scip, &nlhdlrdata) );
2048 nlhdlrdata->isnlhdlrconvex = TRUE;
2049 nlhdlrdata->evalsol = NULL;
2050
2051 SCIP_CALL( SCIPincludeNlhdlrNonlinear(scip, &nlhdlr, CONVEX_NLHDLR_NAME, CONVEX_NLHDLR_DESC,
2052 CONVEX_NLHDLR_DETECTPRIORITY, CONVEX_NLHDLR_ENFOPRIORITY, nlhdlrDetectConvex, nlhdlrEvalAuxConvexConcave, nlhdlrdata) );
2053 assert(nlhdlr != NULL);
2054
2055 SCIP_CALL( SCIPaddBoolParam(scip, "nlhdlr/" CONVEX_NLHDLR_NAME "/detectsum",
2056 "whether to run convexity detection when the root of an expression is a non-quadratic sum",
2057 &nlhdlrdata->detectsum, FALSE, DEFAULT_DETECTSUM, NULL, NULL) );
2058
2059 SCIP_CALL( SCIPaddBoolParam(scip, "nlhdlr/" CONVEX_NLHDLR_NAME "/extendedform",
2060 "whether to create extended formulations instead of looking for maximal convex expressions",
2061 &nlhdlrdata->extendedform, FALSE, DEFAULT_EXTENDEDFORM, NULL, NULL) );
2062
2063 SCIP_CALL( SCIPaddBoolParam(scip, "nlhdlr/" CONVEX_NLHDLR_NAME "/cvxquadratic",
2064 "whether to use convexity check on quadratics",
2065 &nlhdlrdata->cvxquadratic, TRUE, DEFAULT_CVXQUADRATIC_CONVEX, NULL, NULL) );
2066
2067 SCIP_CALL( SCIPaddBoolParam(scip, "nlhdlr/" CONVEX_NLHDLR_NAME "/cvxsignomial",
2068 "whether to use convexity check on signomials",
2069 &nlhdlrdata->cvxsignomial, TRUE, DEFAULT_CVXSIGNOMIAL, NULL, NULL) );
2070
2071 SCIP_CALL( SCIPaddBoolParam(scip, "nlhdlr/" CONVEX_NLHDLR_NAME "/cvxprodcomp",
2072 "whether to use convexity check on product composition f(h)*h",
2073 &nlhdlrdata->cvxprodcomp, TRUE, DEFAULT_CVXPRODCOMP, NULL, NULL) );
2074
2075 SCIP_CALL( SCIPaddBoolParam(scip, "nlhdlr/" CONVEX_NLHDLR_NAME "/handletrivial",
2076 "whether to also handle trivial convex expressions",
2077 &nlhdlrdata->handletrivial, TRUE, DEFAULT_HANDLETRIVIAL, NULL, NULL) );
2078
2079 SCIPnlhdlrSetFreeHdlrData(nlhdlr, nlhdlrfreeHdlrDataConvexConcave);
2080 SCIPnlhdlrSetCopyHdlr(nlhdlr, nlhdlrCopyhdlrConvex);
2081 SCIPnlhdlrSetFreeExprData(nlhdlr, nlhdlrfreeExprDataConvexConcave);
2082 SCIPnlhdlrSetSepa(nlhdlr, nlhdlrInitSepaConvex, NULL, nlhdlrEstimateConvex, NULL);
2083 SCIPnlhdlrSetInitExit(nlhdlr, NULL, nlhdlrExitConvex);
2084
2085 return SCIP_OKAY;
2086 }
2087
2088 /*
2089 * Callback methods of concave nonlinear handler
2090 */
2091
2092 /** deinitialization of problem-specific data */
2093 static
2094 SCIP_DECL_NLHDLREXIT(nlhdlrExitConcave)
2095 {
2096 SCIP_NLHDLRDATA* nlhdlrdata;
2097
2098 nlhdlrdata = SCIPnlhdlrGetData(nlhdlr);
2099 assert(nlhdlrdata != NULL);
2100
2101 if( nlhdlrdata->evalsol != NULL )
2102 {
2103 SCIP_CALL( SCIPfreeSol(scip, &nlhdlrdata->evalsol) );
2104 }
2105
2106 return SCIP_OKAY;
2107 }
2108
2109 /** checks whether expression (or -expression) is concave, possibly after introducing auxiliary variables */
2110 static
2111 SCIP_DECL_NLHDLRDETECT(nlhdlrDetectConcave)
2112 { /*lint --e{715}*/
2113 SCIP_NLHDLRDATA* nlhdlrdata;
2114 SCIP_EXPR* nlexpr = NULL;
2115 SCIP_HASHMAP* nlexpr2origexpr;
2116 int nleafs = 0;
2117
2118 assert(scip != NULL);
2119 assert(nlhdlr != NULL);
2120 assert(expr != NULL);
2121 assert(enforcing != NULL);
2122 assert(participating != NULL);
2123 assert(nlhdlrexprdata != NULL);
2124
2125 /* we currently do not participate if only activity computation is required */
2126 if( (*enforcing & SCIP_NLHDLR_METHOD_SEPABOTH) == SCIP_NLHDLR_METHOD_SEPABOTH )
2127 return SCIP_OKAY;
2128
2129 /* ignore pure constants and variables */
2130 if( SCIPexprGetNChildren(expr) == 0 )
2131 return SCIP_OKAY;
2132
2133 nlhdlrdata = SCIPnlhdlrGetData(nlhdlr);
2134 assert(nlhdlrdata != NULL);
2135 assert(!nlhdlrdata->isnlhdlrconvex);
2136
2137 SCIPdebugMsg(scip, "nlhdlr_concave detect for expr %p\n", (void*)expr);
2138
2139 /* initialize mapping from copied expression to original one
2140 * 20 is not a bad estimate for the size of concave subexpressions that we can usually discover
2141 * when expressions will be allowed to store "user"data, we could get rid of this hashmap (TODO)
2142 */
2143 SCIP_CALL( SCIPhashmapCreate(&nlexpr2origexpr, SCIPblkmem(scip), 20) );
2144
2145 if( (*enforcing & SCIP_NLHDLR_METHOD_SEPABELOW) == 0 ) /* if no separation below yet */
2146 {
2147 SCIP_CALL( constructExpr(scip, nlhdlrdata, &nlexpr, nlexpr2origexpr, &nleafs, expr,
2148 SCIP_EXPRCURV_CONCAVE, NULL, FALSE, NULL) );
2149
2150 if( nlexpr != NULL && nleafs > SCIP_MAXVERTEXPOLYDIM )
2151 {
2152 SCIPdebugMsg(scip, "Too many variables (%d) in constructed expression. Will not be able to estimate. Rejecting.\n", nleafs);
2153 SCIP_CALL( SCIPreleaseExpr(scip, &nlexpr) );
2154 }
2155
2156 if( nlexpr != NULL )
2157 {
2158 assert(SCIPexprGetNChildren(nlexpr) > 0); /* should not be trivial */
2159
2160 *participating |= SCIP_NLHDLR_METHOD_SEPABELOW;
2161
2162 SCIPdebugMsg(scip, "detected expr %p to be concave -> can enforce expr <= auxvar\n", (void*)expr);
2163 }
2164 else
2165 {
2166 SCIP_CALL( SCIPhashmapRemoveAll(nlexpr2origexpr) );
2167 }
2168 }
2169
2170 if( (*enforcing & SCIP_NLHDLR_METHOD_SEPAABOVE) == 0 && nlexpr == NULL ) /* if no separation above and not concave */
2171 {
2172 SCIP_CALL( constructExpr(scip, nlhdlrdata, &nlexpr, nlexpr2origexpr, &nleafs, expr,
2173 SCIP_EXPRCURV_CONVEX, NULL, FALSE, NULL) );
2174
2175 if( nlexpr != NULL && nleafs > SCIP_MAXVERTEXPOLYDIM )
2176 {
2177 SCIPdebugMsg(scip, "Too many variables (%d) in constructed expression. Will not be able to estimate. Rejecting.\n", nleafs);
2178 SCIP_CALL( SCIPreleaseExpr(scip, &nlexpr) );
2179 }
2180
2181 if( nlexpr != NULL )
2182 {
2183 assert(SCIPexprGetNChildren(nlexpr) > 0); /* should not be trivial */
2184
2185 *participating |= SCIP_NLHDLR_METHOD_SEPAABOVE;
2186
2187 SCIPdebugMsg(scip, "detected expr %p to be convex -> can enforce expr >= auxvar\n", (void*)expr);
2188 }
2189 }
2190
2191 /* everything we participate in we also enforce (at the moment) */
2192 *enforcing |= *participating;
2193
2194 assert(*participating || nlexpr == NULL);
2195 if( !*participating )
2196 {
2197 SCIPhashmapFree(&nlexpr2origexpr);
2198 return SCIP_OKAY;
2199 }
2200
2201 /* create the expression data of the nonlinear handler
2202 * notify conshdlr about expr for which we will require auxiliary variables and use activity
2203 */
2204 SCIP_CALL( createNlhdlrExprData(scip, nlhdlrdata, nlhdlrexprdata, expr, nlexpr, nlexpr2origexpr, nleafs, *participating) );
2205
2206 return SCIP_OKAY;
2207 }
2208
2209 /** init sepa callback that initializes LP */
2210 static
2211 SCIP_DECL_NLHDLRINITSEPA(nlhdlrInitSepaConcave)
2212 {
2213 SCIP_EXPR* nlexpr;
2214 SCIP_EXPRCURV curvature;
2215 SCIP_Bool success;
2216 SCIP_ROWPREP* rowprep = NULL;
2217 SCIP_ROW* row;
2218
2219 assert(scip != NULL);
2220 assert(expr != NULL);
2221 assert(nlhdlrexprdata != NULL);
2222
2223 nlexpr = nlhdlrexprdata->nlexpr;
2224 assert(nlexpr != NULL);
2225 assert(SCIPhashmapGetImage(nlhdlrexprdata->nlexpr2origexpr, (void*)nlexpr) == expr);
2226
2227 /* setup nlhdlrexprdata->leafexprs */
2228 SCIP_CALL( collectLeafs(scip, nlhdlrexprdata) );
2229
2230 curvature = SCIPexprGetCurvature(nlexpr);
2231 assert(curvature == SCIP_EXPRCURV_CONVEX || curvature == SCIP_EXPRCURV_CONCAVE);
2232 /* we can only be estimating on non-convex side */
2233 if( curvature == SCIP_EXPRCURV_CONCAVE )
2234 overestimate = FALSE;
2235 else if( curvature == SCIP_EXPRCURV_CONVEX )
2236 underestimate = FALSE;
2237 if( !overestimate && !underestimate )
2238 return SCIP_OKAY;
2239
2240 /* compute estimator and store in rowprep */
2241 SCIP_CALL( SCIPcreateRowprep(scip, &rowprep, overestimate ? SCIP_SIDETYPE_LEFT : SCIP_SIDETYPE_RIGHT, TRUE) );
2242 SCIP_CALL( estimateVertexPolyhedral(scip, conshdlr, nlhdlr, nlhdlrexprdata, NULL, TRUE, overestimate,
2243 overestimate ? SCIPinfinity(scip) : -SCIPinfinity(scip), rowprep, &success) );
2244 if( !success )
2245 {
2246 SCIPdebugMsg(scip, "failed to compute facet of convex hull\n");
2247 goto TERMINATE;
2248 }
2249
2250 /* add auxiliary variable */
2251 SCIP_CALL( SCIPaddRowprepTerm(scip, rowprep, SCIPgetExprAuxVarNonlinear(expr), -1.0) );
2252
2253 /* straighten out numerics */
2254 SCIP_CALL( SCIPcleanupRowprep2(scip, rowprep, NULL, SCIPgetHugeValue(scip), &success) );
2255 if( !success )
2256 {
2257 SCIPdebugMsg(scip, "failed to cleanup rowprep numerics\n");
2258 goto TERMINATE;
2259 }
2260
2261 (void) SCIPsnprintf(SCIProwprepGetName(rowprep), SCIP_MAXSTRLEN, "%sestimate_concave%p_initsepa",
2262 overestimate ? "over" : "under", (void*)expr);
2263 SCIP_CALL( SCIPgetRowprepRowCons(scip, &row, rowprep, cons) );
2264
2265 #ifdef SCIP_DEBUG
2266 SCIPinfoMessage(scip, NULL, "initsepa computed row: ");
2267 SCIPprintRow(scip, row, NULL);
2268 #endif
2269
2270 SCIP_CALL( SCIPaddRow(scip, row, FALSE, infeasible) );
2271 SCIP_CALL( SCIPreleaseRow(scip, &row) );
2272
2273 TERMINATE:
2274 if( rowprep != NULL )
2275 SCIPfreeRowprep(scip, &rowprep);
2276
2277 return SCIP_OKAY;
2278 }
2279
2280 /** estimator callback */
2281 static
2282 SCIP_DECL_NLHDLRESTIMATE(nlhdlrEstimateConcave)
2283 { /*lint --e{715}*/
2284 SCIP_ROWPREP* rowprep;
2285
2286 assert(scip != NULL);
2287 assert(expr != NULL);
2288 assert(nlhdlrexprdata != NULL);
2289 assert(nlhdlrexprdata->nlexpr != NULL);
2290 assert(rowpreps != NULL);
2291 assert(success != NULL);
2292
2293 assert(SCIPhashmapGetImage(nlhdlrexprdata->nlexpr2origexpr, (void*)nlhdlrexprdata->nlexpr) == expr);
2294
2295 /* we must be called only for the side that we indicated to participate in during DETECT */
2296 assert(SCIPexprGetCurvature(nlhdlrexprdata->nlexpr) == SCIP_EXPRCURV_CONVEX
2297 || SCIPexprGetCurvature(nlhdlrexprdata->nlexpr) == SCIP_EXPRCURV_CONCAVE);
2298 assert(!overestimate || SCIPexprGetCurvature(nlhdlrexprdata->nlexpr) == SCIP_EXPRCURV_CONVEX);
2299 assert( overestimate || SCIPexprGetCurvature(nlhdlrexprdata->nlexpr) == SCIP_EXPRCURV_CONCAVE);
2300
2301 *success = FALSE;
2302 *addedbranchscores = FALSE;
2303
2304 SCIP_CALL( SCIPcreateRowprep(scip, &rowprep, overestimate ? SCIP_SIDETYPE_LEFT : SCIP_SIDETYPE_RIGHT, TRUE) );
2305
2306 SCIP_CALL( estimateVertexPolyhedral(scip, conshdlr, nlhdlr, nlhdlrexprdata, sol, FALSE, overestimate, targetvalue, rowprep, success) );
2307
2308 if( *success )
2309 {
2310 SCIP_CALL( SCIPsetPtrarrayVal(scip, rowpreps, 0, rowprep) );
2311
2312 (void) SCIPsnprintf(SCIProwprepGetName(rowprep), SCIP_MAXSTRLEN, "%sestimate_concave%p_%s%d",
2313 overestimate ? "over" : "under",
2314 (void*)expr,
2315 sol != NULL ? "sol" : "lp",
2316 sol != NULL ? SCIPsolGetIndex(sol) : SCIPgetNLPs(scip));
2317 }
2318 else
2319 {
2320 SCIPfreeRowprep(scip, &rowprep);
2321 }
2322
2323 if( addbranchscores )
2324 {
2325 SCIP_Real violation;
2326
2327 /* check how much is the violation on the side that we estimate */
2328 if( auxvalue == SCIP_INVALID )
2329 {
2330 /* if cannot evaluate, then always branch */
2331 violation = SCIPinfinity(scip);
2332 }
2333 else
2334 {
2335 SCIP_Real auxval;
2336
2337 /* get value of auxiliary variable of this expression */
2338 assert(SCIPgetExprAuxVarNonlinear(expr) != NULL);
2339 auxval = SCIPgetSolVal(scip, sol, SCIPgetExprAuxVarNonlinear(expr));
2340
2341 /* compute the violation
2342 * if we underestimate, then we enforce expr <= auxval, so violation is (positive part of) auxvalue - auxval
2343 * if we overestimate, then we enforce expr >= auxval, so violation is (positive part of) auxval - auxvalue
2344 */
2345 if( !overestimate )
2346 violation = MAX(0.0, auxvalue - auxval);
2347 else
2348 violation = MAX(0.0, auxval - auxvalue);
2349 }
2350 assert(violation >= 0.0);
2351
2352 /* add violation as branching-score to expressions; the core will take care distributing this onto variables */
2353 if( nlhdlrexprdata->nleafs == 1 )
2354 {
2355 SCIP_EXPR* e;
2356 e = (SCIP_EXPR*)SCIPhashmapGetImage(nlhdlrexprdata->nlexpr2origexpr, nlhdlrexprdata->leafexprs[0]);
2357 SCIP_CALL( SCIPaddExprsViolScoreNonlinear(scip, &e, 1, violation, sol, addedbranchscores) );
2358 }
2359 else
2360 {
2361 SCIP_EXPR** exprs;
2362 int c;
2363
2364 /* map leaf expressions back to original expressions
2365 * TODO do this once at end of detect and store in nlhdlrexprdata
2366 */
2367 SCIP_CALL( SCIPallocBufferArray(scip, &exprs, nlhdlrexprdata->nleafs) );
2368 for( c = 0; c < nlhdlrexprdata->nleafs; ++c )
2369 exprs[c] = (SCIP_EXPR*)SCIPhashmapGetImage(nlhdlrexprdata->nlexpr2origexpr, nlhdlrexprdata->leafexprs[c]);
2370
2371 SCIP_CALL( SCIPaddExprsViolScoreNonlinear(scip, exprs, nlhdlrexprdata->nleafs, violation, sol, addedbranchscores) );
2372
2373 SCIPfreeBufferArray(scip, &exprs);
2374 }
2375 }
2376
2377 return SCIP_OKAY;
2378 }
2379
2380 /** includes nonlinear handler in another scip instance */
2381 static
2382 SCIP_DECL_NLHDLRCOPYHDLR(nlhdlrCopyhdlrConcave)
2383 { /*lint --e{715}*/
2384 assert(targetscip != NULL);
2385 assert(sourcenlhdlr != NULL);
2386 assert(strcmp(SCIPnlhdlrGetName(sourcenlhdlr), CONCAVE_NLHDLR_NAME) == 0);
2387
2388 SCIP_CALL( SCIPincludeNlhdlrConcave(targetscip) );
2389
2390 return SCIP_OKAY;
2391 }
2392
2393 /** includes concave nonlinear handler in nonlinear constraint handler */
2394 SCIP_EXPORT
2395 SCIP_RETCODE SCIPincludeNlhdlrConcave(
2396 SCIP* scip /**< SCIP data structure */
2397 )
2398 {
2399 SCIP_NLHDLR* nlhdlr;
2400 SCIP_NLHDLRDATA* nlhdlrdata;
2401
2402 assert(scip != NULL);
2403
2404 SCIP_CALL( SCIPallocBlockMemory(scip, &nlhdlrdata) );
2405 nlhdlrdata->isnlhdlrconvex = FALSE;
2406 nlhdlrdata->evalsol = NULL;
2407
2408 SCIP_CALL( SCIPincludeNlhdlrNonlinear(scip, &nlhdlr, CONCAVE_NLHDLR_NAME, CONCAVE_NLHDLR_DESC,
2409 CONCAVE_NLHDLR_DETECTPRIORITY, CONCAVE_NLHDLR_ENFOPRIORITY, nlhdlrDetectConcave, nlhdlrEvalAuxConvexConcave, nlhdlrdata) );
2410 assert(nlhdlr != NULL);
2411
2412 SCIP_CALL( SCIPaddBoolParam(scip, "nlhdlr/" CONCAVE_NLHDLR_NAME "/detectsum",
2413 "whether to run convexity detection when the root of an expression is a sum",
2414 &nlhdlrdata->detectsum, FALSE, DEFAULT_DETECTSUM, NULL, NULL) );
2415
2416 /* "extended" formulations of a concave expressions can give worse estimators */
2417 nlhdlrdata->extendedform = FALSE;
2418
2419 SCIP_CALL( SCIPaddBoolParam(scip, "nlhdlr/" CONCAVE_NLHDLR_NAME "/cvxquadratic",
2420 "whether to use convexity check on quadratics",
2421 &nlhdlrdata->cvxquadratic, TRUE, DEFAULT_CVXQUADRATIC_CONCAVE, NULL, NULL) );
2422
2423 SCIP_CALL( SCIPaddBoolParam(scip, "nlhdlr/" CONCAVE_NLHDLR_NAME "/cvxsignomial",
2424 "whether to use convexity check on signomials",
2425 &nlhdlrdata->cvxsignomial, TRUE, DEFAULT_CVXSIGNOMIAL, NULL, NULL) );
2426
2427 SCIP_CALL( SCIPaddBoolParam(scip, "nlhdlr/" CONCAVE_NLHDLR_NAME "/cvxprodcomp",
2428 "whether to use convexity check on product composition f(h)*h",
2429 &nlhdlrdata->cvxprodcomp, TRUE, DEFAULT_CVXPRODCOMP, NULL, NULL) );
2430
2431 SCIP_CALL( SCIPaddBoolParam(scip, "nlhdlr/" CONCAVE_NLHDLR_NAME "/handletrivial",
2432 "whether to also handle trivial convex expressions",
2433 &nlhdlrdata->handletrivial, TRUE, DEFAULT_HANDLETRIVIAL, NULL, NULL) );
2434
2435 SCIPnlhdlrSetFreeHdlrData(nlhdlr, nlhdlrfreeHdlrDataConvexConcave);
2436 SCIPnlhdlrSetCopyHdlr(nlhdlr, nlhdlrCopyhdlrConcave);
2437 SCIPnlhdlrSetFreeExprData(nlhdlr, nlhdlrfreeExprDataConvexConcave);
2438 SCIPnlhdlrSetSepa(nlhdlr, nlhdlrInitSepaConcave, NULL, nlhdlrEstimateConcave, NULL);
2439 SCIPnlhdlrSetInitExit(nlhdlr, NULL, nlhdlrExitConcave);
2440
2441 return SCIP_OKAY;
2442 }
2443
2444 /** checks whether a given expression is convex or concave w.r.t. the original variables
2445 *
2446 * This function uses the methods that are used in the detection algorithm of the convex nonlinear handler.
2447 */
2448 SCIP_RETCODE SCIPhasExprCurvature(
2449 SCIP* scip, /**< SCIP data structure */
2450 SCIP_EXPR* expr, /**< expression */
2451 SCIP_EXPRCURV curv, /**< curvature to check for */
2452 SCIP_Bool* success, /**< buffer to store whether expression has curvature curv (w.r.t. original variables) */
2453 SCIP_HASHMAP* assumevarfixed /**< hashmap containing variables that should be assumed to be fixed, or NULL */
2454 )
2455 {
2456 SCIP_NLHDLRDATA nlhdlrdata;
2457 SCIP_EXPR* rootnlexpr;
2458 SCIP_HASHMAP* nlexpr2origexpr;
2459 int nleafs;
2460
2461 assert(expr != NULL);
2462 assert(curv != SCIP_EXPRCURV_UNKNOWN);
2463 assert(success != NULL);
2464
2465 /* create temporary hashmap */
2466 SCIP_CALL( SCIPhashmapCreate(&nlexpr2origexpr, SCIPblkmem(scip), 20) );
2467
2468 /* prepare nonlinear handler data */
2469 nlhdlrdata.isnlhdlrconvex = TRUE;
2470 nlhdlrdata.evalsol = NULL;
2471 nlhdlrdata.detectsum = TRUE;
2472 nlhdlrdata.extendedform = FALSE;
2473 nlhdlrdata.cvxquadratic = TRUE;
2474 nlhdlrdata.cvxsignomial = TRUE;
2475 nlhdlrdata.cvxprodcomp = TRUE;
2476 nlhdlrdata.handletrivial = TRUE;
2477
2478 SCIP_CALL( constructExpr(scip, &nlhdlrdata, &rootnlexpr, nlexpr2origexpr, &nleafs, expr, curv, assumevarfixed, FALSE, success) );
2479
2480 /* free created expression */
2481 if( rootnlexpr != NULL )
2482 {
2483 SCIP_CALL( SCIPreleaseExpr(scip, &rootnlexpr) );
2484 }
2485
2486 /* free hashmap */
2487 SCIPhashmapFree(&nlexpr2origexpr);
2488
2489 return SCIP_OKAY;
2490 }
2491