Scippy

SCIP

Solving Constraint Integer Programs

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