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