Scippy

SCIP

Solving Constraint Integer Programs

expr_erf.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 (c) 2002-2023 Zuse Institute Berlin (ZIB) */
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 expr_erf.c
26  * @brief handler for Gaussian error function expressions
27  * @author Benjamin Mueller
28  */
29 
30 /*---+----1----+----2----+----3----+----4----+----5----+----6----+----7----+----8----+----9----+----0----+----1----+----2*/
31 
32 #include "scip/expr_erf.h"
33 #include "scip/expr_value.h"
34 
35 /* fundamental expression handler properties */
36 #define EXPRHDLR_NAME "erf"
37 #define EXPRHDLR_DESC "Gaussian error function"
38 #define EXPRHDLR_PRECEDENCE 79000
39 #define EXPRHDLR_HASHKEY SCIPcalcFibHash(131071.0)
40 
41 /*
42  * Data structures
43  */
44 
45 /*
46  * Local methods
47  */
48 
49 /** evaluates the Gaussian error function at a given point */
50 static
52  SCIP_Real x /**< point to evaluate */
53  )
54 {
55  SCIP_Real a1 = +0.254829592;
56  SCIP_Real a2 = -0.284496736;
57  SCIP_Real a3 = +1.421413741;
58  SCIP_Real a4 = -1.453152027;
59  SCIP_Real a5 = +1.061405429;
60  SCIP_Real p = +0.3275911;
61  int sign = (x >= 0.0) ? 1 : -1;
62  SCIP_Real t = 1.0 / (1.0 + p * REALABS(x));
63  SCIP_Real y = 1.0 - (((((a5 * t + a4) * t) + a3) * t + a2) * t + a1) * t * exp(-x*x);
64 
65  return sign*y;
66 }
67 
68 /*
69  * Callback methods of expression handler
70  */
71 
72 /** expression handler copy callback */
73 static
75 { /*lint --e{715}*/
77 
78  return SCIP_OKAY;
79 }
80 
81 /** simplifies an erf expression */
82 static
84 { /*lint --e{715}*/
85  SCIP_EXPR* child;
86 
87  assert(scip != NULL);
88  assert(expr != NULL);
89  assert(simplifiedexpr != NULL);
90  assert(SCIPexprGetNChildren(expr) == 1);
91 
92  child = SCIPexprGetChildren(expr)[0];
93  assert(child != NULL);
94 
95  /* check for value expression */
96  if( SCIPisExprValue(scip, child) )
97  {
98  SCIP_CALL( SCIPcreateExprValue(scip, simplifiedexpr, errorf(SCIPgetValueExprValue(child)), ownercreate, ownercreatedata) );
99  }
100  else
101  {
102  *simplifiedexpr = expr;
103 
104  /* we have to capture it, since it must simulate a "normal" simplified call in which a new expression is created */
105  SCIPcaptureExpr(*simplifiedexpr);
106  }
107 
108  return SCIP_OKAY;
109 }
110 
111 /** expression parse callback */
112 static
114 { /*lint --e{715}*/
115  SCIP_EXPR* childexpr;
116 
117  assert(expr != NULL);
118 
119  /* parse child expression from remaining string */
120  SCIP_CALL( SCIPparseExpr(scip, &childexpr, string, endstring, ownercreate, ownercreatedata) );
121  assert(childexpr != NULL);
122 
123  /* create gaussian error function expression */
124  SCIP_CALL( SCIPcreateExprErf(scip, expr, childexpr, ownercreate, ownercreatedata) );
125  assert(*expr != NULL);
126 
127  /* release child expression since it has been captured by the gaussian error function expression */
128  SCIP_CALL( SCIPreleaseExpr(scip, &childexpr) );
129 
130  *success = TRUE;
131 
132  return SCIP_OKAY;
133 }
134 
135 /** expression (point-) evaluation callback */
136 static
138 { /*lint --e{715}*/
139  assert(expr != NULL);
140  assert(SCIPexprGetData(expr) == NULL);
141  assert(SCIPexprGetNChildren(expr) == 1);
142  assert(SCIPexprGetEvalValue(SCIPexprGetChildren(expr)[0]) != SCIP_INVALID); /*lint !e777*/
143 
145 
146  return SCIP_OKAY;
147 }
148 
149 /** expression derivative evaluation callback */
150 static
152 { /*lint --e{715}*/
153  assert(expr != NULL);
154 
155  SCIPerrorMessage("method of erf expression handler not implemented yet\n");
156  SCIPABORT(); /*lint --e{527}*/
157 
158  return SCIP_OKAY;
159 }
160 
161 /** expression interval evaluation callback */
162 static
164 { /*lint --e{715}*/
165  SCIP_INTERVAL childinterval;
166 
167  assert(expr != NULL);
168  assert(SCIPexprGetData(expr) == NULL);
169  assert(SCIPexprGetNChildren(expr) == 1);
170 
171  childinterval = SCIPexprGetActivity(SCIPexprGetChildren(expr)[0]);
172 
173  if( SCIPintervalIsEmpty(SCIP_INTERVAL_INFINITY, childinterval) )
174  SCIPintervalSetEmpty(interval);
175  else
176  {
177  SCIP_Real childinf = SCIPintervalGetInf(childinterval);
178  SCIP_Real childsup = SCIPintervalGetSup(childinterval);
179  SCIP_Real inf = childinf <= -SCIP_INTERVAL_INFINITY ? -1.0 : errorf(childinf);
180  SCIP_Real sup = childsup >= +SCIP_INTERVAL_INFINITY ? +1.0 : errorf(childsup);
181  assert(inf <= sup);
182  SCIPintervalSetBounds(interval, inf, sup);
183  }
184 
185  return SCIP_OKAY;
186 }
187 
188 /** erf hash callback */
189 static
191 { /*lint --e{715}*/
192  assert(scip != NULL);
193  assert(expr != NULL);
194  assert(SCIPexprGetNChildren(expr) == 1);
195  assert(hashkey != NULL);
196  assert(childrenhashes != NULL);
197 
198  *hashkey = EXPRHDLR_HASHKEY;
199  *hashkey ^= childrenhashes[0];
200 
201  return SCIP_OKAY;
202 }
203 
204 /** expression curvature detection callback */
205 static
207 { /*lint --e{715}*/
208  assert(scip != NULL);
209  assert(expr != NULL);
210 
211  /* expression is
212  * - convex if child is convex and child <= 0
213  * - concave if child is concave and child >= 0
214  */
215  if( exprcurvature == SCIP_EXPRCURV_CONVEX )
216  {
217  *success = TRUE;
218  *childcurv = SCIP_EXPRCURV_CONVEX;
219  }
220  else
221  *success = FALSE;
222 
223  return SCIP_OKAY;
224 }
225 
226 /** expression monotonicity detection callback */
227 static
229 { /*lint --e{715}*/
230  assert(scip != NULL);
231  assert(expr != NULL);
232  assert(result != NULL);
233 
234  *result = SCIP_MONOTONE_INC;
235 
236  return SCIP_OKAY;
237 }
238 
239 /** expression integrality detection callback */
240 static
242 { /*lint --e{715}*/
243  assert(scip != NULL);
244  assert(expr != NULL);
245  assert(isintegral != NULL);
246 
247  *isintegral = FALSE;
248 
249  return SCIP_OKAY;
250 }
251 
252 /** creates an erf expression
253  *
254  * @attention The implementation of `erf` expressions is incomplete.
255  * They are not usable for most use cases so far.
256  */
258  SCIP* scip, /**< SCIP data structure */
259  SCIP_EXPR** expr, /**< pointer where to store expression */
260  SCIP_EXPR* child, /**< child expression */
261  SCIP_DECL_EXPR_OWNERCREATE((*ownercreate)), /**< function to call to create ownerdata */
262  void* ownercreatedata /**< data to pass to ownercreate */
263  )
264 {
265  SCIP_EXPRHDLR* exprhdlr;
266 
267  assert(expr != NULL);
268  assert(child != NULL);
269 
270  exprhdlr = SCIPfindExprhdlr(scip, EXPRHDLR_NAME);
271  if( exprhdlr == NULL )
272  {
273  SCIPerrorMessage("could not find %s expression handler -> abort\n", EXPRHDLR_NAME);
274  SCIPABORT();
275  return SCIP_PLUGINNOTFOUND;
276  }
277 
278  /* create expression */
279  SCIP_CALL( SCIPcreateExpr(scip, expr, exprhdlr, NULL, 1, &child, ownercreate, ownercreatedata) );
280 
281  return SCIP_OKAY;
282 }
283 
284 /** indicates whether expression is of erf-type */ /*lint -e{715}*/
286  SCIP* scip, /**< SCIP data structure */
287  SCIP_EXPR* expr /**< expression */
288  )
289 { /*lint --e{715}*/
290  assert(expr != NULL);
291 
292  return strcmp(SCIPexprhdlrGetName(SCIPexprGetHdlr(expr)), EXPRHDLR_NAME) == 0;
293 }
294 
295 /** creates the handler for erf expressions and includes it into SCIP
296  *
297  * @attention The implementation of this expression handler is incomplete.
298  * It is not usable for most use cases so far.
299  */
301  SCIP* scip /**< SCIP data structure */
302  )
303 {
304  SCIP_EXPRHDLR* exprhdlr;
305 
306  /* include expression handler */
308  assert(exprhdlr != NULL);
309 
310  SCIPexprhdlrSetCopyFreeHdlr(exprhdlr, copyhdlrErf, NULL);
311  SCIPexprhdlrSetSimplify(exprhdlr, simplifyErf);
312  SCIPexprhdlrSetParse(exprhdlr, parseErf);
313  SCIPexprhdlrSetIntEval(exprhdlr, intevalErf);
314  SCIPexprhdlrSetHash(exprhdlr, hashErf);
315  SCIPexprhdlrSetDiff(exprhdlr, bwdiffErf, NULL, NULL);
316  SCIPexprhdlrSetCurvature(exprhdlr, curvatureErf);
317  SCIPexprhdlrSetMonotonicity(exprhdlr, monotonicityErf);
318  SCIPexprhdlrSetIntegrality(exprhdlr, integralityErf);
319 
320  return SCIP_OKAY;
321 }
static SCIP_DECL_EXPRCOPYHDLR(copyhdlrErf)
Definition: expr_erf.c:74
int SCIPexprGetNChildren(SCIP_EXPR *expr)
Definition: expr.c:3801
#define EXPRHDLR_HASHKEY
Definition: expr_erf.c:39
void SCIPexprhdlrSetDiff(SCIP_EXPRHDLR *exprhdlr, SCIP_DECL_EXPRBWDIFF((*bwdiff)), SCIP_DECL_EXPRFWDIFF((*fwdiff)), SCIP_DECL_EXPRBWFWDIFF((*bwfwdiff)))
Definition: expr.c:471
const char * SCIPexprhdlrGetName(SCIP_EXPRHDLR *exprhdlr)
Definition: expr.c:532
#define EXPRHDLR_DESC
Definition: expr_erf.c:37
static SCIP_DECL_EXPRBWDIFF(bwdiffErf)
Definition: expr_erf.c:151
void SCIPexprhdlrSetParse(SCIP_EXPRHDLR *exprhdlr, SCIP_DECL_EXPRPARSE((*parse)))
Definition: expr.c:405
#define FALSE
Definition: def.h:96
#define TRUE
Definition: def.h:95
enum SCIP_Retcode SCIP_RETCODE
Definition: type_retcode.h:63
void SCIPexprhdlrSetHash(SCIP_EXPRHDLR *exprhdlr, SCIP_DECL_EXPRHASH((*hash)))
Definition: expr.c:449
SCIP_INTERVAL SCIPexprGetActivity(SCIP_EXPR *expr)
Definition: expr.c:3957
void SCIPexprhdlrSetIntegrality(SCIP_EXPRHDLR *exprhdlr, SCIP_DECL_EXPRINTEGRALITY((*integrality)))
Definition: expr.c:438
void SCIPcaptureExpr(SCIP_EXPR *expr)
Definition: scip_expr.c:1399
SCIP_EXPRHDLR * SCIPfindExprhdlr(SCIP *scip, const char *name)
Definition: scip_expr.c:859
SCIP_VAR ** x
Definition: circlepacking.c:63
SCIP_EXPRDATA * SCIPexprGetData(SCIP_EXPR *expr)
Definition: expr.c:3834
SCIP_EXPR ** SCIPexprGetChildren(SCIP_EXPR *expr)
Definition: expr.c:3811
SCIP_Real SCIPexprGetEvalValue(SCIP_EXPR *expr)
Definition: expr.c:3875
SCIP_RETCODE SCIPcreateExpr(SCIP *scip, SCIP_EXPR **expr, SCIP_EXPRHDLR *exprhdlr, SCIP_EXPRDATA *exprdata, int nchildren, SCIP_EXPR **children, SCIP_DECL_EXPR_OWNERCREATE((*ownercreate)), void *ownercreatedata)
Definition: scip_expr.c:964
#define SCIPerrorMessage
Definition: pub_message.h:64
SCIP_Bool SCIPintervalIsEmpty(SCIP_Real infinity, SCIP_INTERVAL operand)
void SCIPintervalSetEmpty(SCIP_INTERVAL *resultant)
SCIP_Bool SCIPisExprValue(SCIP *scip, SCIP_EXPR *expr)
Definition: scip_expr.c:1432
SCIP_Real SCIPintervalGetInf(SCIP_INTERVAL interval)
#define NULL
Definition: lpi_spx1.cpp:164
SCIP_Real SCIPintervalGetSup(SCIP_INTERVAL interval)
#define REALABS(x)
Definition: def.h:210
#define SCIP_CALL(x)
Definition: def.h:394
static SCIP_Real errorf(SCIP_Real x)
Definition: expr_erf.c:51
static SCIP_DECL_EXPRPARSE(parseErf)
Definition: expr_erf.c:113
SCIP_RETCODE SCIPcreateExprValue(SCIP *scip, SCIP_EXPR **expr, SCIP_Real value, SCIP_DECL_EXPR_OWNERCREATE((*ownercreate)), void *ownercreatedata)
Definition: expr_value.c:270
void SCIPexprhdlrSetCopyFreeHdlr(SCIP_EXPRHDLR *exprhdlr, SCIP_DECL_EXPRCOPYHDLR((*copyhdlr)), SCIP_DECL_EXPRFREEHDLR((*freehdlr)))
Definition: expr.c:368
#define SCIP_INTERVAL_INFINITY
Definition: def.h:208
#define EXPRHDLR_NAME
Definition: expr_erf.c:36
SCIP_RETCODE SCIPcreateExprErf(SCIP *scip, SCIP_EXPR **expr, SCIP_EXPR *child, SCIP_DECL_EXPR_OWNERCREATE((*ownercreate)), void *ownercreatedata)
Definition: expr_erf.c:257
#define SCIP_Bool
Definition: def.h:93
SCIP_Bool SCIPisExprErf(SCIP *scip, SCIP_EXPR *expr)
Definition: expr_erf.c:285
#define SCIP_DECL_EXPR_OWNERCREATE(x)
Definition: type_expr.h:140
void SCIPexprhdlrSetMonotonicity(SCIP_EXPRHDLR *exprhdlr, SCIP_DECL_EXPRMONOTONICITY((*monotonicity)))
Definition: expr.c:427
SCIP_RETCODE SCIPreleaseExpr(SCIP *scip, SCIP_EXPR **expr)
Definition: scip_expr.c:1407
SCIP_RETCODE SCIPincludeExprhdlrErf(SCIP *scip)
Definition: expr_erf.c:300
SCIP_RETCODE SCIPparseExpr(SCIP *scip, SCIP_EXPR **expr, const char *exprstr, const char **finalpos, SCIP_DECL_EXPR_OWNERCREATE((*ownercreate)), void *ownercreatedata)
Definition: scip_expr.c:1370
static SCIP_DECL_EXPRSIMPLIFY(simplifyErf)
Definition: expr_erf.c:83
SCIP_EXPRHDLR * SCIPexprGetHdlr(SCIP_EXPR *expr)
Definition: expr.c:3824
static SCIP_DECL_EXPRCURVATURE(curvatureErf)
Definition: expr_erf.c:206
constant value expression handler
static SCIP_DECL_EXPRMONOTONICITY(monotonicityErf)
Definition: expr_erf.c:228
void SCIPexprhdlrSetCurvature(SCIP_EXPRHDLR *exprhdlr, SCIP_DECL_EXPRCURVATURE((*curvature)))
Definition: expr.c:416
void SCIPintervalSetBounds(SCIP_INTERVAL *resultant, SCIP_Real inf, SCIP_Real sup)
#define SCIP_Real
Definition: def.h:186
handler for Gaussian error function expressions
static SCIP_DECL_EXPREVAL(evalErf)
Definition: expr_erf.c:137
SCIP_VAR ** y
Definition: circlepacking.c:64
static SCIP_DECL_EXPRHASH(hashErf)
Definition: expr_erf.c:190
#define SCIP_INVALID
Definition: def.h:206
SCIP_Real SCIPgetValueExprValue(SCIP_EXPR *expr)
Definition: expr_value.c:294
#define EXPRHDLR_PRECEDENCE
Definition: expr_erf.c:38
static SCIP_DECL_EXPRINTEGRALITY(integralityErf)
Definition: expr_erf.c:241
SCIP_RETCODE SCIPincludeExprhdlr(SCIP *scip, SCIP_EXPRHDLR **exprhdlr, const char *name, const char *desc, unsigned int precedence, SCIP_DECL_EXPREVAL((*eval)), SCIP_EXPRHDLRDATA *data)
Definition: scip_expr.c:814
static SCIP_DECL_EXPRINTEVAL(intevalErf)
Definition: expr_erf.c:163
#define SCIPABORT()
Definition: def.h:366
void SCIPexprhdlrSetSimplify(SCIP_EXPRHDLR *exprhdlr, SCIP_DECL_EXPRSIMPLIFY((*simplify)))
Definition: expr.c:497
void SCIPexprhdlrSetIntEval(SCIP_EXPRHDLR *exprhdlr, SCIP_DECL_EXPRINTEVAL((*inteval)))
Definition: expr.c:486