Scippy

SCIP

Solving Constraint Integer Programs

cons_knapsack.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-2017 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 email to scip@zib.de. */
13 /* */
14 /* * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * */
15 
16 /**@file cons_knapsack.c
17  * @brief Constraint handler for knapsack constraints of the form \f$a^T x \le b\f$, x binary and \f$a \ge 0\f$.
18  * @author Tobias Achterberg
19  * @author Xin Liu
20  * @author Kati Wolter
21  * @author Michael Winkler
22  */
23 
24 /*---+----1----+----2----+----3----+----4----+----5----+----6----+----7----+----8----+----9----+----0----+----1----+----2*/
25 
26 #include <assert.h>
27 #include <string.h>
28 #include <limits.h>
29 #include <stdio.h>
30 #include <ctype.h>
31 
32 #include "scip/cons_knapsack.h"
33 #include "scip/cons_linear.h"
34 #include "scip/cons_logicor.h"
35 #include "scip/cons_setppc.h"
36 #include "scip/pub_misc.h"
37 
38 /* constraint handler properties */
39 #define CONSHDLR_NAME "knapsack"
40 #define CONSHDLR_DESC "knapsack constraint of the form a^T x <= b, x binary and a >= 0"
41 #define CONSHDLR_SEPAPRIORITY +600000 /**< priority of the constraint handler for separation */
42 #define CONSHDLR_ENFOPRIORITY -600000 /**< priority of the constraint handler for constraint enforcing */
43 #define CONSHDLR_CHECKPRIORITY -600000 /**< priority of the constraint handler for checking feasibility */
44 #define CONSHDLR_SEPAFREQ 0 /**< frequency for separating cuts; zero means to separate only in the root node */
45 #define CONSHDLR_PROPFREQ 1 /**< frequency for propagating domains; zero means only preprocessing propagation */
46 #define CONSHDLR_EAGERFREQ 100 /**< frequency for using all instead of only the useful constraints in separation,
47  * propagation and enforcement, -1 for no eager evaluations, 0 for first only */
48 #define CONSHDLR_MAXPREROUNDS -1 /**< maximal number of presolving rounds the constraint handler participates in (-1: no limit) */
49 #define CONSHDLR_DELAYSEPA FALSE /**< should separation method be delayed, if other separators found cuts? */
50 #define CONSHDLR_DELAYPROP FALSE /**< should propagation method be delayed, if other propagators found reductions? */
51 #define CONSHDLR_NEEDSCONS TRUE /**< should the constraint handler be skipped, if no constraints are available? */
52 
53 #define CONSHDLR_PRESOLTIMING SCIP_PRESOLTIMING_ALWAYS
54 #define CONSHDLR_PROP_TIMING SCIP_PROPTIMING_BEFORELP
55 
56 #define EVENTHDLR_NAME "knapsack"
57 #define EVENTHDLR_DESC "bound change event handler for knapsack constraints"
58 #define EVENTTYPE_KNAPSACK SCIP_EVENTTYPE_LBCHANGED \
59  | SCIP_EVENTTYPE_UBTIGHTENED \
60  | SCIP_EVENTTYPE_VARFIXED \
61  | SCIP_EVENTTYPE_VARDELETED \
62  | SCIP_EVENTTYPE_IMPLADDED /**< variable events that should be caught by the event handler */
63 
64 #define LINCONSUPGD_PRIORITY +100000 /**< priority of the constraint handler for upgrading of linear constraints */
65 
66 #define MAX_USECLIQUES_SIZE 1000 /**< maximal number of items in knapsack where clique information is used */
67 #define MAX_ZEROITEMS_SIZE 10000 /**< maximal number of items to store in the zero list in preprocessing */
68 
69 #define KNAPSACKRELAX_MAXDELTA 0.1 /**< maximal allowed rounding distance for scaling in knapsack relaxation */
70 #define KNAPSACKRELAX_MAXDNOM 1000LL /**< maximal allowed denominator in knapsack rational relaxation */
71 #define KNAPSACKRELAX_MAXSCALE 1000.0 /**< maximal allowed scaling factor in knapsack rational relaxation */
72 
73 #define DEFAULT_SEPACARDFREQ 1 /**< multiplier on separation frequency, how often knapsack cuts are separated */
74 #define DEFAULT_MAXROUNDS 5 /**< maximal number of separation rounds per node (-1: unlimited) */
75 #define DEFAULT_MAXROUNDSROOT -1 /**< maximal number of separation rounds in the root node (-1: unlimited) */
76 #define DEFAULT_MAXSEPACUTS 50 /**< maximal number of cuts separated per separation round */
77 #define DEFAULT_MAXSEPACUTSROOT 200 /**< maximal number of cuts separated per separation round in the root node */
78 #define DEFAULT_MAXCARDBOUNDDIST 0.0 /**< maximal relative distance from current node's dual bound to primal bound compared
79  * to best node's dual bound for separating knapsack cuts */
80 #define DEFAULT_DISAGGREGATION TRUE /**< should disaggregation of knapsack constraints be allowed in preprocessing? */
81 #define DEFAULT_SIMPLIFYINEQUALITIES TRUE/**< should presolving try to simplify knapsacks */
82 #define DEFAULT_NEGATEDCLIQUE TRUE /**< should negated clique information be used in solving process */
83 
84 #define MAXABSVBCOEF 1e+5 /**< maximal absolute coefficient in variable bounds used for knapsack relaxation */
85 #define USESUPADDLIFT FALSE /**< should lifted minimal cover inequalities using superadditive up-lifting be separated in addition */
86 
87 #define DEFAULT_PRESOLUSEHASHING TRUE /**< should hash table be used for detecting redundant constraints in advance */
88 #define HASHSIZE_KNAPSACKCONS 500 /**< minimal size of hash table in linear constraint tables */
89 
90 #define DEFAULT_PRESOLPAIRWISE TRUE /**< should pairwise constraint comparison be performed in presolving? */
91 #define NMINCOMPARISONS 200000 /**< number for minimal pairwise presolving comparisons */
92 #define MINGAINPERNMINCOMPARISONS 1e-06 /**< minimal gain per minimal pairwise presolving comparisons to repeat pairwise
93  * comparison round */
94 #define DEFAULT_DUALPRESOLVING TRUE /**< should dual presolving steps be performed? */
95 #define DEFAULT_DETECTCUTOFFBOUND TRUE /**< should presolving try to detect constraints parallel to the objective
96  * function defining an upper bound and prevent these constraints from
97  * entering the LP */
98 #define DEFAULT_DETECTLOWERBOUND TRUE /**< should presolving try to detect constraints parallel to the objective
99  * function defining a lower bound and prevent these constraints from
100  * entering the LP */
101 #define DEFAULT_CLIQUEEXTRACTFACTOR 0.5 /**< lower clique size limit for greedy clique extraction algorithm (relative to largest clique) */
102 
103 #define MAXCOVERSIZEITERLEWI 1000 /**< maximal size for which LEWI are iteratively separated by reducing the feasible set */
104 
105 #define DEFAULT_USEGUBS FALSE /**< should GUB information be used for separation? */
106 #define GUBCONSGROWVALUE 6 /**< memory growing value for GUB constraint array */
107 #define GUBSPLITGNC1GUBS FALSE /**< should GNC1 GUB conss without F vars be split into GOC1 and GR GUB conss? */
108 #define DEFAULT_CLQPARTUPDATEFAC 1.5 /**< factor on the growth of global cliques to decide when to update a previous
109  * (negated) clique partition (used only if updatecliquepartitions is set to TRUE) */
110 #define DEFAULT_UPDATECLIQUEPARTITIONS FALSE /**< should clique partition information be updated when old partition seems outdated? */
111 
112 #define MAXNCLIQUEVARSCOMP 1000000 /**< limit on number of pairwise comparisons in clique partitioning algorithm */
114 /* @todo maybe use event SCIP_EVENTTYPE_VARUNLOCKED to decide for another dual-presolving run on a constraint */
116 /*
117  * Data structures
118  */
119 
120 /** constraint handler data */
121 struct SCIP_ConshdlrData
122 {
123  int* ints1; /**< cleared memory array, all entries are set to zero in initpre, if you use this
124  * you have to clear it at the end, exists only in presolving stage */
125  int* ints2; /**< cleared memory array, all entries are set to zero in initpre, if you use this
126  * you have to clear it at the end, exists only in presolving stage */
127  SCIP_Longint* longints1; /**< cleared memory array, all entries are set to zero in initpre, if you use this
128  * you have to clear it at the end, exists only in presolving stage */
129  SCIP_Longint* longints2; /**< cleared memory array, all entries are set to zero in initpre, if you use this
130  * you have to clear it at the end, exists only in presolving stage */
131  SCIP_Bool* bools1; /**< cleared memory array, all entries are set to zero in initpre, if you use this
132  * you have to clear it at the end, exists only in presolving stage */
133  SCIP_Bool* bools2; /**< cleared memory array, all entries are set to zero in initpre, if you use this
134  * you have to clear it at the end, exists only in presolving stage */
135  SCIP_Bool* bools3; /**< cleared memory array, all entries are set to zero in initpre, if you use this
136  * you have to clear it at the end, exists only in presolving stage */
137  SCIP_Bool* bools4; /**< cleared memory array, all entries are set to zero in initpre, if you use this
138  * you have to clear it at the end, exists only in presolving stage */
139  SCIP_Real* reals1; /**< cleared memory array, all entries are set to zero in consinit, if you use this
140  * you have to clear it at the end */
141  int ints1size; /**< size of ints1 array */
142  int ints2size; /**< size of ints2 array */
143  int longints1size; /**< size of longints1 array */
144  int longints2size; /**< size of longints2 array */
145  int bools1size; /**< size of bools1 array */
146  int bools2size; /**< size of bools2 array */
147  int bools3size; /**< size of bools3 array */
148  int bools4size; /**< size of bools4 array */
149  int reals1size; /**< size of reals1 array */
150  SCIP_EVENTHDLR* eventhdlr; /**< event handler for bound change events */
151  SCIP_Real maxcardbounddist; /**< maximal relative distance from current node's dual bound to primal bound compared
152  * to best node's dual bound for separating knapsack cuts */
153  int sepacardfreq; /**< multiplier on separation frequency, how often knapsack cuts are separated */
154  int maxrounds; /**< maximal number of separation rounds per node (-1: unlimited) */
155  int maxroundsroot; /**< maximal number of separation rounds in the root node (-1: unlimited) */
156  int maxsepacuts; /**< maximal number of cuts separated per separation round */
157  int maxsepacutsroot; /**< maximal number of cuts separated per separation round in the root node */
158  SCIP_Bool disaggregation; /**< should disaggregation of knapsack constraints be allowed in preprocessing? */
159  SCIP_Bool simplifyinequalities;/**< should presolving try to cancel down or delete coefficients in inequalities */
160  SCIP_Bool negatedclique; /**< should negated clique information be used in solving process */
161  SCIP_Bool presolpairwise; /**< should pairwise constraint comparison be performed in presolving? */
162  SCIP_Bool presolusehashing; /**< should hash table be used for detecting redundant constraints in advance */
163  SCIP_Bool dualpresolving; /**< should dual presolving steps be performed? */
164  SCIP_Bool usegubs; /**< should GUB information be used for separation? */
165  SCIP_Bool detectcutoffbound; /**< should presolving try to detect constraints parallel to the objective
166  * function defining an upper bound and prevent these constraints from
167  * entering the LP */
168  SCIP_Bool detectlowerbound; /**< should presolving try to detect constraints parallel to the objective
169  * function defining a lower bound and prevent these constraints from
170  * entering the LP */
171  SCIP_Bool updatecliquepartitions; /**< should clique partition information be updated when old partition seems outdated? */
172  SCIP_Real cliqueextractfactor;/**< lower clique size limit for greedy clique extraction algorithm (relative to largest clique) */
173  SCIP_Real clqpartupdatefac; /**< factor on the growth of global cliques to decide when to update a previous
174  * (negated) clique partition (used only if updatecliquepartitions is set to TRUE) */
175 };
176 
177 
178 /** constraint data for knapsack constraints */
179 struct SCIP_ConsData
180 {
181  SCIP_VAR** vars; /**< variables in knapsack constraint */
182  SCIP_Longint* weights; /**< weights of variables in knapsack constraint */
183  SCIP_EVENTDATA** eventdata; /**< event data for bound change events of the variables */
184  int* cliquepartition; /**< clique indices of the clique partition */
185  int* negcliquepartition; /**< clique indices of the negated clique partition */
186  SCIP_ROW* row; /**< corresponding LP row */
187  int nvars; /**< number of variables in knapsack constraint */
188  int varssize; /**< size of vars, weights, and eventdata arrays */
189  int ncliques; /**< number of cliques in the clique partition */
190  int nnegcliques; /**< number of cliques in the negated clique partition */
191  int ncliqueslastnegpart;/**< number of global cliques the last time a negated clique partition was computed */
192  int ncliqueslastpart; /**< number of global cliques the last time a clique partition was computed */
193  SCIP_Longint capacity; /**< capacity of knapsack */
194  SCIP_Longint weightsum; /**< sum of all weights */
195  SCIP_Longint onesweightsum; /**< sum of weights of variables fixed to one */
196  unsigned int presolvedtiming:5; /**< max level in which the knapsack constraint is already presolved */
197  unsigned int sorted:1; /**< are the knapsack items sorted by weight? */
198  unsigned int cliquepartitioned:1;/**< is the clique partition valid? */
199  unsigned int negcliquepartitioned:1;/**< is the negated clique partition valid? */
200  unsigned int merged:1; /**< are the constraint's equal variables already merged? */
201  unsigned int cliquesadded:1; /**< were the cliques of the knapsack already added to clique table? */
202  unsigned int varsdeleted:1; /**< were variables deleted after last cleanup? */
203  unsigned int existmultaggr:1; /**< does this constraint contain multi-aggregations */
204 };
205 
206 /** event data for bound changes events */
207 struct SCIP_EventData
208 {
209  SCIP_CONS* cons; /**< knapsack constraint to process the bound change for */
210  SCIP_Longint weight; /**< weight of variable */
211  int filterpos; /**< position of event in variable's event filter */
212 };
213 
214 
215 /** data structure to combine two sorting key values */
216 struct sortkeypair
217 {
218  SCIP_Real key1; /**< first sort key value */
219  SCIP_Real key2; /**< second sort key value */
220 };
221 typedef struct sortkeypair SORTKEYPAIR;
222 
223 /** status of GUB constraint */
224 enum GUBVarstatus
225 {
226  GUBVARSTATUS_UNINITIAL = -1, /** unintitialized variable status */
227  GUBVARSTATUS_CAPACITYEXCEEDED = 0, /** variable with weight exceeding the knapsack capacity */
228  GUBVARSTATUS_BELONGSTOSET_R = 1, /** variable in noncovervars R */
229  GUBVARSTATUS_BELONGSTOSET_F = 2, /** variable in noncovervars F */
230  GUBVARSTATUS_BELONGSTOSET_C2 = 3, /** variable in covervars C2 */
231  GUBVARSTATUS_BELONGSTOSET_C1 = 4 /** variable in covervars C1 */
232 };
233 typedef enum GUBVarstatus GUBVARSTATUS;
235 /** status of variable in GUB constraint */
237 {
238  GUBCONSSTATUS_UNINITIAL = -1, /** unintitialized GUB constraint status */
239  GUBCONSSTATUS_BELONGSTOSET_GR = 0, /** all GUB variables are in noncovervars R */
240  GUBCONSSTATUS_BELONGSTOSET_GF = 1, /** all GUB variables are in noncovervars F (and noncovervars R) */
241  GUBCONSSTATUS_BELONGSTOSET_GC2 = 2, /** all GUB variables are in covervars C2 */
242  GUBCONSSTATUS_BELONGSTOSET_GNC1 = 3, /** some GUB variables are in covervars C1, others in noncovervars R or F */
243  GUBCONSSTATUS_BELONGSTOSET_GOC1 = 4 /** all GUB variables are in covervars C1 */
244 };
245 typedef enum GUBConsstatus GUBCONSSTATUS;
247 /** data structure of GUB constraints */
249 {
250  int* gubvars; /**< indices of GUB variables in knapsack constraint */
251  GUBVARSTATUS* gubvarsstatus; /**< status of GUB variables */
252  int ngubvars; /**< number of GUB variables */
253  int gubvarssize; /**< size of gubvars array */
254 };
255 typedef struct SCIP_GUBCons SCIP_GUBCONS;
257 /** data structure of a set of GUB constraints */
259 {
260  SCIP_GUBCONS** gubconss; /**< GUB constraints in GUB set */
261  GUBCONSSTATUS* gubconsstatus; /**< status of GUB constraints */
262  int ngubconss; /**< number of GUB constraints */
263  int nvars; /**< number of variables in knapsack constraint */
264  int* gubconssidx; /**< index of GUB constraint (in gubconss array) of each knapsack variable */
265  int* gubvarsidx; /**< index in GUB constraint (in gubvars array) of each knapsack variable */
266 };
267 typedef struct SCIP_GUBSet SCIP_GUBSET;
269 /*
270  * Local methods
271  */
273 /** comparison method for two sorting key pairs */
274 static
275 SCIP_DECL_SORTPTRCOMP(compSortkeypairs)
276 {
277  SORTKEYPAIR* sortkeypair1 = (SORTKEYPAIR*)elem1;
278  SORTKEYPAIR* sortkeypair2 = (SORTKEYPAIR*)elem2;
279 
280  if( sortkeypair1->key1 < sortkeypair2->key1 )
281  return -1;
282  else if( sortkeypair1->key1 > sortkeypair2->key1 )
283  return +1;
284  else if( sortkeypair1->key2 < sortkeypair2->key2 )
285  return -1;
286  else if( sortkeypair1->key2 > sortkeypair2->key2 )
287  return +1;
288  else
289  return 0;
290 }
291 
292 /** creates event data */
293 static
295  SCIP* scip, /**< SCIP data structure */
296  SCIP_EVENTDATA** eventdata, /**< pointer to store event data */
297  SCIP_CONS* cons, /**< constraint */
298  SCIP_Longint weight /**< weight of variable */
299  )
300 {
301  assert(eventdata != NULL);
303  SCIP_CALL( SCIPallocBlockMemory(scip, eventdata) );
304  (*eventdata)->cons = cons;
305  (*eventdata)->weight = weight;
306 
307  return SCIP_OKAY;
308 }
309 
310 /** frees event data */
311 static
313  SCIP* scip, /**< SCIP data structure */
314  SCIP_EVENTDATA** eventdata /**< pointer to event data */
315  )
316 {
317  assert(eventdata != NULL);
318 
319  SCIPfreeBlockMemory(scip, eventdata);
321  return SCIP_OKAY;
322 }
323 
324 /** sorts items in knapsack with nonincreasing weights */
325 static
326 void sortItems(
327  SCIP_CONSDATA* consdata /**< constraint data */
328  )
329 {
330  assert(consdata != NULL);
331  assert(consdata->nvars == 0 || consdata->vars != NULL);
332  assert(consdata->nvars == 0 || consdata->weights != NULL);
333  assert(consdata->nvars == 0 || consdata->eventdata != NULL);
334  assert(consdata->nvars == 0 || (consdata->cliquepartition != NULL && consdata->negcliquepartition != NULL));
335 
336  if( !consdata->sorted )
337  {
338  int pos;
339  int lastcliquenum;
340  int v;
341 
342  /* sort of five joint arrays of Long/pointer/pointer/ints/ints,
343  * sorted by first array in non-increasing order via sort template */
345  consdata->weights,
346  (void**)consdata->vars,
347  (void**)consdata->eventdata,
348  consdata->cliquepartition,
349  consdata->negcliquepartition,
350  consdata->nvars);
351 
352  v = consdata->nvars - 1;
353  /* sort all items with same weight according to their variable index, used for hash value for fast pairwise comparison of all constraints */
354  while( v >= 0 )
355  {
356  int w = v - 1;
357 
358  while( w >= 0 && consdata->weights[v] == consdata->weights[w] )
359  --w;
360 
361  if( v - w > 1 )
362  {
363  /* sort all corresponding parts of arrays for which the weights are equal by using the variable index */
365  (void**)(&(consdata->vars[w+1])),
366  (void**)(&(consdata->eventdata[w+1])),
367  &(consdata->cliquepartition[w+1]),
368  &(consdata->negcliquepartition[w+1]),
369  SCIPvarComp,
370  v - w);
371  }
372  v = w;
373  }
374 
375  /* we need to make sure that our clique numbers of our normal clique will be in increasing order without gaps */
376  if( consdata->cliquepartitioned )
377  {
378  lastcliquenum = 0;
379 
380  for( pos = 0; pos < consdata->nvars; ++pos )
381  {
382  /* if the clique number in the normal clique at position pos is greater than the last found clique number the
383  * partition is invalid */
384  if( consdata->cliquepartition[pos] > lastcliquenum )
385  {
386  consdata->cliquepartitioned = FALSE;
387  break;
388  }
389  else if( consdata->cliquepartition[pos] == lastcliquenum )
390  ++lastcliquenum;
391  }
392  }
393  /* we need to make sure that our clique numbers of our negated clique will be in increasing order without gaps */
394  if( consdata->negcliquepartitioned )
395  {
396  lastcliquenum = 0;
397 
398  for( pos = 0; pos < consdata->nvars; ++pos )
399  {
400  /* if the clique number in the negated clique at position pos is greater than the last found clique number the
401  * partition is invalid */
402  if( consdata->negcliquepartition[pos] > lastcliquenum )
403  {
404  consdata->negcliquepartitioned = FALSE;
405  break;
406  }
407  else if( consdata->negcliquepartition[pos] == lastcliquenum )
408  ++lastcliquenum;
409  }
410  }
411 
412  consdata->sorted = TRUE;
413  }
414 #ifndef NDEBUG
415  {
416  /* check if the weight array is sorted in a non-increasing way */
417  int i;
418  for( i = 0; i < consdata->nvars-1; ++i )
419  assert(consdata->weights[i] >= consdata->weights[i+1]);
420  }
421 #endif
422 }
423 
424 /** calculates a partition of the variables into cliques */
425 static
427  SCIP* scip, /**< SCIP data structure */
428  SCIP_CONSHDLRDATA* conshdlrdata, /**< knapsack constraint handler data */
429  SCIP_CONSDATA* consdata, /**< constraint data */
430  SCIP_Bool normalclique, /**< Should normal cliquepartition be created? */
431  SCIP_Bool negatedclique /**< Should negated cliquepartition be created? */
432  )
433 {
434  SCIP_Bool ispartitionoutdated;
435  SCIP_Bool isnegpartitionoutdated;
436  assert(consdata != NULL);
437  assert(consdata->nvars == 0 || (consdata->cliquepartition != NULL && consdata->negcliquepartition != NULL));
438 
439  /* rerun eventually if number of global cliques increased considerably since last partition */
440  ispartitionoutdated = (conshdlrdata->updatecliquepartitions && consdata->ncliques > 1
441  && SCIPgetNCliques(scip) >= (int)(conshdlrdata->clqpartupdatefac * consdata->ncliqueslastpart));
442 
443  if( normalclique && ( !consdata->cliquepartitioned || ispartitionoutdated ) )
444  {
445  SCIP_CALL( SCIPcalcCliquePartition(scip, consdata->vars, consdata->nvars, consdata->cliquepartition, &consdata->ncliques) );
446  consdata->cliquepartitioned = TRUE;
447  consdata->ncliqueslastpart = SCIPgetNCliques(scip);
448  }
449 
450  /* rerun eventually if number of global cliques increased considerably since last negated partition */
451  isnegpartitionoutdated = (conshdlrdata->updatecliquepartitions && consdata->nnegcliques > 1
452  && SCIPgetNCliques(scip) >= (int)(conshdlrdata->clqpartupdatefac * consdata->ncliqueslastnegpart));
453 
454  if( negatedclique && (!consdata->negcliquepartitioned || isnegpartitionoutdated) )
455  {
456  SCIP_CALL( SCIPcalcNegatedCliquePartition(scip, consdata->vars, consdata->nvars, consdata->negcliquepartition, &consdata->nnegcliques) );
457  consdata->negcliquepartitioned = TRUE;
458  consdata->ncliqueslastnegpart = SCIPgetNCliques(scip);
459  }
460  assert(!consdata->cliquepartitioned || consdata->ncliques <= consdata->nvars);
461  assert(!consdata->negcliquepartitioned || consdata->nnegcliques <= consdata->nvars);
462 
463  return SCIP_OKAY;
464 }
465 
466 /** installs rounding locks for the given variable in the given knapsack constraint */
467 static
469  SCIP* scip, /**< SCIP data structure */
470  SCIP_CONS* cons, /**< knapsack constraint */
471  SCIP_VAR* var /**< variable of constraint entry */
472  )
473 {
474  /* rounding up may violate the constraint */
475  SCIP_CALL( SCIPlockVarCons(scip, var, cons, FALSE, TRUE) );
477  return SCIP_OKAY;
478 }
479 
480 /** removes rounding locks for the given variable in the given knapsack constraint */
481 static
483  SCIP* scip, /**< SCIP data structure */
484  SCIP_CONS* cons, /**< knapsack constraint */
485  SCIP_VAR* var /**< variable of constraint entry */
486  )
487 {
488  /* rounding up may violate the constraint */
489  SCIP_CALL( SCIPunlockVarCons(scip, var, cons, FALSE, TRUE) );
491  return SCIP_OKAY;
492 }
493 
494 /** catches bound change events for variables in knapsack */
495 static
497  SCIP* scip, /**< SCIP data structure */
498  SCIP_CONS* cons, /**< constraint */
499  SCIP_CONSDATA* consdata, /**< constraint data */
500  SCIP_EVENTHDLR* eventhdlr /**< event handler to call for the event processing */
501  )
502 {
503  int i;
505  assert(cons != NULL);
506  assert(consdata != NULL);
507  assert(consdata->nvars == 0 || consdata->vars != NULL);
508  assert(consdata->nvars == 0 || consdata->weights != NULL);
509  assert(consdata->nvars == 0 || consdata->eventdata != NULL);
510 
511  for( i = 0; i < consdata->nvars; i++)
512  {
513  SCIP_CALL( eventdataCreate(scip, &consdata->eventdata[i], cons, consdata->weights[i]) );
514  SCIP_CALL( SCIPcatchVarEvent(scip, consdata->vars[i], EVENTTYPE_KNAPSACK,
515  eventhdlr, consdata->eventdata[i], &consdata->eventdata[i]->filterpos) );
516  }
517 
518  return SCIP_OKAY;
519 }
520 
521 /** drops bound change events for variables in knapsack */
522 static
524  SCIP* scip, /**< SCIP data structure */
525  SCIP_CONSDATA* consdata, /**< constraint data */
526  SCIP_EVENTHDLR* eventhdlr /**< event handler to call for the event processing */
527  )
528 {
529  int i;
530 
531  assert(consdata != NULL);
532  assert(consdata->nvars == 0 || consdata->vars != NULL);
533  assert(consdata->nvars == 0 || consdata->weights != NULL);
534  assert(consdata->nvars == 0 || consdata->eventdata != NULL);
535 
536  for( i = 0; i < consdata->nvars; i++)
537  {
538  SCIP_CALL( SCIPdropVarEvent(scip, consdata->vars[i], EVENTTYPE_KNAPSACK,
539  eventhdlr, consdata->eventdata[i], consdata->eventdata[i]->filterpos) );
540  SCIP_CALL( eventdataFree(scip, &consdata->eventdata[i]) );
541  }
542 
543  return SCIP_OKAY;
544 }
545 
546 /** ensures, that vars and vals arrays can store at least num entries */
547 static
549  SCIP* scip, /**< SCIP data structure */
550  SCIP_CONSDATA* consdata, /**< knapsack constraint data */
551  int num, /**< minimum number of entries to store */
552  SCIP_Bool transformed /**< is constraint from transformed problem? */
553  )
554 {
555  assert(consdata != NULL);
556  assert(consdata->nvars <= consdata->varssize);
557 
558  if( num > consdata->varssize )
559  {
560  int newsize;
561 
562  newsize = SCIPcalcMemGrowSize(scip, num);
563  SCIP_CALL( SCIPreallocBlockMemoryArray(scip, &consdata->vars, consdata->varssize, newsize) );
564  SCIP_CALL( SCIPreallocBlockMemoryArray(scip, &consdata->weights, consdata->varssize, newsize) );
565  if( transformed )
566  {
567  SCIP_CALL( SCIPreallocBlockMemoryArray(scip, &consdata->eventdata, consdata->varssize, newsize) );
568  SCIP_CALL( SCIPreallocBlockMemoryArray(scip, &consdata->cliquepartition, consdata->varssize, newsize) );
569  SCIP_CALL( SCIPreallocBlockMemoryArray(scip, &consdata->negcliquepartition, consdata->varssize, newsize) );
570  }
571  else
572  {
573  assert(consdata->eventdata == NULL);
574  assert(consdata->cliquepartition == NULL);
575  assert(consdata->negcliquepartition == NULL);
576  }
577  consdata->varssize = newsize;
578  }
579  assert(num <= consdata->varssize);
580 
581  return SCIP_OKAY;
582 }
583 
584 /** updates all weight sums for fixed and unfixed variables */
585 static
586 void updateWeightSums(
587  SCIP_CONSDATA* consdata, /**< knapsack constraint data */
588  SCIP_VAR* var, /**< variable for this weight */
589  SCIP_Longint weightdelta /**< difference between the old and the new weight of the variable */
590  )
591 {
592  assert(consdata != NULL);
593  assert(var != NULL);
595  consdata->weightsum += weightdelta;
596 
597  if( SCIPvarGetLbLocal(var) > 0.5 )
598  consdata->onesweightsum += weightdelta;
599 
600  assert(consdata->weightsum >= 0);
601  assert(consdata->onesweightsum >= 0);
602 }
603 
604 /** creates knapsack constraint data */
605 static
607  SCIP* scip, /**< SCIP data structure */
608  SCIP_CONSDATA** consdata, /**< pointer to store constraint data */
609  int nvars, /**< number of variables in knapsack */
610  SCIP_VAR** vars, /**< variables of knapsack */
611  SCIP_Longint* weights, /**< weights of knapsack items */
612  SCIP_Longint capacity /**< capacity of knapsack */
613  )
614 {
615  int v;
616  SCIP_Longint constant;
617 
618  assert(consdata != NULL);
619 
620  SCIP_CALL( SCIPallocBlockMemory(scip, consdata) );
621 
622  constant = 0L;
623  (*consdata)->vars = NULL;
624  (*consdata)->weights = NULL;
625  (*consdata)->nvars = 0;
626  if( nvars > 0 )
627  {
628  SCIP_VAR** varsbuffer;
629  SCIP_Longint* weightsbuffer;
630  int k;
631 
632  SCIP_CALL( SCIPallocBufferArray(scip, &varsbuffer, nvars) );
633  SCIP_CALL( SCIPallocBufferArray(scip, &weightsbuffer, nvars) );
634 
635  k = 0;
636  for( v = 0; v < nvars; ++v )
637  {
638  assert(vars[v] != NULL);
639  assert(SCIPvarIsBinary(vars[v]));
640 
641  /* all weight have to be non negative */
642  assert( weights[v] >= 0 );
643 
644  if( weights[v] > 0 )
645  {
646  /* treat fixed variables as constants if problem compression is enabled */
647  if( SCIPisConsCompressionEnabled(scip) && SCIPvarGetLbGlobal(vars[v]) > SCIPvarGetUbGlobal(vars[v]) - 0.5 )
648  {
649  /* only if the variable is fixed to 1, we add its weight to the constant */
650  if( SCIPvarGetUbGlobal(vars[v]) > 0.5 )
651  constant += weights[v];
652  }
653  else
654  {
655  varsbuffer[k] = vars[v];
656  weightsbuffer[k] = weights[v];
657  ++k;
658  }
659  }
660  }
661  assert(k >= 0);
662 
663  (*consdata)->nvars = k;
664 
665  /* copy the active variables and weights into the constraint data structure */
666  if( k > 0 )
667  {
668  SCIP_CALL( SCIPduplicateBlockMemoryArray(scip, &(*consdata)->vars, varsbuffer, k) );
669  SCIP_CALL( SCIPduplicateBlockMemoryArray(scip, &(*consdata)->weights, weightsbuffer, k) );
670  }
671 
672  /* free buffer storage */
673  SCIPfreeBufferArray(scip, &weightsbuffer);
674  SCIPfreeBufferArray(scip, &varsbuffer);
675  }
676 
677  /* capacity has to be greater or equal to zero */
678  assert(capacity >= 0);
679  assert(constant >= 0);
680 
681  (*consdata)->varssize = (*consdata)->nvars;
682  (*consdata)->capacity = capacity - constant;
683  (*consdata)->eventdata = NULL;
684  (*consdata)->cliquepartition = NULL;
685  (*consdata)->negcliquepartition = NULL;
686  (*consdata)->row = NULL;
687  (*consdata)->weightsum = 0;
688  (*consdata)->onesweightsum = 0;
689  (*consdata)->ncliques = 0;
690  (*consdata)->nnegcliques = 0;
691  (*consdata)->presolvedtiming = 0;
692  (*consdata)->sorted = FALSE;
693  (*consdata)->cliquepartitioned = FALSE;
694  (*consdata)->negcliquepartitioned = FALSE;
695  (*consdata)->ncliqueslastpart = -1;
696  (*consdata)->ncliqueslastnegpart = -1;
697  (*consdata)->merged = FALSE;
698  (*consdata)->cliquesadded = FALSE;
699  (*consdata)->varsdeleted = FALSE;
700  (*consdata)->existmultaggr = FALSE;
701 
702  /* get transformed variables, if we are in the transformed problem */
703  if( SCIPisTransformed(scip) )
704  {
705  SCIP_CALL( SCIPgetTransformedVars(scip, (*consdata)->nvars, (*consdata)->vars, (*consdata)->vars) );
706 
707  for( v = 0; v < (*consdata)->nvars; v++ )
708  {
709  SCIP_VAR* var = SCIPvarGetProbvar((*consdata)->vars[v]);
710  assert(var != NULL);
711  (*consdata)->existmultaggr = (*consdata)->existmultaggr || (SCIPvarGetStatus(var) == SCIP_VARSTATUS_MULTAGGR);
712  }
713 
714  /* allocate memory for additional data structures */
715  SCIP_CALL( SCIPallocBlockMemoryArray(scip, &(*consdata)->eventdata, (*consdata)->nvars) );
716  SCIP_CALL( SCIPallocBlockMemoryArray(scip, &(*consdata)->cliquepartition, (*consdata)->nvars) );
717  SCIP_CALL( SCIPallocBlockMemoryArray(scip, &(*consdata)->negcliquepartition, (*consdata)->nvars) );
718  }
719 
720  /* calculate sum of weights and capture variables */
721  for( v = 0; v < (*consdata)->nvars; ++v )
722  {
723  /* calculate sum of weights */
724  updateWeightSums(*consdata, (*consdata)->vars[v], (*consdata)->weights[v]);
725 
726  /* capture variables */
727  SCIP_CALL( SCIPcaptureVar(scip, (*consdata)->vars[v]) );
728  }
729  return SCIP_OKAY;
730 }
731 
732 /** frees knapsack constraint data */
733 static
735  SCIP* scip, /**< SCIP data structure */
736  SCIP_CONSDATA** consdata, /**< pointer to the constraint data */
737  SCIP_EVENTHDLR* eventhdlr /**< event handler to call for the event processing */
738  )
739 {
740  assert(consdata != NULL);
741  assert(*consdata != NULL);
743  if( (*consdata)->row != NULL )
744  {
745  SCIP_CALL( SCIPreleaseRow(scip, &(*consdata)->row) );
746  }
747  if( (*consdata)->eventdata != NULL )
748  {
749  SCIP_CALL( dropEvents(scip, *consdata, eventhdlr) );
750  SCIPfreeBlockMemoryArray(scip, &(*consdata)->eventdata, (*consdata)->varssize);
751  }
752  if( (*consdata)->negcliquepartition != NULL )
753  {
754  SCIPfreeBlockMemoryArray(scip, &(*consdata)->negcliquepartition, (*consdata)->varssize);
755  }
756  if( (*consdata)->cliquepartition != NULL )
757  {
758  SCIPfreeBlockMemoryArray(scip, &(*consdata)->cliquepartition, (*consdata)->varssize);
759  }
760  if( (*consdata)->vars != NULL )
761  {
762  int v;
763 
764  /* release variables */
765  for( v = 0; v < (*consdata)->nvars; v++ )
766  {
767  assert((*consdata)->vars[v] != NULL);
768  SCIP_CALL( SCIPreleaseVar(scip, &((*consdata)->vars[v])) );
769  }
770 
771  assert( (*consdata)->weights != NULL );
772  assert( (*consdata)->varssize > 0 );
773  SCIPfreeBlockMemoryArray(scip, &(*consdata)->vars, (*consdata)->varssize);
774  SCIPfreeBlockMemoryArray(scip, &(*consdata)->weights, (*consdata)->varssize);
775  }
776 
777  SCIPfreeBlockMemory(scip, consdata);
778 
779  return SCIP_OKAY;
780 }
781 
782 /** changes a single weight in knapsack constraint data */
783 static
784 void consdataChgWeight(
785  SCIP_CONSDATA* consdata, /**< knapsack constraint data */
786  int item, /**< item number */
787  SCIP_Longint newweight /**< new weight of item */
788  )
789 {
790  SCIP_Longint oldweight;
791  SCIP_Longint weightdiff;
793  assert(consdata != NULL);
794  assert(0 <= item && item < consdata->nvars);
795 
796  oldweight = consdata->weights[item];
797  weightdiff = newweight - oldweight;
798  consdata->weights[item] = newweight;
799 
800 
801  /* update weight sums for all and fixed variables */
802  updateWeightSums(consdata, consdata->vars[item], weightdiff);
803 
804  if( consdata->eventdata != NULL )
805  {
806  assert(consdata->eventdata[item] != NULL);
807  assert(consdata->eventdata[item]->weight == oldweight);
808  consdata->eventdata[item]->weight = newweight;
809  }
810 
811  consdata->presolvedtiming = 0;
812  consdata->sorted = FALSE;
813 
814  /* recalculate cliques extraction after a weight was increased */
815  if( oldweight < newweight )
816  {
817  consdata->cliquesadded = FALSE;
818  }
819 }
820 
821 /** creates LP row corresponding to knapsack constraint */
822 static
824  SCIP* scip, /**< SCIP data structure */
825  SCIP_CONS* cons /**< knapsack constraint */
826  )
827 {
828  SCIP_CONSDATA* consdata;
829  int i;
830 
831  consdata = SCIPconsGetData(cons);
832  assert(consdata != NULL);
833  assert(consdata->row == NULL);
834 
835  SCIP_CALL( SCIPcreateEmptyRowCons(scip, &consdata->row, SCIPconsGetHdlr(cons), SCIPconsGetName(cons),
836  -SCIPinfinity(scip), (SCIP_Real)consdata->capacity,
838 
839  SCIP_CALL( SCIPcacheRowExtensions(scip, consdata->row) );
840  for( i = 0; i < consdata->nvars; ++i )
841  {
842  SCIP_CALL( SCIPaddVarToRow(scip, consdata->row, consdata->vars[i], (SCIP_Real)consdata->weights[i]) );
843  }
844  SCIP_CALL( SCIPflushRowExtensions(scip, consdata->row) );
845 
846  return SCIP_OKAY;
847 }
848 
849 /** adds linear relaxation of knapsack constraint to the LP */
850 static
852  SCIP* scip, /**< SCIP data structure */
853  SCIP_CONS* cons, /**< knapsack constraint */
854  SCIP_SOL* sol, /**< primal CIP solution, NULL for current LP solution */
855  SCIP_Bool* cutoff /**< whether a cutoff has been detected */
856  )
857 {
858  SCIP_CONSDATA* consdata;
860  assert( cutoff != NULL );
861  *cutoff = FALSE;
862 
863  consdata = SCIPconsGetData(cons);
864  assert(consdata != NULL);
865 
866  if( consdata->row == NULL )
867  {
868  SCIP_CALL( createRelaxation(scip, cons) );
869  }
870  assert(consdata->row != NULL);
871 
872  /* insert LP row as cut */
873  if( !SCIProwIsInLP(consdata->row) )
874  {
875  SCIPdebugMsg(scip, "adding relaxation of knapsack constraint <%s> (capacity %" SCIP_LONGINT_FORMAT "): ",
876  SCIPconsGetName(cons), consdata->capacity);
877  SCIPdebug( SCIP_CALL(SCIPprintRow(scip, consdata->row, NULL)) );
878  SCIP_CALL( SCIPaddCut(scip, sol, consdata->row, FALSE, cutoff) );
879  }
880 
881  return SCIP_OKAY;
882 }
883 
884 /** checks knapsack constraint for feasibility of given solution: returns TRUE iff constraint is feasible */
885 static
887  SCIP* scip, /**< SCIP data structure */
888  SCIP_CONS* cons, /**< constraint to check */
889  SCIP_SOL* sol, /**< solution to check, NULL for current solution */
890  SCIP_Bool checklprows, /**< Do constraints represented by rows in the current LP have to be checked? */
891  SCIP_Bool printreason, /**< Should the reason for the violation be printed? */
892  SCIP_Bool* violated /**< pointer to store whether the constraint is violated */
893  )
894 {
895  SCIP_CONSDATA* consdata;
896 
897  assert(violated != NULL);
898 
899  consdata = SCIPconsGetData(cons);
900  assert(consdata != NULL);
901 
902  SCIPdebugMsg(scip, "checking knapsack constraint <%s> for feasibility of solution %p (lprows=%u)\n",
903  SCIPconsGetName(cons), (void*)sol, checklprows);
904 
905  *violated = FALSE;
906 
907  if( checklprows || consdata->row == NULL || !SCIProwIsInLP(consdata->row) )
908  {
909  SCIP_Real sum;
910  SCIP_Longint integralsum;
911  SCIP_Bool ishuge;
912  int v;
913 
914  /* increase age of constraint; age is reset to zero, if a violation was found only in case we are in
915  * enforcement
916  */
917  if( sol == NULL )
918  {
919  SCIP_CALL( SCIPincConsAge(scip, cons) );
920  }
921 
922  sum = 0.0;
923  integralsum = 0;
924  /* we perform an more exact comparison if the capacity does not exceed the huge value */
925  if( SCIPisHugeValue(scip, (SCIP_Real) consdata->capacity) )
926  {
927  ishuge = TRUE;
928 
929  /* sum over all weight times the corresponding solution value */
930  for( v = consdata->nvars - 1; v >= 0; --v )
931  {
932  assert(SCIPvarIsBinary(consdata->vars[v]));
933  sum += consdata->weights[v] * SCIPgetSolVal(scip, sol, consdata->vars[v]);
934  }
935  }
936  else
937  {
938  ishuge = FALSE;
939 
940  /* sum over all weight for which the variable has a solution value of 1 in feastol */
941  for( v = consdata->nvars - 1; v >= 0; --v )
942  {
943  assert(SCIPvarIsBinary(consdata->vars[v]));
944 
945  if( SCIPgetSolVal(scip, sol, consdata->vars[v]) > 0.5 )
946  integralsum += consdata->weights[v];
947  }
948  }
949 
950  if( (!ishuge && integralsum > consdata->capacity) || (ishuge && SCIPisFeasGT(scip, sum, (SCIP_Real)consdata->capacity)) )
951  {
952  *violated = TRUE;
953 
954  /* only reset constraint age if we are in enforcement */
955  if( sol == NULL )
956  {
957  SCIP_CALL( SCIPresetConsAge(scip, cons) );
958  }
959 
960  if( printreason )
961  {
962  SCIP_Real viol = ishuge ? sum : (SCIP_Real)integralsum;
963 
964  viol -= consdata->capacity;
965  assert(viol > 0);
966 
967  SCIP_CALL( SCIPprintCons(scip, cons, NULL) );
968 
969  SCIPinfoMessage(scip, NULL, ";\n");
970  SCIPinfoMessage(scip, NULL, "violation: the capacity is violated by %.15g\n", viol);
971  }
972  }
973  }
974 
975  return SCIP_OKAY;
976 }
977 
978 /* IDX computes the integer index for the optimal solution array */
979 #define IDX(j,d) ((j)*(intcap)+(d))
980 
981 /** solves knapsack problem in maximization form exactly using dynamic programming;
982  * if needed, one can provide arrays to store all selected items and all not selected items
983  *
984  * @note in case you provide the solitems or nonsolitems array you also have to provide the counter part as well
985  */
987  SCIP* scip, /**< SCIP data structure */
988  int nitems, /**< number of available items */
989  SCIP_Longint* weights, /**< item weights */
990  SCIP_Real* profits, /**< item profits */
991  SCIP_Longint capacity, /**< capacity of knapsack */
992  int* items, /**< item numbers */
993  int* solitems, /**< array to store items in solution, or NULL */
994  int* nonsolitems, /**< array to store items not in solution, or NULL */
995  int* nsolitems, /**< pointer to store number of items in solution, or NULL */
996  int* nnonsolitems, /**< pointer to store number of items not in solution, or NULL */
997  SCIP_Real* solval, /**< pointer to store optimal solution value, or NULL */
998  SCIP_Bool* success /**< pointer to store if an error occured during solving
999  * (normally a memory problem) */
1000  )
1001 {
1002  SCIP_RETCODE retcode;
1003  SCIP_Real* tempsort;
1004  SCIP_Real* optvalues;
1005  int intcap;
1006  int d;
1007  int j;
1008  SCIP_Longint weightsum;
1009  int* myitems;
1010  SCIP_Longint* myweights;
1011  int* allcurrminweight;
1012  SCIP_Real* myprofits;
1013  int nmyitems;
1014  SCIP_Longint gcd;
1015  SCIP_Longint minweight;
1016  SCIP_Longint maxweight;
1017  int currminweight;
1018  SCIP_Longint greedycap;
1019  SCIP_Longint greedysolweight;
1020  SCIP_Real greedysolvalue;
1021  SCIP_Bool eqweights;
1022  SCIP_Bool isoptimal;
1023  const size_t maxsize_t = (size_t)(-1);
1024 
1025  assert(weights != NULL);
1026  assert(profits != NULL);
1027  assert(capacity >= 0);
1028  assert(items != NULL);
1029  assert(nitems >= 0);
1030  assert(success != NULL);
1031 
1032  *success = TRUE;
1033 
1034 #ifndef NDEBUG
1035  for( j = nitems - 1; j >= 0; --j )
1036  assert(weights[j] >= 0);
1037 #endif
1038 
1039  SCIPdebugMsg(scip, "Solving knapsack exactly.\n");
1040 
1041  /* initializing solution value */
1042  if( solval != NULL )
1043  *solval = 0.0;
1044 
1045  /* produces optimal solution by following the table */
1046  if( solitems != NULL)
1047  {
1048  assert(items != NULL);
1049  assert(nsolitems != NULL);
1050  assert(nonsolitems != NULL);
1051  assert(nnonsolitems != NULL);
1052 
1053  *nnonsolitems = 0;
1054  *nsolitems = 0;
1055  }
1056 
1057  /* allocate temporary memory */
1058  SCIP_CALL( SCIPallocBufferArray(scip, &myweights, nitems) );
1059  SCIP_CALL( SCIPallocBufferArray(scip, &myprofits, nitems) );
1060  SCIP_CALL( SCIPallocBufferArray(scip, &myitems, nitems) );
1061  nmyitems = 0;
1062  weightsum = 0;
1063  minweight = SCIP_LONGINT_MAX;
1064  maxweight = 0;
1065 
1066  /* remove unnecessary items */
1067  for( j = 0; j < nitems; ++j )
1068  {
1069  assert(0 <= weights[j] && weights[j] < SCIP_LONGINT_MAX);
1070  /* items does not fit */
1071  if( weights[j] > capacity )
1072  {
1073  if( solitems != NULL)
1074  {
1075  nonsolitems[*nnonsolitems] = items[j];
1076  ++(*nnonsolitems);
1077  }
1078  }
1079  /* items we does not want */
1080  else if( profits[j] <= 0.0 )
1081  {
1082  if( solitems != NULL)
1083  {
1084  nonsolitems[*nnonsolitems] = items[j];
1085  ++(*nnonsolitems);
1086  }
1087  }
1088  /* items which always fit */
1089  else if( weights[j] == 0 )
1090  {
1091  if( solitems != NULL)
1092  {
1093  solitems[*nsolitems] = items[j];
1094  ++(*nsolitems);
1095  }
1096  if( solval != NULL )
1097  *solval += profits[j];
1098  }
1099  /* all important items */
1100  else
1101  {
1102  myweights[nmyitems] = weights[j];
1103  myprofits[nmyitems] = profits[j];
1104  myitems[nmyitems] = items[j];
1105 
1106  /* remember smallest item */
1107  if( myweights[nmyitems] < minweight )
1108  minweight = myweights[nmyitems];
1109 
1110  /* remember bigest item */
1111  if( myweights[nmyitems] > maxweight )
1112  maxweight = myweights[nmyitems];
1113 
1114  weightsum += myweights[nmyitems];
1115  ++nmyitems;
1116  }
1117  }
1118 
1119  /* no item is left then goto end */
1120  if( nmyitems == 0 )
1121  {
1122  SCIPdebugMsg(scip, "After preprocessing no items are left.\n");
1123 
1124  goto TERMINATE;
1125  }
1126  /* if all items fit, we also do not need to do the expensive stuff later on */
1127  else if( weightsum > 0 && weightsum <= capacity )
1128  {
1129  SCIPdebugMsg(scip, "After preprocessing all items fit into knapsack.\n");
1130 
1131  for( j = nmyitems - 1; j >= 0; --j )
1132  {
1133  if( solitems != NULL )
1134  {
1135  solitems[*nsolitems] = myitems[j];
1136  ++(*nsolitems);
1137  }
1138  if( solval != NULL )
1139  *solval += myprofits[j];
1140  }
1141 
1142  goto TERMINATE;
1143  }
1144 
1145  assert(minweight > 0);
1146  assert(maxweight > 0);
1147 
1148  if( maxweight > 1 )
1149  {
1150  /* determine greatest common divisor */
1151  gcd = myweights[nmyitems - 1];
1152  for( j = nmyitems - 2; j >= 0 && gcd >= 2; --j )
1153  gcd = SCIPcalcGreComDiv(gcd, myweights[j]);
1154 
1155  SCIPdebugMsg(scip, "Gcd is %" SCIP_LONGINT_FORMAT ".\n", gcd);
1156 
1157  /* divide by greatest common divisor */
1158  if( gcd > 1 )
1159  {
1160  eqweights = TRUE;
1161  for( j = nmyitems - 1; j >= 0; --j )
1162  {
1163  myweights[j] /= gcd;
1164  eqweights = eqweights && (myweights[j] == 1);
1165  }
1166  capacity /= gcd;
1167  minweight /= gcd;
1168  }
1169  else
1170  eqweights = FALSE;
1171  }
1172  else
1173  {
1174  assert(maxweight == 1);
1175  eqweights = TRUE;
1176  }
1177 
1178  assert(minweight <= capacity);
1179 
1180  /* only one item fits, than take the best */
1181  if( minweight > capacity / 2 )
1182  {
1183  int p;
1184 
1185  SCIPdebugMsg(scip, "Only one item fits into knapsack, so take the best.\n");
1186 
1187  p = nmyitems - 1;
1188 
1189  /* find best item */
1190  for( j = nmyitems - 2; j >= 0; --j )
1191  if( myprofits[j] > myprofits[p] )
1192  p = j;
1193 
1194  /* update solution information */
1195  if( solitems != NULL)
1196  {
1197  solitems[*nsolitems] = myitems[p];
1198  ++(*nsolitems);
1199  for( j = nmyitems - 1; j >= 0; --j )
1200  if( j != p )
1201  {
1202  nonsolitems[*nnonsolitems] = myitems[j];
1203  ++(*nnonsolitems);
1204  }
1205  }
1206  /* update solution value */
1207  if( solval != NULL )
1208  *solval += myprofits[p];
1209 
1210  goto TERMINATE;
1211  }
1212 
1213  /* all items have the same weight, than take the best */
1214  if( eqweights )
1215  {
1216  SCIP_Real addval;
1217 
1218  SCIPdebugMsg(scip, "All weights are equal, so take the best.\n");
1219 
1220  SCIPsortDownRealIntLong(myprofits, myitems, myweights, nmyitems);
1221 
1222  addval = 0.0;
1223  /* update solution information */
1224  if( solitems != NULL || solval != NULL )
1225  {
1226  SCIP_Longint i;
1227 
1228  /* if all items would fit we had handled this case before */
1229  assert((SCIP_Longint) nmyitems > capacity);
1230 
1231  /* take the first best items into the solution */
1232  for( i = capacity - 1; i >= 0; --i )
1233  {
1234  if( solitems != NULL)
1235  {
1236  assert(nonsolitems != NULL);
1237  solitems[*nsolitems] = myitems[i];
1238  ++(*nsolitems);
1239  }
1240  addval += myprofits[i];
1241  }
1242 
1243  if( solitems != NULL)
1244  {
1245  assert(nonsolitems != NULL);
1246 
1247  /* the rest are not in the solution */
1248  for( i = nmyitems - 1; i >= capacity; --i )
1249  {
1250  nonsolitems[*nnonsolitems] = myitems[i];
1251  ++(*nnonsolitems);
1252  }
1253  }
1254  }
1255  /* update solution value */
1256  if( solval != NULL )
1257  {
1258  assert(addval > 0.0);
1259  *solval += addval;
1260  }
1261 
1262  goto TERMINATE;
1263  }
1264 
1265  /* in the following table we do not need the first minweight columns */
1266  capacity -= (minweight - 1);
1267 
1268  /* we can only handle integers */
1269  if( capacity >= INT_MAX )
1270  {
1271  SCIPdebugMsg(scip, "Capacity is to big, so we cannot handle it here.\n");
1272 
1273  *success = FALSE;
1274  goto TERMINATE;
1275  }
1276  assert(capacity < INT_MAX);
1277 
1278  intcap = (int)capacity;
1279  assert(intcap >= 0);
1280  assert(nmyitems > 0);
1281  assert(sizeof(size_t) >= sizeof(int)); /* no following conversion should be messed up */
1282 
1283  /* this condition checks if we will try to allocate a correct number of bytes and do not have an overflow, while
1284  * computing the size for the allocation
1285  */
1286  if( intcap < 0 || (intcap > 0 && (((size_t)nmyitems) > (maxsize_t / (size_t)intcap / sizeof(*optvalues)) || ((size_t)nmyitems) * ((size_t)intcap) * sizeof(*optvalues) > ((size_t)INT_MAX) )) ) /*lint !e571*/
1287  {
1288  SCIPdebugMsg(scip, "Too much memory (%lu) would be consumed.\n", (unsigned long) (((size_t)nmyitems) * ((size_t)intcap) * sizeof(*optvalues))); /*lint !e571*/
1289 
1290  *success = FALSE;
1291  goto TERMINATE;
1292  }
1293 
1294  /* allocate temporary memory and check for memory exceeding */
1295  retcode = SCIPallocBufferArray(scip, &optvalues, nmyitems * intcap);
1296  if( retcode == SCIP_NOMEMORY )
1297  {
1298  SCIPdebugMsg(scip, "Did not get enough memory.\n");
1299 
1300  *success = FALSE;
1301  goto TERMINATE;
1302  }
1303  else
1304  {
1305  SCIP_CALL( retcode );
1306  }
1307 
1308  /* sort myitems (plus corresponding arrays myweights and myprofits) such that
1309  * p_1/w_1 >= p_2/w_2 >= ... >= p_n/w_n, this is only use for greedy solution
1310  */
1311  SCIP_CALL( SCIPallocBufferArray(scip, &tempsort, nmyitems) );
1312  for( j = nmyitems - 1; j >= 0; --j )
1313  tempsort[j] = myprofits[j]/((SCIP_Real) myweights[j]);
1314 
1315  SCIPsortDownRealLongRealInt(tempsort, myweights, myprofits, myitems, nmyitems);
1316 
1317  /* initialize values for greedy solution information */
1318  greedysolweight = 0;
1319  greedysolvalue = 0.0;
1320  isoptimal = TRUE;
1321  greedycap = capacity + (minweight - 1);
1322 
1323  SCIPdebugMsg(scip, "Determine greedy solution.\n");
1324 
1325  /* determine greedy solution */
1326  for( j = 0; j < nmyitems; ++j )
1327  {
1328  assert(myweights[j] <= greedycap);
1329 
1330  /* take all fitting items */
1331  if( myweights[j] + greedysolweight <= greedycap )
1332  {
1333  /* update greedy solution weight and value */
1334  greedysolweight += myweights[j];
1335  greedysolvalue += myprofits[j];
1336  continue;
1337  }
1338  else if( greedysolweight < greedycap )
1339  isoptimal = FALSE;
1340  break;
1341  }
1342  assert(greedysolweight > 0);
1343  assert(greedysolvalue > 0.0);
1344 
1345  /* greedy solution is optimal */
1346  if( isoptimal )
1347  {
1348  assert(greedysolweight == greedycap);
1349 
1350  SCIPdebugMsg(scip, "Greedy solution is optimal.\n");
1351 
1352  greedysolweight = 0;
1353 
1354  /* update solution information */
1355  if( solitems != NULL)
1356  {
1357  /* take the first best items into the solution */
1358  for( j = 0; j < nmyitems; ++j )
1359  {
1360  /* take all fitting items */
1361  if( myweights[j] + greedysolweight <= greedycap )
1362  {
1363  solitems[*nsolitems] = myitems[j];
1364  ++(*nsolitems);
1365  greedysolweight += myweights[j];
1366  }
1367  else
1368  {
1369  nonsolitems[*nnonsolitems] = myitems[j];
1370  ++(*nnonsolitems);
1371  }
1372  }
1373  }
1374  /* update solution value */
1375  if( solval != NULL )
1376  {
1377  assert(greedysolvalue > 0.0);
1378  *solval += greedysolvalue;
1379  }
1380 
1381  SCIPfreeBufferArray(scip, &tempsort);
1382  SCIPfreeBufferArray(scip, &optvalues);
1383 
1384  goto TERMINATE;
1385  }
1386 
1387  SCIPdebugMsg(scip, "Start real exact algorithm.\n");
1388 
1389  /* we memorize at each step the current minimal weight to later on know which value in our optvalues matrix is valid;
1390  * all values entries of the j-th row of optvalues is valid if the index is >= allcurrminweight[j], otherwise it is
1391  * invalid, a second possibility would be to clear the whole optvalues, which should be more expensive than storing
1392  * 'nmyitem' values
1393  */
1394  SCIP_CALL( SCIPallocBufferArray(scip, &allcurrminweight, nmyitems) );
1395  assert(myweights[0] - minweight < INT_MAX);
1396  currminweight = (int) (myweights[0] - minweight);
1397  allcurrminweight[0] = currminweight;
1398 
1399  /* fills first row of dynamic programming table with optimal values */
1400  for( d = currminweight; d < intcap; ++d )
1401  optvalues[d] = myprofits[0];
1402  /* fills dynamic programming table with optimal values */
1403  for( j = 1; j < nmyitems; ++j )
1404  {
1405  int intweight;
1406 
1407  /* compute important part of weight, which will be represented in the table */
1408  intweight = (int)(myweights[j] - minweight);
1409  assert(0 <= intweight && intweight < intcap);
1410 
1411  /* copy all nonzeros from row above */
1412  for( d = currminweight; d < intweight && d < intcap; ++d )
1413  optvalues[IDX(j,d)] = optvalues[IDX(j-1,d)];
1414 
1415  /* update corresponding row */
1416  for( d = intweight; d < intcap; ++d )
1417  {
1418  /* if index d is smaller the the current minweight then optvalues[IDX(j-1,d)] is not initialized, i.e. should
1419  * be 0
1420  */
1421  if( d < currminweight )
1422  {
1423  optvalues[IDX(j,d)] = myprofits[j];
1424  }
1425  else
1426  {
1427  SCIP_Real sumprofit;
1428 
1429  if( d - myweights[j] < currminweight )
1430  sumprofit = myprofits[j];
1431  else
1432  sumprofit = optvalues[IDX(j-1,(int)(d-myweights[j]))] + myprofits[j];
1433 
1434  optvalues[IDX(j,d)] = MAX(sumprofit, optvalues[IDX(j-1,d)]);
1435  }
1436  }
1437  /* update currminweight */
1438  if( intweight < currminweight )
1439  currminweight = intweight;
1440 
1441  allcurrminweight[j] = currminweight;
1442  }
1443 
1444  /* update optimal solution by following the table */
1445  if( solitems != NULL)
1446  {
1447  d = intcap - 1;
1448 
1449  SCIPdebugMsg(scip, "Fill the solution vector after solving exactly.\n");
1450 
1451  /* insert all items in (non-) solution vector */
1452  for( j = nmyitems - 1; j > 0; --j )
1453  {
1454  /* if we cannot find any item anymore which is in our solution stop, if the following condition holds this
1455  * means all remaining items does not fit anymore
1456  */
1457  if( d < allcurrminweight[j] )
1458  {
1459  /* we cannot have exceeded our capacity */
1460  assert((SCIP_Longint) d >= -minweight);
1461  break;
1462  }
1463  /* collect solution items, first condition means that no next item can fit anymore, but this does */
1464  if( d < allcurrminweight[j-1] || optvalues[IDX(j,d)] > optvalues[IDX(j-1,d)] )
1465  {
1466  solitems[*nsolitems] = myitems[j];
1467  ++(*nsolitems);
1468 
1469  /* check that we do not have an underflow */
1470  assert(myweights[j] <= (INT_MAX + (SCIP_Longint) d));
1471  d = (int)(d - myweights[j]);
1472  }
1473  /* collect non-solution items */
1474  else
1475  {
1476  nonsolitems[*nnonsolitems] = myitems[j];
1477  ++(*nnonsolitems);
1478  }
1479  }
1480 
1481  /* insert remaining items */
1482  if( d >= allcurrminweight[j] )
1483  {
1484  assert(j == 0);
1485  solitems[*nsolitems] = myitems[j];
1486  ++(*nsolitems);
1487  }
1488  else
1489  {
1490  assert(j >= 0);
1491  assert(d < allcurrminweight[j]);
1492 
1493  for( ; j >= 0; --j )
1494  {
1495  nonsolitems[*nnonsolitems] = myitems[j];
1496  ++(*nnonsolitems);
1497  }
1498  }
1499 
1500  assert(*nsolitems + *nnonsolitems == nitems);
1501  }
1502 
1503  /* update solution value */
1504  if( solval != NULL )
1505  *solval += optvalues[IDX(nmyitems-1,intcap-1)];
1506 
1507  SCIPfreeBufferArray(scip, &allcurrminweight);
1508 
1509  /* free all temporary memory */
1510  SCIPfreeBufferArray(scip, &tempsort);
1511  SCIPfreeBufferArray(scip, &optvalues);
1512 
1513  TERMINATE:
1514  SCIPfreeBufferArray(scip, &myitems);
1515  SCIPfreeBufferArray(scip, &myprofits);
1516  SCIPfreeBufferArray(scip, &myweights);
1517 
1518  return SCIP_OKAY;
1519 }
1520 
1521 /** solves knapsack problem in maximization form approximately by solving the LP-relaxation of the problem using Dantzig's
1522  * method and rounding down the solution; if needed, one can provide arrays to store all selected items and all not
1523  * selected items
1524  */
1526  SCIP* scip, /**< SCIP data structure */
1527  int nitems, /**< number of available items */
1528  SCIP_Longint* weights, /**< item weights */
1529  SCIP_Real* profits, /**< item profits */
1530  SCIP_Longint capacity, /**< capacity of knapsack */
1531  int* items, /**< item numbers */
1532  int* solitems, /**< array to store items in solution, or NULL */
1533  int* nonsolitems, /**< array to store items not in solution, or NULL */
1534  int* nsolitems, /**< pointer to store number of items in solution, or NULL */
1535  int* nnonsolitems, /**< pointer to store number of items not in solution, or NULL */
1536  SCIP_Real* solval /**< pointer to store optimal solution value, or NULL */
1537  )
1538 {
1539  SCIP_Real* tempsort;
1540  SCIP_Longint solitemsweight;
1541  SCIP_Real* realweights;
1542  int j;
1543  int criticalindex;
1544 
1545  assert(weights != NULL);
1546  assert(profits != NULL);
1547  assert(capacity >= 0);
1548  assert(items != NULL);
1549  assert(nitems >= 0);
1550 
1551  if( solitems != NULL )
1552  {
1553  *nsolitems = 0;
1554  *nnonsolitems = 0;
1555  }
1556  if( solval != NULL )
1557  *solval = 0.0;
1558 
1559  /* initialize data for median search */
1560  SCIP_CALL( SCIPallocBufferArray(scip, &tempsort, nitems) );
1561  SCIP_CALL( SCIPallocBufferArray(scip, &realweights, nitems) );
1562  for( j = nitems - 1; j >= 0; --j )
1563  {
1564  tempsort[j] = profits[j]/((SCIP_Real) weights[j]);
1565  realweights[j] = (SCIP_Real)weights[j];
1566 
1567  }
1568 
1569  /* partially sort indices such that all elements that are larger than the break item appear first */
1570  SCIPselectWeightedDownRealLongRealInt(tempsort, weights, profits, items, realweights, (SCIP_Real)capacity, nitems, &criticalindex);
1571 
1572  /* selects items as long as they fit into the knapsack */
1573  solitemsweight = 0;
1574  for( j = 0; j < nitems && solitemsweight + weights[j] <= capacity; ++j )
1575  {
1576  if( solitems != NULL )
1577  {
1578  solitems[*nsolitems] = items[j];
1579  (*nsolitems)++;
1580  }
1581  if( solval != NULL )
1582  (*solval) += profits[j];
1583  solitemsweight += weights[j];
1584  }
1585  for( ; j < nitems && solitems != NULL; j++ )
1586  {
1587  nonsolitems[*nnonsolitems] = items[j];
1588  (*nnonsolitems)++;
1589  }
1590 
1591  SCIPfreeBufferArray(scip, &realweights);
1592  SCIPfreeBufferArray(scip, &tempsort);
1593 
1594  return SCIP_OKAY;
1595 }
1596 
1597 #ifdef SCIP_DEBUG
1598 /** prints all nontrivial GUB constraints and their LP solution values */
1599 static
1600 void GUBsetPrint(
1601  SCIP* scip, /**< SCIP data structure */
1602  SCIP_GUBSET* gubset, /**< GUB set data structure */
1603  SCIP_VAR** vars, /**< variables in knapsack constraint */
1604  SCIP_Real* solvals /**< solution values of variables in knapsack constraint; or NULL */
1605  )
1606 {
1607  int nnontrivialgubconss;
1608  int c;
1609 
1610  nnontrivialgubconss = 0;
1611 
1612  SCIPdebugMsg(scip, " Nontrivial GUBs of current GUB set:\n");
1613 
1614  /* print out all nontrivial GUB constraints, i.e., with more than one variable */
1615  for( c = 0; c < gubset->ngubconss; c++ )
1616  {
1617  SCIP_Real gubsolval;
1618 
1619  assert(gubset->gubconss[c]->ngubvars >= 0);
1620 
1621  /* nontrivial GUB */
1622  if( gubset->gubconss[c]->ngubvars > 1 )
1623  {
1624  int v;
1625 
1626  gubsolval = 0.0;
1627  SCIPdebugMsg(scip, " GUB<%d>:\n", c);
1628 
1629  /* print GUB var */
1630  for( v = 0; v < gubset->gubconss[c]->ngubvars; v++ )
1631  {
1632  int currentvar;
1633 
1634  currentvar = gubset->gubconss[c]->gubvars[v];
1635  if( solvals != NULL )
1636  {
1637  gubsolval += solvals[currentvar];
1638  SCIPdebugMsg(scip, " +<%s>(%4.2f)\n", SCIPvarGetName(vars[currentvar]), solvals[currentvar]);
1639  }
1640  else
1641  {
1642  SCIPdebugMsg(scip, " +<%s>\n", SCIPvarGetName(vars[currentvar]));
1643  }
1644  }
1645 
1646  /* check whether LP solution satisfies the GUB constraint */
1647  if( solvals != NULL )
1648  {
1649  SCIPdebugMsg(scip, " =%4.2f <= 1 %s\n", gubsolval,
1650  SCIPisFeasGT(scip, gubsolval, 1.0) ? "--> violated" : "");
1651  }
1652  else
1653  {
1654  SCIPdebugMsg(scip, " <= 1 %s\n", SCIPisFeasGT(scip, gubsolval, 1.0) ? "--> violated" : "");
1655  }
1656  nnontrivialgubconss++;
1657  }
1658  }
1659 
1660  SCIPdebugMsg(scip, " --> %d/%d nontrivial GUBs\n", nnontrivialgubconss, gubset->ngubconss);
1661 }
1662 #endif
1663 
1664 /** creates an empty GUB constraint */
1665 static
1667  SCIP* scip, /**< SCIP data structure */
1668  SCIP_GUBCONS** gubcons /**< pointer to store GUB constraint data */
1669  )
1670 {
1671  assert(scip != NULL);
1672  assert(gubcons != NULL);
1673 
1674  /* allocate memory for GUB constraint data structures */
1675  SCIP_CALL( SCIPallocBuffer(scip, gubcons) );
1676  (*gubcons)->gubvarssize = GUBCONSGROWVALUE;
1677  SCIP_CALL( SCIPallocBufferArray(scip, &(*gubcons)->gubvars, (*gubcons)->gubvarssize) );
1678  SCIP_CALL( SCIPallocBufferArray(scip, &(*gubcons)->gubvarsstatus, (*gubcons)->gubvarssize) );
1679 
1680  (*gubcons)->ngubvars = 0;
1681 
1682  return SCIP_OKAY;
1683 }
1684 
1685 /** frees GUB constraint */
1686 static
1688  SCIP* scip, /**< SCIP data structure */
1689  SCIP_GUBCONS** gubcons /**< pointer to GUB constraint data structure */
1690  )
1691 {
1692  assert(scip != NULL);
1693  assert(gubcons != NULL);
1694  assert((*gubcons)->gubvars != NULL);
1695  assert((*gubcons)->gubvarsstatus != NULL);
1696 
1697  /* free allocated memory */
1698  SCIPfreeBufferArray(scip, &(*gubcons)->gubvarsstatus);
1699  SCIPfreeBufferArray(scip, &(*gubcons)->gubvars);
1700  SCIPfreeBuffer(scip, gubcons);
1701 
1702  return SCIP_OKAY;
1703 }
1704 
1705 /** adds variable to given GUB constraint */
1706 static
1708  SCIP* scip, /**< SCIP data structure */
1709  SCIP_GUBCONS* gubcons, /**< GUB constraint data */
1710  int var /**< index of given variable in knapsack constraint */
1711  )
1712 {
1713  assert(scip != NULL);
1714  assert(gubcons != NULL);
1715  assert(gubcons->ngubvars >= 0 && gubcons->ngubvars < gubcons->gubvarssize);
1716  assert(gubcons->gubvars != NULL);
1717  assert(gubcons->gubvarsstatus != NULL);
1718  assert(var >= 0);
1719 
1720  /* add variable to GUB constraint */
1721  gubcons->gubvars[gubcons->ngubvars] = var;
1722  gubcons->gubvarsstatus[gubcons->ngubvars] = GUBVARSTATUS_UNINITIAL;
1723  gubcons->ngubvars++;
1724 
1725  /* increase space allocated to GUB constraint if the number of variables reaches the size */
1726  if( gubcons->ngubvars == gubcons->gubvarssize )
1727  {
1728  int newlen;
1729 
1730  newlen = gubcons->gubvarssize + GUBCONSGROWVALUE;
1731  SCIP_CALL( SCIPreallocBufferArray(scip, &gubcons->gubvars, newlen) );
1732  SCIP_CALL( SCIPreallocBufferArray(scip, &gubcons->gubvarsstatus, newlen) );
1733 
1734  gubcons->gubvarssize = newlen;
1735  }
1736 
1737  return SCIP_OKAY;
1738 }
1739 
1740 /** deletes variable from its current GUB constraint */
1741 static
1743  SCIP* scip, /**< SCIP data structure */
1744  SCIP_GUBCONS* gubcons, /**< GUB constraint data */
1745  int var, /**< index of given variable in knapsack constraint */
1746  int gubvarsidx /**< index of the variable in its current GUB constraint */
1747  )
1748 {
1749  assert(scip != NULL);
1750  assert(gubcons != NULL);
1751  assert(var >= 0);
1752  assert(gubvarsidx >= 0 && gubvarsidx < gubcons->ngubvars);
1753  assert(gubcons->ngubvars >= gubvarsidx+1);
1754  assert(gubcons->gubvars[gubvarsidx] == var);
1755 
1756  /* delete variable from GUB by swapping it replacing in by the last variable in the GUB constraint */
1757  gubcons->gubvars[gubvarsidx] = gubcons->gubvars[gubcons->ngubvars-1];
1758  gubcons->gubvarsstatus[gubvarsidx] = gubcons->gubvarsstatus[gubcons->ngubvars-1];
1759  gubcons->ngubvars--;
1760 
1761  /* decrease space allocated for the GUB constraint, if the last GUBCONSGROWVALUE+1 array entries are now empty */
1762  if( gubcons->ngubvars < gubcons->gubvarssize - GUBCONSGROWVALUE && gubcons->ngubvars > 0 )
1763  {
1764  int newlen;
1765 
1766  newlen = gubcons->gubvarssize - GUBCONSGROWVALUE;
1767 
1768  SCIP_CALL( SCIPreallocBufferArray(scip, &gubcons->gubvars, newlen) );
1769  SCIP_CALL( SCIPreallocBufferArray(scip, &gubcons->gubvarsstatus, newlen) );
1770 
1771  gubcons->gubvarssize = newlen;
1772  }
1773 
1774  return SCIP_OKAY;
1775 }
1776 
1777 /** moves variable from current GUB constraint to a different existing (nonempty) GUB constraint */
1778 static
1780  SCIP* scip, /**< SCIP data structure */
1781  SCIP_GUBSET* gubset, /**< GUB set data structure */
1782  SCIP_VAR** vars, /**< variables in knapsack constraint */
1783  int var, /**< index of given variable in knapsack constraint */
1784  int oldgubcons, /**< index of old GUB constraint of given variable */
1785  int newgubcons /**< index of new GUB constraint of given variable */
1786  )
1788  int oldgubvaridx;
1789  int replacevar;
1790  int j;
1791 
1792  assert(scip != NULL);
1793  assert(gubset != NULL);
1794  assert(var >= 0);
1795  assert(oldgubcons >= 0 && oldgubcons < gubset->ngubconss);
1796  assert(newgubcons >= 0 && newgubcons < gubset->ngubconss);
1797  assert(oldgubcons != newgubcons);
1798  assert(gubset->gubconssidx[var] == oldgubcons);
1799  assert(gubset->gubconss[oldgubcons]->ngubvars > 0);
1800  assert(gubset->gubconss[newgubcons]->ngubvars >= 0);
1801 
1802  SCIPdebugMsg(scip, " moving variable<%s> from GUB<%d> to GUB<%d>\n", SCIPvarGetName(vars[var]), oldgubcons, newgubcons);
1803 
1804  oldgubvaridx = gubset->gubvarsidx[var];
1805 
1806  /* delete variable from old GUB constraint by replacing it by the last variable of the GUB constraint */
1807  SCIP_CALL( GUBconsDelVar(scip, gubset->gubconss[oldgubcons], var, oldgubvaridx) );
1808 
1809  /* in GUB set, update stored index of variable in old GUB constraint for the variable used for replacement;
1810  * replacement variable is given by old position of the deleted variable
1811  */
1812  replacevar = gubset->gubconss[oldgubcons]->gubvars[oldgubvaridx];
1813  assert(gubset->gubvarsidx[replacevar] == gubset->gubconss[oldgubcons]->ngubvars);
1814  gubset->gubvarsidx[replacevar] = oldgubvaridx;
1815 
1816  /* add variable to the end of new GUB constraint */
1817  SCIP_CALL( GUBconsAddVar(scip, gubset->gubconss[newgubcons], var) );
1818  assert(gubset->gubconss[newgubcons]->gubvars[gubset->gubconss[newgubcons]->ngubvars-1] == var);
1819 
1820  /* in GUB set, update stored index of GUB of moved variable and stored index of variable in this GUB constraint */
1821  gubset->gubconssidx[var] = newgubcons;
1822  gubset->gubvarsidx[var] = gubset->gubconss[newgubcons]->ngubvars-1;
1823 
1824  /* delete old GUB constraint if it became empty */
1825  if( gubset->gubconss[oldgubcons]->ngubvars == 0 )
1826  {
1827  SCIPdebugMsg(scip, "deleting empty GUB cons<%d> from current GUB set\n", oldgubcons);
1828 #ifdef SCIP_DEBUG
1829  GUBsetPrint(scip, gubset, vars, NULL);
1830 #endif
1831 
1832  /* free old GUB constraint */
1833  SCIP_CALL( GUBconsFree(scip, &gubset->gubconss[oldgubcons]) );
1834 
1835  /* if empty GUB was not the last one in GUB set data structure, replace it by last GUB constraint */
1836  if( oldgubcons != gubset->ngubconss-1 )
1837  {
1838  gubset->gubconss[oldgubcons] = gubset->gubconss[gubset->ngubconss-1];
1839  gubset->gubconsstatus[oldgubcons] = gubset->gubconsstatus[gubset->ngubconss-1];
1840 
1841  /* in GUB set, update stored index of GUB constraint for all variable of the GUB constraint used for replacement;
1842  * replacement GUB is given by old position of the deleted GUB
1843  */
1844  for( j = 0; j < gubset->gubconss[oldgubcons]->ngubvars; j++ )
1845  {
1846  assert(gubset->gubconssidx[gubset->gubconss[oldgubcons]->gubvars[j]] == gubset->ngubconss-1);
1847  gubset->gubconssidx[gubset->gubconss[oldgubcons]->gubvars[j]] = oldgubcons;
1848  }
1849  }
1850 
1851  /* update number of GUB constraints */
1852  gubset->ngubconss--;
1853 
1854  /* variable should be at given new position, unless new GUB constraint replaced empty old GUB constraint
1855  * (because it was at the end of the GUB constraint array)
1856  */
1857  assert(gubset->gubconssidx[var] == newgubcons
1858  || (newgubcons == gubset->ngubconss && gubset->gubconssidx[var] == oldgubcons));
1859  }
1860 #ifndef NDEBUG
1861  else
1862  assert(gubset->gubconssidx[var] == newgubcons);
1863 #endif
1864 
1865  return SCIP_OKAY;
1866 }
1867 
1868 /** swaps two variables in the same GUB constraint */
1869 static
1870 void GUBsetSwapVars(
1871  SCIP* scip, /**< SCIP data structure */
1872  SCIP_GUBSET* gubset, /**< GUB set data structure */
1873  int var1, /**< first variable to be swapped */
1874  int var2 /**< second variable to be swapped */
1875  )
1876 {
1877  int gubcons;
1878  int var1idx;
1879  GUBVARSTATUS var1status;
1880  int var2idx;
1881  GUBVARSTATUS var2status;
1882 
1883  assert(scip != NULL);
1884  assert(gubset != NULL);
1885 
1886  gubcons = gubset->gubconssidx[var1];
1887  assert(gubcons == gubset->gubconssidx[var2]);
1888 
1889  /* nothing to be done if both variables are the same */
1890  if( var1 == var2 )
1891  return;
1892 
1893  /* swap index and status of variables in GUB constraint */
1894  var1idx = gubset->gubvarsidx[var1];
1895  var1status = gubset->gubconss[gubcons]->gubvarsstatus[var1idx];
1896  var2idx = gubset->gubvarsidx[var2];
1897  var2status = gubset->gubconss[gubcons]->gubvarsstatus[var2idx];
1898 
1899  gubset->gubvarsidx[var1] = var2idx;
1900  gubset->gubconss[gubcons]->gubvars[var1idx] = var2;
1901  gubset->gubconss[gubcons]->gubvarsstatus[var1idx] = var2status;
1902 
1903  gubset->gubvarsidx[var2] = var1idx;
1904  gubset->gubconss[gubcons]->gubvars[var2idx] = var1;
1905  gubset->gubconss[gubcons]->gubvarsstatus[var2idx] = var1status;
1906 }
1907 
1908 /** initializes partition of knapsack variables into nonoverlapping trivial GUB constraints (GUB with one variable) */
1909 static
1911  SCIP* scip, /**< SCIP data structure */
1912  SCIP_GUBSET** gubset, /**< pointer to store GUB set data structure */
1913  int nvars, /**< number of variables in the knapsack constraint */
1914  SCIP_Longint* weights, /**< weights of variables in knapsack constraint */
1915  SCIP_Longint capacity /**< capacity of knapsack */
1916  )
1917 {
1918  int i;
1919 
1920  assert(scip != NULL);
1921  assert(gubset != NULL);
1922  assert(nvars > 0);
1923  assert(weights != NULL);
1924  assert(capacity >= 0);
1925 
1926  /* allocate memory for GUB set data structures */
1927  SCIP_CALL( SCIPallocBuffer(scip, gubset) );
1928  SCIP_CALL( SCIPallocBufferArray(scip, &(*gubset)->gubconss, nvars) );
1929  SCIP_CALL( SCIPallocBufferArray(scip, &(*gubset)->gubconsstatus, nvars) );
1930  SCIP_CALL( SCIPallocBufferArray(scip, &(*gubset)->gubconssidx, nvars) );
1931  SCIP_CALL( SCIPallocBufferArray(scip, &(*gubset)->gubvarsidx, nvars) );
1932  (*gubset)->ngubconss = nvars;
1933  (*gubset)->nvars = nvars;
1934 
1935  /* initialize the set of GUB constraints */
1936  for( i = 0; i < nvars; i++ )
1937  {
1938  /* assign each variable to a new (trivial) GUB constraint */
1939  SCIP_CALL( GUBconsCreate(scip, &(*gubset)->gubconss[i]) );
1940  SCIP_CALL( GUBconsAddVar(scip, (*gubset)->gubconss[i], i) );
1941 
1942  /* set status of GUB constraint to initial */
1943  (*gubset)->gubconsstatus[i] = GUBCONSSTATUS_UNINITIAL;
1944 
1945  (*gubset)->gubconssidx[i] = i;
1946  (*gubset)->gubvarsidx[i] = 0;
1947  assert((*gubset)->gubconss[i]->ngubvars == 1);
1948 
1949  /* already updated status of variable in GUB constraint if it exceeds the capacity of the knapsack */
1950  if( weights[i] > capacity )
1951  (*gubset)->gubconss[(*gubset)->gubconssidx[i]]->gubvarsstatus[(*gubset)->gubvarsidx[i]] = GUBVARSTATUS_CAPACITYEXCEEDED;
1952 
1953  }
1954 
1955  return SCIP_OKAY;
1956 }
1957 
1958 /** frees GUB set data structure */
1959 static
1961  SCIP* scip, /**< SCIP data structure */
1962  SCIP_GUBSET** gubset /**< pointer to GUB set data structure */
1963  )
1964 {
1965  int i;
1966 
1967  assert(scip != NULL);
1968  assert(gubset != NULL);
1969  assert((*gubset)->gubconss != NULL);
1970  assert((*gubset)->gubconsstatus != NULL);
1971  assert((*gubset)->gubconssidx != NULL);
1972  assert((*gubset)->gubvarsidx != NULL);
1973 
1974  /* free all GUB constraints */
1975  for( i = (*gubset)->ngubconss-1; i >= 0; --i )
1976  {
1977  assert((*gubset)->gubconss[i] != NULL);
1978  SCIP_CALL( GUBconsFree(scip, &(*gubset)->gubconss[i]) );
1979  }
1980 
1981  /* free allocated memory */
1982  SCIPfreeBufferArray( scip, &(*gubset)->gubvarsidx );
1983  SCIPfreeBufferArray( scip, &(*gubset)->gubconssidx );
1984  SCIPfreeBufferArray( scip, &(*gubset)->gubconsstatus );
1985  SCIPfreeBufferArray( scip, &(*gubset)->gubconss );
1986  SCIPfreeBuffer(scip, gubset);
1987 
1988  return SCIP_OKAY;
1989 }
1990 
1991 #ifndef NDEBUG
1992 /** checks whether GUB set data structure is consistent */
1993 static
1995  SCIP* scip, /**< SCIP data structure */
1996  SCIP_GUBSET* gubset, /**< GUB set data structure */
1997  SCIP_VAR** vars /**< variables in the knapsack constraint */
1998  )
1999 {
2000  int i;
2001  int gubconsidx;
2002  int gubvaridx;
2003  SCIP_VAR* var1;
2004  SCIP_VAR* var2;
2005  SCIP_Bool var1negated;
2006  SCIP_Bool var2negated;
2007 
2008  assert(scip != NULL);
2009  assert(gubset != NULL);
2010 
2011  SCIPdebugMsg(scip, " GUB set consistency check:\n");
2012 
2013  /* checks for all knapsack vars consistency of stored index of associated gubcons and corresponding index in gubvars */
2014  for( i = 0; i < gubset->nvars; i++ )
2015  {
2016  gubconsidx = gubset->gubconssidx[i];
2017  gubvaridx = gubset->gubvarsidx[i];
2018 
2019  if( gubset->gubconss[gubconsidx]->gubvars[gubvaridx] != i )
2020  {
2021  SCIPdebugMsg(scip, " var<%d> should be in GUB<%d> at position<%d>, but stored is var<%d> instead\n", i,
2022  gubconsidx, gubvaridx, gubset->gubconss[gubconsidx]->gubvars[gubvaridx] );
2023  }
2024  assert(gubset->gubconss[gubconsidx]->gubvars[gubvaridx] == i);
2025  }
2026 
2027  /* checks for each GUB whether all pairs of its variables have a common clique */
2028  for( i = 0; i < gubset->ngubconss; i++ )
2029  {
2030  int j;
2031 
2032  for( j = 0; j < gubset->gubconss[i]->ngubvars; j++ )
2033  {
2034  int k;
2035 
2036  /* get corresponding active problem variable */
2037  var1 = vars[gubset->gubconss[i]->gubvars[j]];
2038  var1negated = FALSE;
2039  SCIP_CALL( SCIPvarGetProbvarBinary(&var1, &var1negated) );
2040 
2041  for( k = j+1; k < gubset->gubconss[i]->ngubvars; k++ )
2042  {
2043  /* get corresponding active problem variable */
2044  var2 = vars[gubset->gubconss[i]->gubvars[k]];
2045  var2negated = FALSE;
2046  SCIP_CALL( SCIPvarGetProbvarBinary(&var2, &var2negated) );
2047 
2048  if( !SCIPvarsHaveCommonClique(var1, !var1negated, var2, !var2negated, TRUE) )
2049  {
2050  SCIPdebugMsg(scip, " GUB<%d>: var<%d,%s> and var<%d,%s> do not share a clique\n", i, j,
2051  SCIPvarGetName(vars[gubset->gubconss[i]->gubvars[j]]), k,
2052  SCIPvarGetName(vars[gubset->gubconss[i]->gubvars[k]]));
2053  SCIPdebugMsg(scip, " GUB<%d>: var<%d,%s> and var<%d,%s> do not share a clique\n", i, j,
2054  SCIPvarGetName(var1), k,
2055  SCIPvarGetName(var2));
2056  }
2057 
2058  /* @todo: in case we used also negated cliques for the GUB partition, this assert has to be changed */
2059  assert(SCIPvarsHaveCommonClique(var1, !var1negated, var2, !var2negated, TRUE));
2060  }
2061  }
2062  }
2063  SCIPdebugMsg(scip, " --> successful\n");
2064 
2065  return SCIP_OKAY;
2066 }
2067 #endif
2068 
2069 /** calculates a partition of the given set of binary variables into cliques;
2070  * afterwards the output array contains one value for each variable, such that two variables got the same value iff they
2071  * were assigned to the same clique;
2072  * the first variable is always assigned to clique 0, and a variable can only be assigned to clique i if at least one of
2073  * the preceding variables was assigned to clique i-1;
2074  * note: in contrast to SCIPcalcCliquePartition(), variables with LP value 1 are put into trivial cliques (with one
2075  * variable) and for the remaining variables, a partition with a small number of cliques is constructed
2076  */
2077 
2078 static
2080  SCIP*const scip, /**< SCIP data structure */
2081  SCIP_VAR**const vars, /**< binary variables in the clique from which at most one can be set to 1 */
2082  int const nvars, /**< number of variables in the clique */
2083  int*const cliquepartition, /**< array of length nvars to store the clique partition */
2084  int*const ncliques, /**< pointer to store number of cliques actually contained in the partition */
2085  SCIP_Real* solvals /**< solution values of all given binary variables */
2086  )
2088  SCIP_VAR** tmpvars;
2089  SCIP_VAR** cliquevars;
2090  SCIP_Bool* cliquevalues;
2091  SCIP_Bool* tmpvalues;
2092  int* varseq;
2093  int* sortkeys;
2094  int ncliquevars;
2095  int maxncliquevarscomp;
2096  int nignorevars;
2097  int nvarsused;
2098  int i;
2099 
2100  assert(scip != NULL);
2101  assert(nvars == 0 || vars != NULL);
2102  assert(nvars == 0 || cliquepartition != NULL);
2103  assert(ncliques != NULL);
2104 
2105  if( nvars == 0 )
2106  {
2107  *ncliques = 0;
2108  return SCIP_OKAY;
2109  }
2110 
2111  /* allocate temporary memory for storing the variables of the current clique */
2112  SCIP_CALL( SCIPallocBufferArray(scip, &cliquevars, nvars) );
2113  SCIP_CALL( SCIPallocBufferArray(scip, &cliquevalues, nvars) );
2114  SCIP_CALL( SCIPallocBufferArray(scip, &tmpvalues, nvars) );
2115  SCIP_CALL( SCIPduplicateBufferArray(scip, &tmpvars, vars, nvars) );
2116  SCIP_CALL( SCIPallocBufferArray(scip, &varseq, nvars) );
2117  SCIP_CALL( SCIPallocBufferArray(scip, &sortkeys, nvars) );
2118 
2119  /* initialize the cliquepartition array with -1 */
2120  /* initialize the tmpvalues array */
2121  for( i = nvars - 1; i >= 0; --i )
2122  {
2123  tmpvalues[i] = TRUE;
2124  cliquepartition[i] = -1;
2125  }
2126 
2127  /* get corresponding active problem variables */
2128  SCIP_CALL( SCIPvarsGetProbvarBinary(&tmpvars, &tmpvalues, nvars) );
2129 
2130  /* ignore variables with LP value 1 (will be assigned to trivial GUBs at the end) and sort remaining variables
2131  * by nondecreasing number of cliques the variables are in
2132  */
2133  nignorevars = 0;
2134  nvarsused = 0;
2135  for( i = 0; i < nvars; i++ )
2136  {
2137  if( SCIPisFeasEQ(scip, solvals[i], 1.0) )
2138  {
2139  /* variables with LP value 1 are put to the end of varseq array and will not be sorted */
2140  varseq[nvars-1-nignorevars] = i;
2141  nignorevars++;
2142  }
2143  else
2144  {
2145  /* remaining variables are put to the front of varseq array and will be sorted by their number of cliques */
2146  varseq[nvarsused] = i;
2147  sortkeys[nvarsused] = SCIPvarGetNCliques(tmpvars[i], tmpvalues[i]);
2148  nvarsused++;
2149  }
2150  }
2151  assert(nvarsused + nignorevars == nvars);
2152 
2153  /* sort variables with LP value less than 1 by nondecreasing order of the number of cliques they are in */
2154  SCIPsortIntInt(sortkeys, varseq, nvarsused);
2155 
2156  maxncliquevarscomp = MIN(nvars*nvars, MAXNCLIQUEVARSCOMP);
2157 
2158  /* calculate the clique partition */
2159  *ncliques = 0;
2160  for( i = 0; i < nvars; ++i )
2161  {
2162  if( cliquepartition[varseq[i]] == -1 )
2163  {
2164  int j;
2165 
2166  /* variable starts a new clique */
2167  cliquepartition[varseq[i]] = *ncliques;
2168  cliquevars[0] = tmpvars[varseq[i]];
2169  cliquevalues[0] = tmpvalues[varseq[i]];
2170  ncliquevars = 1;
2171 
2172  /* if variable is not active (multi-aggregated or fixed), it cannot be in any clique and
2173  * if the variable has LP value 1 we do not want it to be in nontrivial cliques
2174  */
2175  if( SCIPvarIsActive(tmpvars[varseq[i]]) && i < nvarsused )
2176  {
2177  /* greedily fill up the clique */
2178  for( j = i + 1; j < nvarsused; ++j )
2179  {
2180  /* if variable is not active (multi-aggregated or fixed), it cannot be in any clique */
2181  if( cliquepartition[varseq[j]] == -1 && SCIPvarIsActive(tmpvars[varseq[j]]) )
2182  {
2183  int k;
2184 
2185  /* check if every variable in the actual clique is in clique with the new variable */
2186  for( k = ncliquevars - 1; k >= 0; --k )
2187  {
2188  if( !SCIPvarsHaveCommonClique(tmpvars[varseq[j]], tmpvalues[varseq[j]], cliquevars[k],
2189  cliquevalues[k], TRUE) )
2190  break;
2191  }
2192 
2193  if( k == -1 )
2194  {
2195  /* put the variable into the same clique */
2196  cliquepartition[varseq[j]] = cliquepartition[varseq[i]];
2197  cliquevars[ncliquevars] = tmpvars[varseq[j]];
2198  cliquevalues[ncliquevars] = tmpvalues[varseq[j]];
2199  ++ncliquevars;
2200  }
2201  }
2202  }
2203  }
2204 
2205  /* this clique is finished */
2206  ++(*ncliques);
2207  }
2208  assert(cliquepartition[varseq[i]] >= 0 && cliquepartition[varseq[i]] < i + 1);
2209 
2210  /* break if we reached the maximal number of comparisons */
2211  if( i * nvars > maxncliquevarscomp )
2212  break;
2213  }
2214  /* if we had too many variables fill up the cliquepartition and put each variable in a separate clique */
2215  for( ; i < nvars; ++i )
2216  {
2217  if( cliquepartition[varseq[i]] == -1 )
2218  {
2219  cliquepartition[varseq[i]] = *ncliques;
2220  ++(*ncliques);
2221  }
2222  }
2223 
2224  /* free temporary memory */
2225  SCIPfreeBufferArray(scip, &sortkeys);
2226  SCIPfreeBufferArray(scip, &varseq);
2227  SCIPfreeBufferArray(scip, &tmpvars);
2228  SCIPfreeBufferArray(scip, &tmpvalues);
2229  SCIPfreeBufferArray(scip, &cliquevalues);
2230  SCIPfreeBufferArray(scip, &cliquevars);
2231 
2232  return SCIP_OKAY;
2233 }
2234 
2235 /** constructs sophisticated partition of knapsack variables into non-overlapping GUBs; current partition uses trivial GUBs */
2236 static
2238  SCIP* scip, /**< SCIP data structure */
2239  SCIP_GUBSET* gubset, /**< GUB set data structure */
2240  SCIP_VAR** vars, /**< variables in the knapsack constraint */
2241  SCIP_Real* solvals /**< solution values of all knapsack variables */
2242  )
2243 {
2244  int* cliquepartition;
2245  int* gubfirstvar;
2246  int ncliques;
2247  int currentgubconsidx;
2248  int newgubconsidx;
2249  int cliqueidx;
2250  int nvars;
2251  int i;
2252 
2253  assert(scip != NULL);
2254  assert(gubset != NULL);
2255  assert(vars != NULL);
2256 
2257  nvars = gubset->nvars;
2258  assert(nvars >= 0);
2259 
2260  /* allocate temporary memory for clique partition */
2261  SCIP_CALL( SCIPallocBufferArray(scip, &cliquepartition, nvars) );
2262 
2263  /* compute sophisticated clique partition */
2264  SCIP_CALL( GUBsetCalcCliquePartition(scip, vars, nvars, cliquepartition, &ncliques, solvals) );
2265 
2266  /* allocate temporary memory for GUB set data structure */
2267  SCIP_CALL( SCIPallocBufferArray(scip, &gubfirstvar, ncliques) );
2268 
2269  /* translate GUB partition into GUB set data structure */
2270  for( i = 0; i < ncliques; i++ )
2271  {
2272  /* initialize first variable for every GUB */
2273  gubfirstvar[i] = -1;
2274  }
2275  /* move every knapsack variable into GUB defined by clique partition */
2276  for( i = 0; i < nvars; i++ )
2277  {
2278  assert(cliquepartition[i] >= 0);
2279 
2280  cliqueidx = cliquepartition[i];
2281  currentgubconsidx = gubset->gubconssidx[i];
2282  assert(gubset->gubconss[currentgubconsidx]->ngubvars == 1 );
2283 
2284  /* variable is first element in GUB constraint defined by clique partition */
2285  if( gubfirstvar[cliqueidx] == -1 )
2286  {
2287  /* corresponding GUB constraint in GUB set data structure was already constructed (as initial trivial GUB);
2288  * note: no assert for gubconssidx, because it can changed due to deleting empty GUBs in GUBsetMoveVar()
2289  */
2290  assert(gubset->gubvarsidx[i] == 0);
2291  assert(gubset->gubconss[gubset->gubconssidx[i]]->gubvars[gubset->gubvarsidx[i]] == i);
2292 
2293  /* remember the first variable found for the current GUB */
2294  gubfirstvar[cliqueidx] = i;
2295  }
2296  /* variable is additional element of GUB constraint defined by clique partition */
2297  else
2298  {
2299  assert(gubfirstvar[cliqueidx] >= 0 && gubfirstvar[cliqueidx] < i);
2300 
2301  /* move variable to GUB constraint defined by clique partition; index of this GUB constraint is given by the
2302  * first variable of this GUB constraint
2303  */
2304  newgubconsidx = gubset->gubconssidx[gubfirstvar[cliqueidx]];
2305  assert(newgubconsidx != currentgubconsidx); /* because initially every variable is in a different GUB */
2306  SCIP_CALL( GUBsetMoveVar(scip, gubset, vars, i, currentgubconsidx, newgubconsidx) );
2307 
2308  assert(gubset->gubconss[gubset->gubconssidx[i]]->gubvars[gubset->gubvarsidx[i]] == i);
2309  }
2310  }
2311 
2312 #ifdef SCIP_DEBUG
2313  /* prints GUB set data structure */
2314  GUBsetPrint(scip, gubset, vars, solvals);
2315 #endif
2316 
2317 #ifndef NDEBUG
2318  /* checks consistency of GUB set data structure */
2319  SCIP_CALL( GUBsetCheck(scip, gubset, vars) );
2320 #endif
2321 
2322  /* free temporary memory */
2323  SCIPfreeBufferArray(scip, &gubfirstvar);
2324  SCIPfreeBufferArray(scip, &cliquepartition);
2325 
2326  return SCIP_OKAY;
2327 }
2328 
2329 /** gets a most violated cover C (\f$\sum_{j \in C} a_j > a_0\f$) for a given knapsack constraint \f$\sum_{j \in N} a_j x_j \leq a_0\f$
2330  * taking into consideration the following fixing: \f$j \in C\f$, if \f$j \in N_1 = \{j \in N : x^*_j = 1\}\f$ and
2331  * \f$j \in N \setminus C\f$, if \f$j \in N_0 = \{j \in N : x^*_j = 0\}\f$, if one exists.
2332  */
2333 static
2335  SCIP* scip, /**< SCIP data structure */
2336  SCIP_VAR** vars, /**< variables in knapsack constraint */
2337  int nvars, /**< number of variables in knapsack constraint */
2338  SCIP_Longint* weights, /**< weights of variables in knapsack constraint */
2339  SCIP_Longint capacity, /**< capacity of knapsack */
2340  SCIP_Real* solvals, /**< solution values of all problem variables */
2341  int* covervars, /**< pointer to store cover variables */
2342  int* noncovervars, /**< pointer to store noncover variables */
2343  int* ncovervars, /**< pointer to store number of cover variables */
2344  int* nnoncovervars, /**< pointer to store number of noncover variables */
2345  SCIP_Longint* coverweight, /**< pointer to store weight of cover */
2346  SCIP_Bool* found, /**< pointer to store whether a cover was found */
2347  SCIP_Bool modtransused, /**< should modified transformed separation problem be used to find cover */
2348  int* ntightened, /**< pointer to store number of variables with tightened upper bound */
2349  SCIP_Bool* fractional /**< pointer to store whether the LP sol for knapsack vars is fractional */
2350  )
2351 {
2352  SCIP_Longint* transweights;
2353  SCIP_Real* transprofits;
2354  SCIP_Longint transcapacity;
2355  SCIP_Longint fixedonesweight;
2356  SCIP_Longint itemsweight;
2357  SCIP_Bool infeasible;
2358  int* fixedones;
2359  int* fixedzeros;
2360  int* items;
2361  int nfixedones;
2362  int nfixedzeros;
2363  int nitems;
2364  int j;
2365 
2366  assert(scip != NULL);
2367  assert(vars != NULL);
2368  assert(nvars > 0);
2369  assert(weights != NULL);
2370  assert(capacity >= 0);
2371  assert(solvals != NULL);
2372  assert(covervars != NULL);
2373  assert(noncovervars != NULL);
2374  assert(ncovervars != NULL);
2375  assert(nnoncovervars != NULL);
2376  assert(coverweight != NULL);
2377  assert(found != NULL);
2378  assert(ntightened != NULL);
2379  assert(fractional != NULL);
2380 
2381  SCIPdebugMsg(scip, " get cover for knapsack constraint\n");
2382 
2383  /* allocates temporary memory */
2384  SCIP_CALL( SCIPallocBufferArray(scip, &transweights, nvars) );
2385  SCIP_CALL( SCIPallocBufferArray(scip, &transprofits, nvars) );
2386  SCIP_CALL( SCIPallocBufferArray(scip, &fixedones, nvars) );
2387  SCIP_CALL( SCIPallocBufferArray(scip, &fixedzeros, nvars) );
2388  SCIP_CALL( SCIPallocBufferArray(scip, &items, nvars) );
2389 
2390  *found = FALSE;
2391  *ncovervars = 0;
2392  *nnoncovervars = 0;
2393  *coverweight = 0;
2394  *fractional = TRUE;
2395 
2396  /* gets the following sets
2397  * N_1 = {j in N : x*_j = 1} (fixedones),
2398  * N_0 = {j in N : x*_j = 0} (fixedzeros) and
2399  * N\(N_0 & N_1) (items),
2400  * where x*_j is the solution value of variable x_j
2401  */
2402  nfixedones = 0;
2403  nfixedzeros = 0;
2404  nitems = 0;
2405  fixedonesweight = 0;
2406  itemsweight = 0;
2407  *ntightened = 0;
2408  for( j = 0; j < nvars; j++ )
2409  {
2410  assert(SCIPvarIsBinary(vars[j]));
2411 
2412  /* tightens upper bound of x_j if weight of x_j is greater than capacity of knapsack */
2413  if( weights[j] > capacity )
2414  {
2415  SCIP_CALL( SCIPtightenVarUb(scip, vars[j], 0.0, FALSE, &infeasible, NULL) );
2416  assert(!infeasible);
2417  (*ntightened)++;
2418  continue;
2419  }
2420 
2421  /* variable x_j has solution value one */
2422  if( SCIPisFeasEQ(scip, solvals[j], 1.0) )
2423  {
2424  fixedones[nfixedones] = j;
2425  nfixedones++;
2426  fixedonesweight += weights[j];
2427  }
2428  /* variable x_j has solution value zero */
2429  else if( SCIPisFeasEQ(scip, solvals[j], 0.0) )
2430  {
2431  fixedzeros[nfixedzeros] = j;
2432  nfixedzeros++;
2433  }
2434  /* variable x_j has fractional solution value */
2435  else
2436  {
2437  assert( SCIPisFeasGT(scip, solvals[j], 0.0) && SCIPisFeasLT(scip, solvals[j], 1.0) );
2438  items[nitems] = j;
2439  nitems++;
2440  itemsweight += weights[j];
2441  }
2442  }
2443  assert(nfixedones + nfixedzeros + nitems == nvars - (*ntightened));
2444 
2445  /* sets whether the LP solution x* for the knapsack variables is fractional; if it is not fractional we stop
2446  * the separation routine
2447  */
2448  assert(nitems >= 0);
2449  if( nitems == 0 )
2450  {
2451  *fractional = FALSE;
2452  goto TERMINATE;
2453  }
2454  assert(*fractional);
2455 
2456  /* transforms the traditional separation problem (under consideration of the following fixing:
2457  * z_j = 1 for all j in N_1, z_j = 0 for all j in N_0)
2458  *
2459  * min sum_{j in N\(N_0 & N_1)} (1 - x*_j) z_j
2460  * sum_{j in N\(N_0 & N_1)} a_j z_j >= (a_0 + 1) - sum_{j in N_1} a_j
2461  * z_j in {0,1}, j in N\(N_0 & N_1)
2462  *
2463  * to a knapsack problem in maximization form by complementing the variables
2464  *
2465  * sum_{j in N\(N_0 & N_1)} (1 - x*_j) -
2466  * max sum_{j in N\(N_0 & N_1)} (1 - x*_j) z_j
2467  * sum_{j in N\(N_0 & N_1)} a_j z_j <= sum_{j in N\N_0} a_j - (a_0 + 1)
2468  * z_j in {0,1}, j in N\(N_0 & N_1)
2469  */
2470 
2471  /* gets weight and profit of variables in transformed knapsack problem */
2472  for( j = 0; j < nitems; j++ )
2473  {
2474  transweights[j] = weights[items[j]];
2475  transprofits[j] = 1.0 - solvals[items[j]];
2476  }
2477  /* gets capacity of transformed knapsack problem */
2478  transcapacity = fixedonesweight + itemsweight - capacity - 1;
2479 
2480  /* if capacity of transformed knapsack problem is less than zero, there is no cover
2481  * (when variables fixed to zero are not used)
2482  */
2483  if( transcapacity < 0 )
2484  {
2485  assert(!(*found));
2486  goto TERMINATE;
2487  }
2488 
2489  if( modtransused )
2490  {
2491  /* transforms the modified separation problem (under consideration of the following fixing:
2492  * z_j = 1 for all j in N_1, z_j = 0 for all j in N_0)
2493  *
2494  * min sum_{j in N\(N_0 & N_1)} (1 - x*_j) a_j z_j
2495  * sum_{j in N\(N_0 & N_1)} a_j z_j >= (a_0 + 1) - sum_{j in N_1} a_j
2496  * z_j in {0,1}, j in N\(N_0 & N_1)
2497  *
2498  * to a knapsack problem in maximization form by complementing the variables
2499  *
2500  * sum_{j in N\(N_0 & N_1)} (1 - x*_j) a_j -
2501  * max sum_{j in N\(N_0 & N_1)} (1 - x*_j) a_j z_j
2502  * sum_{j in N\(N_0 & N_1)} a_j z_j <= sum_{j in N\N_0} a_j - (a_0 + 1)
2503  * z_j in {0,1}, j in N\(N_0 & N_1)
2504  */
2505 
2506  /* gets weight and profit of variables in modified transformed knapsack problem */
2507  for( j = 0; j < nitems; j++ )
2508  {
2509  transprofits[j] *= weights[items[j]];
2510  assert(SCIPisFeasPositive(scip, transprofits[j]));
2511  }
2512  }
2513 
2514  /* solves (modified) transformed knapsack problem approximately by solving the LP-relaxation of the (modified)
2515  * transformed knapsack problem using Dantzig's method and rounding down the solution.
2516  * let z* be the solution, then
2517  * j in C, if z*_j = 0 and
2518  * i in N\C, if z*_j = 1.
2519  */
2520  SCIP_CALL( SCIPsolveKnapsackApproximately(scip, nitems, transweights, transprofits, transcapacity, items,
2521  noncovervars, covervars, nnoncovervars, ncovervars, NULL) );
2522  /*assert(checkSolveKnapsack(scip, nitems, transweights, transprofits, items, weights, solvals, modtransused));*/
2523 
2524  /* constructs cover C (sum_{j in C} a_j > a_0) */
2525  for( j = 0; j < *ncovervars; j++ )
2526  {
2527  (*coverweight) += weights[covervars[j]];
2528  }
2529 
2530  /* adds all variables from N_1 to C */
2531  for( j = 0; j < nfixedones; j++ )
2532  {
2533  covervars[*ncovervars] = fixedones[j];
2534  (*ncovervars)++;
2535  (*coverweight) += weights[fixedones[j]];
2536  }
2537 
2538  /* adds all variables from N_0 to N\C */
2539  for( j = 0; j < nfixedzeros; j++ )
2540  {
2541  noncovervars[*nnoncovervars] = fixedzeros[j];
2542  (*nnoncovervars)++;
2543  }
2544  assert((*ncovervars) + (*nnoncovervars) == nvars - (*ntightened));
2545  assert((*coverweight) > capacity);
2546  *found = TRUE;
2547 
2548  TERMINATE:
2549  /* frees temporary memory */
2550  SCIPfreeBufferArray(scip, &items);
2551  SCIPfreeBufferArray(scip, &fixedzeros);
2552  SCIPfreeBufferArray(scip, &fixedones);
2553  SCIPfreeBufferArray(scip, &transprofits);
2554  SCIPfreeBufferArray(scip, &transweights);
2555 
2556  SCIPdebugMsg(scip, " get cover for knapsack constraint -- end\n");
2557 
2558  return SCIP_OKAY;
2559 }
2560 
2561 #ifndef NDEBUG
2562 /** checks if minweightidx is set correctly
2563  */
2564 static
2566  SCIP_Longint* weights, /**< weights of variables in knapsack constraint */
2567  SCIP_Longint capacity, /**< capacity of knapsack */
2568  int* covervars, /**< pointer to store cover variables */
2569  int ncovervars, /**< pointer to store number of cover variables */
2570  SCIP_Longint coverweight, /**< pointer to store weight of cover */
2571  int minweightidx, /**< index of variable in cover variables with minimum weight */
2572  int j /**< current index in cover variables */
2573  )
2574 {
2575  SCIP_Longint minweight;
2576  int i;
2577 
2578  assert(weights != NULL);
2579  assert(covervars != NULL);
2580  assert(ncovervars > 0);
2581 
2582  minweight = weights[covervars[minweightidx]];
2583 
2584  /* checks if all cover variables before index j have weight greater than minweight */
2585  for( i = 0; i < j; i++ )
2586  {
2587  assert(weights[covervars[i]] > minweight);
2588  if( weights[covervars[i]] <= minweight )
2589  return FALSE;
2590  }
2591 
2592  /* checks if all variables before index j cannot be removed, i.e. i cannot be the next minweightidx */
2593  for( i = 0; i < j; i++ )
2594  {
2595  assert(coverweight - weights[covervars[i]] <= capacity);
2596  if( coverweight - weights[covervars[i]] > capacity )
2597  return FALSE;
2598  }
2599  return TRUE;
2600 }
2601 #endif
2602 
2603 
2604 /** gets partition \f$(C_1,C_2)\f$ of minimal cover \f$C\f$, i.e. \f$C_1 \cup C_2 = C\f$ and \f$C_1 \cap C_2 = \emptyset\f$,
2605  * with \f$C_1\f$ not empty; chooses partition as follows \f$C_2 = \{ j \in C : x^*_j = 1 \}\f$ and \f$C_1 = C \setminus C_2\f$
2606  */
2607 static
2609  SCIP* scip, /**< SCIP data structure */
2610  SCIP_Real* solvals, /**< solution values of all problem variables */
2611  int* covervars, /**< cover variables */
2612  int ncovervars, /**< number of cover variables */
2613  int* varsC1, /**< pointer to store variables in C1 */
2614  int* varsC2, /**< pointer to store variables in C2 */
2615  int* nvarsC1, /**< pointer to store number of variables in C1 */
2616  int* nvarsC2 /**< pointer to store number of variables in C2 */
2617  )
2618 {
2619  int j;
2620 
2621  assert(scip != NULL);
2622  assert(ncovervars >= 0);
2623  assert(solvals != NULL);
2624  assert(covervars != NULL);
2625  assert(varsC1 != NULL);
2626  assert(varsC2 != NULL);
2627  assert(nvarsC1 != NULL);
2628  assert(nvarsC2 != NULL);
2629 
2630  *nvarsC1 = 0;
2631  *nvarsC2 = 0;
2632  for( j = 0; j < ncovervars; j++ )
2633  {
2634  assert(SCIPisFeasGT(scip, solvals[covervars[j]], 0.0));
2635 
2636  /* variable has solution value one */
2637  if( SCIPisGE(scip, solvals[covervars[j]], 1.0) )
2638  {
2639  varsC2[*nvarsC2] = covervars[j];
2640  (*nvarsC2)++;
2641  }
2642  /* variable has solution value less than one */
2643  else
2644  {
2645  assert(SCIPisLT(scip, solvals[covervars[j]], 1.0));
2646  varsC1[*nvarsC1] = covervars[j];
2647  (*nvarsC1)++;
2648  }
2649  }
2650  assert((*nvarsC1) + (*nvarsC2) == ncovervars);
2651 }
2652 
2653 /** changes given partition (C_1,C_2) of minimal cover C, if |C1| = 1, by moving one and two (if possible) variables from
2654  * C2 to C1 if |C1| = 1 and |C1| = 0, respectively.
2655  */
2656 static
2658  SCIP* scip, /**< SCIP data structure */
2659  SCIP_Longint* weights, /**< weights of variables in knapsack constraint */
2660  int* varsC1, /**< pointer to store variables in C1 */
2661  int* varsC2, /**< pointer to store variables in C2 */
2662  int* nvarsC1, /**< pointer to store number of variables in C1 */
2663  int* nvarsC2 /**< pointer to store number of variables in C2 */
2664  )
2666  SCIP_Real* sortkeysC2;
2667  int j;
2668 
2669  assert(*nvarsC1 >= 0 && *nvarsC1 <= 1);
2670  assert(*nvarsC2 > 0);
2671 
2672  /* allocates temporary memory */
2673  SCIP_CALL( SCIPallocBufferArray(scip, &sortkeysC2, *nvarsC2) );
2674 
2675  /* sorts variables in C2 such that a_1 >= .... >= a_|C2| */
2676  for( j = 0; j < *nvarsC2; j++ )
2677  sortkeysC2[j] = (SCIP_Real) weights[varsC2[j]];
2678  SCIPsortDownRealInt(sortkeysC2, varsC2, *nvarsC2);
2679 
2680  /* adds one or two variable from C2 with smallest weight to C1 and removes them from C2 */
2681  assert(*nvarsC2 == 1 || weights[varsC2[(*nvarsC2)-1]] <= weights[varsC2[(*nvarsC2)-2]]);
2682  while( *nvarsC1 < 2 && *nvarsC2 > 0 )
2683  {
2684  varsC1[*nvarsC1] = varsC2[(*nvarsC2)-1];
2685  (*nvarsC1)++;
2686  (*nvarsC2)--;
2687  }
2688 
2689  /* frees temporary memory */
2690  SCIPfreeBufferArray(scip, &sortkeysC2);
2691 
2692  return SCIP_OKAY;
2693 }
2694 
2695 /** changes given partition (C_1,C_2) of feasible set C, if |C1| = 1, by moving one variable from C2 to C1 */
2696 static
2698  SCIP* scip, /**< SCIP data structure */
2699  SCIP_Longint* weights, /**< weights of variables in knapsack constraint */
2700  int* varsC1, /**< pointer to store variables in C1 */
2701  int* varsC2, /**< pointer to store variables in C2 */
2702  int* nvarsC1, /**< pointer to store number of variables in C1 */
2703  int* nvarsC2 /**< pointer to store number of variables in C2 */
2704  )
2706  SCIP_Real* sortkeysC2;
2707  int j;
2708 
2709  assert(*nvarsC1 >= 0 && *nvarsC1 <= 1);
2710  assert(*nvarsC2 > 0);
2711 
2712  /* allocates temporary memory */
2713  SCIP_CALL( SCIPallocBufferArray(scip, &sortkeysC2, *nvarsC2) );
2714 
2715  /* sorts variables in C2 such that a_1 >= .... >= a_|C2| */
2716  for( j = 0; j < *nvarsC2; j++ )
2717  sortkeysC2[j] = (SCIP_Real) weights[varsC2[j]];
2718  SCIPsortDownRealInt(sortkeysC2, varsC2, *nvarsC2);
2719 
2720  /* adds variable from C2 with smallest weight to C1 and removes it from C2 */
2721  assert(*nvarsC2 == 1 || weights[varsC2[(*nvarsC2)-1]] <= weights[varsC2[(*nvarsC2)-2]]);
2722  varsC1[*nvarsC1] = varsC2[(*nvarsC2)-1];
2723  (*nvarsC1)++;
2724  (*nvarsC2)--;
2725 
2726  /* frees temporary memory */
2727  SCIPfreeBufferArray(scip, &sortkeysC2);
2728 
2729  return SCIP_OKAY;
2730 }
2731 
2732 
2733 /** gets partition \f$(F,R)\f$ of \f$N \setminus C\f$ where \f$C\f$ is a minimal cover, i.e. \f$F \cup R = N \setminus C\f$
2734  * and \f$F \cap R = \emptyset\f$; chooses partition as follows \f$R = \{ j \in N \setminus C : x^*_j = 0 \}\f$ and
2735  * \f$F = (N \setminus C) \setminus F\f$
2736  */
2737 static
2739  SCIP* scip, /**< SCIP data structure */
2740  SCIP_Real* solvals, /**< solution values of all problem variables */
2741  int* noncovervars, /**< noncover variables */
2742  int nnoncovervars, /**< number of noncover variables */
2743  int* varsF, /**< pointer to store variables in F */
2744  int* varsR, /**< pointer to store variables in R */
2745  int* nvarsF, /**< pointer to store number of variables in F */
2746  int* nvarsR /**< pointer to store number of variables in R */
2747  )
2748 {
2749  int j;
2750 
2751  assert(scip != NULL);
2752  assert(nnoncovervars >= 0);
2753  assert(solvals != NULL);
2754  assert(noncovervars != NULL);
2755  assert(varsF != NULL);
2756  assert(varsR != NULL);
2757  assert(nvarsF != NULL);
2758  assert(nvarsR != NULL);
2759 
2760  *nvarsF = 0;
2761  *nvarsR = 0;
2762 
2763  for( j = 0; j < nnoncovervars; j++ )
2764  {
2765  /* variable has solution value zero */
2766  if( SCIPisFeasEQ(scip, solvals[noncovervars[j]], 0.0) )
2767  {
2768  varsR[*nvarsR] = noncovervars[j];
2769  (*nvarsR)++;
2770  }
2771  /* variable has solution value greater than zero */
2772  else
2773  {
2774  assert(SCIPisFeasGT(scip, solvals[noncovervars[j]], 0.0));
2775  varsF[*nvarsF] = noncovervars[j];
2776  (*nvarsF)++;
2777  }
2778  }
2779  assert((*nvarsF) + (*nvarsR) == nnoncovervars);
2780 }
2781 
2782 /** sorts variables in F, C_2, and R according to the second level lifting sequence that will be used in the sequential
2783  * lifting procedure
2784  */
2785 static
2787  SCIP* scip, /**< SCIP data structure */
2788  SCIP_Real* solvals, /**< solution values of all problem variables */
2789  SCIP_Longint* weights, /**< weights of variables in knapsack constraint */
2790  int* varsF, /**< pointer to store variables in F */
2791  int* varsC2, /**< pointer to store variables in C2 */
2792  int* varsR, /**< pointer to store variables in R */
2793  int nvarsF, /**< number of variables in F */
2794  int nvarsC2, /**< number of variables in C2 */
2795  int nvarsR /**< number of variables in R */
2796  )
2797 {
2798  SORTKEYPAIR** sortkeypairsF;
2799  SORTKEYPAIR* sortkeypairsFstore;
2800  SCIP_Real* sortkeysC2;
2801  SCIP_Real* sortkeysR;
2802  int j;
2803 
2804  assert(scip != NULL);
2805  assert(solvals != NULL);
2806  assert(weights != NULL);
2807  assert(varsF != NULL);
2808  assert(varsC2 != NULL);
2809  assert(varsR != NULL);
2810  assert(nvarsF >= 0);
2811  assert(nvarsC2 >= 0);
2812  assert(nvarsR >= 0);
2813 
2814  /* allocates temporary memory */
2815  SCIP_CALL( SCIPallocBufferArray(scip, &sortkeypairsF, nvarsF) );
2816  SCIP_CALL( SCIPallocBufferArray(scip, &sortkeypairsFstore, nvarsF) );
2817  SCIP_CALL( SCIPallocBufferArray(scip, &sortkeysC2, nvarsC2) );
2818  SCIP_CALL( SCIPallocBufferArray(scip, &sortkeysR, nvarsR) );
2819 
2820  /* gets sorting key for variables in F corresponding to the following lifting sequence
2821  * sequence 1: non-increasing absolute difference between x*_j and the value the variable is fixed to, i.e.
2822  * x*_1 >= x*_2 >= ... >= x*_|F|
2823  * in case of equality uses
2824  * sequence 4: non-increasing a_j, i.e. a_1 >= a_2 >= ... >= a_|C_2|
2825  */
2826  for( j = 0; j < nvarsF; j++ )
2827  {
2828  sortkeypairsF[j] = &(sortkeypairsFstore[j]);
2829  sortkeypairsF[j]->key1 = solvals[varsF[j]];
2830  sortkeypairsF[j]->key2 = (SCIP_Real) weights[varsF[j]];
2831  }
2832 
2833  /* gets sorting key for variables in C_2 corresponding to the following lifting sequence
2834  * sequence 4: non-increasing a_j, i.e. a_1 >= a_2 >= ... >= a_|C_2|
2835  */
2836  for( j = 0; j < nvarsC2; j++ )
2837  sortkeysC2[j] = (SCIP_Real) weights[varsC2[j]];
2838 
2839  /* gets sorting key for variables in R corresponding to the following lifting sequence
2840  * sequence 4: non-increasing a_j, i.e. a_1 >= a_2 >= ... >= a_|R|
2841  */
2842  for( j = 0; j < nvarsR; j++ )
2843  sortkeysR[j] = (SCIP_Real) weights[varsR[j]];
2844 
2845  /* sorts F, C2 and R */
2846  if( nvarsF > 0 )
2847  {
2848  SCIPsortDownPtrInt((void**)sortkeypairsF, varsF, compSortkeypairs, nvarsF);
2849  }
2850  if( nvarsC2 > 0 )
2851  {
2852  SCIPsortDownRealInt(sortkeysC2, varsC2, nvarsC2);
2853  }
2854  if( nvarsR > 0)
2855  {
2856  SCIPsortDownRealInt(sortkeysR, varsR, nvarsR);
2857  }
2858 
2859  /* frees temporary memory */
2860  SCIPfreeBufferArray(scip, &sortkeysR);
2861  SCIPfreeBufferArray(scip, &sortkeysC2);
2862  SCIPfreeBufferArray(scip, &sortkeypairsFstore);
2863  SCIPfreeBufferArray(scip, &sortkeypairsF);
2864 
2865  return SCIP_OKAY;
2866 }
2867 
2868 /** categorizes GUBs of knapsack GUB partion into GOC1, GNC1, GF, GC2, and GR and computes a lifting sequence of the GUBs
2869  * for the sequential GUB wise lifting procedure
2870  */
2871 static
2873  SCIP* scip, /**< SCIP data structure */
2874  SCIP_GUBSET* gubset, /**< GUB set data structure */
2875  SCIP_Real* solvals, /**< solution values of variables in knapsack constraint */
2876  SCIP_Longint* weights, /**< weights of variables in knapsack constraint */
2877  int* varsC1, /**< variables in C1 */
2878  int* varsC2, /**< variables in C2 */
2879  int* varsF, /**< variables in F */
2880  int* varsR, /**< variables in R */
2881  int nvarsC1, /**< number of variables in C1 */
2882  int nvarsC2, /**< number of variables in C2 */
2883  int nvarsF, /**< number of variables in F */
2884  int nvarsR, /**< number of variables in R */
2885  int* gubconsGC1, /**< pointer to store GUBs in GC1(GNC1+GOC1) */
2886  int* gubconsGC2, /**< pointer to store GUBs in GC2 */
2887  int* gubconsGFC1, /**< pointer to store GUBs in GFC1(GNC1+GF) */
2888  int* gubconsGR, /**< pointer to store GUBs in GR */
2889  int* ngubconsGC1, /**< pointer to store number of GUBs in GC1(GNC1+GOC1) */
2890  int* ngubconsGC2, /**< pointer to store number of GUBs in GC2 */
2891  int* ngubconsGFC1, /**< pointer to store number of GUBs in GFC1(GNC1+GF) */
2892  int* ngubconsGR, /**< pointer to store number of GUBs in GR */
2893  int* ngubconscapexceed, /**< pointer to store number of GUBs with only capacity exceeding variables */
2894  int* maxgubvarssize /**< pointer to store the maximal size of GUB constraints */
2895  )
2896 {
2897 #if 0 /* not required */
2898  SORTKEYPAIR** sortkeypairsF;
2899 #endif
2900  SORTKEYPAIR** sortkeypairsGFC1;
2901  SORTKEYPAIR* sortkeypairsGFC1store;
2902  SCIP_Real* sortkeysC1;
2903  SCIP_Real* sortkeysC2;
2904  SCIP_Real* sortkeysR;
2905  int* nC1varsingubcons;
2906  int var;
2907  int gubconsidx;
2908  int varidx;
2909  int ngubconss;
2910  int ngubconsGOC1;
2911  int targetvar;
2912  int nvarsprocessed;
2913  int i;
2914  int j;
2915 
2916 #if GUBSPLITGNC1GUBS
2917  SCIP_Bool gubconswithF;
2918  int origngubconss;
2919  origngubconss = gubset->ngubconss;
2920 #endif
2921 
2922  assert(scip != NULL);
2923  assert(gubset != NULL);
2924  assert(solvals != NULL);
2925  assert(weights != NULL);
2926  assert(varsC1 != NULL);
2927  assert(varsC2 != NULL);
2928  assert(varsF != NULL);
2929  assert(varsR != NULL);
2930  assert(nvarsC1 > 0);
2931  assert(nvarsC2 >= 0);
2932  assert(nvarsF >= 0);
2933  assert(nvarsR >= 0);
2934  assert(gubconsGC1 != NULL);
2935  assert(gubconsGC2 != NULL);
2936  assert(gubconsGFC1 != NULL);
2937  assert(gubconsGR != NULL);
2938  assert(ngubconsGC1 != NULL);
2939  assert(ngubconsGC2 != NULL);
2940  assert(ngubconsGFC1 != NULL);
2941  assert(ngubconsGR != NULL);
2942  assert(maxgubvarssize != NULL);
2943 
2944  ngubconss = gubset->ngubconss;
2945  nvarsprocessed = 0;
2946  ngubconsGOC1 = 0;
2947 
2948  /* GUBs are categorized into different types according to the variables in volved
2949  * - GOC1: involves variables in C1 only -- no C2, R, F
2950  * - GNC1: involves variables in C1 and F (and R) -- no C2
2951  * - GF: involves variables in F (and R) only -- no C1, C2
2952  * - GC2: involves variables in C2 only -- no C1, R, F
2953  * - GR: involves variables in R only -- no C1, C2, F
2954  * which requires splitting GUBs in case they include variable in F and R.
2955  *
2956  * afterwards all GUBs (except GOC1 GUBs, which we do not need to lift) are sorted by a two level lifting sequence.
2957  * - first ordering level is: GFC1 (GNC1+GF), GC2, and GR.
2958  * - second ordering level is
2959  * GFC1: non-increasing number of variables in F and non-increasing max{x*_k : k in GFC1_j} in case of equality
2960  * GC2: non-increasing max{ a_k : k in GC2_j}; note that |GFC2_j| = 1
2961  * GR: non-increasing max{ a_k : k in GR_j}
2962  *
2963  * in additon, another GUB union, which is helpful for the lifting procedure, is formed
2964  * - GC1: GUBs of category GOC1 and GNC1
2965  * with second ordering level non-decreasing min{ a_k : k in GC1_j };
2966  * note that min{ a_k : k in GC1_j } always comes from the first variable in the GUB
2967  */
2968 
2969  /* allocates temporary memory */
2970  SCIP_CALL( SCIPallocBufferArray(scip, &sortkeysC1, nvarsC1) );
2971 #if 0 /* not required */
2972  SCIP_CALL( SCIPallocBufferArray(scip, &sortkeypairsF, nvarsF) );
2973 #endif
2974  SCIP_CALL( SCIPallocBufferArray(scip, &sortkeysC2, nvarsC2) );
2975  SCIP_CALL( SCIPallocBufferArray(scip, &sortkeysR, nvarsR) );
2976 
2977 
2978  /* to get the GUB lifting sequence, we first sort all variables in F, C2, and R
2979  * - F: non-increasing x*_j and non-increasing a_j in case of equality
2980  * - C2: non-increasing a_j
2981  * - R: non-increasing a_j
2982  * furthermore, sort C1 variables as needed for initializing the minweight table (non-increasing a_j).
2983  */
2984 
2985  /* gets sorting key for variables in C1 corresponding to the following ordering
2986  * non-decreasing a_j, i.e. a_1 <= a_2 <= ... <= a_|C_1|
2987  */
2988  for( j = 0; j < nvarsC1; j++ )
2989  {
2990  /* gets sortkeys */
2991  sortkeysC1[j] = (SCIP_Real) weights[varsC1[j]];
2992 
2993  /* update status of variable in its gub constraint */
2994  gubconsidx = gubset->gubconssidx[varsC1[j]];
2995  varidx = gubset->gubvarsidx[varsC1[j]];
2996  gubset->gubconss[gubconsidx]->gubvarsstatus[varidx] = GUBVARSTATUS_BELONGSTOSET_C1;
2997  }
2998 
2999  /* gets sorting key for variables in F corresponding to the following ordering
3000  * non-increasing x*_j, i.e., x*_1 >= x*_2 >= ... >= x*_|F|, and
3001  * non-increasing a_j, i.e., a_1 >= a_2 >= ... >= a_|F| in case of equality
3002  * and updates status of each variable in F in GUB set data structure
3003  */
3004  for( j = 0; j < nvarsF; j++ )
3005  {
3006 #if 0 /* not required */
3007  /* gets sortkeys */
3008  SCIP_CALL( SCIPallocBuffer(scip, &sortkeypairsF[j]) );
3009  sortkeypairsF[j]->key1 = solvals[varsF[j]];
3010  sortkeypairsF[j]->key2 = (SCIP_Real) weights[varsF[j]];
3011 #endif
3012 
3013  /* update status of variable in its gub constraint */
3014  gubconsidx = gubset->gubconssidx[varsF[j]];
3015  varidx = gubset->gubvarsidx[varsF[j]];
3016  gubset->gubconss[gubconsidx]->gubvarsstatus[varidx] = GUBVARSTATUS_BELONGSTOSET_F;
3017  }
3018 
3019  /* gets sorting key for variables in C2 corresponding to the following ordering
3020  * non-increasing a_j, i.e., a_1 >= a_2 >= ... >= a_|C2|
3021  * and updates status of each variable in F in GUB set data structure
3022  */
3023  for( j = 0; j < nvarsC2; j++ )
3024  {
3025  /* gets sortkeys */
3026  sortkeysC2[j] = (SCIP_Real) weights[varsC2[j]];
3027 
3028  /* update status of variable in its gub constraint */
3029  gubconsidx = gubset->gubconssidx[varsC2[j]];
3030  varidx = gubset->gubvarsidx[varsC2[j]];
3031  gubset->gubconss[gubconsidx]->gubvarsstatus[varidx] = GUBVARSTATUS_BELONGSTOSET_C2;
3032  }
3033 
3034  /* gets sorting key for variables in R corresponding to the following ordering
3035  * non-increasing a_j, i.e., a_1 >= a_2 >= ... >= a_|R|
3036  * and updates status of each variable in F in GUB set data structure
3037  */
3038  for( j = 0; j < nvarsR; j++ )
3039  {
3040  /* gets sortkeys */
3041  sortkeysR[j] = (SCIP_Real) weights[varsR[j]];
3042 
3043  /* update status of variable in its gub constraint */
3044  gubconsidx = gubset->gubconssidx[varsR[j]];
3045  varidx = gubset->gubvarsidx[varsR[j]];
3046  gubset->gubconss[gubconsidx]->gubvarsstatus[varidx] = GUBVARSTATUS_BELONGSTOSET_R;
3047  }
3048 
3049  /* sorts C1, F, C2 and R */
3050  if( nvarsC1 > 0 )
3051  {
3052  SCIPsortRealInt(sortkeysC1, varsC1, nvarsC1);
3053  }
3054 #if 0 /* not required */
3055  if( nvarsF > 0 )
3056  {
3057  SCIPsortDownPtrInt((void**)sortkeypairsF, varsF, compSortkeypairs, nvarsF);
3058  }
3059 #endif
3060  if( nvarsC2 > 0 )
3061  {
3062  SCIPsortDownRealInt(sortkeysC2, varsC2, nvarsC2);
3063  }
3064  if( nvarsR > 0)
3065  {
3066  SCIPsortDownRealInt(sortkeysR, varsR, nvarsR);
3067  }
3068 
3069  /* frees temporary memory */
3070  SCIPfreeBufferArray(scip, &sortkeysR);
3071  SCIPfreeBufferArray(scip, &sortkeysC2);
3072 #if 0 /* not required */
3073  for( j = nvarsF-1; j >= 0; j-- )
3074  SCIPfreeBuffer(scip, &sortkeypairsF[j]);
3075  SCIPfreeBufferArray(scip, &sortkeypairsF);
3076 #endif
3077  SCIPfreeBufferArray(scip, &sortkeysC1);
3078 
3079  /* allocate and initialize temporary memory for sorting GUB constraints */
3080  SCIP_CALL( SCIPallocBufferArray(scip, &sortkeypairsGFC1, ngubconss) );
3081  SCIP_CALL( SCIPallocBufferArray(scip, &sortkeypairsGFC1store, ngubconss) );
3082  SCIP_CALL( SCIPallocBufferArray(scip, &nC1varsingubcons, ngubconss) );
3083  BMSclearMemoryArray(nC1varsingubcons, ngubconss);
3084  for( i = 0; i < ngubconss; i++)
3085  {
3086  sortkeypairsGFC1[i] = &(sortkeypairsGFC1store[i]);
3087  sortkeypairsGFC1[i]->key1 = 0.0;
3088  sortkeypairsGFC1[i]->key2 = 0.0;
3089  }
3090  *ngubconsGC1 = 0;
3091  *ngubconsGC2 = 0;
3092  *ngubconsGFC1 = 0;
3093  *ngubconsGR = 0;
3094  *ngubconscapexceed = 0;
3095  *maxgubvarssize = 0;
3096 
3097 #ifndef NDEBUG
3098  for( i = 0; i < gubset->ngubconss; i++ )
3099  assert(gubset->gubconsstatus[i] == GUBCONSSTATUS_UNINITIAL);
3100 #endif
3101 
3102  /* stores GUBs of group GC1 (GOC1+GNC1) and part of the GUBs of group GFC1 (GNC1 GUBs) and sorts variables in these GUBs
3103  * s.t. C1 variables come first (will automatically be sorted by non-decreasing weight).
3104  * gets sorting keys for GUBs of type GFC1 corresponding to the following ordering
3105  * non-increasing number of variables in F, and
3106  * non-increasing max{x*_k : k in GFC1_j} in case of equality
3107  */
3108  for( i = 0; i < nvarsC1; i++ )
3109  {
3110  int nvarsC1capexceed;
3111 
3112  nvarsC1capexceed = 0;
3113 
3114  var = varsC1[i];
3115  gubconsidx = gubset->gubconssidx[var];
3116  varidx = gubset->gubvarsidx[var];
3117 
3118  assert(gubconsidx >= 0 && gubconsidx < ngubconss);
3119  assert(gubset->gubconss[gubconsidx]->gubvarsstatus[varidx] == GUBVARSTATUS_BELONGSTOSET_C1);
3120 
3121  /* current C1 variable is put to the front of its GUB where C1 part is stored by non-decreasing weigth;
3122  * note that variables in C1 are already sorted by non-decreasing weigth
3123  */
3124  targetvar = gubset->gubconss[gubconsidx]->gubvars[nC1varsingubcons[gubconsidx]];
3125  GUBsetSwapVars(scip, gubset, var, targetvar);
3126  nC1varsingubcons[gubconsidx]++;
3127 
3128  /* the GUB was already handled (status set and stored in its group) by another variable of the GUB */
3129  if( gubset->gubconsstatus[gubconsidx] != GUBCONSSTATUS_UNINITIAL )
3130  {
3131  assert(gubset->gubconsstatus[gubconsidx] == GUBCONSSTATUS_BELONGSTOSET_GOC1
3132  || gubset->gubconsstatus[gubconsidx] == GUBCONSSTATUS_BELONGSTOSET_GNC1);
3133  continue;
3134  }
3135 
3136  /* determine the status of the current GUB constraint, GOC1 or GNC1; GUBs involving R variables are split into
3137  * GOC1/GNC1 and GF, if wanted. also update sorting key if GUB is of type GFC1 (GNC1)
3138  */
3139 #if GUBSPLITGNC1GUBS
3140  gubconswithF = FALSE;
3141 #endif
3142  for( j = 0; j < gubset->gubconss[gubconsidx]->ngubvars; j++ )
3143  {
3144  assert(gubset->gubconss[gubconsidx]->gubvarsstatus[j] != GUBVARSTATUS_BELONGSTOSET_C2);
3145 
3146  /* C1-variable: update number of C1/capacity exceeding variables */
3147  if( gubset->gubconss[gubconsidx]->gubvarsstatus[j] == GUBVARSTATUS_BELONGSTOSET_C1 )
3148  {
3149  nvarsC1capexceed++;
3150  nvarsprocessed++;
3151  }
3152  /* F-variable: update sort key (number of F variables in GUB) of corresponding GFC1-GUB */
3153  else if( gubset->gubconss[gubconsidx]->gubvarsstatus[j] == GUBVARSTATUS_BELONGSTOSET_F )
3154  {
3155 #if GUBSPLITGNC1GUBS
3156  gubconswithF = TRUE;
3157 #endif
3158  sortkeypairsGFC1[*ngubconsGFC1]->key1 += 1.0;
3159 
3160  if( solvals[gubset->gubconss[gubconsidx]->gubvars[j]] > sortkeypairsGFC1[*ngubconsGFC1]->key2 )
3161  sortkeypairsGFC1[*ngubconsGFC1]->key2 = solvals[gubset->gubconss[gubconsidx]->gubvars[j]];
3162  }
3163  else if( gubset->gubconss[gubconsidx]->gubvarsstatus[j] == GUBVARSTATUS_CAPACITYEXCEEDED )
3164  {
3165  nvarsC1capexceed++;
3166  }
3167  else
3168  assert(gubset->gubconss[gubconsidx]->gubvarsstatus[j] == GUBVARSTATUS_BELONGSTOSET_R);
3169  }
3170 
3171  /* update set of GC1 GUBs */
3172  gubconsGC1[*ngubconsGC1] = gubconsidx;
3173  (*ngubconsGC1)++;
3174 
3175  /* update maximum size of all GUB constraints */
3176  if( gubset->gubconss[gubconsidx]->gubvarssize > *maxgubvarssize )
3177  *maxgubvarssize = gubset->gubconss[gubconsidx]->gubvarssize;
3178 
3179  /* set status of GC1-GUB (GOC1 or GNC1) and update set of GFC1 GUBs */
3180  if( nvarsC1capexceed == gubset->gubconss[gubconsidx]->ngubvars )
3181  {
3182  gubset->gubconsstatus[gubconsidx] = GUBCONSSTATUS_BELONGSTOSET_GOC1;
3183  ngubconsGOC1++;
3184  }
3185  else
3186  {
3187 #if GUBSPLITGNC1GUBS
3188  /* only variables in C1 and R -- no in F: GUB will be split into GR and GOC1 GUBs */
3189  if( !gubconswithF )
3190  {
3191  GUBVARSTATUS movevarstatus;
3192 
3193  assert(gubset->ngubconss < gubset->nvars);
3194 
3195  /* create a new GUB for GR part of splitting */
3196  SCIP_CALL( GUBconsCreate(scip, &gubset->gubconss[gubset->ngubconss]) );
3197  gubset->ngubconss++;
3198  ngubconss = gubset->ngubconss;
3199 
3200  /* fill GR with R variables in current GUB */
3201  for( j = gubset->gubconss[gubconsidx]->ngubvars-1; j >= 0; j-- )
3202  {
3203  movevarstatus = gubset->gubconss[gubconsidx]->gubvarsstatus[j];
3204  if( movevarstatus != GUBVARSTATUS_BELONGSTOSET_C1 )
3205  {
3206  assert(movevarstatus == GUBVARSTATUS_BELONGSTOSET_R || movevarstatus == GUBVARSTATUS_CAPACITYEXCEEDED);
3207  SCIP_CALL( GUBsetMoveVar(scip, gubset, vars, gubset->gubconss[gubconsidx]->gubvars[j],
3208  gubconsidx, ngubconss-1) );
3209  gubset->gubconss[ngubconss-1]->gubvarsstatus[gubset->gubconss[ngubconss-1]->ngubvars-1] =
3210  movevarstatus;
3211  }
3212  }
3213 
3214  gubset->gubconsstatus[gubconsidx] = GUBCONSSTATUS_BELONGSTOSET_GOC1;
3215  ngubconsGOC1++;
3216 
3217  gubset->gubconsstatus[ngubconss-1] = GUBCONSSTATUS_BELONGSTOSET_GR;
3218  gubconsGR[*ngubconsGR] = ngubconss-1;
3219  (*ngubconsGR)++;
3220  }
3221  /* variables in C1, F, and maybe R: GNC1 GUB */
3222  else
3223  {
3224  assert(gubconswithF);
3225 
3226  gubset->gubconsstatus[gubconsidx] = GUBCONSSTATUS_BELONGSTOSET_GNC1;
3227  gubconsGFC1[*ngubconsGFC1] = gubconsidx;
3228  (*ngubconsGFC1)++;
3229  }
3230 #else
3231  gubset->gubconsstatus[gubconsidx] = GUBCONSSTATUS_BELONGSTOSET_GNC1;
3232  gubconsGFC1[*ngubconsGFC1] = gubconsidx;
3233  (*ngubconsGFC1)++;
3234 #endif
3235  }
3236  }
3237 
3238  /* stores GUBs of group GC2 (only trivial GUBs); sorting is not required because the C2 variables (which we loop over)
3239  * are already sorted correctly
3240  */
3241  for( i = 0; i < nvarsC2; i++ )
3242  {
3243  var = varsC2[i];
3244  gubconsidx = gubset->gubconssidx[var];
3245  varidx = gubset->gubvarsidx[var];
3246 
3247  assert(gubconsidx >= 0 && gubconsidx < ngubconss);
3248  assert(gubset->gubconss[gubconsidx]->ngubvars == 1);
3249  assert(varidx == 0);
3250  assert(gubset->gubconss[gubconsidx]->gubvarsstatus[varidx] == GUBVARSTATUS_BELONGSTOSET_C2);
3251  assert(gubset->gubconsstatus[gubconsidx] == GUBCONSSTATUS_UNINITIAL);
3252 
3253  /* set status of GC2 GUB */
3254  gubset->gubconsstatus[gubconsidx] = GUBCONSSTATUS_BELONGSTOSET_GC2;
3255 
3256  /* update group of GC2 GUBs */
3257  gubconsGC2[*ngubconsGC2] = gubconsidx;
3258  (*ngubconsGC2)++;
3259 
3260  /* update maximum size of all GUB constraints */
3261  if( gubset->gubconss[gubconsidx]->gubvarssize > *maxgubvarssize )
3262  *maxgubvarssize = gubset->gubconss[gubconsidx]->gubvarssize;
3263 
3264  nvarsprocessed++;
3265  }
3266 
3267  /* stores remaining part of the GUBs of group GFC1 (GF GUBs) and gets GUB sorting keys corresp. to following ordering
3268  * non-increasing number of variables in F, and
3269  * non-increasing max{x*_k : k in GFC1_j} in case of equality
3270  */
3271  for( i = 0; i < nvarsF; i++ )
3272  {
3273  var = varsF[i];
3274  gubconsidx = gubset->gubconssidx[var];
3275  varidx = gubset->gubvarsidx[var];
3276 
3277  assert(gubconsidx >= 0 && gubconsidx < ngubconss);
3278  assert(gubset->gubconss[gubconsidx]->gubvarsstatus[varidx] == GUBVARSTATUS_BELONGSTOSET_F);
3279 
3280  nvarsprocessed++;
3281 
3282  /* the GUB was already handled (status set and stored in its group) by another variable of the GUB */
3283  if( gubset->gubconsstatus[gubconsidx] != GUBCONSSTATUS_UNINITIAL )
3284  {
3285  assert(gubset->gubconsstatus[gubconsidx] == GUBCONSSTATUS_BELONGSTOSET_GF
3286  || gubset->gubconsstatus[gubconsidx] == GUBCONSSTATUS_BELONGSTOSET_GNC1);
3287  continue;
3288  }
3289 
3290  /* set status of GF GUB */
3291  gubset->gubconsstatus[gubconsidx] = GUBCONSSTATUS_BELONGSTOSET_GF;
3292 
3293  /* update sorting key of corresponding GFC1 GUB */
3294  for( j = 0; j < gubset->gubconss[gubconsidx]->ngubvars; j++ )
3295  {
3296  assert(gubset->gubconss[gubconsidx]->gubvarsstatus[j] != GUBVARSTATUS_BELONGSTOSET_C2
3297  && gubset->gubconss[gubconsidx]->gubvarsstatus[j] != GUBVARSTATUS_BELONGSTOSET_C1);
3298 
3299  /* F-variable: update sort key (number of F variables in GUB) of corresponding GFC1-GUB */
3300  if( gubset->gubconss[gubconsidx]->gubvarsstatus[j] == GUBVARSTATUS_BELONGSTOSET_F )
3301  {
3302  sortkeypairsGFC1[*ngubconsGFC1]->key1 += 1.0;
3303 
3304  if( solvals[gubset->gubconss[gubconsidx]->gubvars[j]] > sortkeypairsGFC1[*ngubconsGFC1]->key2 )
3305  sortkeypairsGFC1[*ngubconsGFC1]->key2 = solvals[gubset->gubconss[gubconsidx]->gubvars[j]];
3306  }
3307  }
3308 
3309  /* update set of GFC1 GUBs */
3310  gubconsGFC1[*ngubconsGFC1] = gubconsidx;
3311  (*ngubconsGFC1)++;
3312 
3313  /* update maximum size of all GUB constraints */
3314  if( gubset->gubconss[gubconsidx]->gubvarssize > *maxgubvarssize )
3315  *maxgubvarssize = gubset->gubconss[gubconsidx]->gubvarssize;
3316  }
3317 
3318  /* stores GUBs of group GR; sorting is not required because the R variables (which we loop over) are already sorted
3319  * correctly
3320  */
3321  for( i = 0; i < nvarsR; i++ )
3322  {
3323  var = varsR[i];
3324  gubconsidx = gubset->gubconssidx[var];
3325  varidx = gubset->gubvarsidx[var];
3326 
3327  assert(gubconsidx >= 0 && gubconsidx < ngubconss);
3328  assert(gubset->gubconss[gubconsidx]->gubvarsstatus[varidx] == GUBVARSTATUS_BELONGSTOSET_R);
3329 
3330  nvarsprocessed++;
3331 
3332  /* the GUB was already handled (status set and stored in its group) by another variable of the GUB */
3333  if( gubset->gubconsstatus[gubconsidx] != GUBCONSSTATUS_UNINITIAL )
3334  {
3335  assert(gubset->gubconsstatus[gubconsidx] == GUBCONSSTATUS_BELONGSTOSET_GR
3336  || gubset->gubconsstatus[gubconsidx] == GUBCONSSTATUS_BELONGSTOSET_GF
3337  || gubset->gubconsstatus[gubconsidx] == GUBCONSSTATUS_BELONGSTOSET_GNC1);
3338  continue;
3339  }
3340 
3341  /* set status of GR GUB */
3342  gubset->gubconsstatus[gubconsidx] = GUBCONSSTATUS_BELONGSTOSET_GR;
3343 
3344  /* update set of GR GUBs */
3345  gubconsGR[*ngubconsGR] = gubconsidx;
3346  (*ngubconsGR)++;
3347 
3348  /* update maximum size of all GUB constraints */
3349  if( gubset->gubconss[gubconsidx]->gubvarssize > *maxgubvarssize )
3350  *maxgubvarssize = gubset->gubconss[gubconsidx]->gubvarssize;
3351  }
3352  assert(nvarsprocessed == nvarsC1 + nvarsC2 + nvarsF + nvarsR);
3353 
3354  /* update number of GUBs with only capacity exceeding variables (will not be used for lifting) */
3355  (*ngubconscapexceed) = ngubconss - (ngubconsGOC1 + (*ngubconsGC2) + (*ngubconsGFC1) + (*ngubconsGR));
3356  assert(*ngubconscapexceed >= 0);
3357 #ifndef NDEBUG
3358  {
3359  int check;
3360 
3361  check = 0;
3362 
3363  /* remaining not handled GUBs should only contain capacity exceeding variables */
3364  for( i = 0; i < ngubconss; i++ )
3365  {
3366  if( gubset->gubconsstatus[i] == GUBCONSSTATUS_UNINITIAL )
3367  check++;
3368  }
3369  assert(check == *ngubconscapexceed);
3370  }
3371 #endif
3372 
3373  /* sort GFCI GUBs according to computed sorting keys */
3374  if( (*ngubconsGFC1) > 0 )
3375  {
3376  SCIPsortDownPtrInt((void**)sortkeypairsGFC1, gubconsGFC1, compSortkeypairs, (*ngubconsGFC1));
3377  }
3378 
3379  /* free temporary memory */
3380 #if GUBSPLITGNC1GUBS
3381  ngubconss = origngubconss;
3382 #endif
3383  SCIPfreeBufferArray(scip, &nC1varsingubcons);
3384  SCIPfreeBufferArray(scip, &sortkeypairsGFC1store);
3385  SCIPfreeBufferArray(scip, &sortkeypairsGFC1);
3386 
3387  return SCIP_OKAY;
3388 }
3389 
3390 /** enlarges minweight table to at least the given length */
3391 static
3393  SCIP* scip, /**< SCIP data structure */
3394  SCIP_Longint** minweightsptr, /**< pointer to minweights table */
3395  int* minweightslen, /**< pointer to store number of entries in minweights table (incl. z=0) */
3396  int* minweightssize, /**< pointer to current size of minweights table */
3397  int newlen /**< new length of minweights table */
3398  )
3399 {
3400  int j;
3401 
3402  assert(minweightsptr != NULL);
3403  assert(*minweightsptr != NULL);
3404  assert(minweightslen != NULL);
3405  assert(*minweightslen >= 0);
3406  assert(minweightssize != NULL);
3407  assert(*minweightssize >= 0);
3408 
3409  if( newlen > *minweightssize )
3410  {
3411  int newsize;
3412 
3413  /* reallocate table memory */
3414  newsize = SCIPcalcMemGrowSize(scip, newlen);
3415  SCIP_CALL( SCIPreallocBufferArray(scip, minweightsptr, newsize) );
3416  *minweightssize = newsize;
3417  }
3418  assert(newlen <= *minweightssize);
3419 
3420  /* initialize new elements */
3421  for( j = *minweightslen; j < newlen; ++j )
3422  (*minweightsptr)[j] = SCIP_LONGINT_MAX;
3423  *minweightslen = newlen;
3424 
3425  return SCIP_OKAY;
3426 }
3427 
3428 /** lifts given inequality
3429  * sum_{j in M_1} x_j <= alpha_0
3430  * valid for
3431  * S^0 = { x in {0,1}^|M_1| : sum_{j in M_1} a_j x_j <= a_0 - sum_{j in M_2} a_j }
3432  * to a valid inequality
3433  * sum_{j in M_1} x_j + sum_{j in F} alpha_j x_j + sum_{j in M_2} alpha_j x_j + sum_{j in R} alpha_j x_j
3434  * <= alpha_0 + sum_{j in M_2} alpha_j
3435  * for
3436  * S = { x in {0,1}^|N| : sum_{j in N} a_j x_j <= a_0 };
3437  * uses sequential up-lifting for the variables in F, sequential down-lifting for the variable in M_2, and
3438  * sequential up-lifting for the variables in R; procedure can be used to strengthen minimal cover inequalities and
3439  * extended weight inequalities.
3440  */
3441 static
3443  SCIP* scip, /**< SCIP data structure */
3444  SCIP_VAR** vars, /**< variables in knapsack constraint */
3445  int nvars, /**< number of variables in knapsack constraint */
3446  int ntightened, /**< number of variables with tightened upper bound */
3447  SCIP_Longint* weights, /**< weights of variables in knapsack constraint */
3448  SCIP_Longint capacity, /**< capacity of knapsack */
3449  SCIP_Real* solvals, /**< solution values of all problem variables */
3450  int* varsM1, /**< variables in M_1 */
3451  int* varsM2, /**< variables in M_2 */
3452  int* varsF, /**< variables in F */
3453  int* varsR, /**< variables in R */
3454  int nvarsM1, /**< number of variables in M_1 */
3455  int nvarsM2, /**< number of variables in M_2 */
3456  int nvarsF, /**< number of variables in F */
3457  int nvarsR, /**< number of variables in R */
3458  int alpha0, /**< rights hand side of given valid inequality */
3459  int* liftcoefs, /**< pointer to store lifting coefficient of vars in knapsack constraint */
3460  SCIP_Real* cutact, /**< pointer to store activity of lifted valid inequality */
3461  int* liftrhs /**< pointer to store right hand side of the lifted valid inequality */
3462  )
3463 {
3464  SCIP_Longint* minweights;
3465  SCIP_Real* sortkeys;
3466  SCIP_Longint fixedonesweight;
3467  int minweightssize;
3468  int minweightslen;
3469  int j;
3470  int w;
3471 
3472  assert(scip != NULL);
3473  assert(vars != NULL);
3474  assert(nvars >= 0);
3475  assert(weights != NULL);
3476  assert(capacity >= 0);
3477  assert(solvals != NULL);
3478  assert(varsM1 != NULL);
3479  assert(varsM2 != NULL);
3480  assert(varsF != NULL);
3481  assert(varsR != NULL);
3482  assert(nvarsM1 >= 0 && nvarsM1 <= nvars - ntightened);
3483  assert(nvarsM2 >= 0 && nvarsM2 <= nvars - ntightened);
3484  assert(nvarsF >= 0 && nvarsF <= nvars - ntightened);
3485  assert(nvarsR >= 0 && nvarsR <= nvars - ntightened);
3486  assert(nvarsM1 + nvarsM2 + nvarsF + nvarsR == nvars - ntightened);
3487  assert(alpha0 >= 0);
3488  assert(liftcoefs != NULL);
3489  assert(cutact != NULL);
3490  assert(liftrhs != NULL);
3491 
3492  /* allocates temporary memory */
3493  minweightssize = nvarsM1 + 1;
3494  SCIP_CALL( SCIPallocBufferArray(scip, &minweights, minweightssize) );
3495  SCIP_CALL( SCIPallocBufferArray(scip, &sortkeys, nvarsM1) );
3496 
3497  /* initializes data structures */
3498  BMSclearMemoryArray(liftcoefs, nvars);
3499  *cutact = 0.0;
3500 
3501  /* sets lifting coefficient of variables in M1, sorts variables in M1 such that a_1 <= a_2 <= ... <= a_|M1|
3502  * and calculates activity of the current valid inequality
3503  */
3504  for( j = 0; j < nvarsM1; j++ )
3505  {
3506  assert(liftcoefs[varsM1[j]] == 0);
3507  liftcoefs[varsM1[j]] = 1;
3508  sortkeys[j] = (SCIP_Real) (weights[varsM1[j]]);
3509  (*cutact) += solvals[varsM1[j]];
3510  }
3511 
3512  SCIPsortRealInt(sortkeys, varsM1, nvarsM1);
3513 
3514  /* initializes (i = 1) the minweight table, defined as: minweights_i[w] =
3515  * min sum_{j in M_1} a_j x_j + sum_{k=1}^{i-1} a_{j_k} x_{j_k}
3516  * s.t. sum_{j in M_1} x_j + sum_{k=1}^{i-1} alpha_{j_k} x_{j_k} >= w
3517  * x_j in {0,1} for j in M_1 & {j_i,...,j_i-1},
3518  * for i = 1,...,t with t = |N\M1| and w = 0,...,|M1| + sum_{k=1}^{i-1} alpha_{j_k};
3519  */
3520  minweights[0] = 0;
3521  for( w = 1; w <= nvarsM1; w++ )
3522  minweights[w] = minweights[w-1] + weights[varsM1[w-1]];
3523  minweightslen = nvarsM1 + 1;
3524 
3525  /* gets sum of weights of variables fixed to one, i.e. sum of weights of variables in M_2 */
3526  fixedonesweight = 0;
3527  for( j = 0; j < nvarsM2; j++ )
3528  fixedonesweight += weights[varsM2[j]];
3529  assert(fixedonesweight >= 0);
3530 
3531  /* initializes right hand side of lifted valid inequality */
3532  *liftrhs = alpha0;
3533 
3534  /* sequentially up-lifts all variables in F: */
3535  for( j = 0; j < nvarsF; j++ )
3536  {
3537  SCIP_Longint weight;
3538  int liftvar;
3539  int liftcoef;
3540  int z;
3541 
3542  liftvar = varsF[j];
3543  weight = weights[liftvar];
3544  assert(liftvar >= 0 && liftvar < nvars);
3545  assert(SCIPisFeasGT(scip, solvals[liftvar], 0.0));
3546  assert(weight > 0);
3547 
3548  /* knapsack problem is infeasible:
3549  * sets z = 0
3550  */
3551  if( capacity - fixedonesweight - weight < 0 )
3552  {
3553  z = 0;
3554  }
3555  /* knapsack problem is feasible:
3556  * sets z = max { w : 0 <= w <= liftrhs, minweights_i[w] <= a_0 - fixedonesweight - a_{j_i} } = liftrhs,
3557  * if minweights_i[liftrhs] <= a_0 - fixedonesweight - a_{j_i}
3558  */
3559  else if( minweights[*liftrhs] <= capacity - fixedonesweight - weight )
3560  {
3561  z = *liftrhs;
3562  }
3563  /* knapsack problem is feasible:
3564  * uses binary search to find z = max { w : 0 <= w <= liftrhs, minweights_i[w] <= a_0 - fixedonesweight - a_{j_i} }
3565  */
3566  else
3567  {
3568  int left;
3569  int right;
3570  int middle;
3571 
3572  assert((*liftrhs) + 1 >= minweightslen || minweights[(*liftrhs) + 1] > capacity - fixedonesweight - weight);
3573  left = 0;
3574  right = (*liftrhs) + 1;
3575  while( left < right - 1 )
3576  {
3577  middle = (left + right) / 2;
3578  assert(0 <= middle && middle < minweightslen);
3579  if( minweights[middle] <= capacity - fixedonesweight - weight )
3580  left = middle;
3581  else
3582  right = middle;
3583  }
3584  assert(left == right - 1);
3585  assert(0 <= left && left < minweightslen);
3586  assert(minweights[left] <= capacity - fixedonesweight - weight );
3587  assert(left == minweightslen - 1 || minweights[left+1] > capacity - fixedonesweight - weight);
3588 
3589  /* now z = left */
3590  z = left;
3591  assert(z <= *liftrhs);
3592  }
3593 
3594  /* calculates lifting coefficients alpha_{j_i} = liftrhs - z */
3595  liftcoef = (*liftrhs) - z;
3596  liftcoefs[liftvar] = liftcoef;
3597  assert(liftcoef >= 0 && liftcoef <= (*liftrhs) + 1);
3598 
3599  /* minweight table and activity of current valid inequality will not change, if alpha_{j_i} = 0 */
3600  if( liftcoef == 0 )
3601  continue;
3602 
3603  /* updates activity of current valid inequality */
3604  (*cutact) += liftcoef * solvals[liftvar];
3605 
3606  /* enlarges current minweight table:
3607  * from minweightlen = |M1| + sum_{k=1}^{i-1} alpha_{j_k} + 1 entries
3608  * to |M1| + sum_{k=1}^{i } alpha_{j_k} + 1 entries
3609  * and sets minweights_i[w] = infinity for
3610  * w = |M1| + sum_{k=1}^{i-1} alpha_{j_k} + 1 , ... , |M1| + sum_{k=1}^{i} alpha_{j_k}
3611  */
3612  SCIP_CALL( enlargeMinweights(scip, &minweights, &minweightslen, &minweightssize, minweightslen + liftcoef) );
3613 
3614  /* updates minweight table: minweight_i+1[w] =
3615  * min{ minweights_i[w], a_{j_i}}, if w < alpha_j_i
3616  * min{ minweights_i[w], minweights_i[w - alpha_j_i] + a_j_i}, if w >= alpha_j_i
3617  */
3618  for( w = minweightslen - 1; w >= 0; w-- )
3619  {
3620  SCIP_Longint min;
3621  if( w < liftcoef )
3622  {
3623  min = MIN(minweights[w], weight);
3624  minweights[w] = min;
3625  }
3626  else
3627  {
3628  assert(w >= liftcoef);
3629  min = MIN(minweights[w], minweights[w - liftcoef] + weight);
3630  minweights[w] = min;
3631  }
3632  }
3633  }
3634  assert(minweights[0] == 0);
3635 
3636  /* sequentially down-lifts all variables in M_2: */
3637  for( j = 0; j < nvarsM2; j++ )
3638  {
3639  SCIP_Longint weight;
3640  int liftvar;
3641  int liftcoef;
3642  int left;
3643  int right;
3644  int middle;
3645  int z;
3646 
3647  liftvar = varsM2[j];
3648  weight = weights[liftvar];
3649  assert(SCIPisFeasEQ(scip, solvals[liftvar], 1.0));
3650  assert(liftvar >= 0 && liftvar < nvars);
3651  assert(weight > 0);
3652 
3653  /* uses binary search to find
3654  * z = max { w : 0 <= w <= |M_1| + sum_{k=1}^{i-1} alpha_{j_k}, minweights_[w] <= a_0 - fixedonesweight + a_{j_i}}
3655  */
3656  left = 0;
3657  right = minweightslen;
3658  while( left < right - 1 )
3659  {
3660  middle = (left + right) / 2;
3661  assert(0 <= middle && middle < minweightslen);
3662  if( minweights[middle] <= capacity - fixedonesweight + weight )
3663  left = middle;
3664  else
3665  right = middle;
3666  }
3667  assert(left == right - 1);
3668  assert(0 <= left && left < minweightslen);
3669  assert(minweights[left] <= capacity - fixedonesweight + weight );
3670  assert(left == minweightslen - 1 || minweights[left+1] > capacity - fixedonesweight + weight);
3671 
3672  /* now z = left */
3673  z = left;
3674  assert(z >= *liftrhs);
3675 
3676  /* calculates lifting coefficients alpha_{j_i} = z - liftrhs */
3677  liftcoef = z - (*liftrhs);
3678  liftcoefs[liftvar] = liftcoef;
3679  assert(liftcoef >= 0);
3680 
3681  /* updates sum of weights of variables fixed to one */
3682  fixedonesweight -= weight;
3683 
3684  /* updates right-hand side of current valid inequality */
3685  (*liftrhs) += liftcoef;
3686  assert(*liftrhs >= alpha0);
3687 
3688  /* minweight table and activity of current valid inequality will not change, if alpha_{j_i} = 0 */
3689  if( liftcoef == 0 )
3690  continue;
3691 
3692  /* updates activity of current valid inequality */
3693  (*cutact) += liftcoef * solvals[liftvar];
3694 
3695  /* enlarges current minweight table:
3696  * from minweightlen = |M1| + sum_{k=1}^{i-1} alpha_{j_k} + 1 entries
3697  * to |M1| + sum_{k=1}^{i } alpha_{j_k} + 1 entries
3698  * and sets minweights_i[w] = infinity for
3699  * w = |M1| + sum_{k=1}^{i-1} alpha_{j_k} + 1 , ... , |M1| + sum_{k=1}^{i} alpha_{j_k}
3700  */
3701  SCIP_CALL( enlargeMinweights(scip, &minweights, &minweightslen, &minweightssize, minweightslen + liftcoef) );
3702 
3703  /* updates minweight table: minweight_i+1[w] =
3704  * min{ minweights_i[w], a_{j_i}}, if w < alpha_j_i
3705  * min{ minweights_i[w], minweights_i[w - alpha_j_i] + a_j_i}, if w >= alpha_j_i
3706  */
3707  for( w = minweightslen - 1; w >= 0; w-- )
3708  {
3709  SCIP_Longint min;
3710  if( w < liftcoef )
3711  {
3712  min = MIN(minweights[w], weight);
3713  minweights[w] = min;
3714  }
3715  else
3716  {
3717  assert(w >= liftcoef);
3718  min = MIN(minweights[w], minweights[w - liftcoef] + weight);
3719  minweights[w] = min;
3720  }
3721  }
3722  }
3723  assert(fixedonesweight == 0);
3724  assert(*liftrhs >= alpha0);
3725 
3726  /* sequentially up-lifts all variables in R: */
3727  for( j = 0; j < nvarsR; j++ )
3728  {
3729  SCIP_Longint weight;
3730  int liftvar;
3731  int liftcoef;
3732  int z;
3733 
3734  liftvar = varsR[j];
3735  weight = weights[liftvar];
3736  assert(liftvar >= 0 && liftvar < nvars);
3737  assert(SCIPisFeasEQ(scip, solvals[liftvar], 0.0));
3738  assert(weight > 0);
3739  assert(capacity - weight >= 0);
3740  assert((*liftrhs) + 1 >= minweightslen || minweights[(*liftrhs) + 1] > capacity - weight);
3741 
3742  /* sets z = max { w : 0 <= w <= liftrhs, minweights_i[w] <= a_0 - a_{j_i} } = liftrhs,
3743  * if minweights_i[liftrhs] <= a_0 - a_{j_i}
3744  */
3745  if( minweights[*liftrhs] <= capacity - weight )
3746  {
3747  z = *liftrhs;
3748  }
3749  /* uses binary search to find z = max { w : 0 <= w <= liftrhs, minweights_i[w] <= a_0 - a_{j_i} }
3750  */
3751  else
3752  {
3753  int left;
3754  int right;
3755  int middle;
3756 
3757  left = 0;
3758  right = (*liftrhs) + 1;
3759  while( left < right - 1)
3760  {
3761  middle = (left + right) / 2;
3762  assert(0 <= middle && middle < minweightslen);
3763  if( minweights[middle] <= capacity - weight )
3764  left = middle;
3765  else
3766  right = middle;
3767  }
3768  assert(left == right - 1);
3769  assert(0 <= left && left < minweightslen);
3770  assert(minweights[left] <= capacity - weight );
3771  assert(left == minweightslen - 1 || minweights[left+1] > capacity - weight);
3772 
3773  /* now z = left */
3774  z = left;
3775  assert(z <= *liftrhs);
3776  }
3777 
3778  /* calculates lifting coefficients alpha_{j_i} = liftrhs - z */
3779  liftcoef = (*liftrhs) - z;
3780  liftcoefs[liftvar] = liftcoef;
3781  assert(liftcoef >= 0 && liftcoef <= *liftrhs);
3782 
3783  /* minweight table and activity of current valid inequality will not change, if alpha_{j_i} = 0 */
3784  if( liftcoef == 0 )
3785  continue;
3786 
3787  /* updates activity of current valid inequality */
3788  (*cutact) += liftcoef * solvals[liftvar];
3789 
3790  /* updates minweight table: minweight_i+1[w] =
3791  * min{ minweight_i[w], a_{j_i}}, if w < alpha_j_i
3792  * min{ minweight_i[w], minweight_i[w - alpha_j_i] + a_j_i}, if w >= alpha_j_i
3793  */
3794  for( w = *liftrhs; w >= 0; w-- )
3795  {
3796  SCIP_Longint min;
3797  if( w < liftcoef )
3798  {
3799  min = MIN(minweights[w], weight);
3800  minweights[w] = min;
3801  }
3802  else
3803  {
3804  assert(w >= liftcoef);
3805  min = MIN(minweights[w], minweights[w - liftcoef] + weight);
3806  minweights[w] = min;
3807  }
3808  }
3809  }
3810 
3811  /* frees temporary memory */
3812  SCIPfreeBufferArray(scip, &sortkeys);
3813  SCIPfreeBufferArray(scip, &minweights);
3814 
3815  return SCIP_OKAY;
3816 }
3817 
3818 /** adds two minweight values in a safe way, i.e,, ensures no overflow */
3819 static
3821  SCIP_Longint val1, /**< first value to add */
3822  SCIP_Longint val2 /**< second value to add */
3823  )
3824 {
3825  assert(val1 >= 0);
3826  assert(val2 >= 0);
3827 
3828  if( val1 >= SCIP_LONGINT_MAX || val2 >= SCIP_LONGINT_MAX )
3829  return SCIP_LONGINT_MAX;
3830  else
3831  {
3832  assert(val1 <= SCIP_LONGINT_MAX - val2);
3833  return (val1 + val2);
3834  }
3835 }
3836 
3837 /** computes minweights table for lifting with GUBs by combining unfished and fished tables */
3838 static
3840  SCIP_Longint* minweights, /**< minweight table to compute */
3841  SCIP_Longint* finished, /**< given finished table */
3842  SCIP_Longint* unfinished, /**< given unfinished table */
3843  int minweightslen /**< length of minweight, finished, and unfinished tables */
3844  )
3845 {
3846  int w1;
3847  int w2;
3848 
3849  /* minweights_i[w] = min{finished_i[w1] + unfinished_i[w2] : w1>=0, w2>=0, w1+w2=w};
3850  * note that finished and unfished arrays sorted by non-decreasing weight
3851  */
3852 
3853  /* initialize minweight with w2 = 0 */
3854  w2 = 0;
3855  assert(unfinished[w2] == 0);
3856  for( w1 = 0; w1 < minweightslen; w1++ )
3857  minweights[w1] = finished[w1];
3858 
3859  /* consider w2 = 1, ..., minweightslen-1 */
3860  for( w2 = 1; w2 < minweightslen; w2++ )
3861  {
3862  if( unfinished[w2] >= SCIP_LONGINT_MAX )
3863  break;
3864 
3865  for( w1 = 0; w1 < minweightslen - w2; w1++ )
3866  {
3867  SCIP_Longint temp;
3868 
3869  temp = safeAddMinweightsGUB(finished[w1], unfinished[w2]);
3870  if( temp <= minweights[w1+w2] )
3871  minweights[w1+w2] = temp;
3872  }
3873  }
3874 }
3875 
3876 /** lifts given inequality
3877  * sum_{j in C_1} x_j <= alpha_0
3878  * valid for
3879  * S^0 = { x in {0,1}^|C_1| : sum_{j in C_1} a_j x_j <= a_0 - sum_{j in C_2} a_j;
3880  * sum_{j in Q_i} x_j <= 1, forall i in I }
3881  * to a valid inequality
3882  * sum_{j in C_1} x_j + sum_{j in F} alpha_j x_j + sum_{j in C_2} alpha_j x_j + sum_{j in R} alpha_j x_j
3883  * <= alpha_0 + sum_{j in C_2} alpha_j
3884  * for
3885  * S = { x in {0,1}^|N| : sum_{j in N} a_j x_j <= a_0; sum_{j in Q_i} x_j <= 1, forall i in I };
3886  * uses sequential up-lifting for the variables in GUB constraints in gubconsGFC1,
3887  * sequential down-lifting for the variables in GUB constraints in gubconsGC2, and
3888  * sequential up-lifting for the variabels in GUB constraints in gubconsGR.
3889  */
3890 static
3892  SCIP* scip, /**< SCIP data structure */
3893  SCIP_GUBSET* gubset, /**< GUB set data structure */
3894  SCIP_VAR** vars, /**< variables in knapsack constraint */
3895  int ngubconscapexceed, /**< number of GUBs with only capacity exceeding variables */
3896  SCIP_Longint* weights, /**< weights of variables in knapsack constraint */
3897  SCIP_Longint capacity, /**< capacity of knapsack */
3898  SCIP_Real* solvals, /**< solution values of all knapsack variables */
3899  int* gubconsGC1, /**< GUBs in GC1(GNC1+GOC1) */
3900  int* gubconsGC2, /**< GUBs in GC2 */
3901  int* gubconsGFC1, /**< GUBs in GFC1(GNC1+GF) */
3902  int* gubconsGR, /**< GUBs in GR */
3903  int ngubconsGC1, /**< number of GUBs in GC1(GNC1+GOC1) */
3904  int ngubconsGC2, /**< number of GUBs in GC2 */
3905  int ngubconsGFC1, /**< number of GUBs in GFC1(GNC1+GF) */
3906  int ngubconsGR, /**< number of GUBs in GR */
3907  int alpha0, /**< rights hand side of given valid inequality */
3908  int* liftcoefs, /**< pointer to store lifting coefficient of vars in knapsack constraint */
3909  SCIP_Real* cutact, /**< pointer to store activity of lifted valid inequality */
3910  int* liftrhs, /**< pointer to store right hand side of the lifted valid inequality */
3911  int maxgubvarssize /**< maximal size of GUB constraints */
3912  )
3913 {
3914  SCIP_Longint* minweights;
3915  SCIP_Longint* finished;
3916  SCIP_Longint* unfinished;
3917  int* gubconsGOC1;
3918  int* gubconsGNC1;
3919  int* liftgubvars;
3920  SCIP_Longint fixedonesweight;
3921  SCIP_Longint weight;
3922  SCIP_Longint weightdiff1;
3923  SCIP_Longint weightdiff2;
3924  SCIP_Longint min;
3925  int minweightssize;
3926  int minweightslen;
3927  int nvars;
3928  int varidx;
3929  int liftgubconsidx;
3930  int liftvar;
3931  int sumliftcoef;
3932  int liftcoef;
3933  int ngubconsGOC1;
3934  int ngubconsGNC1;
3935  int left;
3936  int right;
3937  int middle;
3938  int nliftgubvars;
3939  int tmplen;
3940  int tmpsize;
3941  int j;
3942  int k;
3943  int w;
3944  int z;
3945 #ifndef NDEBUG
3946  int ngubconss;
3947  int nliftgubC1;
3948 
3949  assert(gubset != NULL);
3950  ngubconss = gubset->ngubconss;
3951 #else
3952  assert(gubset != NULL);
3953 #endif
3954 
3955  nvars = gubset->nvars;
3956 
3957  assert(scip != NULL);
3958  assert(vars != NULL);
3959  assert(nvars >= 0);
3960  assert(weights != NULL);
3961  assert(capacity >= 0);
3962  assert(solvals != NULL);
3963  assert(gubconsGC1 != NULL);
3964  assert(gubconsGC2 != NULL);
3965  assert(gubconsGFC1 != NULL);
3966  assert(gubconsGR != NULL);
3967  assert(ngubconsGC1 >= 0 && ngubconsGC1 <= ngubconss - ngubconscapexceed);
3968  assert(ngubconsGC2 >= 0 && ngubconsGC2 <= ngubconss - ngubconscapexceed);
3969  assert(ngubconsGFC1 >= 0 && ngubconsGFC1 <= ngubconss - ngubconscapexceed);
3970  assert(ngubconsGR >= 0 && ngubconsGR <= ngubconss - ngubconscapexceed);
3971  assert(alpha0 >= 0);
3972  assert(liftcoefs != NULL);
3973  assert(cutact != NULL);
3974  assert(liftrhs != NULL);
3975 
3976  minweightssize = ngubconsGC1+1;
3977 
3978  /* allocates temporary memory */
3979  SCIP_CALL( SCIPallocBufferArray(scip, &liftgubvars, maxgubvarssize) );
3980  SCIP_CALL( SCIPallocBufferArray(scip, &gubconsGOC1, ngubconsGC1) );
3981  SCIP_CALL( SCIPallocBufferArray(scip, &gubconsGNC1, ngubconsGC1) );
3982  SCIP_CALL( SCIPallocBufferArray(scip, &minweights, minweightssize) );
3983  SCIP_CALL( SCIPallocBufferArray(scip, &finished, minweightssize) );
3984  SCIP_CALL( SCIPallocBufferArray(scip, &unfinished, minweightssize) );
3985 
3986  /* initializes data structures */
3987  BMSclearMemoryArray(liftcoefs, nvars);
3988  *cutact = 0.0;
3989 
3990  /* gets GOC1 and GNC1 GUBs, sets lifting coefficient of variables in C1 and calculates activity of the current
3991  * valid inequality
3992  */
3993  ngubconsGOC1 = 0;
3994  ngubconsGNC1 = 0;
3995  for( j = 0; j < ngubconsGC1; j++ )
3996  {
3997  if( gubset->gubconsstatus[gubconsGC1[j]] == GUBCONSSTATUS_BELONGSTOSET_GOC1 )
3998  {
3999  gubconsGOC1[ngubconsGOC1] = gubconsGC1[j];
4000  ngubconsGOC1++;
4001  }
4002  else
4003  {
4004  assert(gubset->gubconsstatus[gubconsGC1[j]] == GUBCONSSTATUS_BELONGSTOSET_GNC1);
4005  gubconsGNC1[ngubconsGNC1] = gubconsGC1[j];
4006  ngubconsGNC1++;
4007  }
4008  for( k = 0; k < gubset->gubconss[gubconsGC1[j]]->ngubvars
4009  && gubset->gubconss[gubconsGC1[j]]->gubvarsstatus[k] == GUBVARSTATUS_BELONGSTOSET_C1; k++ )
4010  {
4011  varidx = gubset->gubconss[gubconsGC1[j]]->gubvars[k];
4012  assert(varidx >= 0 && varidx < nvars);
4013  assert(liftcoefs[varidx] == 0);
4014 
4015  liftcoefs[varidx] = 1;
4016  (*cutact) += solvals[varidx];
4017  }
4018  assert(k >= 1);
4019  }
4020  assert(ngubconsGOC1 + ngubconsGFC1 + ngubconsGC2 + ngubconsGR == ngubconss - ngubconscapexceed);
4021  assert(ngubconsGOC1 + ngubconsGNC1 == ngubconsGC1);
4022 
4023  /* initialize the minweight tables, defined as: for i = 1,...,m with m = |I| and w = 0,...,|gubconsGC1|;
4024  * - finished_i[w] =
4025  * min sum_{k = 1,2,...,i-1} sum_{j in Q_k} a_j x_j
4026  * s.t. sum_{k = 1,2,...,i-1} sum_{j in Q_k} alpha_j x_j >= w
4027  * sum_{j in Q_k} x_j <= 1
4028  * x_j in {0,1} forall j in Q_k forall k = 1,2,...,i-1,
4029  * - unfinished_i[w] =
4030  * min sum_{k = i+1,...,m} sum_{j in Q_k && j in C1} a_j x_j
4031  * s.t. sum_{k = i+1,...,m} sum_{j in Q_k && j in C1} x_j >= w
4032  * sum_{j in Q_k} x_j <= 1
4033  * x_j in {0,1} forall j in Q_k forall k = 1,2,...,i-1,
4034  * - minweights_i[w] = min{finished_i[w1] + unfinished_i[w2] : w1>=0, w2>=0, w1+w2=w};
4035  */
4036 
4037  /* initialize finished table; note that variables in GOC1 GUBs (includes C1 and capacity exceeding variables)
4038  * are sorted s.t. C1 variables come first and are sorted by non-decreasing weight.
4039  * GUBs in the group GCI are sorted by non-decreasing min{ a_k : k in GC1_j } where min{ a_k : k in GC1_j } always
4040  * comes from the first variable in the GUB
4041  */
4042  assert(ngubconsGOC1 <= ngubconsGC1);
4043  finished[0] = 0;
4044  for( w = 1; w <= ngubconsGOC1; w++ )
4045  {
4046  liftgubconsidx = gubconsGOC1[w-1];
4047 
4048  assert(gubset->gubconsstatus[liftgubconsidx] == GUBCONSSTATUS_BELONGSTOSET_GOC1);
4049  assert(gubset->gubconss[liftgubconsidx]->gubvarsstatus[0] == GUBVARSTATUS_BELONGSTOSET_C1);
4050 
4051  varidx = gubset->gubconss[liftgubconsidx]->gubvars[0];
4052 
4053  assert(varidx >= 0 && varidx < nvars);
4054  assert(liftcoefs[varidx] == 1);
4055 
4056  min = weights[varidx];
4057  finished[w] = finished[w-1] + min;
4058 
4059 #ifndef NDEBUG
4060  for( k = 1; k < gubset->gubconss[liftgubconsidx]->ngubvars
4061  && gubset->gubconss[liftgubconsidx]->gubvarsstatus[k] == GUBVARSTATUS_BELONGSTOSET_C1; k++ )
4062  {
4063  varidx = gubset->gubconss[liftgubconsidx]->gubvars[k];
4064  assert(varidx >= 0 && varidx < nvars);
4065  assert(liftcoefs[varidx] == 1);
4066  assert(weights[varidx] >= min);
4067  }
4068 #endif
4069  }
4070  for( w = ngubconsGOC1+1; w <= ngubconsGC1; w++ )
4071  finished[w] = SCIP_LONGINT_MAX;
4072 
4073  /* initialize unfinished table; note that variables in GNC1 GUBs
4074  * are sorted s.t. C1 variables come first and are sorted by non-decreasing weight.
4075  * GUBs in the group GCI are sorted by non-decreasing min{ a_k : k in GC1_j } where min{ a_k : k in GC1_j } always
4076  * comes from the first variable in the GUB
4077  */
4078  assert(ngubconsGNC1 <= ngubconsGC1);
4079  unfinished[0] = 0;
4080  for( w = 1; w <= ngubconsGNC1; w++ )
4081  {
4082  liftgubconsidx = gubconsGNC1[w-1];
4083 
4084  assert(gubset->gubconsstatus[liftgubconsidx] == GUBCONSSTATUS_BELONGSTOSET_GNC1);
4085  assert(gubset->gubconss[liftgubconsidx]->gubvarsstatus[0] == GUBVARSTATUS_BELONGSTOSET_C1);
4086 
4087  varidx = gubset->gubconss[liftgubconsidx]->gubvars[0];
4088 
4089  assert(varidx >= 0 && varidx < nvars);
4090  assert(liftcoefs[varidx] == 1);
4091 
4092  min = weights[varidx];
4093  unfinished[w] = unfinished[w-1] + min;
4094 
4095 #ifndef NDEBUG
4096  for( k = 1; k < gubset->gubconss[liftgubconsidx]->ngubvars
4097  && gubset->gubconss[liftgubconsidx]->gubvarsstatus[k] == GUBVARSTATUS_BELONGSTOSET_C1; k++ )
4098  {
4099  varidx = gubset->gubconss[liftgubconsidx]->gubvars[k];
4100  assert(varidx >= 0 && varidx < nvars);
4101  assert(liftcoefs[varidx] == 1);
4102  assert(weights[varidx] >= min );
4103  }
4104 #endif
4105  }
4106  for( w = ngubconsGNC1 + 1; w <= ngubconsGC1; w++ )
4107  unfinished[w] = SCIP_LONGINT_MAX;
4108 
4109  /* initialize minweights table; note that variables in GC1 GUBs
4110  * are sorted s.t. C1 variables come first and are sorted by non-decreasing weight.
4111  * we can directly initialize minweights instead of computing it from finished and unfinished (which would be more time
4112  * consuming) because is it has to be build using weights from C1 only.
4113  */
4114  assert(ngubconsGOC1 + ngubconsGNC1 == ngubconsGC1);
4115  minweights[0] = 0;
4116  for( w = 1; w <= ngubconsGC1; w++ )
4117  {
4118  liftgubconsidx = gubconsGC1[w-1];
4119 
4120  assert(gubset->gubconsstatus[liftgubconsidx] == GUBCONSSTATUS_BELONGSTOSET_GOC1
4121  || gubset->gubconsstatus[liftgubconsidx] == GUBCONSSTATUS_BELONGSTOSET_GNC1);
4122  assert(gubset->gubconss[liftgubconsidx]->gubvarsstatus[0] == GUBVARSTATUS_BELONGSTOSET_C1);
4123 
4124  varidx = gubset->gubconss[liftgubconsidx]->gubvars[0];
4125 
4126  assert(varidx >= 0 && varidx < nvars);
4127  assert(liftcoefs[varidx] == 1);
4128 
4129  min = weights[varidx];
4130  minweights[w] = minweights[w-1] + min;
4131 
4132 #ifndef NDEBUG
4133  for( k = 1; k < gubset->gubconss[liftgubconsidx]->ngubvars
4134  && gubset->gubconss[liftgubconsidx]->gubvarsstatus[k] == GUBVARSTATUS_BELONGSTOSET_C1; k++ )
4135  {
4136  varidx = gubset->gubconss[liftgubconsidx]->gubvars[k];
4137  assert(varidx >= 0 && varidx < nvars);
4138  assert(liftcoefs[varidx] == 1);
4139  assert(weights[varidx] >= min);
4140  }
4141 #endif
4142  }
4143  minweightslen = ngubconsGC1 + 1;
4144 
4145  /* gets sum of weights of variables fixed to one, i.e. sum of weights of C2 variables GC2 GUBs */
4146  fixedonesweight = 0;
4147  for( j = 0; j < ngubconsGC2; j++ )
4148  {
4149  varidx = gubset->gubconss[gubconsGC2[j]]->gubvars[0];
4150 
4151  assert(gubset->gubconss[gubconsGC2[j]]->ngubvars == 1);
4152  assert(varidx >= 0 && varidx < nvars);
4153  assert(gubset->gubconss[gubconsGC2[j]]->gubvarsstatus[0] == GUBVARSTATUS_BELONGSTOSET_C2);
4154 
4155  fixedonesweight += weights[varidx];
4156  }
4157  assert(fixedonesweight >= 0);
4158 
4159  /* initializes right hand side of lifted valid inequality */
4160  *liftrhs = alpha0;
4161 
4162  /* sequentially up-lifts all variables in GFC1 GUBs */
4163  for( j = 0; j < ngubconsGFC1; j++ )
4164  {
4165  liftgubconsidx = gubconsGFC1[j];
4166  assert(liftgubconsidx >= 0 && liftgubconsidx < ngubconss);
4167 
4168  /* GNC1 GUB: update unfinished table (remove current GUB, i.e., remove min weight of C1 vars in GUB) and
4169  * compute minweight table via updated unfinished table and aleady upto date finished table;
4170  */
4171  k = 0;
4172  if( gubset->gubconsstatus[liftgubconsidx] == GUBCONSSTATUS_BELONGSTOSET_GNC1 )
4173  {
4174  assert(gubset->gubconsstatus[liftgubconsidx] == GUBCONSSTATUS_BELONGSTOSET_GNC1);
4175  assert(gubset->gubconss[liftgubconsidx]->gubvarsstatus[0] == GUBVARSTATUS_BELONGSTOSET_C1);
4176  assert(ngubconsGNC1 > 0);
4177 
4178  /* get number of C1 variables of current GNC1 GUB and put them into array of variables in GUB that
4179  * are considered for the lifting, i.e., not capacity exceeding
4180  */
4181  for( ; k < gubset->gubconss[liftgubconsidx]->ngubvars
4182  && gubset->gubconss[liftgubconsidx]->gubvarsstatus[k] == GUBVARSTATUS_BELONGSTOSET_C1; k++ )
4183  liftgubvars[k] = gubset->gubconss[liftgubconsidx]->gubvars[k];
4184  assert(k >= 1);
4185 
4186  /* update unfinished table by removing current GNC1 GUB, i.e, remove C1 variable with minimal weight
4187  * unfinished[w] = MAX{unfinished[w], unfinished[w+1] - weight}, "weight" is the minimal weight of current GUB
4188  */
4189  weight = weights[liftgubvars[0]];
4190 
4191  weightdiff2 = unfinished[ngubconsGNC1] - weight;
4192  unfinished[ngubconsGNC1] = SCIP_LONGINT_MAX;
4193  for( w = ngubconsGNC1-1; w >= 1; w-- )
4194  {
4195  weightdiff1 = weightdiff2;
4196  weightdiff2 = unfinished[w] - weight;
4197 
4198  if( unfinished[w] < weightdiff1 )
4199  unfinished[w] = weightdiff1;
4200  else
4201  break;
4202  }
4203  ngubconsGNC1--;
4204 
4205  /* computes minweights table by combining unfished and fished tables */
4206  computeMinweightsGUB(minweights, finished, unfinished, minweightslen);
4207  assert(minweights[0] == 0);
4208  }
4209  /* GF GUB: no update of unfinished table (and minweight table) required because GF GUBs have no C1 variables and
4210  * are therefore not in the unfinished table
4211  */
4212  else
4213  assert(gubset->gubconsstatus[liftgubconsidx] == GUBCONSSTATUS_BELONGSTOSET_GF);
4214 
4215 #ifndef NDEBUG
4216  nliftgubC1 = k;
4217 #endif
4218  nliftgubvars = k;
4219  sumliftcoef = 0;
4220 
4221  /* compute lifting coefficient of F and R variables in GNC1 and GF GUBs (C1 vars have already liftcoef 1) */
4222  for( ; k < gubset->gubconss[liftgubconsidx]->ngubvars; k++ )
4223  {
4224  if( gubset->gubconss[liftgubconsidx]->gubvarsstatus[k] == GUBVARSTATUS_BELONGSTOSET_F
4225  || gubset->gubconss[liftgubconsidx]->gubvarsstatus[k] == GUBVARSTATUS_BELONGSTOSET_R )
4226  {
4227  liftvar = gubset->gubconss[liftgubconsidx]->gubvars[k];
4228  weight = weights[liftvar];
4229  assert(weight > 0);
4230  assert(liftvar >= 0 && liftvar < nvars);
4231  assert(capacity - weight >= 0);
4232 
4233  /* put variable into array of variables in GUB that are considered for the lifting,
4234  * i.e., not capacity exceeding
4235  */
4236  liftgubvars[nliftgubvars] = liftvar;
4237  nliftgubvars++;
4238 
4239  /* knapsack problem is infeasible:
4240  * sets z = 0
4241  */
4242  if( capacity - fixedonesweight - weight < 0 )
4243  {
4244  z = 0;
4245  }
4246  /* knapsack problem is feasible:
4247  * sets z = max { w : 0 <= w <= liftrhs, minweights_i[w] <= a_0 - fixedonesweight - a_{j_i} } = liftrhs,
4248  * if minweights_i[liftrhs] <= a_0 - fixedonesweight - a_{j_i}
4249  */
4250  else if( minweights[*liftrhs] <= capacity - fixedonesweight - weight )
4251  {
4252  z = *liftrhs;
4253  }
4254  /* knapsack problem is feasible:
4255  * binary search to find z = max {w : 0 <= w <= liftrhs, minweights_i[w] <= a_0 - fixedonesweight - a_{j_i}}
4256  */
4257  else
4258  {
4259  assert((*liftrhs) + 1 >= minweightslen || minweights[(*liftrhs) + 1] > capacity - fixedonesweight - weight);
4260  left = 0;
4261  right = (*liftrhs) + 1;
4262  while( left < right - 1 )
4263  {
4264  middle = (left + right) / 2;
4265  assert(0 <= middle && middle < minweightslen);
4266  if( minweights[middle] <= capacity - fixedonesweight - weight )
4267  left = middle;
4268  else
4269  right = middle;
4270  }
4271  assert(left == right - 1);
4272  assert(0 <= left && left < minweightslen);
4273  assert(minweights[left] <= capacity - fixedonesweight - weight);
4274  assert(left == minweightslen - 1 || minweights[left+1] > capacity - fixedonesweight - weight);
4275 
4276  /* now z = left */
4277  z = left;
4278  assert(z <= *liftrhs);
4279  }
4280 
4281  /* calculates lifting coefficients alpha_{j_i} = liftrhs - z */
4282  liftcoef = (*liftrhs) - z;
4283  liftcoefs[liftvar] = liftcoef;
4284  assert(liftcoef >= 0 && liftcoef <= (*liftrhs) + 1);
4285 
4286  /* updates activity of current valid inequality */
4287  (*cutact) += liftcoef * solvals[liftvar];
4288 
4289  /* updates sum of all lifting coefficients in GUB */
4290  sumliftcoef += liftcoefs[liftvar];
4291  }
4292  else
4293  assert(gubset->gubconss[liftgubconsidx]->gubvarsstatus[k] == GUBVARSTATUS_CAPACITYEXCEEDED);
4294  }
4295  /* at least one variable is in F or R (j = number of C1 variables in current GUB) */
4296  assert(nliftgubvars > nliftgubC1);
4297 
4298  /* activity of current valid inequality will not change if (sum of alpha_{j_i} in GUB) = 0
4299  * and finished and minweight table can be updated easily as only C1 variables need to be considered;
4300  * not needed for GF GUBs
4301  */
4302  if( sumliftcoef == 0 )
4303  {
4304  if( gubset->gubconsstatus[liftgubconsidx] == GUBCONSSTATUS_BELONGSTOSET_GNC1 )
4305  {
4306  weight = weights[liftgubvars[0]];
4307  /* update finished table and minweights table by applying special case of
4308  * finished[w] = MIN{finished[w], finished[w-1] + weight}, "weight" is the minimal weight of current GUB
4309  * minweights[w] = MIN{minweights[w], minweights[w-1] + weight}, "weight" is the minimal weight of current GUB
4310  */
4311  for( w = minweightslen-1; w >= 1; w-- )
4312  {
4313  SCIP_Longint tmpval;
4314 
4315  tmpval = safeAddMinweightsGUB(finished[w-1], weight);
4316  finished[w] = MIN(finished[w], tmpval);
4317 
4318  tmpval = safeAddMinweightsGUB(minweights[w-1], weight);
4319  minweights[w] = MIN(minweights[w], tmpval);
4320  }
4321  }
4322  else
4323  assert(gubset->gubconsstatus[liftgubconsidx] == GUBCONSSTATUS_BELONGSTOSET_GF);
4324 
4325  continue;
4326  }
4327 
4328  /* enlarges current minweights tables(finished, unfinished, minweights):
4329  * from minweightlen = |gubconsGC1| + sum_{k=1,2,...,i-1}sum_{j in Q_k} alpha_j + 1 entries
4330  * to |gubconsGC1| + sum_{k=1,2,...,i }sum_{j in Q_k} alpha_j + 1 entries
4331  * and sets minweights_i[w] = infinity for
4332  * w = |gubconsGC1| + sum_{k=1,2,..,i-1}sum_{j in Q_k} alpha_j+1,..,|C1| + sum_{k=1,2,..,i}sum_{j in Q_k} alpha_j
4333  */
4334  tmplen = minweightslen; /* will be updated in enlargeMinweights() */
4335  tmpsize = minweightssize;
4336  SCIP_CALL( enlargeMinweights(scip, &unfinished, &tmplen, &tmpsize, tmplen + sumliftcoef) );
4337  tmplen = minweightslen;
4338  tmpsize = minweightssize;
4339  SCIP_CALL( enlargeMinweights(scip, &finished, &tmplen, &tmpsize, tmplen + sumliftcoef) );
4340  SCIP_CALL( enlargeMinweights(scip, &minweights, &minweightslen, &minweightssize, minweightslen + sumliftcoef) );
4341 
4342  /* update finished table and minweight table;
4343  * note that instead of computing minweight table from updated finished and updated unfinished table again
4344  * (for the lifting coefficient, we had to update unfinished table and compute minweight table), we here
4345  * only need to update the minweight table and the updated finished in the same way (i.e., computing for minweight
4346  * not needed because only finished table changed at this point and the change was "adding" one weight)
4347  *
4348  * update formular for minweight table is: minweight_i+1[w] =
4349  * min{ minweights_i[w], min{ minweights_i[w - alpha_k]^{+} + a_k : k in GUB_j_i } }
4350  * formular for finished table has the same pattern.
4351  */
4352  for( w = minweightslen-1; w >= 0; w-- )
4353  {
4354  SCIP_Longint minminweight;
4355  SCIP_Longint minfinished;
4356 
4357  for( k = 0; k < nliftgubvars; k++ )
4358  {
4359  liftcoef = liftcoefs[liftgubvars[k]];
4360  weight = weights[liftgubvars[k]];
4361 
4362  if( w < liftcoef )
4363  {
4364  minfinished = MIN(finished[w], weight);
4365  minminweight = MIN(minweights[w], weight);
4366 
4367  finished[w] = minfinished;
4368  minweights[w] = minminweight;
4369  }
4370  else
4371  {
4372  SCIP_Longint tmpval;
4373 
4374  assert(w >= liftcoef);
4375 
4376  tmpval = safeAddMinweightsGUB(finished[w-liftcoef], weight);
4377  minfinished = MIN(finished[w], tmpval);
4378 
4379  tmpval = safeAddMinweightsGUB(minweights[w-liftcoef], weight);
4380  minminweight = MIN(minweights[w], tmpval);
4381 
4382  finished[w] = minfinished;
4383  minweights[w] = minminweight;
4384  }
4385  }
4386  }
4387  assert(minweights[0] == 0);
4388  }
4389  assert(ngubconsGNC1 == 0);
4390 
4391  /* note: now the unfinished table no longer exists, i.e., it is "0, MAX, MAX, ..." and minweight equals to finished;
4392  * therefore, only work with minweight table from here on
4393  */
4394 
4395  /* sequentially down-lifts C2 variables contained in trivial GC2 GUBs */
4396  for( j = 0; j < ngubconsGC2; j++ )
4397  {
4398  liftgubconsidx = gubconsGC2[j];
4399 
4400  assert(liftgubconsidx >=0 && liftgubconsidx < ngubconss);
4401  assert(gubset->gubconsstatus[liftgubconsidx] == GUBCONSSTATUS_BELONGSTOSET_GC2);
4402  assert(gubset->gubconss[liftgubconsidx]->ngubvars == 1);
4403  assert(gubset->gubconss[liftgubconsidx]->gubvarsstatus[0] == GUBVARSTATUS_BELONGSTOSET_C2);
4404 
4405  liftvar = gubset->gubconss[liftgubconsidx]->gubvars[0]; /* C2 GUBs contain only one variable */
4406  weight = weights[liftvar];
4407 
4408  assert(liftvar >= 0 && liftvar < nvars);
4409  assert(SCIPisFeasEQ(scip, solvals[liftvar], 1.0));
4410  assert(weight > 0);
4411 
4412  /* uses binary search to find
4413  * z = max { w : 0 <= w <= |C_1| + sum_{k=1}^{i-1} alpha_{j_k}, minweights_[w] <= a_0 - fixedonesweight + a_{j_i}}
4414  */
4415  left = 0;
4416  right = minweightslen;
4417  while( left < right - 1 )
4418  {
4419  middle = (left + right) / 2;
4420  assert(0 <= middle && middle < minweightslen);
4421  if( minweights[middle] <= capacity - fixedonesweight + weight )
4422  left = middle;
4423  else
4424  right = middle;
4425  }
4426  assert(left == right - 1);
4427  assert(0 <= left && left < minweightslen);
4428  assert(minweights[left] <= capacity - fixedonesweight + weight);
4429  assert(left == minweightslen - 1 || minweights[left + 1] > capacity - fixedonesweight + weight);
4430 
4431  /* now z = left */
4432  z = left;
4433  assert(z >= *liftrhs);
4434 
4435  /* calculates lifting coefficients alpha_{j_i} = z - liftrhs */
4436  liftcoef = z - (*liftrhs);
4437  liftcoefs[liftvar] = liftcoef;
4438  assert(liftcoef >= 0);
4439 
4440  /* updates sum of weights of variables fixed to one */
4441  fixedonesweight -= weight;
4442 
4443  /* updates right-hand side of current valid inequality */
4444  (*liftrhs) += liftcoef;
4445  assert(*liftrhs >= alpha0);
4446 
4447  /* minweight table and activity of current valid inequality will not change, if alpha_{j_i} = 0 */
4448  if( liftcoef == 0 )
4449  continue;
4450 
4451  /* updates activity of current valid inequality */
4452  (*cutact) += liftcoef * solvals[liftvar];
4453 
4454  /* enlarges current minweight table:
4455  * from minweightlen = |gubconsGC1| + sum_{k=1,2,...,i-1}sum_{j in Q_k} alpha_j + 1 entries
4456  * to |gubconsGC1| + sum_{k=1,2,...,i }sum_{j in Q_k} alpha_j + 1 entries
4457  * and sets minweights_i[w] = infinity for
4458  * w = |C1| + sum_{k=1,2,...,i-1}sum_{j in Q_k} alpha_j + 1 , ... , |C1| + sum_{k=1,2,...,i}sum_{j in Q_k} alpha_j
4459  */
4460  SCIP_CALL( enlargeMinweights(scip, &minweights, &minweightslen, &minweightssize, minweightslen + liftcoef) );
4461 
4462  /* updates minweight table: minweight_i+1[w] =
4463  * min{ minweights_i[w], a_{j_i}}, if w < alpha_j_i
4464  * min{ minweights_i[w], minweights_i[w - alpha_j_i] + a_j_i}, if w >= alpha_j_i
4465  */
4466  for( w = minweightslen - 1; w >= 0; w-- )
4467  {
4468  if( w < liftcoef )
4469  {
4470  min = MIN(minweights[w], weight);
4471  minweights[w] = min;
4472  }
4473  else
4474  {
4475  SCIP_Longint tmpval;
4476 
4477  assert(w >= liftcoef);
4478 
4479  tmpval = safeAddMinweightsGUB(minweights[w-liftcoef], weight);
4480  min = MIN(minweights[w], tmpval);
4481  minweights[w] = min;
4482  }
4483  }
4484  }
4485  assert(fixedonesweight == 0);
4486  assert(*liftrhs >= alpha0);
4487 
4488  /* sequentially up-lifts variables in GUB constraints in GR GUBs */
4489  for( j = 0; j < ngubconsGR; j++ )
4490  {
4491  liftgubconsidx = gubconsGR[j];
4492 
4493  assert(liftgubconsidx >=0 && liftgubconsidx < ngubconss);
4494  assert(gubset->gubconsstatus[liftgubconsidx] == GUBCONSSTATUS_BELONGSTOSET_GR);
4495 
4496  sumliftcoef = 0;
4497  nliftgubvars = 0;
4498  for( k = 0; k < gubset->gubconss[liftgubconsidx]->ngubvars; k++ )
4499  {
4500  if(gubset->gubconss[liftgubconsidx]->gubvarsstatus[k] == GUBVARSTATUS_BELONGSTOSET_R )
4501  {
4502  liftvar = gubset->gubconss[liftgubconsidx]->gubvars[k];
4503  weight = weights[liftvar];
4504  assert(weight > 0);
4505  assert(liftvar >= 0 && liftvar < nvars);
4506  assert(capacity - weight >= 0);
4507  assert((*liftrhs) + 1 >= minweightslen || minweights[(*liftrhs) + 1] > capacity - weight);
4508 
4509  /* put variable into array of variables in GUB that are considered for the lifting,
4510  * i.e., not capacity exceeding
4511  */
4512  liftgubvars[nliftgubvars] = liftvar;
4513  nliftgubvars++;
4514 
4515  /* sets z = max { w : 0 <= w <= liftrhs, minweights_i[w] <= a_0 - a_{j_i} } = liftrhs,
4516  * if minweights_i[liftrhs] <= a_0 - a_{j_i}
4517  */
4518  if( minweights[*liftrhs] <= capacity - weight )
4519  {
4520  z = *liftrhs;
4521  }
4522  /* uses binary search to find z = max { w : 0 <= w <= liftrhs, minweights_i[w] <= a_0 - a_{j_i} }
4523  */
4524  else
4525  {
4526  left = 0;
4527  right = (*liftrhs) + 1;
4528  while( left < right - 1 )
4529  {
4530  middle = (left + right) / 2;
4531  assert(0 <= middle && middle < minweightslen);
4532  if( minweights[middle] <= capacity - weight )
4533  left = middle;
4534  else
4535  right = middle;
4536  }
4537  assert(left == right - 1);
4538  assert(0 <= left && left < minweightslen);
4539  assert(minweights[left] <= capacity - weight);
4540  assert(left == minweightslen - 1 || minweights[left + 1] > capacity - weight);
4541 
4542  /* now z = left */
4543  z = left;
4544  assert(z <= *liftrhs);
4545  }
4546  /* calculates lifting coefficients alpha_{j_i} = liftrhs - z */
4547  liftcoef = (*liftrhs) - z;
4548  liftcoefs[liftvar] = liftcoef;
4549  assert(liftcoef >= 0 && liftcoef <= (*liftrhs) + 1);
4550 
4551  /* updates activity of current valid inequality */
4552  (*cutact) += liftcoef * solvals[liftvar];
4553 
4554  /* updates sum of all lifting coefficients in GUB */
4555  sumliftcoef += liftcoefs[liftvar];
4556  }
4557  else
4558  assert(gubset->gubconss[liftgubconsidx]->gubvarsstatus[k] == GUBVARSTATUS_CAPACITYEXCEEDED);
4559  }
4560  assert(nliftgubvars >= 1); /* at least one variable is in R */
4561 
4562  /* minweight table and activity of current valid inequality will not change if (sum of alpha_{j_i} in GUB) = 0 */
4563  if( sumliftcoef == 0 )
4564  continue;
4565 
4566  /* updates minweight table: minweight_i+1[w] =
4567  * min{ minweights_i[w], min{ minweights_i[w - alpha_k]^{+} + a_k : k in GUB_j_i } }
4568  */
4569  for( w = *liftrhs; w >= 0; w-- )
4570  {
4571  for( k = 0; k < nliftgubvars; k++ )
4572  {
4573  liftcoef = liftcoefs[liftgubvars[k]];
4574  weight = weights[liftgubvars[k]];
4575 
4576  if( w < liftcoef )
4577  {
4578  min = MIN(minweights[w], weight);
4579  minweights[w] = min;
4580  }
4581  else
4582  {
4583  SCIP_Longint tmpval;
4584 
4585  assert(w >= liftcoef);
4586 
4587  tmpval = safeAddMinweightsGUB(minweights[w-liftcoef], weight);
4588  min = MIN(minweights[w], tmpval);
4589  minweights[w] = min;
4590  }
4591  }
4592  }
4593  assert(minweights[0] == 0);
4594  }
4595 
4596  /* frees temporary memory */
4597  SCIPfreeBufferArray(scip, &minweights);
4598  SCIPfreeBufferArray(scip, &finished);
4599  SCIPfreeBufferArray(scip, &unfinished);
4600  SCIPfreeBufferArray(scip, &liftgubvars);
4601  SCIPfreeBufferArray(scip, &gubconsGOC1 );
4602  SCIPfreeBufferArray(scip, &gubconsGNC1);
4603 
4604  return SCIP_OKAY;
4605 }
4606 
4607 /** lifts given minimal cover inequality
4608  * \f[
4609  * \sum_{j \in C} x_j \leq |C| - 1
4610  * \f]
4611  * valid for
4612  * \f[
4613  * S^0 = \{ x \in {0,1}^{|C|} : \sum_{j \in C} a_j x_j \leq a_0 \}
4614  * \f]
4615  * to a valid inequality
4616  * \f[
4617  * \sum_{j \in C} x_j + \sum_{j \in N \setminus C} \alpha_j x_j \leq |C| - 1
4618  * \f]
4619  * for
4620  * \f[
4621  * S = \{ x \in {0,1}^{|N|} : \sum_{j \in N} a_j x_j \leq a_0 \};
4622  * \f]
4623  * uses superadditive up-lifting for the variables in \f$N \setminus C\f$.
4624  */
4625 static
4627  SCIP* scip, /**< SCIP data structure */
4628  SCIP_VAR** vars, /**< variables in knapsack constraint */
4629  int nvars, /**< number of variables in knapsack constraint */
4630  int ntightened, /**< number of variables with tightened upper bound */
4631  SCIP_Longint* weights, /**< weights of variables in knapsack constraint */
4632  SCIP_Longint capacity, /**< capacity of knapsack */
4633  SCIP_Real* solvals, /**< solution values of all problem variables */
4634  int* covervars, /**< cover variables */
4635  int* noncovervars, /**< noncover variables */
4636  int ncovervars, /**< number of cover variables */
4637  int nnoncovervars, /**< number of noncover variables */
4638  SCIP_Longint coverweight, /**< weight of cover */
4639  SCIP_Real* liftcoefs, /**< pointer to store lifting coefficient of vars in knapsack constraint */
4640  SCIP_Real* cutact /**< pointer to store activity of lifted valid inequality */
4641  )
4642 {
4643  SCIP_Longint* maxweightsums;
4644  SCIP_Longint* intervalends;
4645  SCIP_Longint* rhos;
4646  SCIP_Real* sortkeys;
4647  SCIP_Longint lambda;
4648  int j;
4649  int h;
4650 
4651  assert(scip != NULL);
4652  assert(vars != NULL);
4653  assert(nvars >= 0);
4654  assert(weights != NULL);
4655  assert(capacity >= 0);
4656  assert(solvals != NULL);
4657  assert(covervars != NULL);
4658  assert(noncovervars != NULL);
4659  assert(ncovervars > 0 && ncovervars <= nvars);
4660  assert(nnoncovervars >= 0 && nnoncovervars <= nvars - ntightened);
4661  assert(ncovervars + nnoncovervars == nvars - ntightened);
4662  assert(liftcoefs != NULL);
4663  assert(cutact != NULL);
4664 
4665  /* allocates temporary memory */
4666  SCIP_CALL( SCIPallocBufferArray(scip, &sortkeys, ncovervars) );
4667  SCIP_CALL( SCIPallocBufferArray(scip, &maxweightsums, ncovervars + 1) );
4668  SCIP_CALL( SCIPallocBufferArray(scip, &intervalends, ncovervars) );
4669  SCIP_CALL( SCIPallocBufferArray(scip, &rhos, ncovervars) );
4670 
4671  /* initializes data structures */
4672  BMSclearMemoryArray(liftcoefs, nvars);
4673  *cutact = 0.0;
4674 
4675  /* sets lifting coefficient of variables in C, sorts variables in C such that a_1 >= a_2 >= ... >= a_|C|
4676  * and calculates activity of current valid inequality
4677  */
4678  for( j = 0; j < ncovervars; j++ )
4679  {
4680  assert(liftcoefs[covervars[j]] == 0.0);
4681  liftcoefs[covervars[j]] = 1.0;
4682  sortkeys[j] = (SCIP_Real) weights[covervars[j]];
4683  (*cutact) += solvals[covervars[j]];
4684  }
4685  SCIPsortDownRealInt(sortkeys, covervars, ncovervars);
4686 
4687  /* calculates weight excess of cover C */
4688  lambda = coverweight - capacity;
4689  assert(lambda > 0);
4690 
4691  /* calculates A_h for h = 0,...,|C|, I_h for h = 1,...,|C| and rho_h for h = 1,...,|C| */
4692  maxweightsums[0] = 0;
4693  for( h = 1; h <= ncovervars; h++ )
4694  {
4695  maxweightsums[h] = maxweightsums[h-1] + weights[covervars[h-1]];
4696  intervalends[h-1] = maxweightsums[h] - lambda;
4697  rhos[h-1] = MAX(0, weights[covervars[h-1]] - weights[covervars[0]] + lambda);
4698  }
4699 
4700  /* sorts variables in N\C such that a_{j_1} <= a_{j_2} <= ... <= a_{j_t} */
4701  for( j = 0; j < nnoncovervars; j++ )
4702  sortkeys[j] = (SCIP_Real) (weights[noncovervars[j]]);
4703  SCIPsortRealInt(sortkeys, noncovervars, nnoncovervars);
4704 
4705  /* calculates lifting coefficient for all variables in N\C */
4706  h = 0;
4707  for( j = 0; j < nnoncovervars; j++ )
4708  {
4709  int liftvar;
4710  SCIP_Longint weight;
4711  SCIP_Real liftcoef;
4712 
4713  liftvar = noncovervars[j];
4714  weight = weights[liftvar];
4715 
4716  while( intervalends[h] < weight )
4717  h++;
4718 
4719  if( h == 0 )
4720  liftcoef = h;
4721  else
4722  {
4723  if( weight <= intervalends[h-1] + rhos[h] )
4724  {
4725  SCIP_Real tmp1;
4726  SCIP_Real tmp2;
4727  tmp1 = (SCIP_Real) (intervalends[h-1] + rhos[h] - weight);
4728  tmp2 = (SCIP_Real) rhos[1];
4729  liftcoef = h - ( tmp1 / tmp2 );
4730  }
4731  else
4732  liftcoef = h;
4733  }
4734 
4735  /* sets lifting coefficient */
4736  assert(liftcoefs[liftvar] == 0.0);
4737  liftcoefs[liftvar] = liftcoef;
4738 
4739  /* updates activity of current valid inequality */
4740  (*cutact) += liftcoef * solvals[liftvar];
4741  }
4742 
4743  /* frees temporary memory */
4744  SCIPfreeBufferArray(scip, &rhos);
4745  SCIPfreeBufferArray(scip, &intervalends);
4746  SCIPfreeBufferArray(scip, &maxweightsums);
4747  SCIPfreeBufferArray(scip, &sortkeys);
4748 
4749  return SCIP_OKAY;
4750 }
4751 
4752 
4753 /** separates lifted minimal cover inequalities using sequential up- and down-lifting and GUB information, if wanted, for
4754  * given knapsack problem
4755 */
4756 static
4758  SCIP* scip, /**< SCIP data structure */
4759  SCIP_CONS* cons, /**< originating constraint of the knapsack problem, or NULL */
4760  SCIP_SEPA* sepa, /**< originating separator of the knapsack problem, or NULL */
4761  SCIP_VAR** vars, /**< variables in knapsack constraint */
4762  int nvars, /**< number of variables in knapsack constraint */
4763  int ntightened, /**< number of variables with tightened upper bound */
4764  SCIP_Longint* weights, /**< weights of variables in knapsack constraint */
4765  SCIP_Longint capacity, /**< capacity of knapsack */
4766  SCIP_Real* solvals, /**< solution values of all problem variables */
4767  int* mincovervars, /**< mincover variables */
4768  int* nonmincovervars, /**< nonmincover variables */
4769  int nmincovervars, /**< number of mincover variables */
4770  int nnonmincovervars, /**< number of nonmincover variables */
4771  SCIP_SOL* sol, /**< primal SCIP solution to separate, NULL for current LP solution */
4772  SCIP_GUBSET* gubset, /**< GUB set data structure, NULL if no GUB information should be used */
4773  SCIP_Bool* cutoff, /**< pointer to store whether a cutoff has been detected */
4774  int* ncuts /**< pointer to add up the number of found cuts */
4775  )
4776 {
4777  int* varsC1;
4778  int* varsC2;
4779  int* varsF;
4780  int* varsR;
4781  int nvarsC1;
4782  int nvarsC2;
4783  int nvarsF;
4784  int nvarsR;
4785  SCIP_Real cutact;
4786  int* liftcoefs;
4787  int liftrhs;
4788 
4789  assert( cutoff != NULL );
4790  *cutoff = FALSE;
4791 
4792  /* allocates temporary memory */
4793  SCIP_CALL( SCIPallocBufferArray(scip, &varsC1, nvars) );
4794  SCIP_CALL( SCIPallocBufferArray(scip, &varsC2, nvars) );
4795  SCIP_CALL( SCIPallocBufferArray(scip, &varsF, nvars) );
4796  SCIP_CALL( SCIPallocBufferArray(scip, &varsR, nvars) );
4797  SCIP_CALL( SCIPallocBufferArray(scip, &liftcoefs, nvars) );
4798 
4799  /* gets partition (C_1,C_2) of C, i.e. C_1 & C_2 = C and C_1 cap C_2 = emptyset, with C_1 not empty; chooses partition
4800  * as follows
4801  * C_2 = { j in C : x*_j = 1 } and
4802  * C_1 = C\C_2
4803  */
4804  getPartitionCovervars(scip, solvals, mincovervars, nmincovervars, varsC1, varsC2, &nvarsC1, &nvarsC2);
4805  assert(nvarsC1 + nvarsC2 == nmincovervars);
4806  assert(nmincovervars > 0);
4807  assert(nvarsC1 >= 0); /* nvarsC1 > 0 does not always hold, because relaxed knapsack conss may already be violated */
4808 
4809  /* changes partition (C_1,C_2) of minimal cover C, if |C1| = 1, by moving one variable from C2 to C1 */
4810  if( nvarsC1 < 2 && nvarsC2 > 0)
4811  {
4812  SCIP_CALL( changePartitionCovervars(scip, weights, varsC1, varsC2, &nvarsC1, &nvarsC2) );
4813  assert(nvarsC1 >= 1);
4814  }
4815  assert(nvarsC2 == 0 || nvarsC1 >= 1);
4816 
4817  /* gets partition (F,R) of N\C, i.e. F & R = N\C and F cap R = emptyset; chooses partition as follows
4818  * R = { j in N\C : x*_j = 0 } and
4819  * F = (N\C)\F
4820  */
4821  getPartitionNoncovervars(scip, solvals, nonmincovervars, nnonmincovervars, varsF, varsR, &nvarsF, &nvarsR);
4822  assert(nvarsF + nvarsR == nnonmincovervars);
4823  assert(nvarsC1 + nvarsC2 + nvarsF + nvarsR == nvars - ntightened);
4824 
4825  /* lift cuts without GUB information */
4826  if( gubset == NULL )
4827  {
4828  /* sorts variables in F, C_2, R according to the second level lifting sequence that will be used in the sequential
4829  * lifting procedure
4830  */
4831  SCIP_CALL( getLiftingSequence(scip, solvals, weights, varsF, varsC2, varsR, nvarsF, nvarsC2, nvarsR) );
4832 
4833  /* lifts minimal cover inequality sum_{j in C_1} x_j <= |C_1| - 1 valid for
4834  *
4835  * S^0 = { x in {0,1}^|C_1| : sum_{j in C_1} a_j x_j <= a_0 - sum_{j in C_2} a_j }
4836  *
4837  * to a valid inequality sum_{j in C_1} x_j + sum_{j in N\C_1} alpha_j x_j <= |C_1| - 1 + sum_{j in C_2} alpha_j for
4838  *
4839  * S = { x in {0,1}^|N| : sum_{j in N} a_j x_j <= a_0 },
4840  *
4841  * uses sequential up-lifting for the variables in F, sequential down-lifting for the variable in C_2 and sequential
4842  * up-lifting for the variables in R according to the second level lifting sequence
4843  */
4844  SCIP_CALL( sequentialUpAndDownLifting(scip, vars, nvars, ntightened, weights, capacity, solvals, varsC1, varsC2,
4845  varsF, varsR, nvarsC1, nvarsC2, nvarsF, nvarsR, nvarsC1 - 1, liftcoefs, &cutact, &liftrhs) );
4846  }
4847  /* lift cuts with GUB information */
4848  else
4849  {
4850  int* gubconsGC1;
4851  int* gubconsGC2;
4852  int* gubconsGFC1;
4853  int* gubconsGR;
4854  int ngubconsGC1;
4855  int ngubconsGC2;
4856  int ngubconsGFC1;
4857  int ngubconsGR;
4858  int ngubconss;
4859  int nconstightened;
4860  int maxgubvarssize;
4861 
4862  assert(nvars == gubset->nvars);
4863 
4864  ngubconsGC1 = 0;
4865  ngubconsGC2 = 0;
4866  ngubconsGFC1 = 0;
4867  ngubconsGR = 0;
4868  ngubconss = gubset->ngubconss;
4869  nconstightened = 0;
4870  maxgubvarssize = 0;
4871 
4872  /* allocates temporary memory */
4873  SCIP_CALL( SCIPallocBufferArray(scip, &gubconsGC1, ngubconss) );
4874  SCIP_CALL( SCIPallocBufferArray(scip, &gubconsGC2, ngubconss) );
4875  SCIP_CALL( SCIPallocBufferArray(scip, &gubconsGFC1, ngubconss) );
4876  SCIP_CALL( SCIPallocBufferArray(scip, &gubconsGR, ngubconss) );
4877 
4878  /* categorizies GUBs of knapsack GUB partion into GOC1, GNC1, GF, GC2, and GR and computes a lifting sequence of
4879  * the GUBs for the sequential GUB wise lifting procedure
4880  */
4881  SCIP_CALL( getLiftingSequenceGUB(scip, gubset, solvals, weights, varsC1, varsC2, varsF, varsR, nvarsC1,
4882  nvarsC2, nvarsF, nvarsR, gubconsGC1, gubconsGC2, gubconsGFC1, gubconsGR, &ngubconsGC1, &ngubconsGC2,
4883  &ngubconsGFC1, &ngubconsGR, &nconstightened, &maxgubvarssize) );
4884 
4885  /* lifts minimal cover inequality sum_{j in C_1} x_j <= |C_1| - 1 valid for
4886  *
4887  * S^0 = { x in {0,1}^|C_1| : sum_{j in C_1} a_j x_j <= a_0 - sum_{j in C_2} a_j,
4888  * sum_{j in Q_i} x_j <= 1, forall i in I }
4889  *
4890  * to a valid inequality sum_{j in C_1} x_j + sum_{j in N\C_1} alpha_j x_j <= |C_1| - 1 + sum_{j in C_2} alpha_j for
4891  *
4892  * S = { x in {0,1}^|N| : sum_{j in N} a_j x_j <= a_0, sum_{j in Q_i} x_j <= 1, forall i in I },
4893  *
4894  * uses sequential up-lifting for the variables in GUB constraints in gubconsGFC1,
4895  * sequential down-lifting for the variables in GUB constraints in gubconsGC2, and
4896  * sequential up-lifting for the variabels in GUB constraints in gubconsGR.
4897  */
4898  SCIP_CALL( sequentialUpAndDownLiftingGUB(scip, gubset, vars, nconstightened, weights, capacity, solvals, gubconsGC1,
4899  gubconsGC2, gubconsGFC1, gubconsGR, ngubconsGC1, ngubconsGC2, ngubconsGFC1, ngubconsGR,
4900  MIN(nvarsC1 - 1, ngubconsGC1), liftcoefs, &cutact, &liftrhs, maxgubvarssize) );
4901 
4902  /* frees temporary memory */
4903  SCIPfreeBufferArray(scip, &gubconsGR);
4904  SCIPfreeBufferArray(scip, &gubconsGFC1);
4905  SCIPfreeBufferArray(scip, &gubconsGC2);
4906  SCIPfreeBufferArray(scip, &gubconsGC1);
4907  }
4908 
4909  /* checks, if lifting yielded a violated cut */
4910  if( SCIPisEfficacious(scip, (cutact - liftrhs)/sqrt((SCIP_Real)MAX(liftrhs, 1))) )
4911  {
4912  SCIP_ROW* row;
4913  char name[SCIP_MAXSTRLEN];
4914  int j;
4915 
4916  /* creates LP row */
4917  assert( cons == NULL || sepa == NULL );
4918  if ( cons != NULL )
4919  {
4920  (void) SCIPsnprintf(name, SCIP_MAXSTRLEN, "%s_mcseq%" SCIP_LONGINT_FORMAT "", SCIPconsGetName(cons), SCIPconshdlrGetNCutsFound(SCIPconsGetHdlr(cons)));
4921  SCIP_CALL( SCIPcreateEmptyRowCons(scip, &row, SCIPconsGetHdlr(cons), name, -SCIPinfinity(scip), (SCIP_Real)liftrhs,
4922  cons != NULL ? SCIPconsIsLocal(cons) : FALSE, FALSE,
4923  cons != NULL ? SCIPconsIsRemovable(cons) : TRUE) );
4924  }
4925  else if ( sepa != NULL )
4926  {
4927  (void) SCIPsnprintf(name, SCIP_MAXSTRLEN, "%s_mcseq_%" SCIP_LONGINT_FORMAT "", SCIPsepaGetName(sepa), SCIPsepaGetNCutsFound(sepa));
4928  SCIP_CALL( SCIPcreateEmptyRowSepa(scip, &row, sepa, name, -SCIPinfinity(scip), (SCIP_Real)liftrhs, FALSE, FALSE, TRUE) );
4929  }
4930  else
4931  {
4932  (void) SCIPsnprintf(name, SCIP_MAXSTRLEN, "nn_mcseq_%" SCIP_LONGINT_FORMAT "", *ncuts);
4933  SCIP_CALL( SCIPcreateEmptyRowUnspec(scip, &row, name, -SCIPinfinity(scip), (SCIP_Real)liftrhs, FALSE, FALSE, TRUE) );
4934  }
4935 
4936  /* adds all variables in the knapsack constraint with calculated lifting coefficient to the cut */
4937  SCIP_CALL( SCIPcacheRowExtensions(scip, row) );
4938  assert(nvarsC1 + nvarsC2 + nvarsF + nvarsR == nvars - ntightened);
4939  for( j = 0; j < nvarsC1; j++ )
4940  {
4941  SCIP_CALL( SCIPaddVarToRow(scip, row, vars[varsC1[j]], 1.0) );
4942  }
4943  for( j = 0; j < nvarsC2; j++ )
4944  {
4945  if( liftcoefs[varsC2[j]] > 0 )
4946  {
4947  SCIP_CALL( SCIPaddVarToRow(scip, row, vars[varsC2[j]], (SCIP_Real)liftcoefs[varsC2[j]]) );
4948  }
4949  }
4950  for( j = 0; j < nvarsF; j++ )
4951  {
4952  if( liftcoefs[varsF[j]] > 0 )
4953  {
4954  SCIP_CALL( SCIPaddVarToRow(scip, row, vars[varsF[j]], (SCIP_Real)liftcoefs[varsF[j]]) );
4955  }
4956  }
4957  for( j = 0; j < nvarsR; j++ )
4958  {
4959  if( liftcoefs[varsR[j]] > 0 )
4960  {
4961  SCIP_CALL( SCIPaddVarToRow(scip, row, vars[varsR[j]], (SCIP_Real)liftcoefs[varsR[j]]) );
4962  }
4963  }
4964  SCIP_CALL( SCIPflushRowExtensions(scip, row) );
4965 
4966  /* checks, if cut is violated enough */
4967  if( SCIPisCutEfficacious(scip, sol, row) )
4968  {
4969  if( cons != NULL )
4970  {
4971  SCIP_CALL( SCIPresetConsAge(scip, cons) );
4972  }
4973  SCIP_CALL( SCIPaddCut(scip, sol, row, FALSE, cutoff) );
4974  (*ncuts)++;
4975  }
4976  SCIP_CALL( SCIPreleaseRow(scip, &row) );
4977  }
4978 
4979  /* frees temporary memory */
4980  SCIPfreeBufferArray(scip, &liftcoefs);
4981  SCIPfreeBufferArray(scip, &varsR);
4982  SCIPfreeBufferArray(scip, &varsF);
4983  SCIPfreeBufferArray(scip, &varsC2);
4984  SCIPfreeBufferArray(scip, &varsC1);
4985 
4986  return SCIP_OKAY;
4987 }
4988 
4989 /** separates lifted extended weight inequalities using sequential up- and down-lifting for given knapsack problem */
4990 static
4992  SCIP* scip, /**< SCIP data structure */
4993  SCIP_CONS* cons, /**< constraint that originates the knapsack problem, or NULL */
4994  SCIP_SEPA* sepa, /**< originating separator of the knapsack problem, or NULL */
4995  SCIP_VAR** vars, /**< variables in knapsack constraint */
4996  int nvars, /**< number of variables in knapsack constraint */
4997  int ntightened, /**< number of variables with tightened upper bound */
4998  SCIP_Longint* weights, /**< weights of variables in knapsack constraint */
4999  SCIP_Longint capacity, /**< capacity of knapsack */
5000  SCIP_Real* solvals, /**< solution values of all problem variables */
5001  int* feassetvars, /**< variables in feasible set */
5002  int* nonfeassetvars, /**< variables not in feasible set */
5003  int nfeassetvars, /**< number of variables in feasible set */
5004  int nnonfeassetvars, /**< number of variables not in feasible set */
5005  SCIP_SOL* sol, /**< primal SCIP solution to separate, NULL for current LP solution */
5006  SCIP_Bool* cutoff, /**< whether a cutoff has been detected */
5007  int* ncuts /**< pointer to add up the number of found cuts */
5008  )
5009 {
5010  int* varsT1;
5011  int* varsT2;
5012  int* varsF;
5013  int* varsR;
5014  int* liftcoefs;
5015  SCIP_Real cutact;
5016  int nvarsT1;
5017  int nvarsT2;
5018  int nvarsF;
5019  int nvarsR;
5020  int liftrhs;
5021  int j;
5022 
5023  assert( cutoff != NULL );
5024  *cutoff = FALSE;
5025 
5026  /* allocates temporary memory */
5027  SCIP_CALL( SCIPallocBufferArray(scip, &varsT1, nvars) );
5028  SCIP_CALL( SCIPallocBufferArray(scip, &varsT2, nvars) );
5029  SCIP_CALL( SCIPallocBufferArray(scip, &varsF, nvars) );
5030  SCIP_CALL( SCIPallocBufferArray(scip, &varsR, nvars) );
5031  SCIP_CALL( SCIPallocBufferArray(scip, &liftcoefs, nvars) );
5032 
5033  /* gets partition (T_1,T_2) of T, i.e. T_1 & T_2 = T and T_1 cap T_2 = emptyset, with T_1 not empty; chooses partition
5034  * as follows
5035  * T_2 = { j in T : x*_j = 1 } and
5036  * T_1 = T\T_2
5037  */
5038  getPartitionCovervars(scip, solvals, feassetvars, nfeassetvars, varsT1, varsT2, &nvarsT1, &nvarsT2);
5039  assert(nvarsT1 + nvarsT2 == nfeassetvars);
5040 
5041  /* changes partition (T_1,T_2) of feasible set T, if |T1| = 0, by moving one variable from T2 to T1 */
5042  if( nvarsT1 == 0 && nvarsT2 > 0)
5043  {
5044  SCIP_CALL( changePartitionFeasiblesetvars(scip, weights, varsT1, varsT2, &nvarsT1, &nvarsT2) );
5045  assert(nvarsT1 == 1);
5046  }
5047  assert(nvarsT2 == 0 || nvarsT1 > 0);
5048 
5049  /* gets partition (F,R) of N\T, i.e. F & R = N\T and F cap R = emptyset; chooses partition as follows
5050  * R = { j in N\T : x*_j = 0 } and
5051  * F = (N\T)\F
5052  */
5053  getPartitionNoncovervars(scip, solvals, nonfeassetvars, nnonfeassetvars, varsF, varsR, &nvarsF, &nvarsR);
5054  assert(nvarsF + nvarsR == nnonfeassetvars);
5055  assert(nvarsT1 + nvarsT2 + nvarsF + nvarsR == nvars - ntightened);
5056 
5057  /* sorts variables in F, T_2, and R according to the second level lifting sequence that will be used in the sequential
5058  * lifting procedure (the variable removed last from the initial cover does not have to be lifted first, therefore it
5059  * is included in the sorting routine)
5060  */
5061  SCIP_CALL( getLiftingSequence(scip, solvals, weights, varsF, varsT2, varsR, nvarsF, nvarsT2, nvarsR) );
5062 
5063  /* lifts extended weight inequality sum_{j in T_1} x_j <= |T_1| valid for
5064  *
5065  * S^0 = { x in {0,1}^|T_1| : sum_{j in T_1} a_j x_j <= a_0 - sum_{j in T_2} a_j }
5066  *
5067  * to a valid inequality sum_{j in T_1} x_j + sum_{j in N\T_1} alpha_j x_j <= |T_1| + sum_{j in T_2} alpha_j for
5068  *
5069  * S = { x in {0,1}^|N| : sum_{j in N} a_j x_j <= a_0 },
5070  *
5071  * uses sequential up-lifting for the variables in F, sequential down-lifting for the variable in T_2 and sequential
5072  * up-lifting for the variabels in R according to the second level lifting sequence
5073  */
5074  SCIP_CALL( sequentialUpAndDownLifting(scip, vars, nvars, ntightened, weights, capacity, solvals, varsT1, varsT2, varsF, varsR,
5075  nvarsT1, nvarsT2, nvarsF, nvarsR, nvarsT1, liftcoefs, &cutact, &liftrhs) );
5076 
5077  /* checks, if lifting yielded a violated cut */
5078  if( SCIPisEfficacious(scip, (cutact - liftrhs)/sqrt((SCIP_Real)MAX(liftrhs, 1))) )
5079  {
5080  SCIP_ROW* row;
5081  char name[SCIP_MAXSTRLEN];
5082 
5083  /* creates LP row */
5084  assert( cons == NULL || sepa == NULL );
5085  if( cons != NULL )
5086  {
5087  (void) SCIPsnprintf(name, SCIP_MAXSTRLEN, "%s_ewseq%" SCIP_LONGINT_FORMAT "", SCIPconsGetName(cons), SCIPconshdlrGetNCutsFound(SCIPconsGetHdlr(cons)));
5088  SCIP_CALL( SCIPcreateEmptyRowCons(scip, &row, SCIPconsGetHdlr(cons), name, -SCIPinfinity(scip), (SCIP_Real)liftrhs,
5089  cons != NULL ? SCIPconsIsLocal(cons) : FALSE, FALSE,
5090  cons != NULL ? SCIPconsIsRemovable(cons) : TRUE) );
5091  }
5092  else if ( sepa != NULL )
5093  {
5094  (void) SCIPsnprintf(name, SCIP_MAXSTRLEN, "%s_ewseq_%" SCIP_LONGINT_FORMAT "", SCIPsepaGetName(sepa), SCIPsepaGetNCutsFound(sepa));
5095  SCIP_CALL( SCIPcreateEmptyRowSepa(scip, &row, sepa, name, -SCIPinfinity(scip), (SCIP_Real)liftrhs, FALSE, FALSE, TRUE) );
5096  }
5097  else
5098  {
5099  (void) SCIPsnprintf(name, SCIP_MAXSTRLEN, "nn_ewseq_%" SCIP_LONGINT_FORMAT "", *ncuts);
5100  SCIP_CALL( SCIPcreateEmptyRowUnspec(scip, &row, name, -SCIPinfinity(scip), (SCIP_Real)liftrhs, FALSE, FALSE, TRUE) );
5101  }
5102 
5103  /* adds all variables in the knapsack constraint with calculated lifting coefficient to the cut */
5104  SCIP_CALL( SCIPcacheRowExtensions(scip, row) );
5105  assert(nvarsT1 + nvarsT2 + nvarsF + nvarsR == nvars - ntightened);
5106  for( j = 0; j < nvarsT1; j++ )
5107  {
5108  SCIP_CALL( SCIPaddVarToRow(scip, row, vars[varsT1[j]], 1.0) );
5109  }
5110  for( j = 0; j < nvarsT2; j++ )
5111  {
5112  if( liftcoefs[varsT2[j]] > 0 )
5113  {
5114  SCIP_CALL( SCIPaddVarToRow(scip, row, vars[varsT2[j]], (SCIP_Real)liftcoefs[varsT2[j]]) );
5115  }
5116  }
5117  for( j = 0; j < nvarsF; j++ )
5118  {
5119  if( liftcoefs[varsF[j]] > 0 )
5120  {
5121  SCIP_CALL( SCIPaddVarToRow(scip, row, vars[varsF[j]], (SCIP_Real)liftcoefs[varsF[j]]) );
5122  }
5123  }
5124  for( j = 0; j < nvarsR; j++ )
5125  {
5126  if( liftcoefs[varsR[j]] > 0 )
5127  {
5128  SCIP_CALL( SCIPaddVarToRow(scip, row, vars[varsR[j]], (SCIP_Real)liftcoefs[varsR[j]]) );
5129  }
5130  }
5131  SCIP_CALL( SCIPflushRowExtensions(scip, row) );
5132 
5133  /* checks, if cut is violated enough */
5134  if( SCIPisCutEfficacious(scip, sol, row) )
5135  {
5136  if( cons != NULL )
5137  {
5138  SCIP_CALL( SCIPresetConsAge(scip, cons) );
5139  }
5140  SCIP_CALL( SCIPaddCut(scip, sol, row, FALSE, cutoff) );
5141  (*ncuts)++;
5142  }
5143  SCIP_CALL( SCIPreleaseRow(scip, &row) );
5144  }
5145 
5146  /* frees temporary memory */
5147  SCIPfreeBufferArray(scip, &liftcoefs);
5148  SCIPfreeBufferArray(scip, &varsR);
5149  SCIPfreeBufferArray(scip, &varsF);
5150  SCIPfreeBufferArray(scip, &varsT2);
5151  SCIPfreeBufferArray(scip, &varsT1);
5152 
5153  return SCIP_OKAY;
5154 }
5155 
5156 /** separates lifted minimal cover inequalities using superadditive up-lifting for given knapsack problem */
5157 static
5159  SCIP* scip, /**< SCIP data structure */
5160  SCIP_CONS* cons, /**< constraint that originates the knapsack problem, or NULL */
5161  SCIP_SEPA* sepa, /**< originating separator of the knapsack problem, or NULL */
5162  SCIP_VAR** vars, /**< variables in knapsack constraint */
5163  int nvars, /**< number of variables in knapsack constraint */
5164  int ntightened, /**< number of variables with tightened upper bound */
5165  SCIP_Longint* weights, /**< weights of variables in knapsack constraint */
5166  SCIP_Longint capacity, /**< capacity of knapsack */
5167  SCIP_Real* solvals, /**< solution values of all problem variables */
5168  int* mincovervars, /**< mincover variables */
5169  int* nonmincovervars, /**< nonmincover variables */
5170  int nmincovervars, /**< number of mincover variables */
5171  int nnonmincovervars, /**< number of nonmincover variables */
5172  SCIP_Longint mincoverweight, /**< weight of minimal cover */
5173  SCIP_SOL* sol, /**< primal SCIP solution to separate, NULL for current LP solution */
5174  SCIP_Bool* cutoff, /**< whether a cutoff has been detected */
5175  int* ncuts /**< pointer to add up the number of found cuts */
5176  )
5177 {
5178  SCIP_Real* realliftcoefs;
5179  SCIP_Real cutact;
5180  int liftrhs;
5181 
5182  assert( cutoff != NULL );
5183  *cutoff = FALSE;
5184  cutact = 0.0;
5185 
5186  /* allocates temporary memory */
5187  SCIP_CALL( SCIPallocBufferArray(scip, &realliftcoefs, nvars) );
5188 
5189  /* lifts minimal cover inequality sum_{j in C} x_j <= |C| - 1 valid for
5190  *
5191  * S^0 = { x in {0,1}^|C| : sum_{j in C} a_j x_j <= a_0 }
5192  *
5193  * to a valid inequality sum_{j in C} x_j + sum_{j in N\C} alpha_j x_j <= |C| - 1 for
5194  *
5195  * S = { x in {0,1}^|N| : sum_{j in N} a_j x_j <= a_0 },
5196  *
5197  * uses superadditive up-lifting for the variables in N\C.
5198  */
5199  SCIP_CALL( superadditiveUpLifting(scip, vars, nvars, ntightened, weights, capacity, solvals, mincovervars,
5200  nonmincovervars, nmincovervars, nnonmincovervars, mincoverweight, realliftcoefs, &cutact) );
5201  liftrhs = nmincovervars - 1;
5202 
5203  /* checks, if lifting yielded a violated cut */
5204  if( SCIPisEfficacious(scip, (cutact - liftrhs)/sqrt((SCIP_Real)MAX(liftrhs, 1))) )
5205  {
5206  SCIP_ROW* row;
5207  char name[SCIP_MAXSTRLEN];
5208  int j;
5209 
5210  /* creates LP row */
5211  assert( cons == NULL || sepa == NULL );
5212  if ( cons != NULL )
5213  {
5214  (void) SCIPsnprintf(name, SCIP_MAXSTRLEN, "%s_mcsup%" SCIP_LONGINT_FORMAT "", SCIPconsGetName(cons), SCIPconshdlrGetNCutsFound(SCIPconsGetHdlr(cons)));
5215  SCIP_CALL( SCIPcreateEmptyRowCons(scip, &row, SCIPconsGetHdlr(cons), name, -SCIPinfinity(scip), (SCIP_Real)liftrhs,
5216  cons != NULL ? SCIPconsIsLocal(cons) : FALSE, FALSE,
5217  cons != NULL ? SCIPconsIsRemovable(cons) : TRUE) );
5218  }
5219  else if ( sepa != NULL )
5220  {
5221  (void) SCIPsnprintf(name, SCIP_MAXSTRLEN, "%s_mcsup%" SCIP_LONGINT_FORMAT "", SCIPsepaGetName(sepa), SCIPsepaGetNCutsFound(sepa));
5222  SCIP_CALL( SCIPcreateEmptyRowSepa(scip, &row, sepa, name, -SCIPinfinity(scip), (SCIP_Real)liftrhs, FALSE, FALSE, TRUE) );
5223  }
5224  else
5225  {
5226  (void) SCIPsnprintf(name, SCIP_MAXSTRLEN, "nn_mcsup_%" SCIP_LONGINT_FORMAT "", *ncuts);
5227  SCIP_CALL( SCIPcreateEmptyRowUnspec(scip, &row, name, -SCIPinfinity(scip), (SCIP_Real)liftrhs, FALSE, FALSE, TRUE) );
5228  }
5229 
5230  /* adds all variables in the knapsack constraint with calculated lifting coefficient to the cut */
5231  SCIP_CALL( SCIPcacheRowExtensions(scip, row) );
5232  assert(nmincovervars + nnonmincovervars == nvars - ntightened);
5233  for( j = 0; j < nmincovervars; j++ )
5234  {
5235  SCIP_CALL( SCIPaddVarToRow(scip, row, vars[mincovervars[j]], 1.0) );
5236  }
5237  for( j = 0; j < nnonmincovervars; j++ )
5238  {
5239  assert(SCIPisFeasGE(scip, realliftcoefs[nonmincovervars[j]], 0.0));
5240  if( SCIPisFeasGT(scip, realliftcoefs[nonmincovervars[j]], 0.0) )
5241  {
5242  SCIP_CALL( SCIPaddVarToRow(scip, row, vars[nonmincovervars[j]], realliftcoefs[nonmincovervars[j]]) );
5243  }
5244  }
5245  SCIP_CALL( SCIPflushRowExtensions(scip, row) );
5246 
5247  /* checks, if cut is violated enough */
5248  if( SCIPisCutEfficacious(scip, sol, row) )
5249  {
5250  if( cons != NULL )
5251  {
5252  SCIP_CALL( SCIPresetConsAge(scip, cons) );
5253  }
5254  SCIP_CALL( SCIPaddCut(scip, sol, row, FALSE, cutoff) );
5255  (*ncuts)++;
5256  }
5257  SCIP_CALL( SCIPreleaseRow(scip, &row) );
5258  }
5259 
5260  /* frees temporary memory */
5261  SCIPfreeBufferArray(scip, &realliftcoefs);
5262 
5263  return SCIP_OKAY;
5264 }
5265 
5266 /** converts given cover C to a minimal cover by removing variables in the reverse order in which the variables were chosen
5267  * to be in C, i.e. in the order of non-increasing (1 - x*_j)/a_j, if the transformed separation problem was used to find
5268  * C and in the order of non-increasing (1 - x*_j), if the modified transformed separation problem was used to find C;
5269  * note that all variables with x*_j = 1 will be removed last
5270  */
5271 static
5273  SCIP* scip, /**< SCIP data structure */
5274  SCIP_Longint* weights, /**< weights of variables in knapsack constraint */
5275  SCIP_Longint capacity, /**< capacity of knapsack */
5276  SCIP_Real* solvals, /**< solution values of all problem variables */
5277  int* covervars, /**< pointer to store cover variables */
5278  int* noncovervars, /**< pointer to store noncover variables */
5279  int* ncovervars, /**< pointer to store number of cover variables */
5280  int* nnoncovervars, /**< pointer to store number of noncover variables */
5281  SCIP_Longint* coverweight, /**< pointer to store weight of cover */
5282  SCIP_Bool modtransused /**< TRUE if mod trans sepa prob was used to find cover */
5283  )
5284 {
5285  SORTKEYPAIR** sortkeypairs;
5286  SCIP_Longint minweight;
5287  int nsortkeypairs;
5288  int minweightidx;
5289  int j;
5290  int k;
5291 
5292  assert(scip != NULL);
5293  assert(covervars != NULL);
5294  assert(noncovervars != NULL);
5295  assert(ncovervars != NULL);
5296  assert(*ncovervars > 0);
5297  assert(nnoncovervars != NULL);
5298  assert(*nnoncovervars >= 0);
5299  assert(coverweight != NULL);
5300  assert(*coverweight > 0);
5301  assert(*coverweight > capacity);
5302 
5303  /* allocates temporary memory */
5304  nsortkeypairs = *ncovervars;
5305  SCIP_CALL( SCIPallocBufferArray(scip, &sortkeypairs, nsortkeypairs) );
5306 
5307  /* sorts C in the reverse order in which the variables were chosen to be in the cover, i.e.
5308  * such that (1 - x*_1)/a_1 >= ... >= (1 - x*_|C|)/a_|C|, if trans separation problem was used to find C
5309  * such that (1 - x*_1) >= ... >= (1 - x*_|C|), if modified trans separation problem was used to find C
5310  * note that all variables with x*_j = 1 are in the end of the sorted C, so they will be removed last from C
5311  */
5312  assert(*ncovervars == nsortkeypairs);
5313  if( modtransused )
5314  {
5315  for( j = 0; j < *ncovervars; j++ )
5316  {
5317  SCIP_CALL( SCIPallocBuffer(scip, &(sortkeypairs[j])) ); /*lint !e866 */
5318 
5319  sortkeypairs[j]->key1 = solvals[covervars[j]];
5320  sortkeypairs[j]->key2 = (SCIP_Real) weights[covervars[j]];
5321  }
5322  }
5323  else
5324  {
5325  for( j = 0; j < *ncovervars; j++ )
5326  {
5327  SCIP_CALL( SCIPallocBuffer(scip, &(sortkeypairs[j])) ); /*lint !e866 */
5328 
5329  sortkeypairs[j]->key1 = (solvals[covervars[j]] - 1.0) / ((SCIP_Real) weights[covervars[j]]);
5330  sortkeypairs[j]->key2 = (SCIP_Real) (-weights[covervars[j]]);
5331  }
5332  }
5333  SCIPsortPtrInt((void**)sortkeypairs, covervars, compSortkeypairs, *ncovervars);
5334 
5335  /* gets j' with a_j' = min{ a_j : j in C } */
5336  minweightidx = 0;
5337  minweight = weights[covervars[minweightidx]];
5338  for( j = 1; j < *ncovervars; j++ )
5339  {
5340  if( weights[covervars[j]] <= minweight )
5341  {
5342  minweightidx = j;
5343  minweight = weights[covervars[minweightidx]];
5344  }
5345  }
5346  assert(minweightidx >= 0 && minweightidx < *ncovervars);
5347  assert(minweight > 0 && minweight <= *coverweight);
5348 
5349  j = 0;
5350  /* removes variables from C until the remaining variables form a minimal cover */
5351  while( j < *ncovervars && ((*coverweight) - minweight > capacity) )
5352  {
5353  assert(minweightidx >= j);
5354  assert(checkMinweightidx(weights, capacity, covervars, *ncovervars, *coverweight, minweightidx, j));
5355 
5356  /* if sum_{i in C} a_i - a_j <= a_0, j cannot be removed from C */
5357  if( (*coverweight) - weights[covervars[j]] <= capacity )
5358  {
5359  ++j;
5360  continue;
5361  }
5362 
5363  /* adds j to N\C */
5364  noncovervars[*nnoncovervars] = covervars[j];
5365  (*nnoncovervars)++;
5366 
5367  /* removes j from C */
5368  (*coverweight) -= weights[covervars[j]];
5369  for( k = j; k < (*ncovervars) - 1; k++ )
5370  covervars[k] = covervars[k+1];
5371  (*ncovervars)--;
5372 
5373  /* updates j' with a_j' = min{ a_j : j in C } */
5374  if( j == minweightidx )
5375  {
5376  minweightidx = 0;
5377  minweight = weights[covervars[minweightidx]];
5378  for( k = 1; k < *ncovervars; k++ )
5379  {
5380  if( weights[covervars[k]] <= minweight )
5381  {
5382  minweightidx = k;
5383  minweight = weights[covervars[minweightidx]];
5384  }
5385  }
5386  assert(minweight > 0 && minweight <= *coverweight);
5387  assert(minweightidx >= 0 && minweightidx < *ncovervars);
5388  }
5389  else
5390  {
5391  assert(minweightidx > j);
5392  minweightidx--;
5393  }
5394  /* j needs to stay the same */
5395  }
5396  assert((*coverweight) > capacity);
5397  assert((*coverweight) - minweight <= capacity);
5398 
5399  /* frees temporary memory */
5400  for( j = nsortkeypairs-1; j >= 0; j-- )
5401  SCIPfreeBuffer(scip, &(sortkeypairs[j])); /*lint !e866 */
5402  SCIPfreeBufferArray(scip, &sortkeypairs);
5403 
5404  return SCIP_OKAY;
5405 }
5406 
5407 /** converts given initial cover C_init to a feasible set by removing variables in the reverse order in which
5408  * they were chosen to be in C_init:
5409  * non-increasing (1 - x*_j)/a_j, if transformed separation problem was used to find C_init
5410  * non-increasing (1 - x*_j), if modified transformed separation problem was used to find C_init.
5411  * separates lifted extended weight inequalities using sequential up- and down-lifting for this feasible set
5412  * and all subsequent feasible sets.
5413  */
5414 static
5416  SCIP* scip, /**< SCIP data structure */
5417  SCIP_CONS* cons, /**< constraint that originates the knapsack problem */
5418  SCIP_SEPA* sepa, /**< originating separator of the knapsack problem, or NULL */
5419  SCIP_VAR** vars, /**< variables in knapsack constraint */
5420  int nvars, /**< number of variables in knapsack constraint */
5421  int ntightened, /**< number of variables with tightened upper bound */
5422  SCIP_Longint* weights, /**< weights of variables in knapsack constraint */
5423  SCIP_Longint capacity, /**< capacity of knapsack */
5424  SCIP_Real* solvals, /**< solution values of all problem variables */
5425  int* covervars, /**< pointer to store cover variables */
5426  int* noncovervars, /**< pointer to store noncover variables */
5427  int* ncovervars, /**< pointer to store number of cover variables */
5428  int* nnoncovervars, /**< pointer to store number of noncover variables */
5429  SCIP_Longint* coverweight, /**< pointer to store weight of cover */
5430  SCIP_Bool modtransused, /**< TRUE if mod trans sepa prob was used to find cover */
5431  SCIP_SOL* sol, /**< primal SCIP solution to separate, NULL for current LP solution */
5432  SCIP_Bool* cutoff, /**< whether a cutoff has been detected */
5433  int* ncuts /**< pointer to add up the number of found cuts */
5434  )
5435 {
5436  SCIP_Real* sortkeys;
5437  int j;
5438  int k;
5439 
5440  assert(scip != NULL);
5441  assert(covervars != NULL);
5442  assert(noncovervars != NULL);
5443  assert(ncovervars != NULL);
5444  assert(*ncovervars > 0);
5445  assert(nnoncovervars != NULL);
5446  assert(*nnoncovervars >= 0);
5447  assert(coverweight != NULL);
5448  assert(*coverweight > 0);
5449  assert(*coverweight > capacity);
5450  assert(*ncovervars + *nnoncovervars == nvars - ntightened);
5451  assert(cutoff != NULL);
5452 
5453  *cutoff = FALSE;
5454 
5455  /* allocates temporary memory */
5456  SCIP_CALL( SCIPallocBufferArray(scip, &sortkeys, *ncovervars) );
5457 
5458  /* sorts C in the reverse order in which the variables were chosen to be in the cover, i.e.
5459  * such that (1 - x*_1)/a_1 >= ... >= (1 - x*_|C|)/a_|C|, if trans separation problem was used to find C
5460  * such that (1 - x*_1) >= ... >= (1 - x*_|C|), if modified trans separation problem was used to find C
5461  * note that all variables with x*_j = 1 are in the end of the sorted C, so they will be removed last from C
5462  */
5463  if( modtransused )
5464  {
5465  for( j = 0; j < *ncovervars; j++ )
5466  {
5467  sortkeys[j] = solvals[covervars[j]];
5468  assert(SCIPisFeasGE(scip, sortkeys[j], 0.0));
5469  }
5470  }
5471  else
5472  {
5473  for( j = 0; j < *ncovervars; j++ )
5474  {
5475  sortkeys[j] = (solvals[covervars[j]] - 1.0) / ((SCIP_Real) weights[covervars[j]]);
5476  assert(SCIPisFeasLE(scip, sortkeys[j], 0.0));
5477  }
5478  }
5479  SCIPsortRealInt(sortkeys, covervars, *ncovervars);
5480 
5481  /* removes variables from C_init and separates lifted extended weight inequalities using sequential up- and down-lifting;
5482  * in addition to an extended weight inequality this gives cardinality inequalities */
5483  while( *ncovervars >= 2 )
5484  {
5485  /* adds first element of C_init to N\C_init */
5486  noncovervars[*nnoncovervars] = covervars[0];
5487  (*nnoncovervars)++;
5488 
5489  /* removes first element from C_init */
5490  (*coverweight) -= weights[covervars[0]];
5491  for( k = 0; k < (*ncovervars) - 1; k++ )
5492  covervars[k] = covervars[k+1];
5493  (*ncovervars)--;
5494 
5495  assert(*ncovervars + *nnoncovervars == nvars - ntightened);
5496  if( (*coverweight) <= capacity )
5497  {
5498  SCIP_CALL( separateSequLiftedExtendedWeightInequality(scip, cons, sepa, vars, nvars, ntightened, weights, capacity, solvals,
5499  covervars, noncovervars, *ncovervars, *nnoncovervars, sol, cutoff, ncuts) );
5500  }
5501 
5502  /* stop if cover is too large */
5503  if ( *ncovervars >= MAXCOVERSIZEITERLEWI )
5504  break;
5505  }
5506 
5507  /* frees temporary memory */
5508  SCIPfreeBufferArray(scip, &sortkeys);
5509 
5510  return SCIP_OKAY;
5511 }
5512 
5513 /** separates different classes of valid inequalities for the 0-1 knapsack problem */
5515  SCIP* scip, /**< SCIP data structure */
5516  SCIP_CONS* cons, /**< originating constraint of the knapsack problem, or NULL */
5517  SCIP_SEPA* sepa, /**< originating separator of the knapsack problem, or NULL */
5518  SCIP_VAR** vars, /**< variables in knapsack constraint */
5519  int nvars, /**< number of variables in knapsack constraint */
5520  SCIP_Longint* weights, /**< weights of variables in knapsack constraint */
5521  SCIP_Longint capacity, /**< capacity of knapsack */
5522  SCIP_SOL* sol, /**< primal SCIP solution to separate, NULL for current LP solution */
5523  SCIP_Bool usegubs, /**< should GUB information be used for separation? */
5524  SCIP_Bool* cutoff, /**< pointer to store whether a cutoff has been detected */
5525  int* ncuts /**< pointer to add up the number of found cuts */
5526  )
5527 {
5528  SCIP_Real* solvals;
5529  int* covervars;
5530  int* noncovervars;
5531  SCIP_Bool coverfound;
5532  SCIP_Bool fractional;
5533  SCIP_Bool modtransused;
5534  SCIP_Longint coverweight;
5535  int ncovervars;
5536  int nnoncovervars;
5537  int ntightened;
5538 
5539  assert(scip != NULL);
5540  assert(capacity >= 0);
5541  assert(cutoff != NULL);
5542  assert(ncuts != NULL);
5543 
5544  *cutoff = FALSE;
5545 
5546  if( nvars == 0 )
5547  return SCIP_OKAY;
5548 
5549  assert(vars != NULL);
5550  assert(nvars > 0);
5551  assert(weights != NULL);
5552 
5553  /* increase age of constraint (age is reset to zero, if a cut was found) */
5554  if( cons != NULL )
5555  {
5556  SCIP_CALL( SCIPincConsAge(scip, cons) );
5557  }
5558 
5559  /* allocates temporary memory */
5560  SCIP_CALL( SCIPallocBufferArray(scip, &solvals, nvars) );
5561  SCIP_CALL( SCIPallocBufferArray(scip, &covervars, nvars) );
5562  SCIP_CALL( SCIPallocBufferArray(scip, &noncovervars, nvars) );
5563 
5564  /* gets solution values of all problem variables */
5565  SCIP_CALL( SCIPgetSolVals(scip, sol, nvars, vars, solvals) );
5566 
5567 #ifdef SCIP_DEBUG
5568  {
5569  int i;
5570 
5571  SCIPdebugMsg(scip, "separate cuts for knapsack constraint originated by cons <%s>:\n",
5572  cons == NULL ? "-" : SCIPconsGetName(cons));
5573  for( i = 0; i < nvars; ++i )
5574  {
5575  SCIPdebugMsgPrint(scip, "%+" SCIP_LONGINT_FORMAT "<%s>(%g)", weights[i], SCIPvarGetName(vars[i]), solvals[i]);
5576  }
5577  SCIPdebugMsgPrint(scip, " <= %" SCIP_LONGINT_FORMAT "\n", capacity);
5578  }
5579 #endif
5580 
5581  /* LMCI1 (lifted minimal cover inequalities using sequential up- and down-lifting) using GUB information
5582  */
5583  if( usegubs )
5584  {
5585  SCIP_GUBSET* gubset;
5586 
5587  SCIPdebugMsg(scip, "separate LMCI1-GUB cuts:\n");
5588 
5589  /* initializes partion of knapsack variables into nonoverlapping GUB constraints */
5590  SCIP_CALL( GUBsetCreate(scip, &gubset, nvars, weights, capacity) );
5591 
5592  /* constructs sophisticated partition of knapsack variables into nonoverlapping GUBs */
5593  SCIP_CALL( GUBsetGetCliquePartition(scip, gubset, vars, solvals) );
5594  assert(gubset->ngubconss <= nvars);
5595 
5596  /* gets a most violated initial cover C_init ( sum_{j in C_init} a_j > a_0 ) by using the
5597  * MODIFIED transformed separation problem and taking into account the following fixing:
5598  * j in C_init, if j in N_1 = {j in N : x*_j = 1} and
5599  * j in N\C_init, if j in N_0 = {j in N : x*_j = 0},
5600  * if one exists
5601  */
5602  modtransused = TRUE;
5603  SCIP_CALL( getCover(scip, vars, nvars, weights, capacity, solvals, covervars, noncovervars, &ncovervars,
5604  &nnoncovervars, &coverweight, &coverfound, modtransused, &ntightened, &fractional) );
5605 
5606  assert(!coverfound || !fractional || ncovervars + nnoncovervars == nvars - ntightened);
5607 
5608  /* if x* is not fractional we stop the separation routine */
5609  if( !fractional )
5610  {
5611  SCIPdebugMsg(scip, " LMCI1-GUB terminated by no variable with fractional LP value.\n");
5612 
5613  /* frees memory for GUB set data structure */
5614  SCIP_CALL( GUBsetFree(scip, &gubset) );
5615 
5616  goto TERMINATE;
5617  }
5618 
5619  /* if no cover was found we stop the separation routine for lifted minimal cover inequality */
5620  if( coverfound )
5621  {
5622  /* converts initial cover C_init to a minimal cover C by removing variables in the reverse order in which the
5623  * variables were chosen to be in C_init; note that variables with x*_j = 1 will be removed last
5624  */
5625  SCIP_CALL( makeCoverMinimal(scip, weights, capacity, solvals, covervars, noncovervars, &ncovervars,
5626  &nnoncovervars, &coverweight, modtransused) );
5627 
5628  /* only separate with GUB information if we have at least one nontrivial GUB (with more than one variable) */
5629  if( gubset->ngubconss < nvars )
5630  {
5631  /* separates lifted minimal cover inequalities using sequential up- and down-lifting and GUB information */
5632  SCIP_CALL( separateSequLiftedMinimalCoverInequality(scip, cons, sepa, vars, nvars, ntightened, weights, capacity,
5633  solvals, covervars, noncovervars, ncovervars, nnoncovervars, sol, gubset, cutoff, ncuts) );
5634  }
5635  else
5636  {
5637  /* separates lifted minimal cover inequalities using sequential up- and down-lifting, but do not use trivial
5638  * GUB information
5639  */
5640  SCIP_CALL( separateSequLiftedMinimalCoverInequality(scip, cons, sepa, vars, nvars, ntightened, weights, capacity,
5641  solvals, covervars, noncovervars, ncovervars, nnoncovervars, sol, NULL, cutoff, ncuts) );
5642  }
5643  }
5644 
5645  /* frees memory for GUB set data structure */
5646  SCIP_CALL( GUBsetFree(scip, &gubset) );
5647  }
5648  else
5649  {
5650  /* LMCI1 (lifted minimal cover inequalities using sequential up- and down-lifting)
5651  * (and LMCI2 (lifted minimal cover inequalities using superadditive up-lifting))
5652  */
5653 
5654  /* gets a most violated initial cover C_init ( sum_{j in C_init} a_j > a_0 ) by using the
5655  * MODIFIED transformed separation problem and taking into account the following fixing:
5656  * j in C_init, if j in N_1 = {j in N : x*_j = 1} and
5657  * j in N\C_init, if j in N_0 = {j in N : x*_j = 0},
5658  * if one exists
5659  */
5660  SCIPdebugMsg(scip, "separate LMCI1 cuts:\n");
5661  modtransused = TRUE;
5662  SCIP_CALL( getCover(scip, vars, nvars, weights, capacity, solvals, covervars, noncovervars, &ncovervars,
5663  &nnoncovervars, &coverweight, &coverfound, modtransused, &ntightened, &fractional) );
5664  assert(!coverfound || !fractional || ncovervars + nnoncovervars == nvars - ntightened);
5665 
5666  /* if x* is not fractional we stop the separation routine */
5667  if( !fractional )
5668  goto TERMINATE;
5669 
5670  /* if no cover was found we stop the separation routine for lifted minimal cover inequality */
5671  if( coverfound )
5672  {
5673  /* converts initial cover C_init to a minimal cover C by removing variables in the reverse order in which the
5674  * variables were chosen to be in C_init; note that variables with x*_j = 1 will be removed last
5675  */
5676  SCIP_CALL( makeCoverMinimal(scip, weights, capacity, solvals, covervars, noncovervars, &ncovervars,
5677  &nnoncovervars, &coverweight, modtransused) );
5678 
5679  /* separates lifted minimal cover inequalities using sequential up- and down-lifting */
5680  SCIP_CALL( separateSequLiftedMinimalCoverInequality(scip, cons, sepa, vars, nvars, ntightened, weights, capacity,
5681  solvals, covervars, noncovervars, ncovervars, nnoncovervars, sol, NULL, cutoff, ncuts) );
5682 
5683  if( USESUPADDLIFT ) /*lint !e506 !e774*/
5684  {
5685  SCIPdebugMsg(scip, "separate LMCI2 cuts:\n");
5686  /* separates lifted minimal cover inequalities using superadditive up-lifting */
5687  SCIP_CALL( separateSupLiftedMinimalCoverInequality(scip, cons, sepa, vars, nvars, ntightened, weights, capacity,
5688  solvals, covervars, noncovervars, ncovervars, nnoncovervars, coverweight, sol, cutoff, ncuts) );
5689  }
5690  }
5691  }
5692 
5693  /* LEWI (lifted extended weight inequalities using sequential up- and down-lifting) */
5694  if ( ! (*cutoff) )
5695  {
5696  /* gets a most violated initial cover C_init ( sum_{j in C_init} a_j > a_0 ) by using the
5697  * transformed separation problem and taking into account the following fixing:
5698  * j in C_init, if j in N_1 = {j in N : x*_j = 1} and
5699  * j in N\C_init, if j in N_0 = {j in N : x*_j = 0},
5700  * if one exists
5701  */
5702  SCIPdebugMsg(scip, "separate LEWI cuts:\n");
5703  modtransused = FALSE;
5704  SCIP_CALL( getCover(scip, vars, nvars, weights, capacity, solvals, covervars, noncovervars, &ncovervars,
5705  &nnoncovervars, &coverweight, &coverfound, modtransused, &ntightened, &fractional) );
5706  assert(fractional);
5707  assert(!coverfound || ncovervars + nnoncovervars == nvars - ntightened);
5708 
5709  /* if no cover was found we stop the separation routine */
5710  if( coverfound )
5711  {
5712  /* converts initial cover C_init to a feasible set by removing variables in the reverse order in which
5713  * they were chosen to be in C_init and separates lifted extended weight inequalities using sequential
5714  * up- and down-lifting for this feasible set and all subsequent feasible sets.
5715  */
5716  SCIP_CALL( getFeasibleSet(scip, cons, sepa, vars, nvars, ntightened, weights, capacity, solvals, covervars, noncovervars,
5717  &ncovervars, &nnoncovervars, &coverweight, modtransused, sol, cutoff, ncuts) );
5718  }
5719  }
5720 
5721  TERMINATE:
5722  /* frees temporary memory */
5723  SCIPfreeBufferArray(scip, &noncovervars);
5724  SCIPfreeBufferArray(scip, &covervars);
5725  SCIPfreeBufferArray(scip, &solvals);
5726 
5727  return SCIP_OKAY;
5728 }
5729 
5730 /* relaxes given general linear constraint into a knapsack constraint and separates lifted knapsack cover inequalities */
5732  SCIP* scip, /**< SCIP data structure */
5733  SCIP_CONS* cons, /**< originating constraint of the knapsack problem, or NULL */
5734  SCIP_SEPA* sepa, /**< originating separator of the knapsack problem, or NULL */
5735  int nknapvars, /**< number of variables in the continuous knapsack constraint */
5736  SCIP_VAR** knapvars, /**< variables in the continuous knapsack constraint */
5737  SCIP_Real* knapvals, /**< coefficients of the variables in the continuous knapsack constraint */
5738  SCIP_Real valscale, /**< -1.0 if lhs of row is used as rhs of c. k. constraint, +1.0 otherwise */
5739  SCIP_Real rhs, /**< right hand side of the continuous knapsack constraint */
5740  SCIP_SOL* sol, /**< primal CIP solution, NULL for current LP solution */
5741  SCIP_Bool* cutoff, /**< pointer to store whether a cutoff was found */
5742  int* ncuts /**< pointer to add up the number of found cuts */
5743  )
5744 {
5745  SCIP_VAR** binvars;
5746  SCIP_VAR** consvars;
5747  SCIP_Real* binvals;
5748  SCIP_Longint* consvals;
5749  SCIP_Longint minact;
5750  SCIP_Longint maxact;
5751  SCIP_Real intscalar;
5752  SCIP_Bool success;
5753  int nbinvars;
5754  int nconsvars;
5755  int i;
5756 
5757  int* tmpindices;
5758  int tmp;
5759  SCIP_CONSHDLR* conshdlr;
5760  SCIP_CONSHDLRDATA* conshdlrdata;
5761  SCIP_Bool noknapsackconshdlr;
5762  SCIP_Bool usegubs;
5763 
5764  assert(nknapvars > 0);
5765  assert(knapvars != NULL);
5766  assert(cutoff != NULL);
5767 
5768  tmpindices = NULL;
5769 
5770  SCIPdebugMsg(scip, "separate linear constraint <%s> relaxed to knapsack\n", cons != NULL ? SCIPconsGetName(cons) : "-");
5771  SCIPdebug( if( cons != NULL ) { SCIPdebugPrintCons(scip, cons, NULL); } );
5772 
5773  binvars = SCIPgetVars(scip);
5774 
5775  /* all variables which are of integral type can be potentially of binary type; this can be checked via the method SCIPvarIsBinary(var) */
5776  nbinvars = SCIPgetNVars(scip) - SCIPgetNContVars(scip);
5777 
5778  *cutoff = FALSE;
5779 
5780  if( nbinvars == 0 )
5781  return SCIP_OKAY;
5782 
5783  /* set up data structures */
5784  SCIP_CALL( SCIPallocBufferArray(scip, &consvars, nbinvars) );
5785  SCIP_CALL( SCIPallocBufferArray(scip, &consvals, nbinvars) );
5786 
5787  /* get conshdlrdata to use cleared memory */
5788  conshdlr = SCIPfindConshdlr(scip, CONSHDLR_NAME);
5789  if( conshdlr == NULL )
5790  {
5791  noknapsackconshdlr = TRUE;
5792  usegubs = DEFAULT_USEGUBS;
5793 
5794  SCIP_CALL( SCIPallocBufferArray(scip, &binvals, nbinvars) );
5795  BMSclearMemoryArray(binvals, nbinvars);
5796  }
5797  else
5798  {
5799  noknapsackconshdlr = FALSE;
5800  conshdlrdata = SCIPconshdlrGetData(conshdlr);
5801  assert(conshdlrdata != NULL);
5802  usegubs = conshdlrdata->usegubs;
5803 
5804  SCIP_CALL( SCIPallocBufferArray(scip, &tmpindices, nknapvars) );
5805 
5806  /* increase array size to avoid an endless loop in the next block; this might happen if continuous variables
5807  * change their types to SCIP_VARTYPE_BINARY during presolving
5808  */
5809  if( conshdlrdata->reals1size == 0 )
5810  {
5811  SCIP_CALL( SCIPreallocBlockMemoryArray(scip, &conshdlrdata->reals1, conshdlrdata->reals1size, 1) );
5812  conshdlrdata->reals1size = 1;
5813  conshdlrdata->reals1[0] = 0.0;
5814  }
5815 
5816  assert(conshdlrdata->reals1size > 0);
5817 
5818  /* next if condition should normally not be true, because it means that presolving has created more binary
5819  * variables than binary + integer variables existed at the constraint initialization method, but for example if you would
5820  * transform all integers into their binary representation then it maybe happens
5821  */
5822  if( conshdlrdata->reals1size < nbinvars )
5823  {
5824  int oldsize = conshdlrdata->reals1size;
5825 
5826  conshdlrdata->reals1size = nbinvars;
5827  SCIP_CALL( SCIPreallocBlockMemoryArray(scip, &conshdlrdata->reals1, oldsize, conshdlrdata->reals1size) );
5828  BMSclearMemoryArray(&(conshdlrdata->reals1[oldsize]), conshdlrdata->reals1size - oldsize); /*lint !e866 */
5829  }
5830  binvals = conshdlrdata->reals1;
5831 
5832  /* check for cleared array, all entries have to be zero */
5833 #ifndef NDEBUG
5834  for( tmp = nbinvars - 1; tmp >= 0; --tmp )
5835  {
5836  assert(binvals[tmp] == 0);
5837  }
5838 #endif
5839  }
5840 
5841  tmp = 0;
5842 
5843  /* relax continuous knapsack constraint:
5844  * 1. make all variables binary:
5845  * if x_j is continuous or integer variable substitute:
5846  * - a_j < 0: x_j = lb or x_j = b*z + d with variable lower bound b*z + d with binary variable z
5847  * - a_j > 0: x_j = ub or x_j = b*z + d with variable upper bound b*z + d with binary variable z
5848  * 2. convert coefficients of all variables to positive integers:
5849  * - scale all coefficients a_j to a~_j integral
5850  * - substitute x~_j = 1 - x_j if a~_j < 0
5851  */
5852 
5853  /* replace integer and continuous variables with binary variables */
5854  for( i = 0; i < nknapvars; i++ )
5855  {
5856  SCIP_VAR* var;
5857 
5858  var = knapvars[i];
5859 
5860  if( SCIPvarIsBinary(var) && SCIPvarIsActive(var) )
5861  {
5862  assert(0 <= SCIPvarGetProbindex(var) && SCIPvarGetProbindex(var) < nbinvars);
5863  binvals[SCIPvarGetProbindex(var)] += valscale * knapvals[i];
5864  if( !noknapsackconshdlr )
5865  {
5866  assert(tmpindices != NULL);
5867 
5868  tmpindices[tmp] = SCIPvarGetProbindex(var);
5869  ++tmp;
5870  }
5871  SCIPdebugMsg(scip, " -> binary variable %+.15g<%s>(%.15g)\n", valscale * knapvals[i], SCIPvarGetName(var), SCIPgetSolVal(scip, sol, var));
5872  }
5873  else if( valscale * knapvals[i] > 0.0 )
5874  {
5875  SCIP_VAR** zvlb;
5876  SCIP_Real* bvlb;
5877  SCIP_Real* dvlb;
5878  SCIP_Real bestlbsol;
5879  int bestlbtype;
5880  int nvlb;
5881  int j;
5882 
5883  /* a_j > 0: substitution with lb or vlb */
5884  nvlb = SCIPvarGetNVlbs(var);
5885  zvlb = SCIPvarGetVlbVars(var);
5886  bvlb = SCIPvarGetVlbCoefs(var);
5887  dvlb = SCIPvarGetVlbConstants(var);
5888 
5889  /* search for lb or vlb with maximal bound value */
5890  bestlbsol = SCIPvarGetLbGlobal(var);
5891  bestlbtype = -1;
5892  for( j = 0; j < nvlb; j++ )
5893  {
5894  /* use only numerical stable vlb with binary variable z */
5895  if( SCIPvarIsBinary(zvlb[j]) && SCIPvarIsActive(zvlb[j]) && REALABS(bvlb[j]) <= MAXABSVBCOEF )
5896  {
5897  SCIP_Real vlbsol;
5898 
5899  if( (bvlb[j] >= 0.0 && SCIPisGT(scip, bvlb[j] * SCIPvarGetLbLocal(zvlb[j]) + dvlb[j], SCIPvarGetUbLocal(var))) ||
5900  (bvlb[j] <= 0.0 && SCIPisGT(scip, bvlb[j] * SCIPvarGetUbLocal(zvlb[j]) + dvlb[j], SCIPvarGetUbLocal(var))) )
5901  {
5902  *cutoff = TRUE;
5903  SCIPdebugMsg(scip, "variable bound <%s>[%g,%g] >= %g<%s>[%g,%g] + %g implies local cutoff\n",
5905  bvlb[j], SCIPvarGetName(zvlb[j]), SCIPvarGetLbLocal(zvlb[j]), SCIPvarGetUbLocal(zvlb[j]), dvlb[j]);
5906  goto TERMINATE;
5907  }
5908 
5909  assert(0 <= SCIPvarGetProbindex(zvlb[j]) && SCIPvarGetProbindex(zvlb[j]) < nbinvars);
5910  vlbsol = bvlb[j] * SCIPgetSolVal(scip, sol, zvlb[j]) + dvlb[j];
5911  if( SCIPisGE(scip, vlbsol, bestlbsol) )
5912  {
5913  bestlbsol = vlbsol;
5914  bestlbtype = j;
5915  }
5916  }
5917  }
5918 
5919  /* if no lb or vlb with binary variable was found, we have to abort */
5920  if( SCIPisInfinity(scip, -bestlbsol) )
5921  goto TERMINATE;
5922 
5923  if( bestlbtype == -1 )
5924  {
5925  rhs -= valscale * knapvals[i] * bestlbsol;
5926  SCIPdebugMsg(scip, " -> non-binary variable %+.15g<%s>(%.15g) replaced with lower bound %.15g (rhs=%.15g)\n",
5927  valscale * knapvals[i], SCIPvarGetName(var), SCIPgetSolVal(scip, sol, var), SCIPvarGetLbGlobal(var), rhs);
5928  }
5929  else
5930  {
5931  assert(0 <= SCIPvarGetProbindex(zvlb[bestlbtype]) && SCIPvarGetProbindex(zvlb[bestlbtype]) < nbinvars);
5932  rhs -= valscale * knapvals[i] * dvlb[bestlbtype];
5933  binvals[SCIPvarGetProbindex(zvlb[bestlbtype])] += valscale * knapvals[i] * bvlb[bestlbtype];
5934 
5935  if( SCIPisInfinity(scip, REALABS(binvals[SCIPvarGetProbindex(zvlb[bestlbtype])])) )
5936  goto TERMINATE;
5937 
5938  if( !noknapsackconshdlr )
5939  {
5940  assert(tmpindices != NULL);
5941 
5942  tmpindices[tmp] = SCIPvarGetProbindex(zvlb[bestlbtype]);
5943  ++tmp;
5944  }
5945  SCIPdebugMsg(scip, " -> non-binary variable %+.15g<%s>(%.15g) replaced with variable lower bound %+.15g<%s>(%.15g) %+.15g (rhs=%.15g)\n",
5946  valscale * knapvals[i], SCIPvarGetName(var), SCIPgetSolVal(scip, sol, var),
5947  bvlb[bestlbtype], SCIPvarGetName(zvlb[bestlbtype]),
5948  SCIPgetSolVal(scip, sol, zvlb[bestlbtype]), dvlb[bestlbtype], rhs);
5949  }
5950  }
5951  else
5952  {
5953  SCIP_VAR** zvub;
5954  SCIP_Real* bvub;
5955  SCIP_Real* dvub;
5956  SCIP_Real bestubsol;
5957  int bestubtype;
5958  int nvub;
5959  int j;
5960 
5961  assert(valscale * knapvals[i] < 0.0);
5962 
5963  /* a_j < 0: substitution with ub or vub */
5964  nvub = SCIPvarGetNVubs(var);
5965  zvub = SCIPvarGetVubVars(var);
5966  bvub = SCIPvarGetVubCoefs(var);
5967  dvub = SCIPvarGetVubConstants(var);
5968 
5969  /* search for ub or vub with minimal bound value */
5970  bestubsol = SCIPvarGetUbGlobal(var);
5971  bestubtype = -1;
5972  for( j = 0; j < nvub; j++ )
5973  {
5974  /* use only numerical stable vub with active binary variable z */
5975  if( SCIPvarIsBinary(zvub[j]) && SCIPvarIsActive(zvub[j]) && REALABS(bvub[j]) <= MAXABSVBCOEF )
5976  {
5977  SCIP_Real vubsol;
5978 
5979  if( (bvub[j] >= 0.0 && SCIPisLT(scip, bvub[j] * SCIPvarGetUbLocal(zvub[j]) + dvub[j], SCIPvarGetLbLocal(var))) ||
5980  (bvub[j] <= 0.0 && SCIPisLT(scip, bvub[j] * SCIPvarGetLbLocal(zvub[j]) + dvub[j], SCIPvarGetLbLocal(var))) )
5981  {
5982  *cutoff = TRUE;
5983  SCIPdebugMsg(scip, "variable bound <%s>[%g,%g] <= %g<%s>[%g,%g] + %g implies local cutoff\n",
5985  bvub[j], SCIPvarGetName(zvub[j]), SCIPvarGetLbLocal(zvub[j]), SCIPvarGetUbLocal(zvub[j]), dvub[j]);
5986  goto TERMINATE;
5987  }
5988 
5989  assert(0 <= SCIPvarGetProbindex(zvub[j]) && SCIPvarGetProbindex(zvub[j]) < nbinvars);
5990  vubsol = bvub[j] * SCIPgetSolVal(scip, sol, zvub[j]) + dvub[j];
5991  if( SCIPisLE(scip, vubsol, bestubsol) )
5992  {
5993  bestubsol = vubsol;
5994  bestubtype = j;
5995  }
5996  }
5997  }
5998 
5999  /* if no ub or vub with binary variable was found, we have to abort */
6000  if( SCIPisInfinity(scip, bestubsol) )
6001  goto TERMINATE;
6002 
6003  if( bestubtype == -1 )
6004  {
6005  rhs -= valscale * knapvals[i] * bestubsol;
6006  SCIPdebugMsg(scip, " -> non-binary variable %+.15g<%s>(%.15g) replaced with upper bound %.15g (rhs=%.15g)\n",
6007  valscale * knapvals[i], SCIPvarGetName(var), SCIPgetSolVal(scip, sol, var), SCIPvarGetUbGlobal(var), rhs);
6008  }
6009  else
6010  {
6011  assert(0 <= SCIPvarGetProbindex(zvub[bestubtype]) && SCIPvarGetProbindex(zvub[bestubtype]) < nbinvars);
6012  rhs -= valscale * knapvals[i] * dvub[bestubtype];
6013  binvals[SCIPvarGetProbindex(zvub[bestubtype])] += valscale * knapvals[i] * bvub[bestubtype];
6014 
6015  if( SCIPisInfinity(scip, REALABS(binvals[SCIPvarGetProbindex(zvub[bestubtype])])) )
6016  goto TERMINATE;
6017 
6018  if( !noknapsackconshdlr )
6019  {
6020  assert(tmpindices != NULL);
6021 
6022  tmpindices[tmp] = SCIPvarGetProbindex(zvub[bestubtype]);
6023  ++tmp;
6024  }
6025  SCIPdebugMsg(scip, " -> non-binary variable %+.15g<%s>(%.15g) replaced with variable upper bound %+.15g<%s>(%.15g) %+.15g (rhs=%.15g)\n",
6026  valscale * knapvals[i], SCIPvarGetName(var), SCIPgetSolVal(scip, sol, var),
6027  bvub[bestubtype], SCIPvarGetName(zvub[bestubtype]),
6028  SCIPgetSolVal(scip, sol, zvub[bestubtype]), dvub[bestubtype], rhs);
6029  }
6030  }
6031  }
6032 
6033  /* convert coefficients of all (now binary) variables to positive integers:
6034  * - make all coefficients integral
6035  * - make all coefficients positive (substitute negated variable)
6036  */
6037  nconsvars = 0;
6038 
6039  /* calculate scalar which makes all coefficients integral in relative allowed difference in between
6040  * -SCIPepsilon(scip) and KNAPSACKRELAX_MAXDELTA
6041  */
6043  KNAPSACKRELAX_MAXDNOM, KNAPSACKRELAX_MAXSCALE, &intscalar, &success) );
6044  SCIPdebugMsg(scip, " -> intscalar = %.15g\n", intscalar);
6045 
6046  /* if coefficients cannot be made integral, we have to use a scalar of 1.0 and only round fractional coefficients down */
6047  if( !success )
6048  intscalar = 1.0;
6049 
6050  /* make all coefficients integral and positive:
6051  * - scale a~_j = a_j * intscalar
6052  * - substitute x~_j = 1 - x_j if a~_j < 0
6053  */
6054  rhs = rhs*intscalar;
6055 
6056  SCIPdebugMsg(scip, " -> rhs = %.15g\n", rhs);
6057  minact = 0;
6058  maxact = 0;
6059  for( i = 0; i < nbinvars; i++ )
6060  {
6061&