# SCIP

Solving Constraint Integer Programs

sepa_zerohalf.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 /* */
7 /* fuer Informationstechnik Berlin */
8 /* */
10 /* */
12 /* along with SCIP; see the file COPYING. If not email to scip@zib.de. */
13 /* */
14 /* * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * */
15
16 /**@file sepa_zerohalf.c
17  * @brief {0,1/2}-cuts separator
18  * @author Manuel Kutschka
19  * @author Kati Wolter
20  *
21  * {0,1/2}-Chvátal-Gomory cuts separator. It solves the following separation problem:
22  * Consider an integer program
23  * \f[
24  * \min \{ c^T x : Ax \leq b, x \geq 0, x \mbox{ integer} \}
25  * \f]
26  * and a fractional solution \f$x^*\f$ of its LP relaxation. Find a weightvector \f$u\f$ whose entries \f$u_i\f$ are either 0 or
27  * \f$\frac{1}{2}\f$ such that the following inequality is valid for all integral solutions and violated by \f$x^*\f$:
28  * \f[
29  * \lfloor(u^T A) x \rfloor \leq \lfloor u^T b\rfloor
30  * \f]
31  * or (if exact methods are used) give a proof that no such inequality exists.
32  *
33  * References:
34  * - Alberto Caprara, Matteo Fischetti. {0,1/2}-Chvatal-Gomory cuts. Math. Programming, Volume 74, p221--235, 1996.
35  * - Arie M. C. A. Koster, Adrian Zymolka and Manuel Kutschka. \n
36  * Algorithms to separate {0,1/2}-Chvatal-Gomory cuts.
37  * Algorithms - ESA 2007: 15th Annual European Symposium, Eilat, Israel, October 8-10, 2007, \n
38  * Proceedings. Lecture Notes in Computer Science, Volume 4698, p. 693--704, 2007.
39  * - Arie M. C. A. Koster, Adrian Zymolka and Manuel Kutschka. \n
40  * Algorithms to separate {0,1/2}-Chvatal-Gomory cuts (Extended Version). \n
41  * ZIB Report 07-10, Zuse Institute Berlin, 2007. http://www.zib.de/Publications/Reports/ZR-07-10.pdf
42  * - Manuel Kutschka. Algorithmen zur Separierung von {0,1/2}-Schnitten. Diplomarbeit. Technische Universitaet Berlin, 2007.
43  */
44
45 /*---+----1----+----2----+----3----+----4----+----5----+----6----+----7----+----8----+----9----+----0----+----1----+----2*/
46
47 #include "string.h"
48 #include "scip/sepa_zerohalf.h"
49 #include "scip/cons_linear.h"
50 #include "scip/scipdefplugins.h"
51
52 /* prints short statistics (callback, preprocessing, adding cuts) */
53 /* #define ZEROHALF__PRINT_STATISTICS */ /**< print statistics */
54
55 #define SEPA_NAME "zerohalf"
56 #define SEPA_DESC "{0,1/2}-cuts separator"
57 #define SEPA_PRIORITY -6000
58 #define SEPA_FREQ -1
59 #define SEPA_MAXBOUNDDIST 0.0
60 #define SEPA_USESSUBSCIP TRUE
61 #define SEPA_DELAY FALSE
62
63 #define DEFAULT_MAXROUNDS 5 /**< maximal number of zerohalf separation rounds per node (-1: unlimited) */
64 #define DEFAULT_MAXROUNDSROOT 10 /**< maximal number of zerohalf separation rounds in the root node (-1: unlimited) */
65 #define DEFAULT_MAXSEPACUTS 50 /**< maximal number of {0,1/2}-cuts separated per separation round */
66 #define DEFAULT_MAXSEPACUTSROOT 500 /**< maximal number of {0,1/2}-cuts separated per separation round in root node */
67 #define DEFAULT_DYNAMICCUTS TRUE /**< should generated cuts be removed from the LP if they are no longer tight? */
68 #define DEFAULT_DECOMPOSEPROBLEM FALSE /**< should problem be decomposed into subproblems (if possible)? */
69 #define DEFAULT_MAXDEPTH -1 /**< separating cuts only if depth <= maxdepth (-1: unlimited) */
70 #define DEFAULT_MINVIOLATION 0.30 /**< minimal violation of a {0,1/2}-cut to be separated */
71 #define DEFAULT_FORCECUTSTOLP FALSE /**< should the cuts be forced to enter the LP? (bypassing SCIPefficacy criteria) */
72 #define DEFAULT_FORCECUTSTOSEPASTORE FALSE /**< should the cuts be forced to enter SCIP's sepastore?
73  * (bypassing SCIPefficicacy criteria, if no other cut is found) */
74 #define DEFAULT_MAXCUTS 100 /**< maximal number of {0,1/2}-cuts determined per separation round
75  * (this includes separated but inefficacious cuts) */
76 #define DEFAULT_MAXCUTSROOT 1000 /**< maximal number of {0,1/2}-cuts determined per separation round
77  * in the root node (this includes separated but inefficacious cuts) */
78 #define DEFAULT_SUBSCIPOBJECTIVE 'v' /**< auxiliary IP objective function type */
79 #define DEFAULT_RELAXCONTVARS FALSE /**< should continuous variables be relaxed by adding variable bounds? */
80 #define DEFAULT_SCALEFRACCOEFFS TRUE /**< should rows be scaled to make fractional coefficients integer? */
81 #define DEFAULT_SUBSCIPSETTINGS "-" /**< optional settings file of the auxiliary IP (-: none) */
82 #define DEFAULT_SUBSCIPSOLLIMIT -1 /**< limits/solutions setting of the auxiliary IP */
83 #define DEFAULT_SUBSCIPUSEALLSOLS TRUE /**< should all (proper) solutions of the auxiliary IP be used to generate
84  * cuts instead of using only the best? */
85 #define DEFAULT_PPDELTA 0.500 /**< value of delta parameter used in preprocessing method 'd' */
86 #define DEFAULT_SUBSCIPOBJPEN 0.001 /**< penalty factor used with objective function 'p' of auxiliary IP */
87
88 #define DEFAULT_PPMETHODS "CXGXIM" /**< preprocessing methods and ordering */
89 #define DEFAULT_SEPAMETHODS "2g" /**< preprocessing methods and ordering */
90 #define DEFAULT_MAXNCALLS -1LL /**< maximal number of calls (-1: unlimited) */
91 #define DEFAULT_IGNOREPREVIOUSZHCUTS FALSE /**< should zerohalf cuts found in previous callbacks be ignored? */
92 #define DEFAULT_ONLYORIGROWS FALSE /**< should only original LP rows be considered (i.e. ignore previously added LP rows)? */
93 #define DEFAULT_USEZHCUTPOOL TRUE /**< should zerohalf cuts be filtered using a cutpool */
94 #define DEFAULT_DELAYEDCUTS TRUE /**< should cuts be added to the delayed cut pool? */
95 #define DEFAULT_MAXTESTDELTA 10 /**< maximal number of different deltas to try for cmir (-1: unlimited, 0: delta=1) */
96 #define DEFAULT_TRYNEGSCALING TRUE /**< should negative values also be tested in scaling for cmir? */
98 /* cut pool management */
99 #define ORTHOFUNC 'e'
100 #define MINORTHO 0.5
101
102 /* SCIPcalcMIR parameters */
103 #define NNONZOFFSET 500
104 #define BOUNDSWITCH 0.50
105 #define USEVBDS TRUE
106 #define ALLOWLOCAL TRUE
107 #define FIXINTEGRALRHS TRUE
108 #define BOUNDSFORTRANS NULL
109 #define BOUNDTYPESFORTRANS NULL
110 #define MAXWEIGHTRANGE 10000.00
111 #define MINFRAC 0.05
112 #define MAXFRAC 1.00
114 /* SCIPcalcRowIntegralScalar parameters */
115 #define MAXDNOM 1000
116 #define MAXSCALE 1000.0
117
118 /* should variable bounds be used for substituting continuous variables */
119 #define USEVARBOUNDS TRUE
121
122 /* --------------------------------------------------------------------------------------------------------------------
123  * definitions of enums and some related strings
124  * -------------------------------------------------------------------------------------------------------------------- */
125
126 /** preprocessing methods, usable within the ppmethods parameter */
128  {
130  DELETEZEROROWS = 'Z',
132  DELETECOLSINGLETONS = 's',
139  PPCOLUMNS = 'C',
140  PPROWS = 'R'
141  };
142 #if 0 /* currently not used */
143 typedef enum preprocessingmethods PREPROCESSINGMETHODS;
144 #endif
145
146 /** separation methods, usable within the sepamethods parameter */
147 enum sepamethods
148  {
149  STOPIFCUTWASFOUND = '!',
150  SOLVEAUXSCIP = 's',
152  ENUMHEURNMAX1 = 'e',
154  GAUSSHEUR = 'g',
156  };
157 #if 0 /* currently not used */
158 typedef enum sepamethods SEPAMETHODS;
159 #endif
160
161 /** statistics: "origin" of separated cut */
162 enum cutseparatedby
163  {
165  };
167
169
170 /* --------------------------------------------------------------------------------------------------------------------
171  * auxiliary (inline) functions
172  * -------------------------------------------------------------------------------------------------------------------- */
173
174 #define ISEVEN(scip, value) (SCIPisEQ((scip) , SCIPfloor(scip , (value) / 2) , (value) / 2)) /**< is value even? */
175 #define ISODD(scip, value) (!(ISEVEN((scip), (value)))) /**< is value odd? */
176 #define XOR(bool1, bool2) ((bool1) ^ (bool2))
177 #define DIV(value1, powerof2value) (((unsigned int)(value1)) / ((unsigned int)(powerof2value))) /**< integer division using a power of 2 as divisor */
178 #define MOD(value1, powerof2value) (((unsigned int)(value1)) % ((unsigned int)(powerof2value))) /**< remainder of integer division using a power of 2
179  as divisor */
180 #ifndef BMSmoveMemoryArray
181 /** moves array at source with size num to ptr */
182 #define BMSmoveMemoryArray(ptr, source, num) \
183  { \
184  size_t size__ = (num) * sizeof(*(ptr)); \
185  if( size__ > 0 ) \
186  { \
187  assert((void*)(ptr) != NULL); \
188  assert((void*)(source) != NULL); \
189  memmove((void*)(ptr), (void*)(source), size__); \
190  } \
191  }
192 #endif
193
194 #ifdef ZEROHALF__PRINT_STATISTICS
195
196 #define ZEROHALFstatistics(x) x /**< execute if ZEROHALF__PRINT_STATISTICS is defined */
197 #define ZEROHALFstatisticsMessage printf("#### ") ; printf /**< print statistics message */
198 #define ZEROHALFcreateNewTimer(timervar) SCIP_CALL(SCIPcreateClock(scip, &timervar)) /**< create new timer */
199 #define ZEROHALFcreateTimer(timervar) SCIP_CALL(SCIPcreateClock(scip, &timervar)) /**< recreate existing timer */
200 #define ZEROHALFfreeTimer(timervar) SCIP_CALL(SCIPfreeClock(scip, &timervar)) /**< free timer */
201 #define ZEROHALFresetTimer(timervar) SCIP_CALL(SCIPresetClock(scip, timervar)) /**< reset timer */
202 #define ZEROHALFstartTimer(timervar) SCIP_CALL(SCIPstartClock(scip, timervar)) /**< start timer */
203 #define ZEROHALFstopTimer(timervar) SCIP_CALL(SCIPstopClock(scip, timervar)) /**< stop timer */
204 #define ZEROHALFevalTimer(timervar) (SCIPgetClockTime(scip, timervar)) /**< evaluate timer (get time) */
205
206 #else
207
208 #if 0 /* currently not used */
209 #define ZEROHALFstatistics(x) /**/ /**< nothing */
210 #define ZEROHALFstatisticsMessage while( FALSE ) printf /**< nothing */
211 #define ZEROHALFcreateNewTimer(timervar) /**/ /**< nothing */
212 #define ZEROHALFcreateTimer(timervar) /**/ /**< nothing */
213 #define ZEROHALFfreeTimer(timervar) /**/ /**< nothing */
214 #define ZEROHALFresetTimer(timervar) /**/ /**< nothing */
215 #define ZEROHALFstartTimer(timervar) /**/ /**< nothing */
216 #define ZEROHALFstopTimer(timervar) /**/ /**< nothing */
217 #define ZEROHALFevalTimer(timervar) (0.0) /**< nothing */
218 #endif
219
220 #endif
221
222
223 /* --------------------------------------------------------------------------------------------------------------------
224  * symbolic constants for statistical analysis
225  * -------------------------------------------------------------------------------------------------------------------- */
226
227 /* row / columns */
228 #define IRRELEVANT -1 /**< row or column is irrelevant */
229
230 /* row */
231 #define ZERO_ROW -100 /**< row has no nonzero entries */
232 #define IDENT_TO_ROW_WITH_SMALLER_SLACK -101 /**< row is identical to another row but has a larger slack value */
233 #define SLACK_GREATER_THAN_MAXSLACK -102 /**< row has a slack value > maxslack */
234 #define DEFINES_VIOLATED_ZEROHALF_CUT -103 /**< row defines a violated zerohalf cut */
235 #ifdef WITHDECOMPOSE
236 #define ROW_IN_SUBPROB_WITHOUT_ODD_RHS -104 /**< row is part of a subproblem without rows with odd rhs */
237 #endif
238 #define NONEXISTENT_ROW -105 /**< row does not exist (lhs is -infinity or rhs is infinity) */
239 #define NO_RELEVANT_COLUMNS -106 /**< row does not contain relevant columns */
240 #define SLACK_GREATER_THAN_MSL_MINUS_SODD -107 /**< row has even rhs and the sum of its slack value and the minimum
241  slack value of a odd-rhs-row exceeds maxslack */
242 #define LARGE_COL_EXISTS -108 /**< row contains a column which rounding penalty exceeds maxslack
243  and the sum of this row's slack and the minimum slack of another
244  row with the proper columns exceeds maxslack as well */
246 /* column */
247 #define ZERO_COLUMN -200 /**< column has no nonzero entries */
248 #define IDENT_TO_ANOTHER_COLUMN -201 /**< column is identical to another column */
249 #define ZERO_LP_SOL -202 /**< column corresponds to a variable whose LP solution is zero */
250 #define LP_SOL_EQUALS_EVEN_LB -203 /**< column corresponds to a variable whose LP solution equals its even lb */
251 #define LP_SOL_EQUALS_ODD_LB -204 /**< column corresponds to a variable whose LP solution equals its odd lb */
252 #define LP_SOL_EQUALS_EVEN_UB -205 /**< column corresponds to a variable whose LP solution equals its even ub */
253 #define LP_SOL_EQUALS_ODD_UB -206 /**< column corresponds to a variable whose LP solution equals its odd ub */
254 #define SINGLETON_COLUMN -207 /**< column has only one nonzero entry */
255 #define CONTINUOUS_VARIABLE -208 /**< column corresponds to a non-integer variable */
256 #define SMALL_FRACSOL_HEUR -209 /**< column has been omitted (see preprocessColumnsWithSmallFracsol) */
257 #define ALL_MATRIX_ROWS_DELETED -210 /**< all rows (of the current subproblem) have been deleted */
258 #ifdef WITHDECOMPOSE
259 #define COLUMN_IN_SUBPROB_WITHOUT_ODD_RHS -211 /**< column is part of a subproblem without a row with odd rhs value */
260 #endif
263 /* --------------------------------------------------------------------------------------------------------------------
264  * bit array data structures and functions
265  * -------------------------------------------------------------------------------------------------------------------- */
267
268 #define BITARRAYBASETYPE unsigned int /**< base type used for the bitarray data structures */
270
271 /** size of BITARRAYBASETYPE */
272 static const unsigned int Zerohalf_bitarraybasetypesize = sizeof(BITARRAYBASETYPE);
273
274 /** number of bits per BITARRAYBASETYPE */
275 static const unsigned int Zerohalf_bitarraybasetypesize_nbits = sizeof(BITARRAYBASETYPE) << 3;
276
278 #define BITARRAY BITARRAYBASETYPE*
279
280 /** get the bit mask where the pos-th bit is set */
281 #define BITMASK(pos) ((unsigned int)(1u << (pos)))
282
283 /** set the pos-th bit of var */
284 #define BITSET(var, pos) (var) |= BITMASK(pos)
285
286 /** is the pos-th bit of var set? */
287 #define BITISSET(var, pos) (var & BITMASK(pos))
288
289 /** set the pos-th bit of bitarray barray */
290 #define BITARRAYBITSET(barray, pos) BITSET(barray[DIV((pos),Zerohalf_bitarraybasetypesize_nbits)], \
291  MOD(pos,Zerohalf_bitarraybasetypesize_nbits))
292
293 /** is the pos-th bit of bitarray barray set? */
294 #define BITARRAYBITISSET(barray, pos) BITISSET(barray[DIV(pos,Zerohalf_bitarraybasetypesize_nbits)], \
295  MOD(pos,Zerohalf_bitarraybasetypesize_nbits))
297 /** clear bitarray */
298 #define BITARRAYCLEAR(barray, barraysize) BMSclearMemoryArray(barray,barraysize)
300 /** calculates the number of array elements (w.r.t. the bitarray base type) required to create the bitarray */
301 #define GETREQUIREDBITARRAYSIZE(nvalstostore) \
302  ((((unsigned int)(nvalstostore)) % (Zerohalf_bitarraybasetypesize_nbits) == 0) \
303  ? (((unsigned int)(nvalstostore)) / (Zerohalf_bitarraybasetypesize_nbits)) \
304  : ((((unsigned int)(nvalstostore)) / (Zerohalf_bitarraybasetypesize_nbits)) + 1))
305
306 /** get the corresponding array element of a bitarray position */
307 #define GETBITARRAYINDEX(pos) DIV((pos),Zerohalf_bitarraybasetypesize_nbits)
308
309 /** get the bitmask to mask all bits except the pos-th bit of an array element */
311
312 /** apply operation op for all array elements of bitarray barray1 and barray2 */
313 #define BITARRAYSFOREACH(barray1, barray2, size, op) \
314  { \
315  int idx__; \
316  for( idx__ = 0 ; idx__ < (size) ; ++idx__) \
317  { \
318  barray2[idx__] op barray1[idx__]; \
319  } \
320  }
321
322 /** barray2 = barray1 XOR barray2 */
323 #define BITARRAYSXOR(barray1, barray2, size) BITARRAYSFOREACH(barray1,barray2,size,^=)
324
325 /** are barray1 and barray2 equal? */
326 #define BITARRAYSAREEQUAL(barray1, barray2, size) \
327  (memcmp((void*)(barray1), (void*)(barray2), (size_t)((size) * (Zerohalf_bitarraybasetypesize))) == 0)
328
329 #if 0 /* currently not used */
330 /** clear the pos-th bit of var */
331 #define BITCLEAR(var, pos) (var) &= ~BITMASK(pos)
333 /** flip the pos-th bit of var */
334 #define BITFLIP(var, pos) (var) ^= BITMASK(pos)
336 /** clear the pos-th bit of bitarray barray */
337 #define BITARRAYBITCLEAR(barray, pos) BITCLEAR(barray[DIV((pos),Zerohalf_bitarraybasetypesize_nbits)], \
338  MOD(pos,Zerohalf_bitarraybasetypesize_nbits))
339
340 /** flip the pos-th bit of bitarray barray */
341 #define BITARRAYBITFLIP(barray, pos) BITFLIP(barray[DIV((pos),Zerohalf_bitarraybasetypesize_nbits)], \
342  MOD(pos,Zerohalf_bitarraybasetypesize_nbits))
343
344 /** barray2 = barray1 AND barray2 */
345 #define BITARRAYSAND(barray1, barray2, size) BITARRAYSFOREACH(barray1,barray2,size,&=)
346
347 /** barray2 = barray1 OR barray2 */
348 #define BITARRAYSOR(barray1, barray2, size) BITARRAYSFOREACH(barray1,barray2,size,|=)
349
350 /** barray2 = NOT barray1 */
351 #define BITARRAYSNOT(barray1, barray2, size) BITARRAYSFOREACH(barray1,barray2,size,= ~)
352 #endif
353
354
355 /* --------------------------------------------------------------------------------------------------------------------
356  * data structures
357  * -------------------------------------------------------------------------------------------------------------------- */
358
359
360 /** parameters */
362 {
363  int maxrounds; /**< maximal number of {0,1/2} separation rounds per node (-1: unlimited) */
364  int maxroundsroot; /**< maximal number of {0,1/2} separation rounds in the root node (-1: unlimited) */
365  int maxsepacuts; /**< maximal number of {0,1/2} cuts separated per separation round */
366  int maxsepacutsroot; /**< maximal number of {0,1/2} cuts separated per separation round in root node */
367  int maxdepth; /**< separating cuts only if depth <= maxdepth (-1: unlimited) */
368  SCIP_Bool dynamiccuts; /**< should generated cuts be removed from the LP if they are no longer tight? */
369  SCIP_Bool decomposeproblem; /**< should problem be decomposed into subproblems (if possible)? */
370
371  SCIP_Real minviolation; /**< minimal violation of a {0,1/2}-cut to be separated */
372  SCIP_Bool forcecutstolp; /**< should the cuts be forced to enter the LP? */
373  SCIP_Bool forcecutstosepastore; /**< should the cuts be forced to enter SCIP's sepastore? */
374  int maxcuts; /**< maximal number of {0,1/2}-cuts determined per separation round
375  (this includes separated but inefficacious cuts) */
376  int maxcutsroot; /**< maximal number of {0,1/2}-cuts determined per separation round
377  in the root node (this includes separated but inefficacious cuts) */
378
379  char* ppmethods; /**< preprocessing methods */
380  char* sepamethods; /**< separation methods */
381  int nppmethods; /**< length of ppmethods string */
382  int nsepamethods; /**< length of sepamethods string */
383
384  int subscipsollimit; /**< value of auxiliary IP / subscip "limits/sol" */
385  SCIP_Bool subscipuseallsols; /**< should all known feasible solution of the auxiliary IP be considered? */
386  SCIP_Bool relaxcontvars; /**< should continuous vars be relaxed by adding varbounds? */
387  SCIP_Bool scalefraccoeffs; /**< should fractional coeffs be scaled to become integer? */
388  char* subscipsettings; /**< optional settings file for the auxiliary IP / subscip */
389  char subscipobjective; /**< type of objective function of the auxiliary IP / subscip */
390  SCIP_Real ppdelta; /**< value of delta parameter used in preprocessing method 'd' */
391  SCIP_Real subscipobjpen; /**< penalty factor used with objective function 'p' of auxiliary IP */
392  SCIP_Longint maxncalls; /**< maximal number of callbacks */
393  SCIP_Bool ignoreprevzhcuts; /**< should zerohalf cuts found within previous callbacks considered as well? */
394  SCIP_Bool onlyorigrows; /**< should only original LP rows be considered (i.e. ignore previously added LP rows)? */
395  SCIP_Bool usezhcutpool; /**< should zerohalf cuts be filtered using a cutpool? */
396  SCIP_Bool delayedcuts; /**< should cuts be added to the delayed cut pool? */
397
398  SCIP_Real maxslack; /**< initial: 1.0 - 2.0 * minviolation */
399  int norigrows; /**< number of original LP rows */
400  int* origrows; /**< set of SCIP_ROW->index of all original LP rows */
401
402  int maxnnonz; /**< maximal number of nonzeros allowed in a zerohalf cut */
403  int maxtestdelta; /**< maximal number of different deltas to try for cmir (-1: unlimited, 0: delta=1) */
404  SCIP_Bool trynegscaling; /**< should negative values also be tested in scaling for cmir? */
405
406  /* statistics */
407  int totalncutsfound; /**< total number of separated zerohalf cuts, including inefficious ones */
408  int totalnsepacuts; /**< total number of separated zerohalf cuts */
409  SCIP_CLOCK** pptimers; /**< timers of preprocessing methods */
410  SCIP_CLOCK** sepatimers; /**< timers of separation algorithms */
411  SCIP_CLOCK* dtimer; /**< timer of decomposition method */
412  int* nsepacutsalgo; /**< number zerohalf cuts separated by a specific separation algorithm,
413  * including inefficious cuts */
414  int* nzerohalfcutsalgo; /**< number zerohalf cuts separated by a specific separation algorithm */
415 };
416
417
418 /** sub data of the LP or a sub-LP obtained by problem decomposition */
419 struct Zerohalf_SubLPData
420 {
421  int* rrows; /**< relevant rows (indices of elements in rows array of ZEROHALF_LPDATA) */
422  SCIP_Real* rrowsrhs; /**< rhs value of relevant rows; could also be lhs value of the orig row */
423  SCIP_Real* rrowsslack; /**< slack value of relevant rows */
424  int nrrows; /**< number of relevant rows */
425  int rrowssize; /**< size of rrows, rrowsrhs, rrowsslack arrays */
426
427  int* rcols; /**< relevant columns (indices of elements in cols array of ZEROHALF_LPDATA) */
428  SCIP_Real* rcolslbslack; /**< slack value of lower bound constraint: x*_j - lb_j */
429  SCIP_Real* rcolsubslack; /**< slack value of upper bound constraint: ub_j - x*_j */
430  int nrcols; /**< number of relevant columns */
431  int rcolssize; /**< size of rcols, rcolslbslack, rcolsubslack arrays */
432 };
433 typedef struct Zerohalf_SubLPData ZEROHALF_SUBLPDATA;
434
435
436 /** LP data */
437 struct Zerohalf_LPData
438 {
439  SCIP_VAR** vars; /**< LP variables */
440  SCIP_ROW** rows; /**< LP rows */
441  SCIP_COL** cols; /**< LP columns */
442  int nvars; /**< number of LP variables */
443  int nrows; /**< number of LP rows */
444  int ncols; /**< number of LP columns */
445  int nvarbounds; /**< number of variable bounds (-x_j <= -lb_j, x_j <= ub_j) */
446
447  ZEROHALF_SUBLPDATA** subproblems; /**< decomposed subproblems (subset of the variables, rows and columns above) */
448  int nsubproblems; /**< number of subproblems */
449
450  SCIP_Real* intscalarsleftrow; /**< array of scalars that would make left half-rows (-a^Tx <= -lhs) rows integral (0.0 if scalar has not been calculated) */
451  SCIP_Real* intscalarsrightrow; /**< array of scalars that would make right half-rows (a^Tx <= rhs) rows integral (0.0 if scalar has not been calculated) */
452
453  /* row related index sets */
454  int* subproblemsindexofrow; /**< is rows index relevant? value <0: not relevant,
455  value >=0: index of subproblem containing the row */
456  int* rrowsindexofleftrow; /**< maps rows index of lhs <= a^Tx <= rhs to rrows index of -a^Tx <= -lhs */
457  int* rrowsindexofrightrow; /**< maps rows index of lhs <= a^Tx <= rhs to rrows index of a^Tx <= rhs */
458
459  /* col related index sets */
460  int* subproblemsindexofcol; /**< is cols index relevant? value <0: not relevant
461  * value >=0: index of subproblem containing the column */
462  int* rcolsindexofcol; /**< maps cols index to rcols index */
463
464  int* bestlbidxofcol; /**< maps cols index of a continuous variable to the index of its
465  * best lower bound (-2: undetermined, -1: lb, >=0: index of vlb)*/
466  int* bestubidxofcol; /**< maps cols index of a continuous variable to the index of its
467  * best upper bound (-2: undetermined, -1: ub, >=0: index of vub)*/
468
469  /* statistics */
470  int ndelvarbounds; /**< number of deleted variable bounds by basic preprocessing */
471
472 };
473 typedef struct Zerohalf_LPData ZEROHALF_LPDATA;
474
475
476 /** data structure to store data of the auxiliary IP:
477  * (1) minimize violation:
478  * min z := s^T v + x^T y
479  * s.t. (b (mod 2))^T v - 2q = 1 "odd rhs"
480  * (A (mod 2))^T v - y - -2r = 0 "column sum"
481  *
482  * v \\in {0,1}^m
483  * y \\in {0,1}^n
484  * r \\in Z^n_+
485  * q \\in Z_+
486  * (2) minimize (weighted) number of aggregated rows
487  * min z := w^T v
488  * s.t. s^T v + x^T y < 1 "feasibility"
489  * (b (mod 2))^T v - 2q = 1 "odd rhs"
490  * (A (mod 2))^T v - y - -2r = 0 "column sum"
491  *
492  * v \\in {0,1}^m
493  * y \\in {0,1}^n
494  * r \\in Z^n_+
495  * q \\in Z_+
496  *
497  * with w \\in \\Rset^m
498  */
499 struct Zerohalf_AuxIPData
500 {
501  SCIP* subscip; /**< pointer to (sub)SCIP data structure containing the auxiliary IP */
502  int m; /**< number of rows */
503  int n; /**< number of cols */
504
505  SCIP_VAR** v; /**< decision variable: 1 iff row is selected for generating a violated zerohalf cut */
506  SCIP_VAR** y; /**< auxiliary variable used for calculating the rounding down penalties in the objective */
507  SCIP_VAR** r; /**< auxiliary variable used for modelling mod 2 calculus */
508  SCIP_VAR* q; /**< auxiliary variable used for modelling mod 2 calculus */
509
510  SCIP_CONS* feasipcons; /**< feasibility constraint */
511  SCIP_CONS* oddrhscons; /**< odd rhs constraint */
512  SCIP_CONS** columnsumcons; /**< column sum constraints */
513
514  SCIP_Real timelimit; /**< value of "limits/time" of subscip */
515  SCIP_Real memorylimit; /**< value of "limits/memory" of subscip */
516  SCIP_Real objectivelimit; /**< value of objective limit of subscip */
517  int nodelimit; /**< value of "limits/nodes" of subscip */
518 };
519 typedef struct Zerohalf_AuxIPData ZEROHALF_AUXIPDATA;
520
521
522 /** data structure to store data (mod 2) densely*/
523 struct Zerohalf_Mod2Data
524 {
525  ZEROHALF_SUBLPDATA* relatedsubproblem; /**< pointer to corresponding subproblem data structure */
526
527  BITARRAY* rows; /**< dense mod 2 rows */
528  BITARRAY* rowaggregations; /**< Zerohalf_Mod2Data->rows index set, storing the actual row aggregations */
529  SCIP_Bool* rhs; /**< TRUE iff corresponding relatedsubproblem->rrowsrhs is odd */
530
531  SCIP_Real* slacks; /**< slack value (equal to corresponding relatedsubproblem->rrowsslack value) */
532  SCIP_Real* fracsol; /**< LP solution value of variable of SCIP_COL* */
533
534  int nrows; /**< number of Zerohalf_Mod2Data->rows */
535  int nrcols; /**< number of relevant columns */
536
537  int rowsbitarraysize; /**< size (w.r.t. bitarray base type) of Zerohalf_Mod2Data->rows array */
538  int rowaggregationsbitarraysize; /** size (w.r.t. bitarray base type) of Zerohalf_Mod2Data->rowaggregations array */
539  int nvarbounds; /**< number of variable bounds that are part of Zerohalf_Mod2Data->rows */
540
541  int* rowsind; /**< index set of subset of Zerohalf_Mod2Data->rows */
542  int* colsind; /**< index set of subset of relatedsubproblem->rcols */
543
544  int nrowsind; /**< number of rowsind elements */
545  int ncolsind; /**< number of colsind elements */
546
547  /* statistics */
548  int* rowstatistics; /**< stores if and why a row was removed in preprocessing */
549  int* colstatistics; /**< stores if and why a column was removed in preprocessing */
550 };
551 typedef struct Zerohalf_Mod2Data ZEROHALF_MOD2DATA;
552
553
554 /** data structure to store a violated zerohalf cut and related data */
555 struct Zerohalf_CutData
556 {
557  ZEROHALF_SUBLPDATA* relatedsubproblem; /**< pointer to corresponding subproblem data structure */
558  ZEROHALF_MOD2DATA* relatedmod2data; /**< pointer to corresponding mod 2 data structure */
559
560  SCIP_ROW* cut; /**< pointer to SCIP_ROW storing the cut */
561
562  SCIP_Bool success; /**< was SCIPcalcMIR successful? */
563  SCIP_Bool isfeasviolated; /**< is zerohalf cut violated w.r.t. the feasibility tolerance? */
564  SCIP_Bool islocal; /**< is zerohalf cut only locally valid? */
565
566  SCIP_Real activity; /**< activity of the zerohalf cut */
567  SCIP_Real rhs; /**< rhs of the zerohalf cut */
568  SCIP_Real norm; /**< norm of the nonzero elements of the zerohalf cut */
569  SCIP_Real efficacy; /**< efficacy of the zerohalf cut */
570  SCIP_Real violation; /**< violation of the zerohalf cut */
571  int cutrank; /**< rank of cut */
572  int nnonz; /**< number of nonzero coefficients of the zerohalf cut */
573
574  /* statistics */
575  int nrowsincut; /**< number of LP rows combined into the zerohalf cut */
576  int nrrowsincut; /**< number of preprocessed/aggregated LP rows combined into the zerohalf cut */
577  CUTSEPARATEDBY separatedby; /**< flag to store the method that has separated the zerohalf cut */
579
580 };
581 typedef struct Zerohalf_CutData ZEROHALF_CUTDATA;
582
583
584
585 /** auxiliary graph node data structure */
586 struct Zerohalf_AuxGraph_Node;
587 typedef struct Zerohalf_AuxGraph_Node ZEROHALF_AUXGRAPH_NODE;
588
589 struct Zerohalf_AuxGraph_Node
590 {
591  ZEROHALF_AUXGRAPH_NODE** neighbors; /**< node adjacency list */
592  SCIP_Real* edgeweights; /**< weights of outgoing edges */
593  int* relatedrows; /**< label mapping outgoing edges to mod 2 rows */
594  int nneighbors; /**< number of adjacent nodes */
595  SCIP_Real distance; /**< actual distant from start node (used by Dijkstra)*/
596  ZEROHALF_AUXGRAPH_NODE* previous; /**< previous node in shortest-path-tree (used by Dijkstra) */
597 };
598
599
600 /** auxiliary graph data structure */
601 struct Zerohalf_AuxGraph
602 {
603  ZEROHALF_AUXGRAPH_NODE** nodes; /**< list of all original nodes */
604  ZEROHALF_AUXGRAPH_NODE** nodecopies; /**< list of all copies of original nodes */
605  int nnodes; /**< number of original nodes (equals number of copies) */
606 };
607 typedef struct Zerohalf_AuxGraph ZEROHALF_AUXGRAPH;
608
609
610 /* --------------------------------------------------------------------------------------------------------------------
611  * local methods: create / free data structure
612  * -------------------------------------------------------------------------------------------------------------------- */
613
614
615 /** creates and initializes sub LP data structures */
616 static
618  SCIP* scip, /**< SCIP data structure */
619  ZEROHALF_SUBLPDATA** subproblem /**< pointer to store pointer to created data structure */
620  )
621 {
622  assert(scip != NULL);
623  assert(subproblem != NULL);
624
625  SCIP_CALL( SCIPallocBlockMemory(scip, subproblem) );
626  (*subproblem)->rrows = NULL;
627  (*subproblem)->rrowsrhs = NULL;
628  (*subproblem)->rrowsslack = NULL;
629  (*subproblem)->nrrows = 0;
630  (*subproblem)->rrowssize = 0;
631
632  (*subproblem)->rcols = NULL;
633  (*subproblem)->rcolslbslack = NULL;
634  (*subproblem)->rcolsubslack = NULL;
635  (*subproblem)->nrcols = 0;
636  (*subproblem)->rcolssize = 0;
637
638  return SCIP_OKAY;
639 }
640
641
642 /** frees sub LP data structures */
643 static
645  SCIP* scip, /**< SCIP data structure */
646  ZEROHALF_SUBLPDATA** subproblem /**< pointer to pointer of data structure */
647  )
648 {
649  assert(scip != NULL);
650  assert(subproblem != NULL);
651  assert(*subproblem != NULL);
652
653  SCIPfreeBlockMemoryArrayNull(scip, &((*subproblem)->rrows), (*subproblem)->rrowssize);
654  SCIPfreeBlockMemoryArrayNull(scip, &((*subproblem)->rrowsrhs), (*subproblem)->rrowssize);
655  SCIPfreeBlockMemoryArrayNull(scip, &((*subproblem)->rrowsslack), (*subproblem)->rrowssize);
656
657  SCIPfreeBlockMemoryArrayNull(scip, &((*subproblem)->rcols), (*subproblem)->rcolssize);
658  SCIPfreeBlockMemoryArrayNull(scip, &((*subproblem)->rcolslbslack), (*subproblem)->rcolssize);
659  SCIPfreeBlockMemoryArrayNull(scip, &((*subproblem)->rcolsubslack), (*subproblem)->rcolssize);
660
661  SCIPfreeBlockMemory(scip, subproblem);
662 }
663
664
665 /** creates and initializes LP data structures */
666 static
668  SCIP* scip, /**< SCIP data structure */
669  ZEROHALF_LPDATA** lpdata /**< pointer to store pointer to created data structure */
670  )
671 {
672  assert(scip != NULL);
673  assert(lpdata != NULL);
674
675  SCIP_CALL(SCIPallocBlockMemory(scip, lpdata));
676  (*lpdata)->vars = NULL;
677  (*lpdata)->rows = NULL;
678  (*lpdata)->cols = NULL;
679
680  (*lpdata)->subproblems = NULL;
681  (*lpdata)->nsubproblems = 0;
682
683  (*lpdata)->intscalarsleftrow = NULL;
684  (*lpdata)->intscalarsrightrow = NULL;
685
686  (*lpdata)->subproblemsindexofrow = NULL;
687  (*lpdata)->rrowsindexofleftrow = NULL;
688  (*lpdata)->rrowsindexofrightrow = NULL;
689
690  (*lpdata)->subproblemsindexofcol = NULL;
691  (*lpdata)->rcolsindexofcol = NULL;
692
693  (*lpdata)->bestlbidxofcol = NULL;
694  (*lpdata)->bestubidxofcol = NULL;
695
696  return SCIP_OKAY;
697 }
698
699
700 /** frees LP data structures */
701 static
703  SCIP* scip, /**< SCIP data structure */
704  ZEROHALF_LPDATA** lpdata /**< pointer to pointer of data structure */
705  )
706 {
707  int sp;
708
709  assert(scip != NULL);
710  assert(lpdata != NULL);
711  assert(*lpdata != NULL);
712
713  /* ! Do not free (*lpdata)->vars, (*lpdata)->rows and (*lpdata)->cols ! */
714  assert(((*lpdata)->nsubproblems == 0 && (*lpdata)->subproblems == NULL)
715  || ((*lpdata)->nsubproblems > 0 && (*lpdata)->subproblems != NULL));
716
717  if( (*lpdata)->subproblems != NULL )
718  {
719  for( sp = 0; sp < (*lpdata)->nsubproblems; ++sp )
720  {
721  if( (*lpdata)->subproblems[sp] != NULL )
722  ZerohalfSubLPDataFree(scip, &((*lpdata)->subproblems[sp]));
723  }
724  SCIPfreeBlockMemoryArray(scip, &((*lpdata)->subproblems), (*lpdata)->nsubproblems);
725  }
726
727  SCIPfreeBlockMemoryArrayNull(scip, &((*lpdata)->intscalarsleftrow), (*lpdata)->nrows);
728  SCIPfreeBlockMemoryArrayNull(scip, &((*lpdata)->intscalarsrightrow), (*lpdata)->nrows);
729  SCIPfreeBlockMemoryArrayNull(scip, &((*lpdata)->subproblemsindexofrow), (*lpdata)->nrows);
730  SCIPfreeBlockMemoryArrayNull(scip, &((*lpdata)->rrowsindexofleftrow), (*lpdata)->nrows);
731  SCIPfreeBlockMemoryArrayNull(scip, &((*lpdata)->rrowsindexofrightrow), (*lpdata)->nrows);
732  SCIPfreeBlockMemoryArrayNull(scip, &((*lpdata)->subproblemsindexofcol), (*lpdata)->ncols);
733  SCIPfreeBlockMemoryArrayNull(scip, &((*lpdata)->rcolsindexofcol), (*lpdata)->ncols);
734  SCIPfreeBlockMemoryArrayNull(scip, &((*lpdata)->bestlbidxofcol), (*lpdata)->ncols);
735  SCIPfreeBlockMemoryArrayNull(scip, &((*lpdata)->bestubidxofcol), (*lpdata)->ncols);
736  SCIPfreeBlockMemory(scip, lpdata);
737
738  return SCIP_OKAY;
739 }
740
741
742 /** creates and initializes mod 2 data structures */
743 static
745  SCIP* scip, /**< SCIP data structure */
746  ZEROHALF_MOD2DATA** mod2data /**< pointer to store pointer to created data structure */
747  )
748 {
749  assert(scip != NULL);
750  assert(mod2data != NULL);
751
752  SCIP_CALL(SCIPallocBlockMemory(scip, mod2data));
754  (*mod2data)->relatedsubproblem = NULL;
755
756  (*mod2data)->rows = NULL;
757  (*mod2data)->rowaggregations = NULL;
758  (*mod2data)->rhs = NULL;
759
760  (*mod2data)->slacks = NULL;
761  (*mod2data)->fracsol = NULL;
762
763  (*mod2data)->nrows = 0;
764  (*mod2data)->nrcols = 0;
765
766  (*mod2data)->rowstatistics = NULL;
767  (*mod2data)->colstatistics = NULL;
768
769  (*mod2data)->rowsind = NULL;
770  (*mod2data)->colsind = NULL;
771
772  return SCIP_OKAY;
773 }
774
775
776 /** frees data structures */
777 static
779  SCIP* scip, /**< SCIP data structure */
780  ZEROHALF_MOD2DATA** mod2data /**< pointer to pointer of data structure */
781  )
782 {
783  int i;
784
785  assert(scip != NULL);
786  assert(mod2data != NULL);
787  assert(*mod2data != NULL);
788
789  if( (*mod2data)->rows != NULL )
790  {
791  for( i = 0 ; i < (*mod2data)->nrows ; ++i)
792  {
793  if( (*mod2data)->rows[i] != NULL )
794  {
795  SCIPfreeBlockMemoryArray(scip, &((*mod2data)->rows[i]), (*mod2data)->rowsbitarraysize);
796  }
797  }
798  SCIPfreeBlockMemoryArray(scip, &((*mod2data)->rows), (*mod2data)->nrows);
799  }
800
801  if( (*mod2data)->rowaggregations != NULL )
802  {
803  for( i = 0 ; i < (*mod2data)->nrows ; ++i)
804  {
805  if( (*mod2data)->rowaggregations[i] != NULL )
806  {
807  SCIPfreeBlockMemoryArray(scip, &((*mod2data)->rowaggregations[i]), (*mod2data)->rowaggregationsbitarraysize);
808  }
809  }
810  SCIPfreeBlockMemoryArray(scip, &((*mod2data)->rowaggregations), (*mod2data)->nrows);
811  }
812
813  SCIPfreeBlockMemoryArrayNull(scip, &((*mod2data)->rhs), (*mod2data)->nrows);
814  SCIPfreeBlockMemoryArrayNull(scip, &((*mod2data)->slacks), (*mod2data)->nrows);
815  SCIPfreeBlockMemoryArrayNull(scip, &((*mod2data)->fracsol), (*mod2data)->nrcols);
816  SCIPfreeBlockMemoryArrayNull(scip, &((*mod2data)->rowstatistics), (*mod2data)->nrows);
817  SCIPfreeBlockMemoryArrayNull(scip, &((*mod2data)->colstatistics), (*mod2data)->nrcols);
818  SCIPfreeBlockMemoryArrayNull(scip, &((*mod2data)->rowsind), (*mod2data)->nrows);
819  SCIPfreeBlockMemoryArray(scip, &((*mod2data)->colsind), (*mod2data)->nrcols);
820  SCIPfreeBlockMemory(scip, mod2data);
821
822  return SCIP_OKAY;
823 }
824
825
826 /** creates and initializes auxiliary IP data structures */
827 static
829  SCIP* scip, /**< SCIP data structure */
830  ZEROHALF_AUXIPDATA** auxipdata /**< pointer to store pointer to created data structure */
831  )
832 {
833  assert(scip != NULL);
834  assert(auxipdata != NULL);
835
836  SCIP_CALL(SCIPallocBlockMemory(scip, auxipdata));
838  (*auxipdata)->subscip = NULL;
839
840  (*auxipdata)->v = NULL;
841  (*auxipdata)->y = NULL;
842  (*auxipdata)->r = NULL;
843  (*auxipdata)->q = NULL;
844
845  (*auxipdata)->feasipcons = NULL;
846  (*auxipdata)->oddrhscons = NULL;
847  (*auxipdata)->columnsumcons = NULL;
848
849  return SCIP_OKAY;
850 }
851
852
853 /** frees auxiliary IP data structures */
854 static
856  SCIP* scip, /**< SCIP data structure */
857  ZEROHALF_AUXIPDATA** auxipdata /**< pointer to pointer of data structure */
858  )
859 {
860  assert(scip != NULL);
861  assert(auxipdata != NULL);
862  assert(*auxipdata != NULL);
863
864  if( (*auxipdata)->subscip != NULL )
865  {
866  SCIP_CALL( SCIPfree(&((*auxipdata)->subscip)) );
867  }
868
869  SCIPfreeBlockMemoryArrayNull(scip, &((*auxipdata)->columnsumcons), (*auxipdata)->n);
870  SCIPfreeBlockMemoryArrayNull(scip, &((*auxipdata)->r), (*auxipdata)->n);
871  SCIPfreeBlockMemoryArrayNull(scip, &((*auxipdata)->y), (*auxipdata)->n);
872  SCIPfreeBlockMemoryArrayNull(scip, &((*auxipdata)->v), (*auxipdata)->m);
873
874  SCIPfreeBlockMemory(scip, auxipdata);
875
876  return SCIP_OKAY;
877 }
878
879
880 /** creates and initializes cut data structures */
881 static
883  SCIP* scip, /**< SCIP data structure */
884  ZEROHALF_CUTDATA** cutdata, /**< pointer to pointer of data structure */
885  ZEROHALF_SUBLPDATA* relatedsubproblem, /**< pointer to corresponding subproblem */
886  ZEROHALF_MOD2DATA* relatedmod2data, /**< pointer to corresponding mod 2 data */
887  int nrrowsincut, /**< number of preprocessed / mod 2 rows in cut */
888  int nrowsincut, /**< number of original / LP rows in cut */
889  CUTSEPARATEDBY separatedby /**< flag storing the method that separated this cut */
890  )
891 {
892  assert(scip != NULL);
893  assert(cutdata != NULL);
894  assert(nrrowsincut >= 0);
895  assert(nrowsincut >= 0);
896
897  SCIP_CALL(SCIPallocBlockMemory(scip, cutdata));
898
899  (*cutdata)->relatedsubproblem = relatedsubproblem;
900  (*cutdata)->relatedmod2data = relatedmod2data;
901
902  (*cutdata)->cut = NULL;
903
904  (*cutdata)->success = FALSE;
905  (*cutdata)->isfeasviolated = FALSE;
906  (*cutdata)->islocal = TRUE;
907
908  (*cutdata)->activity = 0.0;
909  (*cutdata)->rhs = 0.0;
910  (*cutdata)->norm = 0.0;
911  (*cutdata)->efficacy = 0.0;
912  (*cutdata)->violation = 0.0;
913  (*cutdata)->cutrank = 0;
914  (*cutdata)->nnonz = 0;
915
916  (*cutdata)->nrrowsincut = nrrowsincut;
917  (*cutdata)->nrowsincut = nrowsincut;
918  (*cutdata)->separatedby = separatedby;
920
921  return SCIP_OKAY;
922 }
923
924
925 /** frees cut data structures */
926 static
928  SCIP* scip, /**< SCIP data structure */
929  ZEROHALF_CUTDATA** cutdata /**< pointer to pointer of data structure */
930  )
931 {
932  assert(scip != NULL);
933  assert(cutdata != NULL);
934  assert(*cutdata != NULL);
935
936  if( (*cutdata)->cut != NULL )
937  {
938  SCIP_CALL(SCIPreleaseRow(scip, &((*cutdata)->cut)));
939  }
940  SCIPfreeBlockMemory(scip, cutdata);
941
942  return SCIP_OKAY;
943 }
944
945
946 /** creates and initializes auxiliary graph node data structures */
947 static
949  SCIP* scip, /**< SCIP data structure */
950  ZEROHALF_AUXGRAPH_NODE** node /**< pointer to store pointer to created data structure */
951  )
952 {
953  assert(scip != NULL);
954  assert(node != NULL);
955
956  SCIP_CALL(SCIPallocBlockMemory(scip, node));
958  (*node)->neighbors = NULL;
959  (*node)->edgeweights = NULL;
960  (*node)->relatedrows = NULL;
961  (*node)->nneighbors = 0;
962
963  (*node)->distance = -1.0;
964  (*node)->previous = NULL;
965
966  return SCIP_OKAY;
967 }
968
969
970 /** frees auxiliary graph node data structures */
971 static
973  SCIP* scip, /**< SCIP data structure */
974  ZEROHALF_AUXGRAPH* graph, /**< graph */
975  ZEROHALF_AUXGRAPH_NODE** node /**< pointer to pointer of data structure */
976  )
977 {
978  int maxnumberofneighbors;
979
980  assert(scip != NULL);
981  assert(node != NULL);
982
983  maxnumberofneighbors = 2 * graph->nnodes - 2;
984
985  SCIPfreeBlockMemoryArrayNull(scip, &((*node)->relatedrows), maxnumberofneighbors);
986  SCIPfreeBlockMemoryArrayNull(scip, &((*node)->edgeweights), maxnumberofneighbors);
987  SCIPfreeBlockMemoryArrayNull(scip, &((*node)->neighbors), maxnumberofneighbors);
988
989  SCIPfreeBlockMemory(scip, node);
990  (*node) = NULL;
991
992  return SCIP_OKAY;
993 }
994
995
996 /** creates and initializes auxiliary graph data structures */
997 static
999  SCIP* scip, /**< SCIP data structure */
1000  ZEROHALF_AUXGRAPH** auxgraph /**< pointer to store pointer to created data structure */
1001  )
1002 {
1003  assert(scip != NULL);
1004  assert(auxgraph != NULL);
1005
1006  SCIP_CALL(SCIPallocBlockMemory(scip, auxgraph));
1008  (*auxgraph)->nodes = NULL;
1009  (*auxgraph)->nodecopies = NULL;
1010
1011  (*auxgraph)->nnodes = 0;
1012
1013  return SCIP_OKAY;
1014 }
1015
1016
1017 /** frees auxiliary graph data structures */
1018 static
1020  SCIP* scip, /**< SCIP data structure */
1021  ZEROHALF_AUXGRAPH** auxgraph /**< pointer to pointer of data structure */
1022  )
1023 {
1024  int n;
1025
1026  assert(scip != NULL);
1027  assert(auxgraph != NULL);
1029  if( (*auxgraph)->nodes != NULL )
1030  {
1031  assert((*auxgraph)->nnodes > 0);
1032  for( n = 0; n < (*auxgraph)->nnodes; ++n )
1033  {
1034  if( (*auxgraph)->nodes[n] != NULL )
1035  {
1036  SCIP_CALL( ZerohalfAuxGraphNodeFree(scip, *auxgraph, &((*auxgraph)->nodes[n])) );
1037  }
1038  }
1039  SCIPfreeBlockMemoryArray(scip, (&(*auxgraph)->nodes), (*auxgraph)->nnodes);
1040  }
1041
1042  if( (*auxgraph)->nodecopies != NULL )
1043  {
1044  assert((*auxgraph)->nnodes > 0);
1045  for( n = 0; n < (*auxgraph)->nnodes; ++n )
1046  {
1047  if( (*auxgraph)->nodecopies[n] != NULL )
1048  {
1049  SCIP_CALL( ZerohalfAuxGraphNodeFree(scip, *auxgraph, &((*auxgraph)->nodecopies[n])) );
1050  }
1051  }
1052  SCIPfreeBlockMemoryArray(scip, (&(*auxgraph)->nodecopies), (*auxgraph)->nnodes);
1053  }
1054
1055  SCIPfreeBlockMemory(scip, auxgraph);
1056
1057  return SCIP_OKAY;
1058 }
1059
1060
1061 /* --------------------------------------------------------------------------------------------------------------------
1062  * local methods: debug
1063  * -------------------------------------------------------------------------------------------------------------------- */
1064
1065 #ifdef SCIP_DEBUG
1066 /** returns a string containing the name of the symbolic constant (given as int value) */
1067 static
1068 char* getconstantname(
1069  char* buffer, /**< string containing the name */
1070  int value /**< symbolic constant given as int value */
1071  )
1072 {
1073  switch( value )
1074  {
1075  case IRRELEVANT:
1076  (void) SCIPsnprintf(buffer, SCIP_MAXSTRLEN, "IRRELEVANT"); break;
1077  case ZERO_ROW:
1078  (void) SCIPsnprintf(buffer, SCIP_MAXSTRLEN, "ZEROROW"); break;
1080  (void) SCIPsnprintf(buffer, SCIP_MAXSTRLEN, "IDENTROW"); break;
1082  (void) SCIPsnprintf(buffer, SCIP_MAXSTRLEN, "SLACK>MAXSLACK"); break;
1084  (void) SCIPsnprintf(buffer, SCIP_MAXSTRLEN, "ISZHCUT"); break;
1085 #ifdef WITHDECOMPOSE
1086  case ROW_IN_SUBPROB_WITHOUT_ODD_RHS:
1087  (void) SCIPsnprintf(buffer, SCIP_MAXSTRLEN, "NOODDRHSROW"); break;
1088 #endif
1089  case NONEXISTENT_ROW:
1090  (void) SCIPsnprintf(buffer, SCIP_MAXSTRLEN, "DOESNOTEXIST"); break;
1091  case NO_RELEVANT_COLUMNS:
1092  (void) SCIPsnprintf(buffer, SCIP_MAXSTRLEN, "HASNORELEVANTCOLS"); break;
1094  (void) SCIPsnprintf(buffer, SCIP_MAXSTRLEN, "SLACK+SODD>MAXSLACK"); break;
1095  case LARGE_COL_EXISTS:
1096  (void) SCIPsnprintf(buffer, SCIP_MAXSTRLEN, "SCOL+SODD>MAXSLACK"); break;
1097  case ZERO_COLUMN:
1098  (void) SCIPsnprintf(buffer, SCIP_MAXSTRLEN, "ZEROCOLUMN"); break;
1100  (void) SCIPsnprintf(buffer, SCIP_MAXSTRLEN, "IDENTCOLUMN"); break;
1101  case ZERO_LP_SOL:
1102  (void) SCIPsnprintf(buffer, SCIP_MAXSTRLEN, "PRIMSOL=0"); break;
1103  case LP_SOL_EQUALS_EVEN_LB:
1104  (void) SCIPsnprintf(buffer, SCIP_MAXSTRLEN, "PRIMSOL=EVENLB"); break;
1105  case LP_SOL_EQUALS_ODD_LB:
1106  (void) SCIPsnprintf(buffer, SCIP_MAXSTRLEN, "PRIMSOL=ODDLB"); break;
1107  case LP_SOL_EQUALS_EVEN_UB:
1108  (void) SCIPsnprintf(buffer, SCIP_MAXSTRLEN, "PRIMSOL=EVENUB"); break;
1109  case LP_SOL_EQUALS_ODD_UB:
1110  (void) SCIPsnprintf(buffer, SCIP_MAXSTRLEN, "PRIMSOL=ODDUB"); break;
1111  case SINGLETON_COLUMN:
1112  (void) SCIPsnprintf(buffer, SCIP_MAXSTRLEN, "SINGLETONCOLUMN"); break;
1113  case CONTINUOUS_VARIABLE:
1114  (void) SCIPsnprintf(buffer, SCIP_MAXSTRLEN, "CONTCOLUMN"); break;
1115  case SMALL_FRACSOL_HEUR:
1116  (void) SCIPsnprintf(buffer, SCIP_MAXSTRLEN, "SUMFRACSOL<DELTA"); break;
1118  (void) SCIPsnprintf(buffer, SCIP_MAXSTRLEN, "NOROWSLEFT"); break;
1119 #ifdef WITHDECOMPOSE
1120  case COLUMN_IN_SUBPROB_WITHOUT_ODD_RHS:
1121  (void) SCIPsnprintf(buffer, SCIP_MAXSTRLEN, "NOODDRHSCOL"); break;
1122 #endif
1123  default:
1124  SCIPerrorMessage("parameter <%s> unknown\n", value);
1125  SCIPABORT();
1126  }
1127
1128  return buffer;
1129 }
1130 #endif
1131
1132
1133 #ifdef ZEROHALF__PRINT_STATISTICS
1134 /** prints the preprocessing statistics of the basic preprocessing applied while
1135  reading the LP data and building basic data structures */
1136 static
1137 SCIP_RETCODE printPreprocessingStatistics(
1138  SCIP* scip, /**< SCIP data structure */
1139  ZEROHALF_LPDATA* lpdata /**< data of current LP relaxation */
1140  )
1141 {
1142  int i;
1143  int nrelevantrows;
1144  int nirrelevantrows;
1145  int nerrors;
1146  int nrelevantcols;
1147  int nirrelevantcols;
1148  int nrows;
1149  int ncols;
1150  int nnonexistingrows;
1151
1152  assert(lpdata != NULL);
1153
1154  nrelevantrows = 0;
1155  nirrelevantrows = 0;
1156  nnonexistingrows = 0;
1157  nrelevantcols = 0;
1158  nirrelevantcols = 0;
1159  nerrors = 0;
1160
1161  for( i = 0 ; i < lpdata->nrows ; ++i)
1162  {
1163  if( lpdata->rrowsindexofleftrow[i] >= 0 )
1164  nrelevantrows++;
1165  else
1166  {
1167  if( lpdata->rrowsindexofleftrow[i] < -199 )
1168  nerrors++;
1169  else
1170  if( lpdata->rrowsindexofleftrow[i] == NONEXISTENT_ROW )
1171  nnonexistingrows++;
1172  else
1173  nirrelevantrows++;
1174  }
1175  if( lpdata->rrowsindexofrightrow[i] >= 0 )
1176  nrelevantrows++;
1177  else
1178  {
1179  if( lpdata->rrowsindexofrightrow[i] < -199 )
1180  nerrors++;
1181  else
1182  {
1183  if( lpdata->rrowsindexofrightrow[i] == NONEXISTENT_ROW )
1184  nnonexistingrows++;
1185  else
1186  nirrelevantrows++;
1187  }
1188  }
1189  }
1190
1191  for( i = 0 ; i < lpdata->ncols ; ++i)
1192  {
1193  if( lpdata->rcolsindexofcol[i] >= 0 )
1194  nrelevantcols++;
1195  else
1196  {
1197  if( lpdata->rcolsindexofcol[i] == -1 )
1198  nirrelevantcols++;
1199  else
1200  {
1201  if( lpdata->rcolsindexofcol[i] > -200 || lpdata->rcolsindexofcol[i] < -299 )
1202  nerrors++;
1203  else
1204  nirrelevantcols++;
1205  }
1206  }
1207  }
1208
1209  nrows = nrelevantrows + nirrelevantrows - nnonexistingrows;
1210  ncols = nrelevantcols + nirrelevantcols;
1211
1212  ZEROHALFstatisticsMessage("\n");
1213  ZEROHALFstatisticsMessage(" | ----- lp data ----- | --- (reductions) -- | --- problem data -- | -lpdata- | -(red.)- | -probd.- | --------\n");
1214  ZEROHALFstatisticsMessage(" | nrows | ncols | ndelrows | ndelcols | nrrows | nrcols | nvarbnds | ndlvbnds | nvarbnds | \n");
1215  ZEROHALFstatisticsMessage("%15s | %8d | %8d | %8d | %8d | %8d | %8d | %8d | %8d | %8d |\n",
1216  "READING LPDATA", nrows, ncols, nirrelevantrows - nnonexistingrows, nirrelevantcols, nrelevantrows, nrelevantcols,
1217  lpdata->nvarbounds, lpdata->ndelvarbounds, lpdata->nvarbounds - lpdata->ndelvarbounds);
1218
1219  return SCIP_OKAY;
1220 }
1221 #endif
1222
1223
1224 #ifdef SCIP_DEBUG
1225 /** prints the considered subproblem */
1226 static
1227 void debugPrintSubLpData(
1228  SCIP* scip, /**< SCIP data structure */
1229  ZEROHALF_LPDATA* lpdata, /**< data of current LP relaxation */
1230  ZEROHALF_SUBLPDATA* sublpdata /**< considered subproblem */
1231  )
1232 {
1233  int i;
1234  int j;
1235
1236  assert(scip != NULL);
1237  assert(lpdata != NULL);
1238  assert(sublpdata != NULL);
1239
1240  SCIPdebugMsg(scip, "\n debugPrintSubLpData:\n\n");
1241  SCIPdebugMsg(scip, " rrows: (nrrows=%d)\n", sublpdata->nrrows);
1242
1243  for( i = 0 ; i < sublpdata->nrrows ; ++i)
1244  {
1245  SCIPdebugMsg(scip, " %6d: rrows: %6d rhs: %6g slack: %6g name: %s\n",
1246  i, sublpdata->rrows[i], sublpdata->rrowsrhs[i], sublpdata->rrowsslack[i],
1247  SCIProwGetName(lpdata->rows[sublpdata->rrows[i]]));
1248  }
1249  SCIPdebugMsg(scip, "\n rcols: (nrcols=%d)\n", sublpdata->nrcols);
1250  for( j = 0 ; j < sublpdata->nrcols ; ++j)
1251  {
1252  SCIPdebugMsg(scip, " %6d: rcols: %6d lbslack: %6g ubslack: %6g\n",
1253  i, sublpdata->rcols[i], sublpdata->rcolslbslack[i], sublpdata->rcolsubslack[i]);
1254  }
1255 }
1256 #endif
1257
1258
1259 #ifdef SCIP_DEBUG
1260 /** prints mod 2 data structures */
1261 static
1262 void debugPrintMod2Data(
1263  SCIP* scip, /**< SCIP data structure */
1264  ZEROHALF_LPDATA* lpdata, /**< data of current LP relaxation */
1265  ZEROHALF_MOD2DATA* mod2data, /**< considered mod 2 data structure */
1266  SCIP_Bool printaggregations /**< should row aggregation bitarrays be printed? */
1267  )
1268 {
1269  int i;
1270  int j;
1271  int k;
1272
1273  assert(scip != NULL);
1274  assert(lpdata != NULL);
1275  assert(mod2data != NULL);
1276
1277  SCIPdebugMsg(scip, "\n debugPrintMod2Data:\n\n");
1278
1279  SCIPdebugMsg(scip, " nrows = %d, nvarbounds = %d, nrcols = %d, nrowsind = %d, ncolsind = %d\n",
1280  mod2data->nrows, mod2data->nvarbounds, mod2data->relatedsubproblem->nrcols,
1281  mod2data->nrowsind, mod2data->ncolsind);
1282  SCIPdebugMsg(scip, " rowsbitarraysize = %d, rowaggregationsbitarraysize = %d\n",
1283  mod2data->rowsbitarraysize, mod2data->rowaggregationsbitarraysize);
1284
1285
1286  SCIPdebugMsg(scip, "\n fracsol:\n");
1287  for( j = 0 ; j < mod2data->relatedsubproblem->nrcols ; ++j )
1288  {
1289  for( k = 0; k < mod2data->ncolsind; ++k )
1290  {
1291  if( mod2data->colsind[k] == j )
1292  break;
1293  SCIPdebugMsg(scip, " rcols[%6d]: fracsol: %6g colsind: %6d name: %s\n", j, mod2data->fracsol[j],
1294  k < mod2data->ncolsind ? k : -1,
1295  SCIPvarGetName(SCIPcolGetVar(lpdata->cols[mod2data->relatedsubproblem->rcols[j]])));
1296  }
1297  }
1298
1299  SCIPdebugMsg(scip, "\n (A mod 2, b mod 2, [#nonz] (slacks), R, name(rrows), left(-)/right(+):\n");
1300  if( mod2data->nrowsind == 0 )
1301  {
1302  SCIPdebugMsg(scip, " empty\n");
1303  }
1304  for( i = 0 ; i < mod2data->nrowsind ; ++i )
1305  {
1306  int nnonz = 0;
1307  SCIPdebugMsg(scip, " ");
1308  for( j = 0 ; j < mod2data->ncolsind; ++j )
1309  if( BITARRAYBITISSET(mod2data->rows[mod2data->rowsind[i]], mod2data->colsind[j]) ) /*lint !e701*/
1310  {
1311  nnonz++;
1312  SCIPdebugMsgPrint(scip, "1");
1313  }
1314  else
1315  {
1316  SCIPdebugMsgPrint(scip, ".");
1317  }
1318  if( mod2data->rhs[mod2data->rowsind[i]] )
1319  {
1320  SCIPdebugMsgPrint(scip, " 1");
1321  }
1322  else
1323  {
1324  SCIPdebugMsgPrint(scip, " 0");
1325  }
1326  SCIPdebugMsgPrint(scip, " [%4d] ", nnonz);
1327  SCIPdebugMsgPrint(scip, "(%6g) ", mod2data->slacks[mod2data->rowsind[i]]);
1328
1329  if( printaggregations )
1330  {
1331  for( j = 0 ; j < mod2data->relatedsubproblem->nrrows; ++j )
1332  if( BITARRAYBITISSET(mod2data->rowaggregations[mod2data->rowsind[i]], j) ) /*lint !e701*/
1333  {
1334  SCIPdebugMsgPrint(scip, "1");
1335  }
1336  else
1337  {
1338  SCIPdebugMsgPrint(scip, ".");
1339  }
1340  }
1341
1342  if( mod2data->rowsind[i] < mod2data->nrows - mod2data->nvarbounds )
1343  {
1344  SCIPdebugMsgPrint(scip, " %s ", SCIProwGetName(lpdata->rows[mod2data->relatedsubproblem->rrows[mod2data->rowsind[i]]]));
1345  if( lpdata->rrowsindexofleftrow[mod2data->relatedsubproblem->rrows[mod2data->rowsind[i]]] == mod2data->rowsind[i] )
1346  {
1347  SCIPdebugMsgPrint(scip, " -\n");
1348  }
1349  else
1350  {
1351  SCIPdebugMsgPrint(scip, " +\n");
1352  }
1353  }
1354  else
1355  {
1356  SCIPdebugMsgPrint(scip, " varbound(rows[%d])\n", mod2data->rowsind[i]);
1357  }
1358  }
1359  SCIPdebugMsgPrint(scip, "\n");
1360 }
1361 #endif
1362
1363
1364 #ifdef SCIP_DEBUG
1365 /** prints LP data */
1366 static
1367 SCIP_RETCODE debugPrintLPRowsAndCols(
1368  SCIP* scip, /**< SCIP data structure */
1369  ZEROHALF_LPDATA* lpdata /**< data of current LP relaxation */
1370  )
1371 {
1372  int i;
1373  int j;
1374  char temp[SCIP_MAXSTRLEN];
1375
1376  assert(scip != NULL);
1377  assert(lpdata != NULL);
1378
1379  SCIPdebugMsg(scip, "\n\nLP rows:\n");
1380  for( i = 0 ; i < lpdata->nrows ; ++i)
1381  {
1382  SCIPdebugMsg(scip, "\nrow %d (left): %s[%d,%d] %s:\n", i,
1383  (lpdata->subproblemsindexofrow[i] == IRRELEVANT)
1384  || (lpdata->rrowsindexofleftrow[i] < 0) ? "IRRELEVANT" : "RELEVANT",
1385  lpdata->subproblemsindexofrow[i], lpdata->rrowsindexofleftrow[i],
1386  lpdata->rrowsindexofleftrow[i] < 0 ? getconstantname(temp, lpdata->rrowsindexofleftrow[i]) : "");
1387  SCIPdebugMsg(scip, "row %d (right): %s[%d,%d] %s:\n", i,
1388  (lpdata->subproblemsindexofrow[i] == IRRELEVANT)
1389  || (lpdata->rrowsindexofrightrow[i] < 0) ? "IRRELEVANT" : "RELEVANT",
1390  lpdata->subproblemsindexofrow[i], lpdata->rrowsindexofrightrow[i],
1391  lpdata->rrowsindexofrightrow[i] < 0 ? getconstantname(temp, lpdata->rrowsindexofrightrow[i]) : "");
1392  SCIP_CALL( SCIPprintRow(scip, lpdata->rows[i], NULL) );
1393  }
1394
1395  SCIPdebugMsg(scip, "\n\nLP cols:\n");
1396  for( j = 0 ; j < lpdata->ncols ; ++j)
1397  {
1398  SCIPdebugMsg(scip, "\ncol %d: %s[%d,%d] %s:\n", j, (lpdata->subproblemsindexofcol[j] == IRRELEVANT)
1399  || (lpdata->rcolsindexofcol[j] < 0) ? "IRRELEVANT" : "RELEVANT",
1400  lpdata->subproblemsindexofcol[j], lpdata->rcolsindexofcol[j],
1401  lpdata->rcolsindexofcol[j] < 0 ? getconstantname(temp, lpdata->rcolsindexofcol[j]) : "");
1402  SCIPdebugMsg(scip, "%s = %f\n", SCIPvarGetName(SCIPcolGetVar(lpdata->cols[j])), SCIPcolGetPrimsol(lpdata->cols[j]));
1403  }
1404
1405  return SCIP_OKAY;
1406 }
1407 #endif
1408
1409
1410 /* --------------------------------------------------------------------------------------------------------------------
1411  * local methods: SCIPsortInd comparators
1412  * -------------------------------------------------------------------------------------------------------------------- */
1413
1414
1415 /** comparator function for sorting an index array non-decreasingly according to a real array */
1416 static
1417 SCIP_DECL_SORTINDCOMP(compRealNonDecreasing)
1418 {
1419  SCIP_Real* scores;
1420
1421  scores = (SCIP_Real*) dataptr;
1422
1423  if( scores[ind1] < scores[ind2] )
1424  return -1;
1425  else if( scores[ind1] > scores[ind2] )
1426  return +1;
1427  else
1428  return 0;
1429 }
1430
1431
1432 /** comparator function for sorting an index array non-increasingly according to a real array */
1433 static
1434 SCIP_DECL_SORTINDCOMP(compRealNonIncreasing)
1435 {
1436  SCIP_Real* scores;
1437
1438  scores = (SCIP_Real*) dataptr;
1439
1440  if( scores[ind1] < scores[ind2] )
1441  return +1;
1442  else if( scores[ind1] > scores[ind2] )
1443  return -1;
1444  else
1445  return 0;
1446 }
1447
1448
1449 /* --------------------------------------------------------------------------------------------------------------------
1450  * local methods: LP data
1451  * -------------------------------------------------------------------------------------------------------------------- */
1452
1453
1454 /** searches for relevant columns, i.e., columns that cannot be deleted because of basic preprocessing methods */
1455 static
1457  SCIP* scip, /**< SCIP data structure */
1458  ZEROHALF_LPDATA* lpdata /**< data of current LP relaxation */
1459  )
1460 {
1461  SCIP_VAR* var;
1462  SCIP_COL* col;
1463  SCIP_Real primsol;
1464  SCIP_Real lb;
1466  SCIP_Real lbslack;
1467  SCIP_Real ubslack;
1468  ZEROHALF_SUBLPDATA* problem;
1469  int j;
1470 #ifdef ZEROHALF__PRINT_STATISTICS
1471  int tempnvarbnds;
1472 #endif
1473  int nsubproblems;
1474
1475  assert(scip != NULL);
1476  assert(lpdata != NULL);
1477  assert(lpdata->cols != NULL);
1478  assert(lpdata->ncols > 0);
1479  assert(lpdata->nrows > 0);
1480  assert(lpdata->subproblems == NULL);
1481  assert(lpdata->nsubproblems == 0);
1482  assert(lpdata->subproblemsindexofrow == NULL);
1483  assert(lpdata->rrowsindexofleftrow == NULL);
1484  assert(lpdata->rrowsindexofrightrow == NULL);
1485  assert(lpdata->subproblemsindexofcol == NULL);
1486  assert(lpdata->rcolsindexofcol == NULL);
1487  assert(lpdata->bestlbidxofcol == NULL);
1488  assert(lpdata->bestubidxofcol == NULL);
1489
1490  nsubproblems = 1;
1491
1492  /* allocate temporary memory for column data structures */
1493  SCIP_CALL(SCIPallocBlockMemoryArray(scip, &(lpdata->subproblems), nsubproblems)); /* create one "sub"problem */
1494  SCIP_CALL(SCIPallocBlockMemoryArray(scip, &(lpdata->subproblemsindexofcol), lpdata->ncols));
1495  SCIP_CALL(SCIPallocBlockMemoryArray(scip, &(lpdata->rcolsindexofcol), lpdata->ncols));
1496  SCIP_CALL(SCIPallocBlockMemoryArray(scip, &(lpdata->bestlbidxofcol), lpdata->ncols));
1497  SCIP_CALL(SCIPallocBlockMemoryArray(scip, &(lpdata->bestubidxofcol), lpdata->ncols));
1498
1499  SCIP_CALL(ZerohalfSubLPDataCreate(scip, &problem));
1500  SCIP_CALL(SCIPallocBlockMemoryArray(scip, &(problem->rcols), lpdata->ncols));
1501  SCIP_CALL(SCIPallocBlockMemoryArray(scip, &(problem->rcolslbslack), lpdata->ncols));
1502  SCIP_CALL(SCIPallocBlockMemoryArray(scip, &(problem->rcolsubslack), lpdata->ncols));
1503  problem->rcolssize = lpdata->ncols;
1504
1505  /* initialize data */
1506  BMSclearMemoryArray(lpdata->subproblemsindexofcol, lpdata->ncols);
1507  BMSclearMemoryArray(lpdata->rcolsindexofcol, lpdata->ncols);
1508
1509  lpdata->nsubproblems = 1;
1510  lpdata->subproblems[0] = problem;
1511  lpdata->nvarbounds = 0;
1512  lpdata->ndelvarbounds = 0;
1513
1514  /* check all cols */
1515  for( j = 0 ; j < lpdata->ncols ; ++j)
1516  {
1517  /* initialize best lb and best ub (-2: undetermined)*/
1518  lpdata->bestlbidxofcol[j] = -2;
1519  lpdata->bestubidxofcol[j] = -2;
1520
1521  col = lpdata->cols[j];
1522  var = SCIPcolGetVar(col);
1523
1524  /* check if vartype is continuous */
1526  {
1527  primsol = SCIPcolGetPrimsol(col);
1528  if( SCIPisFeasZero(scip, primsol) )
1529  primsol = 0.0;
1530
1531  assert(SCIPisFeasEQ(scip, SCIPgetVarSol(scip, var), primsol));
1532
1533  lb = SCIPcolGetLb(col);
1534  ub = SCIPcolGetUb(col);
1535  lbslack = primsol - lb;
1536  ubslack = ub - primsol;
1537  if( SCIPisFeasZero(scip, lbslack) )
1538  lbslack = 0.0;
1539  if( SCIPisFeasZero(scip, ubslack) )
1540  ubslack = 0.0;
1541  assert(SCIPisLE(scip, lb, ub));
1542  assert(!SCIPisNegative(scip, lbslack));
1543  assert(!SCIPisNegative(scip, ubslack));
1544
1545 #ifdef ZEROHALF__PRINT_STATISTICS
1546  tempnvarbnds = 0;
1547  if( !SCIPisInfinity(scip, (-1) * lb) )
1548  tempnvarbnds++;
1549  if( !SCIPisInfinity(scip, ub) )
1550  tempnvarbnds++;
1551  lpdata->nvarbounds += tempnvarbnds;
1552 #endif
1553
1554  if( SCIPisNegative(scip, lb) )
1555  {
1556  /* column is declared to be irrelevant because its lower bound is negative,
1557  * variable would have to be shifted, complemented or decomposed.
1558  */
1559  /**@todo consider general integers with negative lower bounds and transform to positive representation
1560  * and propagate through corresponding rows. In the current version, redundant inequalities might be
1561  * considered as cut candidates and valid cuts might be missed, but no wrong cuts should be produced
1562  * (due to SCIPcalcMIR) this leads to performance deterioration in the (rare) case of general integers
1563  * with negative bounds.
1564  */
1565  lpdata->subproblemsindexofcol[j] = IRRELEVANT;
1566  lpdata->rcolsindexofcol[j] = IRRELEVANT;
1567  }
1568  else if( !SCIPisZero(scip, primsol) )
1569  {
1570  if( SCIPisZero(scip, lbslack) )
1571  {
1572 #ifdef ZEROHALF__PRINT_STATISTICS
1573  lpdata->ndelvarbounds += tempnvarbnds;
1574 #endif
1575  lpdata->subproblemsindexofcol[j] = IRRELEVANT;
1576  if( ISODD(scip, lb) )
1577  lpdata->rcolsindexofcol[j] = LP_SOL_EQUALS_ODD_LB;
1578  else
1579  lpdata->rcolsindexofcol[j] = LP_SOL_EQUALS_EVEN_LB;
1580  }
1581  else
1582  {
1583  if( SCIPisZero(scip, ubslack) )
1584  {
1585 #ifdef ZEROHALF__PRINT_STATISTICS
1586  lpdata->ndelvarbounds += tempnvarbnds;
1587 #endif
1588  lpdata->subproblemsindexofcol[j] = IRRELEVANT;
1589  if( ISODD(scip, ub) )
1590  lpdata->rcolsindexofcol[j] = LP_SOL_EQUALS_ODD_UB;
1591  else
1592  lpdata->rcolsindexofcol[j] = LP_SOL_EQUALS_EVEN_UB;
1593  }
1594  else
1595  {
1596  /* relevant col was found */
1597  problem->rcols[problem->nrcols] = j;
1598  problem->rcolslbslack[problem->nrcols] = lbslack;
1599  problem->rcolsubslack[problem->nrcols] = ubslack;
1600  lpdata->subproblemsindexofcol[j] = 0;
1601  lpdata->rcolsindexofcol[j] = problem->nrcols;
1602  problem->nrcols++;
1603  }
1604  }
1605  }
1606  else
1607  {
1608 #ifdef ZEROHALF__PRINT_STATISTICS
1609  lpdata->ndelvarbounds += tempnvarbnds;
1610 #endif
1611  lpdata->subproblemsindexofcol[j] = IRRELEVANT;
1612  lpdata->rcolsindexofcol[j] = ZERO_LP_SOL;
1613  }
1614  }
1615  else
1616  {
1617  /* column is irrelevant because vartype is continuous (is handled in getRelevantRows())*/
1618  lpdata->subproblemsindexofcol[j] = IRRELEVANT;
1619  lpdata->rcolsindexofcol[j] = CONTINUOUS_VARIABLE;
1620  }
1621  }
1622
1623  return SCIP_OKAY;
1624 }
1625
1626
1627 /** finds closest lower bound of col and stores it within lpdata;
1628  * the bound can be the lower bound or the best variable lower bound with nonnegative column variable
1629  */
1630 static
1631 void findClosestLb(
1632  SCIP* scip, /**< SCIP data structure */
1633  ZEROHALF_LPDATA* lpdata, /**< data of current LP relaxation */
1634  SCIP_COL* col, /**< column to get closest lower bound */
1635  SCIP_Real* bestlbsol, /**< pointer to store value of closest lower bound */
1636  int* bestlbtype, /**< pointer to store type of closest lower bound */
1637  SCIP_VAR** bestzvlb, /**< pointer to store variable z in closest variable lower bound b*z + d */
1638  SCIP_Real* bestbvlb, /**< pointer to store coefficient b in closest variable lower bound b*z + d */
1639  SCIP_Real* bestdvlb /**< pointer to store constant d in closest variable lower bound b*z + d */
1641  )
1642 {
1643  SCIP_VAR* var;
1644  SCIP_VAR** zvlb;
1645  SCIP_Real* bvlb;
1646  SCIP_Real* dvlb;
1647  int collppos;
1648  int nvlb;
1649  int j;
1650
1651  assert(lpdata != NULL);
1652  assert(bestlbsol != NULL);
1653  assert(bestlbtype != NULL);
1654  assert(bestzvlb != NULL);
1655  assert(bestbvlb != NULL);
1656  assert(bestdvlb != NULL);
1657
1658
1659  collppos = SCIPcolGetLPPos(col);
1660  var = SCIPcolGetVar(col);
1661  *bestlbsol = SCIPcolGetLb(col);
1662  *bestlbtype = lpdata->bestlbidxofcol[collppos];
1663  *bestzvlb = NULL;
1664  *bestbvlb = 0.0;
1665  *bestdvlb = 0.0;
1666
1667  if( *bestlbtype == -1 )
1668  return;
1669
1670  if( USEVARBOUNDS ) /*lint !e774 !e506*/
1671  {
1672  nvlb = SCIPvarGetNVlbs(var);
1673  zvlb = SCIPvarGetVlbVars(var);
1674  bvlb = SCIPvarGetVlbCoefs(var);
1675  dvlb = SCIPvarGetVlbConstants(var);
1676
1677  assert(zvlb != NULL || nvlb == 0);
1678  assert(bvlb != NULL || nvlb == 0);
1679  assert(dvlb != NULL || nvlb == 0);
1680  }
1681
1682  if( *bestlbtype == -2 )
1683  {
1684  if( USEVARBOUNDS ) /*lint !e774 !e506*/
1685  {
1686
1687  /* search for lb or vlb with maximal bound value */
1688  for( j = 0; j < nvlb; j++ )
1689  {
1690  assert(zvlb != NULL);
1691  assert(bvlb != NULL);
1692  assert(dvlb != NULL);
1693  assert(SCIPvarGetType(zvlb[j]) != SCIP_VARTYPE_CONTINUOUS);
1694
1695  /* use only vlb with nonnegative variable z that are column variables and present in the current LP */
1696  if( SCIPvarGetStatus(zvlb[j]) == SCIP_VARSTATUS_COLUMN && SCIPcolIsInLP(SCIPvarGetCol(zvlb[j])) &&
1697  !SCIPisNegative(scip, SCIPcolGetLb(SCIPvarGetCol(zvlb[j]))) )
1698  {
1699  SCIP_Real vlbsol;
1700
1701  vlbsol = bvlb[j] * SCIPcolGetPrimsol(SCIPvarGetCol(zvlb[j])) + dvlb[j];
1702  if( vlbsol > *bestlbsol )
1703  {
1704  *bestlbsol = vlbsol;
1705  *bestlbtype = j;
1706  }
1707  }
1708  }
1709  }
1710
1711  /* if no better var bound could be found, set type to the fixed bound (-1) */
1712  if( *bestlbtype == -2 )
1713  *bestlbtype = -1;
1714
1715  /* store best bound for substitution */
1716  lpdata->bestlbidxofcol[collppos] = *bestlbtype;
1717  }
1718  assert(lpdata->bestlbidxofcol[collppos] > -2);
1719
1720  if( *bestlbtype >= 0 )
1721  {
1722  assert(USEVARBOUNDS); /*lint !e774 !e506*/
1723  assert(*bestlbtype < nvlb);
1724  assert(zvlb != NULL);
1725  assert(bvlb != NULL);
1726  assert(dvlb != NULL);
1727  *bestzvlb = zvlb[*bestlbtype];
1728  *bestbvlb = bvlb[*bestlbtype];
1729  *bestdvlb = dvlb[*bestlbtype];
1730  }
1731 }
1732
1733
1734 /** finds closest upper bound of col and stores it within lpdata;
1735  * the bound can be the upper bound or the best variable upper bound with nonnegative column variable
1736  */
1737 static
1738 void findClosestUb(
1739  SCIP* scip, /**< SCIP data structure */
1740  ZEROHALF_LPDATA* lpdata, /**< data of current LP relaxation */
1741  SCIP_COL* col, /**< column to get closest upper bound */
1742  SCIP_Real* bestubsol, /**< pointer to store value of closest upper bound */
1743  int* bestubtype, /**< pointer to store type of closest upper bound */
1744  SCIP_VAR** bestzvub, /**< pointer to store variable z in closest variable upper bound b*z + d */
1745  SCIP_Real* bestbvub, /**< pointer to store coefficient b in closest variable upper bound b*z + d */
1746  SCIP_Real* bestdvub /**< pointer to store constant d in closest variable upper bound b*z + d */
1748  )
1749 {
1750  SCIP_VAR* var;
1751  SCIP_VAR** zvub;
1752  SCIP_Real* bvub;
1753  SCIP_Real* dvub;
1754  int collppos;
1755  int nvub;
1756  int j;
1757
1758  assert(lpdata != NULL);
1759  assert(bestubsol != NULL);
1760  assert(bestubtype != NULL);
1761  assert(bestzvub != NULL);
1762  assert(bestbvub != NULL);
1763  assert(bestdvub != NULL);
1764
1765  collppos = SCIPcolGetLPPos(col);
1766  var = SCIPcolGetVar(col);
1767  *bestubsol = SCIPcolGetUb(col);
1768  *bestubtype = lpdata->bestubidxofcol[collppos];
1769  *bestzvub = NULL;
1770  *bestbvub = 0.0;
1771  *bestdvub = 0.0;
1772
1773  if( *bestubtype == -1 )
1774  return;
1775
1776  if( USEVARBOUNDS ) /*lint !e774 !e506*/
1777  {
1778  nvub = SCIPvarGetNVubs(var);
1779  zvub = SCIPvarGetVubVars(var);
1780  bvub = SCIPvarGetVubCoefs(var);
1781  dvub = SCIPvarGetVubConstants(var);
1782
1783  assert(zvub != NULL || nvub == 0);
1784  assert(bvub != NULL || nvub == 0);
1785  assert(dvub != NULL || nvub == 0);
1786  }
1787
1788  if( *bestubtype == -2 )
1789  {
1790  if( USEVARBOUNDS ) /*lint !e774 !e506*/
1791  {
1792  /* search for ub or vub with maximal bound value */
1793  for( j = 0; j < nvub; j++ )
1794  {
1795  assert(zvub != NULL);
1796  assert(bvub != NULL);
1797  assert(dvub != NULL);
1798  assert(SCIPvarGetType(zvub[j]) != SCIP_VARTYPE_CONTINUOUS);
1799
1800  /* use only vub with nonnegative variable z that are column variables and present in the current LP */
1801  if( SCIPvarGetStatus(zvub[j]) == SCIP_VARSTATUS_COLUMN && SCIPcolIsInLP(SCIPvarGetCol(zvub[j])) &&
1802  !SCIPisNegative(scip, SCIPcolGetUb(SCIPvarGetCol(zvub[j]))) )
1803  {
1804  SCIP_Real vubsol;
1805
1806  vubsol = bvub[j] * SCIPcolGetPrimsol(SCIPvarGetCol(zvub[j])) + dvub[j];
1807  if( vubsol < *bestubsol )
1808  {
1809  *bestubsol = vubsol;
1810  *bestubtype = j;
1811  }
1812  }
1813  }
1814  }
1815  /* if no better var bound could be found, set type to the fixed bound (-1) */
1816  if( *bestubtype == -2 )
1817  *bestubtype = -1;
1818
1819  /* store best bound for substitution */
1820  lpdata->bestubidxofcol[collppos] = *bestubtype;
1821  }
1822  assert(lpdata->bestubidxofcol[collppos] > -2);
1823
1824  if( *bestubtype >= 0 )
1825  {
1826  assert(USEVARBOUNDS); /*lint !e774 !e506*/
1827  assert(*bestubtype < nvub);
1828  assert(zvub != NULL);
1829  assert(bvub != NULL);
1830  assert(dvub != NULL);
1831  *bestzvub = zvub[*bestubtype];
1832  *bestbvub = bvub[*bestubtype];
1833  *bestdvub = dvub[*bestubtype];
1834  }
1835 }
1836
1837
1838 /** searches for relevant rows, i.e., rows containing relevant columns that cannot be deleted because of basic
1839  * preprocessing methods
1840  */
1841 static
1843  SCIP* scip, /**< SCIP data structure */
1845  ZEROHALF_LPDATA* lpdata /**< data of current LP relaxation */
1846  )
1847 {
1848  SCIP_COL** colscurrentrow;
1849  SCIP_Real* valscurrentrow;
1850  SCIP_Real* densecoeffscurrentleftrow;
1851  SCIP_Real* densecoeffscurrentrightrow;
1852  SCIP_ROW* row;
1853  SCIP_VAR* var;
1854  SCIP_VAR* bestzvbnd;
1855  SCIP_Real bestbndsol;
1856  SCIP_Real bestbvbnd;
1857  SCIP_Real bestdvbnd;
1858  SCIP_Real intscalarleftrow;
1859  SCIP_Real intscalarrightrow;
1860  SCIP_Real act;
1861  SCIP_Real cst;
1862  SCIP_Real lhs;
1863  SCIP_Real lhsslack;
1864  SCIP_Real maxslack;
1865  SCIP_Real rhs;
1866  SCIP_Real rhsslack;
1867  int bestbndtype;
1868  int nnonzcurrentrow;
1869  int c;
1870  int r;
1871  int k;
1872  SCIP_Bool lhsslackislessequalmaxslack;
1873  SCIP_Bool rhsslackislessequalmaxslack;
1874  SCIP_Bool lhsisinfinity;
1875  SCIP_Bool rhsisinfinity;
1876  SCIP_Bool lhsiseven;
1877  SCIP_Bool rhsiseven;
1878  SCIP_Bool success;
1879  ZEROHALF_SUBLPDATA* problem;
1880
1881  int collppos;
1882  SCIP_Bool rowisrelevant;
1883
1884  assert(scip != NULL);
1886  assert(lpdata != NULL);
1887  assert(lpdata->rows != NULL);
1888  assert(lpdata->nrows > 0);
1889  assert(lpdata->subproblemsindexofcol != NULL);
1890  assert(lpdata->rcolsindexofcol != NULL);
1891  assert(lpdata->subproblems != NULL);
1892  assert(lpdata->nsubproblems > 0);
1893
1894  assert(lpdata->subproblemsindexofrow == NULL);
1895  assert(lpdata->rrowsindexofleftrow == NULL);
1896  assert(lpdata->rrowsindexofrightrow == NULL);
1897
1898  k = 0;
1899  problem = lpdata->subproblems[k];
1900
1901  assert(k >= 0);
1902  assert(k <= lpdata->nsubproblems);
1903  assert(problem != NULL);
1904  assert(problem->rcols != NULL);
1905  assert(problem->nrcols > 0);
1906  assert(problem->rcolslbslack != NULL);
1907  assert(problem->rcolsubslack != NULL);
1908
1909  assert(problem->rrows == NULL);
1910  assert(problem->rrowsrhs == NULL);
1911  assert(problem->rrowsslack == NULL);
1912
1913  /* allocate temporary memory for row data structures */
1914  SCIP_CALL(SCIPallocBlockMemoryArray(scip, &(lpdata->subproblemsindexofrow), lpdata->nrows));
1915  SCIP_CALL(SCIPallocBlockMemoryArray(scip, &(lpdata->rrowsindexofleftrow), lpdata->nrows));
1916  SCIP_CALL(SCIPallocBlockMemoryArray(scip, &(lpdata->rrowsindexofrightrow), lpdata->nrows));
1917
1918  SCIP_CALL(SCIPallocBlockMemoryArray(scip, &(problem->rrows), 2 * lpdata->nrows));
1919  SCIP_CALL(SCIPallocBlockMemoryArray(scip, &(problem->rrowsrhs), 2 * lpdata->nrows));
1920  SCIP_CALL(SCIPallocBlockMemoryArray(scip, &(problem->rrowsslack), 2 * lpdata->nrows));
1921  problem->rrowssize = 2 * lpdata->nrows;
1922
1923  /* allocate temporary memory */
1924  SCIP_CALL(SCIPallocBufferArray(scip, &densecoeffscurrentleftrow, lpdata->ncols));
1925  SCIP_CALL(SCIPallocBufferArray(scip, &densecoeffscurrentrightrow, lpdata->ncols));
1926
1927  /* initialize arrays */
1928  BMSclearMemoryArray(lpdata->subproblemsindexofrow, lpdata->nrows);
1929  BMSclearMemoryArray(lpdata->rrowsindexofleftrow, lpdata->nrows);
1930  BMSclearMemoryArray(lpdata->rrowsindexofrightrow, lpdata->nrows);
1931
1932  BMSclearMemoryArray(densecoeffscurrentleftrow, lpdata->ncols);
1933  BMSclearMemoryArray(densecoeffscurrentrightrow, lpdata->ncols);
1934
1936  problem->nrrows = 0;
1937  for( r = 0 ; r < lpdata->nrows ; ++r)
1938  {
1939  row = lpdata->rows[r];
1940
1942  {
1944  const char* rowname = SCIProwGetName(row);
1945
1946  if( strlen(rowname) > 8 )
1947  {
1948  if(rowname[0] == 'z'
1949  && rowname[1] == 'e'
1950  && rowname[2] == 'r'
1951  && rowname[3] == 'o'
1952  && rowname[4] == 'h'
1953  && rowname[5] == 'a'
1954  && rowname[6] == 'l'
1955  && rowname[7] == 'f' )
1956  {
1957  lpdata->subproblemsindexofrow[r] = IRRELEVANT;
1958  lpdata->rrowsindexofleftrow[r] = NONEXISTENT_ROW;
1959  lpdata->rrowsindexofrightrow[r] = NONEXISTENT_ROW;
1960  continue;
1961  }
1962  }
1963  }
1964
1965  /* check if current row is an original LP row (i.e., was not added) if necessary */
1967  {
1968  int left;
1969  int center;
1970  int right;
1971  int rowindex;
1972
1975
1976  left = 0;
1977  center = 1;
1978  right = sepadata->norigrows - 1;
1979  rowindex = SCIProwGetIndex(row);
1980  while( left <= right && center > -1 )
1981  {
1982  center = left + ((right - left) / 2);
1983  if( sepadata->origrows[center] == rowindex )
1984  center = -1;
1985  if( sepadata->origrows[center] > rowindex )
1986  right = center - 1;
1987  else
1988  left = center + 1;
1989  }
1990  if( center > -1 )
1991  {
1992  lpdata->subproblemsindexofrow[r] = IRRELEVANT;
1993  lpdata->rrowsindexofleftrow[r] = NONEXISTENT_ROW;
1994  lpdata->rrowsindexofrightrow[r] = NONEXISTENT_ROW;
1995  continue;
1996  }
1997  }
1998
1999  /* get row data */
2000  colscurrentrow = SCIProwGetCols(row);
2001  nnonzcurrentrow = SCIProwGetNLPNonz(row);
2002  valscurrentrow = SCIProwGetVals(row);
2003
2004  /* clear dense coeffs arrays */
2005  BMSclearMemoryArray(densecoeffscurrentleftrow, lpdata->ncols);
2006  BMSclearMemoryArray(densecoeffscurrentrightrow, lpdata->ncols);
2007
2008  /* calculate dense coeffs arrays */
2009  for( c = 0; c < nnonzcurrentrow; ++c)
2010  {
2011  collppos = SCIPcolGetLPPos(colscurrentrow[c]);
2012  assert(0 <= collppos && collppos < lpdata->ncols);
2013
2014  var = SCIPcolGetVar(colscurrentrow[c]);
2015
2016  /* check if row contains a continuous variable */
2018  {
2019  bestzvbnd = NULL;
2020  bestbndsol = 0.0;
2021  bestbvbnd = 0.0;
2022  bestdvbnd = 0.0;
2023  bestbndtype = -2;
2024
2025  /* Consider rhs of row and relax continuous variables by substituting for:
2026  * - a_j > 0: x_j = lb or x_j = b*z + d with variable lower bound b*z + d with column var z >= 0
2027  * - a_j < 0: x_j = ub or x_j = b*z + d with variable upper bound b*z + d with column var z >= 0
2028  * and
2029  * consider lhs of row and relax continuous variables by substituting for:
2030  * - a_j < 0: x_j = lb or x_j = b*z + d with variable lower bound b*z + d with column var z >= 0
2031  * - a_j > 0: x_j = ub or x_j = b*z + d with variable upper bound b*z + d with column var z >= 0
2032  */
2033
2034  findClosestLb(scip, lpdata, colscurrentrow[c],
2035  &bestbndsol, &bestbndtype, &bestzvbnd, &bestbvbnd, &bestdvbnd );
2036  assert( bestbndtype > -2 && lpdata->bestlbidxofcol[collppos] == bestbndtype);
2037
2038  if( bestbndtype > -1 )
2039  {
2040  int zlppos;
2041
2042  zlppos = SCIPcolGetLPPos(SCIPvarGetCol(bestzvbnd));
2043  assert(0 <= zlppos && zlppos < lpdata->ncols);
2044
2045  if( valscurrentrow[c] > 0 )
2046  densecoeffscurrentrightrow[zlppos] += valscurrentrow[c] * bestbvbnd;
2047  else
2048  densecoeffscurrentleftrow[zlppos] -= valscurrentrow[c] * bestbvbnd;
2049  }
2050
2051  findClosestUb(scip, lpdata, colscurrentrow[c],
2052  &bestbndsol, &bestbndtype, &bestzvbnd, &bestbvbnd, &bestdvbnd );
2053  assert(bestbndtype > -2 && lpdata->bestubidxofcol[collppos] == bestbndtype);
2054
2055  if( bestbndtype > -1 )
2056  {
2057  int zlppos;
2058
2059  zlppos = SCIPcolGetLPPos(SCIPvarGetCol(bestzvbnd));
2060  assert(0 <= zlppos && zlppos < lpdata->ncols);
2061
2062  if( valscurrentrow[c] > 0 )
2063  densecoeffscurrentleftrow[zlppos] -= valscurrentrow[c] * bestbvbnd;
2064  else
2065  densecoeffscurrentrightrow[zlppos] += valscurrentrow[c] * bestbvbnd;
2066  }
2067  }
2068  else
2069  {
2070  densecoeffscurrentleftrow[collppos] -= valscurrentrow[c];
2071  densecoeffscurrentrightrow[collppos] += valscurrentrow[c];
2072  }
2073  }
2074
2075  /* calculate scalar that would make (left|right) row coefficients integral;
2076  * try to avoid unnecessary or expensive scaling calls
2077  */
2078  intscalarleftrow = 1.0;
2079  intscalarrightrow = 1.0;
2080
2081  if( sepadata->scalefraccoeffs && (!SCIPisIntegral(scip, SCIPgetRowMinCoef(scip, row))
2082  || !SCIPisIntegral(scip, SCIPgetRowMaxCoef(scip, row)) || !SCIPisIntegral(scip, SCIProwGetSumNorm(row))) )
2083  {
2084  SCIP_CALL( SCIPcalcIntegralScalar(densecoeffscurrentleftrow, lpdata->ncols,
2085  -SCIPepsilon(scip), SCIPepsilon(scip), (SCIP_Longint) MAXDNOM, MAXSCALE, &intscalarleftrow, &success) );
2086  if( !success )
2087  {
2088  lpdata->rrowsindexofleftrow[r] = NONEXISTENT_ROW;
2089  }
2090
2091  SCIP_CALL( SCIPcalcIntegralScalar(densecoeffscurrentrightrow, lpdata->ncols,
2092  -SCIPepsilon(scip), SCIPepsilon(scip), (SCIP_Longint) MAXDNOM, MAXSCALE, &intscalarrightrow, &success) );
2093  if( !success )
2094  {
2095  lpdata->rrowsindexofrightrow[r] = NONEXISTENT_ROW;
2096  }
2097
2098  if ( lpdata->rrowsindexofleftrow[r] == NONEXISTENT_ROW
2099  && lpdata->rrowsindexofrightrow[r] == NONEXISTENT_ROW )
2100  {
2101  lpdata->subproblemsindexofrow[r] = IRRELEVANT;
2102  continue;
2103  }
2104  }
2105
2106  lpdata->intscalarsleftrow[r] = intscalarleftrow;
2107  lpdata->intscalarsrightrow[r] = intscalarrightrow;
2108
2109  /* calculate lhs/rhs & slacks */
2110  act = SCIPgetRowLPActivity(scip, row);
2111  lhs = SCIProwGetLhs(row);
2112  rhs = SCIProwGetRhs(row);
2113  cst = SCIProwGetConstant(row);
2114
2115  lhsisinfinity = SCIPisInfinity(scip, -lhs);
2116  rhsisinfinity = SCIPisInfinity(scip, rhs);
2117
2118  lhsslack = SCIPisFeasZero(scip, act - lhs) ? 0.0 : act - lhs;
2119  rhsslack = SCIPisFeasZero(scip, rhs - act) ? 0.0 : rhs - act;
2120
2121  lhs = (lhs - cst) * intscalarleftrow;
2122  rhs = (rhs - cst) * intscalarrightrow;
2123
2124  lhsisinfinity = lhsisinfinity || SCIPisInfinity(scip, -lhs);
2125  rhsisinfinity = rhsisinfinity || SCIPisInfinity(scip, rhs);
2126
2127  lhsslack = lhsslack * intscalarleftrow;
2128  rhsslack = rhsslack * intscalarrightrow;
2129
2130  /* check if the slack value of the row is small enough */
2131  if( (!lhsisinfinity && SCIPisLE(scip, lhsslack, maxslack))
2132  || (!rhsisinfinity && SCIPisLE(scip, rhsslack, maxslack)) )
2133  {
2134  colscurrentrow = SCIProwGetCols(row);
2135  nnonzcurrentrow = SCIProwGetNLPNonz(row);
2136  valscurrentrow = SCIProwGetVals(row);
2137
2138  lhsiseven = ISEVEN(scip, lhs);
2139  rhsiseven = ISEVEN(scip, rhs);
2140
2141  rowisrelevant = FALSE;
2142  for( c = 0 ; c < nnonzcurrentrow ; ++c )
2143  {
2144  collppos = SCIPcolGetLPPos(colscurrentrow[c]);
2145  var = SCIPcolGetVar(colscurrentrow[c]);
2146
2147  /* check if row contains a column with primsol = odd bound and update lhs/rhs parity */
2148  if( lpdata->rcolsindexofcol[collppos] == LP_SOL_EQUALS_ODD_LB
2149  || lpdata->rcolsindexofcol[collppos] == LP_SOL_EQUALS_ODD_UB )
2150  {
2151  assert(SCIPvarGetType(SCIPcolGetVar(colscurrentrow[c])) != SCIP_VARTYPE_CONTINUOUS);
2152  lhsiseven = !lhsiseven;
2153  rhsiseven = !rhsiseven;
2154  }
2155
2156  /* check if row contains a continuous variable */
2158  {
2159  if( !rhsisinfinity )
2160  {
2161  if( valscurrentrow[c] * intscalarrightrow > 0 )
2162  {
2163  findClosestLb(scip, lpdata, colscurrentrow[c], &bestbndsol, &bestbndtype, &bestzvbnd, &bestbvbnd, &bestdvbnd );
2164  assert(bestbndtype > -2 && lpdata->bestlbidxofcol[collppos] == bestbndtype);
2165  }
2166  else
2167  {
2168  assert(valscurrentrow[c] * intscalarrightrow < 0);
2169  findClosestUb(scip, lpdata, colscurrentrow[c], &bestbndsol, &bestbndtype, &bestzvbnd, &bestbvbnd, &bestdvbnd );
2170  assert(bestbndtype > -2 && lpdata->bestubidxofcol[collppos] == bestbndtype);
2171  }
2172  assert(bestbndtype == -1 || bestzvbnd != NULL);
2173
2174  if( SCIPisInfinity(scip, -bestbndsol) || SCIPisInfinity(scip, bestbndsol) )
2175  rhsisinfinity = TRUE;
2176  else
2177  {
2178  /**@todo check whether REALABS is really correct */
2179  if ( bestbndtype == -1 )
2180  rhs -= intscalarrightrow * REALABS(valscurrentrow[c]) * bestbndsol;
2181  else
2182  rhs -= intscalarrightrow * REALABS(valscurrentrow[c]) * bestdvbnd;
2183  rhsslack += intscalarrightrow * REALABS(valscurrentrow[c])
2184  * REALABS(SCIPcolGetPrimsol(colscurrentrow[c]) - bestbndsol);
2185  assert(SCIPisGE(scip, rhsslack, 0.0));
2186  rhsisinfinity = rhsisinfinity || SCIPisInfinity(scip, rhs);
2187  }
2188  }
2189
2190  if( !lhsisinfinity )
2191  {
2192  if( valscurrentrow[c] * intscalarleftrow < 0 )
2193  {
2194  findClosestLb(scip, lpdata, colscurrentrow[c], &bestbndsol, &bestbndtype, &bestzvbnd, &bestbvbnd, &bestdvbnd );
2195  assert(bestbndtype > -2 && lpdata->bestlbidxofcol[collppos] == bestbndtype);
2196  }
2197  else
2198  {
2199  assert(valscurrentrow[c] * intscalarleftrow > 0);
2200  findClosestUb(scip, lpdata, colscurrentrow[c], &bestbndsol, &bestbndtype, &bestzvbnd, &bestbvbnd, &bestdvbnd );
2201  assert(bestbndtype > -2 && lpdata->bestubidxofcol[collppos] == bestbndtype);
2202  }
2203  assert(bestbndtype == -1 || bestzvbnd != NULL);
2204
2205  if( SCIPisInfinity(scip, -bestbndsol) || SCIPisInfinity(scip, bestbndsol) )
2206  lhsisinfinity = TRUE;
2207  else
2208  {
2209  /**@todo check whether REALABS is really correct */
2210  if( bestbndtype == -1 )
2211  lhs -= intscalarleftrow * REALABS(valscurrentrow[c]) * bestbndsol;
2212  else
2213  lhs -= intscalarleftrow * REALABS(valscurrentrow[c]) * bestdvbnd;
2214  lhsslack += intscalarleftrow * REALABS(valscurrentrow[c])
2215  * REALABS(SCIPcolGetPrimsol(colscurrentrow[c]) - bestbndsol);
2216  assert(SCIPisGE(scip, lhsslack, 0.0));
2217  lhsisinfinity = lhsisinfinity || SCIPisInfinity(scip, lhs);
2218  }
2219  }
2220
2221  /* if both lhs and rhs became infinity, then the (relaxed) row is not relevant */
2222  if( lhsisinfinity && rhsisinfinity )
2223  {
2224  rowisrelevant = FALSE;
2225  break;
2226  }
2227
2228  rowisrelevant = TRUE;
2229  }
2230
2231  /* check if row contains at least one relevant column (because k == 0 by initialization) */
2232  if( lpdata->subproblemsindexofcol[collppos] == k )
2233  rowisrelevant = TRUE;
2234
2235  /* check if row contains no relevant columns but an odd lhs or rhs value */
2236  if( c == nnonzcurrentrow - 1 && (!lhsiseven || !rhsiseven) )
2237  rowisrelevant = TRUE;
2238
2239  }
2240  assert(SCIPisGE(scip, lhsslack, 0.0));
2241  assert(SCIPisGE(scip, rhsslack, 0.0));
2242
2243
2244  /* process row if it is relevant */
2245  if( rowisrelevant )
2246  {
2247  /* row is relevant because it contains a relevant column */
2248  problem->rrows[problem->nrrows] = r;
2249
2250  lhsslackislessequalmaxslack = SCIPisLE(scip, lhsslack, maxslack);
2251  rhsslackislessequalmaxslack = SCIPisLE(scip, rhsslack, maxslack);
2252
2253  /* note: due to the relaxation of continuous variables with their bounds the coeffs of nonzero variables
2254  * in left row and right row may be different. hence the row with smaller slack cannot be removed without
2255  * checking the coeffs first.
2256  */
2257  if( !lhsisinfinity && lhsslackislessequalmaxslack )
2258  {
2259  /* "-a^T x <= -lhs" */
2260  lpdata->subproblemsindexofrow[r] = k;
2261  lpdata->rrowsindexofleftrow[r] = problem->nrrows;
2262
2263  problem->rrows[problem->nrrows] = r;
2264  /**@todo check whether lhs is correct, or whether this must be -lhs. do we store the
2265  * -ax <= -lhs constraint or the ax >= lhs constraint? Is this handled correctly above while updating lhs?
2266  */
2267  problem->rrowsrhs[problem->nrrows] = lhs;
2268  problem->rrowsslack[problem->nrrows] = lhsslack;
2269  }
2270  else
2271  {
2272  lpdata->rrowsindexofleftrow[r] = lhsisinfinity ? NONEXISTENT_ROW : SLACK_GREATER_THAN_MAXSLACK;
2273  }
2274  /**@todo check the following: if !lhsinfinity AND !rhsinfinity: then only one of them is stored currently,
2275  * because problem->nrrows++ is only called once. if this is intended, why do we allocate 2 * lpdata->nrows
2276  * entries for rrows?
2277  */
2278  if( !rhsisinfinity && rhsslackislessequalmaxslack )
2279  {
2280  /* "a^T x <= rhs" */
2281  lpdata->subproblemsindexofrow[r] = k;
2282  lpdata->rrowsindexofrightrow[r] = problem->nrrows;
2283
2284  problem->rrows[problem->nrrows] = r;
2285  problem->rrowsrhs[problem->nrrows] = rhs;
2286  problem->rrowsslack[problem->nrrows] = rhsslack;
2287  }
2288  else
2289  {
2290  lpdata->rrowsindexofrightrow[r] = rhsisinfinity ? NONEXISTENT_ROW : SLACK_GREATER_THAN_MAXSLACK;
2291  }
2292
2293  /* increase counter only if at least one half row had a sufficiently small slack */
2294  if( lpdata->rrowsindexofleftrow[r] > -1 || lpdata->rrowsindexofrightrow[r] > -1 )
2295  problem->nrrows++;
2296  }
2297  else /* case: !rowisrelevant */
2298  {
2299  /* row does not contain relevant columns */
2300  lpdata->subproblemsindexofrow[r] = IRRELEVANT;
2301  lpdata->rrowsindexofleftrow[r] = lhsisinfinity ? NONEXISTENT_ROW : NO_RELEVANT_COLUMNS;
2302  lpdata->rrowsindexofrightrow[r] = rhsisinfinity ? NONEXISTENT_ROW : NO_RELEVANT_COLUMNS;
2303  }
2304  }
2305  else /* case: lhsslack > maxslack && rhsslack > maxslack */
2306  {
2307  lpdata->subproblemsindexofrow[r] = IRRELEVANT;
2308  lpdata->rrowsindexofleftrow[r] = lhsisinfinity ? NONEXISTENT_ROW : SLACK_GREATER_THAN_MAXSLACK;
2309  lpdata->rrowsindexofrightrow[r] = rhsisinfinity ? NONEXISTENT_ROW : SLACK_GREATER_THAN_MAXSLACK;
2310  }
2311  }
2312
2313  /* free temporary memory */
2314  SCIPfreeBufferArray(scip, &densecoeffscurrentleftrow);
2315  SCIPfreeBufferArray(scip, &densecoeffscurrentrightrow);
2316
2317  return SCIP_OKAY;
2318 }
2319
2320
2321
2322 /* check if mod 2 data structure contains at most two nonzero entries per row */
2323 static
2325  ZEROHALF_MOD2DATA* mod2data /**< considered mod 2 data */
2326  )
2327 {
2328  int r;
2329  int c;
2330  int nentries;
2331
2332  assert(mod2data != NULL);
2334  if( mod2data->nrowsind == 0 )
2335  return TRUE; /*FALSE;*/
2336
2337  if( mod2data->ncolsind <= 2 )
2338  return TRUE;
2339
2340  for( r = 0; r < mod2data->nrowsind ; ++r )
2341  {
2342  nentries = 0;
2343  for( c = 0; c < mod2data->ncolsind ; ++c )
2344  {
2345  if( BITARRAYBITISSET(mod2data->rows[mod2data->rowsind[r]], mod2data->colsind[c]) ) /*lint !e701*/
2346  {
2347  nentries++;
2348  if( nentries > 2 )
2349  return FALSE;
2350  }
2351  }
2352  }
2353
2354  return TRUE;
2355 }
2356
2357
2358
2359 #ifdef ZEROHALF__PRINT_STATISTICS
2360 /* check if mod 2 data structure contains at most two nonzero entries per column */
2361 static
2362 SCIP_Bool hasMatrixMax2EntriesPerColumn(
2363  ZEROHALF_MOD2DATA* mod2data /**< considered mod 2 data */
2364  )
2365 {
2366  int r;
2367  int c;
2368  int nentries;
2369
2370  assert(mod2data != NULL);
2371
2372  if( mod2data->ncolsind == 0 )
2373  return TRUE; /*FALSE;*/
2374
2375  if( mod2data->nrowsind <= 2 )
2376  return TRUE;
2377
2378  for( c = 0; c < mod2data->ncolsind ; ++c )
2379  {
2380  nentries = 0;
2381  for( r = 0; r < mod2data->nrowsind ; ++r )
2382  {
2383  if( BITARRAYBITISSET(mod2data->rows[mod2data->rowsind[r]], mod2data->colsind[c]) ) /*lint !e701*/
2384  {
2385  nentries++;
2386  if( nentries > 2 )
2387  return FALSE;
2388  }
2389  }
2390  }
2391
2392  return TRUE;
2393 }
2394 #endif
2395
2396
2397
2398 /** stores relevant data into bit arrays (mod 2 data structure) */
2399 static
2401  SCIP* scip, /**< SCIP data structure */
2403  ZEROHALF_LPDATA* lpdata, /**< data of current LP relaxation */
2404  int subproblemindex, /**< index of considered subproblem */
2405  ZEROHALF_MOD2DATA* mod2data /**< data (mod 2) */
2406  )
2407 {
2408  ZEROHALF_SUBLPDATA* problem;
2409  SCIP_COL** colscurrentrow;
2410  SCIP_ROW* row;
2411  SCIP_Real* nonzvalscurrentrow;
2412  SCIP_Real maxslack;
2413  SCIP_Real intscalar;
2414  BITARRAY tempcurrentrow;
2416  int nnonzcurrentrow;
2417  int rcolsindex;
2418  int c;
2419  int i;
2420  int j;
2421 #ifdef ZEROHALF__PRINT_STATISTICS
2422  int nirrelevantvarbounds;
2423 #endif
2424  SCIP_Bool tempmod2rhs;
2425  SCIP_Bool ignorerow;
2426  SCIP_Bool fliplhsrhs;
2427  SCIP_Bool isrhsrow;
2428  SCIP_Real* densecoeffscurrentrow;
2429
2430  assert(scip != NULL);
2432  assert(lpdata != NULL);
2433  assert(lpdata->rows != NULL);
2434  assert(lpdata->nrows > 0);
2435  assert(lpdata->cols != NULL);
2436  assert(lpdata->ncols > 0);
2437  assert(lpdata->subproblems != NULL);
2438  assert(lpdata->nsubproblems > 0);
2439  assert(lpdata->subproblemsindexofrow != NULL);
2440  assert(lpdata->rrowsindexofleftrow != NULL);
2441  assert(lpdata->rrowsindexofrightrow != NULL);
2442  assert(lpdata->subproblemsindexofcol != NULL);
2443  assert(lpdata->rcolsindexofcol != NULL);
2444  assert(0 <= subproblemindex);
2445  assert(subproblemindex <= lpdata->nsubproblems);
2446  problem = lpdata->subproblems[subproblemindex];
2447  assert(problem != NULL);
2448  assert(problem->rrows != NULL);
2449  assert(problem->nrrows > 0);
2450  assert(problem->rrowsrhs != NULL);
2451  assert(problem->rrowsslack != NULL);
2452  assert(problem->rcols != NULL);
2453  assert(problem->nrcols > 0);
2454  assert(problem->rcolslbslack != NULL);
2455  assert(problem->rcolsubslack != NULL);
2456  assert(mod2data != NULL);
2457
2458  /* identify varbounds to be added to the matrix */
2459  SCIP_CALL(SCIPallocBufferArray(scip, &varboundstoadd, 2 * problem->nrcols)); /* <0: lb, >0: ub */
2460
2462
2463 #ifdef ZEROHALF__PRINT_STATISTICS
2464  nirrelevantvarbounds = 0;
2465 #endif
2466  mod2data->nvarbounds = 0;
2467  for( c = 0 ; c < problem->nrcols ; c++ )
2468  {
2469  SCIP_Bool lbslackisok;
2470  SCIP_Bool ubslackisok;
2471
2472  lbslackisok = SCIPisLE(scip, problem->rcolslbslack[c], maxslack);
2473  ubslackisok = SCIPisLE(scip, problem->rcolsubslack[c], maxslack);
2474
2475  if( lbslackisok && ubslackisok )
2476  {
2477  SCIP_Real lb;
2478  SCIP_Real ub;
2479
2480  lb = SCIPcolGetLb(lpdata->cols[problem->rcols[c]]);
2481  ub = SCIPcolGetUb(lpdata->cols[problem->rcols[c]]);
2482
2483  if( ISEVEN(scip, lb) != ISEVEN(scip, ub) )
2484  {
2485  varboundstoadd[mod2data->nvarbounds] = (-1) * (c + 1);
2486  mod2data->nvarbounds++;
2487  varboundstoadd[mod2data->nvarbounds] = c + 1;
2488  mod2data->nvarbounds++;
2489  }
2490  else
2491  {
2492  if( SCIPisLE(scip, problem->rcolslbslack[c], problem->rcolsubslack[c]) )
2493  varboundstoadd[mod2data->nvarbounds] = (-1) * (c + 1);
2494  else
2495  varboundstoadd[mod2data->nvarbounds] = c + 1;
2496  mod2data->nvarbounds++;
2497 #ifdef ZEROHALF__PRINT_STATISTICS
2498  nirrelevantvarbounds++;
2499 #endif
2500  }
2501  }
2502  else
2503  {
2504  if( lbslackisok )
2505  {
2506  varboundstoadd[mod2data->nvarbounds] = (-1) * (c + 1);
2507  mod2data->nvarbounds++;
2508  }
2509 #ifdef ZEROHALF__PRINT_STATISTICS
2510  else
2511  nirrelevantvarbounds++;
2512 #endif
2513  if( ubslackisok )
2514  {
2515  varboundstoadd[mod2data->nvarbounds] = c + 1;
2516  mod2data->nvarbounds++;
2517  }
2518 #ifdef ZEROHALF__PRINT_STATISTICS
2519  else
2520  nirrelevantvarbounds++;
2521 #endif
2522  }
2523  }
2524  mod2data->nrows = problem->nrrows + mod2data->nvarbounds;
2525  mod2data->nrcols = problem->nrcols;
2526
2527  /* allocate temporary memory */
2528  SCIP_CALL( SCIPallocBlockMemoryArray(scip, &(mod2data->rows), mod2data->nrows) );
2529  SCIP_CALL( SCIPallocBlockMemoryArray(scip, &(mod2data->rowaggregations), mod2data->nrows) );
2530  SCIP_CALL( SCIPallocBlockMemoryArray(scip, &(mod2data->rhs), mod2data->nrows) );
2531  SCIP_CALL( SCIPallocBlockMemoryArray(scip, &(mod2data->slacks), mod2data->nrows) );
2532  SCIP_CALL( SCIPallocBlockMemoryArray(scip, &(mod2data->fracsol), problem->nrcols) );
2533  SCIP_CALL( SCIPallocBlockMemoryArray(scip, &(mod2data->rowstatistics), mod2data->nrows) );
2534  SCIP_CALL( SCIPallocBlockMemoryArray(scip, &(mod2data->colstatistics), problem->nrcols) );
2535  SCIP_CALL( SCIPallocBlockMemoryArray(scip, &(mod2data->rowsind), mod2data->nrows) );
2536  SCIP_CALL( SCIPallocBlockMemoryArray(scip, &(mod2data->colsind), problem->nrcols) );
2537  SCIP_CALL( SCIPallocBufferArray(scip, &densecoeffscurrentrow, lpdata->ncols) );
2538
2539  /* initialize temporary memory */
2540  mod2data->relatedsubproblem = problem;
2541  BMSclearMemoryArray(mod2data->rows, mod2data->nrows); /* NULL = 0x0 */
2542  BMSclearMemoryArray(mod2data->rowaggregations, mod2data->nrows); /* NULL = 0x0 */
2543  BMScopyMemoryArray(mod2data->slacks, problem->rrowsslack, problem->nrrows);
2544  for( c = 0 ; c < problem->nrcols ; ++c)
2545  mod2data->fracsol[c] = SCIPcolGetPrimsol(lpdata->cols[problem->rcols[c]]);
2546  for( c = 0 ; c < problem->nrcols ; c++)
2547  mod2data->colsind[c] = c;
2548  mod2data->nrowsind = 0;
2549  mod2data->ncolsind = problem->nrcols;
2550  mod2data->rowsbitarraysize = (int) GETREQUIREDBITARRAYSIZE(problem->nrcols);
2551  mod2data->rowaggregationsbitarraysize = (int) GETREQUIREDBITARRAYSIZE(problem->nrrows);
2552  tempcurrentrow = NULL;
2553
2554  /* (i) for all relevant rows */
2555  for( i = 0; i < problem->nrrows; ++i )
2556  {
2557  row = lpdata->rows[problem->rrows[i]];
2558  colscurrentrow = SCIProwGetCols(row);
2559  nonzvalscurrentrow = SCIProwGetVals(row);
2560  nnonzcurrentrow = SCIProwGetNLPNonz(row);
2561  assert(nnonzcurrentrow > 0);
2562  tempcurrentrow = NULL;
2563  fliplhsrhs = FALSE;
2564  ignorerow = FALSE;
2565
2566  /* check if rrows corresponds to a lhs or rhs row in the LP */
2567  if( lpdata->rrowsindexofleftrow[problem->rrows[i]] == i )
2568  isrhsrow = FALSE;
2569  else
2570  {
2571  assert(lpdata->rrowsindexofrightrow[problem->rrows[i]] == i);
2572  isrhsrow = TRUE;
2573  }
2574  intscalar = isrhsrow ? lpdata->intscalarsrightrow[problem->rrows[i]]
2575  : lpdata->intscalarsleftrow[problem->rrows[i]];
2576
2577  /* clear dense coeffs array */
2578  BMSclearMemoryArray(densecoeffscurrentrow, lpdata->ncols);
2579
2580  /* compute dense coeffs array of current row (including intscaling and bound substitutions) */
2581  for( j = 0 ; j < nnonzcurrentrow; ++j )
2582  {
2583  if( SCIPvarGetType(SCIPcolGetVar(colscurrentrow[j])) == SCIP_VARTYPE_CONTINUOUS )
2584  {
2585  SCIP_Bool ispositivecoeff;
2586  SCIP_VAR* bestzvbnd;
2587  SCIP_Real bestbndsol;
2588  SCIP_Real bestbvbnd;
2589  SCIP_Real bestdvbnd;
2590  int bestbndtype;
2591
2592  /* check sign of coefficient */
2593  if( nonzvalscurrentrow[j] * intscalar > 0.0 )
2594  ispositivecoeff = TRUE;
2595  else
2596  ispositivecoeff = FALSE;
2597
2598  /* get appropriate bound */
2599  if( isrhsrow == ispositivecoeff )
2600  findClosestLb(scip, lpdata, colscurrentrow[j], &bestbndsol, &bestbndtype, &bestzvbnd, &bestbvbnd, &bestdvbnd );
2601  else
2602  findClosestUb(scip, lpdata, colscurrentrow[j], &bestbndsol, &bestbndtype, &bestzvbnd, &bestbvbnd, &bestdvbnd );
2603
2604  /* check bound type */
2605  assert(bestbndtype > -2);
2606  if( !USEVARBOUNDS ) /*lint !e774 !e506*/
2607  assert(bestbndtype == -1);
2608
2609  /* normal lb or ub is used; only rhs would have to be adjusted but this has already been done in getRelevantRows */
2610  if( bestbndtype == -1 )
2611  continue;
2612  assert(USEVARBOUNDS && bestbndtype > -1); /*lint !e774 !e506*/
2613
2614  /* variable bound is used: update coefficient of non-continuous variable z that is used in substitution */
2615  densecoeffscurrentrow[SCIPcolGetLPPos(SCIPvarGetCol(bestzvbnd))] += (nonzvalscurrentrow[j] * intscalar * bestbvbnd);
2616  }
2617  else
2618  {
2619  densecoeffscurrentrow[SCIPcolGetLPPos(colscurrentrow[j])] += (nonzvalscurrentrow[j] * intscalar);
2620  }
2621  }
2622
2623  for( j = 0 ; j < lpdata->ncols; ++j )
2624  {
2625  assert(SCIPcolGetLPPos(lpdata->cols[j]) == j);
2626
2627  if( SCIPisZero(scip, densecoeffscurrentrow[j]) )
2628  continue;
2629
2630  if( intscalar == 1.0 && !SCIPisIntegral(scip, densecoeffscurrentrow[j]) )
2631  {
2632  ignorerow = TRUE;
2633  break;
2634  }
2635  else
2637
2638  /* integral coefficient */
2639  /* coefficient is only integral with respect to tolerances; use really integral values */
2640  if( isrhsrow )
2641  densecoeffscurrentrow[j] = SCIPfloor(scip, densecoeffscurrentrow[j]);
2642  else
2643  densecoeffscurrentrow[j] = SCIPceil(scip, densecoeffscurrentrow[j]);
2644
2645  if( ISODD(scip, densecoeffscurrentrow[j]) )
2646  {
2647  rcolsindex = lpdata->rcolsindexofcol[j];
2648  fliplhsrhs = XOR((int) fliplhsrhs,
2649  (int) (rcolsindex == LP_SOL_EQUALS_ODD_LB || rcolsindex == LP_SOL_EQUALS_ODD_UB));
2650  if( rcolsindex >= 0 ) /* relevant column? */
2651  {
2652  if( tempcurrentrow == NULL )
2653  {
2654  SCIP_CALL( SCIPallocBlockMemoryArray(scip, &tempcurrentrow, mod2data->rowsbitarraysize) );
2655  BITARRAYCLEAR(tempcurrentrow, mod2data->rowsbitarraysize);
2656  }
2657  assert(rcolsindex < problem->nrcols);
2658  BITARRAYBITSET(tempcurrentrow, rcolsindex); /*lint !e701*/
2659  assert(BITARRAYBITISSET(tempcurrentrow, rcolsindex)); /*lint !e701*/
2660  }
2661  }
2662  }
2663
2664  /* check if current row should be ignored, continuing with the next one */
2665  if( ignorerow )
2666  {
2667  SCIPfreeBlockMemoryArrayNull(scip, &tempcurrentrow, mod2data->rowsbitarraysize);
2668  continue;
2669  }
2670
2671  /* consider rhs */
2672  if( XOR((int) ISODD(scip, problem->rrowsrhs[i]), (int) fliplhsrhs) )
2673  tempmod2rhs = TRUE;
2674  else
2675  tempmod2rhs = FALSE;
2676
2677  if( tempcurrentrow == NULL && tempmod2rhs )
2678  {
2679  SCIP_CALL( SCIPallocBlockMemoryArray(scip, &tempcurrentrow, mod2data->rowsbitarraysize) );
2680  BITARRAYCLEAR(tempcurrentrow, mod2data->rowsbitarraysize);
2681  }
2682  assert(tempcurrentrow != NULL || !tempmod2rhs);
2683
2684  /* store temporary data in appropriate (mod 2) data structures */
2685  if( tempcurrentrow != NULL )
2686  {
2687  mod2data->rows[i] = tempcurrentrow;
2688  mod2data->rhs[i] = tempmod2rhs;
2689
2690  assert(mod2data->rowaggregationsbitarraysize > 0);
2691  SCIP_CALL( SCIPallocBlockMemoryArray(scip, &(mod2data->rowaggregations[i]), mod2data->rowaggregationsbitarraysize) ); /*lint !e866*/
2692  BITARRAYCLEAR(mod2data->rowaggregations[i], mod2data->rowaggregationsbitarraysize); /*lint !e866*/
2693  BITARRAYBITSET(mod2data->rowaggregations[i], i); /*lint !e701*/
2694
2695  mod2data->rowsind[mod2data->nrowsind] = i;
2696  mod2data->nrowsind++;
2697
2698  tempcurrentrow = NULL;
2699  }
2700  else
2701  {
2702  /* zero row */
2703  lpdata->subproblemsindexofrow[problem->rrows[i]] = IRRELEVANT;
2704  if( lpdata->rrowsindexofleftrow[problem->rrows[i]] >= 0 )
2705  lpdata->rrowsindexofleftrow[problem->rrows[i]] = ZERO_ROW;
2706  if( lpdata->rrowsindexofrightrow[problem->rrows[i]] >= 0 )
2707  lpdata->rrowsindexofrightrow[problem->rrows[i]] = ZERO_ROW;
2708  }
2709  }
2710
2711
2712  /* (ii) for all relevant varbounds */
2713  i = problem->nrrows;
2714  for( j = 0 ; j < mod2data->nvarbounds ; ++j)
2715  {
2716  SCIP_Real bound;
2717
2718  if( varboundstoadd[j] < 0 )
2719  c = (-1) * varboundstoadd[j] - 1;
2720  else
2721  c = varboundstoadd[j] - 1;
2722
2723  assert(mod2data->rowsbitarraysize > 0);
2724  SCIP_CALL(SCIPallocBlockMemoryArray(scip, &(mod2data->rows[i]), mod2data->rowsbitarraysize)); /*lint !e866*/
2725  BITARRAYCLEAR(mod2data->rows[i], mod2data->rowsbitarraysize); /*lint !e866*/
2726  BITARRAYBITSET(mod2data->rows[i], c); /*lint !e701*/
2727  assert(BITARRAYBITISSET(mod2data->rows[i], c)); /*lint !e701*/
2728
2729  SCIP_CALL(SCIPallocBlockMemoryArray(scip, &(mod2data->rowaggregations[i]), mod2data->rowaggregationsbitarraysize)); /*lint !e866*/
2730  BITARRAYCLEAR(mod2data->rowaggregations[i], mod2data->rowaggregationsbitarraysize); /*lint !e866*/
2731
2732  if( varboundstoadd[j] < 0 )
2733  {
2734  bound = SCIPcolGetLb(lpdata->cols[problem->rcols[c]]);
2735  mod2data->rhs[i] = ISODD(scip, bound);
2736  mod2data->slacks[i] = problem->rcolslbslack[c];
2737  }
2738  else
2739  {
2740  bound = SCIPcolGetUb(lpdata->cols[problem->rcols[c]]);
2741  mod2data->rhs[i] = ISODD(scip, bound);
2742  mod2data->slacks[i] = problem->rcolsubslack[c];
2743  }
2744  if( SCIPisFeasZero(scip, mod2data->slacks[i]) )
2745  mod2data->slacks[i] = 0.0;
2746
2747  mod2data->rowsind[mod2data->nrowsind] = i;
2748  mod2data->nrowsind++;
2749  i++;
2750  }
2751
2752  /* free temporary memory */
2753  SCIPfreeBufferArray(scip, &densecoeffscurrentrow);
2755
2756 #ifdef ZEROHALF__PRINT_STATISTICS
2757  ZEROHALFstatisticsMessage("\n");
2758  ZEROHALFstatisticsMessage(" | ------------------------------- subproblem ------------------------------- | ------------------------------\n");
2759  ZEROHALFstatisticsMessage(" | nrrows | nrcols | nvarbnds | ndlvbnds | max2/row | max2/col | A^T ept | \n");
2760  ZEROHALFstatisticsMessage("%15s | %8d | %8d | %8d | %8d | %8s | %8s | %8s |\n",
2761  "SUBPROBLEMDATA", problem->nrrows, problem->nrcols, mod2data->nvarbounds, nirrelevantvarbounds,
2762  hasMatrixMax2EntriesPerRow(mod2data) ? "yes" : "no", hasMatrixMax2EntriesPerColumn(mod2data) ? "yes" : "no", "n/a");
2763 #endif
2764
2765  return SCIP_OKAY;
2766 }
2767
2768
2769 /* --------------------------------------------------------------------------------------------------------------------
2770  * local methods: cut generation
2771  * -------------------------------------------------------------------------------------------------------------------- */
2772
2773
2774 /** stores nonzero elements of dense coefficient vector as sparse vector, and calculates activity and norm */
2775 static
2777  SCIP* scip, /**< SCIP data structure */
2778  int nvars, /**< number of problem variables */
2779  SCIP_VAR** vars, /**< problem variables */
2780  SCIP_Real* cutcoefs, /**< dense coefficient vector */
2781  SCIP_Real* varsolvals, /**< dense variable LP solution vector */
2782  char normtype, /**< type of norm to use for efficacy norm calculation */
2783
2784  SCIP_VAR** cutvars, /**< array to store variables of sparse cut vector */
2785  SCIP_Real* cutvals, /**< array to store coefficients of sparse cut vector */
2786  int* cutlen, /**< pointer to store number of nonzero entries in cut */
2787  SCIP_Real* cutact, /**< pointer to store activity of cut */
2788  SCIP_Real* cutnorm /**< pointer to store norm of cut vector */
2789  )
2790 {
2791  SCIP_Real val;
2792  SCIP_Real absval;
2793  SCIP_Real cutsqrnorm;
2794  SCIP_Real act;
2795  SCIP_Real norm;
2796  int len;
2797  int v;
2798
2799  assert(nvars == 0 || cutcoefs != NULL);
2800  assert(nvars == 0 || varsolvals != NULL);
2801  assert(cutvars != NULL);
2802  assert(cutvals != NULL);
2803  assert(cutlen != NULL);
2804  assert(cutact != NULL);
2805  assert(cutnorm != NULL);
2806
2807  len = 0;
2808  act = 0.0;
2809  norm = 0.0;
2810  switch(normtype)
2811  {
2812  case 'e':
2813  cutsqrnorm = 0.0;
2814  for( v = 0; v < nvars; ++v)
2815  {
2816  val = cutcoefs[v];
2817  if( !SCIPisZero(scip, val) )
2818  {
2819  act += val * varsolvals[v];
2820  cutsqrnorm += SQR(val);
2821  cutvars[len] = vars[v];
2822  cutvals[len] = val;
2823  len++;
2824  }
2825  }
2826  norm = SQRT(cutsqrnorm);
2827  break;
2828  case 'm':
2829  for( v = 0; v < nvars; ++v)
2830  {
2831  val = cutcoefs[v];
2832  if( !SCIPisZero(scip, val) )
2833  {
2834  act += val * varsolvals[v];
2835  absval = REALABS(val);
2836  norm = MAX(norm, absval);
2837  cutvars[len] = vars[v];
2838  cutvals[len] = val;
2839  len++;
2840  }
2841  }
2842  break;
2843  case 's':
2844  for( v = 0; v < nvars; ++v)
2845  {
2846  val = cutcoefs[v];
2847  if( !SCIPisZero(scip, val) )
2848  {
2849  act += val * varsolvals[v];
2850  norm += REALABS(val);
2851  cutvars[len] = vars[v];
2852  cutvals[len] = val;
2853  len++;
2854  }
2855  }
2856  break;
2857  case 'd':
2858  for( v = 0; v < nvars; ++v)
2859  {
2860  val = cutcoefs[v];
2861  if( !SCIPisZero(scip, val) )
2862  {
2863  act += val * varsolvals[v];
2864  norm = 1.0;
2865  cutvars[len] = vars[v];
2866  cutvals[len] = val;
2867  len++;
2868  }
2869  }
2870  break;
2871  default:
2872  SCIPerrorMessage("invalid efficacy norm parameter '%c'\n", normtype);
2873  return SCIP_INVALIDDATA;
2874  }
2875
2876  *cutlen = len;
2877  *cutact = act;
2878  *cutnorm = norm;
2879
2880  return SCIP_OKAY;
2881 }
2882
2883
2884
2885 /** adds a separated zerohalf cut to SCIP if it was successfully created and is efficacious */
2886 static
2888  SCIP* scip, /**< SCIP data structure */
2890  ZEROHALF_CUTDATA* cutdata, /**< separated zerohalf cut */
2891  int* nsepacuts, /**< pointer to store number of separated (efficacious) zerohalf cuts */
2892  SCIP_RESULT* result /**< pointer to store return code */
2893  )
2894 {
2895  assert(scip != NULL);
2897  assert(cutdata != NULL);
2898  assert(result != NULL);
2899
2900  /* check if SCIPcalcMIR was not successful */
2901  if( !cutdata->isfeasviolated || !cutdata->success )
2902  return SCIP_OKAY;
2903
2904  /* check if norm was not calculated correctly */
2905  if( !SCIPisPositive(scip, cutdata->norm) )
2906  {
2907  SCIPerrorMessage("Zerohalf cut norm is NOT positive!\n");
2908  return SCIP_ERROR;
2909  }
2910
2911  /* check if cut is not efficacious */
2913  && !SCIPisEfficacious(scip, cutdata->efficacy) )
2914  {
2915  return SCIP_OKAY;
2916  }
2917
2918  /* add cut (if no cutpool is used otherwise add it at the end of the separation main method)*/
2920  {
2921  SCIP_Bool cutoff;
2923  if ( cutoff )
2924  {
2925  *result = SCIP_CUTOFF;
2926  return SCIP_OKAY;
2927  }
2928  if( !cutdata->islocal )
2929  {
2931  }
2932  }
2933
2935  (*nsepacuts)++;
2936
2937  *result = SCIP_SEPARATED;
2938
2939  SCIPdebug( SCIP_CALL( SCIPprintRow(scip, cutdata->cut, NULL) ) );
2940
2941  return SCIP_OKAY;
2942 }
2943
2944
2945 /* --------------------------------------------------------------------------------------------------------------------
2946  * local methods: preprocessing
2947  * -------------------------------------------------------------------------------------------------------------------- */
2948
2949
2950
2951 /** marks a row as "removed" and stores why it has been removed using a flag */
2952 static
2953 void markRowAsRemoved(
2954  ZEROHALF_MOD2DATA* mod2data, /**< considered mod 2 data */
2955  int r, /**< mod2data->rows index of row that shall be removed */
2956  int flag /**< flag (cause of removal) */
2957  )
2958 {
2959  assert(mod2data != NULL);
2960  assert(mod2data->relatedsubproblem != NULL);
2961  assert(mod2data->nrowsind > 0);
2962  assert(r >= 0);
2963  assert(r < mod2data->nrowsind);
2964
2965  mod2data->rowstatistics[mod2data->rowsind[r]] = flag;
2966 }
2967
2968
2969
2970 /** marks a row as "removed" and stores why it has been removed using a flag. in addition it clears this column's mod 2 data */
2971 static
2973  ZEROHALF_MOD2DATA* mod2data, /**< considered mod 2 data */
2974  int c, /**< mod2data->relatedsubproblem->rcols index of column that shall be removed */
2975  int flag /**< flag (cause of removal) */
2976  )
2977 {
2978  int i;
2979  int rowsbind;
2982  assert(mod2data != NULL);
2983  assert(mod2data->relatedsubproblem != NULL);
2984  assert(mod2data->ncolsind > 0);
2985  assert(c >= 0);
2986  assert(c < mod2data->ncolsind);
2987
2988
2989  /* mark col */
2990  mod2data->colstatistics[mod2data->colsind[c]] = flag;
2991
2992  /* clear col */
2993  rowsbind = (int) GETBITARRAYINDEX(mod2data->colsind[c]);
2995  for( i = 0 ; i < mod2data->nrowsind ; ++i)
2997 }
2998
2999
3000
3001
3002 /** given a subset of mod 2 rows it returns a {0,1/2} weight vector used to
3003  combine the (original) LP rows. Note: original rows a stored as lhs <= a^Tx
3004  <= rhs by SCIP. Positive weights refer to "right half-rows" a^Tx <= rhs and
3005  negative weights to "left half-rows" -a^Tx <= -lhs */
3006 static
3008  SCIP* scip, /**< SCIP data structure */
3010  ZEROHALF_LPDATA* lpdata, /**< data of current LP relaxation */
3011  ZEROHALF_MOD2DATA* mod2data, /**< considered mod 2 data */
3012  BITARRAY rrowsincut, /**< subset of selected mod2data->rows */
3013  SCIP_Real** weights, /**< pointer to store the {-0.5,0,0.5} weights vector */
3014  int* nrowsincut /**< pointer to store the number of combined original LP rows */
3015  )
3016 { /*lint --e{438}*/
3017  ZEROHALF_SUBLPDATA* problem;
3018  int lppos;
3019  int i;
3020  int nnonz;
3021
3022  assert(scip != NULL);
3023  assert(lpdata != NULL);
3024  assert(lpdata->nrows > 0);
3025  assert(mod2data != NULL);
3026  assert(mod2data->relatedsubproblem != NULL);
3027  assert(mod2data->nrowsind > 0);
3028  assert(rrowsincut != NULL);
3029  assert(weights != NULL);
3030  assert(*weights == NULL);
3031  assert(nrowsincut != NULL);
3032
3033  /* allocate temporary memory */
3034  SCIP_CALL( SCIPallocBlockMemoryArray(scip, weights, lpdata->nrows) );
3035
3036  /* initialize */
3037  BMSclearMemoryArray(*weights, lpdata->nrows);
3038  problem = mod2data->relatedsubproblem;
3039
3040  /* determine row weights */
3041  *nrowsincut = 0;
3042  nnonz = 0;
3043  for( i = 0 ; i < problem->nrrows ; ++i)
3044  {
3045  lppos = problem->rrows[i];
3046  assert(0 <= lppos && lppos <= lpdata->nrows);
3047  if( BITARRAYBITISSET(rrowsincut, i) ) /*lint !e701*/
3048  {
3049  assert(lpdata->rrowsindexofleftrow[lppos] == i || lpdata->rrowsindexofrightrow[lppos] == i);
3050
3051  SCIPdebugMsg(scip, " %1s0.5 (int scaling: %16.4f / %16.4f) row[%d] %s\n",
3052  lpdata->rrowsindexofleftrow[lppos] == i ? "-" : "+",
3053  lpdata->intscalarsleftrow[lppos], lpdata->intscalarsrightrow[lppos],
3054  lppos, SCIProwGetName(lpdata->rows[lppos]));
3055
3056  if( lpdata->rrowsindexofleftrow[lppos] == i )
3057  (*weights)[lppos] = lpdata->intscalarsleftrow[lppos] * (-0.5);
3058  else
3059  (*weights)[lppos] = lpdata->intscalarsrightrow[lppos] * 0.5;
3060
3061  nnonz += SCIProwGetNLPNonz(lpdata->rows[lppos]);
3062  (*nrowsincut)++;
3063  }
3064  }
3065
3066  /* check if row aggregation might be too dense */
3067  if( nnonz >= 5 * sepadata->maxnnonz )
3068  {
3069  SCIPfreeBlockMemoryArray(scip, weights, lpdata->nrows);
3070  }
3071
3072  return SCIP_OKAY;
3073 }
3074
3075
3076 /** creates a zerohalf cut from a given weightvector */
3077 static
3079  SCIP* scip, /**< SCIP data structure */
3080  SCIP_SEPA* sepa, /**< separator */
3082  ZEROHALF_LPDATA* lpdata, /**< data of current LP relaxation */
3083  SCIP_Real* weights, /**< weightvector */
3084  char normtype, /**< SCIP normtype */
3085  int nzerohalfcuts, /**< number of zerohalf cuts (used for naming the cut) */
3086  SCIP_Real** varsolvals, /**< pointer to array of LP solution values of variables */
3087  ZEROHALF_CUTDATA* cutdata, /**< pointer to data structure used for storing the cut */
3088  SCIP_Bool* cutoff /**< whether a cutoff has been detected */
3089  )
3090 {
3091  SCIP_Real* cutcoefs;
3092  SCIP_VAR** cutvars;
3093  SCIP_Real* cutvals;
3094  char cutname[SCIP_MAXSTRLEN];
3095
3096  assert(scip != NULL);
3097  assert(lpdata != NULL);
3098  assert(lpdata->nvars > 0);
3099  assert(weights != NULL);
3100  assert(varsolvals != NULL);
3101  assert(cutdata != NULL);
3102  assert(cutdata->relatedsubproblem != NULL);
3103  assert(cutoff != NULL);
3104
3105  *cutoff = FALSE;
3106
3107  /* note: cutdata->relatedmod2data can be NULL if cut was determined
3108  * before mod 2 data structures were created */
3109
3110  /* allocate temporary memory */
3111  SCIP_CALL( SCIPallocBufferArray(scip, &cutcoefs, lpdata->nvars) );
3112
3113  /* calculate MIR */
3114  cutdata->success = FALSE;
3115  if( sepadata->maxtestdelta == 0 )
3116  {
3117  /* generate cut for delta = 1.0 */
3120  weights, -1.0, NULL, -1, -1, NULL, 1.0, NULL, NULL, cutcoefs, &(cutdata->rhs), &(cutdata->activity),
3121  &(cutdata->success), &(cutdata->islocal), &(cutdata->cutrank)) );
3122
3124  {
3127  weights, -1.0, NULL, -1, -1, NULL, -1.0, NULL, NULL, cutcoefs, &(cutdata->rhs), &(cutdata->activity),
3128  &(cutdata->success), &(cutdata->islocal), &(cutdata->cutrank)) );
3129  }
3130  }
3131  else
3132  {
3133  int ncuts;
3134  SCIP_Real bestdelta;
3135  SCIP_Bool bestdeltavalid;
3136
3137  ncuts = 0;
3138
3139  if( *varsolvals == NULL )
3140  {
3141  /* get the solution values for all active variables */
3142  SCIP_CALL( SCIPallocBlockMemoryArray(scip, varsolvals, lpdata->nvars) );
3143  SCIP_CALL( SCIPgetSolVals(scip, NULL, lpdata->nvars, lpdata->vars, *varsolvals) );
3144
3145 #ifndef NDEBUG
3146  /* because later when calling SCIPcutGenerationHeuristicCmir() varsolvals are used, it is needed that the
3147  * corresponding variables have the same order here and there, so we do the same checking and test that all
3148  * variables are ordered by their problem index
3149  */
3150  {
3151  int i;
3152  for(i = lpdata->nvars - 1; i >= 0; --i )
3153  assert(i == SCIPvarGetProbindex(lpdata->vars[i]));
3154  }
3155 #endif
3156  }
3157  assert(*varsolvals != NULL);
3158
3159  /* find best value of delta */
3160  SCIP_CALL( SCIPcutGenerationHeuristicCmir(scip, sepa, NULL, *varsolvals, sepadata->maxtestdelta, weights,
3161  -1.0, NULL, -1, -1, BOUNDSWITCH, USEVBDS, ALLOWLOCAL, FIXINTEGRALRHS, sepadata->maxnnonz, MAXWEIGHTRANGE,
3162  MINFRAC, MAXFRAC, sepadata->trynegscaling, TRUE, "zerohalf", cutoff, &ncuts, &bestdelta, &bestdeltavalid) );
3163  assert(ncuts == 0);
3164
3165  /* best delta corresponds to an efficient cut */
3166  if( !(*cutoff) && bestdeltavalid )
3167  {
3170  weights, -1.0, NULL, -1, -1, NULL, bestdelta, NULL, NULL, cutcoefs, &(cutdata->rhs), &(cutdata->activity),
3171  &(cutdata->success), &(cutdata->islocal), &(cutdata->cutrank)) );
3172  }
3173  }
3174  assert(ALLOWLOCAL || !cutdata->islocal);
3175
3176  cutdata->violation = cutdata->activity - cutdata->rhs;
3177
3178  /* if successful, convert dense cut into sparse row */
3179  if( !(*cutoff) && cutdata->success )
3180  {
3181  cutdata->isfeasviolated = SCIPisFeasGT(scip, cutdata->activity, cutdata->rhs);
3182  SCIPdebugMsg(scip, "Cut is %sfeasviolated: (act: %e, rhs: %e, viol: %e)\n",
3183  cutdata->isfeasviolated ? "" : "not ", cutdata->activity, cutdata->rhs, cutdata->violation);
3184
3185  if( cutdata->isfeasviolated )
3186  {
3187  if( *varsolvals == NULL )
3188  {
3189  /* get the solution values for all active variables */
3190  SCIP_CALL( SCIPallocBlockMemoryArray(scip, varsolvals, lpdata->nvars) );
3191  SCIP_CALL( SCIPgetSolVals(scip, NULL, lpdata->nvars, lpdata->vars, *varsolvals) );
3192  }
3193  assert(*varsolvals != NULL);
3194
3195  /* get temporary memory for storing the cut as sparse row */
3196  SCIP_CALL(SCIPallocBufferArray(scip, &cutvars, lpdata->nvars));
3197  SCIP_CALL(SCIPallocBufferArray(scip, &cutvals, lpdata->nvars));
3198
3199  /* store the cut as sparse row, calculate activity and norm of cut */
3200  SCIP_CALL(storeCutInArrays(scip, lpdata->nvars, lpdata->vars,
3201  cutcoefs, *varsolvals, normtype, cutvars, cutvals,
3202  &(cutdata->nnonz), &(cutdata->activity), &(cutdata->norm)));
3203
3204  /* check cut norm and efficacy */
3205  if( SCIPisPositive(scip, cutdata->norm) )
3206  {
3207  cutdata->efficacy = (cutdata->activity - cutdata->rhs) / cutdata->norm;
3208
3210  || (SCIPisEfficacious(scip, cutdata->efficacy) && cutdata->nnonz < sepadata->maxnnonz) )
3211  {
3212  /* create cut */
3213  (void) SCIPsnprintf(cutname, SCIP_MAXSTRLEN,"zerohalf%d_%d", SCIPgetNLPs(scip), nzerohalfcuts);
3214  SCIP_CALL(SCIPcreateEmptyRowSepa(scip, &(cutdata->cut), sepa, cutname, -SCIPinfinity(scip), cutdata->rhs,
3216  SCIP_CALL(SCIPaddVarsToRow(scip, cutdata->cut, cutdata->nnonz, cutvars, cutvals));
3217  /* set cut rank */
3218  SCIProwChgRank(cutdata->cut, cutdata->cutrank);
3219  }
3220  else
3221  cutdata->success = FALSE;
3222  }
3223  else
3224  cutdata->success = FALSE;
3225
3226  /* free temporary memory */
3227  SCIPfreeBufferArray(scip, &cutvals);
3228  SCIPfreeBufferArray(scip, &cutvars);
3229  }
3230  else
3231  cutdata->success = FALSE;
3232  }
3233
3234  /* free temporary memory */
3235  SCIPfreeBufferArray(scip, &cutcoefs);
3236
3237  return SCIP_OKAY;
3238 }
3239
3240
3241 /** searches for trivial zerohalf cuts, given as (0,..0) row with rhs=1 and slack <= maxslack */
3242 static
3244  SCIP* scip, /**< SCIP data structure */
3245  SCIP_SEPA* sepa, /**< separator */
3247  ZEROHALF_LPDATA* lpdata, /**< data of current LP relaxation */
3248  ZEROHALF_MOD2DATA* mod2data, /**< considered (preprocessed) subproblem mod 2 */
3249  int firstrowsind, /**< first mod2data->rows index to be considered */
3250  int lastrowsind, /**< last mod2data->rows index to be considered */
3251  char normtype, /**< SCIP normtype */
3252  int maxsepacuts, /**< maximal number of zerohalf cuts separated per separation round */
3253  int maxcuts, /**< maximal number of zerohalf cuts found per separation round (incl. ineff. cuts) */
3254  int* nsepacuts, /**< pointer to store current number of separated zerohalf cuts */
3255  int* nzerohalfcuts, /**< pointer to store current number of found zerohalf cuts */
3256  ZEROHALF_CUTDATA** zerohalfcuts, /**< pointer to store a found zerohalf cut */
3257  SCIP_Real** varsolvals, /**< dense variable LP solution vector */
3258  CUTSEPARATEDBY cutseparatedby, /**< flag */
3259  SCIP_RESULT* result /**< pointer to SCIP result value of separation */
3260  )
3261 {
3262  int r;
3263  int r2;
3264  int nrowsremoved;
3265  SCIP_Real maxslack;
3266  BITARRAY zerorow;
3267  SCIP_Bool* removerow;
3268  SCIP_Real* weights;
3269  int nrowsincut;
3270  SCIP_Bool cutoff = FALSE;
3271
3272  assert(scip != NULL);
3273  assert(lpdata != NULL);
3274  assert(mod2data != NULL);
3275  assert(mod2data->relatedsubproblem != NULL);
3276  assert(firstrowsind >= 0);
3277  assert(lastrowsind <= mod2data->nrowsind);
3278  assert(nsepacuts != NULL);
3279  assert(nzerohalfcuts != NULL);
3280  assert(zerohalfcuts != NULL);
3281  assert(*nsepacuts <= *nzerohalfcuts);
3282  assert(varsolvals != NULL);
3283  assert(result != NULL);
3284
3285  /* check if matrix or colind range is empty */
3286  if( mod2data->nrowsind == 0 || lastrowsind - firstrowsind <= 0 )
3287  return SCIP_OKAY;
3288
3289  /* allocate temporary memory */
3290  SCIP_CALL( SCIPallocBufferArray(scip, &removerow, lastrowsind - firstrowsind) );
3291  SCIP_CALL( SCIPallocBufferArray(scip, &zerorow, mod2data->rowsbitarraysize) );
3292
3293  /* initialize */
3294  BMSclearMemoryArray(zerorow, mod2data->rowsbitarraysize);
3295  BMSclearMemoryArray(removerow, lastrowsind - firstrowsind);
3297  nrowsremoved = 0;
3298
3299  /* check all rows */
3300  for( r = 0 ; r < lastrowsind - firstrowsind && *nsepacuts < maxsepacuts && *nzerohalfcuts < maxcuts; ++r )
3301  {
3302  if( mod2data->rhs[mod2data->rowsind[firstrowsind + r]] == TRUE )
3303  {
3304  if( SCIPisLE(scip, mod2data->slacks[mod2data->rowsind[firstrowsind + r]], maxslack ))
3305  {
3306  if( BITARRAYSAREEQUAL(mod2data->rows[mod2data->rowsind[firstrowsind + r]],
3307  zerorow, mod2data->rowsbitarraysize) ) /*lint !e647 check if row is (0 ... 0 , 1) */
3308  {
3309  /* a violated zerohalf cut has been found */
3310  weights = NULL;
3311  SCIP_CALL( getZerohalfWeightvectorFromSelectedRowsBitarray(scip, sepadata, lpdata, mod2data,
3312  mod2data->rowaggregations[mod2data->rowsind[firstrowsind + r]], &weights, &nrowsincut) );
3313  if( weights == NULL )
3314  continue;
3315  assert(nrowsincut > 0);
3316
3317  /* create zerohalf cut */
3318  SCIP_CALL( ZerohalfCutDataCreate(scip, &(zerohalfcuts[*nzerohalfcuts]),
3319  mod2data->relatedsubproblem, mod2data, 1, nrowsincut, cutseparatedby) );
3321  lpdata, weights, normtype, *nzerohalfcuts, varsolvals, zerohalfcuts[*nzerohalfcuts], &cutoff) );
3322
3323  if ( cutoff )
3324  *result = SCIP_CUTOFF;
3325  else
3326  {
3329  (*nzerohalfcuts)++;
3330  }
3331
3332  /* free temporary memory */
3333  SCIPfreeBlockMemoryArray(scip, &weights, lpdata->nrows);
3334
3335  if ( cutoff )
3336  break;
3337
3338  removerow[r] = TRUE;
3339  nrowsremoved++;
3340  }
3341  }
3342  }
3343  }
3344
3345  /* update mod2data->rowsind if necessary */
3346  if( ! cutoff && nrowsremoved > 0 )
3347  {
3348  r2 = firstrowsind;
3349  for( r = firstrowsind ; r < mod2data->nrowsind && r2 < mod2data->nrowsind; ++r)
3350  {
3351  if( r < lastrowsind - firstrowsind )
3352  while( r2 < mod2data->nrowsind && removerow[r2] )
3353  r2++;
3354  if( r < r2 && r2 < mod2data->nrowsind )
3355  mod2data->rowsind[r] = mod2data->rowsind[r2];
3356  r2++;
3357  }
3358  mod2data->nrowsind -= nrowsremoved;
3359  }
3360
3361
3362  /* free temporary memory */
3363  SCIPfreeBufferArray(scip, &zerorow);
3364  SCIPfreeBufferArray(scip, &removerow);
3365
3366  return SCIP_OKAY;
3367 }
3368
3369
3370 /** applies some row reductions */
3371 static
3373  SCIP* scip, /**< SCIP data structure */
3375  ZEROHALF_LPDATA* lpdata, /**< data of current LP relaxation */
3376  ZEROHALF_MOD2DATA* mod2data, /**< considered (preprocessed) subproblem mod 2 */
3377  int firstrowsind, /**< first mod2data->rows index to be considered */
3378  int lastrowsind, /**< last mod2data->rows index to be considered */
3379  SCIP_Bool removezerorows, /**< should zero rows be removed? */
3380  SCIP_Bool removelargeslackrows, /**< should rows with slack > maxslack be removed? */
3381  SCIP_Bool removeidenticalrows /**< should identical rows be removed? */
3382  )
3383 { /*lint --e{647}*/
3384  int r1;
3385  int r2;
3386  SCIP_Bool* rowisprocessed;
3387  SCIP_Bool* removerow;
3388  int nzerorowsremoved;
3389  int nlargeslackrowsremoved;
3390  int nidenticalrowsremoved;
3391  SCIP_Real maxslack;
3392  BITARRAY zerorow;
3393
3394  assert(scip != NULL);
3395  assert(lpdata != NULL);
3396  assert(mod2data != NULL);
3397  assert(mod2data->relatedsubproblem != NULL);
3398  assert(firstrowsind >= 0);
3399  assert(lastrowsind <= mod2data->nrowsind);
3400
3401
3402  /* check if matrix or colind range is empty */
3403  if( mod2data->nrowsind == 0 || lastrowsind - firstrowsind <= 0 )
3404  return SCIP_OKAY;
3405
3406  /* allocate temporary memory */
3407  SCIP_CALL(SCIPallocBufferArray(scip, &rowisprocessed, lastrowsind - firstrowsind));
3408  SCIP_CALL(SCIPallocBufferArray(scip, &removerow, lastrowsind - firstrowsind));
3409  zerorow = NULL;
3410  if( removezerorows )
3411  {
3412  SCIP_CALL(SCIPallocBufferArray(scip, &zerorow, mod2data->rowsbitarraysize));
3413  }
3414  /* initialize */
3415  BMSclearMemoryArray(rowisprocessed, lastrowsind - firstrowsind);
3416  BMSclearMemoryArray(removerow, lastrowsind - firstrowsind);
3417  if( removezerorows )
3418  {
3419  BMSclearMemoryArray(zerorow, mod2data->rowsbitarraysize);
3420  }
3422  nzerorowsremoved = 0;
3423  nlargeslackrowsremoved = 0;
3424  nidenticalrowsremoved = 0;
3425
3426  /* check all pairs of rows */
3427  for( r1 = 0 ; r1 < lastrowsind - firstrowsind ; ++r1)
3428  if( !rowisprocessed[r1] )
3429  {
3430  rowisprocessed[r1] = TRUE;
3431
3432  if( removezerorows && !removerow[r1] )
3433  if( mod2data->rhs[mod2data->rowsind[firstrowsind + r1]] == FALSE )
3434  {
3435  assert(zerorow != NULL);
3436  if( BITARRAYSAREEQUAL(mod2data->rows[mod2data->rowsind[firstrowsind + r1]],
3437  zerorow, mod2data->rowsbitarraysize) )
3438  {
3439  markRowAsRemoved(mod2data, firstrowsind + r1, ZERO_ROW);
3440  removerow[r1] = TRUE;
3441  nzerorowsremoved++;
3442  }
3443  }
3444  if( removelargeslackrows && !removerow[r1] )
3445  if( SCIPisGT(scip, mod2data->slacks[mod2data->rowsind[firstrowsind + r1]], maxslack) )
3446  {
3447  markRowAsRemoved(mod2data, firstrowsind + r1, SLACK_GREATER_THAN_MAXSLACK);
3448  removerow[r1] = TRUE;
3449  nlargeslackrowsremoved++;
3450  }
3451
3452  if( removeidenticalrows && !removerow[r1] )
3453  for( r2 = r1 + 1 ; r2 < lastrowsind - firstrowsind ; ++r2)
3454  if( !rowisprocessed[r2] )
3455  if( mod2data->rhs[mod2data->rowsind[firstrowsind + r1]]
3456  == mod2data->rhs[mod2data->rowsind[firstrowsind + r2]] )
3457  if( BITARRAYSAREEQUAL(mod2data->rows[mod2data->rowsind[firstrowsind + r1]],
3458  mod2data->rows[mod2data->rowsind[firstrowsind + r2]], mod2data->rowsbitarraysize) )
3459  {
3460  if( SCIPisLT(scip, mod2data->slacks[mod2data->rowsind[firstrowsind + r1]],
3461  mod2data->slacks[mod2data->rowsind[firstrowsind + r2]]) )
3462  {
3463  markRowAsRemoved(mod2data, firstrowsind + r2, IDENT_TO_ROW_WITH_SMALLER_SLACK);
3464  removerow[r2] = TRUE;
3465  nidenticalrowsremoved++;
3466  rowisprocessed[r2] = TRUE;
3467  }
3468  else
3469  {
3470  markRowAsRemoved(mod2data, firstrowsind + r1, IDENT_TO_ROW_WITH_SMALLER_SLACK);
3471  removerow[r1] = TRUE;
3472  nidenticalrowsremoved++;
3473  break;
3474  }
3475  }
3476  }
3477
3478  /* update mod2data->rowsind if necessary */
3479  if( nzerorowsremoved + nlargeslackrowsremoved + nidenticalrowsremoved > 0 )
3480  {
3481  r2 = firstrowsind;
3482  for( r1 = firstrowsind ; r1 < mod2data->nrowsind && r2 < mod2data->nrowsind; ++r1)
3483  {
3484  if( r1 < lastrowsind - firstrowsind )
3485  while( r2 < mod2data->nrowsind && removerow[r2] )
3486  r2++;
3487  if( r1 < r2 && r2 < mod2data->nrowsind )
3488  mod2data->rowsind[r1] = mod2data->rowsind[r2];
3489  r2++;
3490  }
3491  mod2data->nrowsind -= (nzerorowsremoved + nlargeslackrowsremoved + nidenticalrowsremoved);
3492  }
3493
3494
3495  /* free temporary memory */
3496  if( removezerorows && zerorow != NULL )
3497  {
3498  SCIPfreeBufferArray(scip, &zerorow);
3499  }
3500  SCIPfreeBufferArray(scip, &removerow);
3501  SCIPfreeBufferArray(scip, &rowisprocessed);
3502
3503  return SCIP_OKAY;
3504 }
3505
3506
3507 /** applies some column reductions */
3508 static
3510  SCIP* scip, /**< SCIP data structure */
3512  ZEROHALF_LPDATA* lpdata, /**< data of current LP relaxation */
3513  ZEROHALF_MOD2DATA* mod2data, /**< considered (preprocessed) subproblem mod 2 */
3514  int firstcolsind, /**< first mod2data->rows index to be considered */
3515  int lastcolsind, /**< last mod2data->rows index to be considered */
3516  SCIP_Bool removezerocols, /**< should zero columns be removed? */
3517  SCIP_Bool removecolsingletons,/**< should column singletons be removed? */
3518  SCIP_Bool checkresultingrows /**< should rows whose slack becomes larger than maxslack be removed? */
3519  )
3520 {
3521  SCIP_Real maxslack;
3522  int maxnnonzentries;
3523  int nzerocolsremoved;
3524  int ncolsingletonsremoved;
3525  int nunprocessedcols;
3526  int nconsideredcols;
3527  int nnonzentries;
3528  SCIP_Bool* removecol;
3529  SCIP_Bool* colisprocessed;
3530  int rowofcolsingleton;
3531  int rowsbind;
3533  int r;
3534  int c;
3535  int j;
3536
3537  assert(scip != NULL);
3539  assert(lpdata != NULL);
3540  assert(mod2data != NULL);
3541  assert(mod2data->relatedsubproblem != NULL);
3542  assert(firstcolsind >= 0);
3543  assert(lastcolsind <= mod2data->ncolsind);
3544  assert(removezerocols || removecolsingletons);
3545
3546
3547  nconsideredcols = lastcolsind - firstcolsind;
3548
3549  /* check if matrix or colind range is empty */
3550  if( mod2data->ncolsind == 0 || mod2data->nrowsind == 0 || nconsideredcols <= 0 )
3551  return SCIP_OKAY;
3552
3553
3554  /* allocate temporary memory */
3555  SCIP_CALL(SCIPallocBufferArray(scip, &colisprocessed, nconsideredcols));
3556  SCIP_CALL(SCIPallocBufferArray(scip, &removecol, nconsideredcols));
3557
3558  /* initialize */
3559  BMSclearMemoryArray(colisprocessed, nconsideredcols);
3560  BMSclearMemoryArray(removecol, nconsideredcols);
3562  nunprocessedcols = nconsideredcols;
3563  nzerocolsremoved = 0;
3564  ncolsingletonsremoved = 0;
3565  rowofcolsingleton = -1;
3566  if( removecolsingletons )
3567  maxnnonzentries = 1;
3568  else
3569  maxnnonzentries = 0;
3570
3571  /* check all columns if they contain exactly one nonzero entry */
3572  while(nunprocessedcols > 0)
3573  {
3574  for( c = 0 ; c < nconsideredcols ; ++c )
3575  if( colisprocessed[c] == FALSE )
3576  break;
3577  assert(firstcolsind + c < mod2data->ncolsind);
3578
3579  nnonzentries = 0;
3580  rowsbind = (int) GETBITARRAYINDEX(mod2data->colsind[firstcolsind + c]);
3582
3583  for( r = 0 ; r < mod2data->nrowsind ; ++r)
3584  if( mod2data->rows[mod2data->rowsind[r]][rowsbind] & rowsbmask )
3585  {
3586  nnonzentries++;
3587  if( nnonzentries > maxnnonzentries )
3588  break;
3589  rowofcolsingleton = r;
3590  }
3591
3592  /* check if a zero column has been found */
3593  if( removezerocols )
3594  if( nnonzentries == 0 )
3595  {
3596  /* remove zero columns */
3597  removecol[c] = TRUE;
3598  nzerocolsremoved++;
3599  markColAsRemovedAndClearCol(mod2data, firstcolsind + c, ZERO_COLUMN);
3600  }
3601
3602  /* check if a column singleton has been found */
3603  if( removecolsingletons && !removecol[c] )
3604  if( nnonzentries == 1 )
3605  {
3606  r = rowofcolsingleton;
3607  removecol[c] = TRUE;
3608
3609  /* update row slack: slack' = slack + fracsol */
3610  mod2data->slacks[mod2data->rowsind[r]] += mod2data->fracsol[mod2data->colsind[firstcolsind + c]];
3611
3612  /* if removing col results in a row with slack > maxslack,
3613  * then the row can be removed as well */
3614  if( checkresultingrows && SCIPisGT(scip, mod2data->slacks[mod2data->rowsind[r]], maxslack) )
3615  {
3617  for( j = 0 ; j < nconsideredcols ; ++j)
3618  if( !removecol[j] && colisprocessed[j] )
3619  if( BITARRAYBITISSET(mod2data->rows[mod2data->rowsind[r]],
3620  mod2data->colsind[firstcolsind + j]) ) /*lint !e701*/
3621  {
3622  colisprocessed[j] = FALSE; /* re-consider col */
3623  nunprocessedcols++;
3624  }
3625
3626  BMSmoveMemoryArray(&((mod2data->rowsind)[r]), &((mod2data->rowsind)[r + 1]),
3627  mod2data->nrowsind - r - 1); /*lint !e866*/
3628
3629  mod2data->nrowsind--;
3630  }
3631
3632  /* remove column singleton */
3633  ncolsingletonsremoved++;
3634  markColAsRemovedAndClearCol(mod2data, firstcolsind + c, SINGLETON_COLUMN);
3635  }
3636
3637  colisprocessed[c] = TRUE;
3638  nunprocessedcols--;
3639
3640  if( nzerocolsremoved + ncolsingletonsremoved == nconsideredcols || mod2data->nrowsind == 0 )
3641  break;
3642  }
3643
3644  /* if all rows have been deleted, remove cols as well */
3645  if( mod2data->nrowsind == 0 )
3646  {
3647  for( c = firstcolsind ; c < mod2data->ncolsind; ++c)
3648  if( !removecol[c] )
3649  {
3650  removecol[c] = TRUE;
3651  markColAsRemovedAndClearCol(mod2data, firstcolsind + c, ZERO_COLUMN);
3652  }
3653  nzerocolsremoved = nconsideredcols - ncolsingletonsremoved;
3654  assert(nzerocolsremoved + ncolsingletonsremoved == nconsideredcols);
3655  }
3656
3657  /* update mod2data->colsind array if necessary*/
3658  if( mod2data->nrowsind == 0 )
3659  mod2data->ncolsind = 0;
3660  else
3661  if( nzerocolsremoved + ncolsingletonsremoved > 0 )
3662  {
3663  j = firstcolsind;
3664  for( c = firstcolsind ; c < mod2data->ncolsind && j < mod2data->ncolsind; ++c)
3665  {
3666  if( c < nconsideredcols )
3667  while( j < mod2data->ncolsind && removecol[j] )
3668  j++;
3669  if( c < j && j < mod2data->ncolsind )
3670  mod2data->colsind[c] = mod2data->colsind[j];
3671  j++;
3672  }
3673  mod2data->ncolsind -= (nzerocolsremoved + ncolsingletonsremoved);
3674  }
3675
3676  /* free temporary memory */
3677  SCIPfreeBufferArray(scip, &removecol);
3678  SCIPfreeBufferArray(scip, &colisprocessed);
3679
3680  return SCIP_OKAY;
3681 }
3682
3683
3684 /** applies modified Gaussian Elimination reduction */
3685 static
3687  SCIP* scip, /**< SCIP data structure */
3689  ZEROHALF_LPDATA* lpdata, /**< data of current LP relaxation */
3690  ZEROHALF_MOD2DATA* mod2data /**< considered (preprocessed) subproblem mod 2 */
3691  )
3692 {
3693  int nslackzerorows;
3694  int pivotrow;
3695  int pivotcol;
3696  int pivot;
3697  int identsubmatrixsize;
3698  int rowsbind;
3700  int r;
3701  int temp;
3702
3703  assert(scip != NULL);
3705  assert(lpdata != NULL);
3706  assert(mod2data != NULL);
3707  assert(mod2data->relatedsubproblem != NULL);
3708
3709  /* check if matrix or colind range is empty */
3710  if( mod2data->ncolsind == 0 || mod2data->nrowsind == 0 )
3711  return SCIP_OKAY;
3712
3713  /* determine number of slack zero rows */
3714  nslackzerorows = 0;
3715  while( nslackzerorows < mod2data->nrowsind
3716  && SCIPisZero(scip, mod2data->slacks[mod2data->rowsind[nslackzerorows]]) )
3717  nslackzerorows++;
3718  /* check if at least one slack zero row exists */
3719  if( nslackzerorows == 0 )
3720  return SCIP_OKAY;
3721
3722
3723  /* sort column indices sets w.r.t. to their primsol values NON-INCREASINGLY */
3724  if( mod2data->ncolsind > 1 )
3725  {
3726  SCIPsortInd( mod2data->colsind , compRealNonIncreasing , (void*) mod2data->fracsol , mod2data->ncolsind );
3727  }
3728
3729  /* sort row indices sets w.r.t. to their slack values NON-DECREASINGLY */
3730  if( mod2data->nrowsind > 1 )
3731  {
3732  SCIPsortInd( mod2data->rowsind , compRealNonDecreasing , (void*) mod2data->slacks , mod2data->nrowsind );
3733  }
3734
3735
3736  identsubmatrixsize = 0;
3737
3738  /* create maximal identity submatrix */
3739  /* determine pivot col */
3740  for( pivotcol = 0 ; pivotcol < mod2data->ncolsind ; ++pivotcol)
3741  {
3742  if( identsubmatrixsize == mod2data->nrowsind )
3743  break;
3744
3745  /* determine pivot row */
3746  rowsbind = (int) GETBITARRAYINDEX(mod2data->colsind[pivotcol]);
3748  for( pivotrow = identsubmatrixsize ; pivotrow < nslackzerorows ; ++pivotrow)
3749  if( mod2data->rows[mod2data->rowsind[pivotrow]][rowsbind] & rowsbmask )
3750  break;
3751  if( pivotrow == nslackzerorows )
3752  continue;
3753
3754  /* Gaussian elimination step */
3755  for( r = 0 ; r < nslackzerorows ; ++r)
3756  {
3757  if( r == pivotrow )
3758  continue;
3759  if( mod2data->rows[mod2data->rowsind[r]][rowsbind] & rowsbmask )
3760  {
3761  /* add pivot row to r-th row */
3762  BITARRAYSXOR(mod2data->rows[mod2data->rowsind[pivotrow]],
3763  mod2data->rows[mod2data->rowsind[r]],mod2data->rowsbitarraysize);
3764  BITARRAYSXOR(mod2data->rowaggregations[mod2data->rowsind[pivotrow]],
3765  mod2data->rowaggregations[mod2data->rowsind[r]],mod2data->rowaggregationsbitarraysize);
3766  mod2data->rhs[mod2data->rowsind[r]] =
3767  XOR(mod2data->rhs[mod2data->rowsind[pivotrow]],mod2data->rhs[mod2data->rowsind[r]]);
3768  /* all rows have slack zero: */
3769  /* // mod2data->slacks[[mod2data->rowsind[r]] += mod2data->slacks[[mod2data->rowsind[pivotrow]]; */
3770  }
3771  }
3772
3773  /* swap index set positions */
3774  temp = mod2data->rowsind[pivotrow];
3775  mod2data->rowsind[pivotrow] = mod2data->rowsind[identsubmatrixsize];
3776  mod2data->rowsind[identsubmatrixsize] = temp;
3777  temp = mod2data->colsind[pivotcol];
3778  mod2data->colsind[pivotcol] = mod2data->colsind[identsubmatrixsize];
3779  mod2data->colsind[identsubmatrixsize] = temp;
3780
3781  identsubmatrixsize++;
3782  }
3783
3784  if( identsubmatrixsize > 0 )
3785  {
3786  /* add rows of identity submatrix properly to rows with positive slack
3787  * to transform each column of the identity submatrix into a column singleton */
3788  for( pivot = 0 ; pivot < identsubmatrixsize ; ++pivot)
3789  {
3790  rowsbind = (int) GETBITARRAYINDEX(mod2data->colsind[pivot]);
3792  for( r = nslackzerorows ; r < mod2data->nrowsind ; ++r)
3793  if( mod2data->rows[mod2data->rowsind[r]][rowsbind] & rowsbmask )
3794  {
3795  /* add pivot row to r-th row */
3796  BITARRAYSXOR(mod2data->rows[mod2data->rowsind[pivot]],
3797  mod2data->rows[mod2data->rowsind[r]],mod2data->rowsbitarraysize);
3798  BITARRAYSXOR(mod2data->rowaggregations[mod2data->rowsind[pivot]],
3799  mod2data->rowaggregations[mod2data->rowsind[r]],mod2data->rowaggregationsbitarraysize);
3800  mod2data->rhs[mod2data->rowsind[r]] =
3801  XOR(mod2data->rhs[mod2data->rowsind[pivot]],mod2data->rhs[mod2data->rowsind[r]]);
3802  /* all identity submatrix rows have slack zero */
3803  /* // mod2data->slacks[[mod2data->rowsind[r]] += mod2data->slacks[[mod2data->rowsind[pivot]]; */
3804  }
3805  }
3806
3807  /* remove generated column singletons */
3809  0, /*identsubmatrixsize*/ mod2data->ncolsind, FALSE, TRUE, TRUE));
3810
3811  /* remove zero rows */
3813  0, mod2data->nrowsind, TRUE, FALSE, FALSE));
3814  }
3815
3816  return SCIP_OKAY;
3817 }
3818
3819
3820 /** decomposes the problem into subproblems which can be considered separately */
3821 static
3823  SCIP* scip, /**< SCIP data structure */
3825  ZEROHALF_LPDATA* lpdata /**< data of current LP relaxation */
3826  )
3827 {
3828
3829 #ifdef WITHDECOMPOSE
3830  /**@todo this is buggy in different ways.
3831  * 1. it might happen that we ignore a variable of the current row and of all other rows.
3832  * thus at the end, the variable will not occur in any subproblem. BUT, currently we do not update
3833  * lpdata->subproblemsindexofcol[lppos] and lpdata->rcolsindexofcol[lppos] accordingly.
3834  * consequently, it might happen that lpdata->rcolsindexofcol[lppos] > problem->nrcols, with
3835  * with problem being the subproblem still associated to our column. therefore, a corresponding assert
3836  * assert(rcolsindex < problem->nrcols) in storeMod2Data() is violated
3837  * [e.g., for IP/atamtuerk/mik/unbounded/mik.250-1-50.3.mps.gz].
3838  * we could recognize whether a variable is never added to a subproblem and update its data structures,
3839  * but I'm not sure whether this will be correct, i.e., whether it is really ok to ignore some variables here.
3840  * 2. in particular, it seems like that this method has not been adapted to that we can deal with
3841  * continuous variables now. (see the below todo of Manuel)
3842  * 3. in case we end up with only one subproblem, we use the old problem. but in this case we do not update
3843  * the problem data and hence it is not consistent with the lpdata anymore where we might have set some
3844  * rows to be irrelevant.
3845  *
3846  * therefore, we will currently do nothing in here.
3847  */
3848  BITARRAY processedrows;
3849  int nprocessedrows;
3850  int processedrowsbitarraysize;
3851  int unprocessedrowidx;
3852
3853  BITARRAY processedcols;
3854  int nprocessedcols;
3855  int processedcolsbitarraysize;
3856
3857  int* queue;
3858  int queuefirst;
3859  int queuelast;
3860
3861  int i;
3862  int j;
3863  int k;
3864
3865  SCIP_COL** colsofrow;
3866  SCIP_ROW** rowsofcol;
3867
3868  SCIP_Real* colvals;
3869  SCIP_Real* rowvals;
3870  int ncolvals;
3871  int nrowvals;
3872  int cidx;
3873  int ridx;
3874  int lppos;
3875
3876  int rrowsidx;
3877  int rcolsidx;
3878
3879  SCIP_Bool fliplhsrhs;
3880
3881  ZEROHALF_SUBLPDATA* problem;
3882  int problemindex;
3883  ZEROHALF_SUBLPDATA* subproblem;
3884
3885  int* rrowsinsubprob;
3886  int* rcolsinsubprob;
3887  SCIP_Bool* rrowsinsubproboddrhs;
3888  int nrrowsinsubprob;
3889  int nrcolsinsubprob;
3890
3891  int totalnrrows;
3892  int totalnrcols;
3893
3894  int nrrowsinitial;
3895  int nrcolsinitial;
3896  int ndelvarbounds;
3897
3898  int rowindex;
3899  int colindex;
3900
3901  SCIP_Real maxslack;
3902
3903  assert(scip != NULL);
3905  assert(lpdata != NULL);
3906  assert(lpdata->subproblems != NULL);
3907  assert(lpdata->nsubproblems > 0);
3908
3909  problemindex = 0;
3910
3911  assert(problemindex >= 0);
3912  assert(problemindex <= lpdata->nsubproblems);
3913
3914  problem = lpdata->subproblems[problemindex];
3915
3916  assert(problem != NULL);
3917  assert(problem->rcols != NULL);
3918  assert(problem->nrcols > 0);
3919  assert(problem->rcolslbslack != NULL);
3920  assert(problem->rcolsubslack != NULL);
3921  assert(problem->rrows != NULL);
3922  assert(problem->nrrows > 0);
3923  assert(problem->rrowsrhs != NULL);
3924  assert(problem->rrowsslack != NULL);
3925
3926  if( sepadata->dtimer == NULL )
3927  {
3929  }
3931
3932  processedrowsbitarraysize = (int) GETREQUIREDBITARRAYSIZE(problem->nrrows);
3933  processedcolsbitarraysize = (int) GETREQUIREDBITARRAYSIZE(problem->nrcols);
3934
3935  SCIPfreeBlockMemoryArray(scip, &(lpdata->subproblems), lpdata->nsubproblems);
3936
3937  /* allocate temporary memory */
3938  SCIP_CALL( SCIPallocBlockMemoryArray(scip, &(lpdata->subproblems), problem->nrrows) );
3939  SCIP_CALL( SCIPallocBufferArray(scip, &processedrows, processedrowsbitarraysize) );
3940  SCIP_CALL( SCIPallocBufferArray(scip, &processedcols, processedcolsbitarraysize) );
3941  SCIP_CALL( SCIPallocBufferArray(scip, &queue, problem->nrrows) );
3942  SCIP_CALL( SCIPallocBlockMemoryArray(scip, &rrowsinsubprob, problem->nrrows)) ;
3943  SCIP_CALL( SCIPallocBlockMemoryArray(scip, &rrowsinsubproboddrhs, problem->nrrows) );
3944  SCIP_CALL( SCIPallocBlockMemoryArray(scip, &rcolsinsubprob, problem->nrcols) );
3945
3946  /* initialize temporary memory */
3947  BMSclearMemoryArray(processedrows, processedrowsbitarraysize);
3948  BMSclearMemoryArray(processedcols, processedcolsbitarraysize);
3949  lpdata->nsubproblems = 0;
3951
3952  nrrowsinitial = problem->nrrows;
3953  nrcolsinitial = problem->nrcols;
3954  totalnrrows = 0;
3955  totalnrcols = 0;
3956  ndelvarbounds = 0;
3957  nprocessedrows = 0;
3958  nprocessedcols = 0;
3959  k = 0;
3960  unprocessedrowidx = 0;
3961
3962  while( nprocessedrows < problem->nrrows )
3963  {
3964  ++k;
3965  nrrowsinsubprob = 0;
3966  nrcolsinsubprob = 0;
3967
3968  for( i = unprocessedrowidx ; i < problem->nrrows ; ++i)
3969  {
3970  if( BITARRAYBITISSET(processedrows, i) ) /*lint !e701*/
3971  unprocessedrowidx++;
3972  else
3973  break;
3974  }
3975  BITARRAYBITSET(processedrows, i); /*lint !e701*/
3976
3977  queue[0] = i;
3978  queuefirst = 0;
3979  queuelast = 1;
3980
3981  while( queuelast > queuefirst )
3982  {
3983  assert(queuelast <= problem->nrrows);
3984
3985  i = queue[queuefirst];
3986  queuefirst++;
3987
3988  rrowsinsubprob[nrrowsinsubprob] = i;
3989  nrrowsinsubprob++;
3990
3991  fliplhsrhs = FALSE;
3992
3993  colsofrow = SCIProwGetCols(lpdata->rows[problem->rrows[i]]);
3994  colvals = SCIProwGetVals(lpdata->rows[problem->rrows[i]]);
3995  ncolvals = SCIProwGetNLPNonz(lpdata->rows[problem->rrows[i]]);
3996
3997  for( cidx = 0 ; cidx < ncolvals ; ++cidx)
3998  {
3999  lppos = SCIPcolGetLPPos(colsofrow[cidx]);
4000  if( lppos == -1 )
4001  continue;
4002  rcolsidx = lpdata->rcolsindexofcol[lppos];
4003
4004  if( lpdata->subproblemsindexofcol[lppos] != problemindex || rcolsidx < 0 )
4005  {
4006  if( ISODD(scip, colvals[cidx]) )
4007  fliplhsrhs = XOR(fliplhsrhs,
4008  (SCIP_Bool) (lpdata->rcolsindexofcol[lppos] == LP_SOL_EQUALS_ODD_LB
4009  || lpdata->rcolsindexofcol[lppos] == LP_SOL_EQUALS_ODD_UB));
4010
4011  /**@todo analogue for continuous variables? */
4012
4013  continue; /* col is not relevant */
4014  }
4015  if( nprocessedcols == problem->nrcols )
4016  continue;
4017
4018  if( BITARRAYBITISSET(processedcols, rcolsidx) ) /*lint !e701*/
4019  continue;
4020
4021  if( ISEVEN(scip, colvals[cidx]) )
4022  continue;
4023
4024  rcolsinsubprob[nrcolsinsubprob] = rcolsidx;
4025  nrcolsinsubprob++;
4026  BITARRAYBITSET(processedcols, rcolsidx); /*lint !e701*/
4027
4028  rowsofcol = SCIPcolGetRows(colsofrow[cidx]);
4029  rowvals = SCIPcolGetVals(colsofrow[cidx]);
4030  nrowvals = SCIPcolGetNNonz(colsofrow[cidx]);
4031  for( ridx = 0 ; ridx < nrowvals ; ++ridx)
4032  {
4033  lppos = SCIProwGetLPPos(rowsofcol[ridx]);
4034  if( lppos == -1 )
4035  continue;
4036  rrowsidx = lpdata->rrowsindexofleftrow[lppos];
4037  if( lpdata->subproblemsindexofrow[lppos] == problemindex
4038  && rrowsidx >= 0 )
4039  {
4040  if( !BITARRAYBITISSET(processedrows, rrowsidx) ) /*lint !e701*/
4041  if( ISODD(scip, rowvals[ridx]) )
4042  {
4043  queue[queuelast] = rrowsidx;
4044  queuelast++;
4045  BITARRAYBITSET(processedrows, rrowsidx); /*lint !e701*/
4046  }
4047  }
4048  rrowsidx = lpdata->rrowsindexofrightrow[lppos];
4049  if( lpdata->subproblemsindexofrow[lppos] == problemindex
4050  && rrowsidx >= 0 )
4051  {
4052  if( !BITARRAYBITISSET(processedrows, rrowsidx) ) /*lint !e701*/
4053  if( ISODD(scip, rowvals[ridx]) )
4054  {
4055  queue[queuelast] = rrowsidx;
4056  queuelast++;
4057  BITARRAYBITSET(processedrows, rrowsidx); /*lint !e701*/
4058  }
4059  }
4060  }
4061  nprocessedcols++;
4062  }
4063  rrowsinsubproboddrhs[nrrowsinsubprob-1] = XOR((SCIP_Bool) ISODD(scip, problem->rrowsrhs[i]),fliplhsrhs);
4064  nprocessedrows++;
4065  }
4066
4067  /* a subproblem consisting only of rows with even rhs values can be ignored.
4068  * note: varbounds have to be considered! */
4069  for( i = 0 ; i < nrrowsinsubprob ; ++i)
4070  if( rrowsinsubproboddrhs[i] )
4071  break;
4072  for( j = 0; j < nrcolsinsubprob; ++j)
4073  {
4074  if( (SCIPisLE(scip, problem->rcolsubslack[rcolsinsubprob[j]], maxslack)
4075  && ISODD(scip, SCIPcolGetUb(lpdata->cols[problem->rcols[rcolsinsubprob[j]]])))
4076  || (SCIPisLE(scip, problem->rcolslbslack[rcolsinsubprob[j]], maxslack)
4077  && ISODD(scip, SCIPcolGetLb(lpdata->cols[problem->rcols[rcolsinsubprob[j]]]))) )
4078  break; /* a relevant odd varbound has been found */
4079  }
4080  if( i == nrrowsinsubprob && j == nrcolsinsubprob )
4081  {
4082  /* no odd rhs value exists */
4083  for( i = 0 ; i < nrrowsinsubprob ; ++i)
4084  {
4085  lppos = SCIProwGetLPPos(lpdata->rows[problem->rrows[rrowsinsubprob[i]]]);
4086  lpdata->subproblemsindexofrow[lppos] = IRRELEVANT;
4087  if( lpdata->rrowsindexofleftrow[lppos] == rrowsinsubprob[i] )
4088  lpdata->rrowsindexofleftrow[lppos] = ROW_IN_SUBPROB_WITHOUT_ODD_RHS;
4089  if( lpdata->rrowsindexofrightrow[lppos] == rrowsinsubprob[i] )
4090  lpdata->rrowsindexofrightrow[lppos] = ROW_IN_SUBPROB_WITHOUT_ODD_RHS;
4091  }
4092  for( j = 0 ; j < nrcolsinsubprob ; ++j)
4093  {
4094  lppos = SCIPcolGetLPPos(lpdata->cols[problem->rcols[rcolsinsubprob[j]]]);
4095  lpdata->subproblemsindexofcol[lppos] = IRRELEVANT;
4096  lpdata->rcolsindexofcol[lppos] = COLUMN_IN_SUBPROB_WITHOUT_ODD_RHS;
4097
4098  /* statistics */
4099  if( !SCIPisInfinity(scip, problem->rcolslbslack[rcolsinsubprob[j]]) )
4100  ndelvarbounds++;
4101  if( !SCIPisInfinity(scip, problem->rcolsubslack[rcolsinsubprob[j]]) )
4102  ndelvarbounds++;
4103  }
4104  continue;
4105  }
4106
4107  /* don't create new "sub"problem if problem can't be decomposed */
4108  if( lpdata->nsubproblems == 0 && nprocessedrows == problem->nrrows )
4109  continue;
4110
4111  /* create new subproblem */
4112  SCIP_CALL(ZerohalfSubLPDataCreate(scip, &subproblem));
4113  SCIP_CALL(SCIPallocBlockMemoryArray(scip, &(subproblem->rrows), nrrowsinsubprob));
4114  SCIP_CALL(SCIPallocBlockMemoryArray(scip, &(subproblem->rrowsrhs), nrrowsinsubprob));
4115  SCIP_CALL(SCIPallocBlockMemoryArray(scip, &(subproblem->rrowsslack), nrrowsinsubprob));
4116  subproblem->nrrows = nrrowsinsubprob;
4117  SCIP_CALL(SCIPallocBlockMemoryArray(scip, &(subproblem->rcols), nrcolsinsubprob));
4118  SCIP_CALL(SCIPallocBlockMemoryArray(scip, &(subproblem->rcolslbslack), nrcolsinsubprob));
4119  SCIP_CALL(SCIPallocBlockMemoryArray(scip, &(subproblem->rcolsubslack), nrcolsinsubprob));
4120  subproblem->nrcols = nrcolsinsubprob;
4121
4122  for( i = 0 ; i < nrrowsinsubprob ; ++i)
4123  {
4124  rowindex = problem->rrows[rrowsinsubprob[i]];
4125
4126  subproblem->rrows[i] = rowindex;
4127  subproblem->rrowsrhs[i] = problem->rrowsrhs[rrowsinsubprob[i]];
4128  subproblem->rrowsslack[i] = problem->rrowsslack[rrowsinsubprob[i]];
4129
4130  if( lpdata->subproblemsindexofrow[rowindex] != IRRELEVANT )
4131  {
4132  assert(lpdata->rrowsindexofleftrow[rowindex] >= 0
4133  || lpdata->rrowsindexofrightrow[rowindex] >= 0);
4134  lpdata->subproblemsindexofrow[rowindex] = lpdata->nsubproblems;
4135  if( lpdata->rrowsindexofleftrow[rowindex] >= 0 )
4136  lpdata->rrowsindexofleftrow[rowindex] = i;
4137  if( lpdata->rrowsindexofrightrow[rowindex] >= 0 )
4138  lpdata->rrowsindexofrightrow[rowindex] = i;
4139  }
4140  }
4141
4142  for( i = 0 ; i < nrcolsinsubprob ; ++i)
4143  {
4144  colindex = problem->rcols[rcolsinsubprob[i]];
4145
4146  subproblem->rcols[i] = colindex;
4147  subproblem->rcolslbslack[i] = problem->rcolslbslack[rcolsinsubprob[i]];
4148  subproblem->rcolsubslack[i] = problem->rcolsubslack[rcolsinsubprob[i]];
4149
4150  if( lpdata->subproblemsindexofcol[colindex] != IRRELEVANT )
4151  {
4152  assert(lpdata->rcolsindexofcol[colindex] >= 0);
4153  lpdata->subproblemsindexofcol[colindex] = lpdata->nsubproblems;
4154  lpdata->rcolsindexofcol[colindex] = i;
4155  }
4156  }
4157
4158  lpdata->subproblems[lpdata->nsubproblems] = subproblem;
4159  lpdata->nsubproblems++;
4160
4161  totalnrrows += subproblem->nrrows;
4162  totalnrcols += subproblem->nrcols;
4163
4164  SCIPdebugMsg(scip, "subproblem %d: %d rrows, %d rcols\n", k, subproblem->nrrows, subproblem->nrcols);
4165  }
4166  if( lpdata->nsubproblems == 0 )
4167  {
4168  /* problem couldn't be decomposed into different subproblems, hence keep the entire problem */
4169  lpdata->subproblems[0] = problem;
4170  lpdata->nsubproblems = 1;
4171  totalnrrows = problem->nrrows;
4172  totalnrcols = problem->nrcols;
4173
4174  }
4175  else
4176  {
4177  ZerohalfSubLPDataFree(scip, &problem);
4178  }
4179
4180  /* free temporary memory */
4181  SCIPfreeBlockMemoryArray(scip, &rcolsinsubprob, problem->nrcols);
4182  SCIPfreeBlockMemoryArray(scip, &rrowsinsubproboddrhs, problem->nrrows);
4183  SCIPfreeBlockMemoryArray(scip, &rrowsinsubprob, problem->nrrows);
4184  SCIPfreeBufferArray(scip, &queue);
4185  SCIPfreeBufferArray(scip, &processedcols);
4186  SCIPfreeBufferArray(scip, &processedrows);
4187
4189  ZEROHALFstatisticsMessage("\n");
4190  ZEROHALFstatisticsMessage(" | --------------------------------- problem \
4191 -------------------------------- | ----- callback ---- | --total-\n");
4192  ZEROHALFstatisticsMessage(" | nrrows | nrcols | ndlrrows | ndlrcols \
4193 | nsubprob | ndelsubp | ndlvbnds | nsepcuts | ncutsfnd | time\n");
4194  ZEROHALFstatisticsMessage("%15s | %8d | %8d | %8d | %8d | %8d | %8d | %8d | %8d | %8d | %8.4f\n",
4195  "DECOMPOSITION", totalnrrows, totalnrcols, nrrowsinitial - totalnrrows, nrcolsinitial - totalnrcols,
4196  lpdata->nsubproblems, k - lpdata->nsubproblems,
4197  ndelvarbounds,
4199
4200 #else
4201  assert(scip != NULL);
4203  assert(lpdata != NULL);
4204 #endif
4205
4206  return SCIP_OKAY;
4207 }
4208
4209
4210 /** removes the largest number of columns such that the sum of the corresponding variables is at most delta */
4211 static
4213  SCIP* scip, /**< SCIP data structure */
4215  ZEROHALF_MOD2DATA* mod2data, /**< considered (preprocessed) subproblem mod 2 */
4216  SCIP_Real delta /**< delta value */
4217  )
4218 {
4219  int ncolsremoved;
4220  int c;
4221  SCIP_Real maxsumfracsols;
4222  SCIP_Real sumfracsols;
4223
4224
4225  assert(scip != NULL);
4227  assert(mod2data != NULL);
4228  assert(delta >= 0.0);
4229  assert(delta <= 1.0);
4230
4231
4232  /* check if matrix contains rows or columns */
4233  if( mod2data->ncolsind == 0 || mod2data->nrowsind == 0 )
4234  return SCIP_OKAY;
4235
4236  /* check if delta is positive */
4237  if( !SCIPisPositive(scip, delta) )
4238  return SCIP_OKAY;
4239
4240
4241  ncolsremoved = 0;
4242  sumfracsols = 0.0;
4243  maxsumfracsols = sepadata->maxslack * delta;
4244
4245  /* sort column indices sets w.r.t. to their primsol values NON-INCREASINGLY */
4246  if( mod2data->ncolsind > 1 )
4247  {
4248  SCIPsortInd( mod2data->colsind , compRealNonIncreasing , (void*) mod2data->fracsol , mod2data->ncolsind );
4249  }
4250
4251  for( c = mod2data->ncolsind - 1 ; c >= 0 ; --c)
4252  {
4253  if( SCIPisGT(scip, sumfracsols + mod2data->fracsol[mod2data->colsind[c]], maxsumfracsols) )
4254  break;
4255
4256  sumfracsols += mod2data->fracsol[mod2data->colsind[c]];
4258  ncolsremoved++;
4259  }
4260
4261  if( ncolsremoved > 0 )
4262  {
4263  mod2data->ncolsind -= ncolsremoved;
4265  }
4266
4267  return SCIP_OKAY;
4268 }
4269
4270
4271 /** removes some rows that cannot be combined because the resulting slack would be larger than maxslack */
4272 static
4274  SCIP* scip, /**< SCIP data structure */
4276  ZEROHALF_LPDATA* lpdata, /**< data of current LP relaxation */
4277  ZEROHALF_MOD2DATA* mod2data, /**< considered (preprocessed) subproblem mod 2 */
4278  SCIP_Bool removelargeslackrows, /**< should rows with slack + minslack > maxslack be removed? */
4279  SCIP_Bool removelargecolrows /**< should rows with "large valued" columns that cannot be negated be removed? */
4280  )
4281 {
4282  int first;
4283  int last;
4284  int temp;
4285  SCIP_Real minslackoddrhsrows;
4286  SCIP_Real minslackrowwithnonz;
4287  int noddrhsrows;
4288  int nlslrowsremoved;
4289  int nlcolrowsremoved;
4290 #ifndef NDEBUG
4291  int r;
4292 #endif
4293  int c;
4294  SCIP_Bool* removerow;
4295  int i;
4296  int j;
4297  int rowsbind;
4299
4300  assert(scip != NULL);
4302  assert(lpdata != NULL);
4303  assert(mod2data != NULL);
4304  assert(removelargeslackrows || removelargecolrows);
4305
4306
4307  /* check if( A mod 2, b mod 2) is empty */
4308  if( mod2data->nrows == 0 || mod2data->nrowsind == 0 )
4309  return SCIP_OKAY;
4310
4311
4312  /* partition rows into odd-rhs-rows and even-rhs-rows */
4313  first = 0;
4314  last = mod2data->nrowsind - 1;
4315  while(first < last)
4316  {
4317  if( !mod2data->rhs[mod2data->rowsind[first]] )
4318  {
4319  temp = mod2data->rowsind[first];
4320  mod2data->rowsind[first] = mod2data->rowsind[last];
4321  mod2data->rowsind[last] = temp;
4322  --last;
4323  }
4324  else
4325  ++first;
4326  }
4327  noddrhsrows = first + (mod2data->rhs[mod2data->rowsind[first]] ? 1 : 0);
4328
4329
4330  /* check if odd rows exists */
4331  if( noddrhsrows == 0 )
4332  return SCIP_OKAY;
4333
4334  /* sort each partition by nondecreasing slacks */
4335  assert(noddrhsrows >= 0);
4336  SCIPsortInd( mod2data->rowsind , compRealNonDecreasing , (void*) mod2data->slacks , noddrhsrows );
4337  if( noddrhsrows < mod2data->nrowsind )
4338  {
4339  SCIPsortInd( mod2data->rowsind + noddrhsrows , compRealNonDecreasing , (void*) mod2data->slacks ,
4340  mod2data->nrowsind - noddrhsrows );
4341  }
4342
4343  minslackoddrhsrows = mod2data->slacks[mod2data->rowsind[0]];
4344  nlslrowsremoved = 0;
4345  nlcolrowsremoved = 0;
4346
4347  if( SCIPisFeasZero(scip, minslackoddrhsrows) )
4348  return SCIP_OKAY;
4349
4350  /* check if a zerohalf cut may be generated */
4351  if( SCIPisGT(scip, minslackoddrhsrows, sepadata->maxslack) )
4352  {
4353  if( removelargeslackrows )
4354  {
4355  for( i = 0 ; i < mod2data->nrowsind ; ++i )
4357  for( i = 0 ; i < mod2data->ncolsind ; ++i )
4359  mod2data->nrowsind = 0;
4360  mod2data->ncolsind = 0;
4361  }
4362  }
4363  else
4364  {
4365  SCIP_CALL(SCIPallocBufferArray(scip, &removerow, mod2data->nrowsind));
4366  BMSclearMemoryArray(removerow, mod2data->nrowsind);
4367
4368  /* remove all rows with even rhs and slack > maxslack - minslackoddrhsrows */
4369  if( removelargeslackrows )
4370  {
4371  for( i = noddrhsrows ; i < mod2data->nrowsind ; ++i)
4372  if( SCIPisGT(scip, minslackoddrhsrows + mod2data->slacks[mod2data->rowsind[i]], sepadata->maxslack) )
4373  break;
4374  nlslrowsremoved += mod2data->nrowsind - i;
4375  while(i < mod2data->nrowsind)
4376  {
4377 #ifndef NDEBUG
4378  r = mod2data->rowsind[i];
4379 #endif
4380  assert(!mod2data->rhs[r]);
4381  assert(SCIPisGT(scip, minslackoddrhsrows + mod2data->slacks[r], sepadata->maxslack));
4383  removerow[i] = TRUE;
4384  i++;
4385  }
4386  }
4387
4388  /* consider cols */
4389  if( removelargecolrows )
4390  {
4391  if( mod2data->ncolsind > 0 )
4392  {
4393  /* sort column indices sets w.r.t. to their primsol values NON-INCREASINGLY */
4394  if( mod2data->ncolsind > 1 )
4395  {
4396  SCIPsortInd( mod2data->colsind , compRealNonIncreasing , (void*) mod2data->fracsol , mod2data->ncolsind );
4397  }
4398
4399  j = 0;
4400  while( j < mod2data->ncolsind && SCIPisGT(scip, mod2data->fracsol[mod2data->colsind[j]] + minslackoddrhsrows, sepadata->maxslack) )
4401  {
4402  c = mod2data->colsind[j];
4403  minslackrowwithnonz = 1.0;
4404  rowsbind = (int) GETBITARRAYINDEX(c);
4406  for( i = 0 ; i < mod2data->nrowsind ; ++i)
4407  if( !removerow[i] )
4408  if( mod2data->rows[mod2data->rowsind[i]][rowsbind] & rowsbmask )
4409  if( SCIPisLT(scip, mod2data->slacks[mod2data->rowsind[i]], minslackrowwithnonz) )
4410  minslackrowwithnonz = mod2data->slacks[mod2data->rowsind[i]];
4411
4412  if( minslackrowwithnonz < 1.0 )
4413  {
4414  for( i = 0 ; i < mod2data->nrowsind ; ++i)
4415  if( !removerow[i] )
4416  if( mod2data->rows[mod2data->rowsind[i]][rowsbind] & rowsbmask )
4417  if( SCIPisGT(scip, minslackrowwithnonz + mod2data->slacks[mod2data->rowsind[i]], sepadata->maxslack) )
4418  {
4419  markRowAsRemoved(mod2data, i, LARGE_COL_EXISTS);
4420  removerow[i] = TRUE;
4421  nlcolrowsremoved++;
4422  }
4423  }
4424  j++;
4425  }
4426  }
4427  }
4428
4429  /* update mod2data->rowsind if necessary */
4430  if( nlslrowsremoved + nlcolrowsremoved > 0 )
4431  {
4432  j = 0;
4433  for( i = 0 ; i < mod2data->nrowsind && j < mod2data->nrowsind; ++i)
4434  {
4435  if( i < mod2data->nrowsind )
4436  while( j < mod2data->nrowsind && removerow[j] )
4437  j++;
4438  if( i < j && j < mod2data->nrowsind )
4439  mod2data->rowsind[i] = mod2data->rowsind[j];
4440  j++;
4441  }
4442  mod2data->nrowsind -= (nlslrowsremoved + nlcolrowsremoved);
4443  }
4444
4445  /* free temporary memory */
4446  SCIPfreeBufferArray(scip, &removerow);
4447  }
4448
4449  return SCIP_OKAY;
4450 }
4451
4452
4453 /** aggregates identical columns into one column whose (artificial) LP solution is the sum of the aggregated columns */
4454 static
4456  SCIP* scip, /**< SCIP data structure */
4457  ZEROHALF_MOD2DATA* mod2data /**< considered (preprocessed) subproblem mod 2 */
4458  )
4459 {
4460  int c1;
4461  int c2;
4462  int r;
4463  int rowsbind1;
4465  int rowsbind2;
4467  int ncolsremoved;
4468  SCIP_Bool* removecol;
4469
4470  assert(scip != NULL);
4471  assert(mod2data != NULL);
4472
4473  /* check if( A mod 2, b mod 2) is empty */
4474  if( mod2data->nrows == 0 || mod2data->nrowsind == 0 || mod2data->ncolsind == 0 )
4475  return SCIP_OKAY;
4476
4477
4478  /* allocate and initialize temporary memory */
4479  SCIP_CALL(SCIPallocBufferArray(scip, &removecol, mod2data->ncolsind));
4480  BMSclearMemoryArray(removecol, mod2data->ncolsind);
4481  ncolsremoved = 0;
4482
4483
4484  /* check each pair of columns */
4485  for( c1 = 0 ; c1 < mod2data->ncolsind - 1 ; ++c1)
4486  {
4487  rowsbind1 = (int) GETBITARRAYINDEX(mod2data->colsind[c1]);
4489  for( c2 = c1 + 1 ; c2 < mod2data->ncolsind ; ++c2)
4490  {
4491  rowsbind2 = (int) GETBITARRAYINDEX(mod2data->colsind[c2]);
4493  for( r = 0 ; r < mod2data->nrowsind ; ++r)
4495  != (mod2data->rows[mod2data->rowsind[r]][rowsbind2] & rowsbmask2) )
4496  break;
4497  if( r == mod2data->nrowsind )
4498  {
4499  /* a pair of identical columns have been found */
4500
4501  mod2data->fracsol[mod2data->colsind[c2]] += mod2data->fracsol[mod2data->colsind[c1]];
4502  removecol[c1] = TRUE;
4503  ncolsremoved++;
4505  }
4506  }
4507  }
4508
4509  /* update mod2data->colsind array if necessary*/
4510  if( ncolsremoved > 0 )
4511  {
4512  c1 = 0;
4513  for( c2 = 0 ; c2 < mod2data->ncolsind && c1 < mod2data->ncolsind; ++c2)
4514  {
4515  if( c2 < mod2data->ncolsind )
4516  while( c1 < mod2data->ncolsind && removecol[c1] )
4517  c1++;
4518  if( c2 < c1 && c1 < mod2data->ncolsind )
4519  mod2data->colsind[c2] = mod2data->colsind[c1];
4520  c1++;
4521  }
4522  mod2data->ncolsind -= ncolsremoved;
4523  }
4524
4525  /* free temporary memory */
4526  SCIPfreeBufferArray(scip, &removecol);
4527
4528  return SCIP_OKAY;
4529 }
4530
4531
4532 /** preprocess subproblem */
4533 static
4535  SCIP* scip, /**< SCIP data structure */
4536  SCIP_SEPA* sepa, /**< separator */
4538  ZEROHALF_LPDATA* lpdata, /**< data of current LP relaxation */
4539  ZEROHALF_MOD2DATA* mod2data, /**< considered (preprocessed) subproblem mod 2 */
4540  char normtype, /**< type of norm to use for efficacy norm calculation */
4541  int maxsepacuts, /**< maximal number of zerohalf cuts separated per separation round */
4542  int maxcuts, /**< maximal number of zerohalf cuts found per separation round (incl. ineff. cuts) */
4543  int* nsepacuts, /**< pointer to store current number of separated zerohalf cuts */
4544  int* nzerohalfcuts, /**< pointer to store current number of found zerohalf cuts */
4545  ZEROHALF_CUTDATA** zerohalfcuts, /**< array to store found zerohalf cuts */
4546  SCIP_Real** varsolvals, /**< dense variable LP solution vector */
4547  SCIP_RESULT* result /**< pointer to SCIP result value of separation */
4548  )
4549 {
4550  int i;
4551 #ifdef ZEROHALF__PRINT_STATISTICS
4552  int ncolsbeforeppm;
4553  int ncolsinitial;
4554  int nrowsbeforeppm;
4555  int nrowsinitial;
4556  int nsepacutsbeforeppm;
4557  int nsepacutsinitial;
4558  int nzerohalfcutsbeforeppm;
4559  int nzerohalfcutsinitial;
4560  SCIP_CLOCK* timer;
4561  SCIP_CLOCK* pptimer;
4562 #endif
4563  char ppname[SCIP_MAXSTRLEN];
4564
4565  assert(scip != NULL);
4567  assert(lpdata != NULL);
4568  assert(mod2data != NULL);
4569  assert(maxsepacuts >= 0);
4570  assert(maxcuts >= 0);
4571  assert(result != NULL);
4572  assert(nsepacuts != NULL);
4573  assert(nzerohalfcuts != NULL);
4574  assert(zerohalfcuts != NULL);
4575  assert(*nsepacuts <= *nzerohalfcuts);
4576
4577  assert(mod2data->relatedsubproblem != NULL);
4578  assert(mod2data->rows != NULL);
4579  assert(mod2data->rowaggregations != NULL);
4580  assert(mod2data->rhs != NULL);
4581  assert(mod2data->slacks != NULL);
4582  assert(mod2data->fracsol != NULL);
4583  assert(mod2data->nrows > 0);
4584  assert(mod2data->rowsind != NULL);
4585  assert(mod2data->colsind != NULL);
4586
4587  if( sepadata->nppmethods == -1 )
4588  {
4592  }
4593
4594  if( sepadata->nppmethods == 0 )
4595  return SCIP_OKAY;
4596
4597  /* statistics */
4598 #ifdef ZEROHALF__PRINT_STATISTICS
4599  if( sepadata->pptimers == NULL )
4600  {
4602  for( i = 0 ; i <= sepadata->nppmethods; ++i)
4603  {
4605  }
4606  }
4608 #endif
4609
4610  if( mod2data->nrowsind == 0 || mod2data->ncolsind == 0 )
4611  return SCIP_OKAY;
4612
4613 #ifdef ZEROHALF__PRINT_STATISTICS
4614  ncolsinitial = mod2data->ncolsind;
4615  nrowsinitial = mod2data->nrowsind;
4616  nsepacutsinitial = *nsepacuts;
4617  nzerohalfcutsinitial = *nzerohalfcuts;
4618 #endif
4619
4620 #ifdef ZEROHALF__PRINT_STATISTICS
4621  ZEROHALFstatisticsMessage("\n");
4622  ZEROHALFstatisticsMessage(" | ------------------------------- subproblem\
4623  ------------------------------- | ----- callback ---- | --total-\n");
4624  ZEROHALFstatisticsMessage(" | nrowsind | ncolsind | ndelrows | ndelcols \
4625 | nsepcuts | ncutsfnd | time | nsepcuts | ncutsfnd | time\n");
4626  ZEROHALFstatisticsMessage("%15s | %8d | %8d | %8d | %8d | %8d | %8d | %8.4f | %8d | %8d | %8.4f\n",
4627  "START PREPROCSS", mod2data->nrowsind, mod2data->ncolsind, 0, 0, 0, 0, 0.0,
4629  ZEROHALFcreateNewTimer(timer);
4630  ZEROHALFstartTimer(timer);
4631
4632  ZEROHALFcreateNewTimer(pptimer);
4633 #endif
4634
4635  for( i = 0 ; i < sepadata->nppmethods ; ++i)
4636  {
4637  /* abort if enough cuts have already been found */
4638  if( *nsepacuts >= maxsepacuts || *nzerohalfcuts >= maxcuts )
4639  break;
4640 #ifndef ZEROHALF__PRINT_STATISTICS
4641  /* abort preprocessing if matrix is empty */
4642  if( mod2data->nrowsind == 0 && mod2data->ncolsind == 0 )
4643  break;
4644 #endif
4645
4646 #ifdef ZEROHALF__PRINT_STATISTICS
4647  /* statistics*/
4648  ZEROHALFstartTimer(pptimer);
4650  ncolsbeforeppm = mod2data->ncolsind;
4651  nrowsbeforeppm = mod2data->nrowsind;
4652  nsepacutsbeforeppm = *nsepacuts;
4653  nzerohalfcutsbeforeppm = *nzerohalfcuts;
4654 #endif
4655
4656  /* apply preprocessing method */
4658  {
4661  strncpy(ppname,"gauss",SCIP_MAXSTRLEN);
4662  break;
4663  case DELETEZEROROWS:
4665  0, mod2data->nrowsind, TRUE, FALSE, FALSE));
4666  strncpy(ppname,"zero rows",SCIP_MAXSTRLEN);
4667  break;
4668  case DELETEZEROCOLS:
4670  0, mod2data->ncolsind, TRUE, FALSE, FALSE));
4671  strncpy(ppname,"zero columns",SCIP_MAXSTRLEN);
4672  break;
4673  case DELETECOLSINGLETONS:
4675  0, mod2data->ncolsind, FALSE, TRUE, TRUE));
4676  strncpy(ppname,"col singletons",SCIP_MAXSTRLEN);
4677  break;
4679  SCIP_CALL(preprocessTrivialZerohalfCuts(scip, sepa, sepadata, lpdata, mod2data,
4680  0, mod2data->nrowsind, normtype, maxsepacuts, maxcuts, nsepacuts,
4681  nzerohalfcuts, zerohalfcuts, varsolvals, PPZEROONEROW, result));
4682  strncpy(ppname,"trivial cuts",SCIP_MAXSTRLEN);
4683  break;
4684  case DELETEIDENTROWS:
4686  0, mod2data->nrowsind, FALSE, FALSE, TRUE));
4687  strncpy(ppname,"identical rows",SCIP_MAXSTRLEN);
4688  break;
4689  case MERGEIDENTCOLS:
4690  SCIP_CALL(preprocessIdenticalColums(scip, mod2data));
4691  strncpy(ppname,"identical cols",SCIP_MAXSTRLEN);
4692  break;
4693  case DELETELARGESLACKROWS:
4695  0, mod2data->nrowsind, FALSE, TRUE, FALSE));
4696  strncpy(ppname,"sl>maxsl rows",SCIP_MAXSTRLEN);
4697  break;
4698  case PPCOLUMNS:
4700  0, mod2data->ncolsind, TRUE, TRUE, TRUE));
4701  strncpy(ppname,"(pp cols)",SCIP_MAXSTRLEN);
4702  break;
4703  case PPROWS:
4705  0, mod2data->nrowsind, TRUE, TRUE, TRUE));
4706  strncpy(ppname,"(pp rows)",SCIP_MAXSTRLEN);
4707  break;
4711  strncpy(ppname,"delta heur",SCIP_MAXSTRLEN);
4712  break;
4713  case DELETEROWSWRTMINSLACK:
4714  SCIP_CALL(preprocessConsiderMinSlack(scip, sepadata, lpdata, mod2data, TRUE, TRUE));
4715  strncpy(ppname,"s_odd",SCIP_MAXSTRLEN);
4716  break;
4717  default:
4718  SCIPerrorMessage("invalid preprocessing method '%c'\n", sepadata->ppmethods[i]);
4719  return SCIP_INVALIDDATA;
4720  }
4721
4722 #ifdef ZEROHALF__PRINT_STATISTICS
4723  /* statistics */
4725  ZEROHALFstopTimer(pptimer);
4726  ZEROHALFstatisticsMessage("%15s | %8d | %8d | %8d | %8d | %8d | %8d | %8.4f | %8d | %8d | %8.4f\n",
4727  ppname, mod2data->nrowsind, mod2data->ncolsind,
4728  nrowsbeforeppm - mod2data->nrowsind, ncolsbeforeppm - mod2data->ncolsind,
4729  *nsepacuts - nsepacutsbeforeppm, *nzerohalfcuts - nzerohalfcutsbeforeppm,
4730  ZEROHALFevalTimer(pptimer), *nsepacuts, *nzerohalfcuts,
4732  ZEROHALFresetTimer(pptimer);
4733 #endif
4734  }
4735
4736 #ifdef ZEROHALF__PRINT_STATISTICS
4737  /* statistics */
4738  ZEROHALFstopTimer(timer);
4739  ZEROHALFstatisticsMessage("%15s | %8d | %8d | %8d | %8d | %8d | %8d | %8.4f | %8d | %8d | %8.4f\n",
4740  "PREPROCESSED", mod2data->nrowsind, mod2data->ncolsind,
4741  nrowsinitial - mod2data->nrowsind, ncolsinitial - mod2data->ncolsind,
4742  *nsepacuts - nsepacutsinitial, *nzerohalfcuts - nzerohalfcutsinitial,
4743  ZEROHALFevalTimer(timer), *nsepacuts, *nzerohalfcuts,
4745  ZEROHALFstatisticsMessage("\n");
4746  ZEROHALFfreeTimer(timer);
4747  ZEROHALFfreeTimer(pptimer);
4749
4750  ZEROHALFstatisticsMessage(" | ------------------------------- subproblem ------------------------------- | ------------------------------\n");
4751  ZEROHALFstatisticsMessage(" | | max2/row | max2/col | A^T ept | \n");
4752  ZEROHALFstatisticsMessage("%15s | | %8s | %8s | %8s |\n",
4753  "SUBPROBSTRUCT",
4754  hasMatrixMax2EntriesPerRow(mod2data) ? "yes" : "no", hasMatrixMax2EntriesPerColumn(mod2data) ? "yes" : "no", "n/a");
4755  ZEROHALFstatisticsMessage("\n");
4756 #endif
4757
4758  return SCIP_OKAY;
4759 }
4760
4761
4762 /* --------------------------------------------------------------------------------------------------------------------
4763  * local methods: separating methods
4764  * -------------------------------------------------------------------------------------------------------------------- */
4765
4766
4767 /** returns the objective weights for the weighted feasibility AuxIP */
4768 static
4770  BITARRAY rowaggregation, /**< row aggregation bitarray */
4771  int nrrows /**< number of relevant rows */
4772  )
4773 {
4774  int i;
4775  int naggregatedrrows;
4776
4777  assert(rowaggregation != NULL);
4778  assert(nrrows > 0);
4779
4780  naggregatedrrows = 0;
4781  for( i = 0 ; i < nrrows ; ++i)
4782  if( BITARRAYBITISSET(rowaggregation, i) ) /*lint !e701*/
4783  naggregatedrrows++;
4784
4785  return (SCIP_Real) naggregatedrrows;
4786 }
4787
4788
4789 /** creates a "subscip" representing the following auxiliary IP (AuxIP):
4790  * min z := s^T v + x^T y
4791  * s.t. (b (mod 2))^T v - 2q = 1
4792  * (A (mod 2))^T v - y - -2r = 0
4793  *
4794  * v \\in {0,1}^nrowsind
4795  * y \\in {0,1}^ncolsind
4796  * r \\in Z^ncolsind_+
4797  * q \\in Z_+
4798  */
4799 #define BRANCHPRIORITY__AVOID_BRANCHING 0
4800 #define BRANCHPRIORITY__PREFER_BRANCHING 0
4801 static
4803  SCIP* scip, /**< SCIP data structure */
4805  ZEROHALF_LPDATA* lpdata, /**< data of current LP relaxation */
4806  ZEROHALF_MOD2DATA* mod2data, /**< considered (preprocessed) subproblem mod 2 */
4807  ZEROHALF_AUXIPDATA* auxipdata, /**< pointer to data structure to store the auxiliary IP data */
4808  SCIP_Bool setnodelimit /**< should a node limit be set? */
4809  )
4810 {
4811  SCIP_VAR** consvars;
4812  SCIP_Real* consvals;
4813  SCIP_Real maxslack;
4814  char consname[SCIP_MAXSTRLEN];
4815  char varname[SCIP_MAXSTRLEN];
4816  int i;
4817  int j;
4818  int maxnconsvars;
4819  int nconsvars;
4820  SCIP_Bool isfeasip;
4821  SCIP_Bool isweighted;
4822  SCIP_Bool ispenalized;
4823  SCIP_Bool settingsfileexists;
4824
4825  int nrrows;
4826
4827  int rowsbind;
4829
4830  SCIP_Bool success;
4831
4832  SCIP_Real feastol;
4833
4834  assert(scip != NULL);
4836  assert(lpdata != NULL);
4837  assert(mod2data != NULL);
4838  assert(auxipdata != NULL);
4839
4840  assert(mod2data->relatedsubproblem != NULL);
4841  assert(mod2data->rows != NULL);
4842  assert(mod2data->rowaggregations != NULL);
4843  assert(mod2data->rhs != NULL);
4844  assert(mod2data->slacks != NULL);
4845  assert(mod2data->fracsol != NULL);
4846  assert(mod2data->rowsind != NULL);
4847  assert(mod2data->colsind != NULL);
4848
4849  assert(auxipdata->subscip == NULL);
4850  assert(auxipdata->v == NULL);
4851  assert(auxipdata->y == NULL);
4852  assert(auxipdata->r == NULL);
4853  assert(auxipdata->q == NULL);
4854  assert(auxipdata->feasipcons == NULL);
4855  assert(auxipdata->oddrhscons == NULL);
4856  assert(auxipdata->columnsumcons == NULL);
4857
4858  auxipdata->m = mod2data->nrowsind;
4859  auxipdata->n = mod2data->ncolsind;
4860
4861  /* alloc temporary memory for subscipdata elements*/
4862  SCIP_CALL(SCIPallocBlockMemoryArray(scip, &(auxipdata->v), auxipdata->m));
4863  SCIP_CALL(SCIPallocBlockMemoryArray(scip, &(auxipdata->y), auxipdata->n));
4864  SCIP_CALL(SCIPallocBlockMemoryArray(scip, &(auxipdata->r), auxipdata->n));
4865  SCIP_CALL(SCIPallocBlockMemoryArray(scip, &(auxipdata->columnsumcons), auxipdata->n));
4866
4867  /* initialize allocated data structures */
4868  BMSclearMemoryArray(auxipdata->v, auxipdata->m); /* NULL = 0x0 */
4869  BMSclearMemoryArray(auxipdata->y, auxipdata->n); /* NULL = 0x0 */
4870  BMSclearMemoryArray(auxipdata->r, auxipdata->n); /* NULL = 0x0 */
4871  BMSclearMemoryArray(auxipdata->columnsumcons, auxipdata->n); /* NULL = 0x0 */
4872
4874  nrrows = mod2data->relatedsubproblem->nrrows;
4875
4876  /* determine subscip limits */
4877  SCIP_CALL( SCIPgetRealParam(scip, "limits/time", &auxipdata->timelimit) );
4878  if( !SCIPisInfinity(scip, auxipdata->timelimit) )
4879  auxipdata->timelimit -= SCIPgetSolvingTime(scip);
4880
4881  /* substract the memory already used by the main SCIP and the estimated memory usage of external software */
4882  SCIP_CALL( SCIPgetRealParam(scip, "limits/memory", &auxipdata->memorylimit) );
4883  if( !SCIPisInfinity(scip, auxipdata->memorylimit) )
4884  {
4885  auxipdata->memorylimit -= SCIPgetMemUsed(scip)/1048576.0;
4886  auxipdata->memorylimit -= SCIPgetMemExternEstim(scip)/1048576.0;
4887  }
4888
4889  if( setnodelimit == TRUE )
4890  auxipdata->nodelimit = 3000;
4891  else
4892  auxipdata->nodelimit = -1;
4893
4894  feastol = SCIPfeastol(scip);
4895  auxipdata->objectivelimit = MIN(1.0, maxslack + feastol);
4896
4897  /* abort if not enough memory available */
4898  if( auxipdata->memorylimit <= 2.0*SCIPgetMemExternEstim(scip)/1048576.0 )
4899  return SCIP_OKAY;
4900
4901  /* abort if not enough time available */
4902  if( auxipdata->timelimit <= 0.0 )
4903  return SCIP_OKAY;
4904
4905  /* alloc further temporary memory */
4906  maxnconsvars = auxipdata->m + auxipdata->n + 2;
4907  SCIP_CALL(SCIPallocBufferArray(scip, &consvals, maxnconsvars));
4908  SCIP_CALL(SCIPallocBufferArray(scip, &consvars, maxnconsvars));
4909
4910  /* create and initialize framework */
4911
4912  SCIP_CALL( SCIPcreate(&(auxipdata->subscip)) );
4913  success = FALSE;
4914 #ifndef NDEBUG
4915  SCIP_CALL( SCIPcopyPlugins(scip, auxipdata->subscip, TRUE, TRUE, TRUE, TRUE, TRUE, TRUE,
4916  TRUE, TRUE, TRUE, TRUE, TRUE, TRUE, TRUE, TRUE, FALSE, TRUE, &success) );
4917 #else
4918  SCIP_CALL( SCIPcopyPlugins(scip, auxipdata->subscip, FALSE, FALSE, TRUE, TRUE, TRUE, TRUE,
4919  TRUE, TRUE, TRUE, TRUE, TRUE, TRUE, TRUE, FALSE, FALSE, TRUE, &success) );
4920 #endif
4921  SCIPdebugMsg(scip, "Copying the plugins was %s successful.\n", success ? "" : "not");
4922
4923  SCIP_CALL( SCIPcreateProb(auxipdata->subscip, "sepa_zerohalf auxiliary IP (AuxIP)",
4924  NULL, NULL , NULL , NULL , NULL , NULL , NULL) );
4925
4926  settingsfileexists = TRUE;
4927  if( strlen(sepadata->subscipsettings) == 0 )
4928  settingsfileexists = FALSE;
4930  settingsfileexists = FALSE;
4931
4932  if( settingsfileexists )
4933  {
4934  /* read subscip settings file */
4936  }
4937  else
4938  {
4939  /* do not abort subscip on CTRL-C */
4940  SCIP_CALL(SCIPsetBoolParam(auxipdata->subscip, "misc/catchctrlc", FALSE));
4941
4942  /* disable output to console */
4943 #ifdef SCIP_DEBUG
4944  SCIP_CALL(SCIPsetIntParam(auxipdata->subscip, "display/verblevel", 4));
4945  SCIP_CALL(SCIPsetIntParam(auxipdata->subscip, "display/freq", 1));
4946 #else
4947  SCIP_CALL(SCIPsetIntParam(auxipdata->subscip, "display/verblevel", 0));
4948  SCIP_CALL(SCIPsetIntParam(auxipdata->subscip, "display/freq", 1000));
4949 #endif
4950  SCIP_CALL(SCIPsetIntParam(auxipdata->subscip, "display/nsols/active", 2));
4951
4952  /* forbid recursive call of heuristics solving subMIPs */
4953  SCIP_CALL(SCIPsetIntParam(auxipdata->subscip, "heuristics/rins/freq", -1));
4954  SCIP_CALL(SCIPsetIntParam(auxipdata->subscip, "heuristics/rens/freq", -1));
4955  SCIP_CALL(SCIPsetIntParam(auxipdata->subscip, "heuristics/localbranching/freq", -1));
4956  SCIP_CALL(SCIPsetIntParam(auxipdata->subscip, "heuristics/crossover/freq", -1));
4957
4958  /* disable cut separation in subscip */
4959  SCIP_CALL(SCIPsetIntParam(auxipdata->subscip, "separating/zerohalf/freq", -1));
4960  /* SCIP_CALL(SCIPsetIntParam(auxipdata->subscip, "separating/maxrounds", 0)); */
4961  /* SCIP_CALL(SCIPsetIntParam(auxipdata->subscip, "separating/maxroundsroot", 0)); */
4962  /* SCIP_CALL(SCIPsetIntParam(auxipdata->subscip, "separating/maxcuts", 0)); */
4963  /* SCIP_CALL(SCIPsetIntParam(auxipdata->subscip, "separating/maxcutsroot", 0)); */
4964
4965  /* use pseudo cost branching without strong branching */
4966  /* SCIP_CALL(SCIPsetIntParam(auxipdata->subscip, "branching/pscost/priority", INT_MAX/4)); */
4967
4968  /* disable expensive presolving */
4969  /* SCIP_CALL(SCIPsetIntParam(auxipdata->subscip, "presolving/probing/maxrounds", 0)); */
4970  /* SCIP_CALL(SCIPsetIntParam(auxipdata->subscip, "constraints/linear/maxpresolpairrounds", 0)); */
4971  /* SCIP_CALL(SCIPsetRealParam(auxipdata->subscip, "constraints/linear/maxaggrnormscale", 0.0)); */
4972
4973  /* disable conflict analysis */
4974  /* SCIP_CALL(SCIPsetBoolParam(auxipdata->subscip, "conflict/useprop", FALSE)); */
4975  /* SCIP_CALL(SCIPsetBoolParam(auxipdata->subscip, "conflict/useinflp", FALSE)); */
4976  /* SCIP_CALL(SCIPsetBoolParam(auxipdata->subscip, "conflict/useboundlp", FALSE)); */
4977  /* SCIP_CALL(SCIPsetBoolParam(auxipdata->subscip, "conflict/usesb", FALSE)); */
4978  /* SCIP_CALL(SCIPsetBoolParam(auxipdata->subscip, "conflict/usepseudo", FALSE)); */
4979
4980  SCIP_CALL(SCIPsetBoolParam(auxipdata->subscip, "branching/preferbinary", TRUE));
4981  SCIP_CALL(SCIPsetIntParam(auxipdata->subscip, "heuristics/shifting/freq", 3));
4982  SCIP_CALL(SCIPsetIntParam(auxipdata->subscip, "heuristics/simplerounding/freq", 1));
4983  SCIP_CALL(SCIPsetIntParam(auxipdata->subscip, "heuristics/rounding/freq", 1));
4984  SCIP_CALL(SCIPsetIntParam(auxipdata->subscip, "heuristics/oneopt/freq", 1));
4985
4986  /* SCIP_CALL(SCIPsetIntParam(auxipdata->subscip, "heuristics/pscostdiving/freq", 1)); */
4987  /* SCIP_CALL(SCIPsetIntParam(auxipdata->subscip, "heuristics/feaspump/freq", 3)); */
4988
4989  /* SCIP_CALL(SCIPsetIntParam(auxipdata->subscip, "heuristics/coefdiving/freq", -1)); */
4990  /* SCIP_CALL(SCIPsetIntParam(auxipdata->subscip, "heuristics/fracdiving/freq", -1)); */
4991  /* SCIP_CALL(SCIPsetIntParam(auxipdata->subscip, "heuristics/guideddiving/freq", -1)); */
4992  /* SCIP_CALL(SCIPsetIntParam(auxipdata->subscip, "heuristics/linesearchdiving/freq", -1)); */
4993  /* SCIP_CALL(SCIPsetIntParam(auxipdata->subscip, "heuristics/objpscostdiving/freq", -1)); */
4994  /* SCIP_CALL(SCIPsetIntParam(auxipdata->subscip, "heuristics/rootsoldiving/freq", -1)); */
4995  /* SCIP_CALL(SCIPsetIntParam(auxipdata->subscip, "heuristics/veclendiving/freq", -1)); */
4996  }
4997
4998  /* get type of auxiliary IP objective function */
4999  isfeasip = (sepadata->subscipobjective == 'v' ? FALSE : TRUE);
5000  isweighted = (sepadata->subscipobjective == 'w' ? TRUE : FALSE);
5001  ispenalized = (sepadata->subscipobjective == 'p' ? TRUE : FALSE);
5002
5003  /* set limits of subscip */
5004  SCIP_CALL( SCIPsetLongintParam(auxipdata->subscip, "limits/nodes", (SCIP_Longint) auxipdata->nodelimit) );
5005  SCIP_CALL( SCIPsetRealParam(auxipdata->subscip, "limits/time", auxipdata->timelimit) );
5006  SCIP_CALL( SCIPsetRealParam(auxipdata->subscip, "limits/memory", auxipdata->memorylimit) );
5007
5008  if( !isfeasip )
5009  {
5010  SCIP_CALL( SCIPsetObjlimit(auxipdata->subscip, auxipdata->objectivelimit) );
5011  }
5012  SCIP_CALL( SCIPsetIntParam(auxipdata->subscip, "limits/solutions", sepadata->subscipsollimit) );
5013
5014
5015  /* create variables and set objective */
5016  /* q */
5017  SCIP_CALL( SCIPcreateVar(auxipdata->subscip, &(auxipdata->q), "q", 0.0, SCIPinfinity(auxipdata->subscip),
5020  SCIP_CALL( SCIPchgVarBranchPriority(auxipdata->subscip, auxipdata->q, BRANCHPRIORITY__AVOID_BRANCHING) );
5021  /* r */
5022  for( j = 0 ; j < auxipdata->n ; ++j)
5023  {
5024  (void) SCIPsnprintf(varname, SCIP_MAXSTRLEN, "r_colsind%d(rcols%d)", j, mod2data->colsind[j]);
5025  SCIP_CALL( SCIPcreateVar(auxipdata->subscip, &(auxipdata->r[j]), varname, 0.0, SCIPinfinity(auxipdata->subscip),
5028  SCIP_CALL( SCIPchgVarBranchPriority(auxipdata->subscip, auxipdata->q, BRANCHPRIORITY__AVOID_BRANCHING) );
5029  }
5030  /* v */
5031  for( i = 0 ; i < auxipdata->m ; ++i)
5032  {
5033  SCIP_Real objcoef;
5034  assert(mod2data->rows[mod2data->rowsind[i]] != NULL);
5035
5036  (void) SCIPsnprintf(varname, SCIP_MAXSTRLEN, "v_rowsind%d(rrows%d)", i, mod2data->rowsind[i]);
5037  if( isfeasip )
5038  {
5039  if( isweighted )
5040  {
5041  objcoef = calcObjWeight(mod2data->rowaggregations[mod2data->rowsind[i]], nrrows);
5042  }
5043  else if( ispenalized )
5044  {
5045  objcoef = mod2data->slacks[mod2data->rowsind[i]]
5047  * calcObjWeight(mod2data->rowaggregations[mod2data->rowsind[i]], nrrows));
5048  }
5049  else
5050  objcoef = 1.0;
5051  }
5052  else
5053  objcoef = mod2data->slacks[mod2data->rowsind[i]];
5054  assert(!SCIPisFeasNegative(scip, objcoef));
5055  SCIP_CALL( SCIPcreateVar(auxipdata->subscip, &(auxipdata->v[i]), varname, 0.0, 1.0, objcoef,
5058  SCIP_CALL(SCIPchgVarBranchPriority(auxipdata->subscip, auxipdata->q, BRANCHPRIORITY__PREFER_BRANCHING));
5059  }
5060  /* y */
5061  for( j = 0 ; j < auxipdata->n ; ++j)
5062  {
5063  SCIP_Real objcoef;
5064  (void) SCIPsnprintf(varname, SCIP_MAXSTRLEN, "y_%d(%d)", j, mod2data->colsind[j]);
5065  if( isfeasip )
5066  objcoef = 0.0;
5067  else
5068  objcoef = mod2data->fracsol[mod2data->colsind[j]];
5069  SCIP_CALL( SCIPcreateVar(auxipdata->subscip, &(auxipdata->y[j]), varname, 0.0, 1.0, objcoef,
5072  }
5073
5074  /* create constraints */
5075  /* "feasibility constraint" */
5076  if( isfeasip )
5077  {
5078  nconsvars = 0;
5079  for( i = 0 ; i < auxipdata->m ; ++i)
5080  {
5081  consvals[nconsvars] = mod2data->slacks[mod2data->rowsind[i]];
5082  consvars[nconsvars] = auxipdata->v[i];
5083  nconsvars++;
5084  }
5085  for( j = 0 ; j < auxipdata->n ; ++j)
5086  {
5087  consvals[nconsvars] = mod2data->fracsol[mod2data->colsind[j]];
5088  consvars[nconsvars] = auxipdata->y[j];
5089  nconsvars++;
5090  }
5091  SCIP_CALL( SCIPcreateConsLinear(auxipdata->subscip, &(auxipdata->oddrhscons), "feas",
5092  nconsvars, consvars, consvals, 0.0, auxipdata->objectivelimit,
5093  TRUE, TRUE, TRUE, TRUE, TRUE, FALSE, FALSE, FALSE, FALSE, FALSE) );
5095  }
5096  /* "odd rhs" */
5097  nconsvars = 0;
5098  for( i = 0 ; i < auxipdata->m ; ++i)
5099  if( mod2data->rhs[mod2data->rowsind[i]] == TRUE )
5100  {
5101  consvals[nconsvars] = 1.0;
5102  consvars[nconsvars] = auxipdata->v[i];
5103  nconsvars++;
5104  }
5105  consvals[nconsvars] = -2.0;
5106  consvars[nconsvars] = auxipdata->q;
5107  nconsvars++;
5108  SCIP_CALL( SCIPcreateConsLinear(auxipdata->subscip, &(auxipdata->oddrhscons), "odd_rhs",
5109  nconsvars, consvars, consvals, 1.0, 1.0,
5110  TRUE, TRUE, TRUE, TRUE, TRUE, FALSE, FALSE, FALSE, FALSE, FALSE) );
5112  /* "column sum" */
5113  for( j = 0 ; j < auxipdata->n ; ++j)
5114  {
5115  nconsvars = 0;
5116
5117  rowsbind = (int) GETBITARRAYINDEX(mod2data->colsind[j]);
5119  for( i = 0 ; i < auxipdata->m ; ++i) {
5120  if( mod2data->rows[mod2data->rowsind[i]][rowsbind] & rowsbmask )
5121  {
5122  consvals[nconsvars] = 1.0;
5123  consvars[nconsvars] = auxipdata->v[i];
5124  nconsvars++;
5125  }
5126  }
5127  consvals[nconsvars] = -1.0;
5128  consvars[nconsvars] = auxipdata->y[j];
5129  nconsvars++;
5130  consvals[nconsvars] = -2.0;
5131  consvars[nconsvars] = auxipdata->r[j];
5132  nconsvars++;
5133
5134  (void) SCIPsnprintf(consname, SCIP_MAXSTRLEN, "col_%d(%d)_sum", j, mod2data->colsind[j]);
5135  SCIP_CALL( SCIPcreateConsLinear(auxipdata->subscip, &(auxipdata->columnsumcons[j]) , consname, nconsvars, consvars, consvals, 0.0, 0.0,
5136  TRUE, TRUE, TRUE, TRUE, TRUE, FALSE, FALSE, FALSE, FALSE, FALSE) );
5138  }
5139
5140  /* free temporary memory */
5141  SCIPfreeBufferArray(scip, &consvars);
5142  SCIPfreeBufferArray(scip, &consvals);
5143
5144  SCIPdebug( SCIP_CALL( SCIPprintOrigProblem(auxipdata->subscip, NULL, NULL, TRUE) ) );
5145
5146  return SCIP_OKAY;
5147 }
5148
5149
5150 /** solves the auxiliary IP given as subscip */
5151 static
5153  SCIP* scip, /**< SCIP data structure */
5155  ZEROHALF_MOD2DATA* mod2data, /**< considered (preprocessed) subproblem mod 2 */
5156  ZEROHALF_AUXIPDATA* auxipdata, /**< auxiliary IP data */
5157  SCIP_SOL*** sols, /**< pointer to store array of solutions */
5158  int* nsols /**< pointer to store number of solutions */
5159  )
5160 {
5161  SCIP_RETCODE retcode;
5162  SCIP_STAGE subscipstage;
5163  SCIP_Real maxslack;
5164  int i;
5165  int j;
5166
5168  assert(mod2data != NULL);
5169  assert(auxipdata != NULL);
5170  assert(auxipdata->subscip != NULL);
5171  assert(sols != NULL);
5172  assert(*sols == NULL);
5173  assert(*nsols == 0);
5174
5175  /* solve AuxIP */
5176  retcode = SCIPsolve(auxipdata->subscip);
5177
5178  /* errors in solving the subproblem should not kill the overall solving process;
5179  * hence, the return code is caught and a warning is printed, only in debug mode, SCIP will stop. */
5180  if ( retcode != SCIP_OKAY )
5181  {
5182 #ifndef NDEBUG
5183  SCIP_CALL( retcode );
5184 #endif
5185  SCIPwarningMessage(scip, "Error while solving subproblem in zerohalf separator; sub-SCIP terminated with code <%d>\n", retcode);
5186  *nsols = 0;
5187  return SCIP_OKAY;
5188  }
5189
5190  /* print statistic */
5191  SCIPdebug( SCIP_CALL( SCIPprintStatistics(auxipdata->subscip, NULL) ) );
5192
5194
5195  /* check if solving was successful and get solutions */
5196  subscipstage = SCIPgetStage(auxipdata->subscip);
5197  if( subscipstage == SCIP_STAGE_SOLVING || subscipstage == SCIP_STAGE_SOLVED )
5198  *nsols = SCIPgetNSols(auxipdata->subscip);
5199  else
5200  *nsols = 0;
5201  if( *nsols > 0 )
5202  {
5203  *sols = SCIPgetSols(auxipdata->subscip);
5204  /* check if only the best solution should be used */
5206  *nsols = 1;
5207  }
5208
5209  /* check if proper a proper solution was found */
5210  if( sepadata->subscipobjective == 'v' )
5211  {
5212  for( i = 0 ; i < *nsols ; ++i)
5213  if( SCIPisGT(scip, SCIPgetSolOrigObj(auxipdata->subscip , (*sols)[i]), maxslack) )
5214  break;
5215  *nsols = i;
5216  }
5217  else
5218  {
5219  SCIP_Real z;
5220  SCIP_Real* viols;
5221  SCIP_SOL** propersols;
5222
5223  SCIP_Bool swapped;
5224  SCIP_Real tempviol;
5225  SCIP_SOL* tempsol;
5226
5227  int npropersols;
5228
5229  SCIP_CALL(SCIPallocBufferArray(scip, &viols, *nsols));
5230  SCIP_CALL(SCIPallocBlockMemoryArray(scip, &propersols, *nsols));
5231  npropersols = 0;
5232  for( i = 0 ; i < *nsols ; ++i)
5233  {
5234  z = 0.0;
5235  for( j = 0 ; j < auxipdata->m ; ++j)
5236  z += mod2data->slacks[mod2data->rowsind[j]]
5237  * SCIPgetSolVal(auxipdata->subscip, (*sols)[i], auxipdata->v[j]);
5238  for( j = 0 ; j < auxipdata->n ; ++j)
5239  z += mod2data->fracsol[mod2data->colsind[j]]
5240  * SCIPgetSolVal(auxipdata->subscip, (*sols)[i], auxipdata->y[j]);
5241
5242  if( SCIPisLE(scip, z, maxslack) )
5243  {
5244  /* proper sol has been found */
5245  propersols[npropersols] = (*sols)[i];
5246  viols[npropersols] = z;
5247  npropersols++;
5248  }
5249  }
5250
5251
5252  swapped = TRUE;
5253  for( i = 1 ; i < npropersols && swapped; ++i)
5254  {
5255  swapped = FALSE;
5256  for( j = 0 ; j < npropersols - i ; ++j)
5257  {
5258  if( viols[j] > viols[j+1] )
5259  {
5260  tempviol = viols[j+1];
5261  viols[j+1] = viols[j];
5262  viols[j] = tempviol;
5263  tempsol = propersols[j+1];
5264  propersols[j+1] = propersols[j];
5265  propersols[j] = tempsol;
5266  swapped = TRUE;
5267  }
5268  }
5269  }
5270
5271  *sols = propersols;
5272  *nsols = npropersols;
5273  SCIPfreeBufferArray(scip, &viols);
5274
5275  }
5276
5277  return SCIP_OKAY;
5278 }
5279
5280
5281
5282
5283 /** determines the weightvector for a single row */
5284 static
5286  SCIP* scip, /**< SCIP data structure */
5288  ZEROHALF_LPDATA* lpdata, /**< data of current LP relaxation */
5289  int rowsindex, /**< lpdata->rows index */
5290  int rrowsindex, /**< "subproblem"->rrows index */
5291  SCIP_Real** weights /**< pointer to store the weight vector */
5292  )
5293 { /*lint --e{438}*/
5295  assert(scip != NULL);
5296  assert(lpdata != NULL);
5297  assert(lpdata->nrows > 0);
5298  assert(0 <= rowsindex);
5299  assert(rowsindex < lpdata->nrows);
5300  assert(lpdata->rrowsindexofleftrow[rowsindex] == rrowsindex
5301  || lpdata->rrowsindexofrightrow[rowsindex] == rrowsindex);
5302  assert(weights != NULL);
5303  assert(*weights == NULL);
5304
5305  /* allocate temporary memory */
5306  SCIP_CALL(SCIPallocBlockMemoryArray(scip, weights, lpdata->nrows));
5307
5308  /* initialize */
5309  BMSclearMemoryArray(*weights, lpdata->nrows);
5310
5311  /* determine row weights */
5312  if( lpdata->rrowsindexofleftrow[rowsindex] == rrowsindex )
5313  (*weights)[rowsindex] = lpdata->intscalarsleftrow[rowsindex] * (-0.5);
5314  else
5315  (*weights)[rowsindex] = lpdata->intscalarsrightrow[rowsindex] * 0.5;
5316
5317  if( SCIProwGetNLPNonz(lpdata->rows[rowsindex]) >= sepadata->maxnnonz )
5318  {
5319  SCIPfreeBlockMemoryArray(scip, weights, lpdata->nrows);
5320  }
5321
5322  return SCIP_OKAY;
5323 }
5324
5325
5326 /** gets the subset of rows that should be combined to a violated zerohalf cut */
5327 static
5329  SCIP* scip, /**< SCIP data structure */
5330  ZEROHALF_MOD2DATA* mod2data, /**< considered (preprocessed) subproblem mod 2 */
5331  ZEROHALF_AUXIPDATA* auxipdata, /**< auxiliary IP data */
5332  SCIP_SOL* solution, /**< considered solution */
5333  BITARRAY* rrowsincut, /**< pointer to store the subset of rows as bitarray (length: number of relevant rows) */
5334  int* nrrowsincut /**< number of combined relevant rows */
5335  )
5336 {
5337  int i;
5338
5339  assert(scip != NULL);
5340  assert(mod2data != NULL);
5341  assert(mod2data->nrowsind > 0);
5342  assert(auxipdata != NULL);
5343  assert(auxipdata->subscip != NULL);
5344  assert(auxipdata->v != NULL);
5345  assert(solution != NULL);
5346  assert(rrowsincut != NULL);
5347  assert(*rrowsincut == NULL);
5348  assert(nrrowsincut != NULL);
5349
5350
5351  /* allocate and initialize temporary memory for calculating the symmetric difference */
5352  SCIP_CALL(SCIPallocBlockMemoryArray(scip, rrowsincut, mod2data->rowaggregationsbitarraysize));
5353  BITARRAYCLEAR(*rrowsincut, mod2data->rowaggregationsbitarraysize);
5354
5355  *nrrowsincut = 0;
5356
5357  /* calculate symmetric difference of rrowsincut and specific rowaggregations */
5358  for( i = 0 ; i < mod2data->nrowsind ; ++i)
5359  if( auxipdata->v[i] != NULL )
5360  if( !SCIPisZero(scip, SCIPgetSolVal(auxipdata->subscip, solution, auxipdata->v[i])) )
5361  {
5362  BITARRAYSXOR(mod2data->rowaggregations[mod2data->rowsind[i]], (*rrowsincut) ,
5363  mod2data->rowaggregationsbitarraysize);
5364  (*nrrowsincut)++;
5365  }
5366
5367  return SCIP_OKAY;
5368 }
5369
5370
5371 /** separates violated zerohalf cuts by solving an auxiliary IP. (exact method; exponential time) */
5372 static
5374  SCIP* scip, /**< SCIP data structure */
5375  SCIP_SEPA* sepa, /**< separator */
5377  ZEROHALF_LPDATA* lpdata, /**< data of current LP relaxation */
5378  ZEROHALF_MOD2DATA* mod2data, /**< considered (preprocessed) subproblem mod 2 */
5379  char normtype, /**< SCIP normtype */
5380  int maxsepacuts, /**< maximal number of zerohalf cuts separated per separation round */
5381  int maxcuts, /**< maximal number of zerohalf cuts found per separation round (incl. ineff. cuts) */
5382  SCIP_Bool setnodelimit, /**< should a node limit be set for solving the auxiliary IP? */
5383  int* nsepacuts, /**< pointer to store current number of separated zerohalf cuts */
5384  int* nzerohalfcuts, /**< pointer to store current number of found zerohalf cuts */
5385  ZEROHALF_CUTDATA** zerohalfcuts, /**< array to store found zerohalf cuts */
5386  SCIP_Real** varsolvals, /**< dense variable LP solution vector */
5387  SCIP_RESULT* result /**< pointer to SCIP result value of separation */
5388  )
5389 {
5390  ZEROHALF_AUXIPDATA* auxipdata;
5391  SCIP_SOL** sols;
5392  int nsols;
5393  int s;
5394  BITARRAY rrowsincut;
5395  int nrrowsincut;
5396  SCIP_Real* weights;
5397  int nrowsincut;
5398  SCIP_Bool cutoff = FALSE;
5399
5400  assert(scip != NULL);
5402  assert(lpdata != NULL);
5403  assert(mod2data != NULL);
5404  assert(maxsepacuts >= 0);
5405  assert(maxcuts >= 0);
5406  assert(nsepacuts != NULL);
5407  assert(nzerohalfcuts != NULL);
5408  assert(zerohalfcuts != NULL);
5409  assert(*nsepacuts <= *nzerohalfcuts);
5410  assert(varsolvals != NULL);
5411  assert(result != NULL);
5412
5413  assert(mod2data->relatedsubproblem != NULL);
5414  assert(mod2data->rows != NULL);
5415  assert(mod2data->rowaggregations != NULL);
5416  assert(mod2data->rhs != NULL);
5417  assert(mod2data->slacks != NULL);
5418  assert(mod2data->fracsol != NULL);
5419  assert(mod2data->rowsind != NULL);
5420  assert(mod2data->colsind != NULL);
5421
5422
5423  /* check if( A mod 2, b mod 2) is empty */
5424  if( mod2data->nrows == 0 || mod2data->nrowsind == 0 )
5425  return SCIP_OKAY;
5426
5427  /* check if enough cuts have been found */
5428  if( *nsepacuts >= maxsepacuts || *nzerohalfcuts >= maxcuts )
5429  return SCIP_OKAY;
5430
5431  /* allocate temporary memory for subscip data structure */
5432  SCIP_CALL(ZerohalfAuxIPDataCreate(scip, &auxipdata));
5433
5434  /* create subscip */
5435  SCIP_CALL(createSubscip(scip, sepadata, lpdata, mod2data, auxipdata, setnodelimit));
5436
5437  /* abort if subscip was not created */
5438  if( auxipdata->subscip == NULL )
5439  {
5440  SCIP_CALL(ZerohalfAuxIPDataFree(scip, &auxipdata));
5441  return SCIP_OKAY;
5442  }
5443
5444  /* solve subscip and get solutions yielding a zerohalf cut with violation >= minviolation */
5445  sols = NULL;
5446  nsols = 0;
5447  SCIP_CALL( solveSubscip(scip, sepadata, mod2data, auxipdata, &sols, &nsols) );
5448
5449  /* process solutions */
5450  for( s = 0; s < nsols ; ++s)
5451  {
5452  /* check if enough cuts have been found */
5453  if( *nsepacuts >= maxsepacuts || *nzerohalfcuts >= maxcuts )
5454  break;
5455
5456  /* determine rrows of the related subproblem that have to be combined */
5457  rrowsincut = NULL;
5458  SCIP_CALL( getBitarrayOfSelectedRows(scip, mod2data, auxipdata, sols[s], &rrowsincut, &nrrowsincut) );
5459  assert(nrrowsincut > 0);
5460
5461  /* calculate rows zerohalf weightvector */
5462  weights = NULL;
5463  SCIP_CALL( getZerohalfWeightvectorFromSelectedRowsBitarray(scip, sepadata, lpdata, mod2data, rrowsincut, &weights, &nrowsincut) );
5464  if ( weights == NULL )
5465  {
5466  SCIPfreeBlockMemoryArray(scip, &rrowsincut, mod2data->rowaggregationsbitarraysize);
5467  continue;
5468  }
5469  assert(nrowsincut > 0);
5470
5471 #ifdef SCIP_DEBUG
5472  SCIP_CALL( debugPrintLPRowsAndCols(scip, lpdata) );
5473  SCIPdebugMsg(scip, "\n");
5474  debugPrintSubLpData(scip, lpdata, mod2data->relatedsubproblem);
5475  debugPrintMod2Data(scip, lpdata, mod2data, TRUE);
5476  SCIPdebugMsg(scip, "\n");
5477  SCIP_CALL( SCIPprintOrigProblem(auxipdata->subscip, NULL, NULL, TRUE) );
5478  SCIPdebugMsg(scip, "\n");
5479  SCIP_CALL( SCIPprintBestSol(auxipdata->subscip, NULL , FALSE) );
5480 #endif
5481
5482  /* create zerohalf cut */
5483  SCIP_CALL(ZerohalfCutDataCreate(scip, &(zerohalfcuts[*nzerohalfcuts]),
5484  mod2data->relatedsubproblem, mod2data, nrrowsincut, nrowsincut, AUXIP));
5486  lpdata, weights, normtype, *nzerohalfcuts, varsolvals, zerohalfcuts[*nzerohalfcuts], &cutoff));
5487
5488  if ( cutoff )
5489  *result = SCIP_CUTOFF;
5490  else
5491  {
5494  (*nzerohalfcuts)++;
5495  }
5496
5497  /* free temporary memory */
5498  SCIPfreeBlockMemoryArray(scip, &weights, lpdata->nrows);
5499  SCIPfreeBlockMemoryArray(scip, &rrowsincut, mod2data->rowaggregationsbitarraysize);
5500
5501  if ( cutoff )
5502  break;
5503  }
5504
5505  /* free temporary memory */
5506  if( sepadata->subscipobjective != 'v' )
5507  {
5508  SCIPfreeBlockMemoryArray(scip, &sols, nsols);
5509  }
5510  SCIP_CALL( ZerohalfAuxIPDataFree(scip, &auxipdata) );
5511
5512  return SCIP_OKAY;
5513 }
5514
5515
5516 /** calculates the inner product of mod2data->row and the LP solution */
5517 static
5519  SCIP* scip, /**< SCIP data structure */
5520  ZEROHALF_MOD2DATA* mod2data, /**< considered (preprocessed) subproblem mod 2 */
5521  BITARRAY row, /**< considered mod2data->row */
5522  SCIP_Real maxinnerproduct, /**< calculation is aborted if innerproduct >= maxinnerproduct */
5523  SCIP_Real* innerproduct /**< pointer to store the inner product */
5524  )
5525 {
5526  int c;
5527  int rcolindex;
5528
5529  assert(scip != NULL);
5530  assert(mod2data != NULL);
5531  assert(row != NULL);
5532  assert(maxinnerproduct >= 0);
5533  assert(innerproduct != NULL);
5534
5535  *innerproduct = 0.0;
5536
5537  /* check if( A mod 2, b mod 2) is empty */
5538  if( mod2data->nrows == 0 || mod2data->nrowsind == 0
5539  || mod2data->ncolsind == 0 )
5540  return SCIP_OKAY;
5541
5542  /* calculate the inner product of rows[rowsindex] and fracsol */
5543  for( c = 0 ; c < mod2data->ncolsind ; ++c)
5544  {
5545  rcolindex = mod2data->colsind[c];
5546  if( BITARRAYBITISSET(row, rcolindex) ) /*lint !e701*/
5547  *innerproduct += mod2data->fracsol[rcolindex];
5548  if( SCIPisGT(scip, *innerproduct, maxinnerproduct) )
5549  break;
5550  }
5551
5552  return SCIP_OKAY;
5553 }
5554
5555
5556 /** separate violated zerohalf cuts by enumerating possible row combinations. (heuristic; polynomial time) */
5557 static
5559  SCIP* scip, /**< SCIP data structure */
5560  SCIP_SEPA* sepa, /**< separator */
5562  ZEROHALF_LPDATA* lpdata, /**< data of current LP relaxation */
5563  ZEROHALF_MOD2DATA* mod2data, /**< considered (preprocessed) subproblem mod 2 */
5564  char normtype, /**< SCIP normtype */
5565  int maxsepacuts, /**< maximal number of zerohalf cuts separated per separation round */
5566  int maxcuts, /**< maximal number of zerohalf cuts found per separation round (incl. ineff. cuts) */
5567  int* nsepacuts, /**< pointer to store current number of separated zerohalf cuts */
5568  int* nzerohalfcuts, /**< pointer to store current number of found zerohalf cuts */
5569  ZEROHALF_CUTDATA** zerohalfcuts, /**< array to store found zerohalf cuts */
5570  SCIP_Real** varsolvals, /**< dense variable LP solution vector */
5571  SCIP_RESULT* result, /**< pointer to SCIP result value of separation */
5572  int maxncombinedrows /**< maximal number of combined rows; currently only 1 or 2 is supported */
5573  )
5574 {
5575  int first;
5576  int last;
5577  int temp;
5578  SCIP_Real minslackoddrhsrows;
5579  int ncombinedrows;
5580  int noddrhsrows;
5581  int r1;
5582  int r2;
5583  BITARRAY combinedrow;
5584  BITARRAY rrowsincut;
5585  SCIP_Real roundingdownweakening;
5586  SCIP_Real slack1;
5587  SCIP_Real slack2;
5588  int i;
5589  int j;
5590  SCIP_Real* weights;
5591  int nrowsincut;
5592  SCIP_Bool cutoff = FALSE;
5593
5594  assert(scip != NULL);
5596  assert(lpdata != NULL);
5597  assert(mod2data != NULL);
5598  assert(maxsepacuts >= 0);
5599  assert(maxcuts >= 0);
5600  assert(nsepacuts != NULL);
5601  assert(nzerohalfcuts != NULL);
5602  assert(zerohalfcuts != NULL);
5603  assert(*nsepacuts <= *nzerohalfcuts);
5604  assert(varsolvals != NULL);
5605  assert(result != NULL);
5606  assert(maxncombinedrows >= 1);
5607
5608  assert(mod2data->relatedsubproblem != NULL);
5609  assert(mod2data->rows != NULL);
5610  assert(mod2data->rowaggregations != NULL);
5611  assert(mod2data->rhs != NULL);
5612  assert(mod2data->slacks != NULL);
5613  assert(mod2data->fracsol != NULL);
5614  assert(mod2data->rowsind != NULL);
5615  assert(mod2data->colsind != NULL);
5616
5617
5618  /* check if( A mod 2, b mod 2) is empty */
5619  if( mod2data->nrows == 0 || mod2data->nrowsind == 0 )
5620  return SCIP_OKAY;
5621
5622  /* check if enough cuts have been found */
5623  if( *nsepacuts >= maxsepacuts || *nzerohalfcuts >= maxcuts )
5624  return SCIP_OKAY;
5625
5626
5627  /* partition rows into odd-rhs-rows and even-rhs-rows */
5628  first = 0;
5629  last = mod2data->nrowsind - 1;
5630  while(first < last)
5631  {
5632  if( !mod2data->rhs[mod2data->rowsind[first]] )
5633  {
5634  temp = mod2data->rowsind[first];
5635  mod2data->rowsind[first] = mod2data->rowsind[last];
5636  mod2data->rowsind[last] = temp;
5637  --last;
5638  }
5639  else
5640  ++first;
5641  }
5642  noddrhsrows = first + (mod2data->rhs[mod2data->rowsind[first]] ? 1 : 0);
5643
5644  /* check if odd rows exists */
5645  if( noddrhsrows == 0 )
5646  return SCIP_OKAY;
5647
5648  /* allocate and initialize temporary memory for calculating the symmetric difference */
5649  combinedrow = NULL;
5650  rrowsincut = NULL;
5651  if( maxncombinedrows > 1 )
5652  {
5653  SCIP_CALL(SCIPallocBufferArray(scip, &combinedrow, mod2data->rowsbitarraysize));
5654  SCIP_CALL(SCIPallocBufferArray(scip, &rrowsincut, mod2data->rowaggregationsbitarraysize));
5655  BITARRAYCLEAR(rrowsincut, mod2data->rowaggregationsbitarraysize);
5656  }
5657
5658  /* sort each partition by nondecreasing slacks */
5659  assert(noddrhsrows >= 0);
5660  SCIPsortInd( mod2data->rowsind , compRealNonDecreasing , (void*) mod2data->slacks , noddrhsrows );
5661
5662  if( noddrhsrows < mod2data->nrowsind )
5663  {
5664  SCIPsortInd( mod2data->rowsind + noddrhsrows , compRealNonDecreasing , (void*) mod2data->slacks , mod2data->nrowsind - noddrhsrows );
5665  }
5666
5667  minslackoddrhsrows = mod2data->slacks[mod2data->rowsind[0]];
5668
5669  if( SCIPisLE(scip, minslackoddrhsrows, sepadata->maxslack) )
5670  {
5671  for( ncombinedrows = 1 ; ncombinedrows <= maxncombinedrows ; ++ncombinedrows )
5672  {
5673  switch( ncombinedrows )
5674  {
5675  case 1:
5676  /* check all rows r1 with rhs(r1) odd */
5677  for( i = 0 ; i < noddrhsrows ; ++i)
5678  {
5679  r1 = mod2data->rowsind[i];
5680  assert(mod2data->rhs[r1]);
5681  if( *nzerohalfcuts == maxcuts || *nsepacuts == maxsepacuts )
5682  break;
5683  slack1 = mod2data->slacks[r1];
5684  if( SCIPisGT(scip, slack1, sepadata->maxslack) )
5685  break; /* because rowsind is sorted */
5686  SCIP_CALL(calcInnerProductOfRowAndFracsol(scip, mod2data, mod2data->rows[r1],
5688  if( SCIPisLE(scip, roundingdownweakening + slack1, sepadata->maxslack) )
5689  {
5690  /* a violated zerohalf cut has been found */
5691
5692  /* calculate rows zerohalf weightvector */
5693  weights = NULL;
5695  mod2data, mod2data->rowaggregations[r1], &weights, &nrowsincut));
5696  if( weights == NULL )
5697  continue;
5698  assert(nrowsincut > 0);
5699
5700  /* create zerohalf cut */
5701  SCIP_CALL(ZerohalfCutDataCreate(scip, &(zerohalfcuts[*nzerohalfcuts]),
5702  mod2data->relatedsubproblem, mod2data, 2, nrowsincut, HEURISTICSENUM));
5704  lpdata, weights, normtype, *nzerohalfcuts, varsolvals, zerohalfcuts[*nzerohalfcuts], &cutoff));
5705
5706  if ( cutoff )
5707  *result = SCIP_CUTOFF;
5708  else
5709  {
5712  (*nzerohalfcuts)++;
5713  }
5714
5715  /* free temporary memory */
5716  SCIPfreeBlockMemoryArray(scip, &weights, lpdata->nrows);
5717  }
5718  }
5719  break;
5720  case 2:
5721  assert(combinedrow != NULL);
5722  assert(rrowsincut != NULL);
5723  if( noddrhsrows == mod2data->nrowsind )
5724  break;
5725  if( mod2data->nrowsind < 2 )
5726  break;
5727
5728  /* check all pairs (r1,r2) with rhs(r1) odd and rhs(r1) even */
5729  for( i = 0 ; i < noddrhsrows ; ++i)
5730  {
5731  r1 = mod2data->rowsind[i];
5732  assert(mod2data->rhs[r1]);
5733  if( *nzerohalfcuts == maxcuts || *nsepacuts == maxsepacuts )
5734  break;
5735  slack1 = mod2data->slacks[r1];
5736  if( SCIPisGT(scip, slack1, sepadata->maxslack) )
5737  break; /* because rowsind_odd is sorted */
5738
5739  for( j = noddrhsrows ; j < mod2data->ncolsind ; ++j)
5740  {
5741  r2 = mod2data->rowsind[j];
5742  assert(!mod2data->rhs[r2]);
5743  if( *nzerohalfcuts == maxcuts || *nsepacuts == maxsepacuts )
5744  break;
5745  slack2 = mod2data->slacks[r2];
5746  if( SCIPisGT(scip, slack1 + slack2, sepadata->maxslack) )
5747  break; /* because rowsind_even is sorted */
5748  BMScopyMemoryArray(combinedrow, mod2data->rows[r1], mod2data->rowsbitarraysize);
5749  BITARRAYSXOR(mod2data->rows[r2], combinedrow, mod2data->rowsbitarraysize);
5750  SCIP_CALL(calcInnerProductOfRowAndFracsol(scip, mod2data, combinedrow,
5752  if( SCIPisLE(scip, roundingdownweakening + slack1 + slack2, sepadata->maxslack) )
5753  {
5754  /* a violated zerohalf cut has been found */
5755
5756  /* determine rrows of the related subproblem that have to be combined */
5757  BMScopyMemoryArray(rrowsincut, mod2data->rowaggregations[r1], mod2data->rowaggregationsbitarraysize);
5758  BITARRAYSXOR(mod2data->rowaggregations[r2], rrowsincut, mod2data->rowaggregationsbitarraysize);
5759
5760  /* calculate rows zerohalf weightvector */
5761  weights = NULL;
5763  mod2data, rrowsincut, &weights, &nrowsincut));
5764  if ( weights == NULL )
5765  {
5766  continue;
5767  }
5768  assert(nrowsincut > 0);
5769
5770  /* create zerohalf cut */
5771  SCIP_CALL(ZerohalfCutDataCreate(scip, &(zerohalfcuts[*nzerohalfcuts]),
5772  mod2data->relatedsubproblem, mod2data, 2, nrowsincut, HEURISTICSENUM));
5774  lpdata, weights, normtype, *nzerohalfcuts, varsolvals, zerohalfcuts[*nzerohalfcuts], &cutoff));
5775
5776  if ( cutoff )
5777  *result = SCIP_CUTOFF;
5778  else
5779  {
5782  (*nzerohalfcuts)++;
5783  }
5784
5785  /* free temporary memory */
5786  SCIPfreeBlockMemoryArray(scip, &weights, lpdata->nrows);
5787  }
5788  }
5789  }
5790  break;
5791  default:
5792  SCIPerrorMessage("invalid ncombinedrows '%d'\n", ncombinedrows);
5793  return SCIP_INVALIDDATA;
5794  }
5795
5796  if ( cutoff )
5797  break;
5798  }
5799  }
5800
5801  /* free temporary memory */
5802  SCIPfreeBufferArrayNull(scip, &rrowsincut);
5803  SCIPfreeBufferArrayNull(scip, &combinedrow);
5804
5805  return SCIP_OKAY;
5806 }
5807
5808
5809 #if 0
5810 /** prints a node of the auxiliary graph */
5811 static
5812 void debugPrintAuxGraphNode(
5813  ZEROHALF_AUXGRAPH_NODE* node /**< node to be printed */
5814  )
5815 {
5816  int i;
5817
5818  assert(node != NULL);
5819
5820  SCIPdebugMsg(scip, "\nnode: %p\n", node);
5821  for( i = 0 ; i < node->nneighbors ; ++i)
5822  {
5823  SCIPdebugMsg(scip, " neighbor %4d: %p weight: %6f rrow: %4d\n", i, node->neighbors[i], node->edgeweights[i], node->relatedrows[i]);
5824  }
5825  SCIPdebugMsg(scip, " nneighbors: %d distance: %6f previous: %p\n", node->nneighbors, node->distance, node->previous);
5826 }
5827 #endif
5828
5829
5830 /** adds an edge (and its "copy" w.r.t. the node copies) to the auxiliary graph */
5831 static
5833  SCIP* scip, /**< SCIP data structure */
5834  ZEROHALF_AUXGRAPH* graph, /**< auxiliary graph */
5835  int node1index, /**< start node of edge */
5836  int node2index, /**< end node of edge */
5837  SCIP_Bool isodd, /**< is the rhs value of the corresponding mod2data->row odd? */
5838  SCIP_Real weight, /**< weight of the edge */
5839  int relatedrow /**< corresponding mod2data->row */
5840  )
5842  ZEROHALF_AUXGRAPH_NODE* node1;
5843  ZEROHALF_AUXGRAPH_NODE* node2;
5844  ZEROHALF_AUXGRAPH_NODE* node1copy;
5845  ZEROHALF_AUXGRAPH_NODE* node2copy;
5846  int n1;
5847  int n2;
5848
5849  int maxnumberofneighbors;
5850
5851  assert(scip != NULL);
5852  assert(graph != NULL);
5853  assert(node1index >= 0);
5854  assert(node1index < graph->nnodes);
5855  assert(node2index >= 0);
5856  assert(node2index < graph->nnodes);
5857  assert(!SCIPisNegative(scip, weight));
5858
5859  maxnumberofneighbors = 2 * graph->nnodes - 2;
5860
5861  if( isodd )
5862  {
5863  node1 = graph->nodes[node1index];
5864  node2 = graph->nodecopies[node2index];
5865  node1copy = graph->nodecopies[node1index];
5866  node2copy = graph->nodes[node2index];
5867  }
5868  else
5869  {
5870  node1 = graph->nodes[node1index];
5871  node2 = graph->nodes[node2index];
5872  node1copy = graph->nodecopies[node1index];
5873  node2copy = graph->nodecopies[node2index];
5874  }
5875
5876  if( node1->nneighbors == 0 )
5877  {
5878  assert(maxnumberofneighbors > 0);
5879  SCIP_CALL(SCIPallocBlockMemoryArray(scip, &(node1->neighbors), maxnumberofneighbors));
5880  SCIP_CALL(SCIPallocBlockMemoryArray(scip, &(node1->edgeweights), maxnumberofneighbors));
5881  SCIP_CALL(SCIPallocBlockMemoryArray(scip, &(node1->relatedrows), maxnumberofneighbors));
5882
5883  SCIP_CALL(SCIPallocBlockMemoryArray(scip, &(node1copy->neighbors), maxnumberofneighbors));
5884  SCIP_CALL(SCIPallocBlockMemoryArray(scip, &(node1copy->edgeweights), maxnumberofneighbors));
5885  SCIP_CALL(SCIPallocBlockMemoryArray(scip, &(node1copy->relatedrows), maxnumberofneighbors));
5886  }
5887
5888  if( node2->nneighbors == 0 )
5889  {
5890  assert(maxnumberofneighbors > 0);
5891  SCIP_CALL(SCIPallocBlockMemoryArray(scip, &(node2->neighbors), maxnumberofneighbors));
5892  SCIP_CALL(SCIPallocBlockMemoryArray(scip, &(node2->edgeweights), maxnumberofneighbors));
5893  SCIP_CALL(SCIPallocBlockMemoryArray(scip, &(node2->relatedrows), maxnumberofneighbors));
5894
5895  SCIP_CALL(SCIPallocBlockMemoryArray(scip, &(node2copy->neighbors), maxnumberofneighbors));
5896  SCIP_CALL(SCIPallocBlockMemoryArray(scip, &(node2copy->edgeweights), maxnumberofneighbors));
5897  SCIP_CALL(SCIPallocBlockMemoryArray(scip, &(node2copy->relatedrows), maxnumberofneighbors));
5898  }
5899
5900  n2 = node2->nneighbors;
5901  for( n1 = 0 ; n1 < node1->nneighbors ; ++n1)
5902  if( node1->neighbors[n1] == node2 )
5903  break;
5904  if( n1 < node1->nneighbors)
5905  for( n2 = 0 ; n2 < node2->nneighbors ; ++n2)
5906  if( node2->neighbors[n2] == node1 )
5907  break;
5908  if( n1 < node1->nneighbors )
5909  {
5910  /* if node2 is neighbor of node1, then node1 is neighbor of node2 */
5911  assert(node1->neighbors[n1] == node2);
5912  assert(n2 < node2->nneighbors);
5913  assert(node2->neighbors[n2] == node1);
5914  assert(node1->edgeweights[n1] == node2->edgeweights[n2]); /*lint !e777*/
5915  }
5916
5917  if( n1 == node1->nneighbors || SCIPisLT(scip, weight, node1->edgeweights[n1]) )
5918  {
5919  node1->neighbors[n1] = node2;
5920  node1->edgeweights[n1] = weight;
5921  node1->relatedrows[n1] = relatedrow;