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 visit scipopt.org. */
13 /* */
14 /* * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * */
15
16 /**@file lpi_glop.cpp
17 * @ingroup LPIS
18 * @brief LP interface for Glop
19 * @author Frederic Didier
20 * @author Marc Pfetsch
21 */
22
23 /*--+----1----+----2----+----3----+----4----+----5----+----6----+----7----+----8----+----9----+----0----+----1----+----2*/
24
25 /* turn off some warnings from includes */
26 #pragma GCC diagnostic ignored "-Wsign-compare"
27 #pragma GCC diagnostic ignored "-Wpedantic"
28 #pragma GCC diagnostic ignored "-Wignored-qualifiers"
29 #pragma GCC diagnostic ignored "-Wshadow"
30 #pragma GCC diagnostic ignored "-Wnon-virtual-dtor"
31 #pragma GCC diagnostic ignored "-Wctor-dtor-privacy"
32 #pragma GCC diagnostic ignored "-Woverflow"
33
34 #include "ortools/base/version.h"
35 #include "ortools/glop/lp_solver.h"
36 #include "ortools/glop/revised_simplex.h"
37 #include "ortools/lp_data/lp_print_utils.h"
38 #include "ortools/lp_data/lp_data_utils.h"
39 #include "ortools/lp_data/proto_utils.h"
40 #include "ortools/util/file_util.h"
41 #include "ortools/util/stats.h"
42 #include "ortools/util/time_limit.h"
43
44 #include "ortools/base/logging.h"
45 #include "ortools/base/vlog_is_on.h"
46
47 #include "lpi/lpi.h"
48 #include "scip/pub_message.h"
49
50 #include <assert.h>
51
52 /* turn warnings on again */
53 #pragma GCC diagnostic warning "-Wsign-compare"
54 #pragma GCC diagnostic warning "-Wpedantic"
55 #pragma GCC diagnostic warning "-Wignored-qualifiers"
56 #pragma GCC diagnostic warning "-Wshadow"
57 #pragma GCC diagnostic warning "-Wnon-virtual-dtor"
58 #pragma GCC diagnostic warning "-Wctor-dtor-privacy"
59 #pragma GCC diagnostic warning "-Woverflow"
60
61
62 using operations_research::TimeLimit;
63 using operations_research::glop::BasisState;
64 using operations_research::glop::ColIndex;
65 using operations_research::glop::ColIndexVector;
66 using operations_research::glop::ConstraintStatus;
67 using operations_research::glop::ConstraintStatusColumn;
68 using operations_research::glop::DenseBooleanColumn;
69 using operations_research::glop::DenseBooleanRow;
70 using operations_research::glop::DenseColumn;
71 using operations_research::glop::DenseRow;
72 using operations_research::glop::SparseColumn;
73 using operations_research::glop::ScatteredColumn;
74 using operations_research::glop::ScatteredColumnIterator;
75 using operations_research::glop::SparseMatrix;
76 using operations_research::glop::Fractional;
77 using operations_research::glop::GetProblemStatusString;
78 using operations_research::glop::ProblemStatus;
79 using operations_research::glop::RowIndex;
80 using operations_research::glop::ScatteredRow;
81 using operations_research::glop::ScatteredRowIterator;
82 using operations_research::glop::VariableStatus;
83 using operations_research::glop::VariableStatusRow;
84 using operations_research::MPModelProto;
85
86 /** LP interface */
87 struct SCIP_LPi
88 {
89 operations_research::glop::LinearProgram* linear_program; /**< the linear program */
90 operations_research::glop::LinearProgram* scaled_lp; /**< scaled linear program */
91 operations_research::glop::RevisedSimplex* solver; /**< direct reference to the revised simplex, not passing through lp_solver */
92 operations_research::glop::GlopParameters* parameters; /**< parameters */
93 operations_research::glop::LpScalingHelper* scaler; /**< scaler auxiliary class */
94
95 /* the following is used by SCIPlpiWasSolved() */
96 bool lp_modified_since_last_solve;
97 bool lp_time_limit_was_reached;
98
99 /* store the values of some parameters in order to be able to return them */
100 bool lp_info; /**< whether additional output is turned on */
101 SCIP_PRICING pricing; /**< SCIP pricing setting */
102 bool from_scratch; /**< store whether basis is ignored for next solving call */
103 int numthreads; /**< number of threads used to solve the LP (0 = automatic) */
104 SCIP_Real conditionlimit; /**< maximum condition number of LP basis counted as stable (-1.0: no limit) */
105 bool checkcondition; /**< Should condition number of LP basis be checked for stability? */
106 int timing; /**< type of timer (1 - cpu, 2 - wallclock, 0 - off) */
107
108 /* other data */
109 SCIP_Longint niterations; /**< number of iterations used */
110
111 /* Temporary vectors allocated here for speed. This gain is non-negligible
112 * because in many situations, only a few entries of these vectors are
113 * inspected (hypersparsity) and allocating them is in O(num_rows) or
114 * O(num_cols) instead of O(num_non_zeros) to read/clear them. */
115 ScatteredRow* tmp_row; /**< temporary vector */
116 ScatteredColumn* tmp_column; /**< temporary vector */
117 };
118
119 /** uncomment to turn off scaling */
120 /* #define NOSCALING */
121
122 /** define feasibility check to possibly reoptimize: 0: no check, 1: completely new check, 2: check unscaled variable and activity values */
123 #ifdef NOSCALING
124 #define UNSCALEDFEAS_CHECK 0
125 #else
126 #define UNSCALEDFEAS_CHECK 2
127 #endif
128
129
130 /*
131 * LP Interface Methods
132 */
133
134 static char glopname[128];
135
136 /** gets name and version of LP solver */
SCIPlpiGetSolverName(void)137 const char* SCIPlpiGetSolverName(
138 void
139 )
140 {
141 (void) snprintf(glopname, 100, "Glop %d.%d", operations_research::OrToolsMajorVersion(), operations_research::OrToolsMinorVersion());
142 return glopname;
143 }
144
145 /** gets description of LP solver (developer, webpage, ...) */
SCIPlpiGetSolverDesc(void)146 const char* SCIPlpiGetSolverDesc(
147 void
148 )
149 {
150 return "Glop Linear Solver, developed by Google, part of OR-Tools (developers.google.com/optimization)";
151 }
152
153 /** gets pointer for LP solver - use only with great care */
SCIPlpiGetSolverPointer(SCIP_LPI * lpi)154 void* SCIPlpiGetSolverPointer(
155 SCIP_LPI* lpi /**< pointer to an LP interface structure */
156 )
157 {
158 assert( lpi != NULL );
159
160 SCIPerrorMessage("SCIPlpiGetSolverPointer() has not been implemented yet.\n");
161
162 return NULL;
163 }
164
165 /** pass integrality information to LP solver */
SCIPlpiSetIntegralityInformation(SCIP_LPI * lpi,int ncols,int * intInfo)166 SCIP_RETCODE SCIPlpiSetIntegralityInformation(
167 SCIP_LPI* lpi, /**< pointer to an LP interface structure */
168 int ncols, /**< length of integrality array */
169 int* intInfo /**< integrality array (0: continuous, 1: integer). May be NULL iff ncols is 0. */
170 )
171 {
172 assert( lpi != NULL );
173 assert( lpi->linear_program != NULL );
174 assert( ncols == 0 || ncols == lpi->linear_program->num_variables().value() );
175
176 /* Pass on integrality information (currently not used by Glop). */
177 for (ColIndex col(0); col < ColIndex(ncols); ++col)
178 {
179 assert( intInfo != NULL );
180 int info = intInfo[col.value()];
181 assert( info == 0 || info == 1 );
182 if ( info == 0 )
183 lpi->linear_program->SetVariableType(col, operations_research::glop::LinearProgram::VariableType::CONTINUOUS);
184 else
185 lpi->linear_program->SetVariableType(col, operations_research::glop::LinearProgram::VariableType::INTEGER);
186 }
187
188 return SCIP_OKAY;
189 }
190
191 /** informs about availability of a primal simplex solving method */
SCIPlpiHasPrimalSolve(void)192 SCIP_Bool SCIPlpiHasPrimalSolve(
193 void
194 )
195 {
196 return TRUE;
197 }
198
199 /** informs about availability of a dual simplex solving method */
SCIPlpiHasDualSolve(void)200 SCIP_Bool SCIPlpiHasDualSolve(
201 void
202 )
203 {
204 return TRUE;
205 }
206
207 /** informs about availability of a barrier solving method */
SCIPlpiHasBarrierSolve(void)208 SCIP_Bool SCIPlpiHasBarrierSolve(
209 void
210 )
211 {
212 return FALSE;
213 }
214
215
216 /*
217 * LPI Creation and Destruction Methods
218 */
219
220 /**@name LPI Creation and Destruction Methods */
221 /**@{ */
222
223 /** creates an LP problem object */
SCIPlpiCreate(SCIP_LPI ** lpi,SCIP_MESSAGEHDLR * messagehdlr,const char * name,SCIP_OBJSEN objsen)224 SCIP_RETCODE SCIPlpiCreate(
225 SCIP_LPI** lpi, /**< pointer to an LP interface structure */
226 SCIP_MESSAGEHDLR* messagehdlr, /**< message handler to use for printing messages, or NULL */
227 const char* name, /**< problem name */
228 SCIP_OBJSEN objsen /**< objective sense */
229 )
230 {
231 assert( lpi != NULL );
232 assert( name != NULL );
233
234 /* Initilialize memory. */
235 SCIP_ALLOC(BMSallocMemory(lpi));
236 (*lpi)->linear_program = new operations_research::glop::LinearProgram();
237 (*lpi)->scaled_lp = new operations_research::glop::LinearProgram();
238 (*lpi)->solver = new operations_research::glop::RevisedSimplex();
239 (*lpi)->parameters = new operations_research::glop::GlopParameters();
240 (*lpi)->scaler = new operations_research::glop::LpScalingHelper();
241
242 /* Set problem name and objective direction. */
243 (*lpi)->linear_program->SetName(std::string(name));
244 SCIP_CALL( SCIPlpiChgObjsen(*lpi, objsen) );
245
246 (*lpi)->from_scratch = false;
247 (*lpi)->lp_info = false;
248 (*lpi)->pricing = SCIP_PRICING_LPIDEFAULT;
249 (*lpi)->lp_modified_since_last_solve = true;
250 (*lpi)->lp_time_limit_was_reached = false;
251 (*lpi)->conditionlimit = -1.0;
252 (*lpi)->checkcondition = false;
253 (*lpi)->niterations = 0LL;
254
255 (*lpi)->tmp_row = new ScatteredRow();
256 (*lpi)->tmp_column = new ScatteredColumn();
257
258 #ifdef NOSCALING
259 (*lpi)->parameters->set_use_scaling(false);
260 #endif
261
262 return SCIP_OKAY;
263 }
264
265 /** deletes an LP problem object */
SCIPlpiFree(SCIP_LPI ** lpi)266 SCIP_RETCODE SCIPlpiFree(
267 SCIP_LPI** lpi /**< pointer to an LP interface structure */
268 )
269 {
270 SCIPdebugMessage("SCIPlpiFree\n");
271
272 delete (*lpi)->scaler;
273 delete (*lpi)->parameters;
274 delete (*lpi)->solver;
275 delete (*lpi)->scaled_lp;
276 delete (*lpi)->linear_program;
277
278 delete (*lpi)->tmp_row;
279 delete (*lpi)->tmp_column;
280
281 BMSfreeMemory(lpi);
282
283 return SCIP_OKAY;
284 }
285
286 /**@} */
287
288
289
290
291 /*
292 * Modification Methods
293 */
294
295 /**@name Modification Methods */
296 /**@{ */
297
298 /** copies LP data with column matrix into LP solver */
SCIPlpiLoadColLP(SCIP_LPI * lpi,SCIP_OBJSEN objsen,int ncols,const SCIP_Real * obj,const SCIP_Real * lb,const SCIP_Real * ub,char ** colnames,int nrows,const SCIP_Real * lhs,const SCIP_Real * rhs,char ** rownames,int nnonz,const int * beg,const int * ind,const SCIP_Real * val)299 SCIP_RETCODE SCIPlpiLoadColLP(
300 SCIP_LPI* lpi, /**< LP interface structure */
301 SCIP_OBJSEN objsen, /**< objective sense */
302 int ncols, /**< number of columns */
303 const SCIP_Real* obj, /**< objective function values of columns */
304 const SCIP_Real* lb, /**< lower bounds of columns */
305 const SCIP_Real* ub, /**< upper bounds of columns */
306 char** colnames, /**< column names, or NULL */
307 int nrows, /**< number of rows */
308 const SCIP_Real* lhs, /**< left hand sides of rows */
309 const SCIP_Real* rhs, /**< right hand sides of rows */
310 char** rownames, /**< row names, or NULL */
311 int nnonz, /**< number of nonzero elements in the constraint matrix */
312 const int* beg, /**< start index of each column in ind- and val-array */
313 const int* ind, /**< row indices of constraint matrix entries */
314 const SCIP_Real* val /**< values of constraint matrix entries */
315 )
316 {
317 assert( lpi != NULL );
318 assert( lpi->linear_program != NULL );
319 assert( obj != NULL );
320 assert( lb != NULL );
321 assert( ub != NULL );
322 assert( beg != NULL );
323 assert( ind != NULL );
324 assert( val != NULL );
325
326 lpi->linear_program->Clear();
327 SCIP_CALL( SCIPlpiAddRows(lpi, nrows, lhs, rhs, rownames, 0, NULL, NULL, NULL) );
328 SCIP_CALL( SCIPlpiAddCols(lpi, ncols, obj, lb, ub, colnames, nnonz, beg, ind, val) );
329 SCIP_CALL( SCIPlpiChgObjsen(lpi, objsen) );
330
331 return SCIP_OKAY;
332 }
333
334 /** adds columns to the LP */
SCIPlpiAddCols(SCIP_LPI * lpi,int ncols,const SCIP_Real * obj,const SCIP_Real * lb,const SCIP_Real * ub,char ** colnames,int nnonz,const int * beg,const int * ind,const SCIP_Real * val)335 SCIP_RETCODE SCIPlpiAddCols(
336 SCIP_LPI* lpi, /**< LP interface structure */
337 int ncols, /**< number of columns to be added */
338 const SCIP_Real* obj, /**< objective function values of new columns */
339 const SCIP_Real* lb, /**< lower bounds of new columns */
340 const SCIP_Real* ub, /**< upper bounds of new columns */
341 char** colnames, /**< column names, or NULL */
342 int nnonz, /**< number of nonzero elements to be added to the constraint matrix */
343 const int* beg, /**< start index of each column in ind- and val-array, or NULL if nnonz == 0 */
344 const int* ind, /**< row indices of constraint matrix entries, or NULL if nnonz == 0 */
345 const SCIP_Real* val /**< values of constraint matrix entries, or NULL if nnonz == 0 */
346 )
347 {
348 assert( lpi != NULL );
349 assert( lpi->linear_program != NULL );
350 assert( obj != NULL );
351 assert( lb != NULL );
352 assert( ub != NULL );
353 assert( nnonz >= 0) ;
354 assert( ncols >= 0) ;
355
356 SCIPdebugMessage("adding %d columns with %d nonzeros.\n", ncols, nnonz);
357
358 /* @todo add names */
359 if ( nnonz > 0 )
360 {
361 assert( beg != NULL );
362 assert( ind != NULL );
363 assert( val != NULL );
364 assert( ncols > 0 );
365
366 #ifndef NDEBUG
367 /* perform check that no new rows are added */
368 RowIndex num_rows = lpi->linear_program->num_constraints();
369 for (int j = 0; j < nnonz; ++j)
370 {
371 assert( 0 <= ind[j] && ind[j] < num_rows.value() );
372 assert( val[j] != 0.0 );
373 }
374 #endif
375
376 int nz = 0;
377 for (int i = 0; i < ncols; ++i)
378 {
379 const ColIndex col = lpi->linear_program->CreateNewVariable();
380 lpi->linear_program->SetVariableBounds(col, lb[i], ub[i]);
381 lpi->linear_program->SetObjectiveCoefficient(col, obj[i]);
382 const int end = (nnonz == 0 || i == ncols - 1) ? nnonz : beg[i + 1];
383 while ( nz < end )
384 {
385 lpi->linear_program->SetCoefficient(RowIndex(ind[nz]), col, val[nz]);
386 ++nz;
387 }
388 }
389 assert( nz == nnonz );
390 }
391 else
392 {
393 for (int i = 0; i < ncols; ++i)
394 {
395 const ColIndex col = lpi->linear_program->CreateNewVariable();
396 lpi->linear_program->SetVariableBounds(col, lb[i], ub[i]);
397 lpi->linear_program->SetObjectiveCoefficient(col, obj[i]);
398 }
399 }
400
401 lpi->lp_modified_since_last_solve = true;
402
403 return SCIP_OKAY;
404 }
405
406 /** deletes all columns in the given range from LP */
SCIPlpiDelCols(SCIP_LPI * lpi,int firstcol,int lastcol)407 SCIP_RETCODE SCIPlpiDelCols(
408 SCIP_LPI* lpi, /**< LP interface structure */
409 int firstcol, /**< first column to be deleted */
410 int lastcol /**< last column to be deleted */
411 )
412 {
413 assert( lpi != NULL );
414 assert( lpi->linear_program != NULL );
415 assert( 0 <= firstcol && firstcol <= lastcol && lastcol < lpi->linear_program->num_variables() );
416
417 SCIPdebugMessage("deleting columns %d to %d.\n", firstcol, lastcol);
418
419 const ColIndex num_cols = lpi->linear_program->num_variables();
420 DenseBooleanRow columns_to_delete(num_cols, false);
421 for (int i = firstcol; i <= lastcol; ++i)
422 columns_to_delete[ColIndex(i)] = true;
423
424 lpi->linear_program->DeleteColumns(columns_to_delete);
425 lpi->lp_modified_since_last_solve = true;
426
427 return SCIP_OKAY;
428 }
429
430 /** deletes columns from SCIP_LP; the new position of a column must not be greater that its old position */
SCIPlpiDelColset(SCIP_LPI * lpi,int * dstat)431 SCIP_RETCODE SCIPlpiDelColset(
432 SCIP_LPI* lpi, /**< LP interface structure */
433 int* dstat /**< deletion status of columns
434 * input: 1 if column should be deleted, 0 if not
435 * output: new position of column, -1 if column was deleted */
436 )
437 {
438 assert( lpi != NULL );
439 assert( lpi->linear_program != NULL );
440 assert( dstat != NULL );
441
442 const ColIndex num_cols = lpi->linear_program->num_variables();
443 DenseBooleanRow columns_to_delete(num_cols, false);
444 int new_index = 0;
445 int num_deleted_columns = 0;
446 for (ColIndex col(0); col < num_cols; ++col)
447 {
448 int i = col.value();
449 if ( dstat[i] == 1 )
450 {
451 columns_to_delete[col] = true;
452 dstat[i] = -1;
453 ++num_deleted_columns;
454 }
455 else
456 dstat[i] = new_index++;
457 }
458 SCIPdebugMessage("SCIPlpiDelColset: deleting %d columns.\n", num_deleted_columns);
459 lpi->linear_program->DeleteColumns(columns_to_delete);
460 lpi->lp_modified_since_last_solve = true;
461
462 return SCIP_OKAY;
463 }
464
465 /** adds rows to the LP */
SCIPlpiAddRows(SCIP_LPI * lpi,int nrows,const SCIP_Real * lhs,const SCIP_Real * rhs,char ** rownames,int nnonz,const int * beg,const int * ind,const SCIP_Real * val)466 SCIP_RETCODE SCIPlpiAddRows(
467 SCIP_LPI* lpi, /**< LP interface structure */
468 int nrows, /**< number of rows to be added */
469 const SCIP_Real* lhs, /**< left hand sides of new rows */
470 const SCIP_Real* rhs, /**< right hand sides of new rows */
471 char** rownames, /**< row names, or NULL */
472 int nnonz, /**< number of nonzero elements to be added to the constraint matrix */
473 const int* beg, /**< start index of each row in ind- and val-array, or NULL if nnonz == 0 */
474 const int* ind, /**< column indices of constraint matrix entries, or NULL if nnonz == 0 */
475 const SCIP_Real* val /**< values of constraint matrix entries, or NULL if nnonz == 0 */
476 )
477 {
478 assert( lpi != NULL );
479 assert( lpi->linear_program != NULL );
480 assert( lhs != NULL );
481 assert( rhs != NULL );
482 assert( nnonz >= 0) ;
483 assert( nrows >= 0) ;
484
485 SCIPdebugMessage("adding %d rows with %d nonzeros.\n", nrows, nnonz);
486
487 /* @todo add names */
488 if ( nnonz > 0 )
489 {
490 assert( beg != NULL );
491 assert( ind != NULL );
492 assert( val != NULL );
493 assert( nrows > 0 );
494
495 #ifndef NDEBUG
496 /* perform check that no new columns are added - this is likely to be a mistake */
497 const ColIndex num_cols = lpi->linear_program->num_variables();
498 for (int j = 0; j < nnonz; ++j)
499 {
500 assert( val[j] != 0.0 );
501 assert( 0 <= ind[j] && ind[j] < num_cols.value() );
502 }
503 #endif
504
505 int nz = 0;
506 for (int i = 0; i < nrows; ++i)
507 {
508 const RowIndex row = lpi->linear_program->CreateNewConstraint();
509 lpi->linear_program->SetConstraintBounds(row, lhs[i], rhs[i]);
510 const int end = (nnonz == 0 || i == nrows - 1) ? nnonz : beg[i + 1];
511 while ( nz < end )
512 {
513 lpi->linear_program->SetCoefficient(row, ColIndex(ind[nz]), val[nz]);
514 ++nz;
515 }
516 }
517 assert( nz == nnonz );
518 }
519 else
520 {
521 for (int i = 0; i < nrows; ++i)
522 {
523 const RowIndex row = lpi->linear_program->CreateNewConstraint();
524 lpi->linear_program->SetConstraintBounds(row, lhs[i], rhs[i]);
525 }
526 }
527
528 lpi->lp_modified_since_last_solve = true;
529
530 return SCIP_OKAY;
531 }
532
533 /** delete rows from LP and update the current basis */
534 static
deleteRowsAndUpdateCurrentBasis(SCIP_LPI * lpi,const DenseBooleanColumn & rows_to_delete)535 void deleteRowsAndUpdateCurrentBasis(
536 SCIP_LPI* lpi, /**< LP interface structure */
537 const DenseBooleanColumn& rows_to_delete /**< array to mark rows that should be deleted */
538 )
539 {
540 const RowIndex num_rows = lpi->linear_program->num_constraints();
541 const ColIndex num_cols = lpi->linear_program->num_variables();
542
543 /* try to repair basis status if problem size has not changed before */
544 BasisState state = lpi->solver->GetState();
545 if ( state.statuses.size() == num_cols.value() + num_rows.value() )
546 {
547 /* Shift the status of the non-deleted rows. Note that if the deleted rows where part of the basis (i.e., constraint
548 * not tight), then we should be left with a correct basis afterward. This should be the most common use case in SCIP. */
549 ColIndex new_size = num_cols;
550 for (RowIndex row(0); row < num_rows; ++row)
551 {
552 if ( rows_to_delete[row] )
553 continue;
554 state.statuses[new_size++] = state.statuses[num_cols + RowToColIndex(row)];
555 }
556 state.statuses.resize(new_size);
557 lpi->solver->LoadStateForNextSolve(state);
558 }
559
560 lpi->linear_program->DeleteRows(rows_to_delete);
561 lpi->lp_modified_since_last_solve = true;
562 }
563
564 /** deletes all rows in the given range from LP */
SCIPlpiDelRows(SCIP_LPI * lpi,int firstrow,int lastrow)565 SCIP_RETCODE SCIPlpiDelRows(
566 SCIP_LPI* lpi, /**< LP interface structure */
567 int firstrow, /**< first row to be deleted */
568 int lastrow /**< last row to be deleted */
569 )
570 {
571 assert( lpi != NULL );
572 assert( lpi->linear_program != NULL );
573 assert( 0 <= firstrow && firstrow <= lastrow && lastrow < lpi->linear_program->num_constraints() );
574
575 const RowIndex num_rows = lpi->linear_program->num_constraints();
576 DenseBooleanColumn rows_to_delete(num_rows, false);
577 for (int i = firstrow; i <= lastrow; ++i)
578 rows_to_delete[RowIndex(i)] = true;
579
580 SCIPdebugMessage("deleting rows %d to %d.\n", firstrow, lastrow);
581 deleteRowsAndUpdateCurrentBasis(lpi, rows_to_delete);
582
583 return SCIP_OKAY;
584 }
585
586 /** deletes rows from SCIP_LP; the new position of a row must not be greater that its old position */
SCIPlpiDelRowset(SCIP_LPI * lpi,int * dstat)587 SCIP_RETCODE SCIPlpiDelRowset(
588 SCIP_LPI* lpi, /**< LP interface structure */
589 int* dstat /**< deletion status of rows
590 * input: 1 if row should be deleted, 0 if not
591 * output: new position of row, -1 if row was deleted */
592 )
593 {
594 assert( lpi != NULL );
595 assert( lpi->linear_program != NULL );
596
597 const RowIndex num_rows = lpi->linear_program->num_constraints();
598 DenseBooleanColumn rows_to_delete(num_rows, false);
599 int new_index = 0;
600 int num_deleted_rows = 0;
601 for (RowIndex row(0); row < num_rows; ++row)
602 {
603 int i = row.value();
604 if (dstat[i] == 1)
605 {
606 rows_to_delete[row] = true;
607 dstat[i] = -1;
608 ++num_deleted_rows;
609 }
610 else
611 dstat[i] = new_index++;
612 }
613
614 SCIPdebugMessage("SCIPlpiDelRowset: deleting %d rows.\n", num_deleted_rows);
615 deleteRowsAndUpdateCurrentBasis(lpi, rows_to_delete);
616
617 return SCIP_OKAY;
618 }
619
620 /** clears the whole LP */
SCIPlpiClear(SCIP_LPI * lpi)621 SCIP_RETCODE SCIPlpiClear(
622 SCIP_LPI* lpi /**< LP interface structure */
623 )
624 {
625 assert( lpi != NULL );
626 assert( lpi->linear_program != NULL );
627
628 SCIPdebugMessage("SCIPlpiClear\n");
629
630 lpi->linear_program->Clear();
631 lpi->lp_modified_since_last_solve = true;
632
633 return SCIP_OKAY;
634 }
635
636 /** changes lower and upper bounds of columns */
SCIPlpiChgBounds(SCIP_LPI * lpi,int ncols,const int * ind,const SCIP_Real * lb,const SCIP_Real * ub)637 SCIP_RETCODE SCIPlpiChgBounds(
638 SCIP_LPI* lpi, /**< LP interface structure */
639 int ncols, /**< number of columns to change bounds for */
640 const int* ind, /**< column indices or NULL if ncols is zero */
641 const SCIP_Real* lb, /**< values for the new lower bounds or NULL if ncols is zero */
642 const SCIP_Real* ub /**< values for the new upper bounds or NULL if ncols is zero */
643 )
644 {
645 assert( lpi != NULL );
646 assert( lpi->linear_program != NULL );
647 assert( ncols == 0 || (ind != NULL && lb != NULL && ub != NULL) );
648
649 SCIPdebugMessage("changing %d bounds.\n", ncols);
650 if ( ncols <= 0 )
651 return SCIP_OKAY;
652
653 for (int i = 0; i < ncols; ++i)
654 {
655 SCIPdebugMessage(" col %d: [%g,%g]\n", ind[i], lb[i], ub[i]);
656
657 if ( SCIPlpiIsInfinity(lpi, lb[i]) )
658 {
659 SCIPerrorMessage("LP Error: fixing lower bound for variable %d to infinity.\n", ind[i]);
660 return SCIP_LPERROR;
661 }
662 if ( SCIPlpiIsInfinity(lpi, -ub[i]) )
663 {
664 SCIPerrorMessage("LP Error: fixing upper bound for variable %d to -infinity.\n", ind[i]);
665 return SCIP_LPERROR;
666 }
667
668 lpi->linear_program->SetVariableBounds(ColIndex(ind[i]), lb[i], ub[i]);
669 }
670 lpi->lp_modified_since_last_solve = true;
671
672 return SCIP_OKAY;
673 }
674
675 /** changes left and right hand sides of rows */
SCIPlpiChgSides(SCIP_LPI * lpi,int nrows,const int * ind,const SCIP_Real * lhs,const SCIP_Real * rhs)676 SCIP_RETCODE SCIPlpiChgSides(
677 SCIP_LPI* lpi, /**< LP interface structure */
678 int nrows, /**< number of rows to change sides for */
679 const int* ind, /**< row indices */
680 const SCIP_Real* lhs, /**< new values for left hand sides */
681 const SCIP_Real* rhs /**< new values for right hand sides */
682 )
683 {
684 assert( lpi != NULL );
685 assert( lpi->linear_program != NULL );
686
687 if ( nrows <= 0 )
688 return SCIP_OKAY;
689
690 SCIPdebugMessage("changing %d sides\n", nrows);
691
692 for (int i = 0; i < nrows; ++i)
693 lpi->linear_program->SetConstraintBounds(RowIndex(ind[i]), lhs[i], rhs[i]);
694
695 lpi->lp_modified_since_last_solve = true;
696
697 return SCIP_OKAY;
698 }
699
700 /** changes a single coefficient */
SCIPlpiChgCoef(SCIP_LPI * lpi,int row,int col,SCIP_Real newval)701 SCIP_RETCODE SCIPlpiChgCoef(
702 SCIP_LPI* lpi, /**< LP interface structure */
703 int row, /**< row number of coefficient to change */
704 int col, /**< column number of coefficient to change */
705 SCIP_Real newval /**< new value of coefficient */
706 )
707 {
708 assert( lpi != NULL );
709 assert( lpi->linear_program != NULL );
710 assert( 0 <= row && row < lpi->linear_program->num_constraints().value() );
711 assert( 0 <= col && col < lpi->linear_program->num_variables().value() );
712
713 SCIPdebugMessage("Set coefficient (%d,%d) to %f.\n", row, col, newval);
714 lpi->linear_program->CleanUp();
715 assert( lpi->linear_program->IsCleanedUp() );
716 lpi->linear_program->SetCoefficient(RowIndex(row), ColIndex(col), newval);
717
718 lpi->lp_modified_since_last_solve = true;
719
720 return SCIP_OKAY;
721 }
722
723 /** changes the objective sense */
SCIPlpiChgObjsen(SCIP_LPI * lpi,SCIP_OBJSEN objsen)724 SCIP_RETCODE SCIPlpiChgObjsen(
725 SCIP_LPI* lpi, /**< LP interface structure */
726 SCIP_OBJSEN objsen /**< new objective sense */
727 )
728 {
729 assert( lpi != NULL);
730 assert( lpi->linear_program != NULL);
731
732 SCIPdebugMessage("changing objective sense to %d\n", objsen);
733
734 switch (objsen)
735 {
736 case SCIP_OBJSEN_MAXIMIZE:
737 lpi->linear_program->SetMaximizationProblem(true);
738 break;
739 case SCIP_OBJSEN_MINIMIZE:
740 lpi->linear_program->SetMaximizationProblem(false);
741 break;
742 }
743 lpi->lp_modified_since_last_solve = true;
744
745 return SCIP_OKAY;
746 }
747
748 /** changes objective values of columns in the LP */
SCIPlpiChgObj(SCIP_LPI * lpi,int ncols,const int * ind,const SCIP_Real * obj)749 SCIP_RETCODE SCIPlpiChgObj(
750 SCIP_LPI* lpi, /**< LP interface structure */
751 int ncols, /**< number of columns to change objective value for */
752 const int* ind, /**< column indices to change objective value for */
753 const SCIP_Real* obj /**< new objective values for columns */
754 )
755 {
756 assert( lpi != NULL );
757 assert( lpi->linear_program != NULL );
758 assert( ind != NULL );
759 assert( obj != NULL );
760
761 SCIPdebugMessage("changing %d objective values\n", ncols);
762
763 for (int i = 0; i < ncols; ++i)
764 lpi->linear_program->SetObjectiveCoefficient(ColIndex(ind[i]), obj[i]);
765
766 lpi->lp_modified_since_last_solve = true;
767
768 return SCIP_OKAY;
769 }
770
771 /** multiplies a row with a non-zero scalar; for negative scalars, the row's sense is switched accordingly */
SCIPlpiScaleRow(SCIP_LPI * lpi,int row,SCIP_Real scaleval)772 SCIP_RETCODE SCIPlpiScaleRow(
773 SCIP_LPI* lpi, /**< LP interface structure */
774 int row, /**< row number to scale */
775 SCIP_Real scaleval /**< scaling multiplier */
776 )
777 {
778 SCIP_Real* vals;
779 SCIP_Real lhs;
780 SCIP_Real rhs;
781 int nnonz;
782 int* inds;
783 int beg;
784
785 assert( lpi != NULL );
786 assert( lpi->linear_program != NULL );
787
788 SCIPdebugMessage("Scale row %d by %f.\n", row, scaleval);
789
790 /* alloc memory */
791 ColIndex num_cols = lpi->linear_program->num_variables();
792
793 SCIP_ALLOC( BMSallocMemoryArray(&inds, num_cols.value()) );
794 SCIP_ALLOC( BMSallocMemoryArray(&vals, num_cols.value()) );
795
796 /* get the row */
797 SCIP_CALL( SCIPlpiGetRows(lpi, row, row, &lhs, &rhs, &nnonz, &beg, inds, vals) );
798
799 /* scale row coefficients */
800 for (int j = 0; j < nnonz; ++j)
801 {
802 SCIP_CALL( SCIPlpiChgCoef(lpi, row, inds[j], vals[j] * scaleval) );
803 }
804 BMSfreeMemoryArray(&vals);
805 BMSfreeMemoryArray(&inds);
806
807 /* scale row sides */
808 if ( ! SCIPlpiIsInfinity(lpi, -lhs) )
809 lhs *= scaleval;
810 else if ( scaleval < 0.0 )
811 lhs = SCIPlpiInfinity(lpi);
812
813 if ( ! SCIPlpiIsInfinity(lpi, rhs) )
814 rhs *= scaleval;
815 else if ( scaleval < 0.0 )
816 rhs = -SCIPlpiInfinity(lpi);
817
818 if ( scaleval > 0.0 )
819 {
820 SCIP_CALL( SCIPlpiChgSides(lpi, 1, &row, &lhs, &rhs) );
821 }
822 else
823 {
824 SCIP_CALL( SCIPlpiChgSides(lpi, 1, &row, &rhs, &lhs) );
825 }
826
827 lpi->lp_modified_since_last_solve = true;
828
829 return SCIP_OKAY;
830 }
831
832 /** multiplies a column with a non-zero scalar; the objective value is multiplied with the scalar, and the bounds
833 * are divided by the scalar; for negative scalars, the column's bounds are switched
834 */
SCIPlpiScaleCol(SCIP_LPI * lpi,int col,SCIP_Real scaleval)835 SCIP_RETCODE SCIPlpiScaleCol(
836 SCIP_LPI* lpi, /**< LP interface structure */
837 int col, /**< column number to scale */
838 SCIP_Real scaleval /**< scaling multiplier */
839 )
840 {
841 SCIP_Real* vals;
842 SCIP_Real lb;
843 SCIP_Real ub;
844 SCIP_Real obj;
845 int nnonz;
846 int* inds;
847 int beg;
848
849 assert( lpi != NULL );
850 assert( lpi->linear_program != NULL );
851
852 SCIPdebugMessage("Scale column %d by %f.\n", col, scaleval);
853
854 /* alloc memory */
855 RowIndex num_rows = lpi->linear_program->num_constraints();
856
857 SCIP_ALLOC( BMSallocMemoryArray(&inds, num_rows.value()) );
858 SCIP_ALLOC( BMSallocMemoryArray(&vals, num_rows.value()) );
859
860 /* get the column */
861 SCIP_CALL( SCIPlpiGetCols(lpi, col, col, &lb, &ub, &nnonz, &beg, inds, vals) );
862
863 /* scale column coefficients */
864 for (int j = 0; j < nnonz; ++j)
865 {
866 SCIP_CALL( SCIPlpiChgCoef(lpi, col, inds[j], vals[j] * scaleval) );
867 }
868 BMSfreeMemoryArray(&vals);
869 BMSfreeMemoryArray(&inds);
870
871 /* scale objective value */
872 SCIP_CALL( SCIPlpiGetObj(lpi, col, col, &obj) );
873 obj *= scaleval;
874 SCIP_CALL( SCIPlpiChgObj(lpi, 1, &col, &obj) );
875
876 /* scale bound */
877 if ( ! SCIPlpiIsInfinity(lpi, -lb) )
878 lb *= scaleval;
879 else if ( scaleval < 0.0 )
880 lb = SCIPlpiInfinity(lpi);
881
882 if ( ! SCIPlpiIsInfinity(lpi, ub) )
883 ub *= scaleval;
884 else if ( scaleval < 0.0 )
885 ub = -SCIPlpiInfinity(lpi);
886
887 if ( scaleval > 0.0 )
888 {
889 SCIP_CALL( SCIPlpiChgBounds(lpi, 1, &col, &lb, &ub) );
890 }
891 else
892 {
893 SCIP_CALL( SCIPlpiChgBounds(lpi, 1, &col, &ub, &lb) );
894 }
895
896 return SCIP_OKAY;
897 }
898
899
900 /**@} */
901
902
903
904
905 /*
906 * Data Accessing Methods
907 */
908
909 /**@name Data Accessing Methods */
910 /**@{ */
911
912 /** gets the number of rows in the LP */
SCIPlpiGetNRows(SCIP_LPI * lpi,int * nrows)913 SCIP_RETCODE SCIPlpiGetNRows(
914 SCIP_LPI* lpi, /**< LP interface structure */
915 int* nrows /**< pointer to store the number of rows */
916 )
917 {
918 assert( lpi != NULL );
919 assert( lpi->linear_program != NULL );
920 assert( nrows != NULL );
921
922 SCIPdebugMessage("getting number of rows.\n");
923
924 *nrows = lpi->linear_program->num_constraints().value();
925
926 return SCIP_OKAY;
927 }
928
929 /** gets the number of columns in the LP */
SCIPlpiGetNCols(SCIP_LPI * lpi,int * ncols)930 SCIP_RETCODE SCIPlpiGetNCols(
931 SCIP_LPI* lpi, /**< LP interface structure */
932 int* ncols /**< pointer to store the number of cols */
933 )
934 {
935 assert( lpi != NULL );
936 assert( lpi->linear_program != NULL );
937 assert( ncols != NULL );
938
939 SCIPdebugMessage("getting number of columns.\n");
940
941 *ncols = lpi->linear_program->num_variables().value();
942
943 return SCIP_OKAY;
944 }
945
946 /** gets objective sense of the LP */
SCIPlpiGetObjsen(SCIP_LPI * lpi,SCIP_OBJSEN * objsen)947 SCIP_RETCODE SCIPlpiGetObjsen(
948 SCIP_LPI* lpi, /**< LP interface structure */
949 SCIP_OBJSEN* objsen /**< pointer to store objective sense */
950 )
951 {
952 assert( lpi != NULL );
953 assert( lpi->linear_program != NULL );
954 assert( objsen != NULL );
955
956 SCIPdebugMessage("getting objective sense.\n");
957
958 *objsen = lpi->linear_program->IsMaximizationProblem() ? SCIP_OBJSEN_MAXIMIZE : SCIP_OBJSEN_MINIMIZE;
959
960 return SCIP_OKAY;
961 }
962
963 /** gets the number of nonzero elements in the LP constraint matrix */
SCIPlpiGetNNonz(SCIP_LPI * lpi,int * nnonz)964 SCIP_RETCODE SCIPlpiGetNNonz(
965 SCIP_LPI* lpi, /**< LP interface structure */
966 int* nnonz /**< pointer to store the number of nonzeros */
967 )
968 {
969 assert( lpi != NULL );
970 assert( lpi->linear_program != NULL );
971 assert( nnonz != NULL );
972
973 SCIPdebugMessage("getting number of non-zeros.\n");
974
975 *nnonz = (int) lpi->linear_program->num_entries().value();
976
977 return SCIP_OKAY;
978 }
979
980 /** gets columns from LP problem object; the arrays have to be large enough to store all values
981 * Either both, lb and ub, have to be NULL, or both have to be non-NULL,
982 * either nnonz, beg, ind, and val have to be NULL, or all of them have to be non-NULL.
983 */
SCIPlpiGetCols(SCIP_LPI * lpi,int firstcol,int lastcol,SCIP_Real * lb,SCIP_Real * ub,int * nnonz,int * beg,int * ind,SCIP_Real * val)984 SCIP_RETCODE SCIPlpiGetCols(
985 SCIP_LPI* lpi, /**< LP interface structure */
986 int firstcol, /**< first column to get from LP */
987 int lastcol, /**< last column to get from LP */
988 SCIP_Real* lb, /**< buffer to store the lower bound vector, or NULL */
989 SCIP_Real* ub, /**< buffer to store the upper bound vector, or NULL */
990 int* nnonz, /**< pointer to store the number of nonzero elements returned, or NULL */
991 int* beg, /**< buffer to store start index of each column in ind- and val-array, or NULL */
992 int* ind, /**< buffer to store row indices of constraint matrix entries, or NULL */
993 SCIP_Real* val /**< buffer to store values of constraint matrix entries, or NULL */
994 )
995 {
996 assert( lpi != NULL );
997 assert( lpi->linear_program != NULL );
998 assert( 0 <= firstcol && firstcol <= lastcol && lastcol < lpi->linear_program->num_variables() );
999 assert( (lb != NULL && ub != NULL) || (lb == NULL && ub == NULL) );
1000 assert( (nnonz != NULL && beg != NULL && ind != NULL && val != NULL) || (nnonz == NULL && beg == NULL && ind == NULL && val == NULL) );
1001
1002 const DenseRow& tmplb = lpi->linear_program->variable_lower_bounds();
1003 const DenseRow& tmpub = lpi->linear_program->variable_upper_bounds();
1004
1005 if ( nnonz != NULL )
1006 {
1007 assert( beg != NULL );
1008 assert( ind != NULL );
1009 assert( val != NULL );
1010
1011 *nnonz = 0;
1012 int index = 0;
1013 for (ColIndex col(firstcol); col <= ColIndex(lastcol); ++col, ++index)
1014 {
1015 if ( lb != NULL )
1016 lb[index] = tmplb[col];
1017 if ( ub != NULL )
1018 ub[index] = tmpub[col];
1019
1020 beg[index] = *nnonz;
1021 const SparseColumn& column = lpi->linear_program->GetSparseColumn(col);
1022 for (const SparseColumn::Entry& entry : column)
1023 {
1024 const RowIndex row = entry.row();
1025 ind[*nnonz] = row.value();
1026 val[*nnonz] = entry.coefficient();
1027 ++(*nnonz);
1028 }
1029 }
1030 }
1031 else
1032 {
1033 int index = 0;
1034 for (ColIndex col(firstcol); col <= ColIndex(lastcol); ++col, ++index)
1035 {
1036 if ( lb != NULL )
1037 lb[index] = tmplb[col];
1038 if ( ub != NULL )
1039 ub[index] = tmpub[col];
1040 }
1041 }
1042
1043 return SCIP_OKAY;
1044 }
1045
1046 /** gets rows from LP problem object; the arrays have to be large enough to store all values.
1047 * Either both, lhs and rhs, have to be NULL, or both have to be non-NULL,
1048 * either nnonz, beg, ind, and val have to be NULL, or all of them have to be non-NULL.
1049 */
SCIPlpiGetRows(SCIP_LPI * lpi,int firstrow,int lastrow,SCIP_Real * lhs,SCIP_Real * rhs,int * nnonz,int * beg,int * ind,SCIP_Real * val)1050 SCIP_RETCODE SCIPlpiGetRows(
1051 SCIP_LPI* lpi, /**< LP interface structure */
1052 int firstrow, /**< first row to get from LP */
1053 int lastrow, /**< last row to get from LP */
1054 SCIP_Real* lhs, /**< buffer to store left hand side vector, or NULL */
1055 SCIP_Real* rhs, /**< buffer to store right hand side vector, or NULL */
1056 int* nnonz, /**< pointer to store the number of nonzero elements returned, or NULL */
1057 int* beg, /**< buffer to store start index of each row in ind- and val-array, or NULL */
1058 int* ind, /**< buffer to store column indices of constraint matrix entries, or NULL */
1059 SCIP_Real* val /**< buffer to store values of constraint matrix entries, or NULL */
1060 )
1061 {
1062 assert( lpi != NULL );
1063 assert( lpi->linear_program != NULL );
1064 assert( 0 <= firstrow && firstrow <= lastrow && lastrow < lpi->linear_program->num_constraints() );
1065 assert( (lhs == NULL && rhs == NULL) || (rhs != NULL && lhs != NULL) );
1066 assert( (nnonz != NULL && beg != NULL && ind != NULL && val != NULL) || (nnonz == NULL && beg == NULL && ind == NULL && val == NULL) );
1067
1068 const DenseColumn& tmplhs = lpi->linear_program->constraint_lower_bounds();
1069 const DenseColumn& tmprhs = lpi->linear_program->constraint_upper_bounds();
1070
1071 if ( nnonz != NULL )
1072 {
1073 assert( beg != NULL );
1074 assert( ind != NULL );
1075 assert( val != NULL );
1076
1077 const SparseMatrix& matrixtrans = lpi->linear_program->GetTransposeSparseMatrix();
1078
1079 *nnonz = 0;
1080 int index = 0;
1081 for (RowIndex row(firstrow); row <= RowIndex(lastrow); ++row, ++index)
1082 {
1083 if ( lhs != NULL )
1084 lhs[index] = tmplhs[row];
1085 if ( rhs != NULL )
1086 rhs[index] = tmprhs[row];
1087
1088 beg[index] = *nnonz;
1089 const SparseColumn& column = matrixtrans.column(ColIndex(row.value()));
1090 for (const SparseColumn::Entry& entry : column)
1091 {
1092 const RowIndex rowidx = entry.row();
1093 ind[*nnonz] = rowidx.value();
1094 val[*nnonz] = entry.coefficient();
1095 ++(*nnonz);
1096 }
1097 }
1098 }
1099 else
1100 {
1101 int index = 0;
1102 for (RowIndex row(firstrow); row <= RowIndex(lastrow); ++row, ++index)
1103 {
1104 if ( lhs != NULL )
1105 lhs[index] = tmplhs[row];
1106 if ( rhs != NULL )
1107 rhs[index] = tmprhs[row];
1108 }
1109 }
1110
1111 return SCIP_OKAY;
1112 }
1113
1114 /** gets column names */
SCIPlpiGetColNames(SCIP_LPI * lpi,int firstcol,int lastcol,char ** colnames,char * namestorage,int namestoragesize,int * storageleft)1115 SCIP_RETCODE SCIPlpiGetColNames(
1116 SCIP_LPI* lpi, /**< LP interface structure */
1117 int firstcol, /**< first column to get name from LP */
1118 int lastcol, /**< last column to get name from LP */
1119 char** colnames, /**< pointers to column names (of size at least lastcol-firstcol+1) or NULL if namestoragesize is zero */
1120 char* namestorage, /**< storage for col names or NULL if namestoragesize is zero */
1121 int namestoragesize, /**< size of namestorage (if 0, storageleft returns the storage needed) */
1122 int* storageleft /**< amount of storage left (if < 0 the namestorage was not big enough) or NULL if namestoragesize is zero */
1123 )
1124 {
1125 assert( lpi != NULL );
1126 assert( lpi->linear_program != NULL );
1127 assert( colnames != NULL || namestoragesize == 0 );
1128 assert( namestorage != NULL || namestoragesize == 0 );
1129 assert( namestoragesize >= 0 );
1130 assert( storageleft != NULL );
1131 assert( 0 <= firstcol && firstcol <= lastcol && lastcol < lpi->linear_program->num_variables() );
1132
1133 SCIPerrorMessage("SCIPlpiGetColNames() has not been implemented yet.\n");
1134
1135 return SCIP_NOTIMPLEMENTED;
1136 }
1137
1138 /** gets row names */
SCIPlpiGetRowNames(SCIP_LPI * lpi,int firstrow,int lastrow,char ** rownames,char * namestorage,int namestoragesize,int * storageleft)1139 SCIP_RETCODE SCIPlpiGetRowNames(
1140 SCIP_LPI* lpi, /**< LP interface structure */
1141 int firstrow, /**< first row to get name from LP */
1142 int lastrow, /**< last row to get name from LP */
1143 char** rownames, /**< pointers to row names (of size at least lastrow-firstrow+1) or NULL if namestoragesize is zero */
1144 char* namestorage, /**< storage for row names or NULL if namestoragesize is zero */
1145 int namestoragesize, /**< size of namestorage (if 0, -storageleft returns the storage needed) */
1146 int* storageleft /**< amount of storage left (if < 0 the namestorage was not big enough) or NULL if namestoragesize is zero */
1147 )
1148 {
1149 assert( lpi != NULL );
1150 assert( lpi->linear_program != NULL );
1151 assert( rownames != NULL || namestoragesize == 0 );
1152 assert( namestorage != NULL || namestoragesize == 0 );
1153 assert( namestoragesize >= 0 );
1154 assert( storageleft != NULL );
1155 assert( 0 <= firstrow && firstrow <= lastrow && lastrow < lpi->linear_program->num_constraints() );
1156
1157 SCIPerrorMessage("SCIPlpiGetRowNames() has not been implemented yet.\n");
1158
1159 return SCIP_NOTIMPLEMENTED;
1160 }
1161
1162 /** gets objective coefficients from LP problem object */
SCIPlpiGetObj(SCIP_LPI * lpi,int firstcol,int lastcol,SCIP_Real * vals)1163 SCIP_RETCODE SCIPlpiGetObj(
1164 SCIP_LPI* lpi, /**< LP interface structure */
1165 int firstcol, /**< first column to get objective coefficient for */
1166 int lastcol, /**< last column to get objective coefficient for */
1167 SCIP_Real* vals /**< array to store objective coefficients */
1168 )
1169 {
1170 assert( lpi != NULL );
1171 assert( lpi->linear_program != NULL );
1172 assert( firstcol <= lastcol );
1173 assert( vals != NULL );
1174
1175 SCIPdebugMessage("getting objective values %d to %d\n", firstcol, lastcol);
1176
1177 int index = 0;
1178 for (ColIndex col(firstcol); col <= ColIndex(lastcol); ++col)
1179 {
1180 vals[index] = lpi->linear_program->objective_coefficients()[col];
1181 ++index;
1182 }
1183
1184 return SCIP_OKAY;
1185 }
1186
1187 /** gets current bounds from LP problem object */
SCIPlpiGetBounds(SCIP_LPI * lpi,int firstcol,int lastcol,SCIP_Real * lbs,SCIP_Real * ubs)1188 SCIP_RETCODE SCIPlpiGetBounds(
1189 SCIP_LPI* lpi, /**< LP interface structure */
1190 int firstcol, /**< first column to get bounds for */
1191 int lastcol, /**< last column to get bounds for */
1192 SCIP_Real* lbs, /**< array to store lower bound values, or NULL */
1193 SCIP_Real* ubs /**< array to store upper bound values, or NULL */
1194 )
1195 {
1196 assert( lpi != NULL );
1197 assert( lpi->linear_program != NULL );
1198 assert( firstcol <= lastcol );
1199
1200 SCIPdebugMessage("getting bounds %d to %d\n", firstcol, lastcol);
1201
1202 int index = 0;
1203 for (ColIndex col(firstcol); col <= ColIndex(lastcol); ++col)
1204 {
1205 if ( lbs != NULL )
1206 lbs[index] = lpi->linear_program->variable_lower_bounds()[col];
1207
1208 if ( ubs != NULL )
1209 ubs[index] = lpi->linear_program->variable_upper_bounds()[col];
1210
1211 ++index;
1212 }
1213
1214 return SCIP_OKAY;
1215 }
1216
1217 /** gets current row sides from LP problem object */
SCIPlpiGetSides(SCIP_LPI * lpi,int firstrow,int lastrow,SCIP_Real * lhss,SCIP_Real * rhss)1218 SCIP_RETCODE SCIPlpiGetSides(
1219 SCIP_LPI* lpi, /**< LP interface structure */
1220 int firstrow, /**< first row to get sides for */
1221 int lastrow, /**< last row to get sides for */
1222 SCIP_Real* lhss, /**< array to store left hand side values, or NULL */
1223 SCIP_Real* rhss /**< array to store right hand side values, or NULL */
1224 )
1225 {
1226 assert( lpi != NULL );
1227 assert( lpi->linear_program != NULL );
1228 assert( firstrow <= lastrow );
1229
1230 SCIPdebugMessage("getting row sides %d to %d\n", firstrow, lastrow);
1231
1232 int index = 0;
1233 for (RowIndex row(firstrow); row <= RowIndex(lastrow); ++row)
1234 {
1235 if ( lhss != NULL )
1236 lhss[index] = lpi->linear_program->constraint_lower_bounds()[row];
1237
1238 if ( rhss != NULL )
1239 rhss[index] = lpi->linear_program->constraint_upper_bounds()[row];
1240
1241 ++index;
1242 }
1243
1244 return SCIP_OKAY;
1245 }
1246
1247 /** gets a single coefficient */
SCIPlpiGetCoef(SCIP_LPI * lpi,int row,int col,SCIP_Real * val)1248 SCIP_RETCODE SCIPlpiGetCoef(
1249 SCIP_LPI* lpi, /**< LP interface structure */
1250 int row, /**< row number of coefficient */
1251 int col, /**< column number of coefficient */
1252 SCIP_Real* val /**< pointer to store the value of the coefficient */
1253 )
1254 {
1255 assert( lpi != NULL );
1256 assert( lpi->linear_program != NULL );
1257 assert( val != NULL );
1258
1259 /* quite slow method: possibly needs linear time if matrix is not sorted */
1260 const SparseMatrix& matrix = lpi->linear_program->GetSparseMatrix();
1261 *val = matrix.LookUpValue(RowIndex(row), ColIndex(col));
1262
1263 return SCIP_OKAY;
1264 }
1265
1266 /**@} */
1267
1268
1269
1270
1271 /*
1272 * Solving Methods
1273 */
1274
1275 /**@name Solving Methods */
1276 /**@{ */
1277
1278 /** update scaled linear program */
1279 static
updateScaledLP(SCIP_LPI * lpi)1280 void updateScaledLP(
1281 SCIP_LPI* lpi /**< LP interface structure */
1282 )
1283 {
1284 if ( ! lpi->lp_modified_since_last_solve )
1285 return;
1286
1287 lpi->scaled_lp->PopulateFromLinearProgram(*lpi->linear_program);
1288 lpi->scaled_lp->AddSlackVariablesWhereNecessary(false);
1289
1290 /* @todo: Avoid doing a copy if there is no scaling. */
1291 /* @todo: Avoid rescaling if not much changed. */
1292 if ( lpi->parameters->use_scaling() )
1293 lpi->scaler->Scale(lpi->scaled_lp);
1294 else
1295 lpi->scaler->Clear();
1296 }
1297
1298 /** check primal feasibility */
1299 static
checkUnscaledPrimalFeasibility(SCIP_LPI * lpi)1300 bool checkUnscaledPrimalFeasibility(
1301 SCIP_LPI* lpi /**< LP interface structure */
1302 )
1303 {
1304 assert( lpi != NULL );
1305 assert( lpi->solver != NULL );
1306 assert( lpi->linear_program != NULL );
1307
1308 #if UNSCALEDFEAS_CHECK == 1
1309 /* get unscaled solution */
1310 const ColIndex num_cols = lpi->linear_program->num_variables();
1311 DenseRow unscaledsol(num_cols);
1312 for (ColIndex col = ColIndex(0); col < num_cols; ++col)
1313 unscaledsol[col] = lpi->scaler->UnscaleVariableValue(col, lpi->solver->GetVariableValue(col));
1314
1315 /* if the solution is not feasible w.r.t. absolute tolerances, try to fix it in the unscaled problem */
1316 const double feastol = lpi->parameters->primal_feasibility_tolerance();
1317 return lpi->linear_program->SolutionIsLPFeasible(unscaledsol, feastol);
1318
1319 #elif UNSCALEDFEAS_CHECK == 2
1320 const double feastol = lpi->parameters->primal_feasibility_tolerance();
1321
1322 /* check bounds of unscaled solution */
1323 const ColIndex num_cols = lpi->linear_program->num_variables();
1324 for (ColIndex col = ColIndex(0); col < num_cols; ++col)
1325 {
1326 const Fractional val = lpi->scaler->UnscaleVariableValue(col, lpi->solver->GetVariableValue(col));
1327 const Fractional lb = lpi->linear_program->variable_lower_bounds()[col];
1328 if ( val < lb - feastol )
1329 return false;
1330 const Fractional ub = lpi->linear_program->variable_upper_bounds()[col];
1331 if ( val > ub + feastol )
1332 return false;
1333 }
1334
1335 /* check activities of unscaled solution */
1336 const RowIndex num_rows = lpi->linear_program->num_constraints();
1337 for (RowIndex row(0); row < num_rows; ++row)
1338 {
1339 const Fractional val = lpi->scaler->UnscaleConstraintActivity(row, lpi->solver->GetConstraintActivity(row));
1340 const Fractional lhs = lpi->linear_program->constraint_lower_bounds()[row];
1341 if ( val < lhs - feastol )
1342 return false;
1343 const Fractional rhs = lpi->linear_program->constraint_upper_bounds()[row];
1344 if ( val > rhs + feastol )
1345 return false;
1346 }
1347 #endif
1348
1349 return true;
1350 }
1351
1352 /** common function between the two LPI Solve() functions */
1353 static
SolveInternal(SCIP_LPI * lpi,bool recursive,std::unique_ptr<TimeLimit> & time_limit)1354 SCIP_RETCODE SolveInternal(
1355 SCIP_LPI* lpi, /**< LP interface structure */
1356 bool recursive, /**< Is this a recursive call? */
1357 std::unique_ptr<TimeLimit>& time_limit /**< time limit */
1358 )
1359 {
1360 assert( lpi != NULL );
1361 assert( lpi->solver != NULL );
1362 assert( lpi->parameters != NULL );
1363
1364 updateScaledLP(lpi);
1365
1366 lpi->solver->SetParameters(*(lpi->parameters));
1367 lpi->lp_time_limit_was_reached = false;
1368
1369 /* possibly ignore warm start information for next solve */
1370 if ( lpi->from_scratch )
1371 lpi->solver->ClearStateForNextSolve();
1372
1373 if ( ! lpi->solver->Solve(*(lpi->scaled_lp), time_limit.get()).ok() )
1374 {
1375 return SCIP_LPERROR;
1376 }
1377 lpi->lp_time_limit_was_reached = time_limit->LimitReached();
1378 if ( recursive )
1379 lpi->niterations += (SCIP_Longint) lpi->solver->GetNumberOfIterations();
1380 else
1381 lpi->niterations = (SCIP_Longint) lpi->solver->GetNumberOfIterations();
1382
1383 SCIPdebugMessage("status=%s obj=%f iter=%ld.\n", GetProblemStatusString(lpi->solver->GetProblemStatus()).c_str(),
1384 lpi->solver->GetObjectiveValue(), lpi->solver->GetNumberOfIterations());
1385
1386 const ProblemStatus status = lpi->solver->GetProblemStatus();
1387 if ( (status == ProblemStatus::PRIMAL_FEASIBLE || status == ProblemStatus::OPTIMAL) && lpi->parameters->use_scaling() )
1388 {
1389 if ( ! checkUnscaledPrimalFeasibility(lpi) )
1390 {
1391 SCIPdebugMessage("Solution not feasible w.r.t. absolute tolerance %g -> reoptimize.\n", lpi->parameters->primal_feasibility_tolerance());
1392
1393 /* Re-solve without scaling to try to fix the infeasibility. */
1394 lpi->parameters->set_use_scaling(false);
1395 lpi->lp_modified_since_last_solve = true;
1396 SolveInternal(lpi, true, time_limit); /* inherit time limit, so used time is not reset; do not change iteration limit for resolve */
1397 lpi->parameters->set_use_scaling(true);
1398
1399 #ifdef SCIP_DEBUG
1400 if ( ! checkUnscaledPrimalFeasibility(lpi) )
1401 SCIPdebugMessage("Solution still not feasible after turning off scaling.\n");
1402 #endif
1403 }
1404 }
1405
1406 lpi->lp_modified_since_last_solve = false;
1407
1408 return SCIP_OKAY;
1409 }
1410
1411 /** calls primal simplex to solve the LP */
SCIPlpiSolvePrimal(SCIP_LPI * lpi)1412 SCIP_RETCODE SCIPlpiSolvePrimal(
1413 SCIP_LPI* lpi /**< LP interface structure */
1414 )
1415 {
1416 assert( lpi != NULL );
1417 assert( lpi->solver != NULL );
1418 assert( lpi->linear_program != NULL );
1419 assert( lpi->parameters != NULL );
1420
1421 SCIPdebugMessage("SCIPlpiSolvePrimal: %d rows, %d cols.\n", lpi->linear_program->num_constraints().value(), lpi->linear_program->num_variables().value());
1422 std::unique_ptr<TimeLimit> time_limit = TimeLimit::FromParameters(*lpi->parameters);
1423 lpi->niterations = 0;
1424
1425 lpi->parameters->set_use_dual_simplex(false);
1426 return SolveInternal(lpi, false, time_limit);
1427 }
1428
1429 /** calls dual simplex to solve the LP */
SCIPlpiSolveDual(SCIP_LPI * lpi)1430 SCIP_RETCODE SCIPlpiSolveDual(
1431 SCIP_LPI* lpi /**< LP interface structure */
1432 )
1433 {
1434 assert( lpi != NULL );
1435 assert( lpi->solver != NULL );
1436 assert( lpi->linear_program != NULL );
1437 assert( lpi->parameters != NULL );
1438
1439 SCIPdebugMessage("SCIPlpiSolveDual: %d rows, %d cols.\n", lpi->linear_program->num_constraints().value(), lpi->linear_program->num_variables().value());
1440 std::unique_ptr<TimeLimit> time_limit = TimeLimit::FromParameters(*lpi->parameters);
1441 lpi->niterations = 0;
1442
1443 lpi->parameters->set_use_dual_simplex(true);
1444 return SolveInternal(lpi, false, time_limit);
1445 }
1446
1447 /** calls barrier or interior point algorithm to solve the LP with crossover to simplex basis */
SCIPlpiSolveBarrier(SCIP_LPI * lpi,SCIP_Bool crossover)1448 SCIP_RETCODE SCIPlpiSolveBarrier(
1449 SCIP_LPI* lpi, /**< LP interface structure */
1450 SCIP_Bool crossover /**< perform crossover */
1451 )
1452 {
1453 assert( lpi != NULL );
1454 assert( lpi->solver != NULL );
1455 assert( lpi->linear_program != NULL );
1456 assert( lpi->parameters != NULL );
1457
1458 SCIPerrorMessage("SCIPlpiSolveBarrier - Not supported.\n");
1459
1460 return SCIP_NOTIMPLEMENTED;
1461 }
1462
1463 /** start strong branching */
SCIPlpiStartStrongbranch(SCIP_LPI * lpi)1464 SCIP_RETCODE SCIPlpiStartStrongbranch(
1465 SCIP_LPI* lpi /**< LP interface structure */
1466 )
1467 { /*lint --e{715}*/
1468 assert( lpi != NULL );
1469 assert( lpi->linear_program != NULL );
1470 assert( lpi->solver != NULL );
1471
1472 updateScaledLP(lpi);
1473
1474 /* @todo Save state and do all the branching from there. */
1475 return SCIP_OKAY;
1476 }
1477
1478 /** end strong branching */
SCIPlpiEndStrongbranch(SCIP_LPI * lpi)1479 SCIP_RETCODE SCIPlpiEndStrongbranch(
1480 SCIP_LPI* lpi /**< LP interface structure */
1481 )
1482 { /*lint --e{715}*/
1483 assert( lpi != NULL );
1484 assert( lpi->linear_program != NULL );
1485 assert( lpi->solver != NULL );
1486
1487 /* @todo Restore the saved state. */
1488 return SCIP_OKAY;
1489 }
1490
1491 /** determine whether the dual bound is valid */
1492 static
IsDualBoundValid(ProblemStatus status)1493 bool IsDualBoundValid(
1494 ProblemStatus status /**< status to be checked */
1495 )
1496 {
1497 return status == ProblemStatus::OPTIMAL || status == ProblemStatus::DUAL_FEASIBLE || status == ProblemStatus::DUAL_UNBOUNDED;
1498 }
1499
1500 /** performs strong branching iterations */
1501 static
strongbranch(SCIP_LPI * lpi,int col_index,SCIP_Real psol,int itlim,SCIP_Real * down,SCIP_Real * up,SCIP_Bool * downvalid,SCIP_Bool * upvalid,int * iter)1502 SCIP_RETCODE strongbranch(
1503 SCIP_LPI* lpi, /**< LP interface structure */
1504 int col_index, /**< column to apply strong branching on */
1505 SCIP_Real psol, /**< fractional current primal solution value of column */
1506 int itlim, /**< iteration limit for strong branchings */
1507 SCIP_Real* down, /**< stores dual bound after branching column down */
1508 SCIP_Real* up, /**< stores dual bound after branching column up */
1509 SCIP_Bool* downvalid, /**< stores whether the returned down value is a valid dual bound;
1510 * otherwise, it can only be used as an estimate value */
1511 SCIP_Bool* upvalid, /**< stores whether the returned up value is a valid dual bound;
1512 * otherwise, it can only be used as an estimate value */
1513 int* iter /**< stores total number of strong branching iterations, or -1; may be NULL */
1514 )
1515 {
1516 assert( lpi != NULL );
1517 assert( lpi->scaled_lp != NULL );
1518 assert( down != NULL );
1519 assert( up != NULL );
1520 assert( downvalid != NULL );
1521 assert( upvalid != NULL );
1522
1523 SCIPdebugMessage("calling strongbranching on variable %d (%d iterations)\n", col_index, itlim);
1524
1525 /* We work on the scaled problem. */
1526 const ColIndex col(col_index);
1527 const Fractional lb = lpi->scaled_lp->variable_lower_bounds()[col];
1528 const Fractional ub = lpi->scaled_lp->variable_upper_bounds()[col];
1529 const double value = psol * lpi->scaler->VariableScalingFactor(col);
1530
1531 /* Configure solver. */
1532
1533 /* @todo Use the iteration limit once glop supports incrementality. */
1534 int num_iterations = 0;
1535 lpi->parameters->set_use_dual_simplex(true);
1536 lpi->solver->SetParameters(*(lpi->parameters));
1537 const Fractional eps = lpi->parameters->primal_feasibility_tolerance();
1538
1539 std::unique_ptr<TimeLimit> time_limit = TimeLimit::FromParameters(*lpi->parameters);
1540
1541 /* Down branch. */
1542 const Fractional newub = EPSCEIL(value - 1.0, eps);
1543 if ( newub >= lb - 0.5 )
1544 {
1545 lpi->scaled_lp->SetVariableBounds(col, lb, newub);
1546
1547 if ( lpi->solver->Solve(*(lpi->scaled_lp), time_limit.get()).ok() )
1548 {
1549 num_iterations += (int) lpi->solver->GetNumberOfIterations();
1550 *down = lpi->solver->GetObjectiveValue();
1551 *downvalid = IsDualBoundValid(lpi->solver->GetProblemStatus()) ? TRUE : FALSE;
1552
1553 SCIPdebugMessage("down: itlim=%d col=%d [%f,%f] obj=%f status=%d iter=%ld.\n", itlim, col_index, lb, EPSCEIL(value - 1.0, eps),
1554 lpi->solver->GetObjectiveValue(), (int) lpi->solver->GetProblemStatus(), lpi->solver->GetNumberOfIterations());
1555 }
1556 else
1557 {
1558 SCIPerrorMessage("error during solve");
1559 *down = 0.0;
1560 *downvalid = FALSE;
1561 }
1562 }
1563 else
1564 {
1565 if ( lpi->linear_program->IsMaximizationProblem() )
1566 *down = lpi->parameters->objective_lower_limit();
1567 else
1568 *down = lpi->parameters->objective_upper_limit();
1569 *downvalid = TRUE;
1570 }
1571
1572 /* Up branch. */
1573 const Fractional newlb = EPSFLOOR(value + 1.0, eps);
1574 if ( newlb <= ub + 0.5 )
1575 {
1576 lpi->scaled_lp->SetVariableBounds(col, newlb, ub);
1577
1578 if ( lpi->solver->Solve(*(lpi->scaled_lp), time_limit.get()).ok() )
1579 {
1580 num_iterations += (int) lpi->solver->GetNumberOfIterations();
1581 *up = lpi->solver->GetObjectiveValue();
1582 *upvalid = IsDualBoundValid(lpi->solver->GetProblemStatus()) ? TRUE : FALSE;
1583
1584 SCIPdebugMessage("up: itlim=%d col=%d [%f,%f] obj=%f status=%d iter=%ld.\n", itlim, col_index, EPSFLOOR(value + 1.0, eps), ub,
1585 lpi->solver->GetObjectiveValue(), (int) lpi->solver->GetProblemStatus(), lpi->solver->GetNumberOfIterations());
1586 }
1587 else
1588 {
1589 SCIPerrorMessage("error during solve");
1590 *up = 0.0;
1591 *upvalid = FALSE;
1592 }
1593 }
1594 else
1595 {
1596 if (lpi->linear_program->IsMaximizationProblem())
1597 *up = lpi->parameters->objective_lower_limit();
1598 else
1599 *up = lpi->parameters->objective_upper_limit();
1600 *upvalid = TRUE;
1601 }
1602
1603 /* Restore bound. */
1604 lpi->scaled_lp->SetVariableBounds(col, lb, ub);
1605 if ( iter != NULL )
1606 *iter = num_iterations;
1607
1608 return SCIP_OKAY;
1609 }
1610
1611 /** performs strong branching iterations on one @b fractional candidate */
SCIPlpiStrongbranchFrac(SCIP_LPI * lpi,int col_index,SCIP_Real psol,int itlim,SCIP_Real * down,SCIP_Real * up,SCIP_Bool * downvalid,SCIP_Bool * upvalid,int * iter)1612 SCIP_RETCODE SCIPlpiStrongbranchFrac(
1613 SCIP_LPI* lpi, /**< LP interface structure */
1614 int col_index, /**< column to apply strong branching on */
1615 SCIP_Real psol, /**< fractional current primal solution value of column */
1616 int itlim, /**< iteration limit for strong branchings */
1617 SCIP_Real* down, /**< stores dual bound after branching column down */
1618 SCIP_Real* up, /**< stores dual bound after branching column up */
1619 SCIP_Bool* downvalid, /**< stores whether the returned down value is a valid dual bound;
1620 * otherwise, it can only be used as an estimate value */
1621 SCIP_Bool* upvalid, /**< stores whether the returned up value is a valid dual bound;
1622 * otherwise, it can only be used as an estimate value */
1623 int* iter /**< stores total number of strong branching iterations, or -1; may be NULL */
1624 )
1625 {
1626 assert( lpi != NULL );
1627 assert( lpi->scaled_lp != NULL );
1628 assert( down != NULL );
1629 assert( up != NULL );
1630 assert( downvalid != NULL );
1631 assert( upvalid != NULL );
1632
1633 SCIPdebugMessage("calling strongbranching on fractional variable %d (%d iterations)\n", col_index, itlim);
1634
1635 SCIP_CALL( strongbranch(lpi, col_index, psol, itlim, down, up, downvalid, upvalid, iter) );
1636
1637 return SCIP_OKAY;
1638 }
1639
1640 /** performs strong branching iterations on given @b fractional candidates */
SCIPlpiStrongbranchesFrac(SCIP_LPI * lpi,int * cols,int ncols,SCIP_Real * psols,int itlim,SCIP_Real * down,SCIP_Real * up,SCIP_Bool * downvalid,SCIP_Bool * upvalid,int * iter)1641 SCIP_RETCODE SCIPlpiStrongbranchesFrac(
1642 SCIP_LPI* lpi, /**< LP interface structure */
1643 int* cols, /**< columns to apply strong branching on */
1644 int ncols, /**< number of columns */
1645 SCIP_Real* psols, /**< fractional current primal solution values of columns */
1646 int itlim, /**< iteration limit for strong branchings */
1647 SCIP_Real* down, /**< stores dual bounds after branching columns down */
1648 SCIP_Real* up, /**< stores dual bounds after branching columns up */
1649 SCIP_Bool* downvalid, /**< stores whether the returned down values are valid dual bounds;
1650 * otherwise, they can only be used as an estimate values */
1651 SCIP_Bool* upvalid, /**< stores whether the returned up values are a valid dual bounds;
1652 * otherwise, they can only be used as an estimate values */
1653 int* iter /**< stores total number of strong branching iterations, or -1; may be NULL */
1654 )
1655 {
1656 assert( lpi != NULL );
1657 assert( lpi->linear_program != NULL );
1658 assert( cols != NULL );
1659 assert( psols != NULL );
1660 assert( down != NULL) ;
1661 assert( up != NULL );
1662 assert( downvalid != NULL );
1663 assert( upvalid != NULL );
1664
1665 SCIPerrorMessage("SCIPlpiStrongbranchesFrac - not implemented.\n");
1666
1667 return SCIP_NOTIMPLEMENTED;
1668 }
1669
1670 /** performs strong branching iterations on one candidate with @b integral value */
SCIPlpiStrongbranchInt(SCIP_LPI * lpi,int col,SCIP_Real psol,int itlim,SCIP_Real * down,SCIP_Real * up,SCIP_Bool * downvalid,SCIP_Bool * upvalid,int * iter)1671 SCIP_RETCODE SCIPlpiStrongbranchInt(
1672 SCIP_LPI* lpi, /**< LP interface structure */
1673 int col, /**< column to apply strong branching on */
1674 SCIP_Real psol, /**< current integral primal solution value of column */
1675 int itlim, /**< iteration limit for strong branchings */
1676 SCIP_Real* down, /**< stores dual bound after branching column down */
1677 SCIP_Real* up, /**< stores dual bound after branching column up */
1678 SCIP_Bool* downvalid, /**< stores whether the returned down value is a valid dual bound;
1679 * otherwise, it can only be used as an estimate value */
1680 SCIP_Bool* upvalid, /**< stores whether the returned up value is a valid dual bound;
1681 * otherwise, it can only be used as an estimate value */
1682 int* iter /**< stores total number of strong branching iterations, or -1; may be NULL */
1683 )
1684 {
1685 assert( lpi != NULL );
1686 assert( lpi->linear_program != NULL );
1687 assert( down != NULL );
1688 assert( up != NULL );
1689 assert( downvalid != NULL );
1690 assert( upvalid != NULL );
1691
1692 SCIP_CALL( strongbranch(lpi, col, psol, itlim, down, up, downvalid, upvalid, iter) );
1693
1694 return SCIP_OKAY;
1695 }
1696
1697 /** performs strong branching iterations on given candidates with @b integral values */
SCIPlpiStrongbranchesInt(SCIP_LPI * lpi,int * cols,int ncols,SCIP_Real * psols,int itlim,SCIP_Real * down,SCIP_Real * up,SCIP_Bool * downvalid,SCIP_Bool * upvalid,int * iter)1698 SCIP_RETCODE SCIPlpiStrongbranchesInt(
1699 SCIP_LPI* lpi, /**< LP interface structure */
1700 int* cols, /**< columns to apply strong branching on */
1701 int ncols, /**< number of columns */
1702 SCIP_Real* psols, /**< current integral primal solution values of columns */
1703 int itlim, /**< iteration limit for strong branchings */
1704 SCIP_Real* down, /**< stores dual bounds after branching columns down */
1705 SCIP_Real* up, /**< stores dual bounds after branching columns up */
1706 SCIP_Bool* downvalid, /**< stores whether the returned down values are valid dual bounds;
1707 * otherwise, they can only be used as an estimate values */
1708 SCIP_Bool* upvalid, /**< stores whether the returned up values are a valid dual bounds;
1709 * otherwise, they can only be used as an estimate values */
1710 int* iter /**< stores total number of strong branching iterations, or -1; may be NULL */
1711 )
1712 {
1713 assert( lpi != NULL );
1714 assert( lpi->linear_program != NULL );
1715 assert( cols != NULL );
1716 assert( psols != NULL );
1717 assert( down != NULL) ;
1718 assert( up != NULL );
1719 assert( downvalid != NULL );
1720 assert( upvalid != NULL );
1721
1722 SCIPerrorMessage("SCIPlpiStrongbranchesInt - not implemented.\n");
1723
1724 return SCIP_NOTIMPLEMENTED;
1725 }
1726 /**@} */
1727
1728
1729
1730
1731 /*
1732 * Solution Information Methods
1733 */
1734
1735 /**@name Solution Information Methods */
1736 /**@{ */
1737
1738 /** returns whether a solve method was called after the last modification of the LP */
SCIPlpiWasSolved(SCIP_LPI * lpi)1739 SCIP_Bool SCIPlpiWasSolved(
1740 SCIP_LPI* lpi /**< LP interface structure */
1741 )
1742 {
1743 assert( lpi != NULL );
1744
1745 /* @todo Track this to avoid uneeded resolving. */
1746 return ( ! lpi->lp_modified_since_last_solve );
1747 }
1748
1749 /** gets information about primal and dual feasibility of the current LP solution
1750 *
1751 * The feasibility information is with respect to the last solving call and it is only relevant if SCIPlpiWasSolved()
1752 * returns true. If the LP is changed, this information might be invalidated.
1753 *
1754 * Note that @a primalfeasible and @a dualfeasible should only return true if the solver has proved the respective LP to
1755 * be feasible. Thus, the return values should be equal to the values of SCIPlpiIsPrimalFeasible() and
1756 * SCIPlpiIsDualFeasible(), respectively. Note that if feasibility cannot be proved, they should return false (even if
1757 * the problem might actually be feasible).
1758 */
SCIPlpiGetSolFeasibility(SCIP_LPI * lpi,SCIP_Bool * primalfeasible,SCIP_Bool * dualfeasible)1759 SCIP_RETCODE SCIPlpiGetSolFeasibility(
1760 SCIP_LPI* lpi, /**< LP interface structure */
1761 SCIP_Bool* primalfeasible, /**< pointer to store primal feasibility status */
1762 SCIP_Bool* dualfeasible /**< pointer to store dual feasibility status */
1763 )
1764 {
1765 assert( lpi != NULL );
1766 assert( lpi->solver != NULL );
1767 assert( primalfeasible != NULL );
1768 assert( dualfeasible != NULL );
1769
1770 const ProblemStatus status = lpi->solver->GetProblemStatus();
1771 *primalfeasible = (status == ProblemStatus::OPTIMAL || status == ProblemStatus::PRIMAL_FEASIBLE);
1772 *dualfeasible = (status == ProblemStatus::OPTIMAL || status == ProblemStatus::DUAL_FEASIBLE);
1773
1774 SCIPdebugMessage("SCIPlpiGetSolFeasibility primal:%d dual:%d\n", *primalfeasible, *dualfeasible);
1775
1776 return SCIP_OKAY;
1777 }
1778
1779 /** returns TRUE iff LP is proven to have a primal unbounded ray (but not necessary a primal feasible point);
1780 * this does not necessarily mean, that the solver knows and can return the primal ray
1781 */
SCIPlpiExistsPrimalRay(SCIP_LPI * lpi)1782 SCIP_Bool SCIPlpiExistsPrimalRay(
1783 SCIP_LPI* lpi /**< LP interface structure */
1784 )
1785 {
1786 assert( lpi != NULL );
1787 assert( lpi->solver != NULL );
1788
1789 return lpi->solver->GetProblemStatus() == ProblemStatus::PRIMAL_UNBOUNDED;
1790 }
1791
1792 /** returns TRUE iff LP is proven to have a primal unbounded ray (but not necessary a primal feasible point),
1793 * and the solver knows and can return the primal ray
1794 */
SCIPlpiHasPrimalRay(SCIP_LPI * lpi)1795 SCIP_Bool SCIPlpiHasPrimalRay(
1796 SCIP_LPI* lpi /**< LP interface structure */
1797 )
1798 {
1799 assert( lpi != NULL );
1800 assert( lpi->solver != NULL );
1801
1802 return lpi->solver->GetProblemStatus() == ProblemStatus::PRIMAL_UNBOUNDED;
1803 }
1804
1805 /** returns TRUE iff LP is proven to be primal unbounded */
SCIPlpiIsPrimalUnbounded(SCIP_LPI * lpi)1806 SCIP_Bool SCIPlpiIsPrimalUnbounded(
1807 SCIP_LPI* lpi /**< LP interface structure */
1808 )
1809 {
1810 assert( lpi != NULL );
1811 assert( lpi->solver != NULL );
1812
1813 return lpi->solver->GetProblemStatus() == ProblemStatus::PRIMAL_UNBOUNDED;
1814 }
1815
1816 /** returns TRUE iff LP is proven to be primal infeasible */
SCIPlpiIsPrimalInfeasible(SCIP_LPI * lpi)1817 SCIP_Bool SCIPlpiIsPrimalInfeasible(
1818 SCIP_LPI* lpi /**< LP interface structure */
1819 )
1820 {
1821 assert( lpi != NULL );
1822 assert( lpi->solver != NULL );
1823
1824 const ProblemStatus status = lpi->solver->GetProblemStatus();
1825
1826 return status == ProblemStatus::DUAL_UNBOUNDED || status == ProblemStatus::PRIMAL_INFEASIBLE;
1827 }
1828
1829 /** returns TRUE iff LP is proven to be primal feasible */
SCIPlpiIsPrimalFeasible(SCIP_LPI * lpi)1830 SCIP_Bool SCIPlpiIsPrimalFeasible(
1831 SCIP_LPI* lpi /**< LP interface structure */
1832 )
1833 {
1834 assert( lpi != NULL );
1835 assert( lpi->solver != NULL );
1836
1837 const ProblemStatus status = lpi->solver->GetProblemStatus();
1838
1839 return status == ProblemStatus::PRIMAL_FEASIBLE || status == ProblemStatus::OPTIMAL;
1840 }
1841
1842 /** returns TRUE iff LP is proven to have a dual unbounded ray (but not necessary a dual feasible point);
1843 * this does not necessarily mean, that the solver knows and can return the dual ray
1844 */
SCIPlpiExistsDualRay(SCIP_LPI * lpi)1845 SCIP_Bool SCIPlpiExistsDualRay(
1846 SCIP_LPI* lpi /**< LP interface structure */
1847 )
1848 {
1849 assert( lpi != NULL );
1850 assert( lpi->solver != NULL );
1851
1852 const ProblemStatus status = lpi->solver->GetProblemStatus();
1853
1854 return status == ProblemStatus::DUAL_UNBOUNDED;
1855 }
1856
1857 /** returns TRUE iff LP is proven to have a dual unbounded ray (but not necessary a dual feasible point),
1858 * and the solver knows and can return the dual ray
1859 */
SCIPlpiHasDualRay(SCIP_LPI * lpi)1860 SCIP_Bool SCIPlpiHasDualRay(
1861 SCIP_LPI* lpi /**< LP interface structure */
1862 )
1863 {
1864 assert( lpi != NULL );
1865 assert( lpi->solver != NULL );
1866
1867 const ProblemStatus status = lpi->solver->GetProblemStatus();
1868
1869 return status == ProblemStatus::DUAL_UNBOUNDED;
1870 }
1871
1872 /** returns TRUE iff LP is proven to be dual unbounded */
SCIPlpiIsDualUnbounded(SCIP_LPI * lpi)1873 SCIP_Bool SCIPlpiIsDualUnbounded(
1874 SCIP_LPI* lpi /**< LP interface structure */
1875 )
1876 {
1877 assert( lpi != NULL );
1878 assert( lpi->solver != NULL );
1879
1880 const ProblemStatus status = lpi->solver->GetProblemStatus();
1881 return status == ProblemStatus::DUAL_UNBOUNDED;
1882 }
1883
1884 /** returns TRUE iff LP is proven to be dual infeasible */
SCIPlpiIsDualInfeasible(SCIP_LPI * lpi)1885 SCIP_Bool SCIPlpiIsDualInfeasible(
1886 SCIP_LPI* lpi /**< LP interface structure */
1887 )
1888 {
1889 assert( lpi != NULL );
1890 assert( lpi->solver != NULL );
1891
1892 const ProblemStatus status = lpi->solver->GetProblemStatus();
1893 return status == ProblemStatus::PRIMAL_UNBOUNDED || status == ProblemStatus::DUAL_INFEASIBLE;
1894 }
1895
1896 /** returns TRUE iff LP is proven to be dual feasible */
SCIPlpiIsDualFeasible(SCIP_LPI * lpi)1897 SCIP_Bool SCIPlpiIsDualFeasible(
1898 SCIP_LPI* lpi /**< LP interface structure */
1899 )
1900 {
1901 assert( lpi != NULL );
1902 assert( lpi->solver != NULL );
1903
1904 const ProblemStatus status = lpi->solver->GetProblemStatus();
1905
1906 return status == ProblemStatus::DUAL_FEASIBLE || status == ProblemStatus::OPTIMAL;
1907 }
1908
1909 /** returns TRUE iff LP was solved to optimality */
SCIPlpiIsOptimal(SCIP_LPI * lpi)1910 SCIP_Bool SCIPlpiIsOptimal(
1911 SCIP_LPI* lpi /**< LP interface structure */
1912 )
1913 {
1914 assert( lpi != NULL );
1915 assert( lpi->solver != NULL );
1916
1917 return lpi->solver->GetProblemStatus() == ProblemStatus::OPTIMAL;
1918 }
1919
1920 /** returns TRUE iff current LP solution is stable
1921 *
1922 * This function should return true if the solution is reliable, i.e., feasible and optimal (or proven
1923 * infeasible/unbounded) with respect to the original problem. The optimality status might be with respect to a scaled
1924 * version of the problem, but the solution might not be feasible to the unscaled original problem; in this case,
1925 * SCIPlpiIsStable() should return false.
1926 */
SCIPlpiIsStable(SCIP_LPI * lpi)1927 SCIP_Bool SCIPlpiIsStable(
1928 SCIP_LPI* lpi /**< LP interface structure */
1929 )
1930 {
1931 assert( lpi != NULL );
1932 assert( lpi->solver != NULL );
1933
1934 /* For correctness, we need to report "unstable" if Glop was not able to prove optimality because of numerical
1935 * issues. Currently Glop still reports primal/dual feasible if at the end, one status is within the tolerance but not
1936 * the other. */
1937 const ProblemStatus status = lpi->solver->GetProblemStatus();
1938 if ( (status == ProblemStatus::PRIMAL_FEASIBLE || status == ProblemStatus::DUAL_FEASIBLE) &&
1939 ! SCIPlpiIsObjlimExc(lpi) && ! SCIPlpiIsIterlimExc(lpi) && ! SCIPlpiIsTimelimExc(lpi))
1940 {
1941 SCIPdebugMessage("OPTIMAL not reached and no limit: unstable.\n");
1942 return FALSE;
1943 }
1944
1945 if ( status == ProblemStatus::ABNORMAL || status == ProblemStatus::INVALID_PROBLEM || status == ProblemStatus::IMPRECISE )
1946 return FALSE;
1947
1948 /* If we have a regular basis and the condition limit is set, we compute (an upper bound on) the condition number of
1949 * the basis; everything above the specified threshold is then counted as instable. */
1950 if ( lpi->checkcondition && (SCIPlpiIsOptimal(lpi) || SCIPlpiIsObjlimExc(lpi)) )
1951 {
1952 SCIP_RETCODE retcode;
1953 SCIP_Real kappa;
1954
1955 retcode = SCIPlpiGetRealSolQuality(lpi, SCIP_LPSOLQUALITY_ESTIMCONDITION, &kappa);
1956 if ( retcode != SCIP_OKAY )
1957 SCIPABORT();
1958 assert( kappa != SCIP_INVALID ); /*lint !e777*/
1959
1960 if ( kappa > lpi->conditionlimit )
1961 return FALSE;
1962 }
1963
1964 return TRUE;
1965 }
1966
1967 /** returns TRUE iff the objective limit was reached */
SCIPlpiIsObjlimExc(SCIP_LPI * lpi)1968 SCIP_Bool SCIPlpiIsObjlimExc(
1969 SCIP_LPI* lpi /**< LP interface structure */
1970 )
1971 {
1972 assert( lpi != NULL );
1973 assert( lpi->solver != NULL );
1974
1975 return lpi->solver->objective_limit_reached();
1976 }
1977
1978 /** returns TRUE iff the iteration limit was reached */
SCIPlpiIsIterlimExc(SCIP_LPI * lpi)1979 SCIP_Bool SCIPlpiIsIterlimExc(
1980 SCIP_LPI* lpi /**< LP interface structure */
1981 )
1982 {
1983 assert( lpi != NULL );
1984 assert( lpi->solver != NULL );
1985 assert( lpi->niterations >= (int) lpi->solver->GetNumberOfIterations() );
1986
1987 int maxiter = (int) lpi->parameters->max_number_of_iterations();
1988 return maxiter >= 0 && lpi->niterations >= maxiter;
1989 }
1990
1991 /** returns TRUE iff the time limit was reached */
SCIPlpiIsTimelimExc(SCIP_LPI * lpi)1992 SCIP_Bool SCIPlpiIsTimelimExc(
1993 SCIP_LPI* lpi /**< LP interface structure */
1994 )
1995 {
1996 assert( lpi != NULL );
1997 assert( lpi->solver != NULL );
1998
1999 return lpi->lp_time_limit_was_reached;
2000 }
2001
2002 /** returns the internal solution status of the solver */
SCIPlpiGetInternalStatus(SCIP_LPI * lpi)2003 int SCIPlpiGetInternalStatus(
2004 SCIP_LPI* lpi /**< LP interface structure */
2005 )
2006 {
2007 assert( lpi != NULL );
2008 assert( lpi->solver != NULL );
2009
2010 return static_cast<int>(lpi->solver->GetProblemStatus());
2011 }
2012
2013 /** tries to reset the internal status of the LP solver in order to ignore an instability of the last solving call */
SCIPlpiIgnoreInstability(SCIP_LPI * lpi,SCIP_Bool * success)2014 SCIP_RETCODE SCIPlpiIgnoreInstability(
2015 SCIP_LPI* lpi, /**< LP interface structure */
2016 SCIP_Bool* success /**< pointer to store, whether the instability could be ignored */
2017 )
2018 {
2019 assert( lpi != NULL );
2020 assert( lpi->solver != NULL );
2021 assert( success != NULL );
2022
2023 *success = FALSE;
2024
2025 return SCIP_OKAY;
2026 }
2027
2028 /** gets objective value of solution */
SCIPlpiGetObjval(SCIP_LPI * lpi,SCIP_Real * objval)2029 SCIP_RETCODE SCIPlpiGetObjval(
2030 SCIP_LPI* lpi, /**< LP interface structure */
2031 SCIP_Real* objval /**< stores the objective value */
2032 )
2033 {
2034 assert( lpi != NULL );
2035 assert( lpi->solver != NULL );
2036 assert( objval != NULL );
2037
2038 *objval = lpi->solver->GetObjectiveValue();
2039
2040 return SCIP_OKAY;
2041 }
2042
2043 /** gets primal and dual solution vectors for feasible LPs
2044 *
2045 * Before calling this function, the caller must ensure that the LP has been solved to optimality, i.e., that
2046 * SCIPlpiIsOptimal() returns true.
2047 */
SCIPlpiGetSol(SCIP_LPI * lpi,SCIP_Real * objval,SCIP_Real * primsol,SCIP_Real * dualsol,SCIP_Real * activity,SCIP_Real * redcost)2048 SCIP_RETCODE SCIPlpiGetSol(
2049 SCIP_LPI* lpi, /**< LP interface structure */
2050 SCIP_Real* objval, /**< stores the objective value, may be NULL if not needed */
2051 SCIP_Real* primsol, /**< primal solution vector, may be NULL if not needed */
2052 SCIP_Real* dualsol, /**< dual solution vector, may be NULL if not needed */
2053 SCIP_Real* activity, /**< row activity vector, may be NULL if not needed */
2054 SCIP_Real* redcost /**< reduced cost vector, may be NULL if not needed */
2055 )
2056 {
2057 assert( lpi != NULL );
2058 assert( lpi->solver != NULL );
2059
2060 SCIPdebugMessage("SCIPlpiGetSol\n");
2061 if ( objval != NULL )
2062 *objval = lpi->solver->GetObjectiveValue();
2063
2064 const ColIndex num_cols = lpi->linear_program->num_variables();
2065 for (ColIndex col(0); col < num_cols; ++col)
2066 {
2067 int i = col.value();
2068
2069 if ( primsol != NULL )
2070 primsol[i] = lpi->scaler->UnscaleVariableValue(col, lpi->solver->GetVariableValue(col));
2071
2072 if ( redcost != NULL )
2073 redcost[i] = lpi->scaler->UnscaleReducedCost(col, lpi->solver->GetReducedCost(col));
2074 }
2075
2076 const RowIndex num_rows = lpi->linear_program->num_constraints();
2077 for (RowIndex row(0); row < num_rows; ++row)
2078 {
2079 int j = row.value();
2080
2081 if ( dualsol != NULL )
2082 dualsol[j] = lpi->scaler->UnscaleDualValue(row, lpi->solver->GetDualValue(row));
2083
2084 if ( activity != NULL )
2085 activity[j] = lpi->scaler->UnscaleConstraintActivity(row, lpi->solver->GetConstraintActivity(row));
2086 }
2087
2088 return SCIP_OKAY;
2089 }
2090
2091 /** gets primal ray for unbounded LPs */
SCIPlpiGetPrimalRay(SCIP_LPI * lpi,SCIP_Real * ray)2092 SCIP_RETCODE SCIPlpiGetPrimalRay(
2093 SCIP_LPI* lpi, /**< LP interface structure */
2094 SCIP_Real* ray /**< primal ray */
2095 )
2096 {
2097 assert( lpi != NULL );
2098 assert( lpi->solver != NULL );
2099 assert( ray != NULL );
2100
2101 SCIPdebugMessage("SCIPlpiGetPrimalRay\n");
2102
2103 const ColIndex num_cols = lpi->linear_program->num_variables();
2104 const DenseRow& primal_ray = lpi->solver->GetPrimalRay();
2105 for (ColIndex col(0); col < num_cols; ++col)
2106 ray[col.value()] = lpi->scaler->UnscaleVariableValue(col, primal_ray[col]);
2107
2108 return SCIP_OKAY;
2109 }
2110
2111 /** gets dual Farkas proof for infeasibility */
SCIPlpiGetDualfarkas(SCIP_LPI * lpi,SCIP_Real * dualfarkas)2112 SCIP_RETCODE SCIPlpiGetDualfarkas(
2113 SCIP_LPI* lpi, /**< LP interface structure */
2114 SCIP_Real* dualfarkas /**< dual Farkas row multipliers */
2115 )
2116 {
2117 assert( lpi != NULL );
2118 assert( lpi->solver != NULL );
2119 assert( dualfarkas != NULL );
2120
2121 SCIPdebugMessage("SCIPlpiGetDualfarkas\n");
2122
2123 const RowIndex num_rows = lpi->linear_program->num_constraints();
2124 const DenseColumn& dual_ray = lpi->solver->GetDualRay();
2125 for (RowIndex row(0); row < num_rows; ++row)
2126 dualfarkas[row.value()] = -lpi->scaler->UnscaleDualValue(row, dual_ray[row]); /* reverse sign */
2127
2128 return SCIP_OKAY;
2129 }
2130
2131 /** gets the number of LP iterations of the last solve call */
SCIPlpiGetIterations(SCIP_LPI * lpi,int * iterations)2132 SCIP_RETCODE SCIPlpiGetIterations(
2133 SCIP_LPI* lpi, /**< LP interface structure */
2134 int* iterations /**< pointer to store the number of iterations of the last solve call */
2135 )
2136 {
2137 assert( lpi != NULL );
2138 assert( lpi->solver != NULL );
2139 assert( iterations != NULL );
2140
2141 *iterations = (int) lpi->niterations;
2142
2143 return SCIP_OKAY;
2144 }
2145
2146 /** gets information about the quality of an LP solution
2147 *
2148 * Such information is usually only available, if also a (maybe not optimal) solution is available.
2149 * The LPI should return SCIP_INVALID for @p quality, if the requested quantity is not available.
2150 */
SCIPlpiGetRealSolQuality(SCIP_LPI * lpi,SCIP_LPSOLQUALITY qualityindicator,SCIP_Real * quality)2151 SCIP_RETCODE SCIPlpiGetRealSolQuality(
2152 SCIP_LPI* lpi, /**< LP interface structure */
2153 SCIP_LPSOLQUALITY qualityindicator, /**< indicates which quality should be returned */
2154 SCIP_Real* quality /**< pointer to store quality number */
2155 )
2156 {
2157 assert( lpi != NULL );
2158 assert( lpi->solver != NULL );
2159 assert( quality != NULL );
2160
2161 SCIPdebugMessage("Requesting solution quality: quality %d\n", qualityindicator);
2162
2163 switch ( qualityindicator )
2164 {
2165 case SCIP_LPSOLQUALITY_ESTIMCONDITION:
2166 *quality = lpi->solver->GetBasisFactorization().ComputeInfinityNormConditionNumber();
2167 break;
2168
2169 case SCIP_LPSOLQUALITY_EXACTCONDITION:
2170 *quality = lpi->solver->GetBasisFactorization().ComputeInfinityNormConditionNumberUpperBound();
2171 break;
2172
2173 default:
2174 SCIPerrorMessage("Solution quality %d unknown.\n", qualityindicator);
2175 return SCIP_INVALIDDATA;
2176 }
2177
2178 return SCIP_OKAY;
2179 }
2180
2181 /**@} */
2182
2183
2184
2185
2186 /*
2187 * LP Basis Methods
2188 */
2189
2190 /**@name LP Basis Methods */
2191 /**@{ */
2192
2193 /** convert Glop variable basis status to SCIP status */
2194 static
ConvertGlopVariableStatus(VariableStatus status,Fractional rc)2195 SCIP_BASESTAT ConvertGlopVariableStatus(
2196 VariableStatus status, /**< variable status */
2197 Fractional rc /**< reduced cost of variable */
2198 )
2199 {
2200 switch ( status )
2201 {
2202 case VariableStatus::BASIC:
2203 return SCIP_BASESTAT_BASIC;
2204 case VariableStatus::AT_UPPER_BOUND:
2205 return SCIP_BASESTAT_UPPER;
2206 case VariableStatus::AT_LOWER_BOUND:
2207 return SCIP_BASESTAT_LOWER;
2208 case VariableStatus::FREE:
2209 return SCIP_BASESTAT_ZERO;
2210 case VariableStatus::FIXED_VALUE:
2211 return rc > 0.0 ? SCIP_BASESTAT_LOWER : SCIP_BASESTAT_UPPER;
2212 default:
2213 SCIPerrorMessage("invalid Glop basis status.\n");
2214 abort();
2215 }
2216 }
2217
2218 /** convert Glop constraint basis status to SCIP status */
2219 static
ConvertGlopConstraintStatus(ConstraintStatus status,Fractional dual)2220 SCIP_BASESTAT ConvertGlopConstraintStatus(
2221 ConstraintStatus status, /**< constraint status */
2222 Fractional dual /**< dual variable value */
2223 )
2224 {
2225 switch ( status )
2226 {
2227 case ConstraintStatus::BASIC:
2228 return SCIP_BASESTAT_BASIC;
2229 case ConstraintStatus::AT_UPPER_BOUND:
2230 return SCIP_BASESTAT_UPPER;
2231 case ConstraintStatus::AT_LOWER_BOUND:
2232 return SCIP_BASESTAT_LOWER;
2233 case ConstraintStatus::FREE:
2234 return SCIP_BASESTAT_ZERO;
2235 case ConstraintStatus::FIXED_VALUE:
2236 return dual > 0.0 ? SCIP_BASESTAT_LOWER : SCIP_BASESTAT_UPPER;
2237 default:
2238 SCIPerrorMessage("invalid Glop basis status.\n");
2239 abort();
2240 }
2241 }
2242
2243 /** Convert SCIP variable status to Glop status */
2244 static
ConvertSCIPVariableStatus(int status)2245 VariableStatus ConvertSCIPVariableStatus(
2246 int status /**< SCIP variable status */
2247 )
2248 {
2249 switch ( status )
2250 {
2251 case SCIP_BASESTAT_BASIC:
2252 return VariableStatus::BASIC;
2253 case SCIP_BASESTAT_UPPER:
2254 return VariableStatus::AT_UPPER_BOUND;
2255 case SCIP_BASESTAT_LOWER:
2256 return VariableStatus::AT_LOWER_BOUND;
2257 case SCIP_BASESTAT_ZERO:
2258 return VariableStatus::FREE;
2259 default:
2260 SCIPerrorMessage("invalid SCIP basis status.\n");
2261 abort();
2262 }
2263 }
2264
2265 /** Convert a SCIP constraint status to its corresponding Glop slack VariableStatus.
2266 *
2267 * Note that we swap the upper/lower bounds.
2268 */
2269 static
ConvertSCIPConstraintStatusToSlackStatus(int status)2270 VariableStatus ConvertSCIPConstraintStatusToSlackStatus(
2271 int status /**< SCIP constraint status */
2272 )
2273 {
2274 switch ( status )
2275 {
2276 case SCIP_BASESTAT_BASIC:
2277 return VariableStatus::BASIC;
2278 case SCIP_BASESTAT_UPPER:
2279 return VariableStatus::AT_LOWER_BOUND;
2280 case SCIP_BASESTAT_LOWER:
2281 return VariableStatus::AT_UPPER_BOUND;
2282 case SCIP_BASESTAT_ZERO:
2283 return VariableStatus::FREE;
2284 default:
2285 SCIPerrorMessage("invalid SCIP basis status.\n");
2286 abort();
2287 }
2288 }
2289
2290 /** gets current basis status for columns and rows; arrays must be large enough to store the basis status */
SCIPlpiGetBase(SCIP_LPI * lpi,int * cstat,int * rstat)2291 SCIP_RETCODE SCIPlpiGetBase(
2292 SCIP_LPI* lpi, /**< LP interface structure */
2293 int* cstat, /**< array to store column basis status, or NULL */
2294 int* rstat /**< array to store row basis status, or NULL */
2295 )
2296 {
2297 SCIPdebugMessage("SCIPlpiGetBase\n");
2298
2299 assert( lpi->solver->GetProblemStatus() == ProblemStatus::OPTIMAL );
2300
2301 if ( cstat != NULL )
2302 {
2303 const ColIndex num_cols = lpi->linear_program->num_variables();
2304 for (ColIndex col(0); col < num_cols; ++col)
2305 {
2306 int i = col.value();
2307 cstat[i] = (int) ConvertGlopVariableStatus(lpi->solver->GetVariableStatus(col), lpi->solver->GetReducedCost(col));
2308 }
2309 }
2310
2311 if ( rstat != NULL )
2312 {
2313 const RowIndex num_rows = lpi->linear_program->num_constraints();
2314 for (RowIndex row(0); row < num_rows; ++row)
2315 {
2316 int i = row.value();
2317 rstat[i] = (int) ConvertGlopConstraintStatus(lpi->solver->GetConstraintStatus(row), lpi->solver->GetDualValue(row));
2318 }
2319 }
2320
2321 return SCIP_OKAY;
2322 }
2323
2324 /** sets current basis status for columns and rows */
SCIPlpiSetBase(SCIP_LPI * lpi,const int * cstat,const int * rstat)2325 SCIP_RETCODE SCIPlpiSetBase(
2326 SCIP_LPI* lpi, /**< LP interface structure */
2327 const int* cstat, /**< array with column basis status */
2328 const int* rstat /**< array with row basis status */
2329 )
2330 {
2331 assert( lpi != NULL );
2332 assert( lpi->linear_program != NULL );
2333 assert( lpi->solver != NULL );
2334
2335 const ColIndex num_cols = lpi->linear_program->num_variables();
2336 const RowIndex num_rows = lpi->linear_program->num_constraints();
2337
2338 assert( cstat != NULL || num_cols == 0 );
2339 assert( rstat != NULL || num_rows == 0 );
2340
2341 SCIPdebugMessage("SCIPlpiSetBase\n");
2342
2343 BasisState state;
2344 state.statuses.reserve(ColIndex(num_cols.value() + num_rows.value()));
2345
2346 for (ColIndex col(0); col < num_cols; ++col)
2347 state.statuses[col] = ConvertSCIPVariableStatus(cstat[col.value()]);
2348
2349 for (RowIndex row(0); row < num_rows; ++row)
2350 state.statuses[num_cols + RowToColIndex(row)] = ConvertSCIPConstraintStatusToSlackStatus(cstat[row.value()]);
2351
2352 lpi->solver->LoadStateForNextSolve(state);
2353
2354 return SCIP_OKAY;
2355 }
2356
2357 /** returns the indices of the basic columns and rows; basic column n gives value n, basic row m gives value -1-m */
SCIPlpiGetBasisInd(SCIP_LPI * lpi,int * bind)2358 SCIP_RETCODE SCIPlpiGetBasisInd(
2359 SCIP_LPI* lpi, /**< LP interface structure */
2360 int* bind /**< pointer to store basis indices ready to keep number of rows entries */
2361 )
2362 {
2363 assert( lpi != NULL );
2364 assert( lpi->linear_program != NULL );
2365 assert( lpi->solver != NULL );
2366 assert( bind != NULL );
2367
2368 SCIPdebugMessage("SCIPlpiGetBasisInd\n");
2369
2370 /* the order is important! */
2371 const ColIndex num_cols = lpi->linear_program->num_variables();
2372 const RowIndex num_rows = lpi->linear_program->num_constraints();
2373 for (RowIndex row(0); row < num_rows; ++row)
2374 {
2375 const ColIndex col = lpi->solver->GetBasis(row);
2376 if (col < num_cols)
2377 bind[row.value()] = col.value();
2378 else
2379 {
2380 assert( col < num_cols.value() + num_rows.value() );
2381 bind[row.value()] = -1 - (col - num_cols).value();
2382 }
2383 }
2384
2385 return SCIP_OKAY;
2386 }
2387
2388 /** get row of inverse basis matrix B^-1
2389 *
2390 * @note The LP interface defines slack variables to have coefficient +1. This means that if, internally, the LP solver
2391 * uses a -1 coefficient, then rows associated with slacks variables whose coefficient is -1, should be negated;
2392 * see also the explanation in lpi.h.
2393 */
SCIPlpiGetBInvRow(SCIP_LPI * lpi,int r,SCIP_Real * coef,int * inds,int * ninds)2394 SCIP_RETCODE SCIPlpiGetBInvRow(
2395 SCIP_LPI* lpi, /**< LP interface structure */
2396 int r, /**< row number */
2397 SCIP_Real* coef, /**< pointer to store the coefficients of the row */
2398 int* inds, /**< array to store the non-zero indices, or NULL */
2399 int* ninds /**< pointer to store the number of non-zero indices, or NULL
2400 * (-1: if we do not store sparsity information) */
2401 )
2402 {
2403 assert( lpi != NULL );
2404 assert( lpi->linear_program != NULL );
2405 assert( lpi->solver != NULL );
2406 assert( lpi->scaler != NULL );
2407 assert( coef != NULL );
2408
2409 lpi->solver->GetBasisFactorization().LeftSolveForUnitRow(ColIndex(r), lpi->tmp_row);
2410 lpi->scaler->UnscaleUnitRowLeftSolve(lpi->solver->GetBasis(RowIndex(r)), lpi->tmp_row);
2411
2412 const ColIndex size = lpi->tmp_row->values.size();
2413 assert( size.value() == lpi->linear_program->num_constraints() );
2414
2415 /* if we want a sparse vector */
2416 if ( ninds != NULL && inds != NULL )
2417 {
2418 *ninds = 0;
2419 /* Vectors in Glop might be stored in dense or sparse format depending on the values. If non_zeros are given, we
2420 * can directly loop over the non_zeros, otherwise we have to collect the nonzeros. */
2421 if ( ! lpi->tmp_row->non_zeros.empty() )
2422 {
2423 ScatteredRowIterator end = lpi->tmp_row->end();
2424 for (ScatteredRowIterator iter = lpi->tmp_row->begin(); iter != end; ++iter)
2425 {
2426 int idx = (*iter).column().value();
2427 assert( 0 <= idx && idx < lpi->linear_program->num_constraints() );
2428 coef[idx] = (*iter).coefficient();
2429 inds[(*ninds)++] = idx;
2430 }
2431 }
2432 else
2433 {
2434 /* use dense access to tmp_row */
2435 const Fractional eps = lpi->parameters->primal_feasibility_tolerance();
2436 for (ColIndex col(0); col < size; ++col)
2437 {
2438 SCIP_Real val = (*lpi->tmp_row)[col];
2439 if ( fabs(val) >= eps )
2440 {
2441 coef[col.value()] = val;
2442 inds[(*ninds)++] = col.value();
2443 }
2444 }
2445 }
2446 return SCIP_OKAY;
2447 }
2448
2449 /* dense version */
2450 for (ColIndex col(0); col < size; ++col)
2451 coef[col.value()] = (*lpi->tmp_row)[col];
2452
2453 if ( ninds != NULL )
2454 *ninds = -1;
2455
2456 return SCIP_OKAY;
2457 }
2458
2459 /** get column of inverse basis matrix B^-1
2460 *
2461 * @note The LP interface defines slack variables to have coefficient +1. This means that if, internally, the LP solver
2462 * uses a -1 coefficient, then rows associated with slacks variables whose coefficient is -1, should be negated;
2463 * see also the explanation in lpi.h.
2464 */
SCIPlpiGetBInvCol(SCIP_LPI * lpi,int c,SCIP_Real * coef,int * inds,int * ninds)2465 SCIP_RETCODE SCIPlpiGetBInvCol(
2466 SCIP_LPI* lpi, /**< LP interface structure */
2467 int c, /**< column number of B^-1; this is NOT the number of the column in the LP;
2468 * you have to call SCIPlpiGetBasisInd() to get the array which links the
2469 * B^-1 column numbers to the row and column numbers of the LP!
2470 * c must be between 0 and nrows-1, since the basis has the size
2471 * nrows * nrows */
2472 SCIP_Real* coef, /**< pointer to store the coefficients of the column */
2473 int* inds, /**< array to store the non-zero indices, or NULL */
2474 int* ninds /**< pointer to store the number of non-zero indices, or NULL
2475 * (-1: if we do not store sparsity information) */
2476 )
2477 {
2478 assert( lpi != NULL );
2479 assert( lpi->linear_program != NULL );
2480 assert( lpi->solver != NULL );
2481 assert( lpi->scaler != NULL );
2482 assert( coef != NULL );
2483
2484 /* we need to loop through the rows to extract the values for column c */
2485 const ColIndex col(c);
2486 const RowIndex num_rows = lpi->linear_program->num_constraints();
2487
2488 /* if we want a sparse vector */
2489 if ( ninds != NULL && inds != NULL )
2490 {
2491 const Fractional eps = lpi->parameters->primal_feasibility_tolerance();
2492
2493 *ninds = 0;
2494 for (int row = 0; row < num_rows; ++row)
2495 {
2496 lpi->solver->GetBasisFactorization().LeftSolveForUnitRow(ColIndex(row), lpi->tmp_row);
2497 lpi->scaler->UnscaleUnitRowLeftSolve(lpi->solver->GetBasis(RowIndex(row)), lpi->tmp_row);
2498
2499 SCIP_Real val = (*lpi->tmp_row)[col];
2500 if ( fabs(val) >= eps )
2501 {
2502 coef[row] = val;
2503 inds[(*ninds)++] = row;
2504 }
2505 }
2506 return SCIP_OKAY;
2507 }
2508
2509 /* dense version */
2510 for (int row = 0; row < num_rows; ++row)
2511 {
2512 lpi->solver->GetBasisFactorization().LeftSolveForUnitRow(ColIndex(row), lpi->tmp_row);
2513 lpi->scaler->UnscaleUnitRowLeftSolve(lpi->solver->GetBasis(RowIndex(row)), lpi->tmp_row);
2514 coef[row] = (*lpi->tmp_row)[col];
2515 }
2516
2517 if ( ninds != NULL )
2518 *ninds = -1;
2519
2520 return SCIP_OKAY;
2521 }
2522
2523 /** get row of inverse basis matrix times constraint matrix B^-1 * A
2524 *
2525 * @note The LP interface defines slack variables to have coefficient +1. This means that if, internally, the LP solver
2526 * uses a -1 coefficient, then rows associated with slacks variables whose coefficient is -1, should be negated;
2527 * see also the explanation in lpi.h.
2528 */
SCIPlpiGetBInvARow(SCIP_LPI * lpi,int r,const SCIP_Real * binvrow,SCIP_Real * coef,int * inds,int * ninds)2529 SCIP_RETCODE SCIPlpiGetBInvARow(
2530 SCIP_LPI* lpi, /**< LP interface structure */
2531 int r, /**< row number */
2532 const SCIP_Real* binvrow, /**< row in (A_B)^-1 from prior call to SCIPlpiGetBInvRow(), or NULL */
2533 SCIP_Real* coef, /**< vector to return coefficients of the row */
2534 int* inds, /**< array to store the non-zero indices, or NULL */
2535 int* ninds /**< pointer to store the number of non-zero indices, or NULL
2536 * (-1: if we do not store sparsity information) */
2537 )
2538 {
2539 assert( lpi != NULL );
2540 assert( lpi->linear_program != NULL );
2541 assert( lpi->solver != NULL );
2542 assert( lpi->scaler != NULL );
2543 assert( coef != NULL );
2544
2545 /* get row of basis inverse, loop through columns and muliply with matrix */
2546 lpi->solver->GetBasisFactorization().LeftSolveForUnitRow(ColIndex(r), lpi->tmp_row);
2547 lpi->scaler->UnscaleUnitRowLeftSolve(lpi->solver->GetBasis(RowIndex(r)), lpi->tmp_row);
2548
2549 const ColIndex num_cols = lpi->linear_program->num_variables();
2550
2551 /* if we want a sparse vector */
2552 if ( ninds != NULL && inds != NULL )
2553 {
2554 const Fractional eps = lpi->parameters->primal_feasibility_tolerance();
2555
2556 *ninds = 0;
2557 for (ColIndex col(0); col < num_cols; ++col)
2558 {
2559 SCIP_Real val = operations_research::glop::ScalarProduct(lpi->tmp_row->values, lpi->linear_program->GetSparseColumn(col));
2560 if ( fabs(val) >= eps )
2561 {
2562 coef[col.value()] = val;
2563 inds[(*ninds)++] = col.value();
2564 }
2565 }
2566 return SCIP_OKAY;
2567 }
2568
2569 /* dense version */
2570 for (ColIndex col(0); col < num_cols; ++col)
2571 coef[col.value()] = operations_research::glop::ScalarProduct(lpi->tmp_row->values, lpi->linear_program->GetSparseColumn(col));
2572
2573 if ( ninds != NULL )
2574 *ninds = -1;
2575
2576 return SCIP_OKAY;
2577 }
2578
2579 /** get column of inverse basis matrix times constraint matrix B^-1 * A
2580 *
2581 * @note The LP interface defines slack variables to have coefficient +1. This means that if, internally, the LP solver
2582 * uses a -1 coefficient, then rows associated with slacks variables whose coefficient is -1, should be negated;
2583 * see also the explanation in lpi.h.
2584 */
SCIPlpiGetBInvACol(SCIP_LPI * lpi,int c,SCIP_Real * coef,int * inds,int * ninds)2585 SCIP_RETCODE SCIPlpiGetBInvACol(
2586 SCIP_LPI* lpi, /**< LP interface structure */
2587 int c, /**< column number */
2588 SCIP_Real* coef, /**< vector to return coefficients of the column */
2589 int* inds, /**< array to store the non-zero indices, or NULL */
2590 int* ninds /**< pointer to store the number of non-zero indices, or NULL
2591 * (-1: if we do not store sparsity information) */
2592 )
2593 {
2594 assert( lpi != NULL );
2595 assert( lpi->linear_program != NULL );
2596 assert( lpi->solver != NULL );
2597 assert( lpi->scaler != NULL );
2598 assert( coef != NULL );
2599
2600 lpi->solver->GetBasisFactorization().RightSolveForProblemColumn(ColIndex(c), lpi->tmp_column);
2601 lpi->scaler->UnscaleColumnRightSolve(lpi->solver->GetBasisVector(), ColIndex(c), lpi->tmp_column);
2602
2603 const RowIndex num_rows = lpi->tmp_column->values.size();
2604
2605 /* if we want a sparse vector */
2606 if ( ninds != NULL && inds != NULL )
2607 {
2608 *ninds = 0;
2609 /* Vectors in Glop might be stored in dense or sparse format depending on the values. If non_zeros are given, we
2610 * can directly loop over the non_zeros, otherwise we have to collect the nonzeros. */
2611 if ( ! lpi->tmp_column->non_zeros.empty() )
2612 {
2613 ScatteredColumnIterator end = lpi->tmp_column->end();
2614 for (ScatteredColumnIterator iter = lpi->tmp_column->begin(); iter != end; ++iter)
2615 {
2616 int idx = (*iter).row().value();
2617 assert( 0 <= idx && idx < num_rows );
2618 coef[idx] = (*iter).coefficient();
2619 inds[(*ninds)++] = idx;
2620 }
2621 }
2622 else
2623 {
2624 /* use dense access to tmp_column */
2625 const Fractional eps = lpi->parameters->primal_feasibility_tolerance();
2626 for (RowIndex row(0); row < num_rows; ++row)
2627 {
2628 SCIP_Real val = (*lpi->tmp_column)[row];
2629 if ( fabs(val) > eps )
2630 {
2631 coef[row.value()] = val;
2632 inds[(*ninds)++] = row.value();
2633 }
2634 }
2635 }
2636 return SCIP_OKAY;
2637 }
2638
2639 /* dense version */
2640 for (RowIndex row(0); row < num_rows; ++row)
2641 coef[row.value()] = (*lpi->tmp_column)[row];
2642
2643 if ( ninds != NULL )
2644 *ninds = -1;
2645
2646 return SCIP_OKAY;
2647 }
2648
2649 /**@} */
2650
2651
2652
2653
2654 /*
2655 * LP State Methods
2656 */
2657
2658 /**@name LP State Methods */
2659 /**@{ */
2660
2661 /* SCIP_LPiState stores basis information and is implemented by the glop BasisState class.
2662 * However, because in type_lpi.h there is
2663 * typedef struct SCIP_LPiState SCIP_LPISTATE;
2664 * We cannot just use a typedef here and SCIP_LPiState needs to be a struct.
2665 */
2666 struct SCIP_LPiState : BasisState {};
2667
2668 /** stores LPi state (like basis information) into lpistate object */
SCIPlpiGetState(SCIP_LPI * lpi,BMS_BLKMEM * blkmem,SCIP_LPISTATE ** lpistate)2669 SCIP_RETCODE SCIPlpiGetState(
2670 SCIP_LPI* lpi, /**< LP interface structure */
2671 BMS_BLKMEM* blkmem, /**< block memory */
2672 SCIP_LPISTATE** lpistate /**< pointer to LPi state information (like basis information) */
2673 )
2674 {
2675 assert( lpi != NULL );
2676 assert( lpi->solver != NULL );
2677 assert( lpistate != NULL );
2678
2679 *lpistate = static_cast<SCIP_LPISTATE*>(new BasisState(lpi->solver->GetState()));
2680
2681 return SCIP_OKAY;
2682 }
2683
2684 /** loads LPi state (like basis information) into solver; note that the LP might have been extended with additional
2685 * columns and rows since the state was stored with SCIPlpiGetState()
2686 */
SCIPlpiSetState(SCIP_LPI * lpi,BMS_BLKMEM * blkmem,const SCIP_LPISTATE * lpistate)2687 SCIP_RETCODE SCIPlpiSetState(
2688 SCIP_LPI* lpi, /**< LP interface structure */
2689 BMS_BLKMEM* blkmem, /**< block memory */
2690 const SCIP_LPISTATE* lpistate /**< LPi state information (like basis information), or NULL */
2691 )
2692 {
2693 assert( lpi != NULL );
2694 assert( lpi->solver != NULL );
2695 assert( lpistate != NULL );
2696
2697 lpi->solver->LoadStateForNextSolve(*lpistate);
2698
2699 return SCIP_OKAY;
2700 }
2701
2702 /** clears current LPi state (like basis information) of the solver */
SCIPlpiClearState(SCIP_LPI * lpi)2703 SCIP_RETCODE SCIPlpiClearState(
2704 SCIP_LPI* lpi /**< LP interface structure */
2705 )
2706 {
2707 assert( lpi != NULL );
2708 assert( lpi->solver != NULL );
2709
2710 lpi->solver->ClearStateForNextSolve();
2711
2712 return SCIP_OKAY;
2713 }
2714
2715 /** frees LPi state information */
SCIPlpiFreeState(SCIP_LPI * lpi,BMS_BLKMEM * blkmem,SCIP_LPISTATE ** lpistate)2716 SCIP_RETCODE SCIPlpiFreeState(
2717 SCIP_LPI* lpi, /**< LP interface structure */
2718 BMS_BLKMEM* blkmem, /**< block memory */
2719 SCIP_LPISTATE** lpistate /**< pointer to LPi state information (like basis information) */
2720 )
2721 {
2722 assert( lpi != NULL );
2723 assert( lpi->solver != NULL );
2724 assert( lpistate != NULL );
2725
2726 delete *lpistate;
2727 *lpistate = NULL;
2728
2729 return SCIP_OKAY;
2730 }
2731
2732 /** checks, whether the given LP state contains simplex basis information */
SCIPlpiHasStateBasis(SCIP_LPI * lpi,SCIP_LPISTATE * lpistate)2733 SCIP_Bool SCIPlpiHasStateBasis(
2734 SCIP_LPI* lpi, /**< LP interface structure */
2735 SCIP_LPISTATE* lpistate /**< LP state information (like basis information), or NULL */
2736 )
2737 {
2738 assert( lpi != NULL );
2739 assert( lpi->solver != NULL );
2740
2741 return lpistate != NULL;
2742 }
2743
2744 /** reads LP state (like basis information from a file */
SCIPlpiReadState(SCIP_LPI * lpi,const char * fname)2745 SCIP_RETCODE SCIPlpiReadState(
2746 SCIP_LPI* lpi, /**< LP interface structure */
2747 const char* fname /**< file name */
2748 )
2749 {
2750 assert( lpi != NULL );
2751 assert( lpi->solver != NULL );
2752
2753 SCIPerrorMessage("SCIPlpiReadState - not implemented.\n");
2754
2755 return SCIP_NOTIMPLEMENTED;
2756 }
2757
2758 /** writes LPi state (i.e. basis information) to a file */
SCIPlpiWriteState(SCIP_LPI * lpi,const char * fname)2759 SCIP_RETCODE SCIPlpiWriteState(
2760 SCIP_LPI* lpi, /**< LP interface structure */
2761 const char* fname /**< file name */
2762 )
2763 {
2764 assert( lpi != NULL );
2765 assert( lpi->solver != NULL );
2766
2767 SCIPerrorMessage("SCIPlpiWriteState - not implemented.\n");
2768
2769 return SCIP_NOTIMPLEMENTED;
2770 }
2771
2772 /**@} */
2773
2774
2775
2776
2777 /*
2778 * LP Pricing Norms Methods
2779 */
2780
2781 /**@name LP Pricing Norms Methods */
2782 /**@{ */
2783
2784 /* SCIP_LPiNorms stores norm information so they are not recomputed from one state to the next. */
2785 /* @todo Implement this. */
2786 struct SCIP_LPiNorms {};
2787
2788 /** stores LPi pricing norms information
2789 *
2790 * @todo store primal norms as well?
2791 */
SCIPlpiGetNorms(SCIP_LPI * lpi,BMS_BLKMEM * blkmem,SCIP_LPINORMS ** lpinorms)2792 SCIP_RETCODE SCIPlpiGetNorms(
2793 SCIP_LPI* lpi, /**< LP interface structure */
2794 BMS_BLKMEM* blkmem, /**< block memory */
2795 SCIP_LPINORMS** lpinorms /**< pointer to LPi pricing norms information */
2796 )
2797 {
2798 assert( lpi != NULL );
2799 assert( blkmem != NULL );
2800 assert( lpi->solver != NULL );
2801 assert( lpinorms != NULL );
2802
2803 return SCIP_OKAY;
2804 }
2805
2806 /** loads LPi pricing norms into solver; note that the LP might have been extended with additional
2807 * columns and rows since the state was stored with SCIPlpiGetNorms()
2808 */
SCIPlpiSetNorms(SCIP_LPI * lpi,BMS_BLKMEM * blkmem,const SCIP_LPINORMS * lpinorms)2809 SCIP_RETCODE SCIPlpiSetNorms(
2810 SCIP_LPI* lpi, /**< LP interface structure */
2811 BMS_BLKMEM* blkmem, /**< block memory */
2812 const SCIP_LPINORMS* lpinorms /**< LPi pricing norms information, or NULL */
2813 )
2814 {
2815 assert( lpi != NULL );
2816 assert( blkmem != NULL );
2817 assert( lpi->solver != NULL );
2818 assert( lpinorms != NULL );
2819
2820 return SCIP_OKAY;
2821 }
2822
2823 /** frees pricing norms information */
SCIPlpiFreeNorms(SCIP_LPI * lpi,BMS_BLKMEM * blkmem,SCIP_LPINORMS ** lpinorms)2824 SCIP_RETCODE SCIPlpiFreeNorms(
2825 SCIP_LPI* lpi, /**< LP interface structure */
2826 BMS_BLKMEM* blkmem, /**< block memory */
2827 SCIP_LPINORMS** lpinorms /**< pointer to LPi pricing norms information, or NULL */
2828 )
2829 {
2830 assert( lpi != NULL );
2831 assert( blkmem != NULL );
2832 assert( lpi->solver != NULL );
2833 assert( lpinorms != NULL );
2834
2835 return SCIP_OKAY;
2836 }
2837
2838 /**@} */
2839
2840
2841
2842
2843 /*
2844 * Parameter Methods
2845 */
2846
2847 /**@name Parameter Methods */
2848 /**@{ */
2849
2850 /** gets integer parameter of LP */
SCIPlpiGetIntpar(SCIP_LPI * lpi,SCIP_LPPARAM type,int * ival)2851 SCIP_RETCODE SCIPlpiGetIntpar(
2852 SCIP_LPI* lpi, /**< LP interface structure */
2853 SCIP_LPPARAM type, /**< parameter number */
2854 int* ival /**< buffer to store the parameter value */
2855 )
2856 {
2857 assert( lpi != NULL );
2858 assert( lpi->parameters != NULL );
2859
2860 /* Not (yet) supported by Glop: SCIP_LPPAR_FASTMIP, SCIP_LPPAR_POLISHING, SCIP_LPPAR_REFACTOR */
2861 switch ( type )
2862 {
2863 case SCIP_LPPAR_FROMSCRATCH:
2864 *ival = (int) lpi->from_scratch;
2865 SCIPdebugMessage("SCIPlpiGetIntpar: SCIP_LPPAR_FROMSCRATCH = %d.\n", *ival);
2866 break;
2867 case SCIP_LPPAR_LPINFO:
2868 *ival = (int) lpi->lp_info;
2869 SCIPdebugMessage("SCIPlpiGetIntpar: SCIP_LPPAR_LPINFO = %d.\n", *ival);
2870 break;
2871 case SCIP_LPPAR_LPITLIM:
2872 *ival = (int) lpi->parameters->max_number_of_iterations();
2873 SCIPdebugMessage("SCIPlpiGetIntpar: SCIP_LPPAR_LPITLIM = %d.\n", *ival);
2874 break;
2875 case SCIP_LPPAR_PRESOLVING:
2876 *ival = lpi->parameters->use_preprocessing();
2877 SCIPdebugMessage("SCIPlpiGetIntpar: SCIP_LPPAR_PRESOLVING = %d.\n", *ival);
2878 break;
2879 case SCIP_LPPAR_PRICING:
2880 *ival = lpi->pricing;
2881 SCIPdebugMessage("SCIPlpiGetIntpar: SCIP_LPPAR_PRICING = %d.\n", *ival);
2882 break;
2883 #ifndef NOSCALING
2884 case SCIP_LPPAR_SCALING:
2885 *ival = lpi->parameters->use_scaling();
2886 SCIPdebugMessage("SCIPlpiGetIntpar: SCIP_LPPAR_SCALING = %d.\n", *ival);
2887 break;
2888 #endif
2889 case SCIP_LPPAR_THREADS:
2890 *ival = lpi->numthreads;
2891 SCIPdebugMessage("SCIPlpiGetIntpar: SCIP_LPPAR_THREADS = %d.\n", *ival);
2892 break;
2893 case SCIP_LPPAR_TIMING:
2894 *ival = lpi->timing;
2895 SCIPdebugMessage("SCIPlpiGetIntpar: SCIP_LPPAR_TIMING = %d.\n", *ival);
2896 break;
2897 case SCIP_LPPAR_RANDOMSEED:
2898 *ival = (int) lpi->parameters->random_seed();
2899 SCIPdebugMessage("SCIPlpiGetIntpar: SCIP_LPPAR_RANDOMSEED = %d.\n", *ival);
2900 break;
2901 default:
2902 return SCIP_PARAMETERUNKNOWN;
2903 }
2904
2905 return SCIP_OKAY;
2906 }
2907
2908 /** sets integer parameter of LP */
SCIPlpiSetIntpar(SCIP_LPI * lpi,SCIP_LPPARAM type,int ival)2909 SCIP_RETCODE SCIPlpiSetIntpar(
2910 SCIP_LPI* lpi, /**< LP interface structure */
2911 SCIP_LPPARAM type, /**< parameter number */
2912 int ival /**< parameter value */
2913 )
2914 {
2915 assert( lpi != NULL );
2916 assert( lpi->parameters != NULL );
2917
2918 switch ( type )
2919 {
2920 case SCIP_LPPAR_FROMSCRATCH:
2921 SCIPdebugMessage("SCIPlpiSetIntpar: SCIP_LPPAR_FROMSCRATCH -> %d.\n", ival);
2922 lpi->from_scratch = (bool) ival;
2923 break;
2924 case SCIP_LPPAR_LPINFO:
2925 SCIPdebugMessage("SCIPlpiSetIntpar: SCIP_LPPAR_LPINFO -> %d.\n", ival);
2926 if ( ival == 0 )
2927 {
2928 (void) google::SetVLOGLevel("*", google::GLOG_INFO);
2929 lpi->lp_info = false;
2930 }
2931 else
2932 {
2933 (void) google::SetVLOGLevel("*", google::GLOG_ERROR);
2934 lpi->lp_info = true;
2935 }
2936 break;
2937 case SCIP_LPPAR_LPITLIM:
2938 SCIPdebugMessage("SCIPlpiSetIntpar: SCIP_LPPAR_LPITLIM -> %d.\n", ival);
2939 lpi->parameters->set_max_number_of_iterations(ival);
2940 break;
2941 case SCIP_LPPAR_PRESOLVING:
2942 SCIPdebugMessage("SCIPlpiSetIntpar: SCIP_LPPAR_PRESOLVING -> %d.\n", ival);
2943 lpi->parameters->set_use_preprocessing(ival);
2944 break;
2945 case SCIP_LPPAR_PRICING:
2946 SCIPdebugMessage("SCIPlpiSetIntpar: SCIP_LPPAR_PRICING -> %d.\n", ival);
2947 lpi->pricing = (SCIP_Pricing)ival;
2948 switch ( lpi->pricing )
2949 {
2950 case SCIP_PRICING_LPIDEFAULT:
2951 case SCIP_PRICING_AUTO:
2952 case SCIP_PRICING_PARTIAL:
2953 case SCIP_PRICING_STEEP:
2954 case SCIP_PRICING_STEEPQSTART:
2955 lpi->parameters->set_feasibility_rule(operations_research::glop::GlopParameters_PricingRule_STEEPEST_EDGE);
2956 break;
2957 case SCIP_PRICING_FULL:
2958 /* Dantzig does not really fit, but use it anyway */
2959 lpi->parameters->set_feasibility_rule(operations_research::glop::GlopParameters_PricingRule_DANTZIG);
2960 break;
2961 case SCIP_PRICING_DEVEX:
2962 lpi->parameters->set_feasibility_rule(operations_research::glop::GlopParameters_PricingRule_DEVEX);
2963 break;
2964 default:
2965 return SCIP_PARAMETERUNKNOWN;
2966 }
2967 break;
2968 #ifndef NOSCALING
2969 case SCIP_LPPAR_SCALING:
2970 SCIPdebugMessage("SCIPlpiSetIntpar: SCIP_LPPAR_SCALING -> %d.\n", ival);
2971 lpi->parameters->set_use_scaling(ival);
2972 break;
2973 #endif
2974 case SCIP_LPPAR_THREADS:
2975 SCIPdebugMessage("SCIPlpiSetIntpar: SCIP_LPPAR_THREADS -> %d.\n", ival);
2976 assert( ival >= 0 );
2977 lpi->numthreads = ival;
2978 if ( ival == 0 )
2979 lpi->parameters->set_num_omp_threads(1);
2980 else
2981 lpi->parameters->set_num_omp_threads(ival);
2982 break;
2983 case SCIP_LPPAR_TIMING:
2984 SCIPdebugMessage("SCIPlpiSetIntpar: SCIP_LPPAR_TIMING -> %d.\n", ival);
2985 assert( 0 <= ival && ival <= 2 );
2986 lpi->timing = ival;
2987 if ( ival == 1 )
2988 absl::SetFlag(&FLAGS_time_limit_use_usertime, true);
2989 else
2990 absl::SetFlag(&FLAGS_time_limit_use_usertime, false);
2991 break;
2992 case SCIP_LPPAR_RANDOMSEED:
2993 SCIPdebugMessage("SCIPlpiSetIntpar: SCIP_LPPAR_RANDOMSEED -> %d.\n", ival);
2994 assert( ival >= 0 );
2995 lpi->parameters->set_random_seed(ival);
2996 break;
2997 default:
2998 return SCIP_PARAMETERUNKNOWN;
2999 }
3000
3001 return SCIP_OKAY;
3002 }
3003
3004 /** gets floating point parameter of LP */
SCIPlpiGetRealpar(SCIP_LPI * lpi,SCIP_LPPARAM type,SCIP_Real * dval)3005 SCIP_RETCODE SCIPlpiGetRealpar(
3006 SCIP_LPI* lpi, /**< LP interface structure */
3007 SCIP_LPPARAM type, /**< parameter number */
3008 SCIP_Real* dval /**< buffer to store the parameter value */
3009 )
3010 {
3011 assert( lpi != NULL );
3012 assert( lpi->parameters != NULL );
3013
3014 /* Not (yet) supported by Glop: SCIP_LPPAR_ROWREPSWITCH, SCIP_LPPAR_BARRIERCONVTOL */
3015 switch ( type )
3016 {
3017 case SCIP_LPPAR_FEASTOL:
3018 *dval = lpi->parameters->primal_feasibility_tolerance();
3019 SCIPdebugMessage("SCIPlpiGetRealpar: SCIP_LPPAR_FEASTOL = %g.\n", *dval);
3020 break;
3021 case SCIP_LPPAR_DUALFEASTOL:
3022 *dval = lpi->parameters->dual_feasibility_tolerance();
3023 SCIPdebugMessage("SCIPlpiGetRealpar: SCIP_LPPAR_DUALFEASTOL = %g.\n", *dval);
3024 break;
3025 case SCIP_LPPAR_OBJLIM:
3026 if (lpi->linear_program->IsMaximizationProblem())
3027 *dval = lpi->parameters->objective_lower_limit();
3028 else
3029 *dval = lpi->parameters->objective_upper_limit();
3030 SCIPdebugMessage("SCIPlpiGetRealpar: SCIP_LPPAR_OBJLIM = %f.\n", *dval);
3031 break;
3032 case SCIP_LPPAR_LPTILIM:
3033 if ( absl::GetFlag(FLAGS_time_limit_use_usertime) )
3034 *dval = lpi->parameters->max_time_in_seconds();
3035 else
3036 *dval = lpi->parameters->max_deterministic_time();
3037 SCIPdebugMessage("SCIPlpiGetRealpar: SCIP_LPPAR_LPTILIM = %f.\n", *dval);
3038 break;
3039 case SCIP_LPPAR_CONDITIONLIMIT:
3040 *dval = lpi->conditionlimit;
3041 break;
3042 #ifdef SCIP_DISABLED_CODE
3043 /* currently do not apply Markowitz parameter, since the default value does not seem suitable for Glop */
3044 case SCIP_LPPAR_MARKOWITZ:
3045 *dval = lpi->parameters->markowitz_singularity_threshold();
3046 SCIPdebugMessage("SCIPlpiGetRealpar: SCIP_LPPAR_MARKOWITZ = %f.\n", *dval);
3047 break;
3048 #endif
3049 default:
3050 return SCIP_PARAMETERUNKNOWN;
3051 }
3052
3053 return SCIP_OKAY;
3054 }
3055
3056 /** sets floating point parameter of LP */
SCIPlpiSetRealpar(SCIP_LPI * lpi,SCIP_LPPARAM type,SCIP_Real dval)3057 SCIP_RETCODE SCIPlpiSetRealpar(
3058 SCIP_LPI* lpi, /**< LP interface structure */
3059 SCIP_LPPARAM type, /**< parameter number */
3060 SCIP_Real dval /**< parameter value */
3061 )
3062 {
3063 assert( lpi != NULL );
3064 assert( lpi->parameters != NULL );
3065
3066 switch( type )
3067 {
3068 case SCIP_LPPAR_FEASTOL:
3069 SCIPdebugMessage("SCIPlpiSetRealpar: SCIP_LPPAR_FEASTOL -> %g.\n", dval);
3070 lpi->parameters->set_primal_feasibility_tolerance(dval);
3071 break;
3072 case SCIP_LPPAR_DUALFEASTOL:
3073 SCIPdebugMessage("SCIPlpiSetRealpar: SCIP_LPPAR_DUALFEASTOL -> %g.\n", dval);
3074 lpi->parameters->set_dual_feasibility_tolerance(dval);
3075 break;
3076 case SCIP_LPPAR_OBJLIM:
3077 SCIPdebugMessage("SCIPlpiSetRealpar: SCIP_LPPAR_OBJLIM -> %f.\n", dval);
3078 if (lpi->linear_program->IsMaximizationProblem())
3079 lpi->parameters->set_objective_lower_limit(dval);
3080 else
3081 lpi->parameters->set_objective_upper_limit(dval);
3082 break;
3083 case SCIP_LPPAR_LPTILIM:
3084 SCIPdebugMessage("SCIPlpiSetRealpar: SCIP_LPPAR_LPTILIM -> %f.\n", dval);
3085 if ( absl::GetFlag(FLAGS_time_limit_use_usertime) )
3086 lpi->parameters->set_max_time_in_seconds(dval);
3087 else
3088 lpi->parameters->set_max_deterministic_time(dval);
3089 break;
3090 case SCIP_LPPAR_CONDITIONLIMIT:
3091 lpi->conditionlimit = dval;
3092 lpi->checkcondition = (dval >= 0.0);
3093 break;
3094 #ifdef SCIP_DISABLED_CODE
3095 /* currently do not apply Markowitz parameter, since the default value does not seem suitable for Glop */
3096 case SCIP_LPPAR_MARKOWITZ:
3097 SCIPdebugMessage("SCIPlpiSetRealpar: SCIP_LPPAR_MARKOWITZ -> %f.\n", dval);
3098 lpi->parameters->set_markowitz_singularity_threshold(dval);
3099 break;
3100 #endif
3101 default:
3102 return SCIP_PARAMETERUNKNOWN;
3103 }
3104
3105 return SCIP_OKAY;
3106 }
3107
3108 /**@} */
3109
3110
3111
3112
3113 /*
3114 * Numerical Methods
3115 */
3116
3117 /**@name Numerical Methods */
3118 /**@{ */
3119
3120 /** returns value treated as infinity in the LP solver */
SCIPlpiInfinity(SCIP_LPI * lpi)3121 SCIP_Real SCIPlpiInfinity(
3122 SCIP_LPI* lpi /**< LP interface structure */
3123 )
3124 {
3125 assert( lpi != NULL );
3126 return std::numeric_limits<SCIP_Real>::infinity();
3127 }
3128
3129 /** checks if given value is treated as infinity in the LP solver */
SCIPlpiIsInfinity(SCIP_LPI * lpi,SCIP_Real val)3130 SCIP_Bool SCIPlpiIsInfinity(
3131 SCIP_LPI* lpi, /**< LP interface structure */
3132 SCIP_Real val /**< value to be checked for infinity */
3133 )
3134 {
3135 assert( lpi != NULL );
3136
3137 return val == std::numeric_limits<SCIP_Real>::infinity();
3138 }
3139
3140 /**@} */
3141
3142
3143
3144
3145 /*
3146 * File Interface Methods
3147 */
3148
3149 /**@name File Interface Methods */
3150 /**@{ */
3151
3152 /** reads LP from a file */
SCIPlpiReadLP(SCIP_LPI * lpi,const char * fname)3153 SCIP_RETCODE SCIPlpiReadLP(
3154 SCIP_LPI* lpi, /**< LP interface structure */
3155 const char* fname /**< file name */
3156 )
3157 {
3158 assert( lpi != NULL );
3159 assert( lpi->linear_program != NULL );
3160 assert( fname != NULL );
3161
3162 const std::string filespec(fname);
3163 MPModelProto proto;
3164 if ( ! ReadFileToProto(filespec, &proto) )
3165 {
3166 SCIPerrorMessage("Could not read <%s>\n", fname);
3167 return SCIP_READERROR;
3168 }
3169 lpi->linear_program->Clear();
3170 MPModelProtoToLinearProgram(proto, lpi->linear_program);
3171
3172 return SCIP_OKAY;
3173 }
3174
3175 /** writes LP to a file */
SCIPlpiWriteLP(SCIP_LPI * lpi,const char * fname)3176 SCIP_RETCODE SCIPlpiWriteLP(
3177 SCIP_LPI* lpi, /**< LP interface structure */
3178 const char* fname /**< file name */
3179 )
3180 {
3181 assert( lpi != NULL );
3182 assert( lpi->linear_program != NULL );
3183 assert( fname != NULL );
3184
3185 MPModelProto proto;
3186 LinearProgramToMPModelProto(*lpi->linear_program, &proto);
3187 const std::string filespec(fname);
3188 if ( ! WriteProtoToFile(filespec, proto, operations_research::ProtoWriteFormat::kProtoText, true) )
3189 {
3190 SCIPerrorMessage("Could not write <%s>\n", fname);
3191 return SCIP_READERROR;
3192 }
3193
3194 return SCIP_OKAY;
3195 }
3196
3197 /**@} */
3198