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