1 /* * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * */
2 /* */
3 /* This file is part of the class library */
4 /* SoPlex --- the Sequential object-oriented simPlex. */
5 /* */
6 /* Copyright (c) 1996-2023 Zuse Institute Berlin (ZIB) */
7 /* */
8 /* Licensed under the Apache License, Version 2.0 (the "License"); */
9 /* you may not use this file except in compliance with the License. */
10 /* You may obtain a copy of the License at */
11 /* */
12 /* http://www.apache.org/licenses/LICENSE-2.0 */
13 /* */
14 /* Unless required by applicable law or agreed to in writing, software */
15 /* distributed under the License is distributed on an "AS IS" BASIS, */
16 /* WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. */
17 /* See the License for the specific language governing permissions and */
18 /* limitations under the License. */
19 /* */
20 /* You should have received a copy of the Apache-2.0 license */
21 /* along with SoPlex; see the file LICENSE. If not email to soplex@zib.de. */
22 /* */
23 /* * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * */
24
25 /**@file soplex.h
26 * @brief Preconfigured SoPlex LP solver
27 */
28
29 #ifndef _SOPLEX_H_
30 #define _SOPLEX_H_
31
32 #include <string.h>
33
34 #include "soplex/spxgithash.h"
35 #include "soplex/spxdefines.h"
36 #include "soplex/basevectors.h"
37 #include "soplex/spxsolver.h"
38 #include "soplex/slufactor.h"
39 #include "soplex/slufactor_rational.h"
40
41 ///@todo try to move to cpp file by forward declaration
42 #include "soplex/spxsimplifier.h"
43 #include "soplex/spxmainsm.h"
44
45 #include "soplex/spxscaler.h"
46 #include "soplex/spxequilisc.h"
47 #include "soplex/spxleastsqsc.h"
48 #include "soplex/spxgeometsc.h"
49
50 #include "soplex/spxstarter.h"
51 #include "soplex/spxweightst.h"
52 #include "soplex/spxsumst.h"
53 #include "soplex/spxvectorst.h"
54
55 #include "soplex/spxpricer.h"
56 #include "soplex/spxautopr.h"
57 #include "soplex/spxdantzigpr.h"
58 #include "soplex/spxparmultpr.h"
59 #include "soplex/spxdevexpr.h"
60 #include "soplex/spxsteeppr.h"
61 #include "soplex/spxsteepexpr.h"
62 #include "soplex/spxhybridpr.h"
63
64 #include "soplex/spxratiotester.h"
65 #include "soplex/spxdefaultrt.h"
66 #include "soplex/spxharrisrt.h"
67 #include "soplex/spxfastrt.h"
68 #include "soplex/spxboundflippingrt.h"
69
70 #include "soplex/solbase.h"
71 #include "soplex/sol.h"
72
73 #include "soplex/spxlpbase.h"
74
75 #include "soplex/spxpapilo.h"
76
77 #ifdef SOPLEX_WITH_GMP
78 #include <gmp.h>
79 #endif
80
81 #ifdef SOPLEX_WITH_BOOST
82 #ifdef SOPLEX_WITH_MPFR
83 // For multiple precision
84 #include <boost/multiprecision/mpfr.hpp>
85 #endif
86 #ifdef SOPLEX_WITH_CPPMPF
87 #include <boost/multiprecision/cpp_dec_float.hpp>
88 #endif
89
90 // An alias for boost multiprecision
91 namespace mpf = boost::multiprecision;
92 #endif
93
94 #define SOPLEX_DEFAULT_RANDOM_SEED 0 // used to suppress output when the seed was not changed
95
96 ///@todo implement automatic rep switch, based on row/col dim
97 ///@todo introduce status codes for SoPlex, especially for rational solving
98
99 ///@todo record and return "best" solutions found during IR (Ambros)
100 ///@todo implement main IR loop for primal and dual feasible case with fail otherwise (Ambros)
101 ///@todo implement statistical info (time, factor time, iters, ...) since last call to solveReal() or solveRational() (Ambros?)
102 ///@todo implement performInfeasibilityIR (Ambros?)
103 ///@todo extend IR loop to infeasible case (Dan?)
104 ///@todo extend IR loop to unbounded case (Dan?)
105
106 ///@todo interface rational reconstruction code for rational vectors
107 ///@todo integrate rational reconstruction into IR loop
108 ///@todo integrate rational SPxSolver and distinguish between original and transformed rational LP
109 ///@todo rational scalers
110 ///@todo rational simplifier
111
112 namespace soplex
113 {
114
115 /**@class SoPlex
116 * @brief Preconfigured SoPlex LP-solver.
117 * @ingroup Algo
118 */
119 template <class R>
120 class SoPlexBase
121 {
122 public:
123
124 ///@name Construction and destruction
125 ///@{
126
127 /// default constructor
128 SoPlexBase();
129
130 /// assignment operator
131 SoPlexBase<R>& operator=(const SoPlexBase<R>& rhs);
132
133 /// copy constructor
134 SoPlexBase(const SoPlexBase<R>& rhs);
135
136 /// destructor
137 virtual ~SoPlexBase();
138
139 ///@}
140
141
142 ///@name Access to the real LP
143 ///@{
144
145 /// returns number of rows
146 int numRows() const;
147 int numRowsReal() const; /* For SCIP compatibility */
148 int numRowsRational() const;
149
150 /// Templated function that
151 /// returns number of columns
152 int numCols() const;
153 int numColsReal() const; /* For SCIP compatibility */
154 int numColsRational() const;
155
156 /// returns number of nonzeros
157 int numNonzeros() const;
158
159 int numNonzerosRational() const;
160
161 /// returns smallest non-zero element in absolute value
162 R minAbsNonzeroReal() const;
163
164 /// returns biggest non-zero element in absolute value
165 R maxAbsNonzeroReal() const;
166
167 /// returns (unscaled) coefficient
168 R coefReal(int row, int col) const;
169
170 /// returns vector of row \p i, ignoring scaling
171 const SVectorBase<R>& rowVectorRealInternal(int i) const;
172
173 /// gets vector of row \p i
174 void getRowVectorReal(int i, DSVectorBase<R>& row) const;
175
176 /// returns right-hand side vector, ignoring scaling
177 const VectorBase<R>& rhsRealInternal() const;
178
179 /// gets right-hand side vector
180 void getRhsReal(VectorBase<R>& rhs) const;
181
182 /// returns right-hand side of row \p i
183 R rhsReal(int i) const;
184
185 /// returns left-hand side vector, ignoring scaling
186 const VectorBase<R>& lhsRealInternal() const;
187
188 /// gets left-hand side vector
189 void getLhsReal(VectorBase<R>& lhs) const;
190
191 /// returns left-hand side of row \p i
192 R lhsReal(int i) const;
193
194 /// returns inequality type of row \p i
195 typename LPRowBase<R>::Type rowTypeReal(int i) const;
196
197 /// returns vector of col \p i, ignoring scaling
198 const SVectorBase<R>& colVectorRealInternal(int i) const;
199
200 /// gets vector of col \p i
201 void getColVectorReal(int i, DSVectorBase<R>& col) const;
202
203 /// returns upper bound vector
204 const VectorBase<R>& upperRealInternal() const;
205
206 /// returns upper bound of column \p i
207 R upperReal(int i) const;
208
209 /// gets upper bound vector
210 void getUpperReal(VectorBase<R>& upper) const;
211
212 /// returns lower bound vector
213 const VectorBase<R>& lowerRealInternal() const;
214
215 /// returns lower bound of column \p i
216 R lowerReal(int i) const;
217
218 /// gets lower bound vector
219 void getLowerReal(VectorBase<R>& lower) const;
220
221 /// gets objective function vector
222 void getObjReal(VectorBase<R>& obj) const;
223
224 /// returns objective value of column \p i
225 R objReal(int i) const;
226
227 /// returns objective function vector after transformation to a maximization problem; since this is how it is stored
228 /// internally, this is generally faster
229 const VectorBase<R>& maxObjRealInternal() const;
230
231 /// returns objective value of column \p i after transformation to a maximization problem; since this is how it is
232 /// stored internally, this is generally faster
233 R maxObjReal(int i) const;
234
235 /// gets number of available dual norms
236 void getNdualNorms(int& nnormsRow, int& nnormsCol) const;
237
238 /// gets steepest edge norms and returns false if they are not available
239 bool getDualNorms(int& nnormsRow, int& nnormsCol, R* norms) const;
240
241 /// sets steepest edge norms and returns false if that's not possible
242 bool setDualNorms(int nnormsRow, int nnormsCol, R* norms);
243
244 /// pass integrality information about the variables to the solver
245 void setIntegralityInformation(int ncols, int* intInfo);
246
247 ///@}
248
249
250 ///@name Access to the rational LP
251 ///@{
252
253 /// returns smallest non-zero element in absolute value
254 Rational minAbsNonzeroRational() const;
255
256 /// returns biggest non-zero element in absolute value
257 Rational maxAbsNonzeroRational() const;
258
259 /// gets row \p i
260 void getRowRational(int i, LPRowRational& lprow) const;
261
262 /// gets rows \p start, ..., \p end.
263 void getRowsRational(int start, int end, LPRowSetRational& lprowset) const;
264
265 /// returns vector of row \p i
266 const SVectorRational& rowVectorRational(int i) const;
267
268 /// returns right-hand side vector
269 const VectorRational& rhsRational() const;
270
271 /// returns right-hand side of row \p i
272 const Rational& rhsRational(int i) const;
273
274 /// returns left-hand side vector
275 const VectorRational& lhsRational() const;
276
277 /// returns left-hand side of row \p i
278 const Rational& lhsRational(int i) const;
279
280 /// returns inequality type of row \p i
281 LPRowRational::Type rowTypeRational(int i) const;
282
283 /// gets column \p i
284 void getColRational(int i, LPColRational& lpcol) const;
285
286 /// gets columns \p start, ..., \p end
287 void getColsRational(int start, int end, LPColSetRational& lpcolset) const;
288
289 /// returns vector of column \p i
290 const SVectorRational& colVectorRational(int i) const;
291
292 /// returns upper bound vector
293 const VectorRational& upperRational() const;
294
295 /// returns upper bound of column \p i
296 const Rational& upperRational(int i) const;
297
298 /// returns lower bound vector
299 const VectorRational& lowerRational() const;
300
301 /// returns lower bound of column \p i
302 const Rational& lowerRational(int i) const;
303
304 /// gets objective function vector
305 void getObjRational(VectorRational& obj) const;
306
307 /// gets objective value of column \p i
308 void getObjRational(int i, Rational& obj) const;
309
310 /// returns objective value of column \p i
311 Rational objRational(int i) const;
312
313 /// returns objective function vector after transformation to a maximization problem; since this is how it is stored
314 /// internally, this is generally faster
315 const VectorRational& maxObjRational() const;
316
317 /// returns objective value of column \p i after transformation to a maximization problem; since this is how it is
318 /// stored internally, this is generally faster
319 const Rational& maxObjRational(int i) const;
320
321 ///@}
322
323
324 ///@name Modification of the real LP
325 ///@{
326
327 /// adds a single row
328 void addRowReal(const LPRowBase<R>& lprow);
329
330 /// adds multiple rows
331 void addRowsReal(const LPRowSetBase<R>& lprowset);
332
333 /// adds a single column
334 void addColReal(const LPColBase<R>& lpcol);
335
336 /// adds multiple columns
337 void addColsReal(const LPColSetBase<R>& lpcolset);
338
339 /// replaces row \p i with \p lprow
340 void changeRowReal(int i, const LPRowBase<R>& lprow);
341
342 /// changes left-hand side vector for constraints to \p lhs
343 void changeLhsReal(const VectorBase<R>& lhs);
344
345 /// changes left-hand side of row \p i to \p lhs
346 void changeLhsReal(int i, const R& lhs);
347
348 /// changes right-hand side vector to \p rhs
349 void changeRhsReal(const VectorBase<R>& rhs);
350
351 /// changes right-hand side of row \p i to \p rhs
352 void changeRhsReal(int i, const R& rhs);
353
354 /// changes left- and right-hand side vectors
355 void changeRangeReal(const VectorBase<R>& lhs, const VectorBase<R>& rhs);
356
357 /// changes left- and right-hand side of row \p i
358 void changeRangeReal(int i, const R& lhs, const R& rhs);
359
360 /// replaces column \p i with \p lpcol
361 void changeColReal(int i, const LPColReal& lpcol);
362
363 /// changes vector of lower bounds to \p lower
364 void changeLowerReal(const VectorBase<R>& lower);
365
366 /// changes lower bound of column i to \p lower
367 void changeLowerReal(int i, const R& lower);
368
369 /// changes vector of upper bounds to \p upper
370 void changeUpperReal(const VectorBase<R>& upper);
371
372 /// changes \p i 'th upper bound to \p upper
373 void changeUpperReal(int i, const R& upper);
374
375 /// changes vectors of column bounds to \p lower and \p upper
376 void changeBoundsReal(const VectorBase<R>& lower, const VectorBase<R>& upper);
377
378 /// changes bounds of column \p i to \p lower and \p upper
379 void changeBoundsReal(int i, const R& lower, const R& upper);
380
381 /// changes objective function vector to \p obj
382 void changeObjReal(const VectorBase<R>& obj);
383
384 /// changes objective coefficient of column i to \p obj
385 void changeObjReal(int i, const R& obj);
386
387 /// changes matrix entry in row \p i and column \p j to \p val
388 void changeElementReal(int i, int j, const R& val);
389
390 /// removes row \p i
391 void removeRowReal(int i);
392
393 /// removes all rows with an index \p i such that \p perm[i] < 0; upon completion, \p perm[i] >= 0 indicates the
394 /// new index where row \p i has been moved to; note that \p perm must point to an array of size at least
395 /// #numRows()
396 void removeRowsReal(int perm[]);
397
398 /// remove all rows with indices in array \p idx of size \p n; an array \p perm of size #numRows() may be passed
399 /// as buffer memory
400 void removeRowsReal(int idx[], int n, int perm[] = 0);
401
402 /// removes rows \p start to \p end including both; an array \p perm of size #numRows() may be passed as buffer
403 /// memory
404 void removeRowRangeReal(int start, int end, int perm[] = 0);
405
406 /// removes column i
407 void removeColReal(int i);
408
409 /// removes all columns with an index \p i such that \p perm[i] < 0; upon completion, \p perm[i] >= 0 indicates the
410 /// new index where column \p i has been moved to; note that \p perm must point to an array of size at least
411 /// #numColsReal()
412 void removeColsReal(int perm[]);
413
414 /// remove all columns with indices in array \p idx of size \p n; an array \p perm of size #numColsReal() may be
415 /// passed as buffer memory
416 void removeColsReal(int idx[], int n, int perm[] = 0);
417
418 /// removes columns \p start to \p end including both; an array \p perm of size #numColsReal() may be passed as
419 /// buffer memory
420 void removeColRangeReal(int start, int end, int perm[] = 0);
421
422 /// clears the LP
423 void clearLPReal();
424
425 /// synchronizes real LP with rational LP, i.e., copies (rounded) rational LP into real LP, if sync mode is manual
426 void syncLPReal();
427
428 ///@}
429
430
431 ///@name Modification of the rational LP
432 ///@{
433
434 /// adds a single row
435 void addRowRational(const LPRowRational& lprow);
436
437 #ifdef SOPLEX_WITH_GMP
438 /// adds a single row (GMP only method)
439 void addRowRational(const mpq_t* lhs, const mpq_t* rowValues, const int* rowIndices,
440 const int rowSize, const mpq_t* rhs);
441
442 /// adds a set of rows (GMP only method)
443 void addRowsRational(const mpq_t* lhs, const mpq_t* rowValues, const int* rowIndices,
444 const int* rowStarts, const int* rowLengths, const int numRows, const int numValues,
445 const mpq_t* rhs);
446 #endif
447
448 /// adds multiple rows
449 void addRowsRational(const LPRowSetRational& lprowset);
450
451 /// adds a single column
452 void addColRational(const LPColRational& lpcol);
453
454 #ifdef SOPLEX_WITH_GMP
455 /// adds a single column (GMP only method)
456 void addColRational(const mpq_t* obj, const mpq_t* lower, const mpq_t* colValues,
457 const int* colIndices, const int colSize, const mpq_t* upper);
458
459 /// adds a set of columns (GMP only method)
460 void addColsRational(const mpq_t* obj, const mpq_t* lower, const mpq_t* colValues,
461 const int* colIndices, const int* colStarts, const int* colLengths, const int numCols,
462 const int numValues, const mpq_t* upper);
463 #endif
464
465 /// adds multiple columns
466 void addColsRational(const LPColSetRational& lpcolset);
467
468 /// replaces row \p i with \p lprow
469 void changeRowRational(int i, const LPRowRational& lprow);
470
471 /// changes left-hand side vector for constraints to \p lhs
472 void changeLhsRational(const VectorRational& lhs);
473
474 /// changes left-hand side of row \p i to \p lhs
475 void changeLhsRational(int i, const Rational& lhs);
476
477 #ifdef SOPLEX_WITH_GMP
478 /// changes left-hand side of row \p i to \p lhs (GMP only method)
479 void changeLhsRational(int i, const mpq_t* lhs);
480 #endif
481
482 /// changes right-hand side vector to \p rhs
483 void changeRhsRational(const VectorRational& rhs);
484
485 #ifdef SOPLEX_WITH_GMP
486 /// changes right-hand side vector to \p rhs (GMP only method)
487 void changeRhsRational(const mpq_t* rhs, int rhsSize);
488 #endif
489
490 /// changes right-hand side of row \p i to \p rhs
491 void changeRhsRational(int i, const Rational& rhs);
492
493 /// changes left- and right-hand side vectors
494 void changeRangeRational(const VectorRational& lhs, const VectorRational& rhs);
495
496 /// changes left- and right-hand side of row \p i
497 void changeRangeRational(int i, const Rational& lhs, const Rational& rhs);
498
499 #ifdef SOPLEX_WITH_GMP
500 /// changes left- and right-hand side of row \p i (GMP only method)
501 void changeRangeRational(int i, const mpq_t* lhs, const mpq_t* rhs);
502 #endif
503
504 /// replaces column \p i with \p lpcol
505 void changeColRational(int i, const LPColRational& lpcol);
506
507 /// changes vector of lower bounds to \p lower
508 void changeLowerRational(const VectorRational& lower);
509
510 /// changes lower bound of column i to \p lower
511 void changeLowerRational(int i, const Rational& lower);
512
513 #ifdef SOPLEX_WITH_GMP
514 /// changes lower bound of column i to \p lower (GMP only method)
515 void changeLowerRational(int i, const mpq_t* lower);
516 #endif
517
518 /// changes vector of upper bounds to \p upper
519 void changeUpperRational(const VectorRational& upper);
520
521 /// changes \p i 'th upper bound to \p upper
522 void changeUpperRational(int i, const Rational& upper);
523
524 #ifdef SOPLEX_WITH_GMP
525 /// changes upper bound of column i to \p upper (GMP only method)
526 void changeUpperRational(int i, const mpq_t* upper);
527 #endif
528
529 /// changes vectors of column bounds to \p lower and \p upper
530 void changeBoundsRational(const VectorRational& lower, const VectorRational& upper);
531
532 /// changes bounds of column \p i to \p lower and \p upper
533 void changeBoundsRational(int i, const Rational& lower, const Rational& upper);
534
535 #ifdef SOPLEX_WITH_GMP
536 /// changes bounds of column \p i to \p lower and \p upper (GMP only method)
537 void changeBoundsRational(int i, const mpq_t* lower, const mpq_t* upper);
538 #endif
539
540 /// changes objective function vector to \p obj
541 void changeObjRational(const VectorRational& obj);
542
543 /// changes objective coefficient of column i to \p obj
544 void changeObjRational(int i, const Rational& obj);
545
546 #ifdef SOPLEX_WITH_GMP
547 /// changes objective coefficient of column i to \p obj (GMP only method)
548 void changeObjRational(int i, const mpq_t* obj);
549 #endif
550
551 /// changes matrix entry in row \p i and column \p j to \p val
552 void changeElementRational(int i, int j, const Rational& val);
553
554 #ifdef SOPLEX_WITH_GMP
555 /// changes matrix entry in row \p i and column \p j to \p val (GMP only method)
556 void changeElementRational(int i, int j, const mpq_t* val);
557 #endif
558
559 /// removes row \p i
560 void removeRowRational(int i);
561
562 /// removes all rows with an index \p i such that \p perm[i] < 0; upon completion, \p perm[i] >= 0 indicates the new
563 /// index where row \p i has been moved to; note that \p perm must point to an array of size at least
564 /// #numRowsRational()
565 void removeRowsRational(int perm[]);
566
567 /// remove all rows with indices in array \p idx of size \p n; an array \p perm of size #numRowsRational() may be
568 /// passed as buffer memory
569 void removeRowsRational(int idx[], int n, int perm[] = 0);
570
571 /// removes rows \p start to \p end including both; an array \p perm of size #numRowsRational() may be passed as
572 /// buffer memory
573 void removeRowRangeRational(int start, int end, int perm[] = 0);
574
575 /// removes column i
576 void removeColRational(int i);
577
578 /// removes all columns with an index \p i such that \p perm[i] < 0; upon completion, \p perm[i] >= 0 indicates the
579 /// new index where column \p i has been moved to; note that \p perm must point to an array of size at least
580 /// #numColsRational()
581 void removeColsRational(int perm[]);
582
583 /// remove all columns with indices in array \p idx of size \p n; an array \p perm of size #numColsRational() may be
584 /// passed as buffer memory
585 void removeColsRational(int idx[], int n, int perm[] = 0);
586
587 /// removes columns \p start to \p end including both; an array \p perm of size #numColsRational() may be passed as
588 /// buffer memory
589 void removeColRangeRational(int start, int end, int perm[] = 0);
590
591 /// clears the LP
592 void clearLPRational();
593
594 /// synchronizes rational LP with real LP, i.e., copies real LP to rational LP, if sync mode is manual
595 void syncLPRational();
596
597 ///@}
598
599 ///@name Solving and general solution query
600 ///@{
601
602 /// optimize the given LP
603 typename SPxSolverBase<R>::Status optimize(volatile bool* interrupt = NULL);
604
605 // old name for backwards compatibility
606 typename SPxSolverBase<R>::Status solve(volatile bool* interrupt = NULL)
607 {
608 return optimize(interrupt);
609 }
610
611 /// returns the current solver status
612 typename SPxSolverBase<R>::Status status() const;
613
614 /// is stored primal solution feasible?
615 bool isPrimalFeasible() const;
616
617 /// is a solution available (not neccessarily feasible)?
618 bool hasSol() const;
619
620 /// deprecated: use #hasSol() instead
621 bool hasPrimal() const
622 {
623 return hasSol();
624 }
625
626 /// deprecated: use #hasSol() instead
627 bool hasDual() const
628 {
629 return hasSol();
630 }
631
632 /// is a primal unbounded ray available?
633 bool hasPrimalRay() const;
634
635 /// is stored dual solution feasible?
636 bool isDualFeasible() const;
637
638 /// is Farkas proof of infeasibility available?
639 bool hasDualFarkas() const;
640
641 /// sets the status to OPTIMAL in case the LP has been solved with unscaled violations
642 bool ignoreUnscaledViolations()
643 {
644 if(_status == SPxSolverBase<R>::OPTIMAL_UNSCALED_VIOLATIONS)
645 {
646 _status = SPxSolverBase<R>::OPTIMAL;
647 return true;
648 }
649 else
650 return false;
651 }
652 ///@}
653
654
655 ///@name Query for the real solution data
656 ///@{
657
658 /// returns the objective value if a primal solution is available
659 R objValueReal();
660
661 /// gets the primal solution vector if available; returns true on success
662 bool getPrimal(VectorBase<R>& vector);
663 bool getPrimalReal(R* p_vector, int size); // For SCIP compatibility
664 bool getPrimalRational(VectorRational& vector);
665
666 /// gets the vector of slack values if available; returns true on success
667 bool getSlacksReal(VectorBase<R>& vector);
668 bool getSlacksReal(R* p_vector, int dim);
669
670 /// gets the primal ray if available; returns true on success
671 bool getPrimalRay(VectorBase<R>& vector);
672 bool getPrimalRayReal(R* vector, int dim); /* For SCIP compatibility */
673 bool getPrimalRayRational(VectorRational& vector);
674
675 /// gets the dual solution vector if available; returns true on success
676 bool getDual(VectorBase<R>& vector);
677 bool getDualReal(R* p_vector, int dim); /* For SCIP compatibility */
678 bool getDualRational(VectorRational& vector);
679
680 /// gets the vector of reduced cost values if available; returns true on success
681 bool getRedCost(VectorBase<R>& vector);
682 bool getRedCostReal(R* vector, int dim); /* For SCIP compatibility */
683 bool getRedCostRational(VectorRational& vector);
684
685 /// gets the Farkas proof if available; returns true on success
686 bool getDualFarkas(VectorBase<R>& vector);
687 bool getDualFarkasReal(R* vector, int dim);
688 bool getDualFarkasRational(VectorRational& vector);
689
690 /// gets violation of bounds; returns true on success
691 bool getBoundViolation(R& maxviol, R& sumviol);
692 bool getBoundViolationRational(Rational& maxviol, Rational& sumviol);
693
694 /// gets violation of constraints; returns true on success
695 bool getRowViolation(R& maxviol, R& sumviol);
696 bool getRowViolationRational(Rational& maxviol, Rational& sumviol);
697
698 /// gets violation of reduced costs; returns true on success
699 bool getRedCostViolation(R& maxviol, R& sumviol);
700 bool getRedCostViolationRational(Rational& maxviol, Rational& sumviol);
701
702 /// gets violation of dual multipliers; returns true on success
703 bool getDualViolation(R& maxviol, R& sumviol);
704 bool getDualViolationRational(Rational& maxviol, Rational& sumviol);
705
706 ///@}
707
708
709 ///@name Query for the rational solution data
710 ///@{
711
712 /// returns the objective value if a primal solution is available
713 Rational objValueRational();
714
715 /// gets the vector of slack values if available; returns true on success
716 bool getSlacksRational(VectorRational& vector);
717
718 #ifdef SOPLEX_WITH_GMP
719 /// gets the primal solution vector if available; returns true on success (GMP only method)
720 bool getPrimalRational(mpq_t* vector, const int size);
721
722 /// gets the vector of slack values if available; returns true on success (GMP only method)
723 bool getSlacksRational(mpq_t* vector, const int size);
724
725 /// gets the primal ray if LP is unbounded; returns true on success (GMP only method)
726 bool getPrimalRayRational(mpq_t* vector, const int size);
727
728 /// gets the dual solution vector if available; returns true on success (GMP only method)
729 bool getDualRational(mpq_t* vector, const int size);
730
731 /// gets the vector of reduced cost values if available; returns true on success (GMP only method)
732 bool getRedCostRational(mpq_t* vector, const int size);
733
734 /// gets the Farkas proof if LP is infeasible; returns true on success (GMP only method)
735 bool getDualFarkasRational(mpq_t* vector, const int size);
736 #endif
737
738 /// get size of primal solution
739 int totalSizePrimalRational(const int base = 2);
740
741 /// get size of dual solution
742 int totalSizeDualRational(const int base = 2);
743
744 /// get size of least common multiple of denominators in primal solution
745 int dlcmSizePrimalRational(const int base = 2);
746
747 /// get size of least common multiple of denominators in dual solution
748 int dlcmSizeDualRational(const int base = 2);
749
750 /// get size of largest denominator in primal solution
751 int dmaxSizePrimalRational(const int base = 2);
752
753 /// get size of largest denominator in dual solution
754 int dmaxSizeDualRational(const int base = 2);
755
756 ///@}
757
758
759 ///@name Access and modification of basis information
760 ///@{
761
762 /// is an advanced starting basis available?
763 bool hasBasis() const;
764
765 /// returns the current basis status
766 typename SPxBasisBase<R>::SPxStatus basisStatus() const;
767
768 /// returns basis status for a single row
769 typename SPxSolverBase<R>::VarStatus basisRowStatus(int row) const;
770
771 /// returns basis status for a single column
772 typename SPxSolverBase<R>::VarStatus basisColStatus(int col) const;
773
774 /// gets current basis via arrays of statuses
775 void getBasis(typename SPxSolverBase<R>::VarStatus rows[],
776 typename SPxSolverBase<R>::VarStatus cols[]) const;
777
778 /// gets the indices of the basic columns and rows; basic column n gives value n, basic row m gives value -1-m
779 void getBasisInd(int* bind) const;
780
781 /** compute one of several matrix metrics based on the diagonal of the LU factorization
782 * type = 0: max/min ratio
783 * type = 1: trace of U (sum of diagonal elements)
784 * type = 2: determinant (product of diagonal elements)
785 */
786 bool getBasisMetric(R& metric, int type = 0);
787
788 /// computes an estimated condition number for the current basis matrix using the power method; returns true on success
789 bool getEstimatedCondition(R& condition);
790
791 /// computes the exact condition number for the current basis matrix using the power method; returns true on success
792 bool getExactCondition(R& condition);
793
794 /// computes row \p r of basis inverse; returns true on success
795 /// @param r which row of the basis inverse is computed
796 /// @param coef values of result vector (not packed but scattered)
797 /// @param inds indices of result vector (NULL if not to be used)
798 /// @param ninds number of nonzeros in result vector
799 /// @param unscale determines whether the result should be unscaled according to the original LP data
800 bool getBasisInverseRowReal(int r, R* coef, int* inds = NULL, int* ninds = NULL,
801 bool unscale = true);
802
803 /// computes column \p c of basis inverse; returns true on success
804 /// @param c which column of the basis inverse is computed
805 /// @param coef values of result vector (not packed but scattered)
806 /// @param inds indices of result vector (NULL if not to be used)
807 /// @param ninds number of nonzeros in result vector
808 /// @param unscale determines whether the result should be unscaled according to the original LP data
809 bool getBasisInverseColReal(int c, R* coef, int* inds = NULL, int* ninds = NULL,
810 bool unscale = true);
811
812 /// computes dense solution of basis matrix B * \p sol = \p rhs; returns true on success
813 bool getBasisInverseTimesVecReal(R* rhs, R* sol, bool unscale = true);
814
815 /// multiply with basis matrix; B * \p vec (inplace)
816 /// @param vec (dense) vector to be multiplied with
817 /// @param unscale determines whether the result should be unscaled according to the original LP data
818 bool multBasis(R* vec, bool unscale = true);
819
820 /// multiply with transpose of basis matrix; \p vec * B^T (inplace)
821 /// @param vec (dense) vector to be multiplied with
822 /// @param unscale determines whether the result should be unscaled according to the original LP data
823 bool multBasisTranspose(R* vec, bool unscale = true);
824
825 /// compute rational basis inverse; returns true on success
826 bool computeBasisInverseRational();
827
828 /// gets an array of indices for the columns of the rational basis matrix; bind[i] >= 0 means that the i-th column of
829 /// the basis matrix contains variable bind[i]; bind[i] < 0 means that the i-th column of the basis matrix contains
830 /// the slack variable for row -bind[i]-1; performs rational factorization if not available; returns true on success
831 bool getBasisIndRational(DataArray<int>& bind);
832
833 /// computes row r of basis inverse; performs rational factorization if not available; returns true on success
834 bool getBasisInverseRowRational(const int r, SSVectorRational& vec);
835
836 /// computes column c of basis inverse; performs rational factorization if not available; returns true on success
837 bool getBasisInverseColRational(const int c, SSVectorRational& vec);
838
839 /// computes solution of basis matrix B * sol = rhs; performs rational factorization if not available; returns true
840 /// on success
841 bool getBasisInverseTimesVecRational(const SVectorRational& rhs, SSVectorRational& sol);
842
843 /// sets starting basis via arrays of statuses
844 void setBasis(const typename SPxSolverBase<R>::VarStatus rows[],
845 const typename SPxSolverBase<R>::VarStatus cols[]);
846
847 /// clears starting basis
848 void clearBasis();
849
850 ///@}
851
852
853 ///@name Statistical information
854 ///@{
855
856 /// number of iterations since last call to solve
857 int numIterations() const;
858
859 /// number of precision boosts since last call to solve
860 int numPrecisionBoosts() const;
861
862 /// number of iterations in higher precision since last call to solve
863 int numIterationsBoosted() const;
864
865 /// time spen in higher precision since last call to solve
866 Real precisionBoostTime() const;
867
868 /// time spent in last call to solve
869 Real solveTime() const;
870
871 /// statistical information in form of a string
872 std::string statisticString() const;
873
874 /// name of starter
875 const char* getStarterName();
876
877 /// name of simplifier
878 const char* getSimplifierName();
879
880 /// name of scaling method
881 const char* getScalerName();
882
883 /// name of currently loaded pricer
884 const char* getPricerName();
885
886 /// name of currently loaded ratiotester
887 const char* getRatiotesterName();
888
889 ///@}
890
891
892 ///@name File I/O
893 ///@{
894
895 /// reads LP file in LP or MPS format according to READMODE parameter; gets row names, column names, and
896 /// integer variables if desired; returns true on success
897 bool readFile(const char* filename, NameSet* rowNames = 0, NameSet* colNames = 0,
898 DIdxSet* intVars = 0);
899
900 /// Templated write function
901 /// Real
902 /// writes real LP to file; LP or MPS format is chosen from the extension in \p filename; if \p rowNames and \p
903 /// colNames are \c NULL, default names are used; if \p intVars is not \c NULL, the variables contained in it are
904 /// marked as integer; returns true on success
905 /// Rational
906 /// writes rational LP to file; LP or MPS format is chosen from the extension in \p filename; if \p rowNames and \p
907 /// colNames are \c NULL, default names are used; if \p intVars is not \c NULL, the variables contained in it are
908 /// marked as integer; returns true on success
909 bool writeFile(const char* filename, const NameSet* rowNames = 0, const NameSet* colNames = 0,
910 const DIdxSet* intvars = 0, const bool unscale = true) const;
911
912 bool writeFileRational(const char* filename, const NameSet* rowNames = 0,
913 const NameSet* colNames = 0, const DIdxSet* intvars = 0) const;
914
915 /* For SCIP compatibility */
916 bool writeFileReal(const char* filename, const NameSet* rowNames = 0, const NameSet* colNames = 0,
917 const DIdxSet* intvars = 0, const bool unscale = true) const;
918
919 /// writes the dual of the real LP to file; LP or MPS format is chosen from the extension in \p filename;
920 /// if \p rowNames and \p colNames are \c NULL, default names are used; if \p intVars is not \c NULL,
921 /// the variables contained in it are marked as integer; returns true on success
922 bool writeDualFileReal(const char* filename, const NameSet* rowNames = 0,
923 const NameSet* colNames = 0, const DIdxSet* intvars = 0) const;
924
925 /// reads basis information from \p filename and returns true on success; if \p rowNames and \p colNames are \c NULL,
926 /// default names are assumed; returns true on success
927 bool readBasisFile(const char* filename, const NameSet* rowNames = 0, const NameSet* colNames = 0);
928
929 /// writes basis information to \p filename; if \p rowNames and \p colNames are \c NULL, default names are used;
930 /// returns true on success
931 bool writeBasisFile(const char* filename, const NameSet* rowNames = 0, const NameSet* colNames = 0,
932 const bool cpxFormat = false) const;
933
934 /// writes internal LP, basis information, and parameter settings; if \p rowNames and \p colNames are \c NULL,
935 /// default names are used
936 void writeStateReal(const char* filename, const NameSet* rowNames = 0, const NameSet* colNames = 0,
937 const bool cpxFormat = false) const;
938
939 /// writes internal LP, basis information, and parameter settings; if \p rowNames and \p colNames are \c NULL,
940 /// default names are used
941 void writeStateRational(const char* filename, const NameSet* rowNames = 0,
942 const NameSet* colNames = 0, const bool cpxFormat = false) const;
943
944 ///@}
945
946
947 ///@name Parameters
948 ///@{
949
950 /// boolean parameters
951 typedef enum
952 {
953 /// should lifting be used to reduce range of nonzero matrix coefficients?
954 LIFTING = 0,
955
956 /// should LP be transformed to equality form before a rational solve?
957 EQTRANS = 1,
958
959 /// should dual infeasibility be tested in order to try to return a dual solution even if primal infeasible?
960 TESTDUALINF = 2,
961
962 /// should a rational factorization be performed after iterative refinement?
963 RATFAC = 3,
964
965 /// should the decomposition based dual simplex be used to solve the LP? Setting this to true forces the solve mode to
966 /// SOLVEMODE_REAL and the basis representation to REPRESENTATION_ROW
967 USEDECOMPDUALSIMPLEX = 4,
968
969 /// should the degeneracy be computed for each basis?
970 COMPUTEDEGEN = 5,
971
972 /// should the dual of the complementary problem be used in the decomposition simplex?
973 USECOMPDUAL = 6,
974
975 /// should row and bound violations be computed explicitly in the update of reduced problem in the decomposition simplex
976 EXPLICITVIOL = 7,
977
978 /// should cycling solutions be accepted during iterative refinement?
979 ACCEPTCYCLING = 8,
980
981 /// apply rational reconstruction after each iterative refinement?
982 RATREC = 9,
983
984 /// round scaling factors for iterative refinement to powers of two?
985 POWERSCALING = 10,
986
987 /// continue iterative refinement with exact basic solution if not optimal?
988 RATFACJUMP = 11,
989
990 /// use bound flipping also for row representation?
991 ROWBOUNDFLIPS = 12,
992
993 /// use persistent scaling?
994 PERSISTENTSCALING = 13,
995
996 /// perturb the entire problem or only the relevant bounds of s single pivot?
997 FULLPERTURBATION = 14,
998
999 /// re-optimize the original problem to get a proof (ray) of infeasibility/unboundedness?
1000 ENSURERAY = 15,
1001
1002 /// try to enforce that the optimal solution is a basic solution
1003 FORCEBASIC = 16,
1004
1005 // enable presolver SingletonCols in PaPILO?
1006 SIMPLIFIER_SINGLETONCOLS = 17,
1007
1008 // enable presolver ConstraintPropagation in PaPILO?
1009 SIMPLIFIER_CONSTRAINTPROPAGATION = 18,
1010
1011 // enable presolver ParallelRowDetection in PaPILO?
1012 SIMPLIFIER_PARALLELROWDETECTION = 19,
1013
1014 // enable presolver ParallelColDetection in PaPILO?
1015 SIMPLIFIER_PARALLELCOLDETECTION = 20,
1016
1017 // enable presolver SingletonStuffing in PaPILO?
1018 SIMPLIFIER_SINGLETONSTUFFING = 21,
1019
1020 // enable presolver DualFix in PaPILO?
1021 SIMPLIFIER_DUALFIX = 22,
1022
1023 // enable presolver FixContinuous in PaPILO?
1024 SIMPLIFIER_FIXCONTINUOUS = 23,
1025
1026 // enable presolver DominatedCols in PaPILO?
1027 SIMPLIFIER_DOMINATEDCOLS = 24,
1028
1029 // enable iterative refinement ?
1030 ITERATIVE_REFINEMENT = 25,
1031
1032 /// adapt tolerances to the multiprecision used
1033 ADAPT_TOLS_TO_MULTIPRECISION = 26,
1034
1035 /// enable precision boosting ?
1036 PRECISION_BOOSTING = 27,
1037
1038 /// boosted solver start from last basis
1039 BOOSTED_WARM_START = 28,
1040
1041 /// try different settings when solve fails
1042 RECOVERY_MECHANISM = 29,
1043
1044 /// store advanced and stable basis met before each simplex iteration, to better warm start
1045 STORE_BASIS_BEFORE_SIMPLEX_PIVOT = 30,
1046
1047 /// number of boolean parameters
1048 BOOLPARAM_COUNT = 31
1049 } BoolParam;
1050
1051 /// integer parameters
1052 typedef enum
1053 {
1054 /// objective sense
1055 OBJSENSE = 0,
1056
1057 /// type of computational form, i.e., column or row representation
1058 REPRESENTATION = 1,
1059
1060 /// type of algorithm, i.e., primal or dual
1061 ALGORITHM = 2,
1062
1063 /// type of LU update
1064 FACTOR_UPDATE_TYPE = 3,
1065
1066 /// maximum number of updates without fresh factorization
1067 FACTOR_UPDATE_MAX = 4,
1068
1069 /// iteration limit (-1 if unlimited)
1070 ITERLIMIT = 5,
1071
1072 /// refinement limit (-1 if unlimited)
1073 REFLIMIT = 6,
1074
1075 /// stalling refinement limit (-1 if unlimited)
1076 STALLREFLIMIT = 7,
1077
1078 /// display frequency
1079 DISPLAYFREQ = 8,
1080
1081 /// verbosity level
1082 VERBOSITY = 9,
1083
1084 /// type of simplifier
1085 SIMPLIFIER = 10,
1086
1087 /// type of scaler
1088 SCALER = 11,
1089
1090 /// type of starter used to create crash basis
1091 STARTER = 12,
1092
1093 /// type of pricer
1094 PRICER = 13,
1095
1096 /// type of ratio test
1097 RATIOTESTER = 14,
1098
1099 /// mode for synchronizing real and rational LP
1100 SYNCMODE = 15,
1101
1102 /// mode for reading LP files
1103 READMODE = 16,
1104
1105 /// mode for iterative refinement strategy
1106 SOLVEMODE = 17,
1107
1108 /// mode for a posteriori feasibility checks
1109 CHECKMODE = 18,
1110
1111 /// type of timer
1112 TIMER = 19,
1113
1114 /// mode for hyper sparse pricing
1115 HYPER_PRICING = 20,
1116
1117 /// minimum number of stalling refinements since last pivot to trigger rational factorization
1118 RATFAC_MINSTALLS = 21,
1119
1120 /// maximum number of conjugate gradient iterations in least square scaling
1121 LEASTSQ_MAXROUNDS = 22,
1122
1123 /// mode for solution polishing
1124 SOLUTION_POLISHING = 23,
1125
1126 /// the number of iterations before the decomposition simplex initialisation is terminated.
1127 DECOMP_ITERLIMIT = 24,
1128
1129 /// the maximum number of rows that are added in each iteration of the decomposition based simplex
1130 DECOMP_MAXADDEDROWS = 25,
1131
1132 /// the iteration frequency at which the decomposition solve output is displayed.
1133 DECOMP_DISPLAYFREQ = 26,
1134
1135 /// the verbosity of the decomposition based simplex
1136 DECOMP_VERBOSITY = 27,
1137
1138 /// print condition number during the solve
1139 PRINTBASISMETRIC = 28,
1140
1141 /// type of timer for statistics
1142 STATTIMER = 29,
1143
1144 // maximum number of digits for the multiprecision type
1145 MULTIPRECISION_LIMIT = 30,
1146
1147 ///@todo precision-boosting find better parameter name
1148 /// after how many simplex pivots do we store the advanced and stable basis, 1 = every iterations
1149 STORE_BASIS_SIMPLEX_FREQ = 31,
1150
1151 /// number of integer parameters
1152 INTPARAM_COUNT = 32
1153 } IntParam;
1154
1155 /// values for parameter OBJSENSE
1156 enum
1157 {
1158 /// minimization
1159 OBJSENSE_MINIMIZE = -1,
1160
1161 /// maximization
1162 OBJSENSE_MAXIMIZE = 1
1163 };
1164
1165 /// values for parameter REPRESENTATION
1166 enum
1167 {
1168 /// automatic choice according to number of rows and columns
1169 REPRESENTATION_AUTO = 0,
1170
1171 /// column representation Ax - s = 0, lower <= x <= upper, lhs <= s <= rhs
1172 REPRESENTATION_COLUMN = 1,
1173
1174 /// row representation (lower,lhs) <= (x,Ax) <= (upper,rhs)
1175 REPRESENTATION_ROW = 2
1176 };
1177
1178 /// values for parameter ALGORITHM
1179 enum
1180 {
1181 /// primal simplex algorithm, i.e., entering for column and leaving for row representation
1182 ALGORITHM_PRIMAL = 0,
1183
1184 /// dual simplex algorithm, i.e., leaving for column and entering for row representation
1185 ALGORITHM_DUAL = 1
1186 };
1187
1188 /// values for parameter FACTOR_UPDATE_TYPE
1189 enum
1190 {
1191 /// product form update
1192 FACTOR_UPDATE_TYPE_ETA = 0,
1193
1194 /// Forrest-Tomlin type update
1195 FACTOR_UPDATE_TYPE_FT = 1
1196 };
1197
1198 /// values for parameter VERBOSITY
1199 enum
1200 {
1201 /// only error output
1202 VERBOSITY_ERROR = 0,
1203
1204 /// only error and warning output
1205 VERBOSITY_WARNING = 1,
1206
1207 /// only error, warning, and debug output
1208 VERBOSITY_DEBUG = 2,
1209
1210 /// standard verbosity level
1211 VERBOSITY_NORMAL = 3,
1212
1213 /// high verbosity level
1214 VERBOSITY_HIGH = 4,
1215
1216 /// full verbosity level
1217 VERBOSITY_FULL = 5
1218 };
1219
1220 /// values for parameter SIMPLIFIER
1221 enum
1222 {
1223 /// disabling presolving
1224 SIMPLIFIER_OFF = 0,
1225
1226 /// using internal presolving methods
1227 SIMPLIFIER_INTERNAL = 3,
1228
1229 /// using the presolve lib papilo
1230 SIMPLIFIER_PAPILO = 2,
1231
1232 /// SoPlex chooses automatically (currently always "internal")
1233 SIMPLIFIER_AUTO = 1
1234 };
1235
1236 /// values for parameter SCALER
1237 enum
1238 {
1239 /// no scaler
1240 SCALER_OFF = 0,
1241
1242 /// equilibrium scaling on rows or columns
1243 SCALER_UNIEQUI = 1,
1244
1245 /// equilibrium scaling on rows and columns
1246 SCALER_BIEQUI = 2,
1247
1248 /// geometric mean scaling on rows and columns, max 1 round
1249 SCALER_GEO1 = 3,
1250
1251 /// geometric mean scaling on rows and columns, max 8 rounds
1252 SCALER_GEO8 = 4,
1253
1254 /// least square scaling
1255 SCALER_LEASTSQ = 5,
1256
1257 /// geometric mean scaling (max 8 rounds) followed by equilibrium scaling (rows and columns)
1258 SCALER_GEOEQUI = 6
1259 };
1260
1261 /// values for parameter STARTER
1262 enum
1263 {
1264 /// slack basis
1265 STARTER_OFF = 0,
1266
1267 /// greedy crash basis weighted by objective, bounds, and sides
1268 STARTER_WEIGHT = 1,
1269
1270 /// crash basis from a greedy solution
1271 STARTER_SUM = 2,
1272
1273 /// generic solution-based crash basis
1274 STARTER_VECTOR = 3
1275 };
1276
1277 /// values for parameter PRICER
1278 enum
1279 {
1280 /// automatic pricer
1281 PRICER_AUTO = 0,
1282
1283 /// Dantzig pricer
1284 PRICER_DANTZIG = 1,
1285
1286 /// partial multiple pricer based on Dantzig pricing
1287 PRICER_PARMULT = 2,
1288
1289 /// devex pricer
1290 PRICER_DEVEX = 3,
1291
1292 /// steepest edge pricer with initialization to unit norms
1293 PRICER_QUICKSTEEP = 4,
1294
1295 /// steepest edge pricer with exact initialization of norms
1296 PRICER_STEEP = 5
1297 };
1298
1299 /// values for parameter RATIOTESTER
1300 enum
1301 {
1302 /// textbook ratio test without stabilization
1303 RATIOTESTER_TEXTBOOK = 0,
1304
1305 /// standard Harris ratio test
1306 RATIOTESTER_HARRIS = 1,
1307
1308 /// modified Harris ratio test
1309 RATIOTESTER_FAST = 2,
1310
1311 /// bound flipping ratio test for long steps in the dual simplex
1312 RATIOTESTER_BOUNDFLIPPING = 3
1313 };
1314
1315 /// values for parameter SYNCMODE
1316 enum
1317 {
1318 /// store only real LP
1319 SYNCMODE_ONLYREAL = 0,
1320
1321 /// automatic sync of real and rational LP
1322 SYNCMODE_AUTO = 1,
1323
1324 /// user sync of real and rational LP
1325 SYNCMODE_MANUAL = 2
1326 };
1327
1328 /// values for parameter READMODE
1329 enum
1330 {
1331 /// standard floating-point parsing
1332 READMODE_REAL = 0,
1333
1334 /// rational parsing
1335 READMODE_RATIONAL = 1
1336 };
1337
1338 /// values for parameter SOLVEMODE
1339 enum
1340 {
1341 /// apply standard floating-point algorithm
1342 SOLVEMODE_REAL = 0,
1343
1344 /// decide depending on tolerances whether to apply iterative refinement
1345 SOLVEMODE_AUTO = 1,
1346
1347 /// force iterative refinement
1348 SOLVEMODE_RATIONAL = 2
1349 };
1350
1351 /// values for parameter CHECKMODE
1352 enum
1353 {
1354 /// floating-point check
1355 CHECKMODE_REAL = 0,
1356
1357 /// decide according to READMODE
1358 CHECKMODE_AUTO = 1,
1359
1360 /// rational check
1361 CHECKMODE_RATIONAL = 2
1362 };
1363
1364 /// values for parameter TIMER
1365 enum
1366 {
1367 /// disable timing
1368 TIMER_OFF = 0,
1369
1370 /// cpu or user time
1371 TIMER_CPU = 1,
1372
1373 /// wallclock time
1374 TIMER_WALLCLOCK = 2
1375 };
1376
1377 /// values for parameter HYPER_PRICING
1378 enum
1379 {
1380 /// never
1381 HYPER_PRICING_OFF = 0,
1382
1383 /// decide according to problem size
1384 HYPER_PRICING_AUTO = 1,
1385
1386 /// always
1387 HYPER_PRICING_ON = 2
1388 };
1389
1390 /// values for parameter SOLUTION_POLISHING
1391 enum
1392 {
1393 /// no solution polishing
1394 POLISHING_OFF = 0,
1395
1396 /// maximize number of basic slack variables, i.e. more variables on bounds
1397 POLISHING_INTEGRALITY = 1,
1398
1399 /// minimize number of basic slack variables, i.e. more variables between bounds
1400 POLISHING_FRACTIONALITY = 2
1401 };
1402
1403 /// real parameters
1404 typedef enum
1405 {
1406 /// primal feasibility tolerance
1407 FEASTOL = 0,
1408
1409 /// dual feasibility tolerance
1410 OPTTOL = 1,
1411
1412 /// general zero tolerance
1413 EPSILON_ZERO = 2,
1414
1415 /// zero tolerance used in factorization
1416 EPSILON_FACTORIZATION = 3,
1417
1418 /// zero tolerance used in update of the factorization
1419 EPSILON_UPDATE = 4,
1420
1421 /// pivot zero tolerance used in factorization
1422 EPSILON_PIVOT = 5,
1423
1424 /// infinity threshold
1425 INFTY = 6,
1426
1427 /// time limit in seconds (INFTY if unlimited)
1428 TIMELIMIT = 7,
1429
1430 /// lower limit on objective value
1431 OBJLIMIT_LOWER = 8,
1432
1433 /// upper limit on objective value
1434 OBJLIMIT_UPPER = 9,
1435
1436 /// working tolerance for feasibility in floating-point solver during iterative refinement
1437 FPFEASTOL = 10,
1438
1439 /// working tolerance for optimality in floating-point solver during iterative refinement
1440 FPOPTTOL = 11,
1441
1442 /// maximum increase of scaling factors between refinements
1443 MAXSCALEINCR = 12,
1444
1445 /// lower threshold in lifting (nonzero matrix coefficients with smaller absolute value will be reformulated)
1446 LIFTMINVAL = 13,
1447
1448 /// upper threshold in lifting (nonzero matrix coefficients with larger absolute value will be reformulated)
1449 LIFTMAXVAL = 14,
1450
1451 /// sparse pricing threshold (\#violations < dimension * SPARSITY_THRESHOLD activates sparse pricing)
1452 SPARSITY_THRESHOLD = 15,
1453
1454 /// threshold on number of rows vs. number of columns for switching from column to row representations in auto mode
1455 REPRESENTATION_SWITCH = 16,
1456
1457 /// geometric frequency at which to apply rational reconstruction
1458 RATREC_FREQ = 17,
1459
1460 /// minimal reduction (sum of removed rows/cols) to continue simplification
1461 MINRED = 18,
1462
1463 /// refactor threshold for nonzeros in last factorized basis matrix compared to updated basis matrix
1464 REFAC_BASIS_NNZ = 19,
1465
1466 /// refactor threshold for fill-in in current factor update compared to fill-in in last factorization
1467 REFAC_UPDATE_FILL = 20,
1468
1469 /// refactor threshold for memory growth in factorization since last refactorization
1470 REFAC_MEM_FACTOR = 21,
1471
1472 /// accuracy of conjugate gradient method in least squares scaling (higher value leads to more iterations)
1473 LEASTSQ_ACRCY = 22,
1474
1475 /// objective offset
1476 OBJ_OFFSET = 23,
1477
1478 /// minimal Markowitz threshold to control sparsity/stability in LU factorization
1479 MIN_MARKOWITZ = 24,
1480
1481 /// minimal modification threshold to apply presolve reductions
1482 SIMPLIFIER_MODIFYROWFAC = 25,
1483
1484 /// factor by which the precision of the floating-point solver is multiplied
1485 PRECISION_BOOSTING_FACTOR = 26,
1486
1487 /// number of real parameters
1488 REALPARAM_COUNT = 27
1489 } RealParam;
1490
1491 #ifdef SOPLEX_WITH_RATIONALPARAM
1492 /// rational parameters
1493 typedef enum
1494 {
1495 /// number of rational parameters
1496 RATIONALPARAM_COUNT = 0
1497 } RationalParam;
1498 #endif
1499
1500 /// class of parameter settings
(1) Event rule_of_three_violation: |
Class "soplex::SoPlexBase<double>::Settings" has a user definition for at least one special function (copy constructor, copy assignment, destructor) but not all. If one of these functions requires a user definition then the others likely do as well. |
(4) Event remediation: |
Add user-definition for a destructor. |
Also see events: |
[copy_ctor][copy_assign] |
1501 class Settings
1502 {
1503 public:
1504 static struct BoolParam
1505 {
1506 /// constructor
1507 BoolParam();
1508 /// array of names for boolean parameters
1509 std::string name[SoPlexBase<R>::BOOLPARAM_COUNT];
1510 /// array of descriptions for boolean parameters
1511 std::string description[SoPlexBase<R>::BOOLPARAM_COUNT];
1512 /// array of default values for boolean parameters
1513 bool defaultValue[SoPlexBase<R>::BOOLPARAM_COUNT];
1514 } boolParam;
1515
1516 static struct IntParam
1517 {
1518 /// constructor
1519 IntParam();
1520 /// array of names for integer parameters
1521 std::string name[SoPlexBase<R>::INTPARAM_COUNT];
1522 /// array of descriptions for integer parameters
1523 std::string description[SoPlexBase<R>::INTPARAM_COUNT];
1524 /// array of default values for integer parameters
1525 int defaultValue[SoPlexBase<R>::INTPARAM_COUNT];
1526 /// array of lower bounds for int parameter values
1527 int lower[SoPlexBase<R>::INTPARAM_COUNT];
1528 /// array of upper bounds for int parameter values
1529 int upper[SoPlexBase<R>::INTPARAM_COUNT];
1530 } intParam;
1531
1532 static struct RealParam
1533 {
1534 /// constructor
1535 RealParam();
1536 /// array of names for real parameters
1537 std::string name[SoPlexBase<R>::REALPARAM_COUNT];
1538 /// array of descriptions for real parameters
1539 std::string description[SoPlexBase<R>::REALPARAM_COUNT];
1540 /// array of default values for real parameters
1541 Real defaultValue[SoPlexBase<R>::REALPARAM_COUNT];
1542 /// array of lower bounds for real parameter values
1543 Real lower[SoPlexBase<R>::REALPARAM_COUNT];
1544 /// array of upper bounds for real parameter values
1545 Real upper[SoPlexBase<R>::REALPARAM_COUNT];
1546 } realParam;
1547
1548 #ifdef SOPLEX_WITH_RATIONALPARAM
1549 static struct RationalParam
1550 {
1551 /// constructor
1552 RationalParam();
1553 /// array of names for rational parameters
1554 std::string name[SoPlexBase<R>::RATIONALPARAM_COUNT];
1555 /// array of descriptions for rational parameters
1556 std::string description[SoPlexBase<R>::RATIONALPARAM_COUNT];
1557 /// array of default values for rational parameters
1558 Rational defaultValue[SoPlexBase<R>::RATIONALPARAM_COUNT];
1559 /// array of lower bounds for rational parameter values
1560 Rational lower[SoPlexBase<R>::RATIONALPARAM_COUNT];
1561 /// array of upper bounds for rational parameter values
1562 Rational upper[SoPlexBase<R>::RATIONALPARAM_COUNT];
1563 } rationalParam;
1564 #endif
1565
1566 /// array of current boolean parameter values
1567 bool _boolParamValues[SoPlexBase<R>::BOOLPARAM_COUNT];
1568
1569 /// array of current integer parameter values
1570 int _intParamValues[SoPlexBase<R>::INTPARAM_COUNT];
1571
1572 /// array of current real parameter values
1573 Real _realParamValues[SoPlexBase<R>::REALPARAM_COUNT];
1574
1575 #ifdef SOPLEX_WITH_RATIONALPARAM
1576 /// array of current rational parameter values
1577 Rational _rationalParamValues[SoPlexBase<R>::RATIONALPARAM_COUNT];
1578 #endif
1579
1580 /// default constructor initializing default settings
1581 Settings();
1582
1583 /// copy constructor
1584 Settings(const Settings& settings);
1585
1586 /// assignment operator
1587 Settings& operator=(const Settings& settings);
1588 };
1589
1590 mutable SPxOut spxout;
1591
1592 /// returns boolean parameter value
1593 bool boolParam(const BoolParam param) const;
1594
1595 /// returns integer parameter value
1596 int intParam(const IntParam param) const;
1597
1598 /// returns real parameter value
1599 Real realParam(const RealParam param) const;
1600
1601 #ifdef SOPLEX_WITH_RATIONALPARAM
1602 /// returns rational parameter value
1603 Rational rationalParam(const RationalParam param) const;
1604 #endif
1605
1606 /// returns current parameter settings
1607 const Settings& settings() const;
1608
1609 /// returns current tolerances
1610 const std::shared_ptr<Tolerances> tolerances() const;
1611
1612 /// sets boolean parameter value; returns true on success
1613 bool setBoolParam(const BoolParam param, const bool value, const bool init = true);
1614
1615 /// sets integer parameter value; returns true on success
1616 bool setIntParam(const IntParam param, const int value, const bool init = true);
1617
1618 /// sets real parameter value; returns true on success
1619 bool setRealParam(const RealParam param, const Real value, const bool init = true);
1620
1621 #ifdef SOPLEX_WITH_RATIONALPARAM
1622 /// sets rational parameter value; returns true on success
1623 bool setRationalParam(const RationalParam param, const Rational value, const bool init = true);
1624 #endif
1625
1626 /// sets parameter settings; returns true on success
1627 bool setSettings(const Settings& newSettings, const bool init = true);
1628
1629 /// resets default parameter settings
1630 void resetSettings(const bool quiet = false, const bool init = true);
1631
1632 /// print non-default parameter values
1633 void printUserSettings();
1634
1635 /// writes settings file; returns true on success
1636 bool saveSettingsFile(const char* filename, const bool onlyChanged = false,
1637 int solvemode = 1) const;
1638
1639 /// reads settings file; returns true on success
1640 bool loadSettingsFile(const char* filename);
1641
1642 /// parses one setting string and returns true on success; note that string is modified
1643 bool parseSettingsString(char* str);
1644
1645 ///@}
1646
1647
1648 ///@name Statistics
1649 ///@{
1650
1651 /// set statistic timers to a certain type
1652 void setTimings(const Timer::TYPE ttype);
1653
1654 /// prints solution statistics
1655 void printSolutionStatistics(std::ostream& os);
1656
1657 /// prints statistics on solving process
1658 void printSolvingStatistics(std::ostream& os);
1659
1660 /// prints short statistics
1661 void printShortStatistics(std::ostream& os);
1662
1663 /// prints complete statistics
1664 void printStatistics(std::ostream& os);
1665
1666 /// prints status
1667
1668 void printStatus(std::ostream& os, typename SPxSolverBase<R>::Status status);
1669
1670 ///@}
1671
1672
1673 ///@name Miscellaneous
1674 ///@{
1675
1676 /// prints version and compilation options
1677 void printVersion() const;
1678
1679 /// checks if real LP and rational LP are in sync; dimensions will always be compared,
1680 /// vector and matrix values only if the respective parameter is set to true.
1681 /// If quiet is set to true the function will only display which vectors are different.
1682 bool areLPsInSync(const bool checkVecVals = true, const bool checkMatVals = false,
1683 const bool quiet = false) const;
1684
1685 /// set the random seeds of the solver instance
1686 void setRandomSeed(unsigned int seed);
1687
1688 /// returns the current random seed of the solver instance
1689 unsigned int randomSeed() const;
1690
1691 ///@}
1692
1693 private:
1694
1695 ///@name Statistics on solving process
1696 ///@{
1697
1698 /// class of statistics
1699 class Statistics;
1700
1701 /// statistics since last call to solveReal() or solveRational()
1702 Statistics* _statistics;
1703
1704 ///@}
1705
1706
1707 ///@name Parameter settings
1708 ///@{
1709
1710 Settings* _currentSettings;
1711
1712 std::shared_ptr<Tolerances> _tolerances;
1713
1714 Rational _rationalPosInfty;
1715 Rational _rationalNegInfty;
1716 Rational _rationalFeastol;
1717 Rational _rationalOpttol;
1718 Rational _rationalMaxscaleincr;
1719
1720 ///@}
1721
1722
1723 ///@name Data for the real LP
1724 ///@{
1725
1726 SPxSolverBase<R> _solver;
1727 SLUFactor<R> _slufactor;
1728 SPxMainSM<R> _simplifierMainSM;
1729 Presol<R> _simplifierPaPILO;
1730 SPxEquiliSC<R> _scalerUniequi;
1731 SPxEquiliSC<R> _scalerBiequi;
1732 SPxGeometSC<R> _scalerGeo1;
1733 SPxGeometSC<R> _scalerGeo8;
1734 SPxGeometSC<R> _scalerGeoequi;
1735 SPxLeastSqSC<R> _scalerLeastsq;
1736 SPxWeightST<R> _starterWeight;
1737 SPxSumST<R> _starterSum;
1738 SPxVectorST<R> _starterVector;
1739 SPxAutoPR<R> _pricerAuto;
1740 SPxDantzigPR<R> _pricerDantzig;
1741 SPxParMultPR<R> _pricerParMult;
1742 SPxDevexPR<R> _pricerDevex;
1743 SPxSteepPR<R> _pricerQuickSteep;
1744 SPxSteepExPR<R> _pricerSteep;
1745 SPxDefaultRT<R> _ratiotesterTextbook;
1746 SPxHarrisRT<R> _ratiotesterHarris;
1747 SPxFastRT<R> _ratiotesterFast;
1748 SPxBoundFlippingRT<R> _ratiotesterBoundFlipping;
1749
1750 SPxLPBase<R>*
1751 _realLP; // the real LP is also used as the original LP for the decomposition dual simplex
1752 SPxLPBase<R>* _decompLP; // used to store the original LP for the decomposition dual simplex
1753 SPxSimplifier<R>* _simplifier;
1754 SPxScaler<R>* _scaler;
1755 SPxStarter<R>* _starter;
1756
1757 #ifdef SOPLEX_WITH_BOOST
1758 #ifdef SOPLEX_WITH_MPFR
1759 //----------------------------- BOOSTED SOLVER -----------------------------
1760 // multiprecision type used for the boosted solver
1761 using BP = number<mpfr_float_backend<0>, et_off>;
1762 #else
1763 #ifdef SOPLEX_WITH_GMP
1764 using BP = number<gmp_float<50>, et_off>;
1765 #else
1766 using BP = number<cpp_dec_float<50>, et_off>;
1767 #endif
1768 #endif
1769 #else
1770 using BP = double;
1771 #endif
1772
1773 // boosted solver object
1774 SPxSolverBase<BP> _boostedSolver;
1775
1776 // ------------- Main attributes for precision boosting
1777
1778 int _initialPrecision = 50; // initial number of digits for multiprecision
1779 bool _boostingLimitReached; // true if BP::default_precision() > max authorized number of digits
1780 bool _switchedToBoosted; // true if _boostedSolver is used instead of _solver to cope with the numerical failure of _solver
1781 // this attribute remembers wether we are testing feasibility (1), unboundedness (2) or neither (0)
1782 // it is used when storing/loading the right basis in precision boosting.
1783 // example: if _certificateMode == 1, it is the basis for the feasibility LP that should be stored/loaded.
1784 int _certificateMode;
1785
1786 // ------------- Buffers for statistics of precision boosting
1787
1788 // ideally these four attributes would be local variables, however the precision boosting loop
1789 // wraps the solve in a way that it is complicated to declare these variables locally.
1790 int _lastStallPrecBoosts; // number of previous stalling precision boosts
1791 bool _factorSolNewBasisPrecBoost; // false if the current basis has already been factorized (no new iterations have been done)
1792 int _nextRatrecPrecBoost; // the iteration during or after which rational reconstruction can be performed
1793 // buffer storing the number of iterations before a given precision boost
1794 // used to detect stalling (_prevIterations < _statistics->iterations)
1795 int _prevIterations;
1796
1797 // ------------- Tolerances Ratios
1798
1799 /// ratios for computing the tolerances for precision boosting
1800 /// ratio denotes the proportion of precision used by the tolerance
1801 /// e.g. ratio = 0.65, precision = 100 digits, new tol = 10^(0.65*100)
1802 Real _tolPrecisionRatio = 0.65;
1803 Real _epsZeroPrecisionRatio = 1.0;
1804 Real _epsFactorPrecisionRatio = 1.25;
1805 Real _epsUpdatePrecisionRatio = 1.0;
1806 Real _epsPivotPrecisionRatio = 0.625;
1807
1808 // ------------- [Boosted] SLUFactor, Pricers, RatioTesters, Scalers, Simplifiers
1809
1810 SLUFactor<BP> _boostedSlufactor;
1811
1812 SPxAutoPR<BP> _boostedPricerAuto;
1813 SPxDantzigPR<BP> _boostedPricerDantzig;
1814 SPxParMultPR<BP> _boostedPricerParMult;
1815 SPxDevexPR<BP> _boostedPricerDevex;
1816 SPxSteepPR<BP> _boostedPricerQuickSteep;
1817 SPxSteepExPR<BP> _boostedPricerSteep;
1818
1819 SPxDefaultRT<BP> _boostedRatiotesterTextbook;
1820 SPxHarrisRT<BP> _boostedRatiotesterHarris;
1821 SPxFastRT<BP> _boostedRatiotesterFast;
1822 SPxBoundFlippingRT<BP> _boostedRatiotesterBoundFlipping;
1823
1824 SPxScaler<BP>* _boostedScaler;
1825 SPxSimplifier<BP>* _boostedSimplifier;
1826
1827 SPxEquiliSC<BP> _boostedScalerUniequi;
1828 SPxEquiliSC<BP> _boostedScalerBiequi;
1829 SPxGeometSC<BP> _boostedScalerGeo1;
1830 SPxGeometSC<BP> _boostedScalerGeo8;
1831 SPxGeometSC<BP> _boostedScalerGeoequi;
1832 SPxLeastSqSC<BP> _boostedScalerLeastsq;
1833
1834 SPxMainSM<BP> _boostedSimplifierMainSM;
1835 Presol<BP> _boostedSimplifierPaPILO;
1836
1837 //--------------------------------------------------------------------------
1838
1839 bool _isRealLPLoaded; // true indicates that the original LP is loaded in the _solver variable, hence all actions
1840 // are performed on the original LP.
1841 bool _isRealLPScaled;
1842 bool _applyPolishing;
1843
1844 VectorBase<R> _manualLower;
1845 VectorBase<R> _manualUpper;
1846 VectorBase<R> _manualLhs;
1847 VectorBase<R> _manualRhs;
1848 VectorBase<R> _manualObj;
1849 SPxLPBase<R> _manualRealLP;
1850
1851 ///@}
1852
1853
1854 ///@name Data for the rational LP
1855 ///@{
1856
1857 SPxLPRational* _rationalLP;
1858 SLUFactorRational _rationalLUSolver;
1859 DataArray<int> _rationalLUSolverBind;
1860
1861 LPColSetRational _slackCols;
1862 VectorRational _unboundedLower;
1863 VectorRational _unboundedUpper;
1864 VectorRational _unboundedLhs;
1865 VectorRational _unboundedRhs;
1866 DSVectorRational _tauColVector;
1867 VectorRational _feasObj;
1868 VectorRational _feasLhs;
1869 VectorRational _feasRhs;
1870 VectorRational _feasLower;
1871 VectorRational _feasUpper;
1872 VectorRational _modLower;
1873 VectorRational _modUpper;
1874 VectorRational _modLhs;
1875 VectorRational _modRhs;
1876 VectorRational _modObj;
1877 DSVectorRational _primalDualDiff;
1878 DataArray< typename SPxSolverBase<R>::VarStatus > _storedBasisStatusRows;
1879 DataArray< typename SPxSolverBase<R>::VarStatus > _storedBasisStatusCols;
1880 Array< UnitVectorRational* > _unitMatrixRational;
1881 bool _storedBasis;
1882 int _beforeLiftRows;
1883 int _beforeLiftCols;
1884
1885 /// type of bounds and sides
1886 typedef enum
1887 {
1888 /// both bounds are infinite
1889 RANGETYPE_FREE = 0,
1890
1891 /// lower bound is finite, upper bound is infinite
1892 RANGETYPE_LOWER = 1,
1893
1894 /// upper bound is finite, lower bound is infinite
1895 RANGETYPE_UPPER = 2,
1896
1897 /// lower and upper bound finite, but different
1898 RANGETYPE_BOXED = 3,
1899
1900 /// lower bound equals upper bound
1901 RANGETYPE_FIXED = 4
1902 } RangeType;
1903
1904 DataArray< RangeType > _colTypes;
1905 DataArray< RangeType > _rowTypes;
1906
1907 ///@}
1908
1909
1910 ///@name Data for the Decomposition Based Dual Simplex
1911 ///@{
1912
1913 /** row violation structure
1914 */
1915 struct RowViolation
1916 {
1917 R violation; /**< the violation of the row */
1918 int idx; /**< index of corresponding row */
1919 };
1920
1921 /** Compare class for row violations
1922 */
1923 struct RowViolationCompare
1924 {
1925 public:
1926 /** constructor
1927 */
1928 RowViolationCompare()
1929 : entry(0)
1930 {
1931 }
1932
1933 const RowViolation* entry;
1934
1935 R operator()(
1936 RowViolation i,
1937 RowViolation j
1938 ) const
1939 {
1940 return i.violation - j.violation;
1941 }
1942 };
1943
1944
1945 typedef enum
1946 {
1947 // is the original problem currently being solved.
1948 DECOMP_ORIG = 0,
1949
1950 // is the reduced problem currently being solved.
1951 DECOMP_RED = 1,
1952
1953 // is the complementary problem currently being solved.
1954 DECOMP_COMP = 2
1955 } decompStatus;
1956
1957 // the expected sign of the dual variables.
1958 enum DualSign
1959 {
1960 IS_POS = 0,
1961 IS_NEG = 1,
1962 IS_FREE = 2
1963 };
1964
1965 SPxSolverBase<R>
1966 _compSolver; // adding a solver to contain the complementary problem. It is too confusing to switch
1967 // the LP for the reduced and complementary problem in the one solver variable. The reduced
1968 // problem will be stored in _solver and the complementary problem will be stored in
1969 // _compSolver.
1970 SLUFactor<R> _compSlufactor;
1971
1972 SPxBasisBase<R>
1973 _decompTransBasis; // the basis required for the transformation to form the reduced problem
1974
1975 VectorBase<R> _transformedObj; // the objective coefficients of the transformed problem
1976 VectorBase<R> _decompFeasVector; // feasibility vector calculated using unshifted bounds.
1977 LPRowSetBase<R>
1978 _transformedRows; // a set of the original rows that have been transformed using the original basis.
1979 SPxColId _compSlackColId; // column id of the primal complementary problem slack column.
1980 SPxRowId _compSlackDualRowId; // row id in the dual of complementary problem related to the slack column.
1981 bool* _decompReducedProbRows; // flag to indicate the inclusion of a row in the reduced problem.
1982 bool* _decompReducedProbCols; // flag to indicate the inclusion of a col in the reduced problem.
1983 int* _decompRowStatus;
1984 int* _decompColStatus;
1985 int* _decompCompProbColIDsIdx; // the index to _decompPrimalColIDs for a given original col.
1986 DataArray < SPxRowId >
1987 _decompReducedProbRowIDs; // the row IDs for the related rows in the reduced problem
1988 DataArray < SPxRowId >
1989 _decompReducedProbColRowIDs;// the row IDs for the related cols in the reduced problem
1990 DataArray < SPxColId >
1991 _decompReducedProbColIDs; // the col IDs for the related cols in the reduced problem
1992 DataArray < SPxRowId > _decompPrimalRowIDs; // the primal row IDs from the original problem
1993 DataArray < SPxColId > _decompPrimalColIDs; // the primal col IDs from the original problem
1994 DataArray < SPxRowId >
1995 _decompElimPrimalRowIDs; // the primal row IDs eliminated in the complementary problem
1996 DataArray < SPxRowId >
1997 _decompDualRowIDs; // the dual row IDs from the complementary problem
1998 DataArray < SPxColId >
1999 _decompDualColIDs; // the dual col IDs from the complementary problem
2000 DataArray < SPxColId > _decompFixedVarDualIDs; // the column ids related to the fixed variables.
2001 DataArray < SPxColId >
2002 _decompVarBoundDualIDs; // the column ids related to the variable bound constraints.
2003
2004 DataArray < SPxColId >
2005 _decompCompPrimalFixedVarIDs; // the column ids related to the fixed variables in the complementary primal.
2006 DataArray < SPxColId >
2007 _decompCompPrimalVarBoundIDs; // the column ids related to the variable bound constraints in the complementary primal.
2008
2009 DataArray < SPxRowId >
2010 _decompCompPrimalRowIDs; // the primal row IDs from the complementary problem
2011 DataArray < SPxColId >
2012 _decompCompPrimalColIDs; // the primal col IDs from the complementary problem
2013
2014 int _nDecompViolBounds; // the number of violated bound constraints
2015 int _nDecompViolRows; // the number of violated rows
2016 int* _decompViolatedBounds; // the violated bounds given by the solution to the IDS reduced problem
2017 int* _decompViolatedRows; // the violated rows given by the solution to the IDS reduced problem
2018
2019
2020 int* _fixedOrigVars; // the original variables that are at their bounds in the reduced problem.
2021 // 1: fixed to upper, -1: fixed to lower, 0: unfixed.
2022 int _nPrimalRows; // the number of original problem rows included in the complementary problem
2023 int _nPrimalCols; // the number of original problem columns included in the complementary problem
2024 int _nElimPrimalRows; // the number of primal rows from the original problem eliminated from the complementary prob
2025 int _nDualRows; // the number of dual rows in the complementary problem. NOTE: _nPrimalRows = _nDualCols
2026 int _nDualCols; // the number of dual columns in the complementary problem. NOTE: _nPrimalRows = _nDualCols
2027 int _nCompPrimalRows; // the number of rows in the complementary primal problem. NOTE: _nPrimalRows = _nCompPrimalRows
2028 int _nCompPrimalCols; // the number of dual columns in the complementary problem. NOTE: _nPrimalCols = _nCompPrimalCols
2029
2030 int _decompDisplayLine; // the count for the display line
2031
2032 NameSet* _rowNames; // the row names from the input file
2033 NameSet* _colNames; // the col names from the input file
2034
2035 // Statistic information
2036 int numIncludedRows; // the number of rows currently included in the reduced problem.
2037 int numDecompIter; // the number of iterations of the decomposition dual simplex algorithm.
2038 int numRedProbIter; // the number of simplex iterations performed in the reduced problem.
2039 int numCompProbIter; // the number of iterations of the complementary problem.
2040
2041 // problem statistics
2042 int numProbRows;
2043 int numProbCols;
2044 int nNonzeros;
2045 R minAbsNonzero;
2046 R maxAbsNonzero;
2047
2048 int origCountLower;
2049 int origCountUpper;
2050 int origCountBoxed;
2051 int origCountFreeCol;
2052
2053 int origCountEqual;
2054 int origCountLhs;
2055 int origCountRhs;
2056 int origCountRanged;
2057 int origCountFreeRow;
2058
2059 decompStatus _currentProb;
2060
2061 ///@}
2062
2063
2064 ///@name Solution data
2065 ///@{
2066
2067 typename SPxSolverBase<R>::Status _status;
2068 int _lastSolveMode;
2069
2070 DataArray<typename SPxSolverBase<R>::VarStatus > _basisStatusRows;
2071 DataArray<typename SPxSolverBase<R>::VarStatus > _basisStatusCols;
2072
2073 // indicates wether an old basis is currently stored for warm start
2074 bool _hasOldBasis;
2075 bool _hasOldFeasBasis; // basis for testing feasibility
2076 bool _hasOldUnbdBasis; // basis for testing unboundedness
2077
2078 // these vectors store the last basis met in precision boosting when not testing feasibility or unboundedness.
2079 DataArray<typename SPxSolverBase<R>::VarStatus > _oldBasisStatusRows;
2080 DataArray<typename SPxSolverBase<R>::VarStatus > _oldBasisStatusCols;
2081
2082 // these vectors store the last basis met when testing feasibility in precision boosting.
2083 DataArray<typename SPxSolverBase<R>::VarStatus > _oldFeasBasisStatusRows;
2084 DataArray<typename SPxSolverBase<R>::VarStatus > _oldFeasBasisStatusCols;
2085
2086 // these vectors store the last basis met when testing unboundedness in precision boosting.
2087 DataArray<typename SPxSolverBase<R>::VarStatus > _oldUnbdBasisStatusRows;
2088 DataArray<typename SPxSolverBase<R>::VarStatus > _oldUnbdBasisStatusCols;
2089
2090 // these vectors don't replace _basisStatusRows and _basisStatusCols
2091 // they aim to overcome the issue of having the enum VarStatus inside SPxSolverBase.
2092 // When calling setBasis or getBasis (from SPxSolverBase class), a specific conversion is needed.
2093 // Function: SPxSolverBase<BP>::setBasis(...)
2094 // Usage: copy _basisStatusRows(Cols) to _tmpBasisStatusRows(Cols) before calling
2095 // mysolver.setBasis(_tmpBasisStatusRows, _tmpBasisStatusCols)
2096 // Function: SPxSolverBase<BP>::getBasis(...)
2097 // Usage: copy _tmpBasisStatusRows(Cols) to _basisStatusRows(Cols) after calling
2098 // mysolver.getBasis(_tmpBasisStatusRows, _tmpBasisStatusCols, _basisStatusRows.size(), _basisStatusCols.size())
2099 DataArray<typename SPxSolverBase<BP>::VarStatus > _tmpBasisStatusRows;
2100 DataArray<typename SPxSolverBase<BP>::VarStatus > _tmpBasisStatusCols;
2101
2102 SolBase<R> _solReal;
2103 SolRational _solRational;
2104 SolRational _workSol;
2105
2106 bool _hasBasis;
2107 bool _hasSolReal;
2108 bool _hasSolRational;
2109
2110 ///@}
2111
2112 ///@name Miscellaneous
2113 ///@{
2114
2115 int _optimizeCalls;
2116 int _unscaleCalls;
2117
2118 Rational _rationalPosone;
2119 Rational _rationalNegone;
2120 Rational _rationalZero;
2121
2122 ///@}
2123
2124 ///@name Constant helper methods
2125 ///@{
2126
2127 /// extends sparse vector to hold newmax entries if and only if it holds no more free entries
2128 void _ensureDSVectorRationalMemory(DSVectorRational& vec, const int newmax) const;
2129
2130 /// creates a permutation for removing rows/columns from an array of indices
2131 void _idxToPerm(int* idx, int idxSize, int* perm, int permSize) const;
2132
2133 /// creates a permutation for removing rows/columns from a range of indices
2134 void _rangeToPerm(int start, int end, int* perm, int permSize) const;
2135
2136 /// checks consistency for the boosted solver
2137 bool _isBoostedConsistent() const;
2138
2139 /// checks consistency
2140 bool _isConsistent() const;
2141
2142 /// should solving process be stopped?
2143 bool _isSolveStopped(bool& stoppedTime, bool& stoppedIter) const;
2144
2145 /// determines RangeType from real bounds
2146 RangeType _rangeTypeReal(const R& lower, const R& upper) const;
2147
2148 /// determines RangeType from rational bounds
2149 RangeType _rangeTypeRational(const Rational& lower, const Rational& upper) const;
2150
2151 /// switches RANGETYPE_LOWER to RANGETYPE_UPPER and vice versa
2152 RangeType _switchRangeType(const RangeType& rangeType) const;
2153
2154 /// checks whether RangeType corresponds to finite lower bound
2155 bool _lowerFinite(const RangeType& rangeType) const;
2156
2157 /// checks whether RangeType corresponds to finite upper bound
2158 bool _upperFinite(const RangeType& rangeType) const;
2159
2160 ///@}
2161
2162
2163 ///@name Non-constant helper methods
2164 ///@{
2165
2166 /// adds a single row to the real LP and adjusts basis
2167 void _addRowReal(const LPRowBase<R>& lprow);
2168
2169 /// adds a single row to the real LP and adjusts basis
2170 void _addRowReal(R lhs, const SVectorBase<R>& lprow, R rhs);
2171
2172 /// adds multiple rows to the real LP and adjusts basis
2173 void _addRowsReal(const LPRowSetBase<R>& lprowset);
2174
2175 /// adds a single column to the real LP and adjusts basis
2176 void _addColReal(const LPColReal& lpcol);
2177
2178 /// adds a single column to the real LP and adjusts basis
2179 void _addColReal(R obj, R lower, const SVectorBase<R>& lpcol, R upper);
2180
2181 /// adds multiple columns to the real LP and adjusts basis
2182 void _addColsReal(const LPColSetReal& lpcolset);
2183
2184 /// replaces row \p i with \p lprow and adjusts basis
2185 void _changeRowReal(int i, const LPRowBase<R>& lprow);
2186
2187 /// changes left-hand side vector for constraints to \p lhs and adjusts basis
2188 void _changeLhsReal(const VectorBase<R>& lhs);
2189
2190 /// changes left-hand side of row \p i to \p lhs and adjusts basis
2191 void _changeLhsReal(int i, const R& lhs);
2192
2193 /// changes right-hand side vector to \p rhs and adjusts basis
2194 void _changeRhsReal(const VectorBase<R>& rhs);
2195
2196 /// changes right-hand side of row \p i to \p rhs and adjusts basis
2197 void _changeRhsReal(int i, const R& rhs);
2198
2199 /// changes left- and right-hand side vectors and adjusts basis
2200 void _changeRangeReal(const VectorBase<R>& lhs, const VectorBase<R>& rhs);
2201
2202 /// changes left- and right-hand side of row \p i and adjusts basis
2203 void _changeRangeReal(int i, const R& lhs, const R& rhs);
2204
2205 /// replaces column \p i with \p lpcol and adjusts basis
2206 void _changeColReal(int i, const LPColReal& lpcol);
2207
2208 /// changes vector of lower bounds to \p lower and adjusts basis
2209 void _changeLowerReal(const VectorBase<R>& lower);
2210
2211 /// changes lower bound of column i to \p lower and adjusts basis
2212 void _changeLowerReal(int i, const R& lower);
2213
2214 /// changes vector of upper bounds to \p upper and adjusts basis
2215 void _changeUpperReal(const VectorBase<R>& upper);
2216
2217 /// changes \p i 'th upper bound to \p upper and adjusts basis
2218 void _changeUpperReal(int i, const R& upper);
2219
2220 /// changes vectors of column bounds to \p lower and \p upper and adjusts basis
2221 void _changeBoundsReal(const VectorBase<R>& lower, const VectorBase<R>& upper);
2222
2223 /// changes bounds of column \p i to \p lower and \p upper and adjusts basis
2224 void _changeBoundsReal(int i, const R& lower, const R& upper);
2225
2226 /// changes matrix entry in row \p i and column \p j to \p val and adjusts basis
2227 void _changeElementReal(int i, int j, const R& val);
2228
2229 /// removes row \p i and adjusts basis
2230 void _removeRowReal(int i);
2231
2232 /// removes all rows with an index \p i such that \p perm[i] < 0; upon completion, \p perm[i] >= 0 indicates the
2233 /// new index where row \p i has been moved to; note that \p perm must point to an array of size at least
2234 /// #numRows()
2235 void _removeRowsReal(int perm[]);
2236
2237 /// remove all rows with indices in array \p idx of size \p n; an array \p perm of size #numRows() may be passed
2238 /// as buffer memory
2239 void _removeRowsReal(int idx[], int n, int perm[]);
2240
2241 /// removes rows \p start to \p end including both; an array \p perm of size #numRows() may be passed as buffer
2242 /// memory
2243 void _removeRowRangeReal(int start, int end, int perm[]);
2244
2245 /// removes column i
2246 void _removeColReal(int i);
2247
2248 /// removes all columns with an index \p i such that \p perm[i] < 0; upon completion, \p perm[i] >= 0 indicates the
2249 /// new index where column \p i has been moved to; note that \p perm must point to an array of size at least
2250 /// #numColsReal()
2251 void _removeColsReal(int perm[]);
2252
2253 /// remove all columns with indices in array \p idx of size \p n; an array \p perm of size #numColsReal() may be
2254 /// passed as buffer memory
2255 void _removeColsReal(int idx[], int n, int perm[]);
2256
2257 /// removes columns \p start to \p end including both; an array \p perm of size #numColsReal() may be passed as
2258 /// buffer memory
2259 void _removeColRangeReal(int start, int end, int perm[]);
2260
2261 /// invalidates solution
2262 void _invalidateSolution();
2263
2264 /// enables simplifier and scaler according to current parameters
2265 void _enableSimplifierAndScaler();
2266
2267 /// disables simplifier and scaler
2268 void _disableSimplifierAndScaler();
2269
2270 /// ensures that the rational LP is available; performs no sync
2271 void _ensureRationalLP();
2272
2273 /// ensures that the real LP and the basis are loaded in the solver; performs no sync
2274 void _ensureRealLPLoaded();
2275
2276 /// call floating-point solver and update statistics on iterations etc.
2277 void _solveBoostedRealLPAndRecordStatistics(volatile bool* interrupt = NULL);
2278
2279 /// call floating-point solver and update statistics on iterations etc.
2280 void _solveRealLPAndRecordStatistics(volatile bool* interrupt = NULL);
2281
2282 /// reads real LP in LP or MPS format from file and returns true on success; gets row names, column names, and
2283 /// integer variables if desired
2284 bool _readFileReal(const char* filename, NameSet* rowNames = 0, NameSet* colNames = 0,
2285 DIdxSet* intVars = 0);
2286
2287 /// reads rational LP in LP or MPS format from file and returns true on success; gets row names, column names, and
2288 /// integer variables if desired
2289 bool _readFileRational(const char* filename, NameSet* rowNames = 0, NameSet* colNames = 0,
2290 DIdxSet* intVars = 0);
2291
2292 /// completes range type arrays after adding columns and/or rows
2293 void _completeRangeTypesRational();
2294
2295 /// recomputes range types from scratch using real LP
2296 void _recomputeRangeTypesReal();
2297
2298 /// recomputes range types from scratch using rational LP
2299 void _recomputeRangeTypesRational();
2300
2301 /// synchronizes real LP with rational LP, i.e., copies (rounded) rational LP into real LP, without looking at the sync mode
2302 void _syncLPReal(bool time = true);
2303
2304 /// synchronizes rational LP with real LP, i.e., copies real LP to rational LP, without looking at the sync mode
2305 void _syncLPRational(bool time = true);
2306
2307 /// synchronizes rational solution with real solution, i.e., copies (rounded) rational solution to real solution
2308 void _syncRealSolution();
2309
2310 /// synchronizes real solution with rational solution, i.e., copies real solution to rational solution
2311 void _syncRationalSolution();
2312
2313 /// returns pointer to a constant unit vector available until destruction of the SoPlexBase class
2314 const UnitVectorRational* _unitVectorRational(const int i);
2315
2316 /// parses one line in a settings file and returns true on success; note that the string is modified
2317 bool _parseSettingsLine(char* line, const int lineNumber);
2318
2319 ///@}
2320
2321
2322 //**@name Private solving methods implemented in solverational.hpp */
2323 ///@{
2324
2325 /// stores floating-point solution of original LP as current rational solution and ensure that solution vectors have right dimension; ensure that solution is aligned with basis
2326 template <typename T>
2327 void _storeRealSolutionAsRational(
2328 SolRational& sol,
2329 VectorBase<T>& primalReal,
2330 VectorBase<T>& dualReal,
2331 int& dualSize);
2332
2333 /// computes violation of bounds during the refinement loop
2334 void _computeBoundsViolation(SolRational& sol, Rational& boundsViolation);
2335
2336 /// computes violation of sides during the refinement loop
2337 void _computeSidesViolation(SolRational& sol, Rational& sideViolation);
2338
2339 /// computes violation of reduced costs during the refinement loop
2340 void _computeReducedCostViolation(
2341 SolRational& sol,
2342 Rational& redCostViolation,
2343 const bool& maximizing);
2344
2345 /// computes dual violation during the refinement loop
2346 void _computeDualViolation(
2347 SolRational& sol,
2348 Rational& dualViolation,
2349 const bool& maximizing);
2350
2351 /// checks termination criteria for refinement loop
2352 bool _isRefinementOver(
2353 bool& primalFeasible,
2354 bool& dualFeasible,
2355 Rational& boundsViolation,
2356 Rational& sideViolation,
2357 Rational& redCostViolation,
2358 Rational& dualViolation,
2359 int minIRRoundsRemaining,
2360 bool& stoppedTime,
2361 bool& stoppedIter,
2362 int numFailedRefinements);
2363
2364 /// checks refinement loop progress
2365 void _checkRefinementProgress(
2366 Rational& boundsViolation,
2367 Rational& sideViolation,
2368 Rational& redCostViolation,
2369 Rational& dualViolation,
2370 Rational& maxViolation,
2371 Rational& bestViolation,
2372 const Rational& violationImprovementFactor,
2373 int& numFailedRefinements);
2374
2375 /// performs rational reconstruction and/or factorizationd
2376 void _ratrecAndOrRatfac(
2377 int& minIRRoundsRemaining,
2378 int& lastStallIterations,
2379 int& numberOfIterations,
2380 bool& factorSolNewBasis,
2381 int& nextRatrec,
2382 const Rational& errorCorrectionFactor,
2383 Rational& errorCorrection,
2384 Rational& maxViolation,
2385 SolRational& sol,
2386 bool& primalFeasible,
2387 bool& dualFeasible,
2388 bool& stoppedTime,
2389 bool& stoppedIter,
2390 bool& error,
2391 bool& breakAfter,
2392 bool& continueAfter);
2393
2394 /// forces value of given nonbasic variable to bound
2395 void _forceNonbasicToBound(
2396 SolRational& sol,
2397 int& c,
2398 const int& maxDimRational,
2399 bool toLower);
2400
2401 /// computes primal scaling factor; limit increase in scaling by tolerance used in floating point solve
2402 void _computePrimalScalingFactor(
2403 Rational& maxScale,
2404 Rational& primalScale,
2405 Rational& boundsViolation,
2406 Rational& sideViolation,
2407 Rational& redCostViolation);
2408
2409 /// computes dual scaling factor; limit increase in scaling by tolerance used in floating point solve
2410 void _computeDualScalingFactor(
2411 Rational& maxScale,
2412 Rational& primalScale,
2413 Rational& dualScale,
2414 Rational& redCostViolation,
2415 Rational& dualViolation);
2416
2417 /// applies scaled bounds
2418 template <typename T>
2419 void _applyScaledBounds(SPxSolverBase<T>& solver, Rational& primalScale);
2420
2421 /// applies scaled sides
2422 template <typename T>
2423 void _applyScaledSides(SPxSolverBase<T>& solver, Rational& primalScale);
2424
2425 /// applies scaled objective function
2426 template <typename T>
2427 void _applyScaledObj(SPxSolverBase<T>& solver, Rational& dualScale, SolRational& sol);
2428
2429 /// evaluates result of solve. Return true if the algorithm must to stopped, false otherwise.
2430 template <typename T>
2431 bool _evaluateResult(
2432 SPxSolverBase<T>& solver,
2433 typename SPxSolverBase<T>::Status result,
2434 bool usingRefinedLP,
2435 SolRational& sol,
2436 VectorBase<T>& dualReal,
2437 bool& infeasible,
2438 bool& unbounded,
2439 bool& stoppedTime,
2440 bool& stoppedIter,
2441 bool& error);
2442
2443 /// corrects primal solution and aligns with basis
2444 template <typename T>
2445 void _correctPrimalSolution(
2446 SolRational& sol,
2447 Rational& primalScale,
2448 int& primalSize,
2449 const int& maxDimRational,
2450 VectorBase<T>& primalReal);
2451
2452 /// updates or recomputes slacks depending on which looks faster
2453 void _updateSlacks(SolRational& sol, int& primalSize);
2454
2455 /// corrects dual solution and aligns with basis
2456 template <typename T>
2457 void _correctDualSolution(
2458 SPxSolverBase<T>& solver,
2459 SolRational& sol,
2460 const bool& maximizing,
2461 VectorBase<T>& dualReal,
2462 Rational& dualScale,
2463 int& dualSize,
2464 const int& maxDimRational);
2465
2466 /// updates or recomputes reduced cost values depending on which looks faster; adding one to the length of the
2467 /// dual vector accounts for the objective function vector
2468 void _updateReducedCosts(SolRational& sol, int& dualSize, const int& numCorrectedPrimals);
2469
2470 ///@todo precision-boosting move some place else
2471 /// converts the given DataArray of VarStatus to boostedPrecision
2472 void _convertDataArrayVarStatusToBoosted(
2473 DataArray< typename SPxSolverBase<R>::VarStatus >& base,
2474 DataArray< typename SPxSolverBase<BP>::VarStatus >& copy);
2475
2476 ///@todo precision-boosting move some place else
2477 /// converts the given DataArray of VarStatus to R precision
2478 void _convertDataArrayVarStatusToRPrecision(
2479 DataArray< typename SPxSolverBase<BP>::VarStatus >& base,
2480 DataArray< typename SPxSolverBase<R>::VarStatus >& copy);
2481
2482 /// disable initial precision solver and switch to boosted solver
2483 void _switchToBoosted();
2484
2485 /// setup boosted solver before launching iteration
2486 void _setupBoostedSolver();
2487
2488 /// increase the multiprecision, return false if maximum precision is reached, true otherwise
2489 bool _boostPrecision();
2490
2491 /// reset the boosted precision to the default value
2492 void _resetBoostedPrecision();
2493
2494 /// setup recovery mecanism using multiprecision, return false if maximum precision reached, true otherwise
2495 bool _setupBoostedSolverAfterRecovery();
2496
2497 /// return true if slack basis has to be loaded for boosted solver
2498 bool _isBoostedStartingFromSlack(bool initialSolve = true);
2499
2500 /// indicate if we are testing feasibility, unboundedness or neither
2501 void _switchToStandardMode();
2502 void _switchToFeasMode();
2503 void _switchToUnbdMode();
2504
2505 /// check if we are testing feasibility, unboundedness or neither
2506 bool _inStandardMode();
2507 bool _inFeasMode();
2508 bool _inUnbdMode();
2509
2510 // stores given basis in old basis attributes: _oldBasisStatusRows, _oldFeasBasisStatusRows, _oldUnbdBasisStatusRows (and ...Cols)
2511 void _storeBasisAsOldBasis(DataArray< typename SPxSolverBase<R>::VarStatus >& rows,
2512 DataArray< typename SPxSolverBase<R>::VarStatus >& cols);
2513
2514 // stores given basis in old basis attributes: _oldBasisStatusRows, _oldFeasBasisStatusRows, _oldUnbdBasisStatusRows (and ...Cols)
2515 void _storeBasisAsOldBasisBoosted(DataArray< typename SPxSolverBase<BP>::VarStatus >& rows,
2516 DataArray< typename SPxSolverBase<BP>::VarStatus >& cols);
2517
2518 // get the last advanced and stable basis stored by the initial solver and store it as old basis, unsimplify basis if simplifier activated
2519 void _storeLastStableBasis(bool vanished);
2520
2521 // get the last advanced and stable basis stored by the boosted solver and store it as old basis, unsimplify basis if simplifier activated
2522 void _storeLastStableBasisBoosted(bool vanished);
2523
2524 // load old basis in solver. The old basis loaded depends on the certificate mode (feasibility, unboundedness, or neither)
2525 bool _loadBasisFromOldBasis(bool boosted);
2526
2527 // update statistics for precision boosting
2528 void _updateBoostingStatistics();
2529
2530 /// solves current problem using multiprecision floating-point solver
2531 /// return false if a new boosted iteration is necessary, true otherwise
2532 void _solveRealForRationalBoostedStable(
2533 SolRational& sol,
2534 bool& primalFeasible,
2535 bool& dualFeasible,
2536 bool& infeasible,
2537 bool& unbounded,
2538 bool& stoppedTime,
2539 bool& stoppedIter,
2540 bool& error,
2541 bool& needNewBoostedIt);
2542
2543 /// solves current problem with iterative refinement and recovery mechanism using boosted solver
2544 void _performOptIRStableBoosted(
2545 SolRational& sol,
2546 bool acceptUnbounded,
2547 bool acceptInfeasible,
2548 int minIRRoundsRemaining,
2549 bool& primalFeasible,
2550 bool& dualFeasible,
2551 bool& infeasible,
2552 bool& unbounded,
2553 bool& stoppedTime,
2554 bool& stoppedIter,
2555 bool& error,
2556 bool& needNewBoostedIt);
2557
2558 /// perform iterative refinement using the right precision
2559 void _performOptIRWrapper(
2560 SolRational& sol,
2561 bool acceptUnbounded,
2562 bool acceptInfeasible,
2563 int minIRRoundsRemaining,
2564 bool& primalFeasible,
2565 bool& dualFeasible,
2566 bool& infeasible,
2567 bool& unbounded,
2568 bool& stoppedTime,
2569 bool& stoppedIter,
2570 bool& error
2571 );
2572
2573 /// solves current problem using double floating-point solver
2574 void _solveRealForRationalStable(
2575 SolRational& sol,
2576 bool& primalFeasible,
2577 bool& dualFeasible,
2578 bool& infeasible,
2579 bool& unbounded,
2580 bool& stoppedTime,
2581 bool& stoppedIter,
2582 bool& error);
2583
2584 /// solves current problem with iterative refinement and recovery mechanism
2585 void _performOptIRStable(SolRational& sol,
2586 bool acceptUnbounded,
2587 bool acceptInfeasible,
2588 int minIRRoundsRemaining,
2589 bool& primalFeasible,
2590 bool& dualFeasible,
2591 bool& infeasible,
2592 bool& unbounded,
2593 bool& stoppedTime,
2594 bool& stoppedIter,
2595 bool& error);
2596
2597 /// performs iterative refinement on the auxiliary problem for testing unboundedness
2598 void _performUnboundedIRStable(SolRational& sol, bool& hasUnboundedRay, bool& stoppedTime,
2599 bool& stoppedIter, bool& error);
2600
2601 /// performs iterative refinement on the auxiliary problem for testing feasibility
2602 void _performFeasIRStable(SolRational& sol, bool& withDualFarkas, bool& stoppedTime,
2603 bool& stoppedIter, bool& error);
2604
2605 /// reduces matrix coefficient in absolute value by the lifting procedure of Thiele et al. 2013
2606 void _lift();
2607
2608 /// undoes lifting
2609 void _project(SolRational& sol);
2610
2611 /// store basis
2612 void _storeBasis();
2613
2614 /// restore basis
2615 void _restoreBasis();
2616
2617 /// stores objective, bounds, and sides of real LP
2618 void _storeLPReal();
2619
2620 /// restores objective, bounds, and sides of real LP
2621 void _restoreLPReal();
2622
2623 /// introduces slack variables to transform inequality constraints into equations for both rational and real LP,
2624 /// which should be in sync
2625 void _transformEquality();
2626
2627 /// undoes transformation to equality form
2628 void _untransformEquality(SolRational& sol);
2629
2630 /// transforms LP to unboundedness problem by moving the objective function to the constraints, changing right-hand
2631 /// side and bounds to zero, and adding an auxiliary variable for the decrease in the objective function
2632 void _transformUnbounded();
2633
2634 /// undoes transformation to unboundedness problem
2635 void _untransformUnbounded(SolRational& sol, bool unbounded);
2636
2637 /// transforms LP to feasibility problem by removing the objective function, shifting variables, and homogenizing the
2638 /// right-hand side
2639 void _transformFeasibility();
2640
2641 /// undoes transformation to feasibility problem
2642 void _untransformFeasibility(SolRational& sol, bool infeasible);
2643
2644 /** computes radius of infeasibility box implied by an approximate Farkas' proof
2645
2646 Given constraints of the form \f$ lhs <= Ax <= rhs \f$, a farkas proof y should satisfy \f$ y^T A = 0 \f$ and
2647 \f$ y_+^T lhs - y_-^T rhs > 0 \f$, where \f$ y_+, y_- \f$ denote the positive and negative parts of \f$ y \f$.
2648 If \f$ y \f$ is approximate, it may not satisfy \f$ y^T A = 0 \f$ exactly, but the proof is still valid as long
2649 as the following holds for all potentially feasible \f$ x \f$:
2650
2651 \f[
2652 y^T Ax < (y_+^T lhs - y_-^T rhs) (*)
2653 \f]
2654
2655 we may therefore calculate \f$ y^T A \f$ and \f$ y_+^T lhs - y_-^T rhs \f$ exactly and check if the upper and lower
2656 bounds on \f$ x \f$ imply that all feasible \f$ x \f$ satisfy (*), and if not then compute bounds on \f$ x \f$ to
2657 guarantee (*). The simplest way to do this is to compute
2658
2659 \f[
2660 B = (y_+^T lhs - y_-^T rhs) / \sum_i(|(y^T A)_i|)
2661 \f]
2662
2663 noting that if every component of \f$ x \f$ has \f$ |x_i| < B \f$, then (*) holds.
2664
2665 \f$ B \f$ can be increased by iteratively including variable bounds smaller than \f$ B \f$. The speed of this
2666 method can be further improved by using interval arithmetic for all computations. For related information see
2667 Sec. 4 of Neumaier and Shcherbina, Mathematical Programming A, 2004.
2668
2669 Set transformed to true if this method is called after _transformFeasibility().
2670 */
2671 void _computeInfeasBox(SolRational& sol, bool transformed);
2672
2673 /// solves real LP during iterative refinement
2674 typename SPxSolverBase<R>::Status _solveRealForRational(bool fromscratch, VectorBase<R>& primal,
2675 VectorBase<R>& dual,
2676 DataArray< typename SPxSolverBase<R>::VarStatus >& basisStatusRows,
2677 DataArray< typename SPxSolverBase<R>::VarStatus >& basisStatusCols);
2678
2679 /// solves real LP with recovery mechanism
2680 typename SPxSolverBase<R>::Status _solveRealStable(bool acceptUnbounded, bool acceptInfeasible,
2681 VectorBase<R>& primal, VectorBase<R>& dual,
2682 DataArray< typename SPxSolverBase<R>::VarStatus >& basisStatusRows,
2683 DataArray< typename SPxSolverBase<R>::VarStatus >& basisStatusCols,
2684 const bool forceNoSimplifier = false);
2685
2686 /// solves real LP during iterative refinement
2687 void _solveRealForRationalBoosted(
2688 VectorBase<BP>& primal, VectorBase<BP>& dual,
2689 DataArray< typename SPxSolverBase<R>::VarStatus >& basisStatusRows,
2690 DataArray< typename SPxSolverBase<R>::VarStatus >& basisStatusCols,
2691 typename SPxSolverBase<BP>::Status& boostedResult, bool initialSolve);
2692
2693 /// computes rational inverse of basis matrix as defined by _rationalLUSolverBind
2694 void _computeBasisInverseRational();
2695
2696 /// factorizes rational basis matrix in column representation
2697 void _factorizeColumnRational(SolRational& sol,
2698 DataArray< typename SPxSolverBase<R>::VarStatus >& basisStatusRows,
2699 DataArray< typename SPxSolverBase<R>::VarStatus >& basisStatusCols, bool& stoppedTime,
2700 bool& stoppedIter, bool& error, bool& optimal);
2701
2702 /// attempts rational reconstruction of primal-dual solution
2703 bool _reconstructSolutionRational(SolRational& sol,
2704 DataArray< typename SPxSolverBase<R>::VarStatus >& basisStatusRows,
2705 DataArray< typename SPxSolverBase<R>::VarStatus >& basisStatusCols,
2706 const Rational& denomBoundSquared);
2707 ///@}
2708
2709 ///@name Private solving methods implemented in solvereal.cpp
2710 ///@{
2711
2712 /// solves the templated LP
2713 void _optimize(volatile bool* interrupt = NULL);
2714
2715 /// temporary fix for Rational
2716 void _optimizeRational(volatile bool* interrupt = NULL);
2717
2718 /// checks result of the solving process and solves again without preprocessing if necessary
2719 void _evaluateSolutionReal(typename SPxSimplifier<R>::Result simplificationStatus);
2720
2721 /// solves real LP with/without preprocessing
2722 void _preprocessAndSolveReal(bool applyPreprocessing, volatile bool* interrupt = NULL);
2723
2724 /// loads original problem into solver and solves again after it has been solved to optimality with preprocessing
2725 void _resolveWithoutPreprocessing(typename SPxSimplifier<R>::Result simplificationStatus);
2726
2727 /// verify computed solution and resolve if necessary
2728 void _verifySolutionReal();
2729
2730 /// verify computed obj stop and resolve if necessary
2731 void _verifyObjLimitReal();
2732
2733 /// stores solution of the real LP; before calling this, the real LP must be loaded in the solver and solved (again)
2734 void _storeSolutionReal(bool verify = true);
2735
2736 /// stores solution from the simplifier because problem vanished in presolving step
2737 void _storeSolutionRealFromPresol();
2738
2739 /// unscales stored solution to remove internal or external scaling of LP
2740 void _unscaleSolutionReal(SPxLPBase<R>& LP, bool persistent = true);
2741
2742 /// load original LP and possibly setup a slack basis
2743 void _loadRealLP(bool initBasis);
2744
2745 /// check scaling of LP
2746 void _checkScaling(SPxLPBase<R>* origLP) const;
2747
2748 /// check correctness of (un)scaled basis matrix operations
2749 void _checkBasisScaling();
2750
2751 /// check whether persistent scaling is supposed to be reapplied again after unscaling
2752 bool _reapplyPersistentScaling() const;
2753
2754 /// solves LP using the decomposition based dual simplex
2755 void _solveDecompositionDualSimplex();
2756
2757 /// creating copies of the original problem that will be manipulated to form the reduced and complementary problems
2758 void _createDecompReducedAndComplementaryProblems();
2759
2760 /// forms the reduced problem
2761 void _formDecompReducedProblem(bool& stop);
2762
2763 /// solves the reduced problem
2764 void _solveDecompReducedProblem();
2765
2766 /// forms the complementary problem
2767 void _formDecompComplementaryProblem();
2768
2769 /// simplifies the problem and solves
2770 void _decompSimplifyAndSolve(SPxSolverBase<R>& solver, SLUFactor<R>& sluFactor, bool fromScratch,
2771 bool applyPreprocessing);
2772
2773 /// loads original problem into solver and solves again after it has been solved to optimality with preprocessing
2774 void _decompResolveWithoutPreprocessing(SPxSolverBase<R>& solver, SLUFactor<R>& sluFactor,
2775 typename SPxSimplifier<R>::Result result);
2776
2777 /// identifies the columns of the row-form basis that correspond to rows with zero dual multipliers.
2778 void _getZeroDualMultiplierIndices(VectorBase<R> feasVector, int* nonposind, int* colsforremoval,
2779 int* nnonposind, bool& stop);
2780
2781 /// retrieves the compatible columns from the constraint matrix
2782 void _getCompatibleColumns(VectorBase<R> feasVector, int* nonposind, int* compatind,
2783 int* rowsforremoval, int* colsforremoval,
2784 int nnonposind, int* ncompatind, bool formRedProb, bool& stop);
2785
2786 /// computes the reduced problem objective coefficients
2787 void _computeReducedProbObjCoeff(bool& stop);
2788
2789 /// computes the compatible bound constraints and adds them to the reduced problem
2790 void _getCompatibleBoundCons(LPRowSetBase<R>& boundcons, int* compatboundcons, int* nonposind,
2791 int* ncompatboundcons,
2792 int nnonposind, bool& stop);
2793
2794 /// computes the rows to remove from the complementary problem
2795 void _getRowsForRemovalComplementaryProblem(int* nonposind, int* bind, int* rowsforremoval,
2796 int* nrowsforremoval,
2797 int nnonposind);
2798
2799 /// removing rows from the complementary problem.
2800 void _deleteAndUpdateRowsComplementaryProblem(SPxRowId rangedRowIds[], int& naddedrows);
2801
2802 /// evaluates the solution of the reduced problem for the DBDS
2803 void _evaluateSolutionDecomp(SPxSolverBase<R>& solver, SLUFactor<R>& sluFactor,
2804 typename SPxSimplifier<R>::Result result);
2805
2806 /// update the reduced problem with additional columns and rows
2807 void _updateDecompReducedProblem(R objVal, VectorBase<R> dualVector, VectorBase<R> redcostVector,
2808 VectorBase<R> compPrimalVector,
2809 VectorBase<R> compDualVector);
2810
2811 /// update the reduced problem with additional columns and rows based upon the violated original bounds and rows
2812 void _updateDecompReducedProblemViol(bool allrows);
2813
2814 /// builds the update rows with those violated in the complmentary problem
2815 void _findViolatedRows(R compObjValue, Array<RowViolation>& violatedrows, int& nviolatedrows);
2816
2817 /// update the dual complementary problem with additional columns and rows
2818 void _updateDecompComplementaryDualProblem(bool origObj);
2819
2820 /// update the primal complementary problem with additional columns and rows
2821 void _updateDecompComplementaryPrimalProblem(bool origObj);
2822
2823 /// checking the optimality of the original problem.
2824 void _checkOriginalProblemOptimality(VectorBase<R> primalVector, bool printViol);
2825
2826 /// updating the slack column coefficients to adjust for equality constraints
2827 void _updateComplementaryDualSlackColCoeff();
2828
2829 /// updating the slack column coefficients to adjust for equality constraints
2830 void _updateComplementaryPrimalSlackColCoeff();
2831
2832 /// removing the dual columns related to the fixed variables
2833 void _removeComplementaryDualFixedPrimalVars(int* currFixedVars);
2834
2835 /// removing the dual columns related to the fixed variables
2836 void _identifyComplementaryDualFixedPrimalVars(int* currFixedVars);
2837
2838 /// updating the dual columns related to the fixed primal variables.
2839 void _updateComplementaryDualFixedPrimalVars(int* currFixedVars);
2840
2841 /// removing the dual columns related to the fixed variables
2842 void _identifyComplementaryPrimalFixedPrimalVars(int* currFixedVars);
2843
2844 /// updating the dual columns related to the fixed primal variables.
2845 void _updateComplementaryPrimalFixedPrimalVars(int* currFixedVars);
2846
2847 /// updating the complementary dual problem with the original objective function
2848 void _setComplementaryDualOriginalObjective();
2849
2850 /// updating the complementary primal problem with the original objective function
2851 void _setComplementaryPrimalOriginalObjective();
2852
2853 /// determining which bound the primal variables will be fixed to.
2854 int getOrigVarFixedDirection(int colNum);
2855
2856 /// checks the dual feasibility of the current basis
2857 bool checkBasisDualFeasibility(VectorBase<R> feasVec);
2858
2859 /// returns the expected sign of the dual variables for the reduced problem
2860 DualSign getExpectedDualVariableSign(int rowNumber);
2861
2862 /// returns the expected sign of the dual variables for the original problem
2863 DualSign getOrigProbDualVariableSign(int rowNumber);
2864
2865 /// prints a display line of the flying table for the DBDS
2866 void printDecompDisplayLine(SPxSolverBase<R>& solver, const SPxOut::Verbosity origVerb, bool force,
2867 bool forceHead);
2868
2869 /// stores the problem statistics of the original problem
2870 void getOriginalProblemStatistics();
2871
2872 /// stores the problem statistics of the original problem
2873 void printOriginalProblemStatistics(std::ostream& os);
2874
2875 /// gets the coefficient of the slack variable in the primal complementary problem
2876 R getCompSlackVarCoeff(int primalRowNum);
2877
2878 /// gets violation of bounds; returns true on success
2879 bool getDecompBoundViolation(R& maxviol, R& sumviol);
2880
2881 /// gets violation of constraints; returns true on success
2882 bool getDecompRowViolation(R& maxviol, R& sumviol);
2883
2884 /// function call to terminate the decomposition simplex
2885 bool decompTerminate(R timeLimit);
2886
2887 /// function to build a basis for the original problem as given by the solution to the reduced problem
2888 void _writeOriginalProblemBasis(const char* filename, NameSet* rowNames, NameSet* colNames,
2889 bool cpxFormat);
2890
2891 /// function to retrieve the original problem row basis status from the reduced and complementary problems
2892 void getOriginalProblemBasisRowStatus(DataArray< int >& degenerateRowNums,
2893 DataArray< typename SPxSolverBase<R>::VarStatus >& degenerateRowStatus, int& nDegenerateRows,
2894 int& nNonBasicRows);
2895
2896 /// function to retrieve the column status for the original problem basis from the reduced and complementary problems
2897 void getOriginalProblemBasisColStatus(int& nNonBasicCols);
2898
2899 ///@}
2900 };
2901
2902 /* Backwards compatibility */
2903 typedef SoPlexBase<Real> SoPlex;
2904 // A header file containing all the general templated functions
2905
2906 } // namespace soplex
2907
2908 // General templated function
2909 #include "soplex/solvedbds.hpp"
2910 #include "soplex.hpp"
2911 #include "soplex/solverational.hpp"
2912 #include "soplex/testsoplex.hpp"
2913 #include "soplex/solvereal.hpp"
2914
2915 #endif // _SOPLEX_H_
2916