Scippy

SCIP

Solving Constraint Integer Programs

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