Scippy

SCIP

Solving Constraint Integer Programs

sepa_gauge.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-2019 Konrad-Zuse-Zentrum */
7 /* fuer Informationstechnik Berlin */
8 /* */
9 /* SCIP is distributed under the terms of the ZIB Academic License. */
10 /* */
11 /* You should have received a copy of the ZIB Academic License */
12 /* along with SCIP; see the file COPYING. If not visit scip.zib.de. */
13 /* */
14 /* * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * */
15 
16 /**@file sepa_gauge.c
17  * @brief gauge separator
18  * @author Felipe Serrano
19  *
20  * @todo should separator only be run when SCIPallColsInLP is true?
21  * @todo add SCIPisStopped(scip) to the condition of time consuming loops
22  * @todo check if it makes sense to implement the copy callback
23  */
24 
25 /*---+----1----+----2----+----3----+----4----+----5----+----6----+----7----+----8----+----9----+----0----+----1----+----2*/
26 
27 #include <assert.h>
28 #include <string.h>
29 
30 #include "blockmemshell/memory.h"
31 #include "nlpi/exprinterpret.h"
32 #include "nlpi/nlpi.h"
33 #include "nlpi/nlpi_ipopt.h"
34 #include "nlpi/nlpioracle.h"
35 #include "nlpi/pub_expr.h"
36 #include "scip/pub_lp.h"
37 #include "scip/pub_message.h"
38 #include "scip/pub_misc.h"
39 #include "scip/pub_nlp.h"
40 #include "scip/pub_sepa.h"
41 #include "scip/scip_cut.h"
42 #include "scip/scip_lp.h"
43 #include "scip/scip_mem.h"
44 #include "scip/scip_message.h"
45 #include "scip/scip_nlp.h"
46 #include "scip/scip_nonlinear.h"
47 #include "scip/scip_numerics.h"
48 #include "scip/scip_param.h"
49 #include "scip/scip_prob.h"
50 #include "scip/scip_sepa.h"
51 #include "scip/scip_sol.h"
52 #include "scip/scip_solvingstats.h"
53 #include "scip/scip_timing.h"
54 #include "scip/sepa_gauge.h"
55 #include <string.h>
56 
57 
58 #define SEPA_NAME "gauge"
59 #define SEPA_DESC "gauge separator"
60 #define SEPA_PRIORITY 0
61 #define SEPA_FREQ -1
62 #define SEPA_MAXBOUNDDIST 1.0
63 #define SEPA_USESSUBSCIP FALSE /**< does the separator use a secondary SCIP instance? */
64 #define SEPA_DELAY FALSE /**< should separation method be delayed, if other separators found cuts? */
65 
66 #define VIOLATIONFAC 100 /* constraints regarded as violated when violation > VIOLATIONFAC*SCIPfeastol */
67 #define MAX_ITER 75 /* maximum number of iterations for the line search */
68 
69 #define DEFAULT_NLPTIMELIMIT 0.0 /**< default time limit of NLP solver; 0.0 for no limit */
70 #define DEFAULT_NLPITERLIM 1000 /**< default NLP iteration limit */
71 
72 #define NLPFEASFAC 1e-1/**< NLP feasibility tolerance = NLPFEASFAC * SCIP's feasibility tolerance */
73 #define NLPVERBOSITY 0 /**< NLP solver verbosity */
74 
75 #define INTERIOROBJVARLB -100 /**< lower bound of the objective variable when computing interior point */
76 /*
77  * Data structures
78  */
79 
80 /** side that makes a nlrow convex */
82 {
83  LHS = 0, /**< left hand side */
84  RHS = 1 /**< right hand side */
85 };
86 typedef enum ConvexSide CONVEXSIDE;
87 
88 /** position of a point */
90 {
91  INTERIOR = 0, /**< point is in the interior of the region */
92  BOUNDARY = 1, /**< point is in the boundary of the region */
93  EXTERIOR = 2 /**< point is in the exterior of the region */
94 };
95 typedef enum Position POSITION;
96 
97 /** separator data */
98 struct SCIP_SepaData
99 {
100  SCIP_NLROW** nlrows; /**< stores convex nlrows */
101  CONVEXSIDE* convexsides; /**< which sides make the nlrows convex */
102  int* nlrowsidx; /**< indices of nlrows that violate the current lp solution */
103  int nnlrowsidx; /**< total number of convex nonlinear nlrows that violate the current lp solution */
104  int nnlrows; /**< total number of convex nonlinear nlrows */
105  int nlrowssize; /**< memory allocated for nlrows, convexsides and nlrowsidx */
106 
107  SCIP_Bool isintsolavailable; /**< do we have an interior point available? */
108  SCIP_Bool skipsepa; /**< whether separator should be skipped */
109  SCIP_SOL* intsol; /**< stores interior point */
110 
111  SCIP_EXPRINT* exprinterpreter; /**< expression interpreter to compute gradients */
112 
113  int ncuts; /**< number of cuts generated */
114 
115  /* parameters */
116  SCIP_Real nlptimelimit; /**< time limit of NLP solver; 0.0 for no limit */
117  int nlpiterlimit; /**< iteration limit of NLP solver; 0 for no limit */
118 };
119 
120 /*
121  * Local methods
122  */
123 
124 /** stores, from the constraints represented by nlrows, the nonlinear convex ones in sepadata */
125 static
127  SCIP* scip, /**< SCIP data structure */
128  SCIP_SEPADATA* sepadata, /**< separator data */
129  SCIP_NLROW** nlrows, /**< nlrows from which to store convex ones */
130  int nnlrows /**< number of nlrows */
131  )
132 {
133  int i;
134 
135  assert(scip != NULL);
136  assert(sepadata != NULL);
137  assert(nlrows != NULL);
138  assert(nnlrows > 0);
139 
140  SCIPdebugMsg(scip, "storing convex nlrows\n");
141 
142  SCIP_CALL( SCIPallocBlockMemoryArray(scip, &(sepadata->nlrows), nnlrows) );
143  SCIP_CALL( SCIPallocBlockMemoryArray(scip, &(sepadata->convexsides), nnlrows) );
144  SCIP_CALL( SCIPallocBlockMemoryArray(scip, &(sepadata->nlrowsidx), nnlrows) );
145  sepadata->nlrowssize = nnlrows;
146 
147  sepadata->nnlrows = 0;
148  for( i = 0; i < nnlrows; ++i )
149  {
150  SCIP_NLROW* nlrow;
151 
152  nlrow = nlrows[i];
153  assert(nlrow != NULL);
154 
155  /* linear case */
157  (SCIPnlrowGetNQuadElems(nlrow) == 0 && SCIPnlrowGetExprtree(nlrow) == NULL) )
158  continue;
159 
160  /* nonlinear case */
162  {
163  sepadata->convexsides[sepadata->nnlrows] = RHS;
164  sepadata->nlrows[sepadata->nnlrows] = nlrow;
165  ++(sepadata->nnlrows);
166  }
167  else if( !SCIPisInfinity(scip, -SCIPnlrowGetLhs(nlrow)) && SCIPnlrowGetCurvature(nlrow) == SCIP_EXPRCURV_CONCAVE )
168  {
169  sepadata->convexsides[sepadata->nnlrows] = LHS;
170  sepadata->nlrows[sepadata->nnlrows] = nlrow;
171  ++(sepadata->nnlrows);
172  }
173  }
174 
175  return SCIP_OKAY;
176 }
177 
178 /** computes an interior point of a convex NLP relaxation; builds the convex relaxation, modifies it to find an interior
179  * point, solves it and frees it; more details in @ref sepa_gauge.h
180  *
181  * @note the method also counts the number of nonlinear convex constraints and if there are < 2, then the convex
182  * relaxation is not interesting and the separator will not run again
183  */
184 static
186  SCIP* scip, /**< SCIP data structure */
187  SCIP_SEPADATA* sepadata /**< separator data */
188  )
189 {
190  SCIP_NLPIORACLE* nlpioracle;
191  SCIP_NLPIPROBLEM* nlpiprob;
192  SCIP_NLPI* nlpi;
193  SCIP_HASHMAP* var2nlpiidx;
194  SCIP_Real timelimit;
195  SCIP_Real objvarlb;
196  SCIP_Real minusone;
197  SCIP_Real one;
198  int nconvexnlrows;
199  int iterlimit;
200  int objvaridx;
201  int nconss;
202  int nvars;
203  int i;
204 
205  assert(scip != NULL);
206  assert(sepadata != NULL);
207  assert(!sepadata->skipsepa);
208 
209  SCIPdebugMsg(scip, "Computing interior point\n");
210 
211  /* create convex relaxation NLP */
212  assert(SCIPgetNNlpis(scip) > 0);
213 
214  nlpi = SCIPgetNlpis(scip)[0];
215  assert(nlpi != NULL);
216 
217  nvars = SCIPgetNVars(scip);
218  SCIP_CALL( SCIPnlpiCreateProblem(nlpi, &nlpiprob, "gauge-interiorpoint-nlp") );
219  SCIP_CALL( SCIPhashmapCreate(&var2nlpiidx, SCIPblkmem(scip), nvars) );
220  SCIP_CALL( SCIPcreateNlpiProb(scip, nlpi, SCIPgetNLPNlRows(scip), SCIPgetNNLPNlRows(scip), nlpiprob, var2nlpiidx,
221  NULL, SCIPgetCutoffbound(scip), FALSE, TRUE) );
222 
223  /* add objective variable; the problem is \min t, s.t. g(x) <= t, l(x) <= 0, where g are nonlinear and l linear */
224  objvaridx = nvars;
225  objvarlb = INTERIOROBJVARLB;
226  one = 1.0;
227  SCIP_CALL( SCIPnlpiAddVars(nlpi, nlpiprob, 1, &objvarlb, NULL, NULL) );
228  SCIP_CALL( SCIPnlpiSetObjective(nlpi, nlpiprob, 1, &objvaridx, &one, 0, NULL, NULL, NULL, 0.0) );
229 
230  /* add objective variables to constraints; for this we need to get nlpi oracle to have access to number of
231  * constraints and which constraints are nonlinear
232  */
233  /* @todo: this code is only valid when using IPOPT and needs to be changed when new NLP solvers get interfaced */
234  assert(strcmp(SCIPnlpiGetName(nlpi), "ipopt") == 0);
235  nlpioracle = (SCIP_NLPIORACLE *)SCIPgetNlpiOracleIpopt(nlpiprob);
236  assert(nlpioracle != NULL);
237  assert(SCIPnlpiOracleGetNVars(nlpioracle) == objvaridx + 1);
238 
239  minusone = -1.0;
240  nconvexnlrows = 0;
241  nconss = SCIPnlpiOracleGetNConstraints(nlpioracle);
242  for( i = 0; i < nconss; i++ )
243  {
244  if( SCIPnlpiOracleGetConstraintDegree(nlpioracle, i) > 1 )
245  {
246  SCIP_CALL( SCIPnlpiChgLinearCoefs(nlpi, nlpiprob, i, 1, &objvaridx, &minusone) );
247  ++nconvexnlrows;
248  }
249  }
251 
252  /* check if convex relaxation is interesting */
253  if( nconvexnlrows < 2 )
254  {
255  SCIPdebugMsg(scip, "convex relaxation is not interesting, only %d nonlinear convex rows; abort\n", nconvexnlrows);
256  sepadata->skipsepa = TRUE;
257  goto CLEANUP;
258  }
259 
260  /* add linear rows */
261  SCIP_CALL( SCIPaddNlpiProbRows(scip, nlpi, nlpiprob, var2nlpiidx, SCIPgetLPRows(scip), SCIPgetNLPRows(scip)) );
262 
263  /* set parameters in nlpi; time and iterations limit, tolerance, verbosity; for time limit, get time limit of scip;
264  * if scip doesn't have much time left, don't run separator. otherwise, timelimit is the minimum between whats left
265  * for scip and the timelimit setting
266  */
267  SCIP_CALL( SCIPgetRealParam(scip, "limits/time", &timelimit) );
268  if( !SCIPisInfinity(scip, timelimit) )
269  {
270  timelimit -= SCIPgetSolvingTime(scip);
271  if( timelimit <= 1.0 )
272  {
273  SCIPdebugMsg(scip, "skip NLP solve; no time left\n");
274  sepadata->skipsepa = TRUE;
275  goto CLEANUP;
276  }
277  }
278  if( sepadata->nlptimelimit > 0.0 )
279  timelimit = MIN(sepadata->nlptimelimit, timelimit);
280  SCIP_CALL( SCIPnlpiSetRealPar(nlpi, nlpiprob, SCIP_NLPPAR_TILIM, timelimit) );
281 
282  iterlimit = sepadata->nlpiterlimit > 0 ? sepadata->nlpiterlimit : INT_MAX;
283  SCIP_CALL( SCIPnlpiSetIntPar(nlpi, nlpiprob, SCIP_NLPPAR_ITLIM, iterlimit) );
285  SCIP_CALL( SCIPnlpiSetRealPar(nlpi, nlpiprob, SCIP_NLPPAR_RELOBJTOL, MAX(SCIPfeastol(scip), SCIPdualfeastol(scip))) ); /*lint !e666*/
287 
288  /* compute interior point */
289  SCIPdebugMsg(scip, "starting interior point computation\n");
290  SCIP_CALL( SCIPnlpiSolve(nlpi, nlpiprob) );
291  SCIPdebugMsg(scip, "finish interior point computation\n");
292 
293 #ifdef SCIP_DEBUG
294  {
295  SCIP_NLPSTATISTICS* nlpstatistics;
296 
297  /* get statistics */
298  SCIP_CALL( SCIPnlpStatisticsCreate(SCIPblkmem(scip), &nlpstatistics) );
299  SCIP_CALL( SCIPnlpiGetStatistics(nlpi, nlpiprob, nlpstatistics) );
300 
301  SCIPdebugMsg(scip, "nlpi took iters %d, time %g searching for an find interior point: solstat %d\n",
302  SCIPnlpStatisticsGetNIterations(nlpstatistics), SCIPnlpStatisticsGetTotalTime(nlpstatistics),
303  SCIPnlpiGetSolstat(nlpi, nlpiprob));
304 
305  SCIPnlpStatisticsFree(SCIPblkmem(scip), &nlpstatistics);
306  }
307 #endif
308 
309  if( SCIPnlpiGetSolstat(nlpi, nlpiprob) <= SCIP_NLPSOLSTAT_FEASIBLE )
310  {
311  SCIP_Real* nlpisol;
312 
313  SCIP_CALL( SCIPnlpiGetSolution(nlpi, nlpiprob, &nlpisol, NULL, NULL, NULL, NULL) );
314 
315  assert(nlpisol != NULL);
316  SCIPdebugMsg(scip, "NLP solved: sol found has objvalue = %g\n", nlpisol[objvaridx]);
317 
318  /* if we found an interior point store it */
319  if( SCIPisFeasNegative(scip, nlpisol[objvaridx]) )
320  {
321  SCIPdebugMsg(scip, "Interior point found!, storing it\n");
322  SCIP_CALL( SCIPcreateSol(scip, &sepadata->intsol, NULL) );
323  for( i = 0; i < nvars; i ++ )
324  {
325  SCIP_VAR* var;
326 
327  var = SCIPgetVars(scip)[i];
328  assert(SCIPhashmapExists(var2nlpiidx, (void*)var) );
329 
330  /* @todo: filter zero? */
331  SCIP_CALL( SCIPsetSolVal(scip, sepadata->intsol, var,
332  nlpisol[SCIPhashmapGetImageInt(var2nlpiidx, (void *)var)]) );
333  }
334 
335  sepadata->isintsolavailable = TRUE;
336  }
337  else
338  {
339  SCIPdebugMsg(scip, "We got a feasible point but not interior (objval: %g)\n", nlpisol[objvaridx]);
340  sepadata->skipsepa = TRUE;
341  }
342  }
343  else
344  {
345  SCIPdebugMsg(scip, "We couldn't get an interior point (stat: %d)\n", SCIPnlpiGetSolstat(nlpi, nlpiprob));
346  sepadata->skipsepa = TRUE;
347  }
348 
349 CLEANUP:
350  /* free memory */
351  SCIPhashmapFree(&var2nlpiidx);
352  SCIP_CALL( SCIPnlpiFreeProblem(nlpi, &nlpiprob) );
353 
354  return SCIP_OKAY;
355 }
356 
357 
358 /** find whether point is in the interior, at the boundary or in the exterior of the region described by the
359  * intersection of nlrows[i] <= rhs if convexsides[i] = RHS or lhs <= nlrows[i] if convexsides[i] = LHS
360  * @note: point corresponds to a convex combination between the lp solution and the interior point
361  */
362 static
364  SCIP* scip, /**< SCIP data structure */
365  SCIP_NLROW** nlrows, /**< nlrows defining the region */
366  int* nlrowsidx, /**< indices of nlrows defining the region */
367  int nnlrowsidx, /**< number of nlrows indices */
368  CONVEXSIDE* convexsides, /**< sides of the nlrows involved in the region */
369  SCIP_SOL* point, /**< point for which we want to know its position */
370  POSITION* position /**< buffer to store position of sol */
371  )
372 {
373  int i;
374 
375  assert(scip != NULL);
376  assert(nlrows != NULL);
377  assert(convexsides != NULL);
378  assert(nnlrowsidx > 0);
379  assert(point != NULL);
380  assert(position != NULL);
381 
382  *position = INTERIOR;
383  for( i = 0; i < nnlrowsidx; i++ )
384  {
385  SCIP_NLROW* nlrow;
386  SCIP_Real activity;
387  CONVEXSIDE convexside;
388 
389  nlrow = nlrows[nlrowsidx[i]];
390  convexside = convexsides[nlrowsidx[i]];
391 
392  /* compute activity of nlrow at point */
393  SCIP_CALL( SCIPgetNlRowSolActivity(scip, nlrow, point, &activity) );
394 
395  if( convexside == RHS )
396  {
397  assert(!SCIPisInfinity(scip, SCIPnlrowGetRhs(nlrow)));
398 
399  /* if nlrow <= rhs is violated, then we are in the exterior */
400  if( SCIPisFeasGT(scip, activity, SCIPnlrowGetRhs(nlrow)) )
401  {
402  *position = EXTERIOR;
403  SCIPdebugMsg(scip, "exterior because cons <%s> has activity %g. rhs: %g\n", SCIPnlrowGetName(nlrow),
404  activity, SCIPnlrowGetRhs(nlrow));
405  SCIPdebug( SCIPprintNlRow(scip, nlrow, NULL) );
406 
407  return SCIP_OKAY;
408  }
409 
410  /* if nlrow(point) == rhs, then we are currently at the boundary */
411  if( SCIPisFeasEQ(scip, activity, SCIPnlrowGetRhs(nlrow)) )
412  *position = BOUNDARY;
413  }
414  else
415  {
416  assert(!SCIPisInfinity(scip, -SCIPnlrowGetLhs(nlrow)));
417  assert(convexside == LHS);
418 
419  /* if lhs <= nlrow is violated, then we are in the exterior */
420  if( SCIPisFeasLT(scip, activity, SCIPnlrowGetLhs(nlrow)) )
421  {
422  *position = EXTERIOR;
423  return SCIP_OKAY;
424  }
425 
426  /* if lhs == nlrow(point), then we are currently at the boundary */
427  if( SCIPisFeasEQ(scip, activity, SCIPnlrowGetLhs(nlrow)) )
428  *position = BOUNDARY;
429  }
430  }
431 
432  return SCIP_OKAY;
433 }
434 
435 
436 /** returns, in convexcomb, the convex combination
437  * \f$ \lambda \f$ endpoint + (1 - \f$ lambda \f$) startpoint = startpoint + \f$ \lambda \f$ (endpoint - tosepasol) */
438 static
440  SCIP* scip, /**< SCIP data structure */
441  SCIP_Real lambda, /**< convex combination multiplier */
442  SCIP_SOL* startpoint, /**< point corresponding to \f$ \lambda = 0 \f$ */
443  SCIP_SOL* endpoint, /**< point corresponding to \f$ \lambda = 1 \f$ */
444  SCIP_SOL* convexcomb /**< solution to store convex combination of intsol and tosepasol */
445  )
446 {
447  SCIP_VAR** vars;
448  int nvars;
449  int i;
450 
451  assert(scip != NULL);
452  assert(startpoint != NULL);
453  assert(endpoint != NULL);
454  assert(convexcomb != NULL);
455 
456  vars = SCIPgetVars(scip);
457  nvars = SCIPgetNVars(scip);
458 
459  for( i = 0; i < nvars; i++ )
460  {
461  SCIP_Real val;
462  SCIP_VAR* var;
463 
464  var = vars[i];
465  val = lambda * SCIPgetSolVal(scip, endpoint, var) + (1.0 - lambda) * SCIPgetSolVal(scip, startpoint, var);
466 
467  if( !SCIPisZero(scip, val) )
468  {
469  SCIP_CALL( SCIPsetSolVal(scip, convexcomb, var, val) );
470  }
471  else
472  {
473  SCIP_CALL( SCIPsetSolVal(scip, convexcomb, var, 0.0) );
474  }
475  }
476 
477  return SCIP_OKAY;
478 }
479 
480 
481 /** performs binary search to find the point belonging to the segment [intsol, tosepasol] that intersects the boundary
482  * of the region described by the intersection of nlrows[i] <= rhs if convexsides[i] = RHS or lhs <= nlrows[i] if not,
483  * for i in nlrowsidx
484  */
485 static
487  SCIP* scip, /**< SCIP data structure */
488  SCIP_NLROW** nlrows, /**< nlrows defining the region */
489  int* nlrowsidx, /**< indices of nlrows defining the region */
490  int nnlrowsidx, /**< number of nlrows indices */
491  CONVEXSIDE* convexsides, /**< sides of the nlrows involved in the region */
492  SCIP_SOL* intsol, /**< point acting as 'interior point' */
493  SCIP_SOL* tosepasol, /**< solution that should be separated */
494  SCIP_SOL* sol, /**< convex combination of intsol and lpsol */
495  POSITION* position /**< buffer to store position of sol */
496  )
497 {
498  SCIP_Real lb;
499  SCIP_Real ub;
500  int i;
501 
502  assert(scip != NULL);
503  assert(nlrows != NULL);
504  assert(nlrowsidx != NULL);
505  assert(convexsides != NULL);
506  assert(intsol != NULL);
507  assert(tosepasol != NULL);
508  assert(sol != NULL);
509  assert(position != NULL);
510 
511  SCIPdebugMsg(scip, "starting binary search\n");
512  lb = 0.0; /* corresponds to intsol */
513  ub = 1.0; /* corresponds to tosepasol */
514  for( i = 0; i < MAX_ITER; i++ )
515  {
516  /* sol = (ub+lb)/2 * lpsol + (1 - (ub+lb)/2) * intsol */
517  SCIP_CALL( buildConvexCombination(scip, (ub + lb)/2.0, intsol, tosepasol, sol) );
518 
519  /* find poisition of point: boundary, interior, exterior */
520  SCIP_CALL( findPointPosition(scip, nlrows, nlrowsidx, nnlrowsidx, convexsides, sol, position) );
521  SCIPdebugMsg(scip, "Position: %d, lambda: %g\n", position, (ub + lb)/2.0);
522 
523  switch( *position )
524  {
525  case BOUNDARY:
526  SCIPdebugMsg(scip, "Done\n");
527  return SCIP_OKAY;
528 
529  case INTERIOR:
530  /* want to be closer to tosepasol */
531  lb = (ub + lb)/2.0;
532  break;
533 
534  case EXTERIOR:
535  /* want to be closer to intsol */
536  ub = (ub + lb)/2.0;
537  break;
538  }
539  }
540  SCIPdebugMsg(scip, "Done\n");
541  return SCIP_OKAY;
542 }
543 
544 
545 /** computes gradient of exprtree at sol */
546 static
548  SCIP* scip, /**< SCIP data structure */
549  SCIP_EXPRINT* exprint, /**< expressions interpreter */
550  SCIP_SOL* sol, /**< point where we compute gradient */
551  SCIP_EXPRTREE* exprtree, /**< exprtree for which we compute the gradient */
552  SCIP_Real* grad /**< buffer to store the gradient */
553  )
554 {
555  SCIP_Real* x;
556  SCIP_Real val;
557  int nvars;
558  int i;
559 
560  assert(scip != NULL);
561  assert(exprint != NULL);
562  assert(sol != NULL);
563  assert(exprtree != NULL);
564  assert(grad != NULL);
565 
566  nvars = SCIPexprtreeGetNVars(exprtree);
567  assert(nvars > 0);
568 
569  SCIP_CALL( SCIPallocBufferArray(scip, &x, nvars) );
570 
571  /* compile expression exprtree, if not done before */
572  if( SCIPexprtreeGetInterpreterData(exprtree) == NULL )
573  {
574  SCIP_CALL( SCIPexprintCompile(exprint, exprtree) );
575  }
576 
577  for( i = 0; i < nvars; ++i )
578  {
579  x[i] = SCIPgetSolVal(scip, sol, SCIPexprtreeGetVars(exprtree)[i]);
580  }
581 
582  SCIP_CALL( SCIPexprintGrad(exprint, exprtree, x, TRUE, &val, grad) );
583 
584  /*SCIPdebug( for( i = 0; i < nvars; ++i ) printf("%e [%s]\n", grad[i], SCIPvarGetName(SCIPexprtreeGetVars(exprtree)[i])) );*/
585 
586  SCIPfreeBufferArray(scip, &x);
587 
588  return SCIP_OKAY;
589 }
590 
591 
592 /** computes gradient cut (linearization) of nlrow at sol */
593 static
595  SCIP* scip, /**< SCIP data structure */
596  SCIP_SOL* sol, /**< point used to construct gradient cut (x_0) */
597  SCIP_EXPRINT* exprint, /**< expression interpreter */
598  SCIP_NLROW* nlrow, /**< constraint */
599  CONVEXSIDE convexside, /**< whether we use rhs or lhs of nlrow */
600  SCIP_ROW* row /**< storage for cut */
601  )
602 {
603  SCIP_Real activity;
604  SCIP_Real gradx0; /* <grad f(x_0), x_0> */
605  int i;
606 
607  assert(scip != NULL);
608  assert(exprint != NULL);
609  assert(nlrow != NULL);
610  assert(row != NULL);
611 
612  gradx0 = 0.0;
613 
614  /* an nlrow has a linear part, quadratic part and expression tree; ideally one would just build the gradient but we
615  * do not know if the different parts share variables or not, so we can't just build the gradient; for this reason
616  * we create the row right away and compute the gradients of each part independently and add them to the row; the
617  * row takes care to add coeffs corresponding to the same variable when they appear in different parts of the nlrow
618  */
619 
620  SCIP_CALL( SCIPcacheRowExtensions(scip, row) );
621 
622  /* linear part */
623  for( i = 0; i < SCIPnlrowGetNLinearVars(nlrow); i++ )
624  {
625  gradx0 += SCIPgetSolVal(scip, sol, SCIPnlrowGetLinearVars(nlrow)[i]) * SCIPnlrowGetLinearCoefs(nlrow)[i];
626  SCIP_CALL( SCIPaddVarToRow(scip, row, SCIPnlrowGetLinearVars(nlrow)[i], SCIPnlrowGetLinearCoefs(nlrow)[i]) );
627  }
628 
629  /* quadratic part */
630  for( i = 0; i < SCIPnlrowGetNQuadElems(nlrow); i++ )
631  {
632  SCIP_VAR* var1;
633  SCIP_VAR* var2;
634  SCIP_Real grad1;
635  SCIP_Real grad2;
636  SCIP_Real solval1;
637  SCIP_Real solval2;
638 
639  var1 = SCIPnlrowGetQuadVars(nlrow)[SCIPnlrowGetQuadElems(nlrow)[i].idx1];
640  var2 = SCIPnlrowGetQuadVars(nlrow)[SCIPnlrowGetQuadElems(nlrow)[i].idx2];
641  solval1 = SCIPgetSolVal(scip, sol, var1);
642  solval2 = SCIPgetSolVal(scip, sol, var2);
643  grad1 = SCIPnlrowGetQuadElems(nlrow)[i].coef * solval2; /* note that solval2 is correct */
644  grad2 = SCIPnlrowGetQuadElems(nlrow)[i].coef * solval1;
645 
646  SCIP_CALL( SCIPaddVarToRow(scip, row, var1, grad1) );
647  SCIP_CALL( SCIPaddVarToRow(scip, row, var2, grad2) );
648 
649  gradx0 += grad1 * solval1 + grad2 * solval2;
650  }
651 
652  /* expression tree part */
653  {
654  SCIP_Real* grad;
655  SCIP_EXPRTREE* tree;
656 
657  tree = SCIPnlrowGetExprtree(nlrow);
658 
659  if( tree != NULL && SCIPexprtreeGetNVars(tree) > 0 )
660  {
661  SCIP_CALL( SCIPallocBufferArray(scip, &grad, SCIPexprtreeGetNVars(tree)) );
662 
663  SCIP_CALL( computeGradient(scip, exprint, sol, tree, grad) );
664 
665  for( i = 0; i < SCIPexprtreeGetNVars(tree); i++ )
666  {
667  gradx0 += grad[i] * SCIPgetSolVal(scip, sol, SCIPexprtreeGetVars(tree)[i]);
668  SCIP_CALL( SCIPaddVarToRow(scip, row, SCIPexprtreeGetVars(tree)[i], grad[i]) );
669  }
670 
671  SCIPfreeBufferArray(scip, &grad);
672  }
673  }
674 
675  SCIP_CALL( SCIPflushRowExtensions(scip, row) );
676 
677 #ifdef CUT_DEBUG
678  SCIPdebugMsg(scip, "gradient: ");
679  SCIPdebug( SCIP_CALL( SCIPprintRow(scip, row, NULL) ) );
680  SCIPdebugMsg(scip, "gradient dot x_0: %g\n", gradx0);
681 #endif
682 
683  /* gradient cut is f(x_0) - <grad f(x_0), x_0> + <grad f(x_0), x> <= rhs or >= lhs */
684  SCIP_CALL( SCIPgetNlRowSolActivity(scip, nlrow, sol, &activity) );
685  if( convexside == RHS )
686  {
687  assert(!SCIPisInfinity(scip, SCIPnlrowGetRhs(nlrow)));
688  SCIP_CALL( SCIPchgRowRhs(scip, row, SCIPnlrowGetRhs(nlrow) - activity + gradx0) );
689  }
690  else
691  {
692  assert(convexside == LHS);
693  assert(!SCIPisInfinity(scip, -SCIPnlrowGetLhs(nlrow)));
694  SCIP_CALL( SCIPchgRowLhs(scip, row, SCIPnlrowGetLhs(nlrow) - activity + gradx0) );
695  }
696 
697 #ifdef CUT_DEBUG
698  SCIPdebugMsg(scip, "gradient cut: ");
699  SCIPdebug( SCIP_CALL( SCIPprintRow(scip, row, NULL) ) );
700 #endif
701 
702  return SCIP_OKAY;
703 }
704 
705 /** tries to generate gradient cuts at the point on the segment [intsol, tosepasol] that intersecs the boundary of the
706  * convex relaxation
707  * -# checks that the relative interior of the segment actually intersects the boundary
708  * (this check is needed since intsol is not necessarily an interior point)
709  * -# finds point on the boundary
710  * -# generates gradient cut at point on the boundary
711  */
712 static
714  SCIP* scip, /**< SCIP data structure */
715  SCIP_SEPA* sepa, /**< the cut separator itself */
716  SCIP_SOL* tosepasol, /**< solution that should be separated */
717  SCIP_RESULT* result /**< pointer to store the result of the separation call */
718  )
719 {
720  SCIP_SEPADATA* sepadata;
721  SCIP_NLROW** nlrows;
722  CONVEXSIDE* convexsides;
723  SCIP_SOL* sol;
724  SCIP_SOL* intsol;
725  POSITION position;
726  int* nlrowsidx;
727  int nnlrowsidx;
728  int i;
729 
730  assert(sepa != NULL);
731 
732  sepadata = SCIPsepaGetData(sepa);
733  assert(sepadata != NULL);
734 
735  intsol = sepadata->intsol;
736  nlrows = sepadata->nlrows;
737  nlrowsidx = sepadata->nlrowsidx;
738  nnlrowsidx = sepadata->nnlrowsidx;
739  convexsides = sepadata->convexsides;
740 
741  assert(intsol != NULL);
742  assert(nlrows != NULL);
743  assert(nlrowsidx != NULL);
744  assert(nnlrowsidx > 0);
745  assert(convexsides != NULL);
746 
747  /* to evaluate the nlrow one needs a solution */
748  SCIP_CALL( SCIPcreateSol(scip, &sol, NULL) );
749 
750  /* don't separate if, under SCIP tolerances, only a slight perturbation of the interior point in the direction of
751  * tosepasol gives a point that is in the exterior */
752  SCIP_CALL( buildConvexCombination(scip, VIOLATIONFAC * SCIPfeastol(scip), intsol, tosepasol, sol) );
753  SCIP_CALL( findPointPosition(scip, nlrows, nlrowsidx, nnlrowsidx, convexsides, sol, &position) );
754 
755  if( position == EXTERIOR )
756  {
757 #ifdef SCIP_DEBUG
758  SCIPdebugMsg(scip, "segment joining intsol and tosepasol seems to be contained in the exterior of the region, can't separate\n");
759  /* move from intsol in the direction of -tosepasol to check if we are really tangent to the region */
760  SCIP_CALL( buildConvexCombination(scip, -1e-3, intsol, tosepasol, sol) );
761  SCIP_CALL( findPointPosition(scip, nlrows, nlrowsidx, nnlrowsidx, convexsides, sol, &position) );
762  if( position == EXTERIOR )
763  {
764  SCIPdebugMsg(scip, "line through intsol and tosepasol is tangent to region; can't separate\n");
765  }
766  SCIP_CALL( findPointPosition(scip, nlrows, nlrowsidx, nnlrowsidx, convexsides, intsol, &position) );
767  printf("Position of intsol is %s\n",
768  position == EXTERIOR ? "exterior" : position == INTERIOR ? "interior": "boundary");
769  SCIP_CALL( findPointPosition(scip, nlrows, nlrowsidx, nnlrowsidx, convexsides, tosepasol, &position) );
770  printf("Position of tosepasol is %s\n",
771  position == EXTERIOR ? "exterior" : position == INTERIOR ? "interior": "boundary");
772 
773  /* slightly move from intsol in the direction of +-tosepasol */
774  SCIP_CALL( buildConvexCombination(scip, 1e-5, intsol, tosepasol, sol) );
775  SCIP_CALL( findPointPosition(scip, nlrows, nlrowsidx, nnlrowsidx, convexsides, sol, &position) );
776  printf("Position of intsol + 0.00001(tosepasol - inisol) is %s\n",
777  position == EXTERIOR ? "exterior" : position == INTERIOR ? "interior": "boundary");
778  SCIPdebug( SCIPprintSol(scip, sol, NULL, FALSE) );
779 
780  SCIP_CALL( buildConvexCombination(scip, -1e-5, intsol, tosepasol, sol) );
781  SCIP_CALL( findPointPosition(scip, nlrows, nlrowsidx, nnlrowsidx, convexsides, sol, &position) );
782  printf("Position of intsol - 0.00001(tosepasol - inisol) is %s\n",
783  position == EXTERIOR ? "exterior" : position == INTERIOR ? "interior": "boundary");
784  SCIPdebug( SCIPprintSol(scip, sol, NULL, FALSE) );
785 #endif
786  *result = SCIP_DIDNOTFIND;
787  goto CLEANUP;
788  }
789 
790  /* find point on boundary */
791  if( position != BOUNDARY )
792  {
793  SCIP_CALL( findBoundaryPoint(scip, nlrows, nlrowsidx, nnlrowsidx, convexsides, intsol, tosepasol, sol,
794  &position) );
795 
796  /* if MAX_ITER weren't enough to find a point in the boundary we don't separate */
797  if( position != BOUNDARY )
798  {
799  SCIPdebugMsg(scip, "couldn't find boundary point, don't separate\n");
800  goto CLEANUP;
801  }
802  }
803 
804  /* generate cuts at sol */
805  for( i = 0; i < nnlrowsidx; i++ )
806  {
807  SCIP_NLROW* nlrow;
808  SCIP_ROW* row;
809  SCIP_Real activity;
810  CONVEXSIDE convexside;
811  char rowname[SCIP_MAXSTRLEN];
812 
813  nlrow = nlrows[nlrowsidx[i]];
814  convexside = convexsides[nlrowsidx[i]];
815 
816  (void) SCIPsnprintf(rowname, SCIP_MAXSTRLEN, "%s_%u", SCIPnlrowGetName(nlrow), ++(sepadata->ncuts));
817 
818  /* only separate nlrows that are tight at the boundary point */
819  SCIP_CALL( SCIPgetNlRowSolActivity(scip, nlrow, sol, &activity) );
820  SCIPdebugMsg(scip, "cons <%s> at boundary point has activity: %g\n", SCIPnlrowGetName(nlrow), activity);
821 
822  if( (convexside == RHS && !SCIPisFeasEQ(scip, activity, SCIPnlrowGetRhs(nlrow)))
823  || (convexside == LHS && !SCIPisFeasEQ(scip, activity, SCIPnlrowGetLhs(nlrow))) )
824  continue;
825 
826  /* cut is globally valid, since we work on nlrows from the NLP built at the root node, which are globally valid */
827  /* @todo: when local nlrows get supported in SCIP, one can think of recomputing the interior point */
828  SCIP_CALL( SCIPcreateEmptyRowSepa(scip, &row, sepa, rowname, -SCIPinfinity(scip), SCIPinfinity(scip),
829  FALSE, FALSE , TRUE) );
830  SCIP_CALL( generateCut(scip, sol, sepadata->exprinterpreter, nlrow, convexside, row) );
831 
832  /* add cut */
833  SCIPdebugMsg(scip, "cut <%s> has efficacy %g\n", SCIProwGetName(row), SCIPgetCutEfficacy(scip, NULL, row));
834  if( SCIPisCutEfficacious(scip, NULL, row) )
835  {
836  SCIP_Bool infeasible;
837 
838  SCIPdebugMsg(scip, "adding cut\n");
839  SCIP_CALL( SCIPaddRow(scip, row, FALSE, &infeasible) );
840 
841  if( infeasible )
842  {
843  *result = SCIP_CUTOFF;
844  SCIP_CALL( SCIPreleaseRow(scip, &row) );
845  break;
846  }
847  else
848  {
849  *result = SCIP_SEPARATED;
850  }
851  }
852 
853  /* release the row */
854  SCIP_CALL( SCIPreleaseRow(scip, &row) );
855  }
856 
857 CLEANUP:
858  SCIP_CALL( SCIPfreeSol(scip, &sol) );
859 
860  return SCIP_OKAY;
861 }
862 
863 /*
864  * Callback methods of separator
865  */
866 
867 
868 /** destructor of separator to free user data (called when SCIP is exiting) */
869 static
870 SCIP_DECL_SEPAFREE(sepaFreeGauge)
871 { /*lint --e{715}*/
872  SCIP_SEPADATA* sepadata;
873 
874  assert(strcmp(SCIPsepaGetName(sepa), SEPA_NAME) == 0);
875 
876  /* free separator data */
877  sepadata = SCIPsepaGetData(sepa);
878  assert(sepadata != NULL);
879 
880  SCIPfreeBlockMemory(scip, &sepadata);
881 
882  SCIPsepaSetData(sepa, NULL);
883 
884  return SCIP_OKAY;
885 }
886 
887 
888 /** solving process deinitialization method of separator (called before branch and bound process data is freed) */
889 static
890 SCIP_DECL_SEPAEXITSOL(sepaExitsolGauge)
891 { /*lint --e{715}*/
892  SCIP_SEPADATA* sepadata;
893 
894  assert(sepa != NULL);
895 
896  sepadata = SCIPsepaGetData(sepa);
897 
898  assert(sepadata != NULL);
899 
900  /* free memory and reset data */
901  if( sepadata->isintsolavailable )
902  {
903  SCIPfreeBlockMemoryArray(scip, &sepadata->nlrowsidx, sepadata->nlrowssize);
904  SCIPfreeBlockMemoryArray(scip, &sepadata->convexsides, sepadata->nlrowssize);
905  SCIPfreeBlockMemoryArray(scip, &sepadata->nlrows, sepadata->nlrowssize);
906  SCIP_CALL( SCIPfreeSol(scip, &sepadata->intsol) );
907  SCIP_CALL( SCIPexprintFree(&sepadata->exprinterpreter) );
908 
909  sepadata->nnlrows = 0;
910  sepadata->nnlrowsidx = 0;
911  sepadata->nlrowssize = 0;
912  sepadata->isintsolavailable = FALSE;
913  }
914  assert(sepadata->nnlrows == 0);
915  assert(sepadata->nnlrowsidx == 0);
916  assert(sepadata->nlrowssize == 0);
917  assert(sepadata->isintsolavailable == FALSE);
918 
919  sepadata->skipsepa = FALSE;
920 
921  return SCIP_OKAY;
922 }
923 
924 
925 /** LP solution separation method of separator */
926 static
927 SCIP_DECL_SEPAEXECLP(sepaExeclpGauge)
928 { /*lint --e{715}*/
929  SCIP_SEPADATA* sepadata;
930  SCIP_SOL* lpsol;
931  int i;
932 
933  assert(scip != NULL);
934  assert(sepa != NULL);
935 
936  sepadata = SCIPsepaGetData(sepa);
937 
938  assert(sepadata != NULL);
939 
940  *result = SCIP_DIDNOTRUN;
941 
942  /* do not run if there is no interesting convex relaxation (with at least two nonlinear convex constraint) */
943  if( sepadata->skipsepa )
944  {
945  SCIPdebugMsg(scip, "not running because convex relaxation is uninteresting\n");
946  return SCIP_OKAY;
947  }
948 
949  /* do not run if SCIP has not constructed an NLP */
950  if( !SCIPisNLPConstructed(scip) )
951  {
952  SCIPdebugMsg(scip, "NLP not constructed, skipping gauge separator\n");
953  return SCIP_OKAY;
954  }
955 
956  /* do not run if SCIP has no way of solving nonlinear problems */
957  if( SCIPgetNNlpis(scip) == 0 )
958  {
959  SCIPdebugMsg(scip, "Skip gauge separator: no nlpi and SCIP can't solve nonlinear problems without a nlpi\n");
960  return SCIP_OKAY;
961  }
962 
963  /* if we don't have an interior point compute one; if we fail to compute one, then separator will not be run again;
964  * otherwise, we also store the convex nlrows in sepadata
965  */
966  if( !sepadata->isintsolavailable )
967  {
968  /* @todo: one could store the convex nonlinear rows inside computeInteriorPoint */
969  SCIP_CALL( computeInteriorPoint(scip, sepadata) );
970  assert(sepadata->skipsepa || sepadata->isintsolavailable);
971 
972  if( sepadata->skipsepa )
973  return SCIP_OKAY;
974 
976  SCIP_CALL( SCIPexprintCreate(SCIPblkmem(scip), &sepadata->exprinterpreter) );
977  }
978 
979 #ifdef SCIP_DISABLED_CODE
980  /* get interior point: try to compute an interior point, otherwise use primal solution, otherwise use NLP solution */
981  /* @todo: - decide order:
982  * - we can also use convex combination of solutions; there is a function SCIPvarGetAvgSol!
983  * - can add an event handler to only update when a new solution has been found
984  */
985  if( !sepadata->isintsolavailable && sepadata->computeintsol )
986  {
987  SCIP_Bool success;
988 
989  success = FALSE;
990  SCIP_CALL( computeInteriorPoint(scip, sepa, &success) );
991 
992  if( success )
993  sepadata->isintsolavailable = TRUE;
994  }
995 
996  if( !sepadata->isintsolavailable )
997  {
998  if( SCIPgetNSols(scip) > 0 )
999  {
1000  SCIPdebugMsg(scip, "Using current primal solution as interior point!\n");
1001  SCIP_CALL( SCIPcreateSolCopy(scip, &sepadata->intsol, SCIPgetBestSol(scip)) );
1002  sepadata->isintsolavailable = TRUE;
1003  }
1004  else if( SCIPhasNLPSolution(scip, NULL) )
1005  {
1006  SCIPdebugMsg(scip, "Using NLP solution as interior point!\n");
1007  SCIP_CALL( SCIPcreateNLPSol(scip, &sepadata->intsol, NULL) );
1008  sepadata->isintsolavailable = TRUE;
1009  }
1010  else
1011  {
1012  SCIPdebugMsg(scip, "We couldn't find an interior point, don't have a feasible nor an NLP solution; skip separator\n");
1013  return SCIP_OKAY;
1014  }
1015  }
1016 #endif
1017 
1018  /* store lp sol (or pseudo sol when lp is not solved) to be able to use it to compute nlrows' activities */
1019  SCIP_CALL( SCIPcreateCurrentSol(scip, &lpsol, NULL) );
1020 
1021  /* store indices of relevant constraints, ie, the ones that violate the lp sol */
1022  sepadata->nnlrowsidx = 0;
1023  for( i = 0; i < sepadata->nnlrows; i++ )
1024  {
1025  SCIP_NLROW* nlrow;
1026  SCIP_Real activity;
1027 
1028  nlrow = sepadata->nlrows[i];
1029 
1030  SCIP_CALL( SCIPgetNlRowSolActivity(scip, nlrow, lpsol, &activity) );
1031 
1032  if( sepadata->convexsides[i] == RHS )
1033  {
1034  assert(!SCIPisInfinity(scip, SCIPnlrowGetRhs(nlrow)));
1035 
1036  if( activity - SCIPnlrowGetRhs(nlrow) < VIOLATIONFAC * SCIPfeastol(scip) )
1037  continue;
1038  }
1039  else
1040  {
1041  assert(sepadata->convexsides[i] == LHS);
1042  assert(!SCIPisInfinity(scip, -SCIPnlrowGetLhs(nlrow)));
1043 
1044  if( SCIPnlrowGetLhs(nlrow) - activity < VIOLATIONFAC * SCIPfeastol(scip) )
1045  continue;
1046  }
1047 
1048  sepadata->nlrowsidx[sepadata->nnlrowsidx] = i;
1049  ++(sepadata->nnlrowsidx);
1050  }
1051 
1052  /* separate only if there are violated nlrows */
1053  SCIPdebugMsg(scip, "there are %d violated nlrows\n", sepadata->nnlrowsidx);
1054  if( sepadata->nnlrowsidx > 0 )
1055  {
1056  SCIP_CALL( separateCuts(scip, sepa, lpsol, result) );
1057  }
1058 
1059  /* free lpsol */
1060  SCIP_CALL( SCIPfreeSol(scip, &lpsol) );
1061 
1062  return SCIP_OKAY;
1063 }
1064 
1065 
1066 /*
1067  * separator specific interface methods
1068  */
1069 
1070 /** creates the gauge separator and includes it in SCIP */
1072  SCIP* scip /**< SCIP data structure */
1073  )
1074 {
1075  SCIP_SEPADATA* sepadata;
1076  SCIP_SEPA* sepa;
1077 
1078  /* create gauge separator data */
1079  SCIP_CALL( SCIPallocBlockMemory(scip, &sepadata) );
1080 
1081  /* this sets all data in sepadata to 0 */
1082  BMSclearMemory(sepadata);
1083 
1084  /* include separator */
1087  sepaExeclpGauge, NULL,
1088  sepadata) );
1089 
1090  assert(sepa != NULL);
1091 
1092  /* set non fundamental callbacks via setter functions */
1093  SCIP_CALL( SCIPsetSepaFree(scip, sepa, sepaFreeGauge) );
1094  SCIP_CALL( SCIPsetSepaExitsol(scip, sepa, sepaExitsolGauge) );
1095 
1096  /* add gauge separator parameters */
1097  SCIP_CALL( SCIPaddIntParam(scip, "separating/" SEPA_NAME "/nlpiterlimit",
1098  "iteration limit of NLP solver; 0 for no limit",
1099  &sepadata->nlpiterlimit, TRUE, DEFAULT_NLPITERLIM, 0, INT_MAX, NULL, NULL) );
1100 
1101  SCIP_CALL( SCIPaddRealParam(scip, "separating/" SEPA_NAME "/nlptimelimit",
1102  "time limit of NLP solver; 0.0 for no limit",
1103  &sepadata->nlptimelimit, TRUE, DEFAULT_NLPTIMELIMIT, 0.0, SCIP_REAL_MAX, NULL, NULL) );
1104 
1105  return SCIP_OKAY;
1106 }
enum SCIP_Result SCIP_RESULT
Definition: type_result.h:52
#define SEPA_DESC
Definition: sepa_gauge.c:59
#define SCIPfreeBlockMemoryArray(scip, ptr, num)
Definition: scip_mem.h:97
SCIP_Real SCIPfeastol(SCIP *scip)
SCIP_RETCODE SCIPfreeSol(SCIP *scip, SCIP_SOL **sol)
Definition: scip_sol.c:976
SCIP_Real SCIPdualfeastol(SCIP *scip)
SCIP_RETCODE SCIPnlpiAddVars(SCIP_NLPI *nlpi, SCIP_NLPIPROBLEM *problem, int nvars, const SCIP_Real *lbs, const SCIP_Real *ubs, const char **varnames)
Definition: nlpi.c:250
#define NULL
Definition: def.h:253
SCIP_EXPRTREE * SCIPnlrowGetExprtree(SCIP_NLROW *nlrow)
Definition: nlp.c:3364
#define SCIPallocBlockMemoryArray(scip, ptr, num)
Definition: scip_mem.h:80
SCIP_RETCODE SCIPsetSolVal(SCIP *scip, SCIP_SOL *sol, SCIP_VAR *var, SCIP_Real val)
Definition: scip_sol.c:1212
public methods for SCIP parameter handling
#define NLPVERBOSITY
Definition: sepa_gauge.c:73
methods to interpret (evaluate) an expression tree "fast"
void * SCIPgetNlpiOracleIpopt(SCIP_NLPIPROBLEM *nlpiproblem)
int SCIPgetNLPRows(SCIP *scip)
Definition: scip_lp.c:561
public methods for memory management
int SCIPnlpiOracleGetNVars(SCIP_NLPIORACLE *oracle)
Definition: nlpioracle.c:2202
SCIP_NLPI ** SCIPgetNlpis(SCIP *scip)
Definition: scip_nlp.c:118
#define SCIP_MAXSTRLEN
Definition: def.h:274
SCIP_RETCODE SCIPprintNlRow(SCIP *scip, SCIP_NLROW *nlrow, FILE *file)
Definition: scip_nlp.c:2004
SCIP_Real SCIPgetSolVal(SCIP *scip, SCIP_SOL *sol, SCIP_VAR *var)
Definition: scip_sol.c:1352
internal methods for NLPI solver interfaces
SCIP_RETCODE SCIPprintRow(SCIP *scip, SCIP_ROW *row, FILE *file)
Definition: scip_lp.c:2031
public methods for timing
methods to store an NLP and request function, gradient, and hessian values
SCIP_RETCODE SCIPnlpiCreateProblem(SCIP_NLPI *nlpi, SCIP_NLPIPROBLEM **problem, const char *name)
Definition: nlpi.c:211
static SCIP_RETCODE computeInteriorPoint(SCIP *scip, SCIP_SEPADATA *sepadata)
Definition: sepa_gauge.c:185
int SCIPgetNVars(SCIP *scip)
Definition: scip_prob.c:1987
SCIP_Real SCIPgetSolvingTime(SCIP *scip)
Definition: scip_timing.c:359
#define FALSE
Definition: def.h:73
gauge separator
static SCIP_DECL_SEPAFREE(sepaFreeGauge)
Definition: sepa_gauge.c:870
SCIP_Bool SCIPisFeasNegative(SCIP *scip, SCIP_Real val)
#define TRUE
Definition: def.h:72
#define SCIPdebug(x)
Definition: pub_message.h:74
enum Position POSITION
Definition: sepa_gauge.c:95
enum SCIP_Retcode SCIP_RETCODE
Definition: type_retcode.h:53
SCIP_EXPORT SCIP_RETCODE SCIPexprintGrad(SCIP_EXPRINT *exprint, SCIP_EXPRTREE *tree, SCIP_Real *varvals, SCIP_Bool new_varvals, SCIP_Real *val, SCIP_Real *gradient)
enum ConvexSide CONVEXSIDE
Definition: sepa_gauge.c:86
#define SCIPfreeBlockMemory(scip, ptr)
Definition: scip_mem.h:95
#define DEFAULT_NLPITERLIM
Definition: sepa_gauge.c:70
SCIP_RETCODE SCIPflushRowExtensions(SCIP *scip, SCIP_ROW *row)
Definition: scip_lp.c:1502
SCIP_RETCODE SCIPnlpiSetRealPar(SCIP_NLPI *nlpi, SCIP_NLPIPROBLEM *problem, SCIP_NLPPARAM type, SCIP_Real dval)
Definition: nlpi.c:671
#define SCIPfreeBufferArray(scip, ptr)
Definition: scip_mem.h:123
#define SCIPallocBlockMemory(scip, ptr)
Definition: scip_mem.h:78
#define SCIPdebugMsg
Definition: scip_message.h:69
SCIP_RETCODE SCIPcreateNLPSol(SCIP *scip, SCIP_SOL **sol, SCIP_HEUR *heur)
Definition: scip_sol.c:389
public methods for separator plugins
SCIP_Real SCIPgetCutoffbound(SCIP *scip)
SCIP_VAR ** x
Definition: circlepacking.c:54
int SCIPnlpiOracleGetNConstraints(SCIP_NLPIORACLE *oracle)
Definition: nlpioracle.c:2212
SCIP_Bool SCIPisInfinity(SCIP *scip, SCIP_Real val)
SCIP_Real SCIPgetCutEfficacy(SCIP *scip, SCIP_SOL *sol, SCIP_ROW *cut)
Definition: scip_cut.c:64
SCIP_ROW ** SCIPgetLPRows(SCIP *scip)
Definition: scip_lp.c:540
public methods for numerical tolerances
SCIP_RETCODE SCIPnlpiSetObjective(SCIP_NLPI *nlpi, SCIP_NLPIPROBLEM *problem, int nlins, const int *lininds, const SCIP_Real *linvals, int nquadelems, const SCIP_QUADELEM *quadelems, const int *exprvaridxs, const SCIP_EXPRTREE *exprtree, const SCIP_Real constant)
Definition: nlpi.c:300
SCIP_RETCODE SCIPcreateSol(SCIP *scip, SCIP_SOL **sol, SCIP_HEUR *heur)
Definition: scip_sol.c:319
#define SEPA_MAXBOUNDDIST
Definition: sepa_gauge.c:62
public methods for expressions, expression trees, expression graphs, and related stuff ...
static SCIP_RETCODE computeGradient(SCIP *scip, SCIP_EXPRINT *exprint, SCIP_SOL *sol, SCIP_EXPRTREE *exprtree, SCIP_Real *grad)
Definition: sepa_gauge.c:547
public methods for querying solving statistics
int SCIPgetNSols(SCIP *scip)
Definition: scip_sol.c:2205
int SCIPhashmapGetImageInt(SCIP_HASHMAP *hashmap, void *origin)
Definition: misc.c:3098
Definition: sepa_gauge.c:84
static SCIP_RETCODE generateCut(SCIP *scip, SCIP_SOL *sol, SCIP_EXPRINT *exprint, SCIP_NLROW *nlrow, CONVEXSIDE convexside, SCIP_ROW *row)
Definition: sepa_gauge.c:594
static SCIP_DECL_SEPAEXITSOL(sepaExitsolGauge)
Definition: sepa_gauge.c:890
SCIP_Bool SCIPhashmapExists(SCIP_HASHMAP *hashmap, void *origin)
Definition: misc.c:3240
SCIP_RETCODE SCIPcreateSolCopy(SCIP *scip, SCIP_SOL **sol, SCIP_SOL *sourcesol)
Definition: scip_sol.c:609
SCIP_Real coef
Definition: type_expr.h:104
SCIP_Bool SCIPhasNLPSolution(SCIP *scip)
Definition: scip_nlp.c:687
SCIP_RETCODE SCIPcacheRowExtensions(SCIP *scip, SCIP_ROW *row)
Definition: scip_lp.c:1479
SCIP_RETCODE SCIPnlpiGetStatistics(SCIP_NLPI *nlpi, SCIP_NLPIPROBLEM *problem, SCIP_NLPSTATISTICS *statistics)
Definition: nlpi.c:556
int SCIPnlpiOracleGetConstraintDegree(SCIP_NLPIORACLE *oracle, int considx)
Definition: nlpioracle.c:2324
SCIP_NLROW ** SCIPgetNLPNlRows(SCIP *scip)
Definition: scip_nlp.c:416
void SCIPnlpStatisticsFree(BMS_BLKMEM *blkmem, SCIP_NLPSTATISTICS **statistics)
Definition: nlpi.c:801
static SCIP_RETCODE storeNonlinearConvexNlrows(SCIP *scip, SCIP_SEPADATA *sepadata, SCIP_NLROW **nlrows, int nnlrows)
Definition: sepa_gauge.c:126
SCIP_RETCODE SCIPnlpiGetSolution(SCIP_NLPI *nlpi, SCIP_NLPIPROBLEM *problem, SCIP_Real **primalvalues, SCIP_Real **consdualvalues, SCIP_Real **varlbdualvalues, SCIP_Real **varubdualvalues, SCIP_Real *objval)
Definition: nlpi.c:537
SCIP_RETCODE SCIPreleaseRow(SCIP *scip, SCIP_ROW **row)
Definition: scip_lp.c:1406
SCIP_Bool SCIPisFeasEQ(SCIP *scip, SCIP_Real val1, SCIP_Real val2)
public methods for nonlinear functions
SCIP_VAR ** SCIPgetVars(SCIP *scip)
Definition: scip_prob.c:1942
BMS_BLKMEM * SCIPblkmem(SCIP *scip)
Definition: scip_mem.c:47
static SCIP_DECL_SEPAEXECLP(sepaExeclpGauge)
Definition: sepa_gauge.c:927
SCIP_Bool SCIPisZero(SCIP *scip, SCIP_Real val)
SCIP_NLPSOLSTAT SCIPnlpiGetSolstat(SCIP_NLPI *nlpi, SCIP_NLPIPROBLEM *problem)
Definition: nlpi.c:511
int SCIPnlrowGetNQuadElems(SCIP_NLROW *nlrow)
Definition: nlp.c:3323
SCIP_EXPORT SCIP_RETCODE SCIPexprintFree(SCIP_EXPRINT **exprint)
#define SEPA_NAME
Definition: sepa_gauge.c:58
SCIP_RETCODE SCIPcreateCurrentSol(SCIP *scip, SCIP_SOL **sol, SCIP_HEUR *heur)
Definition: scip_sol.c:474
int SCIPexprtreeGetNVars(SCIP_EXPRTREE *tree)
Definition: expr.c:8613
SCIP_RETCODE SCIPcreateNlpiProb(SCIP *scip, SCIP_NLPI *nlpi, SCIP_NLROW **nlrows, int nnlrows, SCIP_NLPIPROBLEM *nlpiprob, SCIP_HASHMAP *var2idx, SCIP_Real *nlscore, SCIP_Real cutoffbound, SCIP_Bool setobj, SCIP_Bool onlyconvex)
#define NLPFEASFAC
Definition: sepa_gauge.c:72
#define SCIP_CALL(x)
Definition: def.h:365
int SCIPnlrowGetNLinearVars(SCIP_NLROW *nlrow)
Definition: nlp.c:3246
int SCIPgetNNLPNlRows(SCIP *scip)
Definition: scip_nlp.c:438
SCIP_EXPORT SCIP_SEPADATA * SCIPsepaGetData(SCIP_SEPA *sepa)
Definition: sepa.c:607
SCIP_EXPORT const char * SCIPsepaGetName(SCIP_SEPA *sepa)
Definition: sepa.c:696
SCIP_QUADELEM * SCIPnlrowGetQuadElems(SCIP_NLROW *nlrow)
Definition: nlp.c:3333
SCIP_EXPORT SCIP_RETCODE SCIPexprintCreate(BMS_BLKMEM *blkmem, SCIP_EXPRINT **exprint)
public methods for NLP management
SCIP_RETCODE SCIPgetRealParam(SCIP *scip, const char *name, SCIP_Real *value)
Definition: scip_param.c:297
SCIP_EXPRCURV SCIPnlrowGetCurvature(SCIP_NLROW *nlrow)
Definition: nlp.c:3394
Ipopt NLP interface.
#define SCIPallocBufferArray(scip, ptr, num)
Definition: scip_mem.h:111
SCIP_Real SCIPinfinity(SCIP *scip)
public data structures and miscellaneous methods
SCIP_RETCODE SCIPnlpiOraclePrintProblem(SCIP_NLPIORACLE *oracle, SCIP_MESSAGEHDLR *messagehdlr, FILE *file)
Definition: nlpioracle.c:2939
SCIP_Real SCIPnlrowGetRhs(SCIP_NLROW *nlrow)
Definition: nlp.c:3384
#define SCIP_Bool
Definition: def.h:70
SCIP_RETCODE SCIPhashmapCreate(SCIP_HASHMAP **hashmap, BMS_BLKMEM *blkmem, int mapsize)
Definition: misc.c:2891
SCIP_Bool SCIPisNLPConstructed(SCIP *scip)
Definition: scip_nlp.c:209
SCIP_RETCODE SCIPincludeSepaBasic(SCIP *scip, SCIP_SEPA **sepa, const char *name, const char *desc, int priority, int freq, SCIP_Real maxbounddist, SCIP_Bool usessubscip, SCIP_Bool delay, SCIP_DECL_SEPAEXECLP((*sepaexeclp)), SCIP_DECL_SEPAEXECSOL((*sepaexecsol)), SCIP_SEPADATA *sepadata)
Definition: scip_sepa.c:99
int SCIPgetNNlpis(SCIP *scip)
Definition: scip_nlp.c:131
int SCIPnlpStatisticsGetNIterations(SCIP_NLPSTATISTICS *statistics)
Definition: nlpi.c:816
#define DEFAULT_NLPTIMELIMIT
Definition: sepa_gauge.c:69
SCIP_EXPRINTDATA * SCIPexprtreeGetInterpreterData(SCIP_EXPRTREE *tree)
Definition: expr.c:8658
SCIP_RETCODE SCIPnlpiSetIntPar(SCIP_NLPI *nlpi, SCIP_NLPIPROBLEM *problem, SCIP_NLPPARAM type, int ival)
Definition: nlpi.c:636
SCIP_MESSAGEHDLR * SCIPgetMessagehdlr(SCIP *scip)
Definition: scip_message.c:90
#define MIN(x, y)
Definition: def.h:223
public methods for LP management
#define SEPA_DELAY
Definition: sepa_gauge.c:64
SCIP_RETCODE SCIPincludeSepaGauge(SCIP *scip)
Definition: sepa_gauge.c:1071
public methods for cuts and aggregation rows
SCIP_VAR ** SCIPexprtreeGetVars(SCIP_EXPRTREE *tree)
Definition: nlp.c:102
SCIP_Bool SCIPisCutEfficacious(SCIP *scip, SCIP_SOL *sol, SCIP_ROW *cut)
Definition: scip_cut.c:87
SCIP_RETCODE SCIPnlpiSolve(SCIP_NLPI *nlpi, SCIP_NLPIPROBLEM *problem)
Definition: nlpi.c:497
SCIP_Real SCIPnlrowGetLhs(SCIP_NLROW *nlrow)
Definition: nlp.c:3374
#define SEPA_PRIORITY
Definition: sepa_gauge.c:60
#define BMSclearMemory(ptr)
Definition: memory.h:119
SCIP_Real * SCIPnlrowGetLinearCoefs(SCIP_NLROW *nlrow)
Definition: nlp.c:3266
SCIP_RETCODE SCIPaddRealParam(SCIP *scip, const char *name, const char *desc, SCIP_Real *valueptr, SCIP_Bool isadvanced, SCIP_Real defaultvalue, SCIP_Real minvalue, SCIP_Real maxvalue, SCIP_DECL_PARAMCHGD((*paramchgd)), SCIP_PARAMDATA *paramdata)
Definition: scip_param.c:129
public methods for the LP relaxation, rows and columns
SCIP_RETCODE SCIPnlpiChgLinearCoefs(SCIP_NLPI *nlpi, SCIP_NLPIPROBLEM *problem, const int idx, int nvals, const int *varidxs, const SCIP_Real *vals)
Definition: nlpi.c:395
#define SCIP_REAL_MAX
Definition: def.h:165
public methods for nonlinear relaxations
SCIP_Real SCIPnlpStatisticsGetTotalTime(SCIP_NLPSTATISTICS *statistics)
Definition: nlpi.c:826
SCIP_RETCODE SCIPaddRow(SCIP *scip, SCIP_ROW *row, SCIP_Bool forcecut, SCIP_Bool *infeasible)
Definition: scip_cut.c:220
#define MAX(x, y)
Definition: def.h:222
Position
Definition: sepa_gauge.c:89
SCIP_Bool SCIPisFeasGT(SCIP *scip, SCIP_Real val1, SCIP_Real val2)
public methods for solutions
SCIP_RETCODE SCIPchgRowRhs(SCIP *scip, SCIP_ROW *row, SCIP_Real rhs)
Definition: scip_lp.c:1451
SCIP_RETCODE SCIPchgRowLhs(SCIP *scip, SCIP_ROW *row, SCIP_Real lhs)
Definition: scip_lp.c:1427
Definition: sepa_gauge.c:83
SCIP_Bool SCIPisFeasLT(SCIP *scip, SCIP_Real val1, SCIP_Real val2)
#define SEPA_FREQ
Definition: sepa_gauge.c:61
SCIP_RETCODE SCIPaddVarToRow(SCIP *scip, SCIP_ROW *row, SCIP_VAR *var, SCIP_Real val)
Definition: scip_lp.c:1539
SCIP_RETCODE SCIPnlpiFreeProblem(SCIP_NLPI *nlpi, SCIP_NLPIPROBLEM **problem)
Definition: nlpi.c:224
public methods for message output
int SCIPsnprintf(char *t, int len, const char *s,...)
Definition: misc.c:10263
SCIP_RETCODE SCIPaddIntParam(SCIP *scip, const char *name, const char *desc, int *valueptr, SCIP_Bool isadvanced, int defaultvalue, int minvalue, int maxvalue, SCIP_DECL_PARAMCHGD((*paramchgd)), SCIP_PARAMDATA *paramdata)
Definition: scip_param.c:73
SCIP_RETCODE SCIPprintSol(SCIP *scip, SCIP_SOL *sol, FILE *file, SCIP_Bool printzeros)
Definition: scip_sol.c:1766
const char * SCIPnlrowGetName(SCIP_NLROW *nlrow)
Definition: nlp.c:3413
void SCIPhashmapFree(SCIP_HASHMAP **hashmap)
Definition: misc.c:2925
#define SCIP_Real
Definition: def.h:164
static SCIP_RETCODE buildConvexCombination(SCIP *scip, SCIP_Real lambda, SCIP_SOL *startpoint, SCIP_SOL *endpoint, SCIP_SOL *convexcomb)
Definition: sepa_gauge.c:439
const char * SCIProwGetName(SCIP_ROW *row)
Definition: lp.c:17025
public methods for message handling
static SCIP_RETCODE separateCuts(SCIP *scip, SCIP_SEPA *sepa, SCIP_SOL *tosepasol, SCIP_RESULT *result)
Definition: sepa_gauge.c:713
SCIP_VAR ** SCIPnlrowGetQuadVars(SCIP_NLROW *nlrow)
Definition: nlp.c:3286
SCIP_RETCODE SCIPsetSepaExitsol(SCIP *scip, SCIP_SEPA *sepa, SCIP_DECL_SEPAEXITSOL((*sepaexitsol)))
Definition: scip_sepa.c:221
SCIP_RETCODE SCIPnlpStatisticsCreate(BMS_BLKMEM *blkmem, SCIP_NLPSTATISTICS **statistics)
Definition: nlpi.c:784
SCIP_EXPORT void SCIPsepaSetData(SCIP_SEPA *sepa, SCIP_SEPADATA *sepadata)
Definition: sepa.c:617
#define INTERIOROBJVARLB
Definition: sepa_gauge.c:75
ConvexSide
SCIP_VAR ** SCIPnlrowGetLinearVars(SCIP_NLROW *nlrow)
Definition: nlp.c:3256
public methods for separators
SCIP_RETCODE SCIPsetSepaFree(SCIP *scip, SCIP_SEPA *sepa, SCIP_DECL_SEPAFREE((*sepafree)))
Definition: scip_sepa.c:157
static SCIP_RETCODE findPointPosition(SCIP *scip, SCIP_NLROW **nlrows, int *nlrowsidx, int nnlrowsidx, CONVEXSIDE *convexsides, SCIP_SOL *point, POSITION *position)
Definition: sepa_gauge.c:363
#define SEPA_USESSUBSCIP
Definition: sepa_gauge.c:63
public methods for global and local (sub)problems
static SCIP_RETCODE findBoundaryPoint(SCIP *scip, SCIP_NLROW **nlrows, int *nlrowsidx, int nnlrowsidx, CONVEXSIDE *convexsides, SCIP_SOL *intsol, SCIP_SOL *tosepasol, SCIP_SOL *sol, POSITION *position)
Definition: sepa_gauge.c:486
SCIP_SOL * SCIPgetBestSol(SCIP *scip)
Definition: scip_sol.c:2304
#define VIOLATIONFAC
Definition: sepa_gauge.c:66
enum ConvexSide CONVEXSIDE
SCIP_RETCODE SCIPaddNlpiProbRows(SCIP *scip, SCIP_NLPI *nlpi, SCIP_NLPIPROBLEM *nlpiprob, SCIP_HASHMAP *var2idx, SCIP_ROW **rows, int nrows)
struct SCIP_SepaData SCIP_SEPADATA
Definition: type_sepa.h:38
SCIP_EXPORT SCIP_RETCODE SCIPexprintCompile(SCIP_EXPRINT *exprint, SCIP_EXPRTREE *tree)
SCIP_RETCODE SCIPgetNlRowSolActivity(SCIP *scip, SCIP_NLROW *nlrow, SCIP_SOL *sol, SCIP_Real *activity)
Definition: scip_nlp.c:1911
SCIP_RETCODE SCIPcreateEmptyRowSepa(SCIP *scip, SCIP_ROW **row, SCIP_SEPA *sepa, const char *name, SCIP_Real lhs, SCIP_Real rhs, SCIP_Bool local, SCIP_Bool modifiable, SCIP_Bool removable)
Definition: scip_lp.c:1297
const char * SCIPnlpiGetName(SCIP_NLPI *nlpi)
Definition: nlpi.c:743
#define MAX_ITER
Definition: sepa_gauge.c:67
memory allocation routines