1 /* * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * */
2 /* */
3 /* This file is part of the program and library */
4 /* SCIP --- Solving Constraint Integer Programs */
5 /* */
6 /* Copyright (C) 2002-2021 Konrad-Zuse-Zentrum */
7 /* fuer Informationstechnik Berlin */
8 /* */
9 /* SCIP is distributed under the terms of the ZIB Academic License. */
10 /* */
11 /* You should have received a copy of the ZIB Academic License */
12 /* along with SCIP; see the file COPYING. If not email to scip@zib.de. */
13 /* */
14 /* * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * */
15
16 /**@file matrix.c
17 * @ingroup OTHER_CFILES
18 * @brief methods for MIP matrix data structure
19 * @author Dieter Weninger
20 * @author Gerald Gamrath
21 *
22 * The MIP matrix is organized as sparse data structure in row and
23 * and column major format.
24 *
25 * @todo disregard relaxation-only variables in lock check and don't copy them to the matrix
26 */
27
28 /*---+----1----+----2----+----3----+----4----+----5----+----6----+----7----+----8----+----9----+----0----+----1----+----2*/
29
30 #include "blockmemshell/memory.h"
31 #include "scip/cons_knapsack.h"
32 #include "scip/cons_linear.h"
33 #include "scip/cons_logicor.h"
34 #include "scip/cons_setppc.h"
35 #include "scip/cons_varbound.h"
36 #include "scip/pub_matrix.h"
37 #include "scip/pub_cons.h"
38 #include "scip/pub_message.h"
39 #include "scip/pub_misc_sort.h"
40 #include "scip/pub_var.h"
41 #include "scip/scip_cons.h"
42 #include "scip/scip_general.h"
43 #include "scip/scip_mem.h"
44 #include "scip/scip_message.h"
45 #include "scip/scip_numerics.h"
46 #include "scip/scip_pricer.h"
47 #include "scip/scip_prob.h"
48 #include "scip/scip_var.h"
49 #include "scip/struct_matrix.h"
50 #include <string.h>
51
52 /*
53 * private functions
54 */
55
56 /** transforms given variables, scalars and constant to the corresponding active variables, scalars and constant */
57 static
getActiveVariables(SCIP * scip,SCIP_VAR *** vars,SCIP_Real ** scalars,int * nvars,SCIP_Real * constant)58 SCIP_RETCODE getActiveVariables(
59 SCIP* scip, /**< SCIP instance */
60 SCIP_VAR*** vars, /**< vars array to get active variables for */
61 SCIP_Real** scalars, /**< scalars a_1, ..., a_n in linear sum a_1*x_1 + ... + a_n*x_n + c */
62 int* nvars, /**< pointer to number of variables and values in vars and vals array */
63 SCIP_Real* constant /**< pointer to constant c in linear sum a_1*x_1 + ... + a_n*x_n + c */
64 )
65 {
66 int requiredsize;
67
68 assert(scip != NULL);
69 assert(vars != NULL);
70 assert(scalars != NULL);
71 assert(*vars != NULL);
72 assert(*scalars != NULL);
73 assert(nvars != NULL);
74 assert(constant != NULL);
75
76 SCIP_CALL( SCIPgetProbvarLinearSum(scip, *vars, *scalars, nvars, *nvars, constant, &requiredsize, TRUE) );
77
78 if( requiredsize > *nvars )
79 {
80 SCIP_CALL( SCIPreallocBufferArray(scip, vars, requiredsize) );
81 SCIP_CALL( SCIPreallocBufferArray(scip, scalars, requiredsize) );
82
83 /* call function a second time with enough memory */
84 SCIP_CALL( SCIPgetProbvarLinearSum(scip, *vars, *scalars, nvars, requiredsize, constant, &requiredsize, TRUE) );
85 assert(requiredsize <= *nvars);
86 }
87
88 return SCIP_OKAY;
89 }
90
91 /** add one row to the constraint matrix */
92 static
addRow(SCIP * scip,SCIP_MATRIX * matrix,SCIP_VAR ** vars,SCIP_Real * vals,int nvars,SCIP_Real lhs,SCIP_Real rhs,int maxnnonzsmem,SCIP_Bool * rowadded)93 SCIP_RETCODE addRow(
94 SCIP* scip, /**< SCIP data structure */
95 SCIP_MATRIX* matrix, /**< constraint matrix */
96 SCIP_VAR** vars, /**< variables of this row */
97 SCIP_Real* vals, /**< coefficients of this row */
98 int nvars, /**< number of variables of this row */
99 SCIP_Real lhs, /**< left hand side */
100 SCIP_Real rhs, /**< right hand side */
101 int maxnnonzsmem, /**< maximal number of fillable elements */
102 SCIP_Bool* rowadded /**< flag indicating if constraint was added to matrix */
103 )
104 {
105 int j;
106 int probindex;
107 int rowidx;
108 SCIP_Real factor;
109 SCIP_Bool rangedorequality;
110
111 assert(vars != NULL);
112 assert(vals != NULL);
113
114 rowidx = matrix->nrows;
115 rangedorequality = FALSE;
116
117 if( SCIPisInfinity(scip, -lhs) )
118 {
119 factor = -1.0;
120 matrix->lhs[rowidx] = -rhs;
121 matrix->rhs[rowidx] = SCIPinfinity(scip);
122 matrix->isrhsinfinite[rowidx] = TRUE;
123 }
124 else
125 {
126 factor = 1.0;
127 matrix->lhs[rowidx] = lhs;
128 matrix->rhs[rowidx] = rhs;
129 matrix->isrhsinfinite[rowidx] = SCIPisInfinity(scip, matrix->rhs[rowidx]);
130
131 if( !SCIPisInfinity(scip, rhs) )
132 rangedorequality = TRUE;
133 }
134
135 if(SCIPisInfinity(scip, -matrix->lhs[rowidx]))
136 {
137 /* ignore redundant constraint */
138 *rowadded = FALSE;
139 return SCIP_OKAY;
140 }
141
142 matrix->rowmatbeg[rowidx] = matrix->nnonzs;
143
144 /* = or ranged */
145 if( rangedorequality )
146 {
147 assert(factor > 0);
148
149 for( j = 0; j < nvars; j++ )
150 {
151 assert(maxnnonzsmem > matrix->nnonzs);
152
153 /* ignore variables with very small coefficients */
154 if( SCIPisZero(scip, vals[j]) )
155 continue;
156
157 matrix->rowmatval[matrix->nnonzs] = factor * vals[j];
158 probindex = SCIPvarGetProbindex(vars[j]);
159 assert(matrix->vars[probindex] == vars[j]);
160
161 matrix->nuplocks[probindex]++;
162 matrix->ndownlocks[probindex]++;
163
164 assert(0 <= probindex && probindex < matrix->ncols);
165 matrix->rowmatind[matrix->nnonzs] = probindex;
166
167 (matrix->nnonzs)++;
168 }
169 }
170 /* >= or <= */
171 else
172 {
173 for( j = 0; j < nvars; j++ )
174 {
175 assert(maxnnonzsmem > matrix->nnonzs);
176
177 /* ignore variables with very small coefficients */
178 if( SCIPisZero(scip, vals[j]) )
179 continue;
180
181 /* due to the factor, <= constraints will be transfered to >= */
182 matrix->rowmatval[matrix->nnonzs] = factor * vals[j];
183 probindex = SCIPvarGetProbindex(vars[j]);
184 assert(matrix->vars[probindex] == vars[j]);
185
186 if( matrix->rowmatval[matrix->nnonzs] > 0 )
187 matrix->ndownlocks[probindex]++;
188 else
189 {
190 assert(matrix->rowmatval[matrix->nnonzs] < 0);
191 matrix->nuplocks[probindex]++;
192 }
193
194 assert(0 <= probindex && probindex < matrix->ncols);
195 matrix->rowmatind[matrix->nnonzs] = probindex;
196
197 (matrix->nnonzs)++;
198 }
199 }
200
201 matrix->rowmatcnt[rowidx] = matrix->nnonzs - matrix->rowmatbeg[rowidx];
202
203 ++(matrix->nrows);
204 *rowadded = TRUE;
205
206 return SCIP_OKAY;
207 }
208
209 /** add one constraint to matrix */
210 static
addConstraint(SCIP * scip,SCIP_MATRIX * matrix,SCIP_VAR ** vars,SCIP_Real * vals,int nvars,SCIP_Real lhs,SCIP_Real rhs,int maxnnonzsmem,SCIP_Bool * rowadded)211 SCIP_RETCODE addConstraint(
212 SCIP* scip, /**< current scip instance */
213 SCIP_MATRIX* matrix, /**< constraint matrix */
214 SCIP_VAR** vars, /**< variables of this constraint */
215 SCIP_Real* vals, /**< variable coefficients of this constraint */
216 int nvars, /**< number of variables */
217 SCIP_Real lhs, /**< left hand side */
218 SCIP_Real rhs, /**< right hand side */
219 int maxnnonzsmem, /**< maximal number of fillable elements */
220 SCIP_Bool* rowadded /**< flag indicating of row was added to matrix */
221 )
222 {
223 SCIP_VAR** activevars;
224 SCIP_Real* activevals;
225 SCIP_Real activeconstant;
226 int nactivevars;
227 int v;
228
229 assert(scip != NULL);
230 assert(matrix != NULL);
231 assert(vars != NULL || nvars == 0);
232 assert(SCIPisLE(scip, lhs, rhs));
233 assert(rowadded != NULL);
234
235 *rowadded = FALSE;
236
237 /* constraint is redundant */
238 if( SCIPisInfinity(scip, -lhs) && SCIPisInfinity(scip, rhs) )
239 return SCIP_OKAY;
240
241 /* we do not add empty constraints to the matrix */
242 if( nvars == 0 )
243 return SCIP_OKAY;
244
245 activevars = NULL;
246 activevals = NULL;
247 nactivevars = nvars;
248 activeconstant = 0.0;
249
250 /* duplicate variable and value array */
251 SCIP_CALL( SCIPduplicateBufferArray(scip, &activevars, vars, nactivevars ) );
252 if( vals != NULL )
253 {
254 SCIP_CALL( SCIPduplicateBufferArray(scip, &activevals, vals, nactivevars ) );
255 }
256 else
257 {
258 SCIP_CALL( SCIPallocBufferArray(scip, &activevals, nactivevars) );
259
260 for( v = 0; v < nactivevars; v++ )
261 activevals[v] = 1.0;
262 }
263
264 /* retransform given variables to active variables */
265 SCIP_CALL( getActiveVariables(scip, &activevars, &activevals, &nactivevars, &activeconstant) );
266
267 /* adapt left and right hand side */
268 if( !SCIPisInfinity(scip, -lhs) )
269 lhs -= activeconstant;
270 if( !SCIPisInfinity(scip, rhs) )
271 rhs -= activeconstant;
272
273 /* add single row to matrix */
274 if( nactivevars > 0 )
275 {
276 SCIP_CALL( addRow(scip, matrix, activevars, activevals, nactivevars, lhs, rhs, maxnnonzsmem, rowadded) );
277 }
278
279 /* free buffer arrays */
280 SCIPfreeBufferArray(scip, &activevals);
281 SCIPfreeBufferArray(scip, &activevars);
282
283 return SCIP_OKAY;
284 }
285
286 /** transform row major format into column major format */
287 static
setColumnMajorFormat(SCIP * scip,SCIP_MATRIX * matrix)288 SCIP_RETCODE setColumnMajorFormat(
289 SCIP* scip, /**< current scip instance */
290 SCIP_MATRIX* matrix /**< constraint matrix */
291 )
292 {
293 int colidx;
294 int i;
295 int* rowpnt;
296 int* rowend;
297 SCIP_Real* valpnt;
298 int* fillidx;
299
300 assert(scip != NULL);
301 assert(matrix != NULL);
302 assert(matrix->colmatval != NULL);
303 assert(matrix->colmatind != NULL);
304 assert(matrix->colmatbeg != NULL);
305 assert(matrix->colmatcnt != NULL);
306 assert(matrix->rowmatval != NULL);
307 assert(matrix->rowmatind != NULL);
308 assert(matrix->rowmatbeg != NULL);
309 assert(matrix->rowmatcnt != NULL);
310
311 SCIP_CALL( SCIPallocBufferArray(scip, &fillidx, matrix->ncols) );
312 BMSclearMemoryArray(fillidx, matrix->ncols);
313 BMSclearMemoryArray(matrix->colmatcnt, matrix->ncols);
314
315 for( i = 0; i < matrix->nrows; i++ )
316 {
317 rowpnt = matrix->rowmatind + matrix->rowmatbeg[i];
318 rowend = rowpnt + matrix->rowmatcnt[i];
319 for( ; rowpnt < rowend; rowpnt++ )
320 {
321 colidx = *rowpnt;
322 (matrix->colmatcnt[colidx])++;
323 }
324 }
325
326 matrix->colmatbeg[0] = 0;
327 for( i = 0; i < matrix->ncols-1; i++ )
328 {
329 matrix->colmatbeg[i+1] = matrix->colmatbeg[i] + matrix->colmatcnt[i];
330 }
331
332 for( i = 0; i < matrix->nrows; i++ )
333 {
334 rowpnt = matrix->rowmatind + matrix->rowmatbeg[i];
335 rowend = rowpnt + matrix->rowmatcnt[i];
336 valpnt = matrix->rowmatval + matrix->rowmatbeg[i];
337
338 for( ; rowpnt < rowend; rowpnt++, valpnt++ )
339 {
340 assert(*rowpnt < matrix->ncols);
341 colidx = *rowpnt;
342 matrix->colmatval[matrix->colmatbeg[colidx] + fillidx[colidx]] = *valpnt;
343 matrix->colmatind[matrix->colmatbeg[colidx] + fillidx[colidx]] = i;
344 fillidx[colidx]++;
345 }
346 }
347
348 SCIPfreeBufferArray(scip, &fillidx);
349
350 return SCIP_OKAY;
351 }
352
353 /** calculate min/max activity per row */
354 static
calcActivityBounds(SCIP * scip,SCIP_MATRIX * matrix)355 SCIP_RETCODE calcActivityBounds(
356 SCIP* scip, /**< current scip instance */
357 SCIP_MATRIX* matrix /**< constraint matrix */
358 )
359 {
360 SCIP_Real val;
361 int* rowpnt;
362 int* rowend;
363 SCIP_Real* valpnt;
364 int col;
365 int row;
366
367 assert(scip != NULL);
368 assert(matrix != NULL);
369
370 for( row = 0; row < matrix->nrows; row++ )
371 {
372 matrix->minactivity[row] = 0;
373 matrix->maxactivity[row] = 0;
374 matrix->minactivityneginf[row] = 0;
375 matrix->minactivityposinf[row] = 0;
376 matrix->maxactivityneginf[row] = 0;
377 matrix->maxactivityposinf[row] = 0;
378
379 rowpnt = matrix->rowmatind + matrix->rowmatbeg[row];
380 rowend = rowpnt + matrix->rowmatcnt[row];
381 valpnt = matrix->rowmatval + matrix->rowmatbeg[row];
382
383 for( ; rowpnt < rowend; rowpnt++, valpnt++ )
384 {
385 /* get column index */
386 col = *rowpnt;
387
388 /* get variable coefficient */
389 val = *valpnt;
390 assert(!SCIPisZero(scip, val));
391
392 assert(matrix->ncols > col);
393
394 assert(!SCIPisInfinity(scip, matrix->lb[col]));
395 assert(!SCIPisInfinity(scip, -matrix->ub[col]));
396
397 /* positive coefficient */
398 if( val > 0.0 )
399 {
400 if( SCIPisInfinity(scip, matrix->ub[col]) )
401 matrix->maxactivityposinf[row]++;
402 else
403 matrix->maxactivity[row] += val * matrix->ub[col];
404
405 if( SCIPisInfinity(scip, -matrix->lb[col]) )
406 matrix->minactivityneginf[row]++;
407 else
408 matrix->minactivity[row] += val * matrix->lb[col];
409 }
410 /* negative coefficient */
411 else
412 {
413 if( SCIPisInfinity(scip, -matrix->lb[col]) )
414 matrix->maxactivityneginf[row]++;
415 else
416 matrix->maxactivity[row] += val * matrix->lb[col];
417
418 if( SCIPisInfinity(scip, matrix->ub[col]) )
419 matrix->minactivityposinf[row]++;
420 else
421 matrix->minactivity[row] += val * matrix->ub[col];
422 }
423 }
424
425 /* consider infinite bound contributions for the activities */
426 if( matrix->maxactivityneginf[row] + matrix->maxactivityposinf[row] > 0 )
427 matrix->maxactivity[row] = SCIPinfinity(scip);
428
429 if( matrix->minactivityneginf[row] + matrix->minactivityposinf[row] > 0 )
430 matrix->minactivity[row] = -SCIPinfinity(scip);
431 }
432
433 return SCIP_OKAY;
434 }
435
436 /*
437 * public functions
438 */
439
440 /** initialize matrix by copying all check constraints
441 *
442 * @note Completeness is checked by testing whether all check constraints are from a list of linear constraint handlers
443 * that can be represented.
444 */
SCIPmatrixCreate(SCIP * scip,SCIP_MATRIX ** matrixptr,SCIP_Bool onlyifcomplete,SCIP_Bool * initialized,SCIP_Bool * complete,SCIP_Bool * infeasible,int * naddconss,int * ndelconss,int * nchgcoefs,int * nchgbds,int * nfixedvars)445 SCIP_RETCODE SCIPmatrixCreate(
446 SCIP* scip, /**< current scip instance */
447 SCIP_MATRIX** matrixptr, /**< pointer to constraint matrix object to be initialized */
448 SCIP_Bool onlyifcomplete, /**< should matrix creation be skipped if matrix will not be complete? */
449 SCIP_Bool* initialized, /**< was the initialization successful? */
450 SCIP_Bool* complete, /**< are all constraint represented within the matrix? */
451 SCIP_Bool* infeasible, /**< pointer to return whether problem was detected to be infeasible during matrix creation */
452 int* naddconss, /**< pointer to count number of added (linear) constraints during matrix creation */
453 int* ndelconss, /**< pointer to count number of deleted specialized linear constraints during matrix creation */
454 int* nchgcoefs, /**< pointer to count number of changed coefficients during matrix creation */
455 int* nchgbds, /**< pointer to count number of changed bounds during matrix creation */
456 int* nfixedvars /**< pointer to count number of fixed variables during matrix creation */
457 )
458 {
459 SCIP_MATRIX* matrix;
460 SCIP_CONSHDLR** conshdlrs;
461 const char* conshdlrname;
462 SCIP_Bool stopped;
463 SCIP_VAR** vars;
464 SCIP_VAR* var;
465 SCIP_CONS* cons;
466 int nconshdlrs;
467 int nconss;
468 int nconssall;
469 int nnonzstmp;
470 int nvars;
471 int c;
472 int i;
473 int v;
474 int cnt;
475
476 nnonzstmp = 0;
477
478 assert(scip != NULL);
479 assert(matrixptr != NULL);
480 assert(initialized != NULL);
481 assert(complete != NULL);
482
483 *initialized = FALSE;
484 *complete = FALSE;
485 *infeasible = FALSE;
486
487 /* return if no variables or constraints are present */
488 if( SCIPgetNVars(scip) == 0 || SCIPgetNConss(scip) == 0 )
489 return SCIP_OKAY;
490
491 /* return if pricers are present and the matrix should only be built when complete */
492 if( onlyifcomplete && SCIPgetNActivePricers(scip) != 0 )
493 return SCIP_OKAY;
494
495 /* loop over all constraint handlers and collect the number of checked constraints */
496 nconshdlrs = SCIPgetNConshdlrs(scip);
497 conshdlrs = SCIPgetConshdlrs(scip);
498 nconss = 0;
499 nconssall = 0;
500
501 for( i = 0; i < nconshdlrs; ++i )
502 {
503 int nconshdlrconss;
504
505 nconshdlrconss = SCIPconshdlrGetNCheckConss(conshdlrs[i]);
506
507 if( nconshdlrconss > 0 )
508 {
509 conshdlrname = SCIPconshdlrGetName(conshdlrs[i]);
510
511 if( (strcmp(conshdlrname, "linear") == 0) || (strcmp(conshdlrname, "setppc") == 0)
512 || (strcmp(conshdlrname, "logicor") == 0) || (strcmp(conshdlrname, "knapsack") == 0)
513 || (strcmp(conshdlrname, "varbound") == 0) )
514 {
515 /* increment number of supported constraints */
516 nconss += nconshdlrconss;
517 }
518 /* disabled because some of the presolvers can currently only handle 1-1 row-cons relationships */
519 #ifdef SCIP_DISABLED_CODE
520 else if( strcmp(conshdlrname, "linking") == 0 )
521 {
522 /* the linear representation of linking constraints involves two linear constraints */
523 nconss += 2* nconshdlrconss;
524 }
525 #endif
526 /* increment number of supported and unsupported constraints */
527 nconssall += nconshdlrconss;
528 }
529 }
530
531 /* print warning if we have unsupported constraint types; we only abort the matrix creation process if requested,
532 * because it makes sometimes sense to work on an incomplete matrix as long as the number of interesting variable
533 * uplocks or downlocks of the matrix and scip are the same
534 */
535 if( nconss < nconssall )
536 {
537 SCIPdebugMsg(scip, "Warning: milp matrix not complete!\n");
538 if( onlyifcomplete )
539 return SCIP_OKAY;
540 }
541 else
542 {
543 /* all constraints represented within the matrix */
544 *complete = TRUE;
545 }
546
547 /* do nothing if we have no checked constraints */
548 if( nconss == 0 )
549 return SCIP_OKAY;
550
551 stopped = FALSE;
552
553 /* first, clean up aggregations and fixings in varbound costraints, since this can lead
554 * to boundchanges and the varbound constraint can get downgraded to a linear constraint
555 */
556 SCIP_CALL( SCIPcleanupConssVarbound(scip, TRUE, infeasible, naddconss, ndelconss, nchgbds ) );
557 if( *infeasible )
558 return SCIP_OKAY;
559
560 /* next, clean up aggregations and fixings in setppc costraints, since this can lead
561 * to fixings and the setppc constraint can get downgraded to a linear constraint
562 */
563 SCIP_CALL( SCIPcleanupConssSetppc(scip, TRUE, infeasible, naddconss, ndelconss, nchgcoefs, nfixedvars ) );
564 if( *infeasible )
565 return SCIP_OKAY;
566
567 /* next, clean up aggregations and fixings in logicor costraints, since this cannot lead
568 * to further fixings but the logicor constraint can also get downgraded to a linear constraint
569 */
570 SCIP_CALL( SCIPcleanupConssLogicor(scip, TRUE, naddconss, ndelconss, nchgcoefs) );
571
572 /* finally, clean up aggregations and fixings in knapsack and linear constraints since now no new linaer constraints
573 * can come up due to downgrading and the remaining cleanup methods cannot fix any more variables
574 */
575
576 SCIP_CALL( SCIPcleanupConssKnapsack(scip, TRUE, infeasible) );
577 if( *infeasible )
578 return SCIP_OKAY;
579
580 SCIP_CALL( SCIPcleanupConssLinear(scip, TRUE, infeasible) );
581 if( *infeasible )
582 return SCIP_OKAY;
583
584 vars = SCIPgetVars(scip);
585 nvars = SCIPgetNVars(scip);
586
587 /* approximate number of nonzeros by taking for each variable the number of up- and downlocks;
588 * this counts nonzeros in equalities twice, but can be at most two times as high as the exact number
589 */
590 for( i = nvars - 1; i >= 0; --i )
591 {
592 nnonzstmp += SCIPvarGetNLocksDownType(vars[i], SCIP_LOCKTYPE_MODEL);
593 nnonzstmp += SCIPvarGetNLocksUpType(vars[i], SCIP_LOCKTYPE_MODEL);
594 }
595
596 /* do nothing if we have no entries */
597 if( nnonzstmp == 0 )
598 return SCIP_OKAY;
599
600 /* build the matrix structure */
601 SCIP_CALL( SCIPallocBuffer(scip, matrixptr) );
602 matrix = *matrixptr;
603
604 /* copy vars array and set number of variables */
605 SCIP_CALL( SCIPduplicateBufferArray(scip, &matrix->vars, vars, nvars) );
606 matrix->ncols = nvars;
607
608 matrix->nrows = 0;
609 matrix->nnonzs = 0;
610
611 /* allocate memory */
612 SCIP_CALL( SCIPallocBufferArray(scip, &matrix->colmatval, nnonzstmp) );
613 SCIP_CALL( SCIPallocBufferArray(scip, &matrix->colmatind, nnonzstmp) );
614 SCIP_CALL( SCIPallocBufferArray(scip, &matrix->colmatbeg, matrix->ncols) );
615 SCIP_CALL( SCIPallocBufferArray(scip, &matrix->colmatcnt, matrix->ncols) );
616 SCIP_CALL( SCIPallocBufferArray(scip, &matrix->lb, matrix->ncols) );
617 SCIP_CALL( SCIPallocBufferArray(scip, &matrix->ub, matrix->ncols) );
618 SCIP_CALL( SCIPallocBufferArray(scip, &matrix->nuplocks, matrix->ncols) );
619 SCIP_CALL( SCIPallocBufferArray(scip, &matrix->ndownlocks, matrix->ncols) );
620
621 BMSclearMemoryArray(matrix->nuplocks, matrix->ncols);
622 BMSclearMemoryArray(matrix->ndownlocks, matrix->ncols);
623
624 /* init bounds */
625 for( v = 0; v < matrix->ncols; v++ )
626 {
627 var = matrix->vars[v];
628 assert(var != NULL);
629
630 matrix->lb[v] = SCIPvarGetLbGlobal(var);
631 matrix->ub[v] = SCIPvarGetUbGlobal(var);
632 }
633
634 /* allocate memory */
635 SCIP_CALL( SCIPallocBufferArray(scip, &matrix->rowmatval, nnonzstmp) );
636 SCIP_CALL( SCIPallocBufferArray(scip, &matrix->rowmatind, nnonzstmp) );
637 SCIP_CALL( SCIPallocBufferArray(scip, &matrix->rowmatbeg, nconss) );
638 SCIP_CALL( SCIPallocBufferArray(scip, &matrix->rowmatcnt, nconss) );
639 SCIP_CALL( SCIPallocBufferArray(scip, &matrix->lhs, nconss) );
640 SCIP_CALL( SCIPallocBufferArray(scip, &matrix->rhs, nconss) );
641 SCIP_CALL( SCIPallocBufferArray(scip, &matrix->cons, nconss) );
642 SCIP_CALL( SCIPallocClearMemoryArray(scip, &matrix->isrhsinfinite, nconss) );
643 SCIP_CALL( SCIPallocBufferArray(scip, &matrix->minactivity, nconss) );
644 SCIP_CALL( SCIPallocBufferArray(scip, &matrix->maxactivity, nconss) );
645 SCIP_CALL( SCIPallocBufferArray(scip, &matrix->minactivityneginf, nconss) );
646 SCIP_CALL( SCIPallocBufferArray(scip, &matrix->minactivityposinf, nconss) );
647 SCIP_CALL( SCIPallocBufferArray(scip, &matrix->maxactivityneginf, nconss) );
648 SCIP_CALL( SCIPallocBufferArray(scip, &matrix->maxactivityposinf, nconss) );
649
650 cnt = 0;
651
652 /* loop a second time over constraints handlers and add supported constraints to the matrix */
653 for( i = 0; i < nconshdlrs; ++i )
654 {
655 SCIP_CONS** conshdlrconss;
656 int nconshdlrconss;
657 SCIP_Bool rowadded;
658
659 if( SCIPisStopped(scip) || (onlyifcomplete && !(*complete)) )
660 {
661 stopped = TRUE;
662 break;
663 }
664
665 conshdlrname = SCIPconshdlrGetName(conshdlrs[i]);
666 conshdlrconss = SCIPconshdlrGetCheckConss(conshdlrs[i]);
667 nconshdlrconss = SCIPconshdlrGetNCheckConss(conshdlrs[i]);
668
669 if( strcmp(conshdlrname, "linear") == 0 )
670 {
671 for( c = 0; c < nconshdlrconss && (c % 1000 != 0 || !SCIPisStopped(scip)); ++c )
672 {
673 cons = conshdlrconss[c];
674 assert(SCIPconsIsTransformed(cons));
675
676 /* do not include constraints that can be altered due to column generation */
677 if( SCIPconsIsModifiable(cons) )
678 {
679 *complete = FALSE;
680
681 if( onlyifcomplete )
682 break;
683
684 continue;
685 }
686
687 SCIP_CALL( addConstraint(scip, matrix, SCIPgetVarsLinear(scip, cons),
688 SCIPgetValsLinear(scip, cons), SCIPgetNVarsLinear(scip, cons),
689 SCIPgetLhsLinear(scip, cons), SCIPgetRhsLinear(scip, cons), nnonzstmp, &rowadded) );
690
691 if(rowadded)
692 {
693 assert(cnt < nconss);
694 matrix->cons[cnt] = cons;
695 cnt++;
696 }
697 }
698 }
699 else if( strcmp(conshdlrname, "setppc") == 0 )
700 {
701 for( c = 0; c < nconshdlrconss && (c % 1000 != 0 || !SCIPisStopped(scip)); ++c )
702 {
703 SCIP_Real lhs;
704 SCIP_Real rhs;
705
706 cons = conshdlrconss[c];
707 assert(SCIPconsIsTransformed(cons));
708
709 /* do not include constraints that can be altered due to column generation */
710 if( SCIPconsIsModifiable(cons) )
711 {
712 *complete = FALSE;
713
714 if( onlyifcomplete )
715 break;
716
717 continue;
718 }
719
720 switch( SCIPgetTypeSetppc(scip, cons) )
721 {
722 case SCIP_SETPPCTYPE_PARTITIONING :
723 lhs = 1.0;
724 rhs = 1.0;
725 break;
726 case SCIP_SETPPCTYPE_PACKING :
727 lhs = -SCIPinfinity(scip);
728 rhs = 1.0;
729 break;
730 case SCIP_SETPPCTYPE_COVERING :
731 lhs = 1.0;
732 rhs = SCIPinfinity(scip);
733 break;
734 default:
735 return SCIP_ERROR;
736 }
737
738 SCIP_CALL( addConstraint(scip, matrix, SCIPgetVarsSetppc(scip, cons), NULL,
739 SCIPgetNVarsSetppc(scip, cons), lhs, rhs, nnonzstmp, &rowadded) );
740
741 if(rowadded)
742 {
743 assert(cnt < nconss);
744 matrix->cons[cnt] = cons;
745 cnt++;
746 }
747 }
748 }
749 else if( strcmp(conshdlrname, "logicor") == 0 )
750 {
751 for( c = 0; c < nconshdlrconss && (c % 1000 != 0 || !SCIPisStopped(scip)); ++c )
752 {
753 cons = conshdlrconss[c];
754 assert(SCIPconsIsTransformed(cons));
755
756 /* do not include constraints that can be altered due to column generation */
757 if( SCIPconsIsModifiable(cons) )
758 {
759 *complete = FALSE;
760
761 if( onlyifcomplete )
762 break;
763
764 continue;
765 }
766
767 SCIP_CALL( addConstraint(scip, matrix, SCIPgetVarsLogicor(scip, cons),
768 NULL, SCIPgetNVarsLogicor(scip, cons), 1.0, SCIPinfinity(scip), nnonzstmp, &rowadded) );
769
770 if(rowadded)
771 {
772 assert(cnt < nconss);
773 matrix->cons[cnt] = cons;
774 cnt++;
775 }
776 }
777 }
778 else if( strcmp(conshdlrname, "knapsack") == 0 )
779 {
780 if( nconshdlrconss > 0 )
781 {
782 SCIP_Real* consvals;
783 int valssize;
784
785 valssize = 100;
786 SCIP_CALL( SCIPallocBufferArray(scip, &consvals, valssize) );
787
788 for( c = 0; c < nconshdlrconss && (c % 1000 != 0 || !SCIPisStopped(scip)); ++c )
789 {
790 SCIP_Longint* weights;
791
792 cons = conshdlrconss[c];
793 assert(SCIPconsIsTransformed(cons));
794
795 /* do not include constraints that can be altered due to column generation */
796 if( SCIPconsIsModifiable(cons) )
797 {
798 *complete = FALSE;
799
800 if( onlyifcomplete )
801 break;
802
803 continue;
804 }
805
806 weights = SCIPgetWeightsKnapsack(scip, cons);
807 nvars = SCIPgetNVarsKnapsack(scip, cons);
808
809 if( nvars > valssize )
810 {
811 valssize = (int) (1.5 * nvars);
812 SCIP_CALL( SCIPreallocBufferArray(scip, &consvals, valssize) );
813 }
814
815 for( v = 0; v < nvars; v++ )
816 consvals[v] = (SCIP_Real)weights[v];
817
818 SCIP_CALL( addConstraint(scip, matrix, SCIPgetVarsKnapsack(scip, cons), consvals,
819 SCIPgetNVarsKnapsack(scip, cons), -SCIPinfinity(scip),
820 (SCIP_Real)SCIPgetCapacityKnapsack(scip, cons), nnonzstmp, &rowadded) );
821
822 if(rowadded)
823 {
824 assert(cnt < nconss);
825 matrix->cons[cnt] = cons;
826 cnt++;
827 }
828 }
829
830 SCIPfreeBufferArray(scip, &consvals);
831 }
832 }
833 else if( strcmp(conshdlrname, "varbound") == 0 )
834 {
835 if( nconshdlrconss > 0 )
836 {
837 SCIP_VAR** consvars;
838 SCIP_Real* consvals;
839
840 SCIP_CALL( SCIPallocBufferArray(scip, &consvars, 2) );
841 SCIP_CALL( SCIPallocBufferArray(scip, &consvals, 2) );
842 consvals[0] = 1.0;
843
844 for( c = 0; c < nconshdlrconss && (c % 1000 != 0 || !SCIPisStopped(scip)); ++c )
845 {
846 cons = conshdlrconss[c];
847 assert(SCIPconsIsTransformed(cons));
848
849 /* do not include constraints that can be altered due to column generation */
850 if( SCIPconsIsModifiable(cons) )
851 {
852 *complete = FALSE;
853
854 if( onlyifcomplete )
855 break;
856
857 continue;
858 }
859
860 consvars[0] = SCIPgetVarVarbound(scip, cons);
861 consvars[1] = SCIPgetVbdvarVarbound(scip, cons);
862
863 consvals[1] = SCIPgetVbdcoefVarbound(scip, cons);
864
865 SCIP_CALL( addConstraint(scip, matrix, consvars, consvals, 2, SCIPgetLhsVarbound(scip, cons),
866 SCIPgetRhsVarbound(scip, cons), nnonzstmp, &rowadded) );
867
868 if(rowadded)
869 {
870 assert(cnt < nconss);
871 matrix->cons[cnt] = cons;
872 cnt++;
873 }
874 }
875
876 SCIPfreeBufferArray(scip, &consvals);
877 SCIPfreeBufferArray(scip, &consvars);
878 }
879 }
880 /* the code below is correct. However, it needs to be disabled
881 * because some of the presolvers can currently only handle 1-1 row-cons relationships,
882 * while the linking constraint handler requires a representation as 2 linear constraints.
883 */
884 #ifdef SCIP_DISABLED_CODE
885 else if( strcmp(conshdlrname, "linking") == 0 )
886 {
887 if( nconshdlrconss > 0 )
888 {
889 SCIP_VAR** consvars;
890 SCIP_VAR** curconsvars;
891 SCIP_Real* consvals;
892 int* curconsvals;
893 int valssize;
894 int nconsvars;
895 int j;
896
897 valssize = 100;
898 SCIP_CALL( SCIPallocBufferArray(scip, &consvars, valssize) );
899 SCIP_CALL( SCIPallocBufferArray(scip, &consvals, valssize) );
900
901 for( c = 0; c < nconshdlrconss && (c % 1000 != 0 || !SCIPisStopped(scip)); ++c )
902 {
903 cons = conshdlrconss[c];
904 assert(SCIPconsIsTransformed(cons));
905
906 /* do not include constraints that can be altered due to column generation */
907 if( SCIPconsIsModifiable(cons) )
908 {
909 *complete = FALSE;
910
911 if( onlyifcomplete )
912 break;
913
914 continue;
915 }
916
917 /* get constraint variables and their amount */
918 SCIP_CALL( SCIPgetBinvarsLinking(scip, cons, &curconsvars, &nconsvars) );
919 curconsvals = SCIPgetValsLinking(scip, cons);
920
921 /* SCIPgetBinVarsLinking returns the number of binary variables, but we also need the integer variable */
922 nconsvars++;
923
924 if( nconsvars > valssize )
925 {
926 valssize = (int) (1.5 * nconsvars);
927 SCIP_CALL( SCIPreallocBufferArray(scip, &consvars, valssize) );
928 SCIP_CALL( SCIPreallocBufferArray(scip, &consvals, valssize) );
929 }
930
931 /* copy vars and vals for binary variables */
932 for( j = 0; j < nconsvars - 1; j++ )
933 {
934 consvars[j] = curconsvars[j];
935 consvals[j] = (SCIP_Real) curconsvals[j];
936 }
937
938 /* set final entry of vars and vals to the linking variable and its coefficient, respectively */
939 consvars[nconsvars - 1] = SCIPgetIntvarLinking(scip, cons);
940 consvals[nconsvars - 1] = -1;
941
942 SCIP_CALL( addConstraint(scip, matrix, consvars, consvals, nconsvars, 0.0, 0.0, nnonzstmp, &rowadded) );
943 SCIP_CALL( addConstraint(scip, matrix, consvars, NULL, nconsvars - 1, 1.0, 1.0, nnonzstmp, &rowadded) );
944
945 if(rowadded)
946 {
947 assert(cnt < nconss);
948 matrix->cons[cnt] = cons;
949 matrix->cons[cnt + 1] = cons;
950 cnt += 2;
951 }
952 }
953
954 SCIPfreeBufferArray(scip, &consvals);
955 SCIPfreeBufferArray(scip, &consvars);
956 }
957 }
958 #endif
959 }
960 assert(matrix->nrows == cnt);
961 assert(matrix->nrows <= nconss);
962 assert(matrix->nnonzs <= nnonzstmp);
963
964 if( *complete )
965 {
966 SCIP_Bool lockmismatch = FALSE;
967
968 for( i = 0; i < matrix->ncols; ++i )
969 {
970 if( SCIPmatrixUplockConflict(matrix, i) || SCIPmatrixDownlockConflict(matrix, i) )
971 {
972 lockmismatch = TRUE;
973 break;
974 }
975 }
976
977 if( lockmismatch )
978 {
979 *complete = FALSE;
980 if( onlyifcomplete )
981 stopped = TRUE;
982 }
983 }
984
985 if( !stopped )
986 {
987 /* calculate row activity bounds */
988 SCIP_CALL( calcActivityBounds(scip, matrix) );
989
990 /* transform row major format into column major format */
991 SCIP_CALL( setColumnMajorFormat(scip, matrix) );
992
993 *initialized = TRUE;
994 }
995 else
996 {
997 SCIPfreeBufferArray(scip, &matrix->maxactivityposinf);
998 SCIPfreeBufferArray(scip, &matrix->maxactivityneginf);
999 SCIPfreeBufferArray(scip, &matrix->minactivityposinf);
1000 SCIPfreeBufferArray(scip, &matrix->minactivityneginf);
1001 SCIPfreeBufferArray(scip, &matrix->maxactivity);
1002 SCIPfreeBufferArray(scip, &matrix->minactivity);
1003
1004 SCIPfreeMemoryArray(scip, &matrix->isrhsinfinite);
1005 SCIPfreeBufferArray(scip, &matrix->cons);
1006
1007 SCIPfreeBufferArray(scip, &matrix->rhs);
1008 SCIPfreeBufferArray(scip, &matrix->lhs);
1009 SCIPfreeBufferArray(scip, &matrix->rowmatcnt);
1010 SCIPfreeBufferArray(scip, &matrix->rowmatbeg);
1011 SCIPfreeBufferArray(scip, &matrix->rowmatind);
1012 SCIPfreeBufferArray(scip, &matrix->rowmatval);
1013
1014 SCIPfreeBufferArray(scip, &matrix->ndownlocks);
1015 SCIPfreeBufferArray(scip, &matrix->nuplocks);
1016 SCIPfreeBufferArray(scip, &matrix->ub);
1017 SCIPfreeBufferArray(scip, &matrix->lb);
1018 SCIPfreeBufferArray(scip, &matrix->colmatcnt);
1019 SCIPfreeBufferArray(scip, &matrix->colmatbeg);
1020 SCIPfreeBufferArray(scip, &matrix->colmatind);
1021 SCIPfreeBufferArray(scip, &matrix->colmatval);
1022 SCIPfreeBufferArrayNull(scip, &matrix->vars);
1023
1024 SCIPfreeBuffer(scip, matrixptr);
1025 }
1026
1027 return SCIP_OKAY;
1028 }
1029
1030
1031 /** frees the constraint matrix */
SCIPmatrixFree(SCIP * scip,SCIP_MATRIX ** matrix)1032 void SCIPmatrixFree(
1033 SCIP* scip, /**< current SCIP instance */
1034 SCIP_MATRIX** matrix /**< constraint matrix object */
1035 )
1036 {
1037 assert(scip != NULL);
1038 assert(matrix != NULL);
1039
1040 if( (*matrix) != NULL )
1041 {
1042 assert((*matrix)->colmatval != NULL);
1043 assert((*matrix)->colmatind != NULL);
1044 assert((*matrix)->colmatbeg != NULL);
1045 assert((*matrix)->colmatcnt != NULL);
1046 assert((*matrix)->lb != NULL);
1047 assert((*matrix)->ub != NULL);
1048 assert((*matrix)->nuplocks != NULL);
1049 assert((*matrix)->ndownlocks != NULL);
1050
1051 assert((*matrix)->rowmatval != NULL);
1052 assert((*matrix)->rowmatind != NULL);
1053 assert((*matrix)->rowmatbeg != NULL);
1054 assert((*matrix)->rowmatcnt != NULL);
1055 assert((*matrix)->lhs != NULL);
1056 assert((*matrix)->rhs != NULL);
1057
1058 SCIPfreeBufferArray(scip, &((*matrix)->maxactivityposinf));
1059 SCIPfreeBufferArray(scip, &((*matrix)->maxactivityneginf));
1060 SCIPfreeBufferArray(scip, &((*matrix)->minactivityposinf));
1061 SCIPfreeBufferArray(scip, &((*matrix)->minactivityneginf));
1062 SCIPfreeBufferArray(scip, &((*matrix)->maxactivity));
1063 SCIPfreeBufferArray(scip, &((*matrix)->minactivity));
1064
1065 SCIPfreeMemoryArray(scip, &((*matrix)->isrhsinfinite));
1066 SCIPfreeBufferArray(scip, &((*matrix)->cons));
1067
1068 SCIPfreeBufferArray(scip, &((*matrix)->rhs));
1069 SCIPfreeBufferArray(scip, &((*matrix)->lhs));
1070 SCIPfreeBufferArray(scip, &((*matrix)->rowmatcnt));
1071 SCIPfreeBufferArray(scip, &((*matrix)->rowmatbeg));
1072 SCIPfreeBufferArray(scip, &((*matrix)->rowmatind));
1073 SCIPfreeBufferArray(scip, &((*matrix)->rowmatval));
1074
1075 SCIPfreeBufferArray(scip, &((*matrix)->ndownlocks));
1076 SCIPfreeBufferArray(scip, &((*matrix)->nuplocks));
1077 SCIPfreeBufferArray(scip, &((*matrix)->ub));
1078 SCIPfreeBufferArray(scip, &((*matrix)->lb));
1079 SCIPfreeBufferArray(scip, &((*matrix)->colmatcnt));
1080 SCIPfreeBufferArray(scip, &((*matrix)->colmatbeg));
1081 SCIPfreeBufferArray(scip, &((*matrix)->colmatind));
1082 SCIPfreeBufferArray(scip, &((*matrix)->colmatval));
1083
1084 (*matrix)->nrows = 0;
1085 (*matrix)->ncols = 0;
1086 (*matrix)->nnonzs = 0;
1087
1088 SCIPfreeBufferArrayNull(scip, &((*matrix)->vars));
1089
1090 SCIPfreeBuffer(scip, matrix);
1091 }
1092 }
1093
1094 /** print one row of the matrix */
SCIPmatrixPrintRow(SCIP * scip,SCIP_MATRIX * matrix,int row)1095 void SCIPmatrixPrintRow(
1096 SCIP* scip, /**< current SCIP instance */
1097 SCIP_MATRIX* matrix, /**< constraint matrix object */
1098 int row /**< row index */
1099 )
1100 {
1101 int* rowpnt;
1102 int* rowend;
1103 int col;
1104 SCIP_Real val;
1105 SCIP_Real* valpnt;
1106
1107 SCIP_UNUSED(scip);
1108
1109 rowpnt = matrix->rowmatind + matrix->rowmatbeg[row];
1110 rowend = rowpnt + matrix->rowmatcnt[row];
1111 valpnt = matrix->rowmatval + matrix->rowmatbeg[row];
1112
1113 printf("### %s: %.15g <=", SCIPconsGetName(matrix->cons[row]), matrix->lhs[row]);
1114 for(; (rowpnt < rowend); rowpnt++, valpnt++)
1115 {
1116 col = *rowpnt;
1117 val = *valpnt;
1118 if( val < 0 )
1119 printf(" %.15g %s [%.15g,%.15g]", val, SCIPvarGetName(matrix->vars[col]),
1120 SCIPvarGetLbGlobal(matrix->vars[col]), SCIPvarGetUbGlobal(matrix->vars[col]));
1121 else
1122 printf(" +%.15g %s [%.15g,%.15g]", val, SCIPvarGetName(matrix->vars[col]),
1123 SCIPvarGetLbGlobal(matrix->vars[col]), SCIPvarGetUbGlobal(matrix->vars[col]));
1124 }
1125 printf(" <= %.15g ###\n", matrix->rhs[row]);
1126 }
1127
1128 /** removes the bounds of a column and updates the activities accordingly */
SCIPmatrixRemoveColumnBounds(SCIP * scip,SCIP_MATRIX * matrix,int col)1129 void SCIPmatrixRemoveColumnBounds(
1130 SCIP* scip, /**< current scip instance */
1131 SCIP_MATRIX* matrix, /**< constraint matrix */
1132 int col /**< column variable to remove bounds from */
1133 )
1134 {
1135 int colmatend = matrix->colmatbeg[col] + matrix->colmatcnt[col];
1136 int i;
1137
1138 for( i = matrix->colmatbeg[col]; i != colmatend; ++i )
1139 {
1140 int row = matrix->colmatind[i];
1141 SCIP_Real val = matrix->colmatval[i];
1142
1143 /* set lower bound to -infinity if necessary */
1144 if( !SCIPisInfinity(scip, -matrix->lb[col]) )
1145 {
1146 if( val > 0.0 )
1147 matrix->minactivityneginf[row]++;
1148 else
1149 matrix->maxactivityneginf[row]++;
1150 }
1151
1152 /* set upper bound to infinity if necessary */
1153 if( !SCIPisInfinity(scip, matrix->ub[col]) )
1154 {
1155 if( val > 0.0 )
1156 matrix->maxactivityposinf[row]++;
1157 else
1158 matrix->minactivityposinf[row]++;
1159 }
1160
1161 assert(matrix->maxactivityneginf[row] + matrix->maxactivityposinf[row] > 0);
1162 assert(matrix->minactivityneginf[row] + matrix->minactivityposinf[row] > 0);
1163
1164 /* mark the activities of the rows to be infinite */
1165 matrix->maxactivity[row] = SCIPinfinity(scip);
1166 matrix->minactivity[row] = -SCIPinfinity(scip);
1167 }
1168
1169 matrix->lb[col] = -SCIPinfinity(scip);
1170 matrix->ub[col] = SCIPinfinity(scip);
1171 }
1172
1173 /** detect parallel rows of matrix. rhs/lhs are ignored. */
SCIPmatrixGetParallelRows(SCIP * scip,SCIP_MATRIX * matrix,SCIP_Real * scale,int * pclass)1174 SCIP_RETCODE SCIPmatrixGetParallelRows(
1175 SCIP* scip, /**< SCIP instance */
1176 SCIP_MATRIX* matrix, /**< matrix containing the constraints */
1177 SCIP_Real* scale, /**< scale factors of rows */
1178 int* pclass /**< parallel row classes */
1179 )
1180 {
1181 SCIP_Real* valpnt;
1182 SCIP_Real* values;
1183 int* classsizes;
1184 int* pcset;
1185 int* colpnt;
1186 int* colend;
1187 int* rowindices;
1188 int* pcs;
1189 SCIP_Real startval;
1190 SCIP_Real aij;
1191 int startpc;
1192 int startk;
1193 int startt;
1194 int pcsetfill;
1195 int rowidx;
1196 int k;
1197 int t;
1198 int m;
1199 int i;
1200 int c;
1201 int newpclass;
1202 int pc;
1203
1204 assert(scip != NULL);
1205 assert(matrix != NULL);
1206 assert(pclass != NULL);
1207
1208 SCIP_CALL( SCIPallocBufferArray(scip, &classsizes, matrix->nrows) );
1209 SCIP_CALL( SCIPallocBufferArray(scip, &pcset, matrix->nrows) );
1210 SCIP_CALL( SCIPallocBufferArray(scip, &values, matrix->nrows) );
1211 SCIP_CALL( SCIPallocBufferArray(scip, &rowindices, matrix->nrows) );
1212 SCIP_CALL( SCIPallocBufferArray(scip, &pcs, matrix->nrows) );
1213
1214 /* init */
1215 BMSclearMemoryArray(scale, matrix->nrows);
1216 BMSclearMemoryArray(pclass, matrix->nrows);
1217 BMSclearMemoryArray(classsizes, matrix->nrows);
1218 classsizes[0] = matrix->nrows;
1219 pcsetfill = 0;
1220 for( t = 1; t < matrix->nrows; ++t )
1221 pcset[pcsetfill++] = t;
1222
1223 /* loop over all columns */
1224 for( c = 0; c < matrix->ncols; ++c )
1225 {
1226 if( matrix->colmatcnt[c] == 0 )
1227 continue;
1228
1229 colpnt = matrix->colmatind + matrix->colmatbeg[c];
1230 colend = colpnt + matrix->colmatcnt[c];
1231 valpnt = matrix->colmatval + matrix->colmatbeg[c];
1232
1233 i = 0;
1234 for( ; (colpnt < colend); colpnt++, valpnt++ )
1235 {
1236 aij = *valpnt;
1237 rowidx = *colpnt;
1238
1239 if( scale[rowidx] == 0.0 )
1240 scale[rowidx] = aij;
1241 assert(scale[rowidx] != 0.0);
1242
1243 rowindices[i] = rowidx;
1244 values[i] = aij / scale[rowidx];
1245 pc = pclass[rowidx];
1246 assert(pc < matrix->nrows);
1247
1248 /* update class sizes and pclass set */
1249 assert(classsizes[pc] > 0);
1250 classsizes[pc]--;
1251 if( classsizes[pc] == 0 )
1252 {
1253 assert(pcsetfill < matrix->nrows);
1254 pcset[pcsetfill++] = pc;
1255 }
1256 pcs[i] = pc;
1257
1258 i++;
1259 }
1260
1261 /* sort on the pclass values */
1262 if( i > 1 )
1263 {
1264 SCIPsortIntIntReal(pcs, rowindices, values, i);
1265 }
1266
1267 k = 0;
1268 while( TRUE ) /*lint !e716*/
1269 {
1270 assert(k < i);
1271 startpc = pcs[k];
1272 startk = k;
1273
1274 /* find pclass-sets */
1275 while( k < i && pcs[k] == startpc )
1276 k++;
1277
1278 /* sort on the A values which have equal pclass values */
1279 if( k - startk > 1 )
1280 SCIPsortRealInt(&(values[startk]), &(rowindices[startk]), k - startk);
1281
1282 t = 0;
1283 while( TRUE ) /*lint !e716*/
1284 {
1285 assert(startk + t < i);
1286 startval = values[startk + t];
1287 startt = t;
1288
1289 /* find A-sets */
1290 while( t < k - startk && SCIPisEQ(scip, startval, values[startk + t]) )
1291 t++;
1292
1293 /* get new pclass */
1294 newpclass = pcset[0];
1295 assert(pcsetfill > 0);
1296 pcset[0] = pcset[--pcsetfill];
1297
1298 /* renumbering */
1299 for( m = startk + startt; m < startk + t; m++ )
1300 {
1301 assert(m < i);
1302 assert(rowindices[m] < matrix->nrows);
1303 assert(newpclass < matrix->nrows);
1304
1305 pclass[rowindices[m]] = newpclass;
1306 classsizes[newpclass]++;
1307 }
1308
1309 if( t == k - startk )
1310 break;
1311 }
1312
1313 if( k == matrix->colmatcnt[c] )
1314 break;
1315 }
1316 }
1317
1318 SCIPfreeBufferArray(scip, &pcs);
1319 SCIPfreeBufferArray(scip, &rowindices);
1320 SCIPfreeBufferArray(scip, &values);
1321 SCIPfreeBufferArray(scip, &pcset);
1322 SCIPfreeBufferArray(scip, &classsizes);
1323
1324 return SCIP_OKAY;
1325 }
1326
1327 /** detect parallel rows of matrix.
1328 * obj coefficients are ignored.
1329 */
SCIPmatrixGetParallelCols(SCIP * scip,SCIP_MATRIX * matrix,SCIP_Real * scale,int * pclass,SCIP_Bool * varineq)1330 SCIP_RETCODE SCIPmatrixGetParallelCols(
1331 SCIP* scip, /**< SCIP instance */
1332 SCIP_MATRIX* matrix, /**< matrix containing the constraints */
1333 SCIP_Real* scale, /**< scale factors of cols */
1334 int* pclass, /**< parallel column classes */
1335 SCIP_Bool* varineq /**< indicating if variable is within an equation */
1336 )
1337 {
1338 SCIP_Real* valpnt;
1339 SCIP_Real* values;
1340 int* classsizes;
1341 int* pcset;
1342 int* rowpnt;
1343 int* rowend;
1344 int* colindices;
1345 int* pcs;
1346 SCIP_Real startval;
1347 SCIP_Real aij;
1348 int startpc;
1349 int startk;
1350 int startt;
1351 int pcsetfill;
1352 int colidx;
1353 int k;
1354 int t;
1355 int m;
1356 int i;
1357 int r;
1358 int newpclass;
1359 int pc;
1360
1361 assert(scip != NULL);
1362 assert(matrix != NULL);
1363 assert(pclass != NULL);
1364 assert(varineq != NULL);
1365
1366 SCIP_CALL( SCIPallocBufferArray(scip, &classsizes, matrix->ncols) );
1367 SCIP_CALL( SCIPallocBufferArray(scip, &pcset, matrix->ncols) );
1368 SCIP_CALL( SCIPallocBufferArray(scip, &values, matrix->ncols) );
1369 SCIP_CALL( SCIPallocBufferArray(scip, &colindices, matrix->ncols) );
1370 SCIP_CALL( SCIPallocBufferArray(scip, &pcs, matrix->ncols) );
1371
1372 /* init */
1373 BMSclearMemoryArray(scale, matrix->ncols);
1374 BMSclearMemoryArray(pclass, matrix->ncols);
1375 BMSclearMemoryArray(classsizes, matrix->ncols);
1376 classsizes[0] = matrix->ncols;
1377 pcsetfill = 0;
1378 for( t = 1; t < matrix->ncols; ++t )
1379 pcset[pcsetfill++] = t;
1380
1381 /* loop over all rows */
1382 for( r = 0; r < matrix->nrows; ++r )
1383 {
1384 /* we consider only equations or ranged rows */
1385 if( !matrix->isrhsinfinite[r] )
1386 {
1387 rowpnt = matrix->rowmatind + matrix->rowmatbeg[r];
1388 rowend = rowpnt + matrix->rowmatcnt[r];
1389 valpnt = matrix->rowmatval + matrix->rowmatbeg[r];
1390
1391 i = 0;
1392 for( ; (rowpnt < rowend); rowpnt++, valpnt++ )
1393 {
1394 aij = *valpnt;
1395 colidx = *rowpnt;
1396
1397 /* remember variable was part of an equation or ranged row */
1398 varineq[colidx] = TRUE;
1399
1400 if( scale[colidx] == 0.0 )
1401 scale[colidx] = aij;
1402 assert(scale[colidx] != 0.0);
1403
1404 colindices[i] = colidx;
1405 values[i] = aij / scale[colidx];
1406 pc = pclass[colidx];
1407 assert(pc < matrix->ncols);
1408
1409 /* update class sizes and pclass set */
1410 assert(classsizes[pc] > 0);
1411 classsizes[pc]--;
1412 if( classsizes[pc] == 0 )
1413 {
1414 assert(pcsetfill < matrix->ncols);
1415 pcset[pcsetfill++] = pc;
1416 }
1417 pcs[i] = pc;
1418
1419 i++;
1420 }
1421
1422 /* sort on the pclass values */
1423 if( i > 1 )
1424 {
1425 SCIPsortIntIntReal(pcs, colindices, values, i);
1426 }
1427
1428 k = 0;
1429 while( TRUE ) /*lint !e716*/
1430 {
1431 assert(k < i);
1432 startpc = pcs[k];
1433 startk = k;
1434
1435 /* find pclass-sets */
1436 while( k < i && pcs[k] == startpc )
1437 k++;
1438
1439 /* sort on the A values which have equal pclass values */
1440 if( k - startk > 1 )
1441 SCIPsortRealInt(&(values[startk]), &(colindices[startk]), k - startk);
1442
1443 t = 0;
1444 while( TRUE ) /*lint !e716*/
1445 {
1446 assert(startk + t < i);
1447 startval = values[startk + t];
1448 startt = t;
1449
1450 /* find A-sets */
1451 while( t < k - startk && SCIPisEQ(scip, startval, values[startk + t]) )
1452 t++;
1453
1454 /* get new pclass */
1455 newpclass = pcset[0];
1456 assert(pcsetfill > 0);
1457 pcset[0] = pcset[--pcsetfill];
1458
1459 /* renumbering */
1460 for( m = startk + startt; m < startk + t; m++ )
1461 {
1462 assert(m < i);
1463 assert(colindices[m] < matrix->ncols);
1464 assert(newpclass < matrix->ncols);
1465
1466 pclass[colindices[m]] = newpclass;
1467 classsizes[newpclass]++;
1468 }
1469
1470 if( t == k - startk )
1471 break;
1472 }
1473
1474 if( k == matrix->rowmatcnt[r] )
1475 break;
1476 }
1477 }
1478 }
1479
1480 SCIPfreeBufferArray(scip, &pcs);
1481 SCIPfreeBufferArray(scip, &colindices);
1482 SCIPfreeBufferArray(scip, &values);
1483 SCIPfreeBufferArray(scip, &pcset);
1484 SCIPfreeBufferArray(scip, &classsizes);
1485
1486 return SCIP_OKAY;
1487 }
1488
1489
1490 /*
1491 * access functions implemented as defines
1492 */
1493
1494 /* In debug mode, the following methods are implemented as function calls to ensure
1495 * type validity.
1496 * In optimized mode, the methods are implemented as defines to improve performance.
1497 * However, we want to have them in the library anyways, so we have to undef the defines.
1498 */
1499
1500 #undef SCIPmatrixGetColValPtr
1501 #undef SCIPmatrixGetColIdxPtr
1502 #undef SCIPmatrixGetColNNonzs
1503 #undef SCIPmatrixGetNColumns
1504 #undef SCIPmatrixGetColUb
1505 #undef SCIPmatrixGetColLb
1506 #undef SCIPmatrixGetColNUplocks
1507 #undef SCIPmatrixGetColNDownlocks
1508 #undef SCIPmatrixGetVar
1509 #undef SCIPmatrixGetColName
1510 #undef SCIPmatrixGetRowValPtr
1511 #undef SCIPmatrixGetRowIdxPtr
1512 #undef SCIPmatrixGetRowNNonzs
1513 #undef SCIPmatrixGetRowName
1514 #undef SCIPmatrixGetNRows
1515 #undef SCIPmatrixGetRowLhs
1516 #undef SCIPmatrixGetRowRhs
1517 #undef SCIPmatrixIsRowRhsInfinity
1518 #undef SCIPmatrixGetNNonzs
1519 #undef SCIPmatrixGetRowMinActivity
1520 #undef SCIPmatrixGetRowMaxActivity
1521 #undef SCIPmatrixGetRowNMinActNegInf
1522 #undef SCIPmatrixGetRowNMinActPosInf
1523 #undef SCIPmatrixGetRowNMaxActNegInf
1524 #undef SCIPmatrixGetRowNMaxActPosInf
1525 #undef SCIPmatrixGetCons
1526
1527 /** get column based start pointer of values */
SCIPmatrixGetColValPtr(SCIP_MATRIX * matrix,int col)1528 SCIP_Real* SCIPmatrixGetColValPtr(
1529 SCIP_MATRIX* matrix, /**< matrix instance */
1530 int col /**< column index */
1531 )
1532 {
1533 assert(matrix != NULL);
1534 assert(0 <= col && col < matrix->ncols);
1535
1536 return matrix->colmatval + matrix->colmatbeg[col];
1537 }
1538
1539 /** get column based start pointer of row indices */
SCIPmatrixGetColIdxPtr(SCIP_MATRIX * matrix,int col)1540 int* SCIPmatrixGetColIdxPtr(
1541 SCIP_MATRIX* matrix, /**< matrix instance */
1542 int col /**< column index */
1543 )
1544 {
1545 assert(matrix != NULL);
1546 assert(0 <= col && col < matrix->ncols);
1547
1548 return matrix->colmatind + matrix->colmatbeg[col];
1549 }
1550
1551 /** get the number of non-zero entries of this column */
SCIPmatrixGetColNNonzs(SCIP_MATRIX * matrix,int col)1552 int SCIPmatrixGetColNNonzs(
1553 SCIP_MATRIX* matrix, /**< matrix instance */
1554 int col /**< column index */
1555 )
1556 {
1557 assert(matrix != NULL);
1558 assert(0 <= col && col < matrix->ncols);
1559
1560 return matrix->colmatcnt[col];
1561 }
1562
1563 /** get number of columns of the matrix */
SCIPmatrixGetNColumns(SCIP_MATRIX * matrix)1564 int SCIPmatrixGetNColumns(
1565 SCIP_MATRIX* matrix /**< matrix instance */
1566 )
1567 {
1568 assert(matrix != NULL);
1569
1570 return matrix->ncols;
1571 }
1572
1573 /** get upper bound of column */
SCIPmatrixGetColUb(SCIP_MATRIX * matrix,int col)1574 SCIP_Real SCIPmatrixGetColUb(
1575 SCIP_MATRIX* matrix, /**< matrix instance */
1576 int col /**< column index */
1577 )
1578 {
1579 assert(matrix != NULL);
1580
1581 return matrix->ub[col];
1582 }
1583
1584 /** get lower bound of column */
SCIPmatrixGetColLb(SCIP_MATRIX * matrix,int col)1585 SCIP_Real SCIPmatrixGetColLb(
1586 SCIP_MATRIX* matrix, /**< matrix instance */
1587 int col /**< column index */
1588 )
1589 {
1590 assert(matrix != NULL);
1591
1592 return matrix->lb[col];
1593 }
1594
1595 /** get number of uplocks of column */
SCIPmatrixGetColNUplocks(SCIP_MATRIX * matrix,int col)1596 int SCIPmatrixGetColNUplocks(
1597 SCIP_MATRIX* matrix, /**< matrix instance */
1598 int col /**< column index */
1599 )
1600 {
1601 assert(matrix != NULL);
1602 assert(0 <= col && col < matrix->ncols);
1603
1604 return matrix->nuplocks[col];
1605 }
1606
1607 /** get number of downlocks of column */
SCIPmatrixGetColNDownlocks(SCIP_MATRIX * matrix,int col)1608 int SCIPmatrixGetColNDownlocks(
1609 SCIP_MATRIX* matrix, /**< matrix instance */
1610 int col /**< column index */
1611 )
1612 {
1613 assert(matrix != NULL);
1614 assert(0 <= col && col < matrix->ncols);
1615
1616 return matrix->ndownlocks[col];
1617 }
1618
1619 /** get variable pointer of column */
SCIPmatrixGetVar(SCIP_MATRIX * matrix,int col)1620 SCIP_VAR* SCIPmatrixGetVar(
1621 SCIP_MATRIX* matrix, /**< matrix instance */
1622 int col /**< column index */
1623 )
1624 {
1625 assert(matrix != NULL);
1626 assert(0 <= col && col < matrix->ncols);
1627
1628 return matrix->vars[col];
1629 }
1630
1631 /** get name of column/variable */
SCIPmatrixGetColName(SCIP_MATRIX * matrix,int col)1632 const char* SCIPmatrixGetColName(
1633 SCIP_MATRIX* matrix, /**< matrix instance */
1634 int col /**< column index */
1635 )
1636 {
1637 assert(matrix != NULL);
1638 assert(0 <= col && col < matrix->ncols);
1639
1640 return SCIPvarGetName(matrix->vars[col]);
1641 }
1642
1643 /** get row based start pointer of values */
SCIPmatrixGetRowValPtr(SCIP_MATRIX * matrix,int row)1644 SCIP_Real* SCIPmatrixGetRowValPtr(
1645 SCIP_MATRIX* matrix, /**< matrix instance */
1646 int row /**< row index */
1647 )
1648 {
1649 assert(matrix != NULL);
1650 assert(0 <= row && row < matrix->nrows);
1651
1652 return matrix->rowmatval + matrix->rowmatbeg[row];
1653 }
1654
1655 /** get row based start pointer of column indices */
SCIPmatrixGetRowIdxPtr(SCIP_MATRIX * matrix,int row)1656 int* SCIPmatrixGetRowIdxPtr(
1657 SCIP_MATRIX* matrix, /**< matrix instance */
1658 int row /**< row index */
1659 )
1660 {
1661 assert(matrix != NULL);
1662 assert(0 <= row && row < matrix->nrows);
1663
1664 return matrix->rowmatind + matrix->rowmatbeg[row];
1665 }
1666
1667 /** get number of non-zeros of this row */
SCIPmatrixGetRowNNonzs(SCIP_MATRIX * matrix,int row)1668 int SCIPmatrixGetRowNNonzs(
1669 SCIP_MATRIX* matrix, /**< matrix instance */
1670 int row /**< row index */
1671 )
1672 {
1673 assert(matrix != NULL);
1674 assert(0 <= row && row < matrix->nrows);
1675
1676 return matrix->rowmatcnt[row];
1677 }
1678
1679 /** get name of row */
SCIPmatrixGetRowName(SCIP_MATRIX * matrix,int row)1680 const char* SCIPmatrixGetRowName(
1681 SCIP_MATRIX* matrix, /**< matrix instance */
1682 int row /**< row index */
1683 )
1684 {
1685 assert(matrix != NULL);
1686 assert(0 <= row && row < matrix->nrows);
1687
1688 return SCIPconsGetName(matrix->cons[row]);
1689 }
1690
1691 /** get number of rows of the matrix */
SCIPmatrixGetNRows(SCIP_MATRIX * matrix)1692 int SCIPmatrixGetNRows(
1693 SCIP_MATRIX* matrix /**< matrix instance */
1694 )
1695 {
1696 assert(matrix != NULL);
1697
1698 return matrix->nrows;
1699 }
1700
1701 /** get left-hand-side of row */
SCIPmatrixGetRowLhs(SCIP_MATRIX * matrix,int row)1702 SCIP_Real SCIPmatrixGetRowLhs(
1703 SCIP_MATRIX* matrix, /**< matrix instance */
1704 int row /**< row index */
1705 )
1706 {
1707 assert(matrix != NULL);
1708 assert(0 <= row && row < matrix->nrows);
1709
1710 return matrix->lhs[row];
1711 }
1712
1713 /** get right-hand-side of row */
SCIPmatrixGetRowRhs(SCIP_MATRIX * matrix,int row)1714 SCIP_Real SCIPmatrixGetRowRhs(
1715 SCIP_MATRIX* matrix, /**< matrix instance */
1716 int row /**< row index */
1717 )
1718 {
1719 assert(matrix != NULL);
1720 assert(0 <= row && row < matrix->nrows);
1721
1722 return matrix->rhs[row];
1723 }
1724
1725 /** flag indicating if right-hand-side of row is infinity */
SCIPmatrixIsRowRhsInfinity(SCIP_MATRIX * matrix,int row)1726 SCIP_Bool SCIPmatrixIsRowRhsInfinity(
1727 SCIP_MATRIX* matrix, /**< matrix instance */
1728 int row /**< row index */
1729 )
1730 {
1731 assert(matrix != NULL);
1732 assert(0 <= row && row < matrix->nrows);
1733
1734 return matrix->isrhsinfinite[row];
1735 }
1736
1737 /** get number of non-zeros of matrix */
SCIPmatrixGetNNonzs(SCIP_MATRIX * matrix)1738 int SCIPmatrixGetNNonzs(
1739 SCIP_MATRIX* matrix /**< matrix instance */
1740 )
1741 {
1742 assert(matrix != NULL);
1743
1744 return matrix->nnonzs;
1745 }
1746
1747 /** get minimal activity of row */
SCIPmatrixGetRowMinActivity(SCIP_MATRIX * matrix,int row)1748 SCIP_Real SCIPmatrixGetRowMinActivity(
1749 SCIP_MATRIX* matrix, /**< matrix instance */
1750 int row /**< row index */
1751 )
1752 {
1753 assert(matrix != NULL);
1754 assert(0 <= row && row < matrix->nrows);
1755
1756 return matrix->minactivity[row];
1757 }
1758
1759 /** get maximal activity of row */
SCIPmatrixGetRowMaxActivity(SCIP_MATRIX * matrix,int row)1760 SCIP_Real SCIPmatrixGetRowMaxActivity(
1761 SCIP_MATRIX* matrix, /**< matrix instance */
1762 int row /**< row index */
1763 )
1764 {
1765 assert(matrix != NULL);
1766 assert(0 <= row && row < matrix->nrows);
1767
1768 return matrix->maxactivity[row];
1769 }
1770
1771 /** get number of negative infinities present within minimal activity */
SCIPmatrixGetRowNMinActNegInf(SCIP_MATRIX * matrix,int row)1772 int SCIPmatrixGetRowNMinActNegInf(
1773 SCIP_MATRIX* matrix, /**< matrix instance */
1774 int row /**< row index */
1775 )
1776 {
1777 assert(matrix != NULL);
1778 assert(0 <= row && row < matrix->nrows);
1779
1780 return matrix->minactivityneginf[row];
1781 }
1782
1783 /** get number of positive infinities present within minimal activity */
SCIPmatrixGetRowNMinActPosInf(SCIP_MATRIX * matrix,int row)1784 int SCIPmatrixGetRowNMinActPosInf(
1785 SCIP_MATRIX* matrix, /**< matrix instance */
1786 int row /**< row index */
1787 )
1788 {
1789 assert(matrix != NULL);
1790 assert(0 <= row && row < matrix->nrows);
1791
1792 return matrix->minactivityposinf[row];
1793 }
1794
1795 /** get number of negative infinities present within maximal activity */
SCIPmatrixGetRowNMaxActNegInf(SCIP_MATRIX * matrix,int row)1796 int SCIPmatrixGetRowNMaxActNegInf(
1797 SCIP_MATRIX* matrix, /**< matrix instance */
1798 int row /**< row index */
1799 )
1800 {
1801 assert(matrix != NULL);
1802 assert(0 <= row && row < matrix->nrows);
1803
1804 return matrix->maxactivityneginf[row];
1805 }
1806
1807 /** get number of positive infinities present within maximal activity */
SCIPmatrixGetRowNMaxActPosInf(SCIP_MATRIX * matrix,int row)1808 int SCIPmatrixGetRowNMaxActPosInf(
1809 SCIP_MATRIX* matrix, /**< matrix instance */
1810 int row /**< row index */
1811 )
1812 {
1813 assert(matrix != NULL);
1814 assert(0 <= row && row < matrix->nrows);
1815
1816 return matrix->maxactivityposinf[row];
1817 }
1818
1819 /** get constraint pointer for constraint representing row */
SCIPmatrixGetCons(SCIP_MATRIX * matrix,int row)1820 SCIP_CONS* SCIPmatrixGetCons(
1821 SCIP_MATRIX* matrix, /**< matrix instance */
1822 int row /**< row index */
1823 )
1824 {
1825 assert(matrix != NULL);
1826 assert(0 <= row && row < matrix->nrows);
1827
1828 return matrix->cons[row];
1829 }
1830
1831 /** get if conflicting uplocks of a specific variable present */
SCIPmatrixUplockConflict(SCIP_MATRIX * matrix,int col)1832 SCIP_Bool SCIPmatrixUplockConflict(
1833 SCIP_MATRIX* matrix, /**< matrix instance */
1834 int col /**< column index */
1835 )
1836 {
1837 assert(matrix != NULL);
1838 assert(0 <= col && col < matrix->ncols);
1839
1840 return (SCIPvarGetNLocksUpType(matrix->vars[col], SCIP_LOCKTYPE_MODEL) != matrix->nuplocks[col]);
1841 }
1842
1843 /** get if conflicting downlocks of a specific variable present */
SCIPmatrixDownlockConflict(SCIP_MATRIX * matrix,int col)1844 SCIP_Bool SCIPmatrixDownlockConflict(
1845 SCIP_MATRIX* matrix, /**< matrix instance */
1846 int col /**< column index */
1847 )
1848 {
1849 assert(matrix != NULL);
1850 assert(0 <= col && col < matrix->ncols);
1851
1852 return (SCIPvarGetNLocksDownType(matrix->vars[col], SCIP_LOCKTYPE_MODEL) != matrix->ndownlocks[col]);
1853 }
1854