SoPlex Documentation
Loading...
Searching...
No Matches
spxpapilo.h
Go to the documentation of this file.
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
26#ifndef SOPLEX_WITH_PAPILO
27
28namespace soplex
29{
30
31template <class R> class Presol : public SPxSimplifier<R>
32{
33
34public:
35
36 //------------------------------------
37 ///@name Constructors / destructors
38 ///@{
39 /// default constructor.
41 { ; };
42
43 /// copy constructor.
44 Presol(const Presol& old) : SPxSimplifier<R>(old) { ; }
45
46 /// assignment operator
47 Presol& operator=(const Presol& rhs)
48 {
49 return *this;
50 }
51
52 /// destructor.
53 virtual ~Presol()
54 {
55 ;
56 }
57
59 {
60 return new Presol(*this);
61 }
62
63 virtual typename SPxSimplifier<R>::Result simplify(SPxLPBase<R>& lp, R eps, R delta,
65 {
67 }
68
72 {
73 assert(false);
75 };
76
77 virtual void unsimplify(const VectorBase<R>&, const VectorBase<R>&,
78 const VectorBase<R>&, const VectorBase<R>&,
79 const typename SPxSolverBase<R>::VarStatus[],
80 const typename SPxSolverBase<R>::VarStatus[],
81 bool isOptimal)
82 {
83 assert(false);
84 };
85
86 /// returns result status of the simplification
87 virtual typename SPxSimplifier<R>::Result result() const
88 {
89 assert(false);
91 }
92
93 /// specifies whether an optimal solution has already been unsimplified.
94 virtual bool isUnsimplified() const
95 {
96 assert(false);
97 return false;
98 }
99
100 /// returns a reference to the unsimplified primal solution.
102 {
103 assert(false);
104 static const VectorBase<R>& emptyVector = VectorBase<R>();
105 return emptyVector;
106 }
107
108 /// returns a reference to the unsimplified dual solution.
110 {
111 assert(false);
112 static const VectorBase<R>& emptyVector = VectorBase<R>();
113 return emptyVector;
114 }
115
116 /// returns a reference to the unsimplified slack values.
118 {
119 assert(false);
120 static const VectorBase<R>& emptyVector = VectorBase<R>();
121 return emptyVector;
122 }
123
124 /// returns a reference to the unsimplified reduced costs.
126 {
127 assert(false);
128 static const VectorBase<R>& emptyVector = VectorBase<R>();
129 return emptyVector;
130 }
131
132 /// gets basis status for a single row.
134 {
135 assert(false);
137 }
138
139 /// gets basis status for a single column.
141 {
142 assert(false);
144 }
145
146 /// get optimal basis.
147 virtual void getBasis(typename SPxSolverBase<R>::VarStatus rows[],
148 typename SPxSolverBase<R>::VarStatus cols[], const int rowsSize = -1,
149 const int colsSize = -1) const
150 {
151 assert(false);
152 }
153
154};
155
156}
157
158
159#else
160
161#include <memory>
162
163#include "papilo/core/Presolve.hpp"
164#include "papilo/core/ProblemBuilder.hpp"
165#include "papilo/Config.hpp"
166
167#include "soplex/spxsimplifier.h"
168
169#include "soplex/spxdefines.h"
170#include "soplex/spxsimplifier.h"
171#include "soplex/array.h"
172#include "soplex/exceptions.h"
173#include "soplex/spxdefines.h"
174
175
176namespace soplex
177{
178
179template<class R>
180class Presol : public SPxSimplifier<R>
181{
182private:
183
184#ifdef SOPLEX_DEBUG
185 const papilo::VerbosityLevel verbosityLevel = papilo::VerbosityLevel::kInfo;
186#else
187 const papilo::VerbosityLevel verbosityLevel = papilo::VerbosityLevel::kQuiet;
188#endif
189
190 VectorBase<R> m_prim; ///< unsimplified primal solution VectorBase<R>.
191 VectorBase<R> m_slack; ///< unsimplified slack VectorBase<R>.
192 VectorBase<R> m_dual; ///< unsimplified dual solution VectorBase<R>.
193 VectorBase<R> m_redCost; ///< unsimplified reduced cost VectorBase<R>.
194 DataArray<typename SPxSolverBase<R>::VarStatus> m_cBasisStat; ///< basis status of columns.
195 DataArray<typename SPxSolverBase<R>::VarStatus> m_rBasisStat; ///< basis status of rows.
196
197
198 papilo::PostsolveStorage<R>
199 postsolveStorage; ///< stored postsolve to recalculate the original solution
200 bool noChanges = false; ///< did PaPILO reduce the problem?
201
202 bool postsolved; ///< was the solution already postsolve?
203 bool vanished = false;
204 R m_epsilon; ///< epsilon zero.
205 R m_feastol; ///< primal feasibility tolerance.
206 R m_opttol; ///< dual feasibility tolerance.
207 R modifyRowsFac; ///<
208 DataArray<int> m_stat; ///< preprocessing history.
209 typename SPxLPBase<R>::SPxSense m_thesense; ///< optimization sense.
210
211 // TODO: the following parameters were ignored? Maybe I don't exactly know what they suppose to be
212 bool m_keepbounds; ///< keep some bounds (for boundflipping)
213 typename SPxSimplifier<R>::Result m_result; ///< result of the simplification.
214
215protected:
216
217 R epsZero() const
218 {
219 return m_epsilon;
220 }
221
222 R feastol() const
223 {
224 return m_feastol;
225 }
226
227 R opttol() const
228 {
229 return m_opttol;
230 }
231
232public:
233
234 //------------------------------------
235 ///@name Constructors / destructors
236 ///@{
237 /// default constructor.
238 explicit Presol(Timer::TYPE ttype = Timer::USER_TIME)
239 : SPxSimplifier<R>("PaPILO", ttype), postsolved(false), m_epsilon(DEFAULT_EPS_ZERO),
240 m_feastol(DEFAULT_BND_VIOL), m_opttol(DEFAULT_BND_VIOL), modifyRowsFac(1.0),
241 m_thesense(SPxLPBase<R>::MAXIMIZE),
242 m_keepbounds(false), m_result(this->OKAY)
243 { ; };
244
245 /// copy constructor.
246 Presol(const Presol& old)
247 : SPxSimplifier<R>(old), m_prim(old.m_prim), m_slack(old.m_slack), m_dual(old.m_dual),
248 m_redCost(old.m_redCost), m_cBasisStat(old.m_cBasisStat), m_rBasisStat(old.m_rBasisStat),
249 postsolveStorage(old.postsolveStorage), postsolved(old.postsolved), m_epsilon(old.m_epsilon),
250 m_feastol(old.m_feastol), m_opttol(old.m_opttol),
251 modifyRowsFac(old.modifyRowsFac), m_thesense(old.m_thesense),
252 m_keepbounds(old.m_keepbounds), m_result(old.m_result)
253 {
254 ;
255 }
256
257 /// assignment operator
258 Presol& operator=(const Presol& rhs)
259 {
260 if(this != &rhs)
261 {
263 m_prim = rhs.m_prim;
264 m_slack = rhs.m_slack;
265 m_dual = rhs.m_dual;
266 m_redCost = rhs.m_redCost;
267 m_cBasisStat = rhs.m_cBasisStat;
268 m_rBasisStat = rhs.m_rBasisStat;
269 postsolved = rhs.postsolved;
270 m_epsilon = rhs.m_epsilon;
271 m_feastol = rhs.m_feastol;
272 m_opttol = rhs.m_opttol;
273 m_thesense = rhs.m_thesense;
274 m_keepbounds = rhs.m_keepbounds;
275 m_result = rhs.m_result;
276 postsolveStorage = rhs.postsolveStorage;
277 modifyRowsFac = rhs.modifyRowsFac;
278 }
279
280 return *this;
281 }
282
283 /// destructor.
284 virtual ~Presol()
285 {
286 ;
287 }
288
289 SPxSimplifier<R>* clone() const
290 {
291 return new Presol(*this);
292 }
293
294 void
295 setModifyConsFrac(R value)
296 {
297 modifyRowsFac = value;
298 }
299
300 virtual typename SPxSimplifier<R>::Result simplify(SPxLPBase<R>& lp, R eps, R delta,
301 Real remainingTime)
302 {
303 return simplify(lp, eps, delta, delta, remainingTime, false, 0);
304 }
305
306 virtual typename SPxSimplifier<R>::Result simplify(SPxLPBase<R>& lp, R eps, R ftol, R otol,
307 Real remainingTime,
308 bool keepbounds, uint32_t seed);
309
310 virtual void unsimplify(const VectorBase<R>&, const VectorBase<R>&, const VectorBase<R>&,
311 const VectorBase<R>&,
312 const typename SPxSolverBase<R>::VarStatus[],
313 const typename SPxSolverBase<R>::VarStatus[],
314 bool isOptimal);
315
316 /// returns result status of the simplification
317 virtual typename SPxSimplifier<R>::Result result() const
318 {
319 return m_result;
320 }
321
322 /// specifies whether an optimal solution has already been unsimplified.
323 virtual bool isUnsimplified() const
324 {
325 return postsolved;
326 }
327
328 /// returns a reference to the unsimplified primal solution.
329 virtual const VectorBase<R>& unsimplifiedPrimal()
330 {
331 assert(postsolved);
332 return m_prim;
333 }
334
335 /// returns a reference to the unsimplified dual solution.
336 virtual const VectorBase<R>& unsimplifiedDual()
337 {
338 assert(postsolved);
339 return m_dual;
340 }
341
342 /// returns a reference to the unsimplified slack values.
343 virtual const VectorBase<R>& unsimplifiedSlacks()
344 {
345 assert(postsolved);
346 return m_slack;
347 }
348
349 /// returns a reference to the unsimplified reduced costs.
350 virtual const VectorBase<R>& unsimplifiedRedCost()
351 {
352 assert(postsolved);
353 return m_redCost;
354 }
355
356 /// gets basis status for a single row.
357 virtual typename SPxSolverBase<R>::VarStatus getBasisRowStatus(int i) const
358 {
359 assert(postsolved);
360 return m_rBasisStat[i];
361 }
362
363 /// gets basis status for a single column.
364 virtual typename SPxSolverBase<R>::VarStatus getBasisColStatus(int j) const
365 {
366 assert(postsolved);
367 return m_cBasisStat[j];
368 }
369
370 /// get optimal basis.
371 virtual void getBasis(typename SPxSolverBase<R>::VarStatus rows[],
372 typename SPxSolverBase<R>::VarStatus cols[], const int rowsSize = -1,
373 const int colsSize = -1) const
374 {
375 assert(postsolved);
376 assert(rowsSize < 0 || rowsSize >= m_rBasisStat.size());
377 assert(colsSize < 0 || colsSize >= m_cBasisStat.size());
378
379 for(int i = 0; i < m_rBasisStat.size(); ++i)
380 rows[i] = m_rBasisStat[i];
381
382 for(int j = 0; j < m_cBasisStat.size(); ++j)
383 cols[j] = m_cBasisStat[j];
384 }
385
386private:
387
388 void initLocalVariables(const SPxLPBase <R>& lp);
389
390 void configurePapilo(papilo::Presolve<R>& presolve, R feasTolerance, R epsilon, uint32_t seed,
391 Real remainingTime) const;
392
393 void applyPresolveResultsToColumns(SPxLPBase <R>& lp, const papilo::Problem<R>& problem,
394 const papilo::PresolveResult<R>& res) const;
395
396 void applyPresolveResultsToRows(SPxLPBase <R>& lp, const papilo::Problem<R>& problem,
397 const papilo::PresolveResult<R>& res) const;
398
399 papilo::VarBasisStatus
400 convertToPapiloStatus(typename SPxSolverBase<R>::VarStatus status) const;
401
403 convertToSoplexStatus(papilo::VarBasisStatus status) const ;
404};
405
406template <class R>
407void Presol<R>::unsimplify(const VectorBase<R>& x, const VectorBase<R>& y,
408 const VectorBase<R>& s, const VectorBase<R>& r,
409 const typename SPxSolverBase<R>::VarStatus rows[],
410 const typename SPxSolverBase<R>::VarStatus cols[],
411 bool isOptimal)
412{
413
414 MSG_INFO1((*this->spxout), (*this->spxout)
415 << " --- unsimplifying solution and basis"
416 << std::endl;
417 )
418
419 assert(x.dim() <= m_prim.dim());
420 assert(y.dim() <= m_dual.dim());
421 assert(x.dim() == r.dim());
422 assert(y.dim() == s.dim());
423
424 //if presolving made no changes then copy the reduced solution to the original
425 if(noChanges)
426 {
427 for(int j = 0; j < x.dim(); ++j)
428 {
429 m_prim[j] = x[j];
430 m_redCost[j] = r[j];
431 m_cBasisStat[j] = cols[j];
432 }
433
434 for(int i = 0; i < y.dim(); ++i)
435 {
436 m_dual[i] = y[i];
437 m_slack[i] = s[i];
438 m_rBasisStat[i] = rows[i];
439 }
440
441 postsolved = true;
442 return;
443 }
444
445 int nColsReduced = (int)postsolveStorage.origcol_mapping.size();
446 int nRowsReduced = (int)postsolveStorage.origrow_mapping.size();
447 assert(x.dim() == (int)postsolveStorage.origcol_mapping.size() || vanished);
448 assert(y.dim() == (int)postsolveStorage.origrow_mapping.size() || vanished);
449
450 papilo::Solution<R> originalSolution{};
451 papilo::Solution<R> reducedSolution{};
452 reducedSolution.type = papilo::SolutionType::kPrimalDual;
453 reducedSolution.basisAvailabe = true;
454
455 reducedSolution.primal.clear();
456 reducedSolution.reducedCosts.clear();
457 reducedSolution.varBasisStatus.clear();
458 reducedSolution.dual.clear();
459 reducedSolution.rowBasisStatus.clear();
460
461 reducedSolution.primal.resize(nColsReduced);
462 reducedSolution.reducedCosts.resize(nColsReduced);
463 reducedSolution.varBasisStatus.resize(nColsReduced);
464 reducedSolution.dual.resize(nRowsReduced);
465 reducedSolution.rowBasisStatus.resize(nRowsReduced);
466
467 postsolved = true;
468
469 // NOTE: for maximization problems, we have to switch signs of dual and
470 // reduced cost values, since simplifier assumes minimization problem
471 R switch_sign = m_thesense == SPxLPBase<R>::MAXIMIZE ? -1 : 1;
472
473 // assign values of variables in reduced LP
474 for(int j = 0; j < nColsReduced; ++j)
475 {
476 reducedSolution.primal[j] = isZero(x[j], this->epsZero()) ? 0.0 : x[j];
477 reducedSolution.reducedCosts[j] =
478 isZero(r[j], this->epsZero()) ? 0.0 : switch_sign * r[j];
479 reducedSolution.varBasisStatus[j] = convertToPapiloStatus(cols[j]);
480 }
481
482 for(int i = 0; i < nRowsReduced; ++i)
483 {
484 reducedSolution.dual[i] = isZero(y[i], this->epsZero()) ? 0.0 : switch_sign * y[i];
485 reducedSolution.rowBasisStatus[i] = convertToPapiloStatus(rows[i]);
486 }
487
488 papilo::Num<R> num {};
489 num.setEpsilon(m_epsilon);
490 num.setFeasTol(m_feastol);
491 /* since PaPILO verbosity is quiet it's irrelevant what the messenger is */
492 papilo::Message msg{};
493 msg.setVerbosityLevel(verbosityLevel);
494
495 papilo::Postsolve<R> postsolve {msg, num};
497
498 if(status == PostsolveStatus::kFailed && isOptimal)
499 {
500 MSG_ERROR(std::cerr << "PaPILO did not pass validation" << std::endl;)
501 assert(false);
502 }
503
504 for(int j = 0; j < (int)postsolveStorage.nColsOriginal; ++j)
505 {
506 m_prim[j] = originalSolution.primal[j];
507 m_redCost[j] = switch_sign * originalSolution.reducedCosts[j];
508 m_cBasisStat[j] = convertToSoplexStatus(originalSolution.varBasisStatus[j]);
509 }
510
511 for(int i = 0; i < (int)postsolveStorage.nRowsOriginal; ++i)
512 {
513 m_dual[i] = switch_sign * originalSolution.dual[i];
514 m_slack[i] = originalSolution.slack[i];
515 m_rBasisStat[i] = convertToSoplexStatus(originalSolution.rowBasisStatus[i]);
516 }
517
518}
519
520template <class R>
521papilo::VarBasisStatus
522Presol<R>::convertToPapiloStatus(const typename SPxSolverBase<R>::VarStatus status) const
523{
524 switch(status)
525 {
527 return papilo::VarBasisStatus::ON_UPPER;
528
530 return papilo::VarBasisStatus::ON_LOWER;
531
533 return papilo::VarBasisStatus::FIXED;
534
536 return papilo::VarBasisStatus::BASIC;
537
539 return papilo::VarBasisStatus::UNDEFINED;
540
542 return papilo::VarBasisStatus::ZERO;
543 }
544
545 return papilo::VarBasisStatus::UNDEFINED;
546}
547
548template <class R>
549typename SPxSolverBase<R>::VarStatus
550Presol<R>::convertToSoplexStatus(papilo::VarBasisStatus status) const
551{
552 switch(status)
553 {
554 case papilo::VarBasisStatus::ON_UPPER:
556
557 case papilo::VarBasisStatus::ON_LOWER:
559
560 case papilo::VarBasisStatus::ZERO:
562
563 case papilo::VarBasisStatus::FIXED:
565
566 case papilo::VarBasisStatus::UNDEFINED:
568
569 case papilo::VarBasisStatus::BASIC:
571 }
572
574}
575
576
577template<class R>
578papilo::Problem<R> buildProblem(SPxLPBase<R>& lp)
579{
580 papilo::ProblemBuilder<R> builder;
581
582 /* build problem from matrix */
583 int nnz = lp.nNzos();
584 int ncols = lp.nCols();
585 int nrows = lp.nRows();
586 builder.reserve(nnz, nrows, ncols);
587
588 /* set up columns */
589 builder.setNumCols(ncols);
590
591 R switch_sign = lp.spxSense() == SPxLPBase<R>::MAXIMIZE ? -1 : 1;
592
593 for(int i = 0; i < ncols; ++i)
594 {
595 R lowerbound = lp.lower(i);
596 R upperbound = lp.upper(i);
597 R objective = lp.obj(i);
598 builder.setColLb(i, lowerbound);
599 builder.setColUb(i, upperbound);
600 builder.setColLbInf(i, lowerbound <= -R(infinity));
601 builder.setColUbInf(i, upperbound >= R(infinity));
602
603 builder.setColIntegral(i, false);
604 builder.setObj(i, objective * switch_sign);
605 }
606
607 /* set up rows */
608 builder.setNumRows(nrows);
609
610 for(int i = 0; i < nrows; ++i)
611 {
612 const SVectorBase<R> rowVector = lp.rowVector(i);
613 int rowlength = rowVector.size();
614 int* indices = new int[rowlength];
615 R* rowValues = new R[rowlength];
616
617 for(int j = 0; j < rowlength; j++)
618 {
619 const Nonzero<R> element = rowVector.element(j);
620 indices[j] = element.idx;
621 rowValues[j] = element.val;
622 }
623
624 builder.addRowEntries(i, rowlength, indices, rowValues);
625
626 R lhs = lp.lhs(i);
627 R rhs = lp.rhs(i);
628 builder.setRowLhs(i, lhs);
629 builder.setRowRhs(i, rhs);
630 builder.setRowLhsInf(i, lhs <= -R(infinity));
631 builder.setRowRhsInf(i, rhs >= R(infinity));
632 }
633
634 return builder.build();
635}
636
637
638template<class R>
639typename SPxSimplifier<R>::Result
640Presol<R>::simplify(SPxLPBase<R>& lp, R eps, R ftol, R otol,
641 Real remainingTime, bool keepbounds, uint32_t seed)
642{
643
644 //TODO: how to use the keepbounds parameter?
645 m_keepbounds = keepbounds;
646
647 if(m_keepbounds)
648 MSG_WARNING((*this->spxout), (*this->spxout) << "==== PaPILO doesn't handle parameter keepbounds" <<
649 std::endl;)
650
652
653 papilo::Problem<R> problem = buildProblem(lp);
654 papilo::Presolve<R> presolve;
655
657 MSG_INFO1((*this->spxout), (*this->spxout)
658 << " --- starting PaPILO" << std::endl;
659 )
660
661 papilo::PresolveResult<R> res = presolve.apply(problem);
662
663 switch(res.status)
664 {
665 case papilo::PresolveStatus::kInfeasible:
667 MSG_INFO1((*this->spxout), (*this->spxout)
668 << " --- presolving detected infeasibility" << std::endl;
669 )
671
672 case papilo::PresolveStatus::kUnbndOrInfeas:
673 case papilo::PresolveStatus::kUnbounded:
675 MSG_INFO1((*this->spxout), (*this->spxout) <<
676 "==== Presolving detected unboundedness of the problem" << std::endl;
677 )
679
680 case papilo::PresolveStatus::kUnchanged:
681 // since Soplex has no state unchanged store the value in a new variable
682 noChanges = true;
683 MSG_INFO1((*this->spxout), (*this->spxout)
684 << "==== Presolving found nothing " << std::endl;
685 )
687
688 case papilo::PresolveStatus::kReduced:
689 break;
690 }
691
692
693 int newNonzeros = problem.getConstraintMatrix().getNnz();
694
695 if(newNonzeros == 0 || ((problem.getNRows() <= modifyRowsFac * lp.nRows() ||
696 newNonzeros <= modifyRowsFac * lp.nNzos())))
697 {
698 MSG_INFO1((*this->spxout), (*this->spxout)
699 << " --- presolved problem has " << problem.getNRows() <<
700 " rows, "
701 << problem.getNCols() << " cols and "
702 << newNonzeros << " non-zeros and "
703 << presolve.getStatistics().nboundchgs << " boundchanges and "
704 << presolve.getStatistics().nsidechgs << " sidechanges"
705 << std::endl;
706 )
707 postsolveStorage = res.postsolve;
708
709 // remove all constraints and variables
710 for(int j = lp.nCols() - 1; j >= 0; j--)
711 lp.removeCol(j);
712
713 for(int i = lp.nRows() - 1; i >= 0; i--)
714 lp.removeRow(i);
715
718 assert(newNonzeros == lp.nNzos());
719 }
720 else
721 {
722 noChanges = true;
723 MSG_INFO1((*this->spxout),
724 (*this->spxout)
725
726 << " --- presolve results smaller than the modifyconsfac"
727 << std::endl;
728 )
729 }
730
731 if(newNonzeros == 0)
732 {
733 vanished = true;
735 }
736
737 return m_result;
738}
739
740template<class R>
741void Presol<R>::initLocalVariables(const SPxLPBase <R>& lp)
742{
743 m_result = SPxSimplifier<R>::OKAY;
744
745 m_thesense = lp.spxSense();
746 postsolved = false;
747
748 m_prim.reDim(lp.nCols());
749 m_slack.reDim(lp.nRows());
750 m_dual.reDim(lp.nRows());
751 m_redCost.reDim(lp.nCols());
752 m_cBasisStat.reSize(lp.nCols());
753 m_rBasisStat.reSize(lp.nRows());
754
755 this->m_timeUsed->reset();
756 this->m_timeUsed->start();
757}
758
759template<class R>
760void Presol<R>::configurePapilo(papilo::Presolve<R>& presolve, R feasTolerance, R epsilon,
761 uint32_t seed, Real remainingTime) const
762{
763 /* communicate the SOPLEX parameters to the presolve libary */
764
765 /* communicate the random seed */
766 presolve.getPresolveOptions().randomseed = (unsigned int) seed;
767
768 /* set number of threads to be used for presolve */
769 /* TODO: set threads for PaPILO? Can Soplex be run with multiple threads?*/
770 // presolve.getPresolveOptions().threads = data->threads;
771
772 presolve.getPresolveOptions().tlim = remainingTime;
773 presolve.getPresolveOptions().feastol = double(feasTolerance);
774 presolve.getPresolveOptions().epsilon = double(epsilon);
775 presolve.getPresolveOptions().detectlindep = 0;
776 presolve.getPresolveOptions().componentsmaxint = -1;
777 presolve.getPresolveOptions().calculate_basis_for_dual = true;
778
779 presolve.setVerbosityLevel(verbosityLevel);
780
781 /* enable lp presolvers with dual postsolve*/
782 using uptr = std::unique_ptr<papilo::PresolveMethod<R>>;
783
784 /* fast presolvers*/
785 presolve.addPresolveMethod(uptr(new papilo::SingletonCols<R>()));
786 presolve.addPresolveMethod(uptr(new papilo::ConstraintPropagation<R>()));
787
788 /* medium presolver */
789 presolve.addPresolveMethod(uptr(new papilo::ParallelRowDetection<R>()));
790 presolve.addPresolveMethod(uptr(new papilo::ParallelColDetection<R>()));
791 presolve.addPresolveMethod(uptr(new papilo::SingletonStuffing<R>()));
792 presolve.addPresolveMethod(uptr(new papilo::DualFix<R>()));
793 presolve.addPresolveMethod(uptr(new papilo::FixContinuous<R>()));
794
795 /* exhaustive presolvers*/
796 presolve.addPresolveMethod(uptr(new papilo::DominatedCols<R>()));
797
798 /**
799 * TODO: PaPILO doesn't support dualpostsolve for those presolvers
800 * presolve.addPresolveMethod(uptr(new papilo::SimpleSubstitution<R>()));
801 * presolve.addPresolveMethod(uptr(new papilo::DualInfer<R>()));
802 * presolve.addPresolveMethod(uptr(new papilo::Substitution<R>()));
803 * presolve.addPresolveMethod(uptr(new papilo::Sparsify<R>()));
804 * presolve.getPresolveOptions().removeslackvars = false;
805 * presolve.getPresolveOptions().maxfillinpersubstitution
806 * =data->maxfillinpersubstitution;
807 * presolve.getPresolveOptions().maxshiftperrow = data->maxshiftperrow;
808 */
809
810
811}
812
813template<class R>
814void Presol<R>::applyPresolveResultsToColumns(SPxLPBase <R>& lp, const papilo::Problem<R>& problem,
815 const papilo::PresolveResult<R>& res) const
816{
817
818 const papilo::Objective<R>& objective = problem.getObjective();
819 const papilo::Vec<R>& upperBounds = problem.getUpperBounds();
820 const papilo::Vec<R>& lowerBounds = problem.getLowerBounds();
821 const papilo::Vec<papilo::ColFlags>& colFlags = problem.getColFlags();
822
823 R switch_sign = lp.spxSense() == SPxLPBase<R>::MAXIMIZE ? -1 : 1;
824
825 for(int col = 0; col < problem.getNCols(); col++)
826 {
828 R lb = lowerBounds[col];
829
830 if(colFlags[col].test(papilo::ColFlag::kLbInf))
831 lb = -R(infinity);
832
833 R ub = upperBounds[col];
834
835 if(colFlags[col].test(papilo::ColFlag::kUbInf))
836 ub = R(infinity);
837
839 lp.addCol(column);
840 assert(lp.lower(col) == lb);
841 assert(lp.upper(col) == ub);
842 }
843
844 lp.changeObjOffset(objective.offset);
845
846 assert(problem.getNCols() == lp.nCols());
847}
848
849template<class R>
850void Presol<R>::applyPresolveResultsToRows(SPxLPBase <R>& lp, const papilo::Problem<R>& problem,
851 const papilo::PresolveResult<R>& res) const
852{
853 int size = res.postsolve.origrow_mapping.size();
854
855 //add the adjusted constraints
856 for(int row = 0; row < size; row++)
857 {
858 R rhs = problem.getConstraintMatrix().getRightHandSides()[row];
859
860 if(problem.getRowFlags()[row].test(papilo::RowFlag::kRhsInf))
861 rhs = R(infinity);
862
863 R lhs = problem.getConstraintMatrix().getLeftHandSides()[row];
864
865 if(problem.getRowFlags()[row].test(papilo::RowFlag::kLhsInf))
866 lhs = -R(infinity);
867
868 const papilo::SparseVectorView<R> papiloRowVector =
869 problem.getConstraintMatrix().getRowCoefficients(row);
870 const int* indices = papiloRowVector.getIndices();
871 const R* values = papiloRowVector.getValues();
872
873 int length = papiloRowVector.getLength();
875
876 for(int i = 0; i < length; i++)
877 {
878 soplexRowVector.add(indices[i], values[i]);
879 }
880
882 lp.addRow(lpRowBase);
883 assert(lp.lhs(row) == lhs);
884 assert(lp.rhs(row) == rhs);
885 }
886
887 assert(problem.getNRows() == lp.nRows());
888}
889
890} // namespace soplex
891
892#endif
Save arrays of arbitrary types.
Safe arrays of data objects.
Definition dataarray.h:75
DataArray(const DataArray &old)
copy constructor
Definition dataarray.h:328
int size() const
return nr. of elements.
Definition dataarray.h:227
virtual SPxSimplifier< R >::Result simplify(SPxLPBase< R > &lp, R eps, R ftol, R otol, Real remainingTime, bool keepbounds, uint32_t seed)
Definition spxpapilo.h:69
virtual void unsimplify(const VectorBase< R > &, const VectorBase< R > &, const VectorBase< R > &, const VectorBase< R > &, const typename SPxSolverBase< R >::VarStatus[], const typename SPxSolverBase< R >::VarStatus[], bool isOptimal)
Definition spxpapilo.h:77
Presol(Timer::TYPE ttype=Timer::USER_TIME)
Definition spxpapilo.h:40
Presol & operator=(const Presol &rhs)
assignment operator
Definition spxpapilo.h:47
virtual const VectorBase< R > & unsimplifiedSlacks()
returns a reference to the unsimplified slack values.
Definition spxpapilo.h:117
virtual const VectorBase< R > & unsimplifiedRedCost()
returns a reference to the unsimplified reduced costs.
Definition spxpapilo.h:125
SPxSimplifier< R > * clone() const
Definition spxpapilo.h:58
virtual ~Presol()
destructor.
Definition spxpapilo.h:53
virtual SPxSolverBase< R >::VarStatus getBasisColStatus(int j) const
gets basis status for a single column.
Definition spxpapilo.h:140
virtual SPxSimplifier< R >::Result simplify(SPxLPBase< R > &lp, R eps, R delta, Real remainingTime)
Definition spxpapilo.h:63
virtual SPxSimplifier< R >::Result result() const
returns result status of the simplification
Definition spxpapilo.h:87
virtual const VectorBase< R > & unsimplifiedDual()
returns a reference to the unsimplified dual solution.
Definition spxpapilo.h:109
virtual SPxSolverBase< R >::VarStatus getBasisRowStatus(int i) const
gets basis status for a single row.
Definition spxpapilo.h:133
virtual void getBasis(typename SPxSolverBase< R >::VarStatus rows[], typename SPxSolverBase< R >::VarStatus cols[], const int rowsSize=-1, const int colsSize=-1) const
get optimal basis.
Definition spxpapilo.h:147
Presol(const Presol &old)
copy constructor.
Definition spxpapilo.h:44
virtual const VectorBase< R > & unsimplifiedPrimal()
returns a reference to the unsimplified primal solution.
Definition spxpapilo.h:101
virtual bool isUnsimplified() const
specifies whether an optimal solution has already been unsimplified.
Definition spxpapilo.h:94
SPxSense
Optimization sense.
Definition spxlpbase.h:125
LP simplification abstract base class.
Result
Result of the simplification.
@ INFEASIBLE
primal infeasibility was detected
@ UNBOUNDED
primal unboundedness was detected
@ OKAY
simplification could be done
@ VANISHED
the problem was so much simplified that it vanished
SPxSimplifier & operator=(const SPxSimplifier &rhs)
assignment operator
SPxSimplifier(const char *p_name, Timer::TYPE ttype=Timer::USER_TIME)
constructor
Sequential object-oriented SimPlex.
Definition spxsolver.h:104
@ BASIC
variable is basic.
Definition spxsolver.h:213
@ ON_LOWER
variable set to its lower bound.
Definition spxsolver.h:210
@ ON_UPPER
variable set to its upper bound.
Definition spxsolver.h:209
@ UNDEFINED
nothing known about basis status (possibly due to a singular basis in transformed problem)
Definition spxsolver.h:214
@ FIXED
variable fixed to identical bounds.
Definition spxsolver.h:211
@ ZERO
free variable fixed to zero.
Definition spxsolver.h:212
TYPE
types of timers
Definition timer.h:109
Exception classes for SoPlex.
Everything should be within this namespace.
THREADLOCAL const Real infinity
double Real
Definition spxdefines.h:266
Debugging, floating point type and parameter definitions.
#define MSG_WARNING(spxout, x)
Prints out message x if the verbosity level is at least SPxOut::WARNING.
Definition spxdefines.h:164
#define DEFAULT_BND_VIOL
default allowed bound violation
Definition spxdefines.h:274
#define DEFAULT_EPS_ZERO
default allowed additive zero: 1.0 + EPS_ZERO == 1.0
Definition spxdefines.h:278
#define MSG_INFO1(spxout, x)
Prints out message x if the verbosity level is at least SPxOut::INFO1.
Definition spxdefines.h:166
#define MSG_ERROR(x)
Prints out message x if the verbosity level is at least SPxOut::ERROR.
Definition spxdefines.h:162
LP simplification base class.