Scippy

SCIP

Solving Constraint Integer Programs

sepa_eccuts.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-2022 Konrad-Zuse-Zentrum */
7 /* fuer Informationstechnik Berlin */
8 /* */
9 /* SCIP is distributed under the terms of the ZIB Academic License. */
10 /* */
11 /* You should have received a copy of the ZIB Academic License */
12 /* along with SCIP; see the file COPYING. If not visit scipopt.org. */
13 /* */
14 /* * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * */
15 
16 /**@file sepa_eccuts.c
17  * @ingroup DEFPLUGINS_SEPA
18  * @brief edge concave cut separator
19  * @author Benjamin Mueller
20  */
21 
22 /**@todo only count number of fixed variables in the edge concave terms */
23 /**@todo only add nonlinear row aggregations where at least ...% of the variables (bilinear terms?) are in edge concave
24  * terms */
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 "scip/scipdefplugins.h"
31 #include "scip/sepa_eccuts.h"
32 #include "scip/cons_xor.h"
33 #include "scip/nlp.h"
34 #include "tclique/tclique.h"
35 
36 #define SEPA_NAME "eccuts"
37 #define SEPA_DESC "separator for edge-concave functions"
38 #define SEPA_PRIORITY -13000
39 #define SEPA_FREQ -1
40 #define SEPA_MAXBOUNDDIST 1.0
41 #define SEPA_USESSUBSCIP FALSE /**< does the separator use a secondary SCIP instance? */
42 #define SEPA_DELAY FALSE /**< should separation method be delayed, if other separators found cuts? */
43 
44 #define CLIQUE_MAXFIRSTNODEWEIGHT 1000 /**< maximum weight of branching nodes in level 0; 0 if not used for cliques
45  * with at least one fractional node) */
46 #define CLIQUE_MINWEIGHT 0 /**< lower bound for weight of generated cliques */
47 #define CLIQUE_MAXNTREENODES 10000 /**< maximal number of nodes of b&b tree */
48 #define CLIQUE_BACKTRACKFREQ 10000 /**< frequency to backtrack to first level of tree (0: no premature backtracking) */
49 
50 #define DEFAULT_DYNAMICCUTS TRUE /**< should generated cuts be removed from the LP if they are no longer tight? */
51 #define DEFAULT_MAXROUNDS 10 /**< maximal number of separation rounds per node (-1: unlimited) */
52 #define DEFAULT_MAXROUNDSROOT 250 /**< maximal number of separation rounds in the root node (-1: unlimited) */
53 #define DEFAULT_MAXDEPTH -1 /**< maximal depth at which the separator is applied */
54 #define DEFAULT_MAXSEPACUTS 10 /**< maximal number of e.c. cuts separated per separation round */
55 #define DEFAULT_MAXSEPACUTSROOT 50 /**< maximal number of e.c. cuts separated per separation round in root node */
56 #define DEFAULT_CUTMAXRANGE 1e+7 /**< maximal coefficient range of a cut (maximal coefficient divided by minimal
57  * coefficient) in order to be added to LP relaxation */
58 #define DEFAULT_MINVIOLATION 0.3 /**< minimal violation of an e.c. cut to be separated */
59 #define DEFAULT_MINAGGRSIZE 3 /**< search for e.c. aggregation of at least this size (has to be >= 3) */
60 #define DEFAULT_MAXAGGRSIZE 4 /**< search for e.c. aggregation of at most this size (has to be >= minaggrsize) */
61 #define DEFAULT_MAXBILINTERMS 500 /**< maximum number of bilinear terms allowed to be in a quadratic constraint */
62 #define DEFAULT_MAXSTALLROUNDS 5 /**< maximum number of unsuccessful rounds in the e.c. aggregation search */
63 
64 #define SUBSCIP_NODELIMIT 100LL /**< node limit to solve the sub-SCIP */
65 
66 #define ADJUSTFACETTOL 1e-6 /**< adjust resulting facets in checkRikun() up to a violation of this value */
67 #define USEDUALSIMPLEX TRUE /**< use dual or primal simplex algorithm? */
68 
69 /** first values for \f$2^n\f$ */
70 static const int poweroftwo[] = { 1, 2, 4, 8, 16, 32, 64, 128, 256, 512, 1024, 2048, 4096, 8192 };
71 
72 /*
73  * Data structures
74  */
75 
76 /** data to store a single edge-concave aggregations; an edge-concave aggregation of a quadratic constraint is a subset
77  * of nonconvex bilinear terms
78  */
79 struct EcAggr
80 {
81  SCIP_VAR** vars; /**< variables */
82  int nvars; /**< number of variables */
83  int varsize; /**< size of vars array */
84 
85  SCIP_Real* termcoefs; /**< coefficients of bilinear terms */
86  int* termvars1; /**< index of the first variable of each bilinear term */
87  int* termvars2; /**< index of the second variable of each bilinear term*/
88  int nterms; /**< number of bilinear terms in the aggregation */
89  int termsize; /**< size of term{coefs,vars1,vars2} arrays */
90 };
91 typedef struct EcAggr SCIP_ECAGGR;
92 
93 /** data to store all edge-concave aggregations and the remaining part of a nonlinear row of the form g(x) <= rhs */
94 struct NlrowAggr
95 {
96  SCIP_NLROW* nlrow; /**< nonlinear row aggregation */
97  SCIP_Bool rhsaggr; /**< consider nonlinear row aggregation for g(x) <= rhs (TRUE) or
98  * g(x) >= lhs (FALSE) */
99 
100  SCIP_ECAGGR** ecaggr; /**< array with all edge-concave aggregations */
101  int necaggr; /**< number of edge-concave aggregation */
102 
103  SCIP_VAR** linvars; /**< linear variables */
104  SCIP_Real* lincoefs; /**< linear coefficients */
105  int nlinvars; /**< number of linear variables */
106  int linvarssize; /**< size of linvars array */
107 
108  SCIP_VAR** quadvars; /**< quadratic variables */
109  int* quadvar2aggr; /**< stores in which edge-concave aggregation the i-th quadratic variable
110  * is contained (< 0: in no edge-concave aggregation) */
111  int nquadvars; /**< number of quadratic variables */
112  int quadvarssize; /**< size of quadvars array */
113 
114  SCIP_VAR** remtermvars1; /**< first quadratic variable of remaining bilinear terms */
115  SCIP_VAR** remtermvars2; /**< second quadratic variable of remaining bilinear terms */
116  SCIP_Real* remtermcoefs; /**< coefficients for each remaining bilinear term */
117  int nremterms; /**< number of remaining bilinear terms */
118  int remtermsize; /**< size of remterm* arrays */
119 
120  SCIP_Real rhs; /**< rhs of the nonlinear row */
121  SCIP_Real constant; /**< constant part of the nonlinear row */
122 };
123 typedef struct NlrowAggr SCIP_NLROWAGGR;
124 
125 /** separator data */
126 struct SCIP_SepaData
127 {
128  SCIP_NLROWAGGR** nlrowaggrs; /**< array containing all nonlinear row aggregations */
129  int nnlrowaggrs; /**< number of nonlinear row aggregations */
130  int nlrowaggrssize; /**< size of nlrowaggrs array */
131  SCIP_Bool searchedforaggr; /**< flag if we already searched for nlrow aggregation candidates */
132  int minaggrsize; /**< only search for e.c. aggregations of at least this size (has to be >= 3) */
133  int maxaggrsize; /**< only search for e.c. aggregations of at most this size (has to be >= minaggrsize) */
134  int maxecsize; /**< largest edge concave aggregation size */
135  int maxbilinterms; /**< maximum number of bilinear terms allowed to be in a quadratic constraint */
136  int maxstallrounds; /**< maximum number of unsuccessful rounds in the e.c. aggregation search */
137 
138  SCIP_LPI* lpi; /**< LP interface to solve the LPs to compute the facets of the convex envelopes */
139  int lpisize; /**< maximum size of e.c. aggregations which can be handled by the LP interface */
140 
141  SCIP_Real cutmaxrange; /**< maximal coef range of a cut (maximal coefficient divided by minimal
142  * coefficient) in order to be added to LP relaxation */
143  SCIP_Bool dynamiccuts; /**< should generated cuts be removed from the LP if they are no longer tight? */
144  SCIP_Real minviolation; /**< minimal violation of an e.c. cut to be separated */
145 
146  int maxrounds; /**< maximal number of separation rounds per node (-1: unlimited) */
147  int maxroundsroot; /**< maximal number of separation rounds in the root node (-1: unlimited) */
148  int maxdepth; /**< maximal depth at which the separator is applied */
149  int maxsepacuts; /**< maximal number of e.c. cuts separated per separation round */
150  int maxsepacutsroot; /**< maximal number of e.c. cuts separated per separation round in root node */
151 
152 #ifdef SCIP_STATISTIC
153  SCIP_Real aggrsearchtime; /**< total time spent for searching edge concave aggregations */
154  int nlhsnlrowaggrs; /**< number of found nonlinear row aggregations for SCIP_NLROWs of the form g(x) <= rhs */
155  int nrhsnlrowaggrs; /**< number of found nonlinear row aggregations for SCIP_NLROWs of the form g(x) >= lhs */
156 #endif
157 };
158 
159 
160 /*
161  * Local methods
162  */
163 
164 /** creates an empty edge-concave aggregation (without bilinear terms) */
165 static
166 SCIP_RETCODE ecaggrCreateEmpty(
167  SCIP* scip, /**< SCIP data structure */
168  SCIP_ECAGGR** ecaggr, /**< pointer to store the edge-concave aggregation */
169  int nquadvars, /**< number of quadratic variables */
170  int nquadterms /**< number of bilinear terms */
171  )
172 {
173  assert(scip != NULL);
174  assert(ecaggr != NULL);
175  assert(nquadvars > 0);
176  assert(nquadterms >= nquadvars);
177 
178  SCIP_CALL( SCIPallocBlockMemory(scip, ecaggr) );
179 
180  (*ecaggr)->nvars = 0;
181  (*ecaggr)->nterms = 0;
182  (*ecaggr)->varsize = nquadvars;
183  (*ecaggr)->termsize = nquadterms;
184 
185  /* allocate enough memory for the quadratic variables and bilinear terms */
186  SCIP_CALL( SCIPallocBlockMemoryArray(scip, &(*ecaggr)->vars, nquadvars) );
187  SCIP_CALL( SCIPallocBlockMemoryArray(scip, &(*ecaggr)->termcoefs, nquadterms) );
188  SCIP_CALL( SCIPallocBlockMemoryArray(scip, &(*ecaggr)->termvars1, nquadterms) );
189  SCIP_CALL( SCIPallocBlockMemoryArray(scip, &(*ecaggr)->termvars2, nquadterms) );
190 
191  return SCIP_OKAY;
192 }
193 
194 /** frees an edge-concave aggregation */
195 static
196 SCIP_RETCODE ecaggrFree(
197  SCIP* scip, /**< SCIP data structure */
198  SCIP_ECAGGR** ecaggr /**< pointer to store the edge-concave aggregation */
199  )
200 {
201  assert(scip != NULL);
202  assert(ecaggr != NULL);
203 
204  SCIPfreeBlockMemoryArray(scip, &((*ecaggr)->termcoefs), (*ecaggr)->termsize);
205  SCIPfreeBlockMemoryArray(scip, &((*ecaggr)->termvars1), (*ecaggr)->termsize);
206  SCIPfreeBlockMemoryArray(scip, &((*ecaggr)->termvars2), (*ecaggr)->termsize);
207  SCIPfreeBlockMemoryArray(scip, &((*ecaggr)->vars), (*ecaggr)->varsize);
208 
209  SCIPfreeBlockMemory(scip, ecaggr);
210  *ecaggr = NULL;
211 
212  return SCIP_OKAY;
213 }
214 
215 /** adds a quadratic variable to an edge-concave aggregation */
216 static
217 SCIP_RETCODE ecaggrAddQuadvar(
218  SCIP_ECAGGR* ecaggr, /**< pointer to store the edge-concave aggregation */
219  SCIP_VAR* x /**< first variable */
220  )
221 {
222  ecaggr->vars[ ecaggr->nvars++ ] = x;
223  return SCIP_OKAY;
224 }
225 
226 /** adds a bilinear term to an edge-concave aggregation */
227 static
228 SCIP_RETCODE ecaggrAddBilinTerm(
229  SCIP* scip, /**< SCIP data structure */
230  SCIP_ECAGGR* ecaggr, /**< pointer to store the edge-concave aggregation */
231  SCIP_VAR* x, /**< first variable */
232  SCIP_VAR* y, /**< second variable */
233  SCIP_Real coef /**< bilinear coefficient */
234  )
235 {
236  int idx1;
237  int idx2;
238  int i;
239 
240  assert(x != NULL);
241  assert(y != NULL);
242  assert(ecaggr->nterms + 1 <= ((ecaggr->nvars + 1) * ecaggr->nvars) / 2);
243  assert(!SCIPisZero(scip, coef));
244 
245  idx1 = -1;
246  idx2 = -1;
247 
248  /* search for the quadratic variables in the e.c. aggregation */
249  for( i = 0; i < ecaggr->nvars && (idx1 == -1 || idx2 == -1); ++i )
250  {
251  if( ecaggr->vars[i] == x )
252  idx1 = i;
253  if( ecaggr->vars[i] == y )
254  idx2 = i;
255  }
256 
257  assert(idx1 != -1 && idx2 != -1);
258 
259  ecaggr->termcoefs[ ecaggr->nterms ] = coef;
260  ecaggr->termvars1[ ecaggr->nterms ] = idx1;
261  ecaggr->termvars2[ ecaggr->nterms ] = idx2;
262  ++(ecaggr->nterms);
263 
264  return SCIP_OKAY;
265 }
266 
267 #ifdef SCIP_DEBUG
268 /** prints an edge-concave aggregation */
269 static
270 void ecaggrPrint(
271  SCIP* scip, /**< SCIP data structure */
272  SCIP_ECAGGR* ecaggr /**< pointer to store the edge-concave aggregation */
273  )
274 {
275  int i;
276 
277  assert(scip != NULL);
278  assert(ecaggr != NULL);
279 
280  SCIPdebugMsg(scip, " nvars = %d nterms = %d\n", ecaggr->nvars, ecaggr->nterms);
281  SCIPdebugMsg(scip, " vars: ");
282  for( i = 0; i < ecaggr->nvars; ++i )
283  SCIPdebugMsgPrint(scip, "%s ", SCIPvarGetName(ecaggr->vars[i]));
284  SCIPdebugMsgPrint(scip, "\n");
285 
286  SCIPdebugMsg(scip, " terms: ");
287  for( i = 0; i < ecaggr->nterms; ++i )
288  {
289  SCIP_VAR* x;
290  SCIP_VAR* y;
291 
292  x = ecaggr->vars[ ecaggr->termvars1[i] ];
293  y = ecaggr->vars[ ecaggr->termvars2[i] ];
294  SCIPdebugMsgPrint(scip, "%e %s * %s ", ecaggr->termcoefs[i], SCIPvarGetName(x), SCIPvarGetName(y) );
295  }
296  SCIPdebugMsgPrint(scip, "\n");
297 }
298 #endif
299 
300 /** stores linear terms in a given nonlinear row aggregation */
301 static
302 SCIP_RETCODE nlrowaggrStoreLinearTerms(
303  SCIP* scip, /**< SCIP data structure */
304  SCIP_NLROWAGGR* nlrowaggr, /**< nonlinear row aggregation */
305  SCIP_VAR** linvars, /**< linear variables */
306  SCIP_Real* lincoefs, /**< linear coefficients */
307  int nlinvars /**< number of linear variables */
308  )
309 {
310  assert(scip != NULL);
311  assert(nlrowaggr != NULL);
312  assert(linvars != NULL || nlinvars == 0);
313  assert(lincoefs != NULL || nlinvars == 0);
314  assert(nlinvars >= 0);
315 
316  nlrowaggr->nlinvars = 0;
317  nlrowaggr->linvarssize = 0;
318  nlrowaggr->linvars = NULL;
319  nlrowaggr->lincoefs = NULL;
320 
321  if( nlinvars == 0 )
322  return SCIP_OKAY;
323 
324  SCIP_CALL( SCIPduplicateBlockMemoryArray(scip, &nlrowaggr->linvars, linvars, nlinvars) );
325  SCIP_CALL( SCIPduplicateBlockMemoryArray(scip, &nlrowaggr->lincoefs, lincoefs, nlinvars) );
326  nlrowaggr->nlinvars = nlinvars;
327  nlrowaggr->linvarssize = nlinvars;
328 
329  /* if we have a nlrow of the form g(x) >= lhs, multiply every coefficient by -1 */
330  if( !nlrowaggr->rhsaggr )
331  {
332  int i;
333 
334  for( i = 0; i < nlrowaggr->nlinvars; ++i )
335  nlrowaggr->lincoefs[i] *= -1.0;
336  }
337 
338  return SCIP_OKAY;
339 }
340 
341 /** adds linear term to a given nonlinear row aggregation */
342 static
343 SCIP_RETCODE nlrowaggrAddLinearTerm(
344  SCIP* scip, /**< SCIP data structure */
345  SCIP_NLROWAGGR* nlrowaggr, /**< nonlinear row aggregation */
346  SCIP_VAR* linvar, /**< linear variable */
347  SCIP_Real lincoef /**< coefficient */
348  )
349 {
350  assert(scip != NULL);
351  assert(nlrowaggr != NULL);
352  assert(linvar != NULL);
353 
354  if( nlrowaggr->nlinvars == nlrowaggr->linvarssize )
355  {
356  int newsize = SCIPcalcMemGrowSize(scip, nlrowaggr->linvarssize+1);
357  SCIP_CALL( SCIPreallocBlockMemoryArray(scip, &nlrowaggr->linvars, nlrowaggr->linvarssize, newsize) );
358  SCIP_CALL( SCIPreallocBlockMemoryArray(scip, &nlrowaggr->lincoefs, nlrowaggr->linvarssize, newsize) );
359  nlrowaggr->linvarssize = newsize;
360  }
361  assert(nlrowaggr->linvarssize > nlrowaggr->nlinvars);
362 
363  /* if we have a nlrow of the form g(x) >= lhs, multiply coefficient by -1 */
364  if( !nlrowaggr->rhsaggr )
365  lincoef = -lincoef;
366 
367  nlrowaggr->linvars[nlrowaggr->nlinvars] = linvar;
368  nlrowaggr->lincoefs[nlrowaggr->nlinvars] = lincoef;
369  ++nlrowaggr->nlinvars;
370 
371  return SCIP_OKAY;
372 }
373 
374 /** adds quadratic variable to a given nonlinear row aggregation */
375 static
376 SCIP_RETCODE nlrowaggrAddQuadraticVar(
377  SCIP* scip, /**< SCIP data structure */
378  SCIP_NLROWAGGR* nlrowaggr, /**< nonlinear row aggregation */
379  SCIP_VAR* quadvar /**< quadratic variable */
380  )
381 {
382  assert(scip != NULL);
383  assert(nlrowaggr != NULL);
384  assert(quadvar != NULL);
385 
386  SCIP_CALL( SCIPensureBlockMemoryArray(scip, &nlrowaggr->quadvars, &nlrowaggr->quadvarssize, nlrowaggr->nquadvars+1) );
387  assert(nlrowaggr->quadvarssize > nlrowaggr->nquadvars);
388  nlrowaggr->quadvars[nlrowaggr->nquadvars] = quadvar;
389  ++nlrowaggr->nquadvars;
390 
391  return SCIP_OKAY;
392 }
393 
394 /** adds a remaining bilinear term to a given nonlinear row aggregation */
395 static
396 SCIP_RETCODE nlrowaggrAddRemBilinTerm(
397  SCIP_NLROWAGGR* nlrowaggr, /**< nonlinear row aggregation */
398  SCIP_VAR* x, /**< first variable */
399  SCIP_VAR* y, /**< second variable */
400  SCIP_Real coef /**< bilinear coefficient */
401  )
402 {
403  assert(nlrowaggr != NULL);
404  assert(x != NULL);
405  assert(y != NULL);
406  assert(coef != 0.0);
407  assert(nlrowaggr->remtermcoefs != NULL);
408  assert(nlrowaggr->remtermvars1 != NULL);
409  assert(nlrowaggr->remtermvars2 != NULL);
410 
411  nlrowaggr->remtermcoefs[ nlrowaggr->nremterms ] = coef;
412  nlrowaggr->remtermvars1[ nlrowaggr->nremterms ] = x;
413  nlrowaggr->remtermvars2[ nlrowaggr->nremterms ] = y;
414  ++(nlrowaggr->nremterms);
415 
416  return SCIP_OKAY;
417 }
418 
419 /** creates a nonlinear row aggregation */
420 static
421 SCIP_RETCODE nlrowaggrCreate(
422  SCIP* scip, /**< SCIP data structure */
423  SCIP_NLROW* nlrow, /**< nonlinear row */
424  SCIP_NLROWAGGR** nlrowaggr, /**< pointer to store the nonlinear row aggregation */
425  int* quadvar2aggr, /**< mapping between quadratic variables and edge-concave aggregation
426  * stores a negative value if the quadratic variables does not belong
427  * to any aggregation */
428  int nfound, /**< number of edge-concave aggregations */
429  SCIP_Bool rhsaggr /**< consider nonlinear row aggregation for g(x) <= rhs (TRUE) or
430  * lhs <= g(x) (FALSE) */
431  )
432 {
433  SCIP_EXPR* expr;
434  int* aggrnvars; /* count the number of variables in each e.c. aggregations */
435  int* aggrnterms; /* count the number of bilinear terms in each e.c. aggregations */
436  int nquadvars;
437  int nremterms;
438  int i;
439 
440  assert(scip != NULL);
441  assert(nlrow != NULL);
442  assert(nlrowaggr != NULL);
443  assert(quadvar2aggr != NULL);
444  assert(nfound > 0);
445 
446  expr = SCIPnlrowGetExpr(nlrow);
447  SCIPexprGetQuadraticData(expr, NULL, NULL, NULL, NULL, &nquadvars, NULL, NULL, NULL);
448  nremterms = 0;
449 
450  SCIP_CALL( SCIPallocClearBufferArray(scip, &aggrnvars, nfound) );
451  SCIP_CALL( SCIPallocClearBufferArray(scip, &aggrnterms, nfound) );
452 
453  /* create an empty nonlinear row aggregation */
454  SCIP_CALL( SCIPallocBlockMemory(scip, nlrowaggr) );
455  (*nlrowaggr)->nlrow = nlrow;
456  (*nlrowaggr)->rhsaggr = rhsaggr;
457  (*nlrowaggr)->rhs = rhsaggr ? SCIPnlrowGetRhs(nlrow) : -SCIPnlrowGetLhs(nlrow);
458  (*nlrowaggr)->constant = rhsaggr ? SCIPnlrowGetConstant(nlrow) : -SCIPnlrowGetConstant(nlrow);
459 
460  (*nlrowaggr)->quadvars = NULL;
461  (*nlrowaggr)->nquadvars = 0;
462  (*nlrowaggr)->quadvarssize = 0;
463  (*nlrowaggr)->quadvar2aggr = NULL;
464  (*nlrowaggr)->remtermcoefs = NULL;
465  (*nlrowaggr)->remtermvars1 = NULL;
466  (*nlrowaggr)->remtermvars2 = NULL;
467  (*nlrowaggr)->nremterms = 0;
468 
469  /* copy quadvar2aggr array */
470  SCIP_CALL( SCIPduplicateBlockMemoryArray(scip, &(*nlrowaggr)->quadvar2aggr, quadvar2aggr, nquadvars) );
471 
472  /* store all linear terms */
473  SCIP_CALL( nlrowaggrStoreLinearTerms(scip, *nlrowaggr, SCIPnlrowGetLinearVars(nlrow), SCIPnlrowGetLinearCoefs(nlrow),
474  SCIPnlrowGetNLinearVars(nlrow)) );
475 
476  /* store all quadratic variables and additional linear terms */
477  /* count the number of variables in each e.c. aggregation */
478  /* count the number of square and bilinear terms in each e.c. aggregation */
479  for( i = 0; i < nquadvars; ++i )
480  {
481  SCIP_EXPR* qterm;
482  SCIP_Real lincoef;
483  SCIP_Real sqrcoef;
484  int idx1;
485  int nadjbilin;
486  int* adjbilin;
487  int j;
488 
489  SCIPexprGetQuadraticQuadTerm(expr, i, &qterm, &lincoef, &sqrcoef, &nadjbilin, &adjbilin, NULL);
490  assert(SCIPisExprVar(scip, qterm));
491 
492  SCIP_CALL( nlrowaggrAddQuadraticVar(scip, *nlrowaggr, SCIPgetVarExprVar(qterm)) );
493 
494  if( lincoef != 0.0 )
495  {
496  SCIP_CALL( nlrowaggrAddLinearTerm(scip, *nlrowaggr, SCIPgetVarExprVar(qterm), lincoef) );
497  }
498 
499  if( quadvar2aggr[i] >= 0)
500  ++aggrnvars[ quadvar2aggr[i] ];
501 
502  idx1 = quadvar2aggr[i];
503  if( rhsaggr )
504  sqrcoef = -sqrcoef;
505 
506  /* variable has to belong to an e.c. aggregation; square term has to be concave */
507  if( idx1 >= 0 && SCIPisNegative(scip, sqrcoef) )
508  ++aggrnterms[idx1];
509  else
510  ++nremterms;
511 
512  for( j = 0; j < nadjbilin; ++j )
513  {
514  SCIP_EXPR* qterm1;
515  int pos2;
516  int idx2;
517 
518  SCIPexprGetQuadraticBilinTerm(expr, adjbilin[j], &qterm1, NULL, NULL, &pos2, NULL);
519 
520  /* only handle qterm1 == qterm here; the other case will be handled when its turn for qterm2 to be qterm */
521  if( qterm1 != qterm )
522  continue;
523 
524  idx2 = quadvar2aggr[pos2];
525 
526  /* variables have to belong to the same e.c. aggregation; bilinear term has to be concave */
527  if( idx1 >= 0 && idx2 >= 0 && idx1 == idx2 )
528  ++aggrnterms[idx1];
529  else
530  ++nremterms;
531  }
532  }
533  assert((*nlrowaggr)->nquadvars == nquadvars);
534 
535  /* create all edge-concave aggregations (empty) and remaining terms */
536  SCIP_CALL( SCIPallocBlockMemoryArray(scip, &(*nlrowaggr)->ecaggr, nfound) );
537  if( nremterms > 0 )
538  {
539  SCIP_CALL( SCIPallocBlockMemoryArray(scip, &(*nlrowaggr)->remtermcoefs, nremterms) );
540  SCIP_CALL( SCIPallocBlockMemoryArray(scip, &(*nlrowaggr)->remtermvars1, nremterms) );
541  SCIP_CALL( SCIPallocBlockMemoryArray(scip, &(*nlrowaggr)->remtermvars2, nremterms) );
542  (*nlrowaggr)->remtermsize = nremterms;
543  }
544  (*nlrowaggr)->necaggr = nfound;
545 
546  for( i = 0; i < nfound; ++i )
547  {
548  SCIP_CALL( ecaggrCreateEmpty(scip, &(*nlrowaggr)->ecaggr[i], aggrnvars[i], aggrnterms[i]) );
549  }
550 
551  /* add quadratic variables to the edge-concave aggregations */
552  for( i = 0; i < nquadvars; ++i )
553  {
554  int idx;
555 
556  idx = quadvar2aggr[i];
557 
558  if( idx >= 0)
559  {
560  SCIP_EXPR* qterm;
561 
562  SCIPdebugMsg(scip, "add quadvar %d to aggr. %d\n", i, idx);
563 
564  SCIPexprGetQuadraticQuadTerm(expr, i, &qterm, NULL, NULL, NULL, NULL, NULL);
565  assert(SCIPisExprVar(scip, qterm));
566 
567  SCIP_CALL( ecaggrAddQuadvar((*nlrowaggr)->ecaggr[idx], SCIPgetVarExprVar(qterm)) );
568  }
569  }
570 
571  /* add the bilinear/square terms to the edge-concave aggregations or in the remaining part */
572  for( i = 0; i < nquadvars; ++i )
573  {
574  SCIP_EXPR* qterm;
575  SCIP_VAR* x;
576  SCIP_Real coef;
577  int idx1;
578  int nadjbilin;
579  int* adjbilin;
580  int j;
581 
582  SCIPexprGetQuadraticQuadTerm(expr, i, &qterm, NULL, &coef, &nadjbilin, &adjbilin, NULL);
583 
584  x = SCIPgetVarExprVar(qterm);
585 
586  idx1 = quadvar2aggr[i];
587  if( rhsaggr )
588  coef = -coef;
589 
590  if( idx1 >= 0 && SCIPisNegative(scip, coef) )
591  {
592  SCIP_CALL( ecaggrAddBilinTerm(scip, (*nlrowaggr)->ecaggr[idx1], x, x, coef) );
593  SCIPdebugMsg(scip, "add term %e *%d^2 to aggr. %d\n", coef, i, idx1);
594  }
595  else
596  {
597  SCIP_CALL( nlrowaggrAddRemBilinTerm(*nlrowaggr, x, x, coef) );
598  SCIPdebugMsg(scip, "add term %e *%d^2 to the remaining part\n", coef, idx1);
599  }
600 
601  for( j = 0; j < nadjbilin; ++j )
602  {
603  SCIP_EXPR* qterm1;
604  SCIP_EXPR* qterm2;
605  int pos2;
606  int idx2;
607  SCIP_VAR* y;
608 
609  SCIPexprGetQuadraticBilinTerm(expr, adjbilin[j], &qterm1, &qterm2, &coef, &pos2, NULL);
610 
611  /* only handle qterm1 == qterm here; the other case will be handled when its turn for qterm2 to be qterm */
612  if( qterm1 != qterm )
613  continue;
614 
615  y = SCIPgetVarExprVar(qterm2);
616 
617  idx2 = quadvar2aggr[pos2];
618  if( rhsaggr )
619  coef = -coef;
620 
621  if( idx1 >= 0 && idx2 >= 0 && idx1 == idx2 )
622  {
623  SCIP_CALL( ecaggrAddBilinTerm(scip, (*nlrowaggr)->ecaggr[idx1], x, y, coef) );
624  SCIPdebugMsg(scip, "add term %e *%d*%d to aggr. %d\n", coef, i, pos2, idx1);
625  }
626  else
627  {
628  SCIP_CALL( nlrowaggrAddRemBilinTerm(*nlrowaggr, x, y, coef) );
629  SCIPdebugMsg(scip, "add term %e *%d*%d to the remaining part\n", coef, i, pos2);
630  }
631  }
632  }
633 
634  /* free allocated memory */
635  SCIPfreeBufferArray(scip, &aggrnterms);
636  SCIPfreeBufferArray(scip, &aggrnvars);
637 
638  return SCIP_OKAY;
639 }
640 
641 /** frees a nonlinear row aggregation */
642 static
643 SCIP_RETCODE nlrowaggrFree(
644  SCIP* scip, /**< SCIP data structure */
645  SCIP_NLROWAGGR** nlrowaggr /**< pointer to free the nonlinear row aggregation */
646  )
647 {
648  int i;
649 
650  assert(scip != NULL);
651  assert(nlrowaggr != NULL);
652  assert(*nlrowaggr != NULL);
653  (*nlrowaggr)->nlrow = NULL;
654  assert((*nlrowaggr)->quadvars != NULL);
655  assert((*nlrowaggr)->nquadvars > 0);
656  assert((*nlrowaggr)->nremterms >= 0);
657 
658  /* free remaining part */
659  SCIPfreeBlockMemoryArrayNull(scip, &(*nlrowaggr)->remtermcoefs, (*nlrowaggr)->remtermsize);
660  SCIPfreeBlockMemoryArrayNull(scip, &(*nlrowaggr)->remtermvars1, (*nlrowaggr)->remtermsize);
661  SCIPfreeBlockMemoryArrayNull(scip, &(*nlrowaggr)->remtermvars2, (*nlrowaggr)->remtermsize);
662 
663  /* free quadratic variables */
664  SCIPfreeBlockMemoryArray(scip, &(*nlrowaggr)->quadvars, (*nlrowaggr)->quadvarssize);
665  SCIPfreeBlockMemoryArray(scip, &(*nlrowaggr)->quadvar2aggr, (*nlrowaggr)->nquadvars);
666 
667  /* free linear part */
668  if( (*nlrowaggr)->nlinvars > 0 )
669  {
670  SCIPfreeBlockMemoryArray(scip, &(*nlrowaggr)->linvars, (*nlrowaggr)->linvarssize);
671  SCIPfreeBlockMemoryArray(scip, &(*nlrowaggr)->lincoefs, (*nlrowaggr)->linvarssize);
672  }
673 
674  /* free edge-concave aggregations */
675  for( i = 0; i < (*nlrowaggr)->necaggr; ++i )
676  {
677  SCIP_CALL( ecaggrFree(scip, &(*nlrowaggr)->ecaggr[i]) );
678  }
679  SCIPfreeBlockMemoryArray(scip, &(*nlrowaggr)->ecaggr, (*nlrowaggr)->necaggr);
680 
681  /* free nlrow aggregation */
682  SCIPfreeBlockMemory(scip, nlrowaggr);
683 
684  return SCIP_OKAY;
685 }
686 
687 #ifdef SCIP_DEBUG
688 /** prints a nonlinear row aggregation */
689 static
690 void nlrowaggrPrint(
691  SCIP* scip, /**< SCIP data structure */
692  SCIP_NLROWAGGR* nlrowaggr /**< nonlinear row aggregation */
693  )
694 {
695  int i;
696 
697  SCIPdebugMsg(scip, " nlrowaggr rhs = %e\n", nlrowaggr->rhs);
698  SCIPdebugMsg(scip, " #remaining terms = %d\n", nlrowaggr->nremterms);
699 
700  SCIPdebugMsg(scip, "remaining terms: ");
701  for( i = 0; i < nlrowaggr->nremterms; ++i )
702  SCIPdebugMsgPrint(scip, "%e %s * %s + ", nlrowaggr->remtermcoefs[i], SCIPvarGetName(nlrowaggr->remtermvars1[i]),
703  SCIPvarGetName(nlrowaggr->remtermvars2[i]) );
704  for( i = 0; i < nlrowaggr->nlinvars; ++i )
705  SCIPdebugMsgPrint(scip, "%e %s + ", nlrowaggr->lincoefs[i], SCIPvarGetName(nlrowaggr->linvars[i]) );
706  SCIPdebugMsgPrint(scip, "\n");
707 
708  for( i = 0; i < nlrowaggr->necaggr; ++i )
709  {
710  SCIPdebugMsg(scip, "print e.c. aggr %d\n", i);
711  ecaggrPrint(scip, nlrowaggr->ecaggr[i]);
712  }
713  return;
714 }
715 #endif
716 
717 /** creates separator data */
718 static
719 SCIP_RETCODE sepadataCreate(
720  SCIP* scip, /**< SCIP data structure */
721  SCIP_SEPADATA** sepadata /**< pointer to store separator data */
722  )
723 {
724  assert(scip != NULL);
725  assert(sepadata != NULL);
726 
727  SCIP_CALL( SCIPallocBlockMemory(scip, sepadata) );
728  BMSclearMemory(*sepadata);
729 
730  return SCIP_OKAY;
731 }
732 
733 /** frees all nonlinear row aggregations */
734 static
735 SCIP_RETCODE sepadataFreeNlrows(
736  SCIP* scip, /**< SCIP data structure */
737  SCIP_SEPADATA* sepadata /**< pointer to store separator data */
738  )
739 {
740  assert(scip != NULL);
741  assert(sepadata != NULL);
742 
743  /* free nonlinear row aggregations */
744  if( sepadata->nlrowaggrs != NULL )
745  {
746  int i;
747 
748  for( i = sepadata->nnlrowaggrs - 1; i >= 0; --i )
749  {
750  SCIP_CALL( nlrowaggrFree(scip, &sepadata->nlrowaggrs[i]) );
751  }
752 
753  SCIPfreeBlockMemoryArray(scip, &sepadata->nlrowaggrs, sepadata->nlrowaggrssize);
754 
755  sepadata->nlrowaggrs = NULL;
756  sepadata->nnlrowaggrs = 0;
757  sepadata->nlrowaggrssize = 0;
758  }
759 
760  return SCIP_OKAY;
761 }
762 
763 /** frees separator data */
764 static
765 SCIP_RETCODE sepadataFree(
766  SCIP* scip, /**< SCIP data structure */
767  SCIP_SEPADATA** sepadata /**< pointer to store separator data */
768  )
769 {
770  assert(scip != NULL);
771  assert(sepadata != NULL);
772  assert(*sepadata != NULL);
773 
774  /* free nonlinear row aggregations */
775  SCIP_CALL( sepadataFreeNlrows(scip, *sepadata) );
776 
777  /* free LP interface */
778  if( (*sepadata)->lpi != NULL )
779  {
780  SCIP_CALL( SCIPlpiFree(&((*sepadata)->lpi)) );
781  (*sepadata)->lpisize = 0;
782  }
783 
784  SCIPfreeBlockMemory(scip, sepadata);
785 
786  return SCIP_OKAY;
787 }
788 
789 /** adds a nonlinear row aggregation to the separator data */
790 static
791 SCIP_RETCODE sepadataAddNlrowaggr(
792  SCIP* scip, /**< SCIP data structure */
793  SCIP_SEPADATA* sepadata, /**< separator data */
794  SCIP_NLROWAGGR* nlrowaggr /**< non-linear row aggregation */
795  )
796 {
797  int i;
798 
799  assert(scip != NULL);
800  assert(sepadata != NULL);
801  assert(nlrowaggr != NULL);
802 
803  if( sepadata->nlrowaggrssize == 0 )
804  {
805  SCIP_CALL( SCIPallocBlockMemoryArray(scip, &sepadata->nlrowaggrs, 2) ); /*lint !e506*/
806  sepadata->nlrowaggrssize = 2;
807  }
808  else if( sepadata->nlrowaggrssize < sepadata->nnlrowaggrs + 1 )
809  {
810  SCIP_CALL( SCIPreallocBlockMemoryArray(scip, &sepadata->nlrowaggrs, sepadata->nlrowaggrssize, 2 * sepadata->nlrowaggrssize) ); /*lint !e506 !e647*/
811  sepadata->nlrowaggrssize *= 2;
812  assert(sepadata->nlrowaggrssize >= sepadata->nnlrowaggrs + 1);
813  }
814 
815  sepadata->nlrowaggrs[ sepadata->nnlrowaggrs ] = nlrowaggr;
816  ++(sepadata->nnlrowaggrs);
817 
818  /* update maximum e.c. aggregation size */
819  for( i = 0; i < nlrowaggr->necaggr; ++i )
820  sepadata->maxecsize = MAX(sepadata->maxecsize, nlrowaggr->ecaggr[i]->nvars);
821 
822 #ifdef SCIP_STATISTIC
823  /* update statistics */
824  if( nlrowaggr->rhsaggr )
825  ++(sepadata->nrhsnlrowaggrs);
826  else
827  ++(sepadata->nlhsnlrowaggrs);
828 #endif
829 
830  return SCIP_OKAY;
831 }
832 
833 /** returns min{val-lb,ub-val} / (ub-lb) */
834 static
835 SCIP_Real phi(
836  SCIP* scip, /**< SCIP data structure */
837  SCIP_Real val, /**< solution value */
838  SCIP_Real lb, /**< lower bound */
839  SCIP_Real ub /**< upper bound */
840  )
841 {
842  if( SCIPisFeasEQ(scip, lb, ub) )
843  return 0.0;
844 
845  /* adjust */
846  val = MAX(val, lb);
847  val = MIN(val, ub);
848 
849  return MIN(ub - val, val - lb) / (ub - lb);
850 }
851 
852 /** creates an MIP to search for cycles with an odd number of positive edges in the graph representation of a nonlinear row
853  *
854  * The model uses directed binary arc flow variables.
855  * We introduce for all quadratic elements a forward and backward edge.
856  * If the term is quadratic (e.g., loop in the graph) we fix the corresponding variables to zero.
857  * This leads to an easy mapping between quadratic elements and the variables of the MIP.
858  */
859 static
860 SCIP_RETCODE createMIP(
861  SCIP* scip, /**< SCIP data structure */
862  SCIP* subscip, /**< auxiliary SCIP to search aggregations */
863  SCIP_SEPADATA* sepadata, /**< separator data */
864  SCIP_NLROW* nlrow, /**< nonlinear row */
865  SCIP_Bool rhsaggr, /**< consider nonlinear row aggregation for g(x) <= rhs (TRUE) or
866  * lhs <= g(x) (FALSE) */
867  SCIP_VAR** forwardarcs, /**< array to store all forward arc variables */
868  SCIP_VAR** backwardarcs, /**< array to store all backward arc variables */
869  SCIP_Real* nodeweights, /**< weights for each node of the graph */
870  int* nedges, /**< pointer to store the number of nonexcluded edges in the graph */
871  int* narcs /**< pointer to store the number of created arc variables (number of square and bilinear terms) */
872  )
873 {
874  SCIP_VAR** oddcyclearcs;
875  SCIP_CONS** flowcons;
876  SCIP_CONS* cyclelengthcons;
877  SCIP_CONS* oddcyclecons;
878  char name[SCIP_MAXSTRLEN];
879  SCIP_EXPR* expr;
880  int noddcyclearcs;
881  int nnodes;
882  int nquadexprs;
883  int nbilinexprs;
884  int i;
885  int arcidx;
886 
887  assert(subscip != NULL);
888  assert(forwardarcs != NULL);
889  assert(backwardarcs != NULL);
890  assert(nedges != NULL);
891  assert(sepadata->minaggrsize <= sepadata->maxaggrsize);
892 
893  expr = SCIPnlrowGetExpr(nlrow);
894  SCIPexprGetQuadraticData(expr, NULL, NULL, NULL, NULL, &nquadexprs, &nbilinexprs, NULL, NULL);
895 
896  nnodes = nquadexprs;
897  *nedges = 0;
898  *narcs = 0;
899 
900  assert(nnodes > 0);
901 
902  noddcyclearcs = 0;
903  SCIP_CALL( SCIPallocBufferArray(subscip, &oddcyclearcs, 2*nbilinexprs) );
904 
905  /* create problem with default plug-ins */
906  SCIP_CALL( SCIPcreateProbBasic(subscip, "E.C. aggregation MIP") );
907  SCIP_CALL( SCIPsetObjsense(subscip, SCIP_OBJSENSE_MAXIMIZE) );
908  SCIP_CALL( SCIPincludeDefaultPlugins(subscip) );
909 
910  /* create forward and backward arc variables */
911  for( i = 0; i < nquadexprs; ++i )
912  {
913  SCIP_EXPR* qterm;
914  SCIP_Real coef;
915  int nadjbilin;
916  int* adjbilin;
917  int j;
918 
919  SCIPexprGetQuadraticQuadTerm(expr, i, &qterm, NULL, &coef, &nadjbilin, &adjbilin, NULL);
920 
921  if( !SCIPisZero(scip, coef) )
922  {
923  /* squares (loops) are fixed to zero */
924  SCIPdebugMsg(scip, "edge {%d,%d} = {%s,%s} coeff=%e edgeweight=0\n", i, i,
925  SCIPvarGetName(SCIPgetVarExprVar(qterm)), SCIPvarGetName(SCIPgetVarExprVar(qterm)),
926  coef);
927 
928  (void) SCIPsnprintf(name, SCIP_MAXSTRLEN, "x#%d#%d", i, i);
929  SCIP_CALL( SCIPcreateVarBasic(subscip, &forwardarcs[*narcs], name, 0.0, 0.0, 0.01, SCIP_VARTYPE_BINARY) );
930  SCIP_CALL( SCIPaddVar(subscip, forwardarcs[*narcs]) );
931 
932  SCIP_CALL( SCIPcreateVarBasic(subscip, &backwardarcs[*narcs], name, 0.0, 0.0, 0.01, SCIP_VARTYPE_BINARY) );
933  SCIP_CALL( SCIPaddVar(subscip, backwardarcs[*narcs]) );
934 
935  ++*narcs;
936  }
937 
938  for( j = 0 ; j < nadjbilin; ++j )
939  {
940  SCIP_EXPR* qterm1;
941  SCIP_EXPR* qterm2;
942  int pos2;
943  SCIP_Real edgeweight;
944  SCIP_CONS* noparallelcons;
945 
946  SCIPexprGetQuadraticBilinTerm(expr, adjbilin[j], &qterm1, &qterm2, &coef, &pos2, NULL);
947 
948  /* handle qterm == qterm2 later */
949  if( qterm1 != qterm )
950  continue;
951 
952  edgeweight = nodeweights[i] + nodeweights[pos2];
953  SCIPdebugMsg(scip, "edge {%d,%d} = {%s,%s} coeff=%e edgeweight=%e\n", i, pos2,
954  SCIPvarGetName(SCIPgetVarExprVar(qterm1)), SCIPvarGetName(SCIPgetVarExprVar(qterm2)),
955  coef, edgeweight);
956 
957  (void) SCIPsnprintf(name, SCIP_MAXSTRLEN, "x#%d#%d", i, pos2);
958  SCIP_CALL( SCIPcreateVarBasic(subscip, &forwardarcs[*narcs], name, 0.0, 1.0, 0.01 + edgeweight, SCIP_VARTYPE_BINARY) );
959  SCIP_CALL( SCIPaddVar(subscip, forwardarcs[*narcs]) );
960 
961  (void) SCIPsnprintf(name, SCIP_MAXSTRLEN, "x#%d#%d", i, pos2);
962  SCIP_CALL( SCIPcreateVarBasic(subscip, &backwardarcs[*narcs], name, 0.0, 1.0, 0.01 + edgeweight, SCIP_VARTYPE_BINARY) );
963  SCIP_CALL( SCIPaddVar(subscip, backwardarcs[*narcs]) );
964 
965  ++(*nedges);
966 
967  /* store all arcs which are important for the odd cycle property (no loops) */
968  if( rhsaggr && SCIPisPositive(scip, coef) )
969  {
970  assert(noddcyclearcs < 2*nbilinexprs-1);
971  oddcyclearcs[noddcyclearcs++] = forwardarcs[i];
972  oddcyclearcs[noddcyclearcs++] = backwardarcs[i];
973  }
974 
975  if( !rhsaggr && SCIPisNegative(scip, coef) )
976  {
977  assert(noddcyclearcs < 2*nbilinexprs-1);
978  oddcyclearcs[noddcyclearcs++] = forwardarcs[i];
979  oddcyclearcs[noddcyclearcs++] = backwardarcs[i];
980  }
981 
982  /* add constraints to ensure no parallel edges */
983  (void) SCIPsnprintf(name, SCIP_MAXSTRLEN, "cons_noparalleledges");
984  SCIP_CALL( SCIPcreateConsBasicLinear(subscip, &noparallelcons, name, 0, NULL, NULL, 0.0, 1.0) );
985  SCIP_CALL( SCIPaddCoefLinear(subscip, noparallelcons, forwardarcs[*narcs], 1.0) );
986  SCIP_CALL( SCIPaddCoefLinear(subscip, noparallelcons, backwardarcs[*narcs], 1.0) );
987  SCIP_CALL( SCIPaddCons(subscip, noparallelcons) );
988  SCIP_CALL( SCIPreleaseCons(subscip, &noparallelcons) );
989 
990  ++*narcs;
991  }
992  }
993  assert(*narcs > 0);
994 
995  /* odd cycle property constraint */
996  (void) SCIPsnprintf(name, SCIP_MAXSTRLEN, "cons_oddcycle");
997  SCIP_CALL( SCIPcreateConsBasicXor(subscip, &oddcyclecons, name, TRUE, noddcyclearcs, oddcyclearcs) );
998  SCIP_CALL( SCIPaddCons(subscip, oddcyclecons) );
999  SCIP_CALL( SCIPreleaseCons(subscip, &oddcyclecons) );
1000  SCIPfreeBufferArray(subscip, &oddcyclearcs);
1001 
1002  /* cycle length constraint */
1003  (void) SCIPsnprintf(name, SCIP_MAXSTRLEN, "cons_cyclelength");
1004  SCIP_CALL( SCIPcreateConsBasicLinear(subscip, &cyclelengthcons, name, 0, NULL, NULL,
1005  (SCIP_Real) sepadata->minaggrsize, (SCIP_Real) sepadata->maxaggrsize) );
1006 
1007  for( i = 0; i < *narcs; ++i )
1008  {
1009  SCIP_CALL( SCIPaddCoefLinear(subscip, cyclelengthcons, forwardarcs[i], 1.0) );
1010  SCIP_CALL( SCIPaddCoefLinear(subscip, cyclelengthcons, backwardarcs[i], 1.0) );
1011  }
1012 
1013  SCIP_CALL( SCIPaddCons(subscip, cyclelengthcons) );
1014  SCIP_CALL( SCIPreleaseCons(subscip, &cyclelengthcons) );
1015 
1016  /* create flow conservation constraints */
1017  SCIP_CALL( SCIPallocBufferArray(subscip, &flowcons, nnodes) );
1018 
1019  for( i = 0; i < nnodes; ++i )
1020  {
1021  (void) SCIPsnprintf(name, SCIP_MAXSTRLEN, "cons_flowconservation#%d", i);
1022  SCIP_CALL( SCIPcreateConsBasicLinear(subscip, &flowcons[i], name, 0, NULL, NULL, 0.0, 0.0) );
1023  }
1024 
1025  arcidx = 0;
1026  for( i = 0; i < nquadexprs; ++i )
1027  {
1028  SCIP_EXPR* qterm;
1029  SCIP_Real coef;
1030  int nadjbilin;
1031  int* adjbilin;
1032  int j;
1033 
1034  SCIPexprGetQuadraticQuadTerm(expr, i, &qterm, NULL, &coef, &nadjbilin, &adjbilin, NULL);
1035 
1036  if( !SCIPisZero(scip, coef) )
1037  ++arcidx;
1038 
1039  for( j = 0 ; j < nadjbilin; ++j )
1040  {
1041  SCIP_EXPR* qterm1;
1042  int pos2;
1043 
1044  SCIPexprGetQuadraticBilinTerm(expr, adjbilin[j], &qterm1, NULL, NULL, &pos2, NULL);
1045 
1046  /* handle qterm == qterm2 later */
1047  if( qterm1 != qterm )
1048  continue;
1049 
1050  SCIP_CALL( SCIPaddCoefLinear(subscip, flowcons[i], forwardarcs[arcidx], 1.0) );
1051  SCIP_CALL( SCIPaddCoefLinear(subscip, flowcons[i], backwardarcs[arcidx], -1.0) );
1052 
1053  SCIP_CALL( SCIPaddCoefLinear(subscip, flowcons[pos2], forwardarcs[arcidx], -1.0) );
1054  SCIP_CALL( SCIPaddCoefLinear(subscip, flowcons[pos2], backwardarcs[arcidx], 1.0) );
1055 
1056  ++arcidx;
1057  }
1058  }
1059  assert(arcidx == *narcs);
1060 
1061  for( i = 0; i < nnodes; ++i )
1062  {
1063  SCIP_CALL( SCIPaddCons(subscip, flowcons[i]) );
1064  SCIP_CALL( SCIPreleaseCons(subscip, &flowcons[i]) );
1065  }
1066 
1067  SCIPfreeBufferArray(subscip, &flowcons);
1068 
1069  return SCIP_OKAY;
1070 }
1071 
1072 /** fixed all arc variables (u,v) for which u or v is already in an edge-concave aggregation */
1073 static
1074 SCIP_RETCODE updateMIP(
1075  SCIP* subscip, /**< auxiliary SCIP to search aggregations */
1076  SCIP_NLROW* nlrow, /**< nonlinear row */
1077  SCIP_VAR** forwardarcs, /**< forward arc variables */
1078  SCIP_VAR** backwardarcs, /**< backward arc variables */
1079  int* quadvar2aggr, /**< mapping of quadvars to e.c. aggr. index (< 0: in no aggr.) */
1080  int* nedges /**< pointer to store the number of nonexcluded edges */
1081  )
1082 {
1083  SCIP_EXPR* expr;
1084  int nquadexprs;
1085  int arcidx;
1086  int i;
1087 
1088  assert(subscip != NULL);
1089  assert(nlrow != NULL);
1090  assert(forwardarcs != NULL);
1091  assert(backwardarcs != NULL);
1092  assert(quadvar2aggr != NULL);
1093  assert(nedges != NULL);
1094 
1095  SCIP_CALL( SCIPfreeTransform(subscip) );
1096 
1097  /* recompute the number of edges */
1098  *nedges = 0;
1099 
1100  expr = SCIPnlrowGetExpr(nlrow);
1101  SCIPexprGetQuadraticData(expr, NULL, NULL, NULL, NULL, &nquadexprs, NULL, NULL, NULL);
1102 
1103  /* fix each arc to 0 if at least one of its nodes is contained in an e.c. aggregation */
1104  arcidx = 0;
1105  for( i = 0; i < nquadexprs; ++i )
1106  {
1107  SCIP_EXPR* qterm;
1108  SCIP_Real coef;
1109  int nadjbilin;
1110  int* adjbilin;
1111  int j;
1112 
1113  SCIPexprGetQuadraticQuadTerm(expr, i, &qterm, NULL, &coef, &nadjbilin, &adjbilin, NULL);
1114 
1115  if( !SCIPisZero(subscip, coef) )
1116  {
1117  if( quadvar2aggr[i] != -1 )
1118  {
1119  SCIP_CALL( SCIPchgVarUb(subscip, forwardarcs[arcidx], 0.0) );
1120  SCIP_CALL( SCIPchgVarUb(subscip, backwardarcs[arcidx], 0.0) );
1121  }
1122  ++arcidx;
1123  }
1124 
1125  for( j = 0 ; j < nadjbilin; ++j )
1126  {
1127  SCIP_EXPR* qterm1;
1128  int pos2;
1129 
1130  SCIPexprGetQuadraticBilinTerm(expr, adjbilin[j], &qterm1, NULL, NULL, &pos2, NULL);
1131 
1132  /* handle qterm == qterm2 later */
1133  if( qterm1 != qterm )
1134  continue;
1135 
1136  if( quadvar2aggr[i] != -1 || quadvar2aggr[pos2] != -1 )
1137  {
1138  SCIP_CALL( SCIPchgVarUb(subscip, forwardarcs[arcidx], 0.0) );
1139  SCIP_CALL( SCIPchgVarUb(subscip, backwardarcs[arcidx], 0.0) );
1140  }
1141  else
1142  ++*nedges;
1143 
1144  ++arcidx;
1145  }
1146  }
1147 
1148  return SCIP_OKAY;
1149 }
1150 
1151 /** stores the best edge-concave aggregation found by the MIP model */
1152 static
1153 SCIP_RETCODE storeAggrFromMIP(
1154  SCIP* subscip, /**< auxiliary SCIP to search aggregations */
1155  SCIP_NLROW* nlrow, /**< nonlinear row */
1156  SCIP_VAR** forwardarcs, /**< forward arc variables */
1157  SCIP_VAR** backwardarcs, /**< backward arc variables */
1158  int* quadvar2aggr, /**< mapping of quadvars to e.c. aggr. index (< 0: in no aggr.) */
1159  int nfoundsofar /**< number of e.c. aggregation found so far */
1160  )
1161 {
1162  SCIP_SOL* sol;
1163  SCIP_EXPR* expr;
1164  int nquadexprs;
1165  int arcidx;
1166  int i;
1167 
1168  assert(subscip != NULL);
1169  assert(nlrow != NULL);
1170  assert(forwardarcs != NULL);
1171  assert(backwardarcs != NULL);
1172  assert(quadvar2aggr != NULL);
1173  assert(nfoundsofar >= 0);
1174  assert(SCIPgetStatus(subscip) < SCIP_STATUS_INFEASIBLE);
1175  assert(SCIPgetNSols(subscip) > 0);
1176 
1177  sol = SCIPgetBestSol(subscip);
1178  assert(sol != NULL);
1179 
1180  expr = SCIPnlrowGetExpr(nlrow);
1181  SCIPexprGetQuadraticData(expr, NULL, NULL, NULL, NULL, &nquadexprs, NULL, NULL, NULL);
1182 
1183  /* fix each arc to 0 if at least one of its nodes is contained in an e.c. aggregation */
1184  arcidx = 0;
1185  for( i = 0; i < nquadexprs; ++i )
1186  {
1187  SCIP_EXPR* qterm;
1188  SCIP_Real coef;
1189  int nadjbilin;
1190  int* adjbilin;
1191  int j;
1192 
1193  SCIPexprGetQuadraticQuadTerm(expr, i, &qterm, NULL, &coef, &nadjbilin, &adjbilin, NULL);
1194 
1195  if( !SCIPisZero(subscip, coef) )
1196  {
1197  if( SCIPisGT(subscip, SCIPgetSolVal(subscip, sol, forwardarcs[arcidx]), 0.5) ||
1198  SCIPisGT(subscip, SCIPgetSolVal(subscip, sol, backwardarcs[arcidx]), 0.5) )
1199  {
1200  assert(quadvar2aggr[i] == -1 || quadvar2aggr[i] == nfoundsofar);
1201  quadvar2aggr[i] = nfoundsofar;
1202  }
1203 
1204  ++arcidx;
1205  }
1206 
1207  for( j = 0; j < nadjbilin; ++j )
1208  {
1209  SCIP_EXPR* qterm1;
1210  int pos2;
1211 
1212  SCIPexprGetQuadraticBilinTerm(expr, adjbilin[j], &qterm1, NULL, NULL, &pos2, NULL);
1213 
1214  /* handle qterm == qterm2 later */
1215  if( qterm1 != qterm )
1216  continue;
1217 
1218  if( SCIPisGT(subscip, SCIPgetSolVal(subscip, sol, forwardarcs[arcidx]), 0.5) ||
1219  SCIPisGT(subscip, SCIPgetSolVal(subscip, sol, backwardarcs[arcidx]), 0.5) )
1220  {
1221  assert(quadvar2aggr[i] == -1 || quadvar2aggr[i] == nfoundsofar);
1222  assert(quadvar2aggr[pos2] == -1 || quadvar2aggr[pos2] == nfoundsofar);
1223 
1224  quadvar2aggr[i] = nfoundsofar;
1225  quadvar2aggr[pos2] = nfoundsofar;
1226  }
1227 
1228  ++arcidx;
1229  }
1230  }
1231 
1232  return SCIP_OKAY;
1233 }
1234 
1235 /** searches for edge-concave aggregations with a MIP model based on binary flow variables */
1236 static
1237 SCIP_RETCODE searchEcAggrWithMIP(
1238  SCIP* subscip, /**< SCIP data structure */
1239  SCIP_Real timelimit, /**< time limit to solve the MIP */
1240  int nedges, /**< number of nonexcluded undirected edges */
1241  SCIP_Bool* aggrleft, /**< pointer to store if there might be a left aggregation */
1242  SCIP_Bool* found /**< pointer to store if we have found an aggregation */
1243  )
1244 {
1245  assert(subscip != NULL);
1246  assert(aggrleft != NULL);
1247  assert(found != NULL);
1248  assert(nedges >= 0);
1249 
1250  *aggrleft = TRUE;
1251  *found = FALSE;
1252 
1253  if( SCIPisLE(subscip, timelimit, 0.0) )
1254  return SCIP_OKAY;
1255 
1256  /* set working limits */
1257  SCIP_CALL( SCIPsetRealParam(subscip, "limits/time", timelimit) );
1258  SCIP_CALL( SCIPsetLongintParam(subscip, "limits/totalnodes", SUBSCIP_NODELIMIT) );
1259 
1260  /* set heuristics to aggressive */
1261  SCIP_CALL( SCIPsetHeuristics(subscip, SCIP_PARAMSETTING_AGGRESSIVE, TRUE) );
1262 
1263  /* disable output to console in optimized mode, enable in SCIP's debug mode */
1264 #ifdef SCIP_DEBUG
1265  SCIP_CALL( SCIPsetIntParam(subscip, "display/verblevel", 5) );
1266  SCIP_CALL( SCIPsetIntParam(subscip, "display/freq", 1) );
1267 #else
1268  SCIP_CALL( SCIPsetIntParam(subscip, "display/verblevel", 0) );
1269 #endif
1270 
1271  SCIP_CALL( SCIPsolve(subscip) );
1272 
1273  /* no more aggregation left if the MIP is infeasible */
1274  if( SCIPgetStatus(subscip) >= SCIP_STATUS_INFEASIBLE )
1275  {
1276  *found = FALSE;
1277  *aggrleft = FALSE;
1278  return SCIP_OKAY;
1279  }
1280 
1281  if( SCIPgetNSols(subscip) > 0 )
1282  {
1283  *found = TRUE;
1284  *aggrleft = TRUE;
1285 
1286 #ifdef SCIP_DEBUG
1287  if( SCIPgetNSols(subscip) > 0 )
1288  {
1289  SCIP_CALL( SCIPprintSol(subscip, SCIPgetBestSol(subscip), NULL , FALSE) );
1290  }
1291 #endif
1292  }
1293 
1294  return SCIP_OKAY;
1295 }
1296 
1297 /** creates a tclique graph from a given nonlinear row
1298  *
1299  * SCIP's clique code can only handle integer node weights; all node weights are scaled by a factor of 100; since the
1300  * clique code ignores nodes with weight of zero, we add an offset of 100 to each weight
1301  */
1302 static
1303 SCIP_RETCODE createTcliqueGraph(
1304  SCIP_NLROW* nlrow, /**< nonlinear row */
1305  TCLIQUE_GRAPH** graph, /**< TCLIQUE graph structure */
1306  SCIP_Real* nodeweights /**< weights for each quadratic variable (nodes in the graph) */
1307  )
1308 {
1309  SCIP_EXPR* expr;
1310  int nquadexprs;
1311  int i;
1312 
1313  assert(graph != NULL);
1314  assert(nlrow != NULL);
1315 
1316  /* create the tclique graph */
1317  if( !tcliqueCreate(graph) )
1318  {
1319  SCIPerrorMessage("could not create clique graph\n");
1320  return SCIP_ERROR;
1321  }
1322 
1323  expr = SCIPnlrowGetExpr(nlrow);
1324  SCIPexprGetQuadraticData(expr, NULL, NULL, NULL, NULL, &nquadexprs, NULL, NULL, NULL);
1325 
1326  /* add all nodes to the tclique graph */
1327  for( i = 0; i < nquadexprs; ++i )
1328  {
1329  int nodeweight;
1330 
1331  /* note: clique code can only handle integer weights */
1332  nodeweight = 100 + (int)(100 * nodeweights[i]);
1333  /* SCIPdebugMsg(scip, "%d (%s): nodeweight %d \n", i, SCIPvarGetName(SCIPnlrowGetQuadVars(nlrow)[i]), nodeweight); */
1334 
1335  if( !tcliqueAddNode(*graph, i, nodeweight) )
1336  {
1337  SCIPerrorMessage("could not add node to clique graph\n");
1338  return SCIP_ERROR;
1339  }
1340  }
1341 
1342  /* add all edges */
1343  for( i = 0; i < nquadexprs; ++i )
1344  {
1345  SCIP_EXPR* qterm;
1346  int nadjbilin;
1347  int* adjbilin;
1348  int j;
1349 
1350  SCIPexprGetQuadraticQuadTerm(expr, i, &qterm, NULL, NULL, &nadjbilin, &adjbilin, NULL);
1351 
1352  for( j = 0; j < nadjbilin; ++j )
1353  {
1354  SCIP_EXPR* qterm1;
1355  SCIP_EXPR* qterm2;
1356  int pos2;
1357 
1358  SCIPexprGetQuadraticBilinTerm(expr, adjbilin[j], &qterm1, &qterm2, NULL, &pos2, NULL);
1359 
1360  /* handle qterm == qterm2 later */
1361  if( qterm1 != qterm )
1362  continue;
1363 
1364 #ifdef SCIP_DEBUG_DETAILED
1365  SCIPdebugMessage(" add edge (%d, %d) = (%s,%s) to tclique graph\n",
1366  SCIPvarGetIndex(SCIPgetVarExprVar(qterm1)), SCIPvarGetIndex(SCIPgetVarExprVar(qterm2)),
1367  SCIPvarGetName(SCIPgetVarExprVar(qterm1)), SCIPvarGetName(SCIPgetVarExprVar(qterm2)));
1368 #endif
1369 
1370  if( !tcliqueAddEdge(*graph, i, pos2) )
1371  {
1372  SCIPerrorMessage("could not add edge to clique graph\n");
1373  return SCIP_ERROR;
1374  }
1375  }
1376  }
1377 
1378  /* flush the clique graph */
1379  if( !tcliqueFlush(*graph) )
1380  {
1381  SCIPerrorMessage("could not flush the clique graph\n");
1382  return SCIP_ERROR;
1383  }
1384 
1385  return SCIP_OKAY;
1386 }
1387 
1388 /** searches for edge-concave aggregations by computing cliques in the graph representation of a given nonlinear row
1389  *
1390  * update graph, compute clique, store clique; after computing a clique we heuristically check if the clique contains
1391  * at least one good cycle
1392  */
1393 static
1394 SCIP_RETCODE searchEcAggrWithCliques(
1395  SCIP* scip, /**< SCIP data structure */
1396  TCLIQUE_GRAPH* graph, /**< TCLIQUE graph structure */
1397  SCIP_SEPADATA* sepadata, /**< separator data */
1398  SCIP_NLROW* nlrow, /**< nonlinear row */
1399  int* quadvar2aggr, /**< mapping of quadvars to e.c. aggr. index (< 0: in no aggr.) */
1400  int nfoundsofar, /**< number of e.c. aggregation found so far */
1401  SCIP_Bool rhsaggr, /**< consider nonlinear row aggregation for g(x) <= rhs (TRUE) or
1402  * lhs <= g(x) (FALSE) */
1403  SCIP_Bool* foundaggr, /**< pointer to store if we have found an aggregation */
1404  SCIP_Bool* foundclique /**< pointer to store if we have found a clique */
1405  )
1406 {
1407  SCIP_HASHMAP* cliquemap;
1408  TCLIQUE_STATUS status;
1409  SCIP_EXPR* expr;
1410  int nquadexprs;
1411  int* maxcliquenodes;
1412  int* degrees;
1413  int nmaxcliquenodes;
1414  int maxcliqueweight;
1415  int noddcycleedges;
1416  int ntwodegrees;
1417  int aggrsize;
1418  int i;
1419 
1420  assert(graph != NULL);
1421  assert(nfoundsofar >= 0);
1422  assert(foundaggr != NULL);
1423  assert(foundclique != NULL);
1424 
1425  cliquemap = NULL;
1426  *foundaggr = FALSE;
1427  *foundclique = FALSE;
1428 
1429  expr = SCIPnlrowGetExpr(nlrow);
1430  SCIPexprGetQuadraticData(expr, NULL, NULL, NULL, NULL, &nquadexprs, NULL, NULL, NULL);
1431  assert(nquadexprs == tcliqueGetNNodes(graph));
1432 
1433  /* exclude all nodes which are already in an edge-concave aggregation (no flush is needed) */
1434  for( i = 0; i < nquadexprs; ++i )
1435  {
1436  if( quadvar2aggr[i] != -1 )
1437  {
1438  SCIPdebugMsg(scip, "exclude node %d from clique graph\n", i);
1439  tcliqueChangeWeight(graph, i, 0);
1440  }
1441  }
1442 
1443  SCIP_CALL( SCIPallocBufferArray(scip, &maxcliquenodes, nquadexprs) );
1444 
1445  /* solve clique problem */
1446  tcliqueMaxClique(tcliqueGetNNodes, tcliqueGetWeights, tcliqueIsEdge, tcliqueSelectAdjnodes, graph, NULL, NULL,
1447  maxcliquenodes, &nmaxcliquenodes, &maxcliqueweight, CLIQUE_MAXFIRSTNODEWEIGHT, CLIQUE_MINWEIGHT,
1448  CLIQUE_MAXNTREENODES, CLIQUE_BACKTRACKFREQ, 0, -1, NULL, &status);
1449 
1450  if( status != TCLIQUE_OPTIMAL || nmaxcliquenodes < sepadata->minaggrsize )
1451  goto TERMINATE;
1452 
1453  *foundclique = TRUE;
1454  aggrsize = MIN(sepadata->maxaggrsize, nmaxcliquenodes);
1455  SCIP_CALL( SCIPhashmapCreate(&cliquemap, SCIPblkmem(scip), aggrsize) );
1456 
1457  for( i = 0; i < aggrsize; ++i )
1458  {
1459  SCIP_CALL( SCIPhashmapInsertInt(cliquemap, (void*) (size_t) maxcliquenodes[i], 0) ); /*lint !e571*/
1460  }
1461 
1462  /* count the degree of good cycle edges for each node in the clique */
1463  SCIP_CALL( SCIPallocBufferArray(scip, &degrees, aggrsize) );
1464  BMSclearMemoryArray(degrees, aggrsize);
1465  ntwodegrees = 0;
1466 
1467  /* count the number of positive or negative edges (depending on <= rhs or >= lhs) */
1468  noddcycleedges = 0;
1469  for( i = 0; i < nquadexprs; ++i )
1470  {
1471  SCIP_Bool isoddcycleedge;
1472  SCIP_EXPR* qterm;
1473  SCIP_Real coef;
1474  int nadjbilin;
1475  int* adjbilin;
1476  int j;
1477 
1478  SCIPexprGetQuadraticQuadTerm(expr, i, &qterm, NULL, &coef, &nadjbilin, &adjbilin, NULL);
1479 
1480  isoddcycleedge = (rhsaggr && SCIPisPositive(scip, coef)) || (!rhsaggr && SCIPisNegative(scip, coef));
1481 
1482  if( isoddcycleedge && SCIPhashmapExists(cliquemap, (void*) (size_t) i) )
1483  {
1484  ++noddcycleedges;
1485  ++degrees[i];
1486  }
1487 
1488  for( j = 0; j < nadjbilin; ++j )
1489  {
1490  SCIP_EXPR* qterm1;
1491  SCIP_EXPR* qterm2;
1492  int pos2;
1493 
1494  SCIPexprGetQuadraticBilinTerm(expr, adjbilin[j], &qterm1, &qterm2, &coef, &pos2, NULL);
1495 
1496  /* handle qterm == qterm2 later */
1497  if( qterm1 != qterm )
1498  continue;
1499 
1500  isoddcycleedge = (rhsaggr && SCIPisPositive(scip, coef)) || (!rhsaggr && SCIPisNegative(scip, coef));
1501 
1502  if( isoddcycleedge
1503  && SCIPhashmapExists(cliquemap, (void*) (size_t) i)
1504  && SCIPhashmapExists(cliquemap, (void*) (size_t) pos2) )
1505  {
1506  ++noddcycleedges;
1507  ++degrees[i];
1508  ++degrees[pos2];
1509  }
1510  }
1511  }
1512 
1513  /* count the number of nodes with exactly two incident odd cycle edges */
1514  for( i = 0; i < aggrsize; ++i )
1515  if( degrees[i] == 2 )
1516  ++ntwodegrees;
1517 
1518  /* check cases for which we are sure that there are no good cycles in the clique */
1519  if( noddcycleedges == 0 || (aggrsize == 3 && noddcycleedges == 2) || (aggrsize == 4 && ntwodegrees == 4) )
1520  *foundaggr = FALSE;
1521  else
1522  *foundaggr = TRUE;
1523 
1524  /* add the found clique as an edge-concave aggregation or exclude the nodes from the remaining search */
1525  for( i = 0; i < aggrsize; ++i )
1526  {
1527  quadvar2aggr[ maxcliquenodes[i] ] = *foundaggr ? nfoundsofar : -2;
1528  SCIPdebugMsg(scip, "%s %d\n", *foundaggr ? "aggregate node: " : "exclude node: ", maxcliquenodes[i]);
1529  }
1530 
1531  SCIPfreeBufferArray(scip, &degrees);
1532 
1533 TERMINATE:
1534  if( cliquemap != NULL )
1535  SCIPhashmapFree(&cliquemap);
1536  SCIPfreeBufferArray(scip, &maxcliquenodes);
1537 
1538  return SCIP_OKAY;
1539 }
1540 
1541 /** helper function for searchEcAggr() */
1542 static
1543 SCIP_RETCODE doSeachEcAggr(
1544  SCIP* scip, /**< SCIP data structure */
1545  SCIP* subscip, /**< sub-SCIP data structure */
1546  SCIP_SEPADATA* sepadata, /**< separator data */
1547  SCIP_NLROW* nlrow, /**< nonlinear row */
1548  SCIP_SOL* sol, /**< current solution (might be NULL) */
1549  SCIP_Bool rhsaggr, /**< consider nonlinear row aggregation for g(x) <= rhs (TRUE) or g(x) >= lhs (FALSE) */
1550  int* quadvar2aggr, /**< array to store for each quadratic variable in which edge-concave
1551  * aggregation it is stored (< 0: in no aggregation); size has to be at
1552  * least SCIPnlrowGetNQuadVars(nlrow) */
1553  int* nfound /**< pointer to store the number of found e.c. aggregations */
1554  )
1555 {
1556  TCLIQUE_GRAPH* graph = NULL;
1557  SCIP_EXPR* expr;
1558  SCIP_VAR** forwardarcs;
1559  SCIP_VAR** backwardarcs;
1560  SCIP_Real* nodeweights;
1561  SCIP_Real timelimit;
1562  SCIP_RETCODE retcode;
1563  int nunsucces = 0;
1564  int nedges = 0;
1565  int narcs;
1566  int nquadvars;
1567  int nbilinexprs;
1568  int i;
1569 
1570  assert(subscip != NULL);
1571  assert(quadvar2aggr != NULL);
1572  assert(nfound != NULL);
1573 
1574  expr = SCIPnlrowGetExpr(nlrow);
1575  SCIPexprGetQuadraticData(expr, NULL, NULL, NULL, NULL, &nquadvars, &nbilinexprs, NULL, NULL);
1576 
1577  retcode = SCIP_OKAY;
1578  *nfound = 0;
1579 
1580  /* arrays to store all arc variables of the MIP model; note that we introduce variables even for loops in the graph
1581  * to have an easy mapping from the edges of the graph to the quadratic elements
1582  * nquadvars + nbilinexprs is an upper bound on the actual number of square and bilinear terms
1583  */
1584  SCIP_CALL( SCIPallocBufferArray(scip, &nodeweights, nquadvars) );
1585  SCIP_CALL( SCIPallocBufferArray(scip, &forwardarcs, nquadvars + nbilinexprs) );
1586  SCIP_CALL( SCIPallocBufferArray(scip, &backwardarcs, nquadvars + nbilinexprs) );
1587 
1588  /* initialize mapping from quadvars to e.c. aggregation index (-1: quadvar is in no aggregation); compute node
1589  * weights
1590  */
1591  for( i = 0; i < nquadvars; ++i )
1592  {
1593  SCIP_EXPR* qterm;
1594  SCIP_VAR* var;
1595 
1596  SCIPexprGetQuadraticQuadTerm(expr, i, &qterm, NULL, NULL, NULL, NULL, NULL);
1597  assert(SCIPisExprVar(scip, qterm));
1598  var = SCIPgetVarExprVar(qterm);
1599 
1600  quadvar2aggr[i] = -1;
1601  nodeweights[i] = phi(scip, SCIPgetSolVal(scip, sol, var), SCIPvarGetLbLocal(var), SCIPvarGetUbLocal(var));
1602  SCIPdebugMsg(scip, "%s = %e (%e in [%e, %e])\n", SCIPvarGetName(var), nodeweights[i], SCIPgetSolVal(scip, sol, var),
1603  SCIPvarGetLbLocal(var), SCIPvarGetUbLocal(var));
1604  }
1605 
1606  SCIP_CALL( createMIP(scip, subscip, sepadata, nlrow, rhsaggr, forwardarcs, backwardarcs, nodeweights, &nedges, &narcs) );
1607  assert(nedges >= 0);
1608  assert(narcs > 0);
1609  SCIPdebugMsg(scip, "nedges (without loops) = %d\n", nedges);
1610  SCIPdebugMsg(scip, "narcs (number of quadratic terms) = %d\n", narcs);
1611 
1612  SCIP_CALL( SCIPgetRealParam(scip, "limits/time", &timelimit) );
1613 
1614  /* main loop to search for edge-concave aggregations */
1615  while( !SCIPisStopped(scip) )
1616  {
1617  SCIP_Bool aggrleft;
1618  SCIP_Bool found;
1619 
1620  SCIPdebugMsg(scip, "#remaining edges = %d\n", nedges);
1621 
1622  /* not enough edges left */
1623  if( nedges < sepadata->minaggrsize )
1624  break;
1625 
1626  /* check whether there is enough time left; update the remaining time */
1627  if( !SCIPisInfinity(scip, timelimit) )
1628  {
1629  timelimit -= SCIPgetSolvingTime(scip);
1630  if( timelimit <= 0.0 )
1631  {
1632  SCIPdebugMsg(scip, "skip aggregation search since no time left\n");
1633  goto TERMINATE;
1634  }
1635  }
1636 
1637  /* 1.a - search for edge-concave aggregation with the help of the MIP model */
1638  SCIP_CALL( searchEcAggrWithMIP(subscip, timelimit, nedges, &aggrleft, &found) );
1639 
1640  /* 1.b - there are no more edge-concave aggregations left */
1641  if( !aggrleft )
1642  {
1643  SCIPdebugMsg(scip, "no more aggregation left\n");
1644  break;
1645  }
1646 
1647  if( found )
1648  {
1649  SCIP_CALL( storeAggrFromMIP(subscip, nlrow, forwardarcs, backwardarcs, quadvar2aggr, *nfound) );
1650  ++(*nfound);
1651  nunsucces = 0;
1652  }
1653  /* try to find an edge-concave aggregation by computing cliques */
1654  else
1655  {
1656  SCIP_Bool foundaggr;
1657  SCIP_Bool foundclique;
1658 
1659  ++nunsucces;
1660 
1661  /* create graph if necessary */
1662  if( graph == NULL )
1663  {
1664  SCIP_CALL_TERMINATE( retcode, createTcliqueGraph(nlrow, &graph, nodeweights), TERMINATE );
1665  }
1666 
1667  /* 2.a - search and store a single edge-concave aggregation by computing a clique with a good cycle */
1668  SCIP_CALL_FINALLY( searchEcAggrWithCliques(scip, graph, sepadata, nlrow, quadvar2aggr, *nfound, rhsaggr,
1669  &foundaggr, &foundclique), tcliqueFree(&graph) );
1670 
1671  if( foundaggr )
1672  {
1673  assert(foundclique);
1674  ++(*nfound);
1675  nunsucces = 0;
1676  }
1677  else
1678  ++nunsucces;
1679 
1680  /* 2.b - no clique of at least minaggrsize size found */
1681  if( !foundclique )
1682  {
1683  assert(!foundaggr);
1684  SCIPdebugMsg(scip, "did not find a clique to exclude -> leave aggregation search\n");
1685  break;
1686  }
1687  }
1688 
1689  /* leave the algorithm if we did not find something for maxstallrounds many iterations */
1690  if( nunsucces >= sepadata->maxstallrounds && *nfound == 0 )
1691  {
1692  SCIPdebugMsg(scip, "did not find an e.c. aggregation for %d iterations\n", nunsucces);
1693  break;
1694  }
1695 
1696  /* exclude all edges used in the last aggregation and nodes found in the clique solution */
1697  SCIP_CALL_FINALLY( updateMIP(subscip, nlrow, forwardarcs, backwardarcs, quadvar2aggr, &nedges), tcliqueFree(&graph) );
1698  }
1699 
1700 TERMINATE:
1701 
1702 #ifdef SCIP_DEBUG
1703  SCIPdebugMsg(scip, "aggregations found:\n");
1704  for( i = 0; i < nquadvars; ++i )
1705  {
1706  SCIPdebugMsg(scip, " %d in %d\n", i, quadvar2aggr[i]);
1707  }
1708 #endif
1709 
1710  /* free clique graph */
1711  if( graph != NULL )
1712  tcliqueFree(&graph);
1713 
1714  /* free sub-SCIP */
1715  for( i = 0; i < narcs; ++i )
1716  {
1717  SCIP_CALL( SCIPreleaseVar(subscip, &forwardarcs[i]) );
1718  SCIP_CALL( SCIPreleaseVar(subscip, &backwardarcs[i]) );
1719  }
1720 
1721  SCIPfreeBufferArray(scip, &backwardarcs);
1722  SCIPfreeBufferArray(scip, &forwardarcs);
1723  SCIPfreeBufferArray(scip, &nodeweights);
1724 
1725  return retcode;
1726 }
1727 
1728 /** computes a partitioning into edge-concave aggregations for a given (quadratic) nonlinear row
1729  *
1730  * Each aggregation has to contain a cycle with an odd number of positive weighted edges (good cycles) in the corresponding graph representation.
1731  * For this we use the following algorithm:
1732  * -# use a MIP model based on binary flow variables to compute good cycles and store the implied subgraphs as an e.c. aggr.
1733  * -# if we find a good cycle, store the implied subgraph, delete it from the graph representation and go to 1)
1734  * -# if the MIP model is infeasible (there are no good cycles), STOP
1735  * -# we compute a large clique C if the MIP model fails (because of working limits, etc)
1736  * -# if we find a good cycle in C, store the implied subgraph of C, delete it from the graph representation and go to 1)
1737  * -# if C is not large enough, STOP
1738  */
1739 static
1740 SCIP_RETCODE searchEcAggr(
1741  SCIP* scip, /**< SCIP data structure */
1742  SCIP_SEPADATA* sepadata, /**< separator data */
1743  SCIP_NLROW* nlrow, /**< nonlinear row */
1744  SCIP_SOL* sol, /**< current solution (might be NULL) */
1745  SCIP_Bool rhsaggr, /**< consider nonlinear row aggregation for g(x) <= rhs (TRUE) or g(x) >= lhs (FALSE) */
1746  int* quadvar2aggr, /**< array to store for each quadratic variable in which edge-concave
1747  * aggregation it is stored (< 0: in no aggregation); size has to be at
1748  * least SCIPnlrowGetNQuadVars(nlrow) */
1749  int* nfound /**< pointer to store the number of found e.c. aggregations */
1750  )
1751 {
1752  SCIP* subscip;
1753  SCIP_RETCODE retcode;
1754 
1755  /* create and set up a sub-SCIP */
1756  SCIP_CALL_FINALLY( SCIPcreate(&subscip), (void)SCIPfree(&subscip) );
1757 
1758  retcode = doSeachEcAggr(scip, subscip, sepadata, nlrow, sol, rhsaggr, quadvar2aggr, nfound);
1759 
1760  SCIP_CALL( SCIPfree(&subscip) );
1761  SCIP_CALL( retcode );
1762 
1763  return SCIP_OKAY;
1764 }
1765 
1766 /** returns whether a given nonlinear row can be used to compute edge-concave aggregations for which their convex
1767  * envelope could dominate the termwise bilinear relaxation
1768  *
1769  * This is the case if there exists at least one cycle with
1770  * an odd number of positive edges in the corresponding graph representation of the nonlinear row.
1771  */
1772 static
1773 SCIP_RETCODE isCandidate(
1774  SCIP* scip, /**< SCIP data structure */
1775  SCIP_SEPADATA* sepadata, /**< separator data */
1776  SCIP_NLROW* nlrow, /**< nonlinear row representation of a nonlinear constraint */
1777  SCIP_Bool* rhscandidate, /**< pointer to store if we should compute edge-concave aggregations for
1778  * the <= rhs case */
1779  SCIP_Bool* lhscandidate /**< pointer to store if we should compute edge-concave aggregations for
1780  * the >= lhs case */
1781  )
1782 {
1783  SCIP_EXPR* expr = NULL;
1784  SCIP_Bool takerow = FALSE;
1785  int nquadvars = 0;
1786  int* degrees;
1787  int ninterestingnodes;
1788  int nposedges;
1789  int nnegedges;
1790  int i;
1791 
1792  assert(rhscandidate != NULL);
1793  assert(lhscandidate != NULL);
1794 
1795  *rhscandidate = TRUE;
1796  *lhscandidate = TRUE;
1797 
1798  /* check whether nlrow is in the NLP, is quadratic in variables, and there are enough quadratic variables */
1799  if( SCIPnlrowIsInNLP(nlrow) && SCIPnlrowGetExpr(nlrow) != NULL )
1800  {
1801  expr = SCIPnlrowGetExpr(nlrow);
1802  SCIP_CALL( SCIPcheckExprQuadratic(scip, expr, &takerow) );
1803  }
1804  if( takerow )
1805  takerow = SCIPexprAreQuadraticExprsVariables(expr);
1806  if( takerow )
1807  {
1808  SCIPexprGetQuadraticData(expr, NULL, NULL, NULL, NULL, &nquadvars, NULL, NULL, NULL);
1809  takerow = nquadvars >= sepadata->minaggrsize;
1810  }
1811  if( !takerow )
1812  {
1813  *rhscandidate = FALSE;
1814  *lhscandidate = FALSE;
1815  return SCIP_OKAY;
1816  }
1817 
1818  /* check for infinite rhs or lhs */
1819  if( SCIPisInfinity(scip, REALABS(SCIPnlrowGetRhs(nlrow))) )
1820  *rhscandidate = FALSE;
1821  if( SCIPisInfinity(scip, REALABS(SCIPnlrowGetLhs(nlrow))) )
1822  *lhscandidate = FALSE;
1823 
1824  SCIP_CALL( SCIPallocClearBufferArray(scip, &degrees, nquadvars) );
1825 
1826  ninterestingnodes = 0;
1827  nposedges = 0;
1828  nnegedges = 0;
1829 
1830  for( i = 0; i < nquadvars; ++i )
1831  {
1832  SCIP_EXPR* qterm;
1833  SCIP_VAR* var1;
1834  int nadjbilin;
1835  int* adjbilin;
1836  int j;
1837 
1838  SCIPexprGetQuadraticQuadTerm(expr, i, &qterm, NULL, NULL, &nadjbilin, &adjbilin, NULL);
1839  assert(SCIPisExprVar(scip, qterm));
1840 
1841  var1 = SCIPgetVarExprVar(qterm);
1842 
1843  /* do not consider global fixed variables */
1844  if( SCIPisEQ(scip, SCIPvarGetLbGlobal(var1), SCIPvarGetUbGlobal(var1)) )
1845  continue;
1846 
1847  for( j = 0; j < nadjbilin; ++j )
1848  {
1849  SCIP_EXPR* qterm1;
1850  SCIP_EXPR* qterm2;
1851  SCIP_VAR* var2;
1852  SCIP_Real coef;
1853  int pos2;
1854 
1855  SCIPexprGetQuadraticBilinTerm(expr, adjbilin[j], &qterm1, &qterm2, &coef, &pos2, NULL);
1856 
1857  if( qterm1 != qterm )
1858  continue;
1859 
1860  var2 = SCIPgetVarExprVar(qterm2);
1861 
1862  /* do not consider loops or global fixed variables */
1863  if( SCIPisEQ(scip, SCIPvarGetLbGlobal(var2), SCIPvarGetUbGlobal(var2)) )
1864  continue;
1865 
1866  ++degrees[i];
1867  ++degrees[pos2];
1868 
1869  /* count the number of nodes with a degree of at least 2 */
1870  if( degrees[i] == 2 )
1871  ++ninterestingnodes;
1872  if( degrees[pos2] == 2 )
1873  ++ninterestingnodes;
1874 
1875  nposedges += SCIPisPositive(scip, coef) ? 1 : 0;
1876  nnegedges += SCIPisNegative(scip, coef) ? 1 : 0;
1877  }
1878  }
1879 
1880  SCIPfreeBufferArray(scip, &degrees);
1881 
1882  SCIPdebugMsg(scip, "nlrow contains: %d edges\n", nposedges + nnegedges);
1883 
1884  /* too many edges, too few edges, or to few nodes with degree at least 2 in the graph */
1885  if( nposedges + nnegedges > sepadata->maxbilinterms || nposedges + nnegedges < sepadata->minaggrsize
1886  || ninterestingnodes < sepadata->minaggrsize )
1887  {
1888  *rhscandidate = FALSE;
1889  *lhscandidate = FALSE;
1890  return SCIP_OKAY;
1891  }
1892 
1893  /* check if there are enough positive/negative edges; for a 3-clique there has to be an odd number of those edges */
1894  if( nposedges == 0 || (nposedges + nnegedges == 3 && (nposedges % 2) == 0) )
1895  *rhscandidate = FALSE;
1896  if( nnegedges == 0 || (nposedges + nnegedges == 3 && (nnegedges % 2) == 0) )
1897  *lhscandidate = FALSE;
1898 
1899  return SCIP_OKAY;
1900 }
1901 
1902 /** finds and stores edge-concave aggregations for a given nonlinear row */
1903 static
1904 SCIP_RETCODE findAndStoreEcAggregations(
1905  SCIP* scip, /**< SCIP data structure */
1906  SCIP_SEPADATA* sepadata, /**< separator data */
1907  SCIP_NLROW* nlrow, /**< nonlinear row */
1908  SCIP_SOL* sol /**< current solution (might be NULL) */
1909  )
1910 {
1911  int nquadvars;
1912  int* quadvar2aggr;
1913  SCIP_Bool rhscandidate;
1914  SCIP_Bool lhscandidate;
1915 
1916  assert(scip != NULL);
1917  assert(nlrow != NULL);
1918  assert(sepadata != NULL);
1919 
1920 #ifdef SCIP_DEBUG
1921  SCIPdebugMsg(scip, "search for edge-concave aggregation for the nonlinear row: \n");
1922  SCIP_CALL( SCIPprintNlRow(scip, nlrow, NULL) );
1923 #endif
1924 
1925  /* check obvious conditions for existing cycles with an odd number of positive/negative edges */
1926  SCIP_CALL( isCandidate(scip, sepadata, nlrow, &rhscandidate, &lhscandidate) );
1927  SCIPdebugMsg(scip, "rhs candidate = %u lhs candidate = %u\n", rhscandidate, lhscandidate);
1928 
1929  if( !rhscandidate && !lhscandidate )
1930  return SCIP_OKAY;
1931 
1932  SCIPexprGetQuadraticData(SCIPnlrowGetExpr(nlrow), NULL, NULL, NULL, NULL, &nquadvars, NULL, NULL, NULL);
1933  SCIP_CALL( SCIPallocBufferArray(scip, &quadvar2aggr, nquadvars) ); /*lint !e705*/
1934 
1935  /* search for edge-concave aggregations (consider <= rhs) */
1936  if( rhscandidate )
1937  {
1938  SCIP_NLROWAGGR* nlrowaggr;
1939  int nfound;
1940 
1941  assert(!SCIPisInfinity(scip, REALABS(SCIPnlrowGetRhs(nlrow))));
1942 
1943  SCIPdebugMsg(scip, "consider <= rhs\n");
1944  SCIP_CALL( searchEcAggr(scip, sepadata, nlrow, sol, TRUE, quadvar2aggr, &nfound) );
1945 
1946  if( nfound > 0 )
1947  {
1948  SCIP_CALL( nlrowaggrCreate(scip, nlrow, &nlrowaggr, quadvar2aggr, nfound, TRUE) );
1949  assert(nlrow != NULL);
1950  SCIPdebug(nlrowaggrPrint(scip, nlrowaggr));
1951  SCIP_CALL( sepadataAddNlrowaggr(scip, sepadata, nlrowaggr) );
1952  }
1953  }
1954 
1955  /* search for edge-concave aggregations (consider <= lhs) */
1956  if( lhscandidate )
1957  {
1958  SCIP_NLROWAGGR* nlrowaggr;
1959  int nfound;
1960 
1961  assert(!SCIPisInfinity(scip, REALABS(SCIPnlrowGetLhs(nlrow))));
1962 
1963  SCIPdebugMsg(scip, "consider >= lhs\n");
1964  SCIP_CALL( searchEcAggr(scip, sepadata, nlrow, sol, FALSE, quadvar2aggr, &nfound) );
1965 
1966  if( nfound > 0 )
1967  {
1968  SCIP_CALL( nlrowaggrCreate(scip, nlrow, &nlrowaggr, quadvar2aggr, nfound, FALSE) );
1969  assert(nlrow != NULL);
1970  SCIPdebug(nlrowaggrPrint(scip, nlrowaggr));
1971  SCIP_CALL( sepadataAddNlrowaggr(scip, sepadata, nlrowaggr) );
1972  }
1973  }
1974 
1975  SCIPfreeBufferArray(scip, &quadvar2aggr);
1976  return SCIP_OKAY;
1977 }
1978 
1979 /*
1980  * methods to compute edge-concave cuts
1981  */
1982 
1983 #ifdef SCIP_DEBUG
1984 /** prints a given facet (candidate) */
1985 static
1986 void printFacet(
1987  SCIP* scip, /**< SCIP data structure */
1988  SCIP_VAR** vars, /**< variables contained in the edge-concave aggregation */
1989  int nvars, /**< number of variables contained in the edge-concave aggregation */
1990  SCIP_Real* facet, /**< current facet candidate */
1991  SCIP_Real facetval /**< facet evaluated at the current solution */
1992  )
1993 {
1994  int i;
1995 
1996  SCIPdebugMsg(scip, "print facet (val=%e): ", facetval);
1997  for( i = 0; i < nvars; ++i )
1998  SCIPdebugMsgPrint(scip, "%e %s + ", facet[i], SCIPvarGetName(vars[i]));
1999  SCIPdebugMsgPrint(scip, "%e\n", facet[nvars]);
2000 }
2001 #endif
2002 
2003 /** checks if a facet is really an underestimate for all corners of the domain [l,u]
2004  *
2005  * Because of numerics it can happen that a facet violates a corner of the domain.
2006  * To make the facet valid we subtract the maximum violation from the constant part of the facet.
2007  */
2008 static
2009 SCIP_Bool checkRikun(
2010  SCIP* scip, /**< SCIP data structure */
2011  SCIP_ECAGGR* ecaggr, /**< edge-concave aggregation data */
2012  SCIP_Real* fvals, /**< array containing all corner values of the aggregation */
2013  SCIP_Real* facet /**< current facet candidate (of dimension ecaggr->nvars + 1) */
2014  )
2015 {
2016  SCIP_Real maxviolation;
2017  SCIP_Real val;
2018  unsigned int i;
2019  unsigned int ncorner;
2020  unsigned int prev;
2021 
2022  assert(scip != NULL);
2023  assert(ecaggr != NULL);
2024  assert(fvals != NULL);
2025  assert(facet != NULL);
2026 
2027  ncorner = (unsigned int) poweroftwo[ecaggr->nvars];
2028  maxviolation = 0.0;
2029 
2030  /* check for the origin */
2031  val = facet[ecaggr->nvars];
2032  for( i = 0; i < (unsigned int) ecaggr->nvars; ++i )
2033  val += facet[i] * SCIPvarGetLbLocal(ecaggr->vars[i]);
2034 
2035  /* update maximum violation */
2036  maxviolation = MAX(val - fvals[0], maxviolation);
2037  assert(SCIPisFeasEQ(scip, maxviolation, 0.0));
2038 
2039  prev = 0;
2040  for( i = 1; i < ncorner; ++i )
2041  {
2042  unsigned int gray;
2043  unsigned int diff;
2044  unsigned int pos;
2045 
2046  gray = i ^ (i >> 1);
2047  diff = gray ^ prev;
2048 
2049  /* compute position of unique 1 of diff */
2050  pos = 0;
2051  while( (diff >>= 1) != 0 )
2052  ++pos;
2053 
2054  if( gray > prev )
2055  val += facet[pos] * (SCIPvarGetUbLocal(ecaggr->vars[pos]) - SCIPvarGetLbLocal(ecaggr->vars[pos]));
2056  else
2057  val -= facet[pos] * (SCIPvarGetUbLocal(ecaggr->vars[pos]) - SCIPvarGetLbLocal(ecaggr->vars[pos]));
2058 
2059  /* update maximum violation */
2060  maxviolation = MAX(val - fvals[gray], maxviolation);
2061  assert(SCIPisFeasEQ(scip, maxviolation, 0.0));
2062 
2063  prev = gray;
2064  }
2065 
2066  SCIPdebugMsg(scip, "maximum violation of facet: %2.8e\n", maxviolation);
2067 
2068  /* there seem to be numerical problems if the violation is too large; in this case we reject the facet */
2069  if( maxviolation > ADJUSTFACETTOL )
2070  return FALSE;
2071 
2072  /* adjust constant part of the facet */
2073  facet[ecaggr->nvars] -= maxviolation;
2074 
2075  return TRUE;
2076 }
2077 
2078 /** set up LP interface to solve LPs to compute the facet of the convex envelope */
2079 static
2080 SCIP_RETCODE createLP(
2081  SCIP* scip, /**< SCIP data structure */
2082  SCIP_SEPADATA* sepadata /**< separation data */
2083  )
2084 {
2085  SCIP_Real* obj;
2086  SCIP_Real* lb;
2087  SCIP_Real* ub;
2088  SCIP_Real* val;
2089  int* beg;
2090  int* ind;
2091  int nnonz;
2092  int ncols;
2093  int nrows;
2094  int i;
2095  int k;
2096 
2097  assert(scip != NULL);
2098  assert(sepadata != NULL);
2099  assert(sepadata->nnlrowaggrs > 0);
2100 
2101  /* LP interface has been already created with enough rows/columns*/
2102  if( sepadata->lpi != NULL && sepadata->lpisize >= sepadata->maxecsize )
2103  return SCIP_OKAY;
2104 
2105  /* size of lpi is too small; reconstruct lpi */
2106  if( sepadata->lpi != NULL )
2107  {
2108  SCIP_CALL( SCIPlpiFree(&sepadata->lpi) );
2109  sepadata->lpi = NULL;
2110  }
2111 
2112  assert(sepadata->lpi == NULL);
2113  SCIP_CALL( SCIPlpiCreate(&(sepadata->lpi), SCIPgetMessagehdlr(scip), "e.c. LP", SCIP_OBJSEN_MINIMIZE) );
2114  sepadata->lpisize = sepadata->maxecsize;
2115 
2116  nrows = sepadata->maxecsize + 1;
2117  ncols = poweroftwo[nrows - 1];
2118  nnonz = (ncols * (nrows + 1)) / 2;
2119  k = 0;
2120 
2121  /* allocate necessary memory */
2122  SCIP_CALL( SCIPallocBufferArray(scip, &obj, ncols) );
2123  SCIP_CALL( SCIPallocBufferArray(scip, &lb, ncols) );
2124  SCIP_CALL( SCIPallocBufferArray(scip, &ub, ncols) );
2125  SCIP_CALL( SCIPallocBufferArray(scip, &beg, ncols) );
2126  SCIP_CALL( SCIPallocBufferArray(scip, &val, nnonz) );
2127  SCIP_CALL( SCIPallocBufferArray(scip, &ind, nnonz) );
2128 
2129  /* calculate nonzero entries in the LP; set obj, lb, and ub to zero */
2130  for( i = 0; i < ncols; ++i )
2131  {
2132  int row;
2133  int a;
2134 
2135  obj[i] = 0.0;
2136  lb[i] = 0.0;
2137  ub[i] = 0.0;
2138 
2139  SCIPdebugMsg(scip, "col %i starts at position %d\n", i, k);
2140  beg[i] = k;
2141  row = 0;
2142  a = 1;
2143 
2144  /* iterate through the bit representation of i */
2145  while( a <= i )
2146  {
2147  if( (a & i) != 0 )
2148  {
2149  val[k] = 1.0;
2150  ind[k] = row;
2151 
2152  SCIPdebugMsg(scip, " val[%d][%d] = 1 (position %d)\n", row, i, k);
2153 
2154  ++k;
2155  }
2156 
2157  a <<= 1; /*lint !e701*/
2158  ++row;
2159  assert(poweroftwo[row] == a);
2160  }
2161 
2162  /* put 1 as a coefficient for sum_{i} \lambda_i = 1 row (last row) */
2163  val[k] = 1.0;
2164  ind[k] = nrows - 1;
2165  ++k;
2166  SCIPdebugMsg(scip, " val[%d][%d] = 1 (position %d)\n", nrows - 1, i, k);
2167  }
2168  assert(k == nnonz);
2169 
2170  /*
2171  * add all columns to the LP interface
2172  * CPLEX needs the row to exist before adding columns, so we create the rows with dummy sides
2173  * note that the assert is not needed once somebody fixes the LPI
2174  */
2175  assert(nrows <= ncols);
2176  SCIP_CALL( SCIPlpiAddRows(sepadata->lpi, nrows, obj, obj, NULL, 0, NULL, NULL, NULL) );
2177  SCIP_CALL( SCIPlpiAddCols(sepadata->lpi, ncols, obj, lb, ub, NULL, nnonz, beg, ind, val) );
2178 
2179  /* free allocated memory */
2180  SCIPfreeBufferArray(scip, &ind);
2181  SCIPfreeBufferArray(scip, &val);
2182  SCIPfreeBufferArray(scip, &beg);
2183  SCIPfreeBufferArray(scip, &ub);
2184  SCIPfreeBufferArray(scip, &lb);
2185  SCIPfreeBufferArray(scip, &obj);
2186 
2187  return SCIP_OKAY;
2188 }
2189 
2190 /** evaluates an edge-concave aggregation at a corner of the domain [l,u] */
2191 static
2192 SCIP_Real evalCorner(
2193  SCIP_ECAGGR* ecaggr, /**< edge-concave aggregation data */
2194  int k /**< k-th corner */
2195  )
2196 {
2197  SCIP_Real val;
2198  int i;
2199 
2200  assert(ecaggr != NULL);
2201  assert(k >= 0 && k < poweroftwo[ecaggr->nvars]);
2202 
2203  val = 0.0;
2204 
2205  for( i = 0; i < ecaggr->nterms; ++i )
2206  {
2207  SCIP_Real coef;
2208  SCIP_Real bound1;
2209  SCIP_Real bound2;
2210  int idx1;
2211  int idx2;
2212 
2213  idx1 = ecaggr->termvars1[i];
2214  idx2 = ecaggr->termvars2[i];
2215  coef = ecaggr->termcoefs[i];
2216  assert(idx1 >= 0 && idx1 < ecaggr->nvars);
2217  assert(idx2 >= 0 && idx2 < ecaggr->nvars);
2218 
2219  bound1 = ((poweroftwo[idx1]) & k) == 0 ? SCIPvarGetLbLocal(ecaggr->vars[idx1]) : SCIPvarGetUbLocal(ecaggr->vars[idx1]); /*lint !e661*/
2220  bound2 = ((poweroftwo[idx2]) & k) == 0 ? SCIPvarGetLbLocal(ecaggr->vars[idx2]) : SCIPvarGetUbLocal(ecaggr->vars[idx2]); /*lint !e661*/
2221 
2222  val += coef * bound1 * bound2;
2223  }
2224 
2225  return val;
2226 }
2227 
2228 /** returns (val - lb) / (ub - lb) for a in [lb, ub] */
2229 static
2230 SCIP_Real transformValue(
2231  SCIP* scip, /**< SCIP data structure */
2232  SCIP_Real lb, /**< lower bound */
2233  SCIP_Real ub, /**< upper bound */
2234  SCIP_Real val /**< value in [lb,ub] */
2235  )
2236 {
2237  assert(scip != NULL);
2238  assert(!SCIPisInfinity(scip, -lb));
2239  assert(!SCIPisInfinity(scip, ub));
2240  assert(!SCIPisInfinity(scip, REALABS(val)));
2241  assert(!SCIPisFeasEQ(scip, ub - lb, 0.0)); /* this would mean that a variable has been fixed */
2242 
2243  /* adjust val */
2244  val = MIN(val, ub);
2245  val = MAX(val, lb);
2246 
2247  val = (val - lb) / (ub - lb);
2248  assert(val >= 0.0 && val <= 1.0);
2249 
2250  return val;
2251 }
2252 
2253 /** computes a facet of the convex envelope of an edge concave aggregation
2254  *
2255  * The algorithm solves the following LP:
2256  * \f{align}{
2257  * \min & \sum_i \lambda_i f(v_i)\\
2258  * s.t. & \sum_i \lambda_i v_i = x\\
2259  * & \sum_i \lambda_i = 1\\
2260  * & \lambda \geq 0
2261  * \f}
2262  * where \f$f\f$ is an edge concave function, \f$x\in [l,u]\f$ is a solution of the current relaxation, and \f$v_i\f$ are the vertices of \f$[l,u]\f$.
2263  * The method transforms the problem to the domain \f$[0,1]^n\f$, computes a facet, and transforms this facet to the
2264  * original space. The dual solution of the LP above are the coefficients of the facet.
2265  *
2266  * The complete algorithm works as follows:
2267  * -# compute \f$f(v_i)\f$ for each corner \f$v_i\f$ of \f$[l,u]\f$
2268  * -# set up the described LP for the transformed space
2269  * -# solve the LP and store the resulting facet for the transformed space
2270  * -# transform the facet to original space
2271  * -# adjust and check facet with the algorithm of Rikun et al.
2272  */
2273 static
2274 SCIP_RETCODE computeConvexEnvelopeFacet(
2275  SCIP* scip, /**< SCIP data structure */
2276  SCIP_SEPADATA* sepadata, /**< separation data */
2277  SCIP_SOL* sol, /**< solution (might be NULL) */
2278  SCIP_ECAGGR* ecaggr, /**< edge-concave aggregation data */
2279  SCIP_Real* facet, /**< array to store the coefficients of the resulting facet; size has to be at least (ecaggr->nvars + 1) */
2280  SCIP_Real* facetval, /**< pointer to store the value of the facet evaluated at the current solution */
2281  SCIP_Bool* success /**< pointer to store if we have found a facet */
2282  )
2283 {
2284  SCIP_Real* fvals;
2285  SCIP_Real* side;
2286  SCIP_Real* lb;
2287  SCIP_Real* ub;
2288  SCIP_Real perturbation;
2289  int* inds;
2290  int ncorner;
2291  int ncols;
2292  int nrows;
2293  int i;
2294 
2295  assert(scip != NULL);
2296  assert(sepadata != NULL);
2297  assert(ecaggr != NULL);
2298  assert(facet != NULL);
2299  assert(facetval != NULL);
2300  assert(success != NULL);
2301  assert(ecaggr->nvars <= sepadata->maxecsize);
2302 
2303  *facetval = -SCIPinfinity(scip);
2304  *success = FALSE;
2305 
2306  /* create LP if this has not been done yet */
2307  SCIP_CALL( createLP(scip, sepadata) );
2308 
2309  assert(sepadata->lpi != NULL);
2310  assert(sepadata->lpisize >= ecaggr->nvars);
2311 
2312  SCIP_CALL( SCIPlpiGetNCols(sepadata->lpi, &ncols) );
2313  SCIP_CALL( SCIPlpiGetNRows(sepadata->lpi, &nrows) );
2314  ncorner = poweroftwo[ecaggr->nvars];
2315 
2316  assert(ncorner <= ncols);
2317  assert(ecaggr->nvars + 1 <= nrows);
2318  assert(nrows <= ncols);
2319 
2320  /* allocate necessary memory */
2321  SCIP_CALL( SCIPallocBufferArray(scip, &fvals, ncols) );
2322  SCIP_CALL( SCIPallocBufferArray(scip, &inds, ncols) );
2323  SCIP_CALL( SCIPallocBufferArray(scip, &lb, ncols) );
2324  SCIP_CALL( SCIPallocBufferArray(scip, &ub, ncols) );
2325  SCIP_CALL( SCIPallocBufferArray(scip, &side, ncols) );
2326 
2327  /*
2328  * 1. compute f(v_i) for each corner v_i of [l,u]
2329  * 2. set up the described LP for the transformed space
2330  */
2331  for( i = 0; i < ncols; ++i )
2332  {
2333  fvals[i] = i < ncorner ? evalCorner(ecaggr, i) : 0.0;
2334  inds[i] = i;
2335 
2336  /* update bounds; fix variables to zero which are currently not in the LP */
2337  lb[i] = 0.0;
2338  ub[i] = i < ncorner ? 1.0 : 0.0;
2339  SCIPdebugMsg(scip, "bounds of LP col %d = [%e, %e]; obj = %e\n", i, lb[i], ub[i], fvals[i]);
2340  }
2341 
2342  /* update lhs and rhs */
2343  perturbation = 0.001;
2344  for( i = 0; i < nrows; ++i )
2345  {
2346  /* note that the last row corresponds to sum_{j} \lambda_j = 1 */
2347  if( i < ecaggr->nvars )
2348  {
2349  SCIP_VAR* x;
2350 
2351  x = ecaggr->vars[i];
2352  assert(x != NULL);
2353 
2354  side[i] = transformValue(scip, SCIPvarGetLbLocal(x), SCIPvarGetUbLocal(x), SCIPgetSolVal(scip, sol, x));
2355 
2356  /* perturb point to enforce an LP solution with ecaggr->nvars + 1 nonzero */
2357  side[i] += side[i] > perturbation ? -perturbation : perturbation;
2358  perturbation /= 1.2;
2359  }
2360  else
2361  {
2362  side[i] = (i == nrows - 1) ? 1.0 : 0.0;
2363  }
2364 
2365  SCIPdebugMsg(scip, "LP row %d in [%e, %e]\n", i, side[i], side[i]);
2366  }
2367 
2368  /* update LP */
2369  SCIP_CALL( SCIPlpiChgObj(sepadata->lpi, ncols, inds, fvals) );
2370  SCIP_CALL( SCIPlpiChgBounds(sepadata->lpi, ncols, inds, lb, ub) );
2371  SCIP_CALL( SCIPlpiChgSides(sepadata->lpi, nrows, inds, side, side) );
2372 
2373  /* free memory used to build the LP */
2374  SCIPfreeBufferArray(scip, &side);
2375  SCIPfreeBufferArray(scip, &ub);
2376  SCIPfreeBufferArray(scip, &lb);
2377  SCIPfreeBufferArray(scip, &inds);
2378 
2379  /*
2380  * 3. solve the LP and store the resulting facet for the transformed space
2381  */
2382  if( USEDUALSIMPLEX ) /*lint !e774 !e506*/
2383  {
2384  SCIP_CALL( SCIPlpiSolveDual(sepadata->lpi) );
2385  }
2386  else
2387  {
2388  SCIP_CALL( SCIPlpiSolvePrimal(sepadata->lpi) );
2389  }
2390 
2391  /* the dual solution corresponds to the coefficients of the facet in the transformed problem; note that it might be
2392  * the case that the dual solution has more components than the facet array
2393  */
2394  if( ecaggr->nvars + 1 == ncols )
2395  {
2396  SCIP_CALL( SCIPlpiGetSol(sepadata->lpi, NULL, NULL, facet, NULL, NULL) );
2397  }
2398  else
2399  {
2400  SCIP_Real* dualsol;
2401 
2402  SCIP_CALL( SCIPallocBufferArray(scip, &dualsol, nrows) );
2403 
2404  /* get the dual solution */
2405  SCIP_CALL( SCIPlpiGetSol(sepadata->lpi, NULL, NULL, dualsol, NULL, NULL) );
2406 
2407  for( i = 0; i < ecaggr->nvars; ++i )
2408  facet[i] = dualsol[i];
2409 
2410  /* constant part of the facet is the last component of the dual solution */
2411  facet[ecaggr->nvars] = dualsol[nrows - 1];
2412 
2413  SCIPfreeBufferArray(scip, &dualsol);
2414  }
2415 
2416 #ifdef SCIP_DEBUG
2417  SCIPdebugMsg(scip, "facet for the transformed problem: ");
2418  for( i = 0; i < ecaggr->nvars; ++i )
2419  {
2420  SCIPdebugMsgPrint(scip, "%3.4e * %s + ", facet[i], SCIPvarGetName(ecaggr->vars[i]));
2421  }
2422  SCIPdebugMsgPrint(scip, "%3.4e\n", facet[ecaggr->nvars]);
2423 #endif
2424 
2425  /*
2426  * 4. transform the facet to original space
2427  * we now have the linear underestimator L(x) = beta^T x + beta_0, which needs to be transform to the original space
2428  * the underestimator in the original space, G(x) = alpha^T x + alpha_0, is given by G(x) = L(T(x)), where T(.) is
2429  * the transformation applied in step 2; therefore,
2430  * alpha_i = beta_i/(ub_i - lb_i)
2431  * alpha_0 = beta_0 - sum_i lb_i * beta_i/(ub_i - lb_i)
2432  */
2433 
2434  SCIPdebugMsg(scip, "facet in orig. space: ");
2435  *facetval = 0.0;
2436 
2437  for( i = 0; i < ecaggr->nvars; ++i )
2438  {
2439  SCIP_Real varlb;
2440  SCIP_Real varub;
2441 
2442  varlb = SCIPvarGetLbLocal(ecaggr->vars[i]);
2443  varub = SCIPvarGetUbLocal(ecaggr->vars[i]);
2444  assert(!SCIPisEQ(scip, varlb, varub));
2445 
2446  /* substract (\beta_i * lb_i) / (ub_i - lb_i) from current alpha_0 */
2447  facet[ecaggr->nvars] -= (facet[i] * varlb) / (varub - varlb);
2448 
2449  /* set \alpha_i := \beta_i / (ub_i - lb_i) */
2450  facet[i] = facet[i] / (varub - varlb);
2451  *facetval += facet[i] * SCIPgetSolVal(scip, sol, ecaggr->vars[i]);
2452 
2453  SCIPdebugMsgPrint(scip, "%3.4e * %s + ", facet[i], SCIPvarGetName(ecaggr->vars[i]));
2454  }
2455 
2456  /* add constant part to the facet value */
2457  *facetval += facet[ecaggr->nvars];
2458  SCIPdebugMsgPrint(scip, "%3.4e\n", facet[ecaggr->nvars]);
2459 
2460  /*
2461  * 5. adjust and check facet with the algorithm of Rikun et al.
2462  */
2463 
2464  if( checkRikun(scip, ecaggr, fvals, facet) )
2465  {
2466  SCIPdebugMsg(scip, "facet pass the check of Rikun et al.\n");
2467  *success = TRUE;
2468  }
2469 
2470  /* free allocated memory */
2471  SCIPfreeBufferArray(scip, &fvals);
2472 
2473  return SCIP_OKAY;
2474 }
2475 
2476 /*
2477  * miscellaneous methods
2478  */
2479 
2480 /** method to add a facet of the convex envelope of an edge-concave aggregation to a given cut */
2481 static
2482 SCIP_RETCODE addFacetToCut(
2483  SCIP* scip, /**< SCIP data structure */
2484  SCIP_SOL* sol, /**< current solution (might be NULL) */
2485  SCIP_ROW* cut, /**< current cut (modifiable) */
2486  SCIP_Real* facet, /**< coefficient of the facet (dimension nvars + 1) */
2487  SCIP_VAR** vars, /**< variables of the facet */
2488  int nvars, /**< number of variables in the facet */
2489  SCIP_Real* cutconstant, /**< pointer to update the constant part of the facet */
2490  SCIP_Real* cutactivity, /**< pointer to update the activity of the cut */
2491  SCIP_Bool* success /**< pointer to store if everything went fine */
2492  )
2493 {
2494  int i;
2495 
2496  assert(cut != NULL);
2497  assert(facet != NULL);
2498  assert(vars != NULL);
2499  assert(nvars > 0);
2500  assert(cutconstant != NULL);
2501  assert(cutactivity != NULL);
2502  assert(success != NULL);
2503 
2504  *success = TRUE;
2505 
2506  for( i = 0; i < nvars; ++i )
2507  {
2508  if( SCIPisInfinity(scip, REALABS(facet[i])) )
2509  {
2510  *success = FALSE;
2511  return SCIP_OKAY;
2512  }
2513 
2514  if( !SCIPisZero(scip, facet[i]) )
2515  {
2516  /* add only a constant if the variable has been fixed */
2517  if( SCIPvarGetLbLocal(vars[i]) == SCIPvarGetUbLocal(vars[i]) ) /*lint !e777*/
2518  {
2519  assert(SCIPisFeasEQ(scip, SCIPvarGetLbLocal(vars[i]), SCIPgetSolVal(scip, sol, vars[i])));
2520  *cutconstant += facet[i] * SCIPgetSolVal(scip, sol, vars[i]);
2521  *cutactivity += facet[i] * SCIPgetSolVal(scip, sol, vars[i]);
2522  }
2523  else
2524  {
2525  *cutactivity += facet[i] * SCIPgetSolVal(scip, sol, vars[i]);
2526  SCIP_CALL( SCIPaddVarToRow(scip, cut, vars[i], facet[i]) );
2527  }
2528  }
2529  }
2530 
2531  /* add constant part of the facet */
2532  *cutconstant += facet[nvars];
2533  *cutactivity += facet[nvars];
2534 
2535  return SCIP_OKAY;
2536 }
2537 
2538 /** method to add a linear term to a given cut */
2539 static
2540 SCIP_RETCODE addLinearTermToCut(
2541  SCIP* scip, /**< SCIP data structure */
2542  SCIP_SOL* sol, /**< current solution (might be NULL) */
2543  SCIP_ROW* cut, /**< current cut (modifiable) */
2544  SCIP_VAR* x, /**< linear variable */
2545  SCIP_Real coeff, /**< coefficient */
2546  SCIP_Real* cutconstant, /**< pointer to update the constant part of the facet */
2547  SCIP_Real* cutactivity, /**< pointer to update the activity of the cut */
2548  SCIP_Bool* success /**< pointer to store if everything went fine */
2549  )
2550 {
2551  SCIP_Real activity;
2552 
2553  assert(cut != NULL);
2554  assert(x != NULL);
2555  assert(!SCIPisZero(scip, coeff));
2556  assert(!SCIPisInfinity(scip, coeff));
2557  assert(cutconstant != NULL);
2558  assert(cutactivity != NULL);
2559  assert(success != NULL);
2560 
2561  *success = TRUE;
2562  activity = SCIPgetSolVal(scip, sol, x) * coeff;
2563 
2564  /* do not add a term if the activity is -infinity */
2565  if( SCIPisInfinity(scip, -1.0 * REALABS(activity)) )
2566  {
2567  *success = FALSE;
2568  return SCIP_OKAY;
2569  }
2570 
2571  /* add activity to the constant part if the variable has been fixed */
2572  if( SCIPvarGetLbLocal(x) == SCIPvarGetUbLocal(x) ) /*lint !e777*/
2573  {
2574  assert(SCIPisFeasEQ(scip, SCIPvarGetLbLocal(x), SCIPgetSolVal(scip, sol, x)));
2575  *cutconstant += activity;
2576  SCIPdebugMsg(scip, "add to cut: %e\n", activity);
2577  }
2578  else
2579  {
2580  SCIP_CALL( SCIPaddVarToRow(scip, cut, x, coeff) );
2581  SCIPdebugMsg(scip, "add to cut: %e * %s\n", coeff, SCIPvarGetName(x));
2582  }
2583 
2584  *cutactivity += activity;
2585 
2586  return SCIP_OKAY;
2587 }
2588 
2589 /** method to add an underestimate of a bilinear term to a given cut */
2590 static
2591 SCIP_RETCODE addBilinearTermToCut(
2592  SCIP* scip, /**< SCIP data structure */
2593  SCIP_SOL* sol, /**< current solution (might be NULL) */
2594  SCIP_ROW* cut, /**< current cut (modifiable) */
2595  SCIP_VAR* x, /**< first bilinear variable */
2596  SCIP_VAR* y, /**< seconds bilinear variable */
2597  SCIP_Real coeff, /**< coefficient */
2598  SCIP_Real* cutconstant, /**< pointer to update the constant part of the facet */
2599  SCIP_Real* cutactivity, /**< pointer to update the activity of the cut */
2600  SCIP_Bool* success /**< pointer to store if everything went fine */
2601  )
2602 {
2603  SCIP_Real activity;
2604 
2605  assert(cut != NULL);
2606  assert(x != NULL);
2607  assert(y != NULL);
2608  assert(!SCIPisZero(scip, coeff));
2609  assert(cutconstant != NULL);
2610  assert(cutactivity != NULL);
2611  assert(success != NULL);
2612 
2613  *success = TRUE;
2614  activity = coeff * SCIPgetSolVal(scip, sol, x) * SCIPgetSolVal(scip, sol, y);
2615 
2616  if( SCIPisInfinity(scip, REALABS(coeff)) )
2617  {
2618  *success = FALSE;
2619  return SCIP_OKAY;
2620  }
2621 
2622  /* do not add a term if the activity is -infinity */
2623  if( SCIPisInfinity(scip, -1.0 * REALABS(activity)) )
2624  {
2625  *success = FALSE;
2626  return SCIP_OKAY;
2627  }
2628 
2629  /* quadratic case */
2630  if( x == y )
2631  {
2632  SCIP_Real refpoint;
2633  SCIP_Real lincoef;
2634  SCIP_Real linconst;
2635 
2636  lincoef = 0.0;
2637  linconst = 0.0;
2638  refpoint = SCIPgetSolVal(scip, sol, x);
2639 
2640  /* adjust the reference point */
2641  refpoint = SCIPisLT(scip, refpoint, SCIPvarGetLbLocal(x)) ? SCIPvarGetLbLocal(x) : refpoint;
2642  refpoint = SCIPisGT(scip, refpoint, SCIPvarGetUbLocal(x)) ? SCIPvarGetUbLocal(x) : refpoint;
2643  assert(SCIPisLE(scip, refpoint, SCIPvarGetUbLocal(x)) && SCIPisGE(scip, refpoint, SCIPvarGetLbLocal(x)));
2644 
2645  if( SCIPisPositive(scip, coeff) )
2646  SCIPaddSquareLinearization(scip, coeff, refpoint, SCIPvarIsIntegral(x), &lincoef, &linconst, success);
2647  else
2648  SCIPaddSquareSecant(scip, coeff, SCIPvarGetLbLocal(x), SCIPvarGetUbLocal(x), &lincoef, &linconst, success);
2649 
2650  *cutactivity += lincoef * refpoint + linconst;
2651  *cutconstant += linconst;
2652 
2653  /* add underestimate to cut */
2654  SCIP_CALL( SCIPaddVarToRow(scip, cut, x, lincoef) );
2655 
2656  SCIPdebugMsg(scip, "add to cut: %e * %s + %e\n", lincoef, SCIPvarGetName(x), linconst);
2657  }
2658  /* bilinear case */
2659  else
2660  {
2661  SCIP_Real refpointx;
2662  SCIP_Real refpointy;
2663  SCIP_Real lincoefx;
2664  SCIP_Real lincoefy;
2665  SCIP_Real linconst;
2666 
2667  lincoefx = 0.0;
2668  lincoefy = 0.0;
2669  linconst = 0.0;
2670  refpointx = SCIPgetSolVal(scip, sol, x);
2671  refpointy = SCIPgetSolVal(scip, sol, y);
2672 
2673  /* adjust the reference points */
2674  refpointx = SCIPisLT(scip, refpointx, SCIPvarGetLbLocal(x)) ? SCIPvarGetLbLocal(x) : refpointx;
2675  refpointx = SCIPisGT(scip, refpointx, SCIPvarGetUbLocal(x)) ? SCIPvarGetUbLocal(x) : refpointx;
2676  refpointy = SCIPisLT(scip, refpointy, SCIPvarGetLbLocal(y)) ? SCIPvarGetLbLocal(y) : refpointy;
2677  refpointy = SCIPisGT(scip, refpointy, SCIPvarGetUbLocal(y)) ? SCIPvarGetUbLocal(y) : refpointy;
2678  assert(SCIPisLE(scip, refpointx, SCIPvarGetUbLocal(x)) && SCIPisGE(scip, refpointx, SCIPvarGetLbLocal(x)));
2679  assert(SCIPisLE(scip, refpointy, SCIPvarGetUbLocal(y)) && SCIPisGE(scip, refpointy, SCIPvarGetLbLocal(y)));
2680 
2681  SCIPaddBilinMcCormick(scip, coeff, SCIPvarGetLbLocal(x), SCIPvarGetUbLocal(x), refpointx, SCIPvarGetLbLocal(y),
2682  SCIPvarGetUbLocal(y), refpointy, FALSE, &lincoefx, &lincoefy, &linconst, success);
2683 
2684  *cutactivity += lincoefx * refpointx + lincoefy * refpointy + linconst;
2685  *cutconstant += linconst;
2686 
2687  /* add underestimate to cut */
2688  SCIP_CALL( SCIPaddVarToRow(scip, cut, x, lincoefx) );
2689  SCIP_CALL( SCIPaddVarToRow(scip, cut, y, lincoefy) );
2690 
2691  SCIPdebugMsg(scip, "add to cut: %e * %s + %e * %s + %e\n", lincoefx, SCIPvarGetName(x), lincoefy,
2692  SCIPvarGetName(y), linconst);
2693  }
2694 
2695  return SCIP_OKAY;
2696 }
2697 
2698 /** method to compute and add a cut for a nonlinear row aggregation and a given solution
2699  *
2700  * we compute for each edge concave aggregation one facet;
2701  * the remaining bilinear terms will be underestimated with McCormick, secants or linearizations;
2702  * constant and linear terms will be added to the cut directly
2703  */
2704 static
2705 SCIP_RETCODE computeCut(
2706  SCIP* scip, /**< SCIP data structure */
2707  SCIP_SEPA* sepa, /**< separator */
2708  SCIP_SEPADATA* sepadata, /**< separator data */
2709  SCIP_NLROWAGGR* nlrowaggr, /**< nonlinear row aggregation */
2710  SCIP_SOL* sol, /**< current solution (might be NULL) */
2711  SCIP_Bool* separated, /**< pointer to store if we could separate the current solution */
2712  SCIP_Bool* cutoff /**< pointer to store if the current node gets cut off */
2713  )
2714 {
2715  SCIP_ROW* cut;
2716  SCIP_Real* bestfacet;
2717  SCIP_Real bestfacetval;
2718  SCIP_Real cutconstant;
2719  SCIP_Real cutactivity;
2720  int bestfacetsize;
2721  char cutname[SCIP_MAXSTRLEN];
2722  SCIP_Bool found;
2723  SCIP_Bool islocalcut;
2724  int i;
2725 
2726  assert(separated != NULL);
2727  assert(cutoff != NULL);
2728  assert(nlrowaggr->necaggr > 0);
2729  assert(nlrowaggr->nlrow != NULL);
2730  assert(SCIPnlrowIsInNLP(nlrowaggr->nlrow));
2731 
2732  *separated = FALSE;
2733  *cutoff = FALSE;
2734  /* we use SCIPgetDepth because we add the cut to the global cut pool if cut is globally valid */
2735  islocalcut = SCIPgetDepth(scip) != 0;
2736 
2737  /* create the cut */
2738  (void) SCIPsnprintf(cutname, SCIP_MAXSTRLEN, "ec");
2739  SCIP_CALL( SCIPcreateEmptyRowSepa(scip, &cut, sepa, cutname, -SCIPinfinity(scip), SCIPinfinity(scip), islocalcut, FALSE,
2740  sepadata->dynamiccuts) );
2741  SCIP_CALL( SCIPcacheRowExtensions(scip, cut) );
2742 
2743  /* track rhs and activity of the cut */
2744  cutconstant = nlrowaggr->constant;
2745  cutactivity = 0.0;
2746 
2747  /* allocate necessary memory */
2748  bestfacetsize = sepadata->maxaggrsize + 1;
2749  SCIP_CALL( SCIPallocBufferArray(scip, &bestfacet, bestfacetsize) );
2750 
2751 #ifdef SCIP_DEBUG
2752  SCIP_CALL( SCIPprintNlRow(scip, nlrowaggr->nlrow, NULL) );
2753 
2754  SCIPdebugMsg(scip, "current solution:\n");
2755  for( i = 0; i < SCIPgetNVars(scip); ++i )
2756  {
2757  SCIP_VAR* var = SCIPgetVars(scip)[i];
2758  SCIPdebugMsg(scip, " %s = [%e, %e] solval = %e\n", SCIPvarGetName(var), SCIPvarGetLbLocal(var),
2759  SCIPvarGetUbLocal(var), SCIPgetSolVal(scip, sol, var));
2760  }
2761 #endif
2762 
2763  /* compute a facet for each edge-concave aggregation */
2764  for( i = 0; i < nlrowaggr->necaggr; ++i )
2765  {
2766  SCIP_ECAGGR* ecaggr;
2767  SCIP_Bool success;
2768 
2769  ecaggr = nlrowaggr->ecaggr[i];
2770  assert(ecaggr != NULL);
2771 
2772  /* compute a facet of the convex envelope */
2773  SCIP_CALL( computeConvexEnvelopeFacet(scip, sepadata, sol, ecaggr, bestfacet, &bestfacetval, &found) );
2774 
2775  SCIPdebugMsg(scip, "found facet for edge-concave aggregation %d/%d ? %s\n", i, nlrowaggr->necaggr,
2776  found ? "yes" : "no");
2777 
2778 #ifdef SCIP_DEBUG
2779  if( found )
2780  printFacet(scip, ecaggr->vars, ecaggr->nvars, bestfacet, bestfacetval);
2781 #endif
2782 
2783  /* do not add any cut because we did not found a facet for at least one edge-concave aggregation */
2784  if( !found ) /*lint !e774*/
2785  goto TERMINATE;
2786 
2787  /* add facet to the cut and update the rhs and activity of the cut */
2788  SCIP_CALL( addFacetToCut(scip, sol, cut, bestfacet, ecaggr->vars, ecaggr->nvars, &cutconstant, &cutactivity,
2789  &success) );
2790 
2791  if( !success )
2792  goto TERMINATE;
2793  }
2794 
2795  /* compute an underestimate for each bilinear term which is not in any edge-concave aggregation */
2796  for( i = 0; i < nlrowaggr->nremterms; ++i )
2797  {
2798  SCIP_VAR* x;
2799  SCIP_VAR* y;
2800  SCIP_Bool success;
2801 
2802  x = nlrowaggr->remtermvars1[i];
2803  y = nlrowaggr->remtermvars2[i];
2804  assert(x != NULL);
2805  assert(y != NULL);
2806 
2807  SCIP_CALL( addBilinearTermToCut(scip, sol, cut, x, y, nlrowaggr->remtermcoefs[i], &cutconstant, &cutactivity,
2808  &success) );
2809 
2810  if( !success )
2811  goto TERMINATE;
2812  }
2813 
2814  /* add all linear terms to the cut */
2815  for( i = 0; i < nlrowaggr->nlinvars; ++i )
2816  {
2817  SCIP_VAR* x;
2818  SCIP_Real coef;
2819  SCIP_Bool success;
2820 
2821  x = nlrowaggr->linvars[i];
2822  assert(x != NULL);
2823 
2824  coef = nlrowaggr->lincoefs[i];
2825 
2826  SCIP_CALL( addLinearTermToCut(scip, sol, cut, x, coef, &cutconstant, &cutactivity, &success) );
2827 
2828  if( !success )
2829  goto TERMINATE;
2830  }
2831 
2832  SCIPdebugMsg(scip, "cut activity = %e rhs(nlrow) = %e\n", cutactivity, nlrowaggr->rhs);
2833 
2834  /* set rhs of the cut (substract the constant part of the cut) */
2835  SCIP_CALL( SCIPchgRowRhs(scip, cut, nlrowaggr->rhs - cutconstant) );
2836  SCIP_CALL( SCIPflushRowExtensions(scip, cut) );
2837 
2838  /* check activity of the row; this assert can fail because of numerics */
2839  /* assert(SCIPisFeasEQ(scip, cutactivity - cutconstant, SCIPgetRowSolActivity(scip, cut, sol)) ); */
2840 
2841 #ifdef SCIP_DEBUG
2842  SCIP_CALL( SCIPprintRow(scip, cut, NULL) );
2843 #endif
2844 
2845  SCIPdebugMsg(scip, "EC cut <%s>: act=%f eff=%f rank=%d range=%e\n",
2846  SCIProwGetName(cut), SCIPgetRowSolActivity(scip, cut, sol), SCIPgetCutEfficacy(scip, sol, cut),
2847  SCIProwGetRank(cut), SCIPgetRowMaxCoef(scip, cut) / SCIPgetRowMinCoef(scip, cut) );
2848 
2849  /* try to add the cut has a finite rhs, is efficacious, and does not exceed the maximum cut range */
2850  if( !SCIPisInfinity(scip, nlrowaggr->rhs - cutconstant) && SCIPisCutEfficacious(scip, sol, cut)
2851  && SCIPgetRowMaxCoef(scip, cut) / SCIPgetRowMinCoef(scip, cut) < sepadata->cutmaxrange )
2852  {
2853  /* add the cut if it is separating the given solution by at least minviolation */
2854  if( SCIPisGE(scip, cutactivity - nlrowaggr->rhs, sepadata->minviolation) )
2855  {
2856  SCIP_CALL( SCIPaddRow(scip, cut, FALSE, cutoff) );
2857  *separated = TRUE;
2858  SCIPdebugMsg(scip, "added separating cut\n");
2859  }
2860 
2861  if( !(*cutoff) && !islocalcut )
2862  {
2863  SCIP_CALL( SCIPaddPoolCut(scip, cut) );
2864  SCIPdebugMsg(scip, "added cut to cut pool\n");
2865  }
2866  }
2867 
2868 TERMINATE:
2869  /* free allocated memory */
2870  SCIPfreeBufferArray(scip, &bestfacet);
2871 
2872  /* release the row */
2873  SCIP_CALL( SCIPreleaseRow(scip, &cut) );
2874 
2875  return SCIP_OKAY;
2876 }
2877 
2878 /** returns whether it is possible to compute a cut for a given nonlinear row aggregation */
2879 static
2880 SCIP_Bool isPossibleToComputeCut(
2881  SCIP* scip, /**< SCIP data structure */
2882  SCIP_SOL* sol, /**< current solution (might be NULL) */
2883  SCIP_NLROWAGGR* nlrowaggr /**< nonlinear row aggregation */
2884  )
2885 {
2886  int i;
2887 
2888  assert(scip != NULL);
2889  assert(nlrowaggr != NULL);
2890 
2891  if( !SCIPnlrowIsInNLP(nlrowaggr->nlrow) )
2892  {
2893  SCIPdebugMsg(scip, "nlrow is not in NLP anymore\n");
2894  return FALSE;
2895  }
2896 
2897  for( i = 0; i < nlrowaggr->nquadvars; ++i )
2898  {
2899  SCIP_VAR* var = nlrowaggr->quadvars[i];
2900  assert(var != NULL);
2901 
2902  /* check whether the variable has infinite bounds */
2903  if( SCIPisInfinity(scip, REALABS(SCIPvarGetLbLocal(var))) || SCIPisInfinity(scip, REALABS(SCIPvarGetUbLocal(var)))
2904  || SCIPisInfinity(scip, REALABS(SCIPgetSolVal(scip, sol, var))) )
2905  {
2906  SCIPdebugMsg(scip, "nlrow aggregation contains unbounded variables\n");
2907  return FALSE;
2908  }
2909 
2910  /* check whether the variable has been fixed and is in one edge-concave aggregation */
2911  if( nlrowaggr->quadvar2aggr[i] >= 0 && SCIPisFeasEQ(scip, SCIPvarGetLbLocal(var), SCIPvarGetUbLocal(var)) )
2912  {
2913  SCIPdebugMsg(scip, "nlrow aggregation contains fixed variables in an e.c. aggregation\n");
2914  return FALSE;
2915  }
2916  }
2917 
2918  return TRUE;
2919 }
2920 
2921 /** searches and tries to add edge-concave cuts */
2922 static
2923 SCIP_RETCODE separateCuts(
2924  SCIP* scip, /**< SCIP data structure */
2925  SCIP_SEPA* sepa, /**< separator */
2926  SCIP_SEPADATA* sepadata, /**< separator data */
2927  int depth, /**< current depth */
2928  SCIP_SOL* sol, /**< current solution */
2929  SCIP_RESULT* result /**< pointer to store the result of the separation call */
2930  )
2931 {
2932  int nmaxcuts;
2933  int ncuts;
2934  int i;
2935 
2936  assert(*result == SCIP_DIDNOTRUN);
2937 
2938  SCIPdebugMsg(scip, "separate cuts...\n");
2939 
2940  /* skip if there are no nonlinear row aggregations */
2941  if( sepadata->nnlrowaggrs == 0 )
2942  {
2943  SCIPdebugMsg(scip, "no aggregations exists -> skip call\n");
2944  return SCIP_OKAY;
2945  }
2946 
2947  /* get the maximal number of cuts allowed in a separation round */
2948  nmaxcuts = depth == 0 ? sepadata->maxsepacutsroot : sepadata->maxsepacuts;
2949  ncuts = 0;
2950 
2951  /* try to compute cuts for each nonlinear row independently */
2952  for( i = 0; i < sepadata->nnlrowaggrs && ncuts < nmaxcuts && !SCIPisStopped(scip); ++i )
2953  {
2954  SCIP_NLROWAGGR* nlrowaggr;
2955  SCIP_Bool separated;
2956  SCIP_Bool cutoff;
2957 
2958  nlrowaggr = sepadata->nlrowaggrs[i];
2959  assert(nlrowaggr != NULL);
2960 
2961  /* skip nonlinear aggregations for which it is obviously not possible to compute a cut */
2962  if( !isPossibleToComputeCut(scip, sol, nlrowaggr) )
2963  return SCIP_OKAY;
2964 
2965  *result = (*result == SCIP_DIDNOTRUN) ? SCIP_DIDNOTFIND : *result;
2966 
2967  SCIPdebugMsg(scip, "try to compute a cut for nonlinear row aggregation %d\n", i);
2968 
2969  /* compute and add cut */
2970  SCIP_CALL( computeCut(scip, sepa, sepadata, nlrowaggr, sol, &separated, &cutoff) );
2971  SCIPdebugMsg(scip, "found a cut: %s cutoff: %s\n", separated ? "yes" : "no", cutoff ? "yes" : "no");
2972 
2973  /* stop if the current node gets cut off */
2974  if( cutoff )
2975  {
2976  assert(separated);
2977  *result = SCIP_CUTOFF;
2978  return SCIP_OKAY;
2979  }
2980 
2981  /* do not compute more cuts if we already separated the given solution */
2982  if( separated )
2983  {
2984  assert(!cutoff);
2985  *result = SCIP_SEPARATED;
2986  ++ncuts;
2987  }
2988  }
2989 
2990  return SCIP_OKAY;
2991 }
2992 
2993 /*
2994  * Callback methods of separator
2995  */
2996 
2997 /** copy method for separator plugins (called when SCIP copies plugins) */
2998 static
2999 SCIP_DECL_SEPACOPY(sepaCopyEccuts)
3000 { /*lint --e{715}*/
3001  assert(scip != NULL);
3002  assert(sepa != NULL);
3003  assert(strcmp(SCIPsepaGetName(sepa), SEPA_NAME) == 0);
3004 
3005  /* call inclusion method of constraint handler */
3006  SCIP_CALL( SCIPincludeSepaEccuts(scip) );
3007 
3008  return SCIP_OKAY;
3009 }
3010 
3011 /** destructor of separator to free user data (called when SCIP is exiting) */
3012 static
3013 SCIP_DECL_SEPAFREE(sepaFreeEccuts)
3014 { /*lint --e{715}*/
3015  SCIP_SEPADATA* sepadata;
3016 
3017  sepadata = SCIPsepaGetData(sepa);
3018  assert(sepadata != NULL);
3019 
3020  SCIP_CALL( sepadataFree(scip, &sepadata) );
3021  SCIPsepaSetData(sepa, NULL);
3022 
3023  return SCIP_OKAY;
3024 }
3025 
3026 /** solving process deinitialization method of separator (called before branch and bound process data is freed) */
3027 static
3028 SCIP_DECL_SEPAEXITSOL(sepaExitsolEccuts)
3029 { /*lint --e{715}*/
3030  SCIP_SEPADATA* sepadata;
3031 
3032  sepadata = SCIPsepaGetData(sepa);
3033  assert(sepadata != NULL);
3034 
3035  /* print statistics */
3036 #ifdef SCIP_STATISTIC
3037  SCIPstatisticMessage("rhs-AGGR %d\n", sepadata->nrhsnlrowaggrs);
3038  SCIPstatisticMessage("lhs-AGGR %d\n", sepadata->nlhsnlrowaggrs);
3039  SCIPstatisticMessage("aggr. search time = %f\n", sepadata->aggrsearchtime);
3040 #endif
3041 
3042  /* free nonlinear row aggregations */
3043  SCIP_CALL( sepadataFreeNlrows(scip, sepadata) );
3044 
3045  /* mark that we should search again for nonlinear row aggregations */
3046  sepadata->searchedforaggr = FALSE;
3047 
3048  SCIPdebugMsg(scip, "exitsol\n");
3049 
3050  return SCIP_OKAY;
3051 }
3052 
3053 /** LP solution separation method of separator */
3054 static
3055 SCIP_DECL_SEPAEXECLP(sepaExeclpEccuts)
3056 { /*lint --e{715}*/
3057  SCIP_SEPADATA* sepadata;
3058  int ncalls;
3059 
3060  sepadata = SCIPsepaGetData(sepa);
3061  assert(sepadata != NULL);
3062 
3063  *result = SCIP_DIDNOTRUN;
3064 
3065  if( !allowlocal )
3066  return SCIP_OKAY;
3067 
3068  /* check min- and maximal aggregation size */
3069  if( sepadata->maxaggrsize < sepadata->minaggrsize )
3070  return SCIP_PARAMETERWRONGVAL;
3071 
3072  /* only call separator, if we are not close to terminating */
3073  if( SCIPisStopped(scip) )
3074  return SCIP_OKAY;
3075 
3076  /* skip if the LP is not constructed yet */
3077  if( !SCIPisNLPConstructed(scip) )
3078  {
3079  SCIPdebugMsg(scip, "Skip since NLP is not constructed yet.\n");
3080  return SCIP_OKAY;
3081  }
3082 
3083  /* only call separator up to a maximum depth */
3084  if ( sepadata->maxdepth >= 0 && depth > sepadata->maxdepth )
3085  return SCIP_OKAY;
3086 
3087  /* only call separator a given number of times at each node */
3088  ncalls = SCIPsepaGetNCallsAtNode(sepa);
3089  if ( (depth == 0 && sepadata->maxroundsroot >= 0 && ncalls >= sepadata->maxroundsroot)
3090  || (depth > 0 && sepadata->maxrounds >= 0 && ncalls >= sepadata->maxrounds) )
3091  return SCIP_OKAY;
3092 
3093  /* search for nonlinear row aggregations */
3094  if( !sepadata->searchedforaggr )
3095  {
3096  int i;
3097 
3098  SCIPstatistic( sepadata->aggrsearchtime -= SCIPgetTotalTime(scip) );
3099 
3100  SCIPdebugMsg(scip, "search for nonlinear row aggregations\n");
3101  for( i = 0; i < SCIPgetNNLPNlRows(scip) && !SCIPisStopped(scip); ++i )
3102  {
3103  SCIP_NLROW* nlrow = SCIPgetNLPNlRows(scip)[i];
3104  SCIP_CALL( findAndStoreEcAggregations(scip, sepadata, nlrow, NULL) );
3105  }
3106  sepadata->searchedforaggr = TRUE;
3107 
3108  SCIPstatistic( sepadata->aggrsearchtime += SCIPgetTotalTime(scip) );
3109  }
3110 
3111  /* search for edge-concave cuts */
3112  SCIP_CALL( separateCuts(scip, sepa, sepadata, depth, NULL, result) );
3113 
3114  return SCIP_OKAY;
3115 }
3116 
3117 /*
3118  * separator specific interface methods
3119  */
3120 
3121 /** creates the edge-concave separator and includes it in SCIP
3122  *
3123  * @ingroup SeparatorIncludes
3124  */
3125 SCIP_RETCODE SCIPincludeSepaEccuts(
3126  SCIP* scip /**< SCIP data structure */
3127  )
3128 {
3129  SCIP_SEPADATA* sepadata;
3130  SCIP_SEPA* sepa;
3131 
3132  /* create eccuts separator data */
3133  SCIP_CALL( sepadataCreate(scip, &sepadata) );
3134 
3135  /* include separator */
3136  SCIP_CALL( SCIPincludeSepaBasic(scip, &sepa, SEPA_NAME, SEPA_DESC, SEPA_PRIORITY, SEPA_FREQ, SEPA_MAXBOUNDDIST,
3137  SEPA_USESSUBSCIP, SEPA_DELAY, sepaExeclpEccuts, NULL, sepadata) );
3138 
3139  assert(sepa != NULL);
3140 
3141  /* set non fundamental callbacks via setter functions */
3142  SCIP_CALL( SCIPsetSepaCopy(scip, sepa, sepaCopyEccuts) );
3143  SCIP_CALL( SCIPsetSepaFree(scip, sepa, sepaFreeEccuts) );
3144  SCIP_CALL( SCIPsetSepaExitsol(scip, sepa, sepaExitsolEccuts) );
3145 
3146  /* add eccuts separator parameters */
3147  SCIP_CALL( SCIPaddBoolParam(scip,
3148  "separating/" SEPA_NAME "/dynamiccuts",
3149  "should generated cuts be removed from the LP if they are no longer tight?",
3150  &sepadata->dynamiccuts, FALSE, DEFAULT_DYNAMICCUTS, NULL, NULL) );
3151 
3152  SCIP_CALL( SCIPaddIntParam(scip,
3153  "separating/" SEPA_NAME "/maxrounds",
3154  "maximal number of eccuts separation rounds per node (-1: unlimited)",
3155  &sepadata->maxrounds, FALSE, DEFAULT_MAXROUNDS, -1, INT_MAX, NULL, NULL) );
3156 
3157  SCIP_CALL( SCIPaddIntParam(scip,
3158  "separating/" SEPA_NAME "/maxroundsroot",
3159  "maximal number of eccuts separation rounds in the root node (-1: unlimited)",
3160  &sepadata->maxroundsroot, FALSE, DEFAULT_MAXROUNDSROOT, -1, INT_MAX, NULL, NULL) );
3161 
3162  SCIP_CALL( SCIPaddIntParam(scip,
3163  "separating/" SEPA_NAME "/maxdepth",
3164  "maximal depth at which the separator is applied (-1: unlimited)",
3165  &sepadata->maxdepth, FALSE, DEFAULT_MAXDEPTH, -1, INT_MAX, NULL, NULL) );
3166 
3167  SCIP_CALL( SCIPaddIntParam(scip,
3168  "separating/" SEPA_NAME "/maxsepacuts",
3169  "maximal number of edge-concave cuts separated per separation round",
3170  &sepadata->maxsepacuts, FALSE, DEFAULT_MAXSEPACUTS, 0, INT_MAX, NULL, NULL) );
3171 
3172  SCIP_CALL( SCIPaddIntParam(scip,
3173  "separating/" SEPA_NAME "/maxsepacutsroot",
3174  "maximal number of edge-concave cuts separated per separation round in the root node",
3175  &sepadata->maxsepacutsroot, FALSE, DEFAULT_MAXSEPACUTSROOT, 0, INT_MAX, NULL, NULL) );
3176 
3177  SCIP_CALL( SCIPaddRealParam(scip, "separating/" SEPA_NAME "/cutmaxrange",
3178  "maximal coef. range of a cut (max coef. divided by min coef.) in order to be added to LP relaxation",
3179  &sepadata->cutmaxrange, FALSE, DEFAULT_CUTMAXRANGE, 0.0, SCIPinfinity(scip), NULL, NULL) );
3180 
3181  SCIP_CALL( SCIPaddRealParam(scip, "separating/" SEPA_NAME "/minviolation",
3182  "minimal violation of an edge-concave cut to be separated",
3183  &sepadata->minviolation, FALSE, DEFAULT_MINVIOLATION, 0.0, 0.5, NULL, NULL) );
3184 
3185  SCIP_CALL( SCIPaddIntParam(scip,
3186  "separating/" SEPA_NAME "/minaggrsize",
3187  "search for edge-concave aggregations of at least this size",
3188  &sepadata->minaggrsize, TRUE, DEFAULT_MINAGGRSIZE, 3, 5, NULL, NULL) );
3189 
3190  SCIP_CALL( SCIPaddIntParam(scip,
3191  "separating/" SEPA_NAME "/maxaggrsize",
3192  "search for edge-concave aggregations of at most this size",
3193  &sepadata->maxaggrsize, TRUE, DEFAULT_MAXAGGRSIZE, 3, 5, NULL, NULL) );
3194 
3195  SCIP_CALL( SCIPaddIntParam(scip,
3196  "separating/" SEPA_NAME "/maxbilinterms",
3197  "maximum number of bilinear terms allowed to be in a quadratic constraint",
3198  &sepadata->maxbilinterms, TRUE, DEFAULT_MAXBILINTERMS, 0, INT_MAX, NULL, NULL) );
3199 
3200  SCIP_CALL( SCIPaddIntParam(scip,
3201  "separating/" SEPA_NAME "/maxstallrounds",
3202  "maximum number of unsuccessful rounds in the edge-concave aggregation search",
3203  &sepadata->maxstallrounds, TRUE, DEFAULT_MAXSTALLROUNDS, 0, INT_MAX, NULL, NULL) );
3204 
3205  return SCIP_OKAY;
3206 }
SCIP_RESULT
enum SCIP_Result SCIP_RESULT
Definition: type_result.h:52
SCIPfreeBlockMemoryArray
#define SCIPfreeBlockMemoryArray(scip, ptr, num)
Definition: scip_mem.h:101
USEDUALSIMPLEX
#define USEDUALSIMPLEX
Definition: sepa_eccuts.c:69
NlrowAggr::remtermcoefs
SCIP_Real * remtermcoefs
Definition: sepa_eccuts.c:118
SCIPexprGetQuadraticData
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:4057
SCIPreallocBlockMemoryArray
#define SCIPreallocBlockMemoryArray(scip, ptr, oldnum, newnum)
Definition: scip_mem.h:90
NlrowAggr::quadvar2aggr
int * quadvar2aggr
Definition: sepa_eccuts.c:111
NlrowAggr::rhsaggr
SCIP_Bool rhsaggr
Definition: sepa_eccuts.c:99
SCIPlpiGetNRows
SCIP_RETCODE SCIPlpiGetNRows(SCIP_LPI *lpi, int *nrows)
Definition: lpi_clp.cpp:1408
SCIP_DIDNOTRUN
Definition: type_result.h:33
SCIPgetNNLPNlRows
int SCIPgetNNLPNlRows(SCIP *scip)
Definition: scip_nlp.c:332
SCIPgetSolvingTime
SCIP_Real SCIPgetSolvingTime(SCIP *scip)
Definition: scip_timing.c:369
TCLIQUE_STATUS
enum TCLIQUE_Status TCLIQUE_STATUS
Definition: tclique.h:59
SCIPlpiFree
SCIP_RETCODE SCIPlpiFree(SCIP_LPI **lpi)
Definition: lpi_clp.cpp:634
narcs
#define narcs
Definition: gastrans.c:68
SCIPallocBlockMemoryArray
#define SCIPallocBlockMemoryArray(scip, ptr, num)
Definition: scip_mem.h:84
SCIPisNLPConstructed
SCIP_Bool SCIPisNLPConstructed(SCIP *scip)
Definition: scip_nlp.c:101
tcliqueFree
void tcliqueFree(TCLIQUE_GRAPH **tcliquegraph)
Definition: tclique_graph.c:202
SCIPcacheRowExtensions
SCIP_RETCODE SCIPcacheRowExtensions(SCIP *scip, SCIP_ROW *row)
Definition: scip_lp.c:1626
addFacetToCut
static SCIP_RETCODE addFacetToCut(SCIP *scip, SCIP_SOL *sol, SCIP_ROW *cut, SCIP_Real *facet, SCIP_VAR **vars, int nvars, SCIP_Real *cutconstant, SCIP_Real *cutactivity, SCIP_Bool *success)
Definition: sepa_eccuts.c:2484
SCIPisFeasEQ
SCIP_Bool SCIPisFeasEQ(SCIP *scip, SCIP_Real val1, SCIP_Real val2)
Definition: scip_numerics.c:773
SCIPcreateConsBasicLinear
SCIP_RETCODE SCIPcreateConsBasicLinear(SCIP *scip, SCIP_CONS **cons, const char *name, int nvars, SCIP_VAR **vars, SCIP_Real *vals, SCIP_Real lhs, SCIP_Real rhs)
Definition: cons_linear.c:18010
tcliqueCreate
TCLIQUE_Bool tcliqueCreate(TCLIQUE_GRAPH **tcliquegraph)
Definition: tclique_graph.c:176
NlrowAggr::nremterms
int nremterms
Definition: sepa_eccuts.c:119
SCIP_OBJSENSE_MAXIMIZE
Definition: type_prob.h:38
TCLIQUE_GRAPH
struct TCLIQUE_Graph TCLIQUE_GRAPH
Definition: tclique.h:40
NlrowAggr::linvarssize
int linvarssize
Definition: sepa_eccuts.c:108
SCIPflushRowExtensions
SCIP_RETCODE SCIPflushRowExtensions(SCIP *scip, SCIP_ROW *row)
Definition: scip_lp.c:1649
SCIPallocClearBufferArray
#define SCIPallocClearBufferArray(scip, ptr, num)
Definition: scip_mem.h:117
SCIPlpiGetSol
SCIP_RETCODE SCIPlpiGetSol(SCIP_LPI *lpi, SCIP_Real *objval, SCIP_Real *primsol, SCIP_Real *dualsol, SCIP_Real *activity, SCIP_Real *redcost)
Definition: lpi_clp.cpp:2774
NlrowAggr::rhs
SCIP_Real rhs
Definition: sepa_eccuts.c:122
tcliqueAddNode
TCLIQUE_Bool tcliqueAddNode(TCLIQUE_GRAPH *tcliquegraph, int node, TCLIQUE_WEIGHT weight)
Definition: tclique_graph.c:331
SCIPvarGetLbGlobal
SCIP_Real SCIPvarGetLbGlobal(SCIP_VAR *var)
Definition: var.c:17910
SEPA_DELAY
#define SEPA_DELAY
Definition: sepa_eccuts.c:42
SEPA_DESC
#define SEPA_DESC
Definition: sepa_eccuts.c:37
NlrowAggr::remtermvars1
SCIP_VAR ** remtermvars1
Definition: sepa_eccuts.c:116
SCIPgetRealParam
SCIP_RETCODE SCIPgetRealParam(SCIP *scip, const char *name, SCIP_Real *value)
Definition: scip_param.c:298
SCIP_MAXSTRLEN
#define SCIP_MAXSTRLEN
Definition: def.h:293
SCIPlpiSolvePrimal
SCIP_RETCODE SCIPlpiSolvePrimal(SCIP_LPI *lpi)
Definition: lpi_clp.cpp:1796
SCIP_DECL_SEPACOPY
static SCIP_DECL_SEPACOPY(sepaCopyEccuts)
Definition: sepa_eccuts.c:3001
SCIPcalcMemGrowSize
int SCIPcalcMemGrowSize(SCIP *scip, int num)
Definition: scip_mem.c:130
SCIPaddVarToRow
SCIP_RETCODE SCIPaddVarToRow(SCIP *scip, SCIP_ROW *row, SCIP_VAR *var, SCIP_Real val)
Definition: scip_lp.c:1686
CLIQUE_MINWEIGHT
#define CLIQUE_MINWEIGHT
Definition: sepa_eccuts.c:47
SCIPisPositive
SCIP_Bool SCIPisPositive(SCIP *scip, SCIP_Real val)
Definition: scip_numerics.c:563
SCIPvarGetLbLocal
SCIP_Real SCIPvarGetLbLocal(SCIP_VAR *var)
Definition: var.c:17966
SCIP_CALL_FINALLY
#define SCIP_CALL_FINALLY(x, y)
Definition: def.h:426
SCIP_DECL_SEPAEXITSOL
static SCIP_DECL_SEPAEXITSOL(sepaExitsolEccuts)
Definition: sepa_eccuts.c:3030
SCIPisGE
SCIP_Bool SCIPisGE(SCIP *scip, SCIP_Real val1, SCIP_Real val2)
Definition: scip_numerics.c:490
sepadataFreeNlrows
static SCIP_RETCODE sepadataFreeNlrows(SCIP *scip, SCIP_SEPADATA *sepadata)
Definition: sepa_eccuts.c:737
SEPA_PRIORITY
#define SEPA_PRIORITY
Definition: sepa_eccuts.c:38
EcAggr::termcoefs
SCIP_Real * termcoefs
Definition: sepa_eccuts.c:87
SCIPlpiChgSides
SCIP_RETCODE SCIPlpiChgSides(SCIP_LPI *lpi, int nrows, const int *ind, const SCIP_Real *lhs, const SCIP_Real *rhs)
Definition: lpi_clp.cpp:1158
CLIQUE_MAXNTREENODES
#define CLIQUE_MAXNTREENODES
Definition: sepa_eccuts.c:48
SCIProwGetName
const char * SCIProwGetName(SCIP_ROW *row)
Definition: lp.c:17284
SCIP_Var
Definition: struct_var.h:198
searchEcAggrWithCliques
static SCIP_RETCODE searchEcAggrWithCliques(SCIP *scip, TCLIQUE_GRAPH *graph, SCIP_SEPADATA *sepadata, SCIP_NLROW *nlrow, int *quadvar2aggr, int nfoundsofar, SCIP_Bool rhsaggr, SCIP_Bool *foundaggr, SCIP_Bool *foundclique)
Definition: sepa_eccuts.c:1396
SCIPnlrowGetLinearVars
SCIP_VAR ** SCIPnlrowGetLinearVars(SCIP_NLROW *nlrow)
Definition: nlp.c:1774
SCIPreleaseVar
SCIP_RETCODE SCIPreleaseVar(SCIP *scip, SCIP_VAR **var)
Definition: scip_var.c:1245
EcAggr::varsize
int varsize
Definition: sepa_eccuts.c:85
searchEcAggrWithMIP
static SCIP_RETCODE searchEcAggrWithMIP(SCIP *subscip, SCIP_Real timelimit, int nedges, SCIP_Bool *aggrleft, SCIP_Bool *found)
Definition: sepa_eccuts.c:1239
SCIPsetHeuristics
SCIP_RETCODE SCIPsetHeuristics(SCIP *scip, SCIP_PARAMSETTING paramsetting, SCIP_Bool quiet)
Definition: scip_param.c:897
SCIPcheckExprQuadratic
SCIP_RETCODE SCIPcheckExprQuadratic(SCIP *scip, SCIP_EXPR *expr, SCIP_Bool *isquadratic)
Definition: scip_expr.c:2340
FALSE
#define FALSE
Definition: def.h:87
SCIPhashmapCreate
SCIP_RETCODE SCIPhashmapCreate(SCIP_HASHMAP **hashmap, BMS_BLKMEM *blkmem, int mapsize)
Definition: misc.c:3014
NlrowAggr::nlrow
SCIP_NLROW * nlrow
Definition: sepa_eccuts.c:98
CLIQUE_MAXFIRSTNODEWEIGHT
#define CLIQUE_MAXFIRSTNODEWEIGHT
Definition: sepa_eccuts.c:44
NlrowAggr::ecaggr
SCIP_ECAGGR ** ecaggr
Definition: sepa_eccuts.c:102
SCIPinfinity
SCIP_Real SCIPinfinity(SCIP *scip)
Definition: scip_numerics.c:503
SCIPsnprintf
int SCIPsnprintf(char *t, int len, const char *s,...)
Definition: misc.c:10755
SCIPisNegative
SCIP_Bool SCIPisNegative(SCIP *scip, SCIP_Real val)
Definition: scip_numerics.c:575
createTcliqueGraph
static SCIP_RETCODE createTcliqueGraph(SCIP_NLROW *nlrow, TCLIQUE_GRAPH **graph, SCIP_Real *nodeweights)
Definition: sepa_eccuts.c:1305
TRUE
#define TRUE
Definition: def.h:86
SCIPdebug
#define SCIPdebug(x)
Definition: pub_message.h:84
SCIPsepaGetName
const char * SCIPsepaGetName(SCIP_SEPA *sepa)
Definition: sepa.c:720
computeConvexEnvelopeFacet
static SCIP_RETCODE computeConvexEnvelopeFacet(SCIP *scip, SCIP_SEPADATA *sepadata, SCIP_SOL *sol, SCIP_ECAGGR *ecaggr, SCIP_Real *facet, SCIP_Real *facetval, SCIP_Bool *success)
Definition: sepa_eccuts.c:2276
SCIP_RETCODE
enum SCIP_Retcode SCIP_RETCODE
Definition: type_retcode.h:54
SCIPstatisticMessage
#define SCIPstatisticMessage
Definition: pub_message.h:114
SCIP_VARTYPE_BINARY
Definition: type_var.h:53
SCIPhashmapInsertInt
SCIP_RETCODE SCIPhashmapInsertInt(SCIP_HASHMAP *hashmap, void *origin, int image)
Definition: misc.c:3132
SCIPnlrowGetRhs
SCIP_Real SCIPnlrowGetRhs(SCIP_NLROW *nlrow)
Definition: nlp.c:1814
SCIPcreateVarBasic
SCIP_RETCODE SCIPcreateVarBasic(SCIP *scip, SCIP_VAR **var, const char *name, SCIP_Real lb, SCIP_Real ub, SCIP_Real obj, SCIP_VARTYPE vartype)
Definition: scip_var.c:185
sepadataFree
static SCIP_RETCODE sepadataFree(SCIP *scip, SCIP_SEPADATA **sepadata)
Definition: sepa_eccuts.c:767
NlrowAggr::remtermvars2
SCIP_VAR ** remtermvars2
Definition: sepa_eccuts.c:117
SCIPnlrowGetNLinearVars
int SCIPnlrowGetNLinearVars(SCIP_NLROW *nlrow)
Definition: nlp.c:1764
SCIPlpiGetNCols
SCIP_RETCODE SCIPlpiGetNCols(SCIP_LPI *lpi, int *ncols)
Definition: lpi_clp.cpp:1426
sepa_eccuts.h
edge concave cut separator
SCIPfreeBlockMemory
#define SCIPfreeBlockMemory(scip, ptr)
Definition: scip_mem.h:99
NlrowAggr::quadvars
SCIP_VAR ** quadvars
Definition: sepa_eccuts.c:110
SCIPdebugMessage
#define SCIPdebugMessage
Definition: pub_message.h:87
SCIP_SEPARATED
Definition: type_result.h:40
SCIPgetMessagehdlr
SCIP_MESSAGEHDLR * SCIPgetMessagehdlr(SCIP *scip)
Definition: scip_message.c:79
tclique.h
tclique user interface
SCIPisEQ
SCIP_Bool SCIPisEQ(SCIP *scip, SCIP_Real val1, SCIP_Real val2)
Definition: scip_numerics.c:438
SCIP_Sepa
Definition: struct_sepa.h:37
SCIPfreeBufferArray
#define SCIPfreeBufferArray(scip, ptr)
Definition: scip_mem.h:127
SCIPcreate
SCIP_RETCODE SCIPcreate(SCIP **scip)
Definition: scip_general.c:283
SCIPsetSepaCopy
SCIP_RETCODE SCIPsetSepaCopy(SCIP *scip, SCIP_SEPA *sepa, SCIP_DECL_SEPACOPY((*sepacopy)))
Definition: scip_sepa.c:142
SCIPsetRealParam
SCIP_RETCODE SCIPsetRealParam(SCIP *scip, const char *name, SCIP_Real value)
Definition: scip_param.c:594
SCIPexprAreQuadraticExprsVariables
SCIP_Bool SCIPexprAreQuadraticExprsVariables(SCIP_EXPR *expr)
Definition: expr.c:4185
SCIPdebugMsgPrint
#define SCIPdebugMsgPrint
Definition: scip_message.h:70
SCIPdebugMsg
#define SCIPdebugMsg
Definition: scip_message.h:69
SCIPaddIntParam
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:74
nlrowaggrAddRemBilinTerm
static SCIP_RETCODE nlrowaggrAddRemBilinTerm(SCIP_NLROWAGGR *nlrowaggr, SCIP_VAR *x, SCIP_VAR *y, SCIP_Real coef)
Definition: sepa_eccuts.c:398
x
SCIP_VAR ** x
Definition: circlepacking.c:54
SCIPlpiCreate
SCIP_RETCODE SCIPlpiCreate(SCIP_LPI **lpi, SCIP_MESSAGEHDLR *messagehdlr, const char *name, SCIP_OBJSEN objsen)
Definition: lpi_clp.cpp:522
SCIPaddCoefLinear
SCIP_RETCODE SCIPaddCoefLinear(SCIP *scip, SCIP_CONS *cons, SCIP_VAR *var, SCIP_Real val)
Definition: cons_linear.c:18162
SCIPcreateProbBasic
SCIP_RETCODE SCIPcreateProbBasic(SCIP *scip, const char *name)
Definition: scip_prob.c:171
SCIPgetRowMaxCoef
SCIP_Real SCIPgetRowMaxCoef(SCIP *scip, SCIP_ROW *row)
Definition: scip_lp.c:1907
SCIPlpiAddCols
SCIP_RETCODE SCIPlpiAddCols(SCIP_LPI *lpi, int ncols, const SCIP_Real *obj, const SCIP_Real *lb, const SCIP_Real *ub, char **colnames, int nnonz, const int *beg, const int *ind, const SCIP_Real *val)
Definition: lpi_clp.cpp:749
SCIPsepaGetData
SCIP_SEPADATA * SCIPsepaGetData(SCIP_SEPA *sepa)
Definition: sepa.c:610
ecaggrCreateEmpty
static SCIP_RETCODE ecaggrCreateEmpty(SCIP *scip, SCIP_ECAGGR **ecaggr, int nquadvars, int nquadterms)
Definition: sepa_eccuts.c:168
SCIP_Sol
Definition: struct_sol.h:64
DEFAULT_MAXSTALLROUNDS
#define DEFAULT_MAXSTALLROUNDS
Definition: sepa_eccuts.c:64
SCIPhashmapExists
SCIP_Bool SCIPhashmapExists(SCIP_HASHMAP *hashmap, void *origin)
Definition: misc.c:3363
SCIP_NlRow::rhs
SCIP_Real rhs
Definition: struct_nlp.h:59
SCIPlpiSolveDual
SCIP_RETCODE SCIPlpiSolveDual(SCIP_LPI *lpi)
Definition: lpi_clp.cpp:1871
ecaggrFree
static SCIP_RETCODE ecaggrFree(SCIP *scip, SCIP_ECAGGR **ecaggr)
Definition: sepa_eccuts.c:198
EcAggr::termvars1
int * termvars1
Definition: sepa_eccuts.c:88
SCIPvarGetUbGlobal
SCIP_Real SCIPvarGetUbGlobal(SCIP_VAR *var)
Definition: var.c:17920
isPossibleToComputeCut
static SCIP_Bool isPossibleToComputeCut(SCIP *scip, SCIP_SOL *sol, SCIP_NLROWAGGR *nlrowaggr)
Definition: sepa_eccuts.c:2882
SEPA_NAME
#define SEPA_NAME
Definition: sepa_eccuts.c:36
SCIPcreateConsBasicXor
SCIP_RETCODE SCIPcreateConsBasicXor(SCIP *scip, SCIP_CONS **cons, const char *name, SCIP_Bool rhs, int nvars, SCIP_VAR **vars)
Definition: cons_xor.c:5900
SCIPduplicateBlockMemoryArray
#define SCIPduplicateBlockMemoryArray(scip, ptr, source, num)
Definition: scip_mem.h:96
SCIP_HashMap
Definition: struct_misc.h:128
SCIPisCutEfficacious
SCIP_Bool SCIPisCutEfficacious(SCIP *scip, SCIP_SOL *sol, SCIP_ROW *cut)
Definition: scip_cut.c:108
SCIPsetObjsense
SCIP_RETCODE SCIPsetObjsense(SCIP *scip, SCIP_OBJSENSE objsense)
Definition: scip_prob.c:1241
SCIPlpiAddRows
SCIP_RETCODE SCIPlpiAddRows(SCIP_LPI *lpi, int nrows, const SCIP_Real *lhs, const SCIP_Real *rhs, char **rownames, int nnonz, const int *beg, const int *ind, const SCIP_Real *val)
Definition: lpi_clp.cpp:905
SCIPnlrowIsInNLP
SCIP_Bool SCIPnlrowIsInNLP(SCIP_NLROW *nlrow)
Definition: nlp.c:1863
SCIPsolve
SCIP_RETCODE SCIPsolve(SCIP *scip)
Definition: scip_solve.c:2613
EcAggr::nterms
int nterms
Definition: sepa_eccuts.c:90
tcliqueAddEdge
TCLIQUE_Bool tcliqueAddEdge(TCLIQUE_GRAPH *tcliquegraph, int node1, int node2)
Definition: tclique_graph.c:371
doSeachEcAggr
static SCIP_RETCODE doSeachEcAggr(SCIP *scip, SCIP *subscip, SCIP_SEPADATA *sepadata, SCIP_NLROW *nlrow, SCIP_SOL *sol, SCIP_Bool rhsaggr, int *quadvar2aggr, int *nfound)
Definition: sepa_eccuts.c:1545
SCIPgetRowMinCoef
SCIP_Real SCIPgetRowMinCoef(SCIP *scip, SCIP_ROW *row)
Definition: scip_lp.c:1889
SCIP_DIDNOTFIND
Definition: type_result.h:35
SCIP_Cons
Definition: struct_cons.h:37
SCIPerrorMessage
#define SCIPerrorMessage
Definition: pub_message.h:55
SCIPaddCons
SCIP_RETCODE SCIPaddCons(SCIP *scip, SCIP_CONS *cons)
Definition: scip_prob.c:2769
EcAggr::vars
SCIP_VAR ** vars
Definition: sepa_eccuts.c:83
SCIPisLT
SCIP_Bool SCIPisLT(SCIP *scip, SCIP_Real val1, SCIP_Real val2)
Definition: scip_numerics.c:451
tcliqueChangeWeight
void tcliqueChangeWeight(TCLIQUE_GRAPH *tcliquegraph, int node, TCLIQUE_WEIGHT weight)
Definition: tclique_graph.c:353
SCIPgetVarExprVar
SCIP_VAR * SCIPgetVarExprVar(SCIP_EXPR *expr)
Definition: expr_var.c:407
SCIP_PARAMETERWRONGVAL
Definition: type_retcode.h:48
SCIP_OBJSEN_MINIMIZE
Definition: type_lpi.h:34
SEPA_FREQ
#define SEPA_FREQ
Definition: sepa_eccuts.c:39
SCIPsepaGetNCallsAtNode
int SCIPsepaGetNCallsAtNode(SCIP_SEPA *sepa)
Definition: sepa.c:847
SCIPgetStatus
SCIP_STATUS SCIPgetStatus(SCIP *scip)
Definition: scip_general.c:474
SCIPblkmem
BMS_BLKMEM * SCIPblkmem(SCIP *scip)
Definition: scip_mem.c:48
SCIPchgVarUb
SCIP_RETCODE SCIPchgVarUb(SCIP *scip, SCIP_VAR *var, SCIP_Real newbound)
Definition: scip_var.c:4763
phi
static SCIP_Real phi(SCIP *scip, SCIP_Real val, SCIP_Real lb, SCIP_Real ub)
Definition: sepa_eccuts.c:837
SCIPvarGetName
const char * SCIPvarGetName(SCIP_VAR *var)
Definition: var.c:17251
sepadataAddNlrowaggr
static SCIP_RETCODE sepadataAddNlrowaggr(SCIP *scip, SCIP_SEPADATA *sepadata, SCIP_NLROWAGGR *nlrowaggr)
Definition: sepa_eccuts.c:793
storeAggrFromMIP
static SCIP_RETCODE storeAggrFromMIP(SCIP *subscip, SCIP_NLROW *nlrow, SCIP_VAR **forwardarcs, SCIP_VAR **backwardarcs, int *quadvar2aggr, int nfoundsofar)
Definition: sepa_eccuts.c:1155
nlp.h
internal methods for NLP management
SCIPhashmapFree
void SCIPhashmapFree(SCIP_HASHMAP **hashmap)
Definition: misc.c:3048
SCIPsepaSetData
void SCIPsepaSetData(SCIP_SEPA *sepa, SCIP_SEPADATA *sepadata)
Definition: sepa.c:620
SCIPexprGetQuadraticQuadTerm
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:4104
NULL
#define NULL
Definition: lpi_spx1.cpp:155
REALABS
#define REALABS(x)
Definition: def.h:201
NlrowAggr::constant
SCIP_Real constant
Definition: sepa_eccuts.c:123
addBilinearTermToCut
static SCIP_RETCODE addBilinearTermToCut(SCIP *scip, SCIP_SOL *sol, SCIP_ROW *cut, SCIP_VAR *x, SCIP_VAR *y, SCIP_Real coeff, SCIP_Real *cutconstant, SCIP_Real *cutactivity, SCIP_Bool *success)
Definition: sepa_eccuts.c:2593
SCIP_OKAY
Definition: type_retcode.h:33
tcliqueMaxClique
void tcliqueMaxClique(TCLIQUE_GETNNODES((*getnnodes)), TCLIQUE_GETWEIGHTS((*getweights)), TCLIQUE_ISEDGE((*isedge)), TCLIQUE_SELECTADJNODES((*selectadjnodes)), TCLIQUE_GRAPH *tcliquegraph, TCLIQUE_NEWSOL((*newsol)), TCLIQUE_DATA *tcliquedata, int *maxcliquenodes, int *nmaxcliquenodes, TCLIQUE_WEIGHT *maxcliqueweight, TCLIQUE_WEIGHT maxfirstnodeweight, TCLIQUE_WEIGHT minweight, int maxntreenodes, int backtrackfreq, int maxnzeroextensions, int fixednode, int *ntreenodes, TCLIQUE_STATUS *status)
Definition: tclique_branch.c:1001
DEFAULT_MAXAGGRSIZE
#define DEFAULT_MAXAGGRSIZE
Definition: sepa_eccuts.c:62
DEFAULT_CUTMAXRANGE
#define DEFAULT_CUTMAXRANGE
Definition: sepa_eccuts.c:57
SCIP_CALL
#define SCIP_CALL(x)
Definition: def.h:384
SCIPensureBlockMemoryArray
#define SCIPensureBlockMemoryArray(scip, ptr, arraysizeptr, minsize)
Definition: scip_mem.h:98
SCIPgetNLPNlRows
SCIP_NLROW ** SCIPgetNLPNlRows(SCIP *scip)
Definition: scip_nlp.c:310
EcAggr::termsize
int termsize
Definition: sepa_eccuts.c:91
NlrowAggr::nlinvars
int nlinvars
Definition: sepa_eccuts.c:107
SCIPaddRow
SCIP_RETCODE SCIPaddRow(SCIP *scip, SCIP_ROW *row, SCIP_Bool forcecut, SCIP_Bool *infeasible)
Definition: scip_cut.c:241
NlrowAggr::lincoefs
SCIP_Real * lincoefs
Definition: sepa_eccuts.c:106
DEFAULT_MAXDEPTH
#define DEFAULT_MAXDEPTH
Definition: sepa_eccuts.c:54
SCIP_ERROR
Definition: type_retcode.h:34
DEFAULT_MINVIOLATION
#define DEFAULT_MINVIOLATION
Definition: sepa_eccuts.c:60
SCIPincludeSepaBasic
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:100
DEFAULT_MAXSEPACUTS
#define DEFAULT_MAXSEPACUTS
Definition: sepa_eccuts.c:55
SCIP_Expr
Definition: struct_expr.h:95
SCIPallocBufferArray
#define SCIPallocBufferArray(scip, ptr, num)
Definition: scip_mem.h:115
SCIPfreeTransform
SCIP_RETCODE SCIPfreeTransform(SCIP *scip)
Definition: scip_solve.c:3432
SCIPsetSepaExitsol
SCIP_RETCODE SCIPsetSepaExitsol(SCIP *scip, SCIP_SEPA *sepa, SCIP_DECL_SEPAEXITSOL((*sepaexitsol)))
Definition: scip_sepa.c:222
DEFAULT_MAXSEPACUTSROOT
#define DEFAULT_MAXSEPACUTSROOT
Definition: sepa_eccuts.c:56
SCIP_Bool
#define SCIP_Bool
Definition: def.h:84
SCIPincludeDefaultPlugins
SCIP_RETCODE SCIPincludeDefaultPlugins(SCIP *scip)
Definition: scipdefplugins.c:28
SCIPchgRowRhs
SCIP_RETCODE SCIPchgRowRhs(SCIP *scip, SCIP_ROW *row, SCIP_Real rhs)
Definition: scip_lp.c:1598
SCIPincludeSepaEccuts
SCIP_RETCODE SCIPincludeSepaEccuts(SCIP *scip)
Definition: sepa_eccuts.c:3127
SCIPgetDepth
int SCIPgetDepth(SCIP *scip)
Definition: scip_tree.c:661
NlrowAggr::linvars
SCIP_VAR ** linvars
Definition: sepa_eccuts.c:105
searchEcAggr
static SCIP_RETCODE searchEcAggr(SCIP *scip, SCIP_SEPADATA *sepadata, SCIP_NLROW *nlrow, SCIP_SOL *sol, SCIP_Bool rhsaggr, int *quadvar2aggr, int *nfound)
Definition: sepa_eccuts.c:1742
MAX
#define MAX(x, y)
Definition: tclique_def.h:83
SCIP_Row
Definition: struct_lp.h:192
transformValue
static SCIP_Real transformValue(SCIP *scip, SCIP_Real lb, SCIP_Real ub, SCIP_Real val)
Definition: sepa_eccuts.c:2232
SCIPaddPoolCut
SCIP_RETCODE SCIPaddPoolCut(SCIP *scip, SCIP_ROW *row)
Definition: scip_cut.c:352
SCIPsetIntParam
SCIP_RETCODE SCIPsetIntParam(SCIP *scip, const char *name, int value)
Definition: scip_param.c:478
SUBSCIP_NODELIMIT
#define SUBSCIP_NODELIMIT
Definition: sepa_eccuts.c:66
SCIPcreateEmptyRowSepa
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:1444
DEFAULT_MAXBILINTERMS
#define DEFAULT_MAXBILINTERMS
Definition: sepa_eccuts.c:63
SCIPgetCutEfficacy
SCIP_Real SCIPgetCutEfficacy(SCIP *scip, SCIP_SOL *sol, SCIP_ROW *cut)
Definition: scip_cut.c:85
SCIPaddSquareSecant
void SCIPaddSquareSecant(SCIP *scip, SCIP_Real sqrcoef, SCIP_Real lb, SCIP_Real ub, SCIP_Real *lincoef, SCIP_Real *linconstant, SCIP_Bool *success)
Definition: expr_pow.c:3294
NlrowAggr::nquadvars
int nquadvars
Definition: sepa_eccuts.c:113
SCIPgetNSols
int SCIPgetNSols(SCIP *scip)
Definition: scip_sol.c:2205
nlrowaggrCreate
static SCIP_RETCODE nlrowaggrCreate(SCIP *scip, SCIP_NLROW *nlrow, SCIP_NLROWAGGR **nlrowaggr, int *quadvar2aggr, int nfound, SCIP_Bool rhsaggr)
Definition: sepa_eccuts.c:423
EcAggr::nvars
int nvars
Definition: sepa_eccuts.c:84
SCIPaddBilinMcCormick
void SCIPaddBilinMcCormick(SCIP *scip, SCIP_Real bilincoef, SCIP_Real lbx, SCIP_Real ubx, SCIP_Real refpointx, SCIP_Real lby, SCIP_Real uby, SCIP_Real refpointy, SCIP_Bool overestimate, SCIP_Real *lincoefx, SCIP_Real *lincoefy, SCIP_Real *linconstant, SCIP_Bool *success)
Definition: nlhdlr_bilinear.c:1752
isCandidate
static SCIP_RETCODE isCandidate(SCIP *scip, SCIP_SEPADATA *sepadata, SCIP_NLROW *nlrow, SCIP_Bool *rhscandidate, SCIP_Bool *lhscandidate)
Definition: sepa_eccuts.c:1775
nlrowaggrAddQuadraticVar
static SCIP_RETCODE nlrowaggrAddQuadraticVar(SCIP *scip, SCIP_NLROWAGGR *nlrowaggr, SCIP_VAR *quadvar)
Definition: sepa_eccuts.c:378
SCIPisInfinity
SCIP_Bool SCIPisInfinity(SCIP *scip, SCIP_Real val)
Definition: scip_numerics.c:514
SCIPgetRowSolActivity
SCIP_Real SCIPgetRowSolActivity(SCIP *scip, SCIP_ROW *row, SCIP_SOL *sol)
Definition: scip_lp.c:2129
nlrowaggrStoreLinearTerms
static SCIP_RETCODE nlrowaggrStoreLinearTerms(SCIP *scip, SCIP_NLROWAGGR *nlrowaggr, SCIP_VAR **linvars, SCIP_Real *lincoefs, int nlinvars)
Definition: sepa_eccuts.c:304
BMSclearMemory
#define BMSclearMemory(ptr)
Definition: memory.h:122
SCIProwGetRank
int SCIProwGetRank(SCIP_ROW *row)
Definition: lp.c:17314
SCIPlpiChgBounds
SCIP_RETCODE SCIPlpiChgBounds(SCIP_LPI *lpi, int ncols, const int *ind, const SCIP_Real *lb, const SCIP_Real *ub)
Definition: lpi_clp.cpp:1075
SCIPgetNVars
int SCIPgetNVars(SCIP *scip)
Definition: scip_prob.c:1991
nlrowaggrFree
static SCIP_RETCODE nlrowaggrFree(SCIP *scip, SCIP_NLROWAGGR **nlrowaggr)
Definition: sepa_eccuts.c:645
SCIPnlrowGetLinearCoefs
SCIP_Real * SCIPnlrowGetLinearCoefs(SCIP_NLROW *nlrow)
Definition: nlp.c:1784
separateCuts
static SCIP_RETCODE separateCuts(SCIP *scip, SCIP_SEPA *sepa, SCIP_SEPADATA *sepadata, int depth, SCIP_SOL *sol, SCIP_RESULT *result)
Definition: sepa_eccuts.c:2925
cons_xor.h
Constraint handler for XOR constraints, .
SCIPreleaseRow
SCIP_RETCODE SCIPreleaseRow(SCIP *scip, SCIP_ROW **row)
Definition: scip_lp.c:1553
CLIQUE_BACKTRACKFREQ
#define CLIQUE_BACKTRACKFREQ
Definition: sepa_eccuts.c:49
addLinearTermToCut
static SCIP_RETCODE addLinearTermToCut(SCIP *scip, SCIP_SOL *sol, SCIP_ROW *cut, SCIP_VAR *x, SCIP_Real coeff, SCIP_Real *cutconstant, SCIP_Real *cutactivity, SCIP_Bool *success)
Definition: sepa_eccuts.c:2542
SCIPsetSepaFree
SCIP_RETCODE SCIPsetSepaFree(SCIP *scip, SCIP_SEPA *sepa, SCIP_DECL_SEPAFREE((*sepafree)))
Definition: scip_sepa.c:158
NlrowAggr::quadvarssize
int quadvarssize
Definition: sepa_eccuts.c:114
SCIPgetBestSol
SCIP_SOL * SCIPgetBestSol(SCIP *scip)
Definition: scip_sol.c:2304
SCIPisGT
SCIP_Bool SCIPisGT(SCIP *scip, SCIP_Real val1, SCIP_Real val2)
Definition: scip_numerics.c:477
createLP
static SCIP_RETCODE createLP(SCIP *scip, SCIP_SEPADATA *sepadata)
Definition: sepa_eccuts.c:2082
SCIPaddVar
SCIP_RETCODE SCIPaddVar(SCIP *scip, SCIP_VAR *var)
Definition: scip_prob.c:1667
findAndStoreEcAggregations
static SCIP_RETCODE findAndStoreEcAggregations(SCIP *scip, SCIP_SEPADATA *sepadata, SCIP_NLROW *nlrow, SCIP_SOL *sol)
Definition: sepa_eccuts.c:1906
DEFAULT_MAXROUNDSROOT
#define DEFAULT_MAXROUNDSROOT
Definition: sepa_eccuts.c:53
nlrowaggrAddLinearTerm
static SCIP_RETCODE nlrowaggrAddLinearTerm(SCIP *scip, SCIP_NLROWAGGR *nlrowaggr, SCIP_VAR *linvar, SCIP_Real lincoef)
Definition: sepa_eccuts.c:345
SCIPreleaseCons
SCIP_RETCODE SCIPreleaseCons(SCIP *scip, SCIP_CONS **cons)
Definition: scip_cons.c:1110
a
SCIP_VAR * a
Definition: circlepacking.c:57
SCIPisExprVar
SCIP_Bool SCIPisExprVar(SCIP *scip, SCIP_EXPR *expr)
Definition: scip_expr.c:1421
SCIPexprGetQuadraticBilinTerm
void SCIPexprGetQuadraticBilinTerm(SCIP_EXPR *expr, int termidx, SCIP_EXPR **expr1, SCIP_EXPR **expr2, SCIP_Real *coef, int *pos2, SCIP_EXPR **prodexpr)
Definition: expr.c:4147
SCIPgetVars
SCIP_VAR ** SCIPgetVars(SCIP *scip)
Definition: scip_prob.c:1946
DEFAULT_DYNAMICCUTS
#define DEFAULT_DYNAMICCUTS
Definition: sepa_eccuts.c:51
SCIPstatistic
#define SCIPstatistic(x)
Definition: pub_message.h:111
SCIP_Real
#define SCIP_Real
Definition: def.h:177
SCIP_LPi
Definition: lpi_clp.cpp:95
SCIPisStopped
SCIP_Bool SCIPisStopped(SCIP *scip)
Definition: scip_general.c:694
SCIP_DECL_SEPAEXECLP
static SCIP_DECL_SEPAEXECLP(sepaExeclpEccuts)
Definition: sepa_eccuts.c:3057
SCIP_CALL_TERMINATE
#define SCIP_CALL_TERMINATE(retcode, x, TERM)
Definition: def.h:405
y
SCIP_VAR ** y
Definition: circlepacking.c:55
SCIP_NlRow
Definition: struct_nlp.h:55
SCIP_DECL_SEPAFREE
static SCIP_DECL_SEPAFREE(sepaFreeEccuts)
Definition: sepa_eccuts.c:3015
createMIP
static SCIP_RETCODE createMIP(SCIP *scip, SCIP *subscip, SCIP_SEPADATA *sepadata, SCIP_NLROW *nlrow, SCIP_Bool rhsaggr, SCIP_VAR **forwardarcs, SCIP_VAR **backwardarcs, SCIP_Real *nodeweights, int *nedges, int *narcs)
Definition: sepa_eccuts.c:862
SCIPprintRow
SCIP_RETCODE SCIPprintRow(SCIP *scip, SCIP_ROW *row, FILE *file)
Definition: scip_lp.c:2197
SCIPprintNlRow
SCIP_RETCODE SCIPprintNlRow(SCIP *scip, SCIP_NLROW *nlrow, FILE *file)
Definition: scip_nlp.c:1547
SCIPvarGetIndex
int SCIPvarGetIndex(SCIP_VAR *var)
Definition: var.c:17590
SCIPgetTotalTime
SCIP_Real SCIPgetTotalTime(SCIP *scip)
Definition: scip_timing.c:342
sepadataCreate
static SCIP_RETCODE sepadataCreate(SCIP *scip, SCIP_SEPADATA **sepadata)
Definition: sepa_eccuts.c:721
evalCorner
static SCIP_Real evalCorner(SCIP_ECAGGR *ecaggr, int k)
Definition: sepa_eccuts.c:2194
SCIPlpiChgObj
SCIP_RETCODE SCIPlpiChgObj(SCIP_LPI *lpi, int ncols, const int *ind, const SCIP_Real *obj)
Definition: lpi_clp.cpp:1231
SCIPisZero
SCIP_Bool SCIPisZero(SCIP *scip, SCIP_Real val)
Definition: scip_numerics.c:551
SCIPisLE
SCIP_Bool SCIPisLE(SCIP *scip, SCIP_Real val1, SCIP_Real val2)
Definition: scip_numerics.c:464
tcliqueFlush
TCLIQUE_Bool tcliqueFlush(TCLIQUE_GRAPH *tcliquegraph)
Definition: tclique_graph.c:408
ADJUSTFACETTOL
#define ADJUSTFACETTOL
Definition: sepa_eccuts.c:68
computeCut
static SCIP_RETCODE computeCut(SCIP *scip, SCIP_SEPA *sepa, SCIP_SEPADATA *sepadata, SCIP_NLROWAGGR *nlrowaggr, SCIP_SOL *sol, SCIP_Bool *separated, SCIP_Bool *cutoff)
Definition: sepa_eccuts.c:2707
SCIPvarGetUbLocal
SCIP_Real SCIPvarGetUbLocal(SCIP_VAR *var)
Definition: var.c:17976
nnodes
#define nnodes
Definition: gastrans.c:65
SCIPfreeBlockMemoryArrayNull
#define SCIPfreeBlockMemoryArrayNull(scip, ptr, num)
Definition: scip_mem.h:102
DEFAULT_MINAGGRSIZE
#define DEFAULT_MINAGGRSIZE
Definition: sepa_eccuts.c:61
SCIPaddSquareLinearization
void SCIPaddSquareLinearization(SCIP *scip, SCIP_Real sqrcoef, SCIP_Real refpoint, SCIP_Bool isint, SCIP_Real *lincoef, SCIP_Real *linconstant, SCIP_Bool *success)
Definition: expr_pow.c:3226
BMSclearMemoryArray
#define BMSclearMemoryArray(ptr, num)
Definition: memory.h:123
SCIPallocBlockMemory
SCIPallocBlockMemory(scip, subsol))
SCIPnlrowGetConstant
SCIP_Real SCIPnlrowGetConstant(SCIP_NLROW *nlrow)
Definition: nlp.c:1754
ecaggrAddBilinTerm
static SCIP_RETCODE ecaggrAddBilinTerm(SCIP *scip, SCIP_ECAGGR *ecaggr, SCIP_VAR *x, SCIP_VAR *y, SCIP_Real coef)
Definition: sepa_eccuts.c:230
SCIPnlrowGetLhs
SCIP_Real SCIPnlrowGetLhs(SCIP_NLROW *nlrow)
Definition: nlp.c:1804
poweroftwo
static const int poweroftwo[]
Definition: sepa_eccuts.c:72
SEPA_USESSUBSCIP
#define SEPA_USESSUBSCIP
Definition: sepa_eccuts.c:41
SCIPvarIsIntegral
SCIP_Bool SCIPvarIsIntegral(SCIP_VAR *var)
Definition: var.c:17442
DEFAULT_MAXROUNDS
#define DEFAULT_MAXROUNDS
Definition: sepa_eccuts.c:52
SEPA_MAXBOUNDDIST
#define SEPA_MAXBOUNDDIST
Definition: sepa_eccuts.c:40
EcAggr::termvars2
int * termvars2
Definition: sepa_eccuts.c:89
SCIPgetSolVal
SCIP_Real SCIPgetSolVal(SCIP *scip, SCIP_SOL *sol, SCIP_VAR *var)
Definition: scip_sol.c:1352
scipdefplugins.h
default SCIP plugins
SCIP_STATUS_INFEASIBLE
Definition: type_stat.h:53
SCIPaddRealParam
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:130
ecaggrAddQuadvar
static SCIP_RETCODE ecaggrAddQuadvar(SCIP_ECAGGR *ecaggr, SCIP_VAR *x)
Definition: sepa_eccuts.c:219
updateMIP
static SCIP_RETCODE updateMIP(SCIP *subscip, SCIP_NLROW *nlrow, SCIP_VAR **forwardarcs, SCIP_VAR **backwardarcs, int *quadvar2aggr, int *nedges)
Definition: sepa_eccuts.c:1076
SCIPsetLongintParam
SCIP_RETCODE SCIPsetLongintParam(SCIP *scip, const char *name, SCIP_Longint value)
Definition: scip_param.c:536
SCIP_CUTOFF
Definition: type_result.h:39
SCIP_SEPADATA
struct SCIP_SepaData SCIP_SEPADATA
Definition: type_sepa.h:43
SCIP_PARAMSETTING_AGGRESSIVE
Definition: type_paramset.h:52
SCIPaddBoolParam
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:48
SCIPnlrowGetExpr
SCIP_EXPR * SCIPnlrowGetExpr(SCIP_NLROW *nlrow)
Definition: nlp.c:1794
SCIPfree
SCIP_RETCODE SCIPfree(SCIP **scip)
Definition: scip_general.c:315
checkRikun
static SCIP_Bool checkRikun(SCIP *scip, SCIP_ECAGGR *ecaggr, SCIP_Real *fvals, SCIP_Real *facet)
Definition: sepa_eccuts.c:2011
SCIPprintSol
SCIP_RETCODE SCIPprintSol(SCIP *scip, SCIP_SOL *sol, FILE *file, SCIP_Bool printzeros)
Definition: scip_sol.c:1766