Storm 1.10.0.1
A Modern Probabilistic Model Checker
Loading...
Searching...
No Matches
IterativeMinMaxLinearEquationSolver.cpp
Go to the documentation of this file.
1#include <functional>
2#include <limits>
3#include <type_traits>
4
9
12
29
30namespace storm {
31namespace solver {
32
33template<typename ValueType, typename SolutionType>
35 std::unique_ptr<LinearEquationSolverFactory<SolutionType>>&& linearEquationSolverFactory)
36 : linearEquationSolverFactory(std::move(linearEquationSolverFactory)) {
37 // Intentionally left empty
38}
39
40template<typename ValueType, typename SolutionType>
42 storm::storage::SparseMatrix<ValueType> const& A, std::unique_ptr<LinearEquationSolverFactory<SolutionType>>&& linearEquationSolverFactory)
43 : StandardMinMaxLinearEquationSolver<ValueType, SolutionType>(A), linearEquationSolverFactory(std::move(linearEquationSolverFactory)) {
44 // Intentionally left empty.
45}
46
47template<typename ValueType, typename SolutionType>
49 storm::storage::SparseMatrix<ValueType>&& A, std::unique_ptr<LinearEquationSolverFactory<SolutionType>>&& linearEquationSolverFactory)
50 : StandardMinMaxLinearEquationSolver<ValueType, SolutionType>(std::move(A)), linearEquationSolverFactory(std::move(linearEquationSolverFactory)) {
51 // Intentionally left empty.
52}
53
54template<typename ValueType, typename SolutionType>
55MinMaxMethod IterativeMinMaxLinearEquationSolver<ValueType, SolutionType>::getMethod(Environment const& env, bool isExactMode) const {
56 // Adjust the method if none was specified and we want exact or sound computations.
57 auto method = env.solver().minMax().getMethod();
58
59 if (isExactMode && method != MinMaxMethod::PolicyIteration && method != MinMaxMethod::RationalSearch && method != MinMaxMethod::ViToPi) {
60 if (env.solver().minMax().isMethodSetFromDefault()) {
62 "Selecting 'Policy iteration' as the solution technique to guarantee exact results. If you want to override this, please explicitly specify a "
63 "different method.");
64 method = MinMaxMethod::PolicyIteration;
65 } else {
66 STORM_LOG_WARN("The selected solution method " << toString(method) << " does not guarantee exact results.");
67 }
68 } else if (env.solver().isForceSoundness() && method != MinMaxMethod::SoundValueIteration && method != MinMaxMethod::IntervalIteration &&
69 method != MinMaxMethod::PolicyIteration && method != MinMaxMethod::RationalSearch && method != MinMaxMethod::OptimisticValueIteration &&
70 method != MinMaxMethod::GuessingValueIteration) {
71 if (env.solver().minMax().isMethodSetFromDefault()) {
72 method = MinMaxMethod::OptimisticValueIteration;
74 "Selecting '"
75 << toString(method)
76 << "' as the solution technique to guarantee sound results. If you want to override this, please explicitly specify a different method.");
77 } else {
78 STORM_LOG_WARN("The selected solution method does not guarantee sound results.");
79 }
80 }
81 STORM_LOG_THROW(method == MinMaxMethod::ValueIteration || method == MinMaxMethod::PolicyIteration || method == MinMaxMethod::RationalSearch ||
82 method == MinMaxMethod::SoundValueIteration || method == MinMaxMethod::IntervalIteration ||
83 method == MinMaxMethod::OptimisticValueIteration || method == MinMaxMethod::GuessingValueIteration || method == MinMaxMethod::ViToPi,
84 storm::exceptions::InvalidEnvironmentException, "This solver does not support the selected method '" << toString(method) << "'.");
85 return method;
86}
87
88template<typename ValueType, typename SolutionType>
90 std::vector<SolutionType>& x, std::vector<ValueType> const& b) const {
91 bool result = false;
92 switch (getMethod(env, storm::NumberTraits<ValueType>::IsExact || env.solver().isForceExact())) {
93 case MinMaxMethod::ValueIteration:
94 result = solveEquationsValueIteration(env, dir, x, b);
95 break;
96 case MinMaxMethod::OptimisticValueIteration:
97 result = solveEquationsOptimisticValueIteration(env, dir, x, b);
98 break;
99 case MinMaxMethod::GuessingValueIteration:
100 result = solveEquationsGuessingValueIteration(env, dir, x, b);
101 break;
102 case MinMaxMethod::PolicyIteration:
103 result = solveEquationsPolicyIteration(env, dir, x, b);
104 break;
105 case MinMaxMethod::RationalSearch:
106 result = solveEquationsRationalSearch(env, dir, x, b);
107 break;
108 case MinMaxMethod::IntervalIteration:
109 result = solveEquationsIntervalIteration(env, dir, x, b);
110 break;
111 case MinMaxMethod::SoundValueIteration:
112 result = solveEquationsSoundValueIteration(env, dir, x, b);
113 break;
114 case MinMaxMethod::ViToPi:
115 result = solveEquationsViToPi(env, dir, x, b);
116 break;
117 default:
118 STORM_LOG_THROW(false, storm::exceptions::InvalidEnvironmentException, "This solver does not implement the selected solution method");
119 }
120
121 return result;
122}
123
124template<typename ValueType, typename SolutionType>
126 if (!viOperatorTriv && !viOperatorNontriv) {
127 if (this->A->hasTrivialRowGrouping()) {
128 // The trivial row grouping minmax operator makes sense over intervals.
129 viOperatorTriv = std::make_shared<helper::ValueIterationOperator<ValueType, true, SolutionType>>();
130 viOperatorTriv->setMatrixBackwards(*this->A);
131 if constexpr (!std::is_same_v<ValueType, storm::Interval>) {
132 // It might be that someone is using a minmaxlinearequationsolver with an advanced VI algorithm
133 // but is just passing a DTMC over doubles. In this case we need to populate this VI operator.
134 // It behaves exactly the same as the trivial row grouping operator, but it is currently hardcoded
135 // to be used by, e.g., optimistic value iteration.
136 viOperatorNontriv = std::make_shared<helper::ValueIterationOperator<ValueType, false, SolutionType>>();
137 viOperatorNontriv->setMatrixBackwards(*this->A);
138 }
139 } else {
140 // The nontrivial row grouping minmax operator makes sense for MDPs.
141 viOperatorNontriv = std::make_shared<helper::ValueIterationOperator<ValueType, false, SolutionType>>();
142 viOperatorNontriv->setMatrixBackwards(*this->A);
143 }
144 }
145 if (this->choiceFixedForRowGroup) {
146 // Ignore those rows that are not selected
147 assert(this->initialScheduler);
148 auto callback = [&](uint64_t groupIndex, uint64_t localRowIndex) {
149 return this->choiceFixedForRowGroup->get(groupIndex) && this->initialScheduler->at(groupIndex) != localRowIndex;
150 };
151 if (viOperatorTriv) {
152 viOperatorTriv->setIgnoredRows(true, callback);
153 }
154 if (viOperatorNontriv) {
155 viOperatorNontriv->setIgnoredRows(true, callback);
156 }
157 }
158}
159
160template<typename ValueType, typename SolutionType>
161void IterativeMinMaxLinearEquationSolver<ValueType, SolutionType>::extractScheduler(std::vector<SolutionType>& x, std::vector<ValueType> const& b,
162 OptimizationDirection const& dir, bool updateX, bool robust) const {
163 if (std::is_same_v<ValueType, storm::Interval> && this->A->hasTrivialRowGrouping()) {
164 // Create robust scheduler index if it doesn't exist yet
165 if (!this->robustSchedulerIndex) {
166 this->robustSchedulerIndex = std::vector<uint64_t>(x.size(), 0);
167 } else {
168 this->robustSchedulerIndex->resize(x.size(), 0);
169 }
170 uint64_t numSchedulerChoices = 0;
171 for (uint64_t row = 0; row < this->A->getRowCount(); ++row) {
172 this->robustSchedulerIndex->at(row) = numSchedulerChoices;
173 numSchedulerChoices += this->A->getRow(row).getNumberOfEntries();
174 }
175 // Make sure that storage for scheduler choices is available
176 if (!this->schedulerChoices) {
177 this->schedulerChoices = std::vector<uint64_t>(numSchedulerChoices, 0);
178 } else {
179 this->schedulerChoices->resize(numSchedulerChoices, 0);
180 }
181 } else {
182 // Make sure that storage for scheduler choices is available
183 if (!this->schedulerChoices) {
184 this->schedulerChoices = std::vector<uint64_t>(x.size(), 0);
185 } else {
186 this->schedulerChoices->resize(x.size(), 0);
187 }
188 }
189
190 // Set the correct choices.
191 STORM_LOG_WARN_COND(!viOperatorTriv && !viOperatorNontriv,
192 "Expected VI operator to be initialized for scheduler extraction. Initializing now, but this is inefficient.");
193 if (!viOperatorTriv && !viOperatorNontriv) {
194 setUpViOperator();
195 }
196 if (viOperatorTriv) {
197 if constexpr (std::is_same<ValueType, storm::Interval>() && std::is_same<SolutionType, double>()) {
199 schedHelper.computeScheduler(x, b, dir, *this->schedulerChoices, robust, updateX ? &x : nullptr, this->robustSchedulerIndex);
200 } else {
201 STORM_LOG_ERROR("SchedulerTrackingHelper not implemented for this setting (trivial row grouping but not Interval->double).");
202 }
203 }
204 if (viOperatorNontriv) {
206 schedHelper.computeScheduler(x, b, dir, *this->schedulerChoices, robust, updateX ? &x : nullptr, this->robustSchedulerIndex);
207 }
208}
209
210template<typename ValueType, typename SolutionType>
211bool IterativeMinMaxLinearEquationSolver<ValueType, SolutionType>::solveInducedEquationSystem(
212 Environment const& env, std::unique_ptr<LinearEquationSolver<SolutionType>>& linearEquationSolver, std::vector<uint64_t> const& scheduler,
213 std::vector<SolutionType>& x, std::vector<ValueType>& subB, std::vector<ValueType> const& originalB, OptimizationDirection dir) const {
214 if constexpr (std::is_same_v<ValueType, storm::Interval>) {
215 if constexpr (std::is_same_v<SolutionType, storm::Interval>) {
216 STORM_LOG_THROW(false, storm::exceptions::NotImplementedException,
217 "We did not implement solving induced equation systems for interval-based models outside of robust VI.");
218 // Implementing this requires linear equation systems with different value types and solution types (or some appropriate casting)
219 return false;
220 }
221 STORM_LOG_ASSERT(subB.size() == x.size(), "Sizes of subB and x do not coincide.");
222 STORM_LOG_ASSERT(this->linearEquationSolverFactory != nullptr, "Wrong constructor was called.");
223
224 bool convertToEquationSystem = this->linearEquationSolverFactory->getEquationProblemFormat(env) == LinearEquationSolverProblemFormat::EquationSystem;
225
226 storm::storage::SparseMatrixBuilder<SolutionType> newMatrixBuilder(this->A->getRowCount(), this->A->getColumnCount(), this->A->getEntryCount());
227
228 // Robust VI scheduler is an order, compute the correct values for this order
229 auto schedulerIterator = scheduler.begin();
230 for (uint64_t rowIndex = 0; rowIndex < this->A->getRowCount(); rowIndex++) {
231 auto const& row = this->A->getRow(rowIndex);
232
233 static std::vector<SolutionType> tmp;
234 tmp.clear();
235
236 SolutionType probLeft = storm::utility::one<SolutionType>();
237
238 for (auto const& entry : row) {
239 tmp.push_back(entry.getValue().lower());
240 probLeft -= entry.getValue().lower();
241 }
242
243 auto const& rowIter = row.begin();
244 for (uint64_t i = 0; i < row.getNumberOfEntries(); i++, schedulerIterator++) {
245 if (!utility::isZero(probLeft)) {
246 auto const& entry = rowIter[*schedulerIterator];
247 auto const diameter = entry.getValue().upper() - entry.getValue().lower();
248 auto const value = utility::min(probLeft, diameter);
249 tmp[*schedulerIterator] += value;
250 probLeft -= value;
251 } else {
252 // Intentionally left empty: advance schedulerIterator to end of row
253 }
254 }
255
256 for (uint64_t i = 0; i < row.getNumberOfEntries(); i++) {
257 auto const& entry = rowIter[i];
258 newMatrixBuilder.addNextValue(rowIndex, entry.getColumn(), tmp[i]);
259 }
260 }
261 STORM_LOG_ASSERT(schedulerIterator == scheduler.end(), "Offset issue in scheduler?");
262
263 subB = originalB;
264
265 std::vector<SolutionType> b;
266 for (auto const& entry : subB) {
267 if (dir == OptimizationDirection::Maximize) {
268 b.push_back(entry.upper());
269 } else {
270 b.push_back(entry.lower());
271 }
272 }
273
274 auto const submatrix = newMatrixBuilder.build();
275
276 // Check whether the linear equation solver is already initialized
277 if (!linearEquationSolver) {
278 // Initialize the equation solver
279 linearEquationSolver = this->linearEquationSolverFactory->create(env, std::move(submatrix));
280 linearEquationSolver->setBoundsFromOtherSolver(*this);
281 linearEquationSolver->setCachingEnabled(true);
282 } else {
283 // If the equation solver is already initialized, it suffices to update the matrix
284 linearEquationSolver->setMatrix(std::move(submatrix));
285 }
286 // Solve the equation system for the 'DTMC' and return true upon success
287 return linearEquationSolver->solveEquations(env, x, b);
288 } else {
289 STORM_LOG_ASSERT(subB.size() == x.size(), "Sizes of subB and x do not coincide.");
290 STORM_LOG_ASSERT(this->linearEquationSolverFactory != nullptr, "Wrong constructor was called.");
291
292 // Resolve the nondeterminism according to the given scheduler.
293 bool convertToEquationSystem = this->linearEquationSolverFactory->getEquationProblemFormat(env) == LinearEquationSolverProblemFormat::EquationSystem;
295
296 submatrix = this->A->selectRowsFromRowGroups(scheduler, convertToEquationSystem);
297 if (convertToEquationSystem) {
298 submatrix.convertToEquationSystem();
299 }
300 storm::utility::vector::selectVectorValues<ValueType>(subB, scheduler, this->A->getRowGroupIndices(), originalB);
301
302 // Check whether the linear equation solver is already initialized
303 if (!linearEquationSolver) {
304 // Initialize the equation solver
305 linearEquationSolver = this->linearEquationSolverFactory->create(env, std::move(submatrix));
306 linearEquationSolver->setBoundsFromOtherSolver(*this);
307 linearEquationSolver->setCachingEnabled(true);
308 } else {
309 // If the equation solver is already initialized, it suffices to update the matrix
310 linearEquationSolver->setMatrix(std::move(submatrix));
311 }
312 // Solve the equation system for the 'DTMC' and return true upon success
313 return linearEquationSolver->solveEquations(env, x, subB);
314 }
315}
316
317template<typename ValueType, typename SolutionType>
318bool IterativeMinMaxLinearEquationSolver<ValueType, SolutionType>::solveEquationsPolicyIteration(Environment const& env, OptimizationDirection dir,
319 std::vector<SolutionType>& x,
320 std::vector<ValueType> const& b) const {
321 std::vector<storm::storage::sparse::state_type> scheduler =
322 this->hasInitialScheduler() ? this->getInitialScheduler() : std::vector<storm::storage::sparse::state_type>(this->A->getRowGroupCount());
323 return performPolicyIteration(env, dir, x, b, std::move(scheduler));
324}
325
326template<typename ValueType, typename SolutionType>
327bool IterativeMinMaxLinearEquationSolver<ValueType, SolutionType>::performPolicyIteration(
328 Environment const& env, OptimizationDirection dir, std::vector<SolutionType>& x, std::vector<ValueType> const& b,
329 std::vector<storm::storage::sparse::state_type>&& initialPolicy) const {
330 if constexpr (std::is_same_v<ValueType, storm::Interval>) {
331 STORM_LOG_THROW(false, storm::exceptions::NotImplementedException, "We did not implement policy iteration for interval-based models.");
332 return false;
333 } else {
334 std::vector<storm::storage::sparse::state_type> scheduler = std::move(initialPolicy);
335 // Get a vector for storing the right-hand side of the inner equation system.
336 if (!auxiliaryRowGroupVector) {
337 auxiliaryRowGroupVector = std::make_unique<std::vector<ValueType>>(this->A->getRowGroupCount());
338 }
339 std::vector<ValueType>& subB = *auxiliaryRowGroupVector;
340
341 // The solver that we will use throughout the procedure.
342 std::unique_ptr<storm::solver::LinearEquationSolver<ValueType>> solver;
343 // The linear equation solver should be at least as precise as this solver
344 std::unique_ptr<storm::Environment> environmentOfSolverStorage;
345 auto precOfSolver = env.solver().getPrecisionOfLinearEquationSolver(env.solver().getLinearEquationSolverType());
347 bool changePrecision = precOfSolver.first && precOfSolver.first.get() > env.solver().minMax().getPrecision();
348 bool changeRelative = precOfSolver.second && !precOfSolver.second.get() && env.solver().minMax().getRelativeTerminationCriterion();
349 if (changePrecision || changeRelative) {
350 environmentOfSolverStorage = std::make_unique<storm::Environment>(env);
351 boost::optional<storm::RationalNumber> newPrecision;
352 boost::optional<bool> newRelative;
353 if (changePrecision) {
354 newPrecision = env.solver().minMax().getPrecision();
355 }
356 if (changeRelative) {
357 newRelative = true;
358 }
359 environmentOfSolverStorage->solver().setLinearEquationSolverPrecision(newPrecision, newRelative);
360 }
361 }
362 storm::Environment const& environmentOfSolver = environmentOfSolverStorage ? *environmentOfSolverStorage : env;
363
365 uint64_t iterations = 0;
366 this->startMeasureProgress();
367 do {
368 // Solve the equation system for the 'DTMC'.
369 solveInducedEquationSystem(environmentOfSolver, solver, scheduler, x, subB, b, dir);
370
371 // Go through the multiplication result and see whether we can improve any of the choices.
372 bool schedulerImproved = false;
373 // Group refers to the state number
374 for (uint_fast64_t group = 0; group < this->A->getRowGroupCount(); ++group) {
375 if (!this->choiceFixedForRowGroup || !this->choiceFixedForRowGroup.get()[group]) {
376 // Only update when the choice is not fixed
377 uint_fast64_t currentChoice = scheduler[group];
378 for (uint_fast64_t choice = this->A->getRowGroupIndices()[group]; choice < this->A->getRowGroupIndices()[group + 1]; ++choice) {
379 // If the choice is the currently selected one, we can skip it.
380 if (choice - this->A->getRowGroupIndices()[group] == currentChoice) {
381 continue;
382 }
383
384 // Create the value of the choice.
385 ValueType choiceValue = storm::utility::zero<ValueType>();
386 for (auto const& entry : this->A->getRow(choice)) {
387 choiceValue += entry.getValue() * x[entry.getColumn()];
388 }
389 choiceValue += b[choice];
390
391 // If the value is strictly better than the solution of the inner system, we need to improve the scheduler.
392 // TODO: If the underlying solver is not precise, this might run forever (i.e. when a state has two choices where the (exact) values are
393 // equal). only changing the scheduler if the values are not equal (modulo precision) would make this unsound.
394 if (valueImproved(dir, x[group], choiceValue)) {
395 schedulerImproved = true;
396 scheduler[group] = choice - this->A->getRowGroupIndices()[group];
397 x[group] = std::move(choiceValue);
398 }
399 }
400 }
401 }
402
403 // If the scheduler did not improve, we are done.
404 if (!schedulerImproved) {
406 }
407
408 // Update environment variables.
409 ++iterations;
410 status =
411 this->updateStatus(status, x, dir == storm::OptimizationDirection::Minimize ? SolverGuarantee::GreaterOrEqual : SolverGuarantee::LessOrEqual,
412 iterations, env.solver().minMax().getMaximalNumberOfIterations());
413
414 // Potentially show progress.
415 this->showProgressIterative(iterations);
416 } while (status == SolverStatus::InProgress);
417
418 STORM_LOG_INFO("Number of iterations: " << iterations);
419 this->reportStatus(status, iterations);
420
421 // If requested, we store the scheduler for retrieval.
422 if (this->isTrackSchedulerSet()) {
423 this->schedulerChoices = std::move(scheduler);
424 }
425
426 if (!this->isCachingEnabled()) {
427 clearCache();
428 }
429
430 return status == SolverStatus::Converged || status == SolverStatus::TerminatedEarly;
431 }
432}
433
434template<typename ValueType, typename SolutionType>
435bool IterativeMinMaxLinearEquationSolver<ValueType, SolutionType>::valueImproved(OptimizationDirection dir, ValueType const& value1,
436 ValueType const& value2) const {
437 if (dir == OptimizationDirection::Minimize) {
438 return value2 < value1;
439 } else {
440 return value2 > value1;
441 }
442}
443
444template<typename ValueType, typename SolutionType>
446 Environment const& env, boost::optional<storm::solver::OptimizationDirection> const& direction, bool const& hasInitialScheduler) const {
447 auto method = getMethod(env, storm::NumberTraits<ValueType>::IsExact || env.solver().isForceExact());
448
449 // Check whether a linear equation solver is needed and potentially start with its requirements
450 bool needsLinEqSolver = false;
451 needsLinEqSolver |= method == MinMaxMethod::PolicyIteration;
452 needsLinEqSolver |= method == MinMaxMethod::ValueIteration && (this->hasInitialScheduler() || hasInitialScheduler);
453 needsLinEqSolver |= method == MinMaxMethod::ViToPi;
454
456 if constexpr (std::is_same_v<ValueType, storm::Interval>) {
457 STORM_LOG_ASSERT(!needsLinEqSolver, "Intervals should not require a linear equation solver.");
458 // nothing to be done;
459 } else if (needsLinEqSolver) {
460 requirements = MinMaxLinearEquationSolverRequirements(this->linearEquationSolverFactory->getRequirements(env));
461 } else {
462 // nothing to be done.
463 }
464
465 if (method == MinMaxMethod::ValueIteration) {
466 if (!this->hasUniqueSolution()) { // Traditional value iteration has no requirements if the solution is unique.
467 // Computing a scheduler is only possible if the solution is unique
468 if (env.solver().minMax().isForceRequireUnique() || this->isTrackSchedulerSet()) {
469 requirements.requireUniqueSolution();
470 } else {
471 // As we want the smallest (largest) solution for maximizing (minimizing) equation systems, we have to approach the solution from below (above).
472 if (!direction || direction.get() == OptimizationDirection::Maximize) {
473 requirements.requireLowerBounds();
474 }
475 if (!direction || direction.get() == OptimizationDirection::Minimize) {
476 requirements.requireUpperBounds();
477 }
478 }
479 }
480 } else if (method == MinMaxMethod::OptimisticValueIteration) {
481 // OptimisticValueIteration always requires lower bounds and a unique solution.
482 if (!this->hasUniqueSolution()) {
483 requirements.requireUniqueSolution();
484 }
485 requirements.requireLowerBounds();
486
487 } else if (method == MinMaxMethod::GuessingValueIteration) {
488 // Guessing value iteration requires a unique solution and lower+upper bounds
489 if (!this->hasUniqueSolution()) {
490 requirements.requireUniqueSolution();
491 }
492 requirements.requireBounds();
493 } else if (method == MinMaxMethod::IntervalIteration) {
494 // Interval iteration requires a unique solution and lower+upper bounds
495 if (!this->hasUniqueSolution()) {
496 requirements.requireUniqueSolution();
497 }
498 requirements.requireBounds();
499 } else if (method == MinMaxMethod::RationalSearch) {
500 // Rational search needs to approach the solution from below.
501 requirements.requireLowerBounds();
502 // The solution needs to be unique in case of minimizing or in cases where we want a scheduler.
503 if (!this->hasUniqueSolution()) {
504 // RationalSearch guesses and verifies a fixpoint and terminates once a fixpoint is found. To ensure that the guessed fixpoint is the
505 // correct one, we enforce uniqueness.
506 requirements.requireUniqueSolution();
507 }
508 } else if (method == MinMaxMethod::PolicyIteration) {
509 // The initial scheduler shall not select an end component
510 if (!this->hasUniqueSolution() && env.solver().minMax().isForceRequireUnique()) {
511 requirements.requireUniqueSolution();
512 }
513 if (!this->hasNoEndComponents() && !this->hasInitialScheduler()) {
514 requirements.requireValidInitialScheduler();
515 }
516 } else if (method == MinMaxMethod::SoundValueIteration) {
517 if (!this->hasUniqueSolution()) {
518 requirements.requireUniqueSolution();
519 }
520 requirements.requireBounds(false);
521 } else if (method == MinMaxMethod::ViToPi) {
522 // Since we want to use value iteration to extract an initial scheduler, the solution has to be unique.
523 if (!this->hasUniqueSolution()) {
524 requirements.requireUniqueSolution();
525 }
526 } else {
527 STORM_LOG_THROW(false, storm::exceptions::InvalidEnvironmentException, "Unsupported technique for iterative MinMax linear equation solver.");
528 }
529 return requirements;
530}
531
532template<typename ValueType, typename SolutionType>
533ValueType computeMaxAbsDiff(std::vector<ValueType> const& allValues, storm::storage::BitVector const& relevantValues, std::vector<ValueType> const& oldValues) {
534 ValueType result = storm::utility::zero<ValueType>();
535 auto oldValueIt = oldValues.begin();
536 for (auto value : relevantValues) {
537 result = storm::utility::max<ValueType>(result, storm::utility::abs<ValueType>(allValues[value] - *oldValueIt));
538 ++oldValueIt;
539 }
540 return result;
541}
542
543template<typename ValueType, typename SolutionType>
544ValueType computeMaxAbsDiff(std::vector<ValueType> const& allOldValues, std::vector<ValueType> const& allNewValues,
545 storm::storage::BitVector const& relevantValues) {
546 ValueType result = storm::utility::zero<ValueType>();
547 for (auto value : relevantValues) {
548 result = storm::utility::max<ValueType>(result, storm::utility::abs<ValueType>(allNewValues[value] - allOldValues[value]));
549 }
550 return result;
551}
552
553template<typename ValueType, typename SolutionType>
554bool IterativeMinMaxLinearEquationSolver<ValueType, SolutionType>::solveEquationsOptimisticValueIteration(Environment const& env, OptimizationDirection dir,
555 std::vector<SolutionType>& x,
556 std::vector<ValueType> const& b) const {
557 if constexpr (std::is_same_v<ValueType, storm::Interval>) {
558 STORM_LOG_THROW(false, storm::exceptions::NotImplementedException, "We did not implement optimistic value iteration for interval-based models.");
559 return false;
560 } else {
562 // If all entries are zero, OVI might run in an endless loop. However, the result is easy in this case.
563 x.assign(x.size(), storm::utility::zero<SolutionType>());
564 if (this->isTrackSchedulerSet()) {
565 this->schedulerChoices = std::vector<uint_fast64_t>(x.size(), 0);
566 }
567 return true;
568 }
569
570 setUpViOperator();
571
572 helper::OptimisticValueIterationHelper<ValueType, false> oviHelper(viOperatorNontriv);
573 auto prec = storm::utility::convertNumber<ValueType>(env.solver().minMax().getPrecision());
574 std::optional<ValueType> lowerBound, upperBound;
575 if (this->hasLowerBound()) {
576 lowerBound = this->getLowerBound(true);
577 }
578 if (this->hasUpperBound()) {
579 upperBound = this->getUpperBound(true);
580 }
581 uint64_t numIterations{0};
582 auto oviCallback = [&](SolverStatus const& current, std::vector<ValueType> const& v) {
583 this->showProgressIterative(numIterations);
584 return this->updateStatus(current, v, SolverGuarantee::LessOrEqual, numIterations, env.solver().minMax().getMaximalNumberOfIterations());
585 };
586 this->createLowerBoundsVector(x);
587 std::optional<ValueType> guessingFactor;
588 if (env.solver().ovi().getUpperBoundGuessingFactor()) {
589 guessingFactor = storm::utility::convertNumber<ValueType>(*env.solver().ovi().getUpperBoundGuessingFactor());
590 }
591 this->startMeasureProgress();
592 auto status = oviHelper.OVI(x, b, numIterations, env.solver().minMax().getRelativeTerminationCriterion(), prec, dir, guessingFactor, lowerBound,
593 upperBound, oviCallback);
594 this->reportStatus(status, numIterations);
595
596 // If requested, we store the scheduler for retrieval.
597 if (this->isTrackSchedulerSet()) {
598 this->extractScheduler(x, b, dir, this->isUncertaintyRobust());
599 }
600
601 if (!this->isCachingEnabled()) {
602 clearCache();
603 }
604
605 return status == SolverStatus::Converged || status == SolverStatus::TerminatedEarly;
606 }
607}
608
609template<typename ValueType, typename SolutionType>
610bool IterativeMinMaxLinearEquationSolver<ValueType, SolutionType>::solveEquationsGuessingValueIteration(Environment const& env, OptimizationDirection dir,
611 std::vector<SolutionType>& x,
612 std::vector<ValueType> const& b) const {
613 if constexpr (std::is_same_v<ValueType, storm::Interval>) {
614 STORM_LOG_THROW(false, storm::exceptions::NotImplementedException, "We did not implement guessing value iteration for interval-based models.");
615 return false;
616 } else {
617 setUpViOperator();
618
619 auto& lowerX = x;
620 auto upperX = std::make_unique<std::vector<SolutionType>>(x.size());
621
623
624 uint64_t numIterations{0};
625
626 auto gviCallback = [&](helper::GVIData<ValueType> const& data) {
627 this->showProgressIterative(numIterations);
628 bool terminateEarly = this->hasCustomTerminationCondition() && this->getTerminationCondition().terminateNow(data.x, SolverGuarantee::LessOrEqual) &&
629 this->getTerminationCondition().terminateNow(data.y, SolverGuarantee::GreaterOrEqual);
630 return this->updateStatus(data.status, terminateEarly, numIterations, env.solver().minMax().getMaximalNumberOfIterations());
631 };
632
633 this->createLowerBoundsVector(lowerX);
634 this->createUpperBoundsVector(*upperX);
635
636 this->startMeasureProgress();
637 auto statusIters = helper.solveEquations(lowerX, *upperX, b, numIterations,
638 storm::utility::convertNumber<ValueType>(env.solver().minMax().getPrecision()), dir, gviCallback);
639 auto two = storm::utility::convertNumber<ValueType>(2.0);
640 storm::utility::vector::applyPointwise<ValueType, ValueType, ValueType>(
641 lowerX, *upperX, x, [&two](ValueType const& a, ValueType const& b) -> ValueType { return (a + b) / two; });
642
643 this->reportStatus(statusIters, numIterations);
644
645 // If requested, we store the scheduler for retrieval.
646 if (this->isTrackSchedulerSet()) {
647 this->extractScheduler(x, b, dir, this->isUncertaintyRobust());
648 }
649
650 if (!this->isCachingEnabled()) {
651 clearCache();
652 }
653
654 return statusIters == SolverStatus::Converged || statusIters == SolverStatus::TerminatedEarly;
655 }
656}
657
658template<typename ValueType, typename SolutionType>
659bool IterativeMinMaxLinearEquationSolver<ValueType, SolutionType>::solveEquationsValueIteration(Environment const& env, OptimizationDirection dir,
660 std::vector<SolutionType>& x,
661 std::vector<ValueType> const& b) const {
662 setUpViOperator();
663 // By default, we can not provide any guarantee
665
666 if (this->hasInitialScheduler()) {
667 if (!auxiliaryRowGroupVector) {
668 auxiliaryRowGroupVector = std::make_unique<std::vector<ValueType>>(this->A->getRowGroupCount());
669 }
670 // Solve the equation system induced by the initial scheduler.
671 std::unique_ptr<storm::solver::LinearEquationSolver<SolutionType>> linEqSolver;
672 // The linear equation solver should be at least as precise as this solver
673 std::unique_ptr<storm::Environment> environmentOfSolverStorage;
674 auto precOfSolver = env.solver().getPrecisionOfLinearEquationSolver(env.solver().getLinearEquationSolverType());
676 bool changePrecision = precOfSolver.first && precOfSolver.first.get() > env.solver().minMax().getPrecision();
677 bool changeRelative = precOfSolver.second && !precOfSolver.second.get() && env.solver().minMax().getRelativeTerminationCriterion();
678 if (changePrecision || changeRelative) {
679 environmentOfSolverStorage = std::make_unique<storm::Environment>(env);
680 boost::optional<storm::RationalNumber> newPrecision;
681 boost::optional<bool> newRelative;
682 if (changePrecision) {
683 newPrecision = env.solver().minMax().getPrecision();
684 }
685 if (changeRelative) {
686 newRelative = true;
687 }
688 environmentOfSolverStorage->solver().setLinearEquationSolverPrecision(newPrecision, newRelative);
689 }
690 }
691 storm::Environment const& environmentOfSolver = environmentOfSolverStorage ? *environmentOfSolverStorage : env;
692
693 solveInducedEquationSystem(environmentOfSolver, linEqSolver, this->getInitialScheduler(), x, *auxiliaryRowGroupVector, b, dir);
694 // If we were given an initial scheduler and are maximizing (minimizing), our current solution becomes
695 // always less-or-equal (greater-or-equal) than the actual solution.
697 } else if (!this->hasUniqueSolution()) {
698 if (maximize(dir)) {
699 this->createLowerBoundsVector(x);
701 } else {
702 this->createUpperBoundsVector(x);
704 }
705 } else if (this->hasCustomTerminationCondition()) {
706 if (this->getTerminationCondition().requiresGuarantee(SolverGuarantee::LessOrEqual) && this->hasLowerBound()) {
707 this->createLowerBoundsVector(x);
709 } else if (this->getTerminationCondition().requiresGuarantee(SolverGuarantee::GreaterOrEqual) && this->hasUpperBound()) {
710 this->createUpperBoundsVector(x);
712 }
713 }
714
715 uint64_t numIterations{0};
716 auto viCallback = [&](SolverStatus const& current) {
717 this->showProgressIterative(numIterations);
718 return this->updateStatus(current, x, guarantee, numIterations, env.solver().minMax().getMaximalNumberOfIterations());
719 };
720 this->startMeasureProgress();
721 // This code duplication is necessary because the helper class is different for the two cases.
722 if (this->A->hasTrivialRowGrouping()) {
724
725 auto status = viHelper.VI(x, b, numIterations, env.solver().minMax().getRelativeTerminationCriterion(),
726 storm::utility::convertNumber<SolutionType>(env.solver().minMax().getPrecision()), dir, viCallback,
727 env.solver().minMax().getMultiplicationStyle(), this->isUncertaintyRobust());
728 this->reportStatus(status, numIterations);
729
730 // If requested, we store the scheduler for retrieval.
731 if (this->isTrackSchedulerSet()) {
732 this->extractScheduler(x, b, dir, this->isUncertaintyRobust());
733 }
734
735 if (!this->isCachingEnabled()) {
736 clearCache();
737 }
738
739 return status == SolverStatus::Converged || status == SolverStatus::TerminatedEarly;
740 } else {
742
743 auto status = viHelper.VI(x, b, numIterations, env.solver().minMax().getRelativeTerminationCriterion(),
744 storm::utility::convertNumber<SolutionType>(env.solver().minMax().getPrecision()), dir, viCallback,
745 env.solver().minMax().getMultiplicationStyle(), this->isUncertaintyRobust());
746 this->reportStatus(status, numIterations);
747
748 // If requested, we store the scheduler for retrieval.
749 if (this->isTrackSchedulerSet()) {
750 this->extractScheduler(x, b, dir, this->isUncertaintyRobust());
751 }
752
753 if (!this->isCachingEnabled()) {
754 clearCache();
755 }
756
757 return status == SolverStatus::Converged || status == SolverStatus::TerminatedEarly;
758 }
759}
760
761template<typename ValueType, typename SolutionType>
762void preserveOldRelevantValues(std::vector<ValueType> const& allValues, storm::storage::BitVector const& relevantValues, std::vector<ValueType>& oldValues) {
763 storm::utility::vector::selectVectorValues(oldValues, relevantValues, allValues);
764}
765
772template<typename ValueType, typename SolutionType>
773bool IterativeMinMaxLinearEquationSolver<ValueType, SolutionType>::solveEquationsIntervalIteration(Environment const& env, OptimizationDirection dir,
774 std::vector<SolutionType>& x,
775 std::vector<ValueType> const& b) const {
776 if constexpr (std::is_same_v<ValueType, storm::Interval>) {
777 STORM_LOG_THROW(false, storm::exceptions::NotImplementedException, "We did not implement intervaliteration for interval-based models");
778 return false;
779 } else {
780 setUpViOperator();
781 helper::IntervalIterationHelper<ValueType, false> iiHelper(viOperatorNontriv);
782 auto prec = storm::utility::convertNumber<ValueType>(env.solver().minMax().getPrecision());
783 auto lowerBoundsCallback = [&](std::vector<SolutionType>& vector) { this->createLowerBoundsVector(vector); };
784 auto upperBoundsCallback = [&](std::vector<SolutionType>& vector) { this->createUpperBoundsVector(vector); };
785
786 uint64_t numIterations{0};
787 auto iiCallback = [&](helper::IIData<ValueType> const& data) {
788 this->showProgressIterative(numIterations);
789 bool terminateEarly = this->hasCustomTerminationCondition() && this->getTerminationCondition().terminateNow(data.x, SolverGuarantee::LessOrEqual) &&
790 this->getTerminationCondition().terminateNow(data.y, SolverGuarantee::GreaterOrEqual);
791 return this->updateStatus(data.status, terminateEarly, numIterations, env.solver().minMax().getMaximalNumberOfIterations());
792 };
793 std::optional<storm::storage::BitVector> optionalRelevantValues;
794 if (this->hasRelevantValues()) {
795 optionalRelevantValues = this->getRelevantValues();
796 }
797 this->startMeasureProgress();
798 auto status = iiHelper.II(x, b, numIterations, env.solver().minMax().getRelativeTerminationCriterion(), prec, lowerBoundsCallback, upperBoundsCallback,
799 dir, iiCallback, optionalRelevantValues);
800 this->reportStatus(status, numIterations);
801
802 // If requested, we store the scheduler for retrieval.
803 if (this->isTrackSchedulerSet()) {
804 this->extractScheduler(x, b, dir, this->isUncertaintyRobust());
805 }
806
807 if (!this->isCachingEnabled()) {
808 clearCache();
809 }
810
811 return status == SolverStatus::Converged || status == SolverStatus::TerminatedEarly;
812 }
813}
814
815template<typename ValueType, typename SolutionType>
816bool IterativeMinMaxLinearEquationSolver<ValueType, SolutionType>::solveEquationsSoundValueIteration(Environment const& env, OptimizationDirection dir,
817 std::vector<SolutionType>& x,
818 std::vector<ValueType> const& b) const {
819 if constexpr (std::is_same_v<ValueType, storm::Interval>) {
820 STORM_LOG_THROW(false, storm::exceptions::NotImplementedException, "SoundVI does not handle interval-based models");
821 return false;
822 } else {
823 // Prepare the solution vectors and the helper.
824 assert(x.size() == this->A->getRowGroupCount());
825
826 std::optional<ValueType> lowerBound, upperBound;
827 if (this->hasLowerBound()) {
828 lowerBound = this->getLowerBound(true);
829 }
830 if (this->hasUpperBound()) {
831 upperBound = this->getUpperBound(true);
832 }
833
834 setUpViOperator();
835
836 auto precision = storm::utility::convertNumber<ValueType>(env.solver().minMax().getPrecision());
837 uint64_t numIterations{0};
838 auto sviCallback = [&](typename helper::SoundValueIterationHelper<ValueType, false>::SVIData const& current) {
839 this->showProgressIterative(numIterations);
840 return this->updateStatus(current.status,
841 this->hasCustomTerminationCondition() && current.checkCustomTerminationCondition(this->getTerminationCondition()),
842 numIterations, env.solver().minMax().getMaximalNumberOfIterations());
843 };
844 this->startMeasureProgress();
845 helper::SoundValueIterationHelper<ValueType, false> sviHelper(viOperatorNontriv);
846 std::optional<storm::storage::BitVector> optionalRelevantValues;
847 if (this->hasRelevantValues()) {
848 optionalRelevantValues = this->getRelevantValues();
849 }
850 auto status = sviHelper.SVI(x, b, numIterations, env.solver().minMax().getRelativeTerminationCriterion(), precision, dir, lowerBound, upperBound,
851 sviCallback, optionalRelevantValues);
852
853 // If requested, we store the scheduler for retrieval.
854 if (this->isTrackSchedulerSet()) {
855 this->extractScheduler(x, b, dir, this->isUncertaintyRobust());
856 }
857
858 this->reportStatus(status, numIterations);
859
860 if (!this->isCachingEnabled()) {
861 clearCache();
862 }
863
864 return status == SolverStatus::Converged || status == SolverStatus::TerminatedEarly;
865 }
866}
867
868template<typename ValueType, typename SolutionType>
869bool IterativeMinMaxLinearEquationSolver<ValueType, SolutionType>::solveEquationsViToPi(Environment const& env, OptimizationDirection dir,
870 std::vector<SolutionType>& x, std::vector<ValueType> const& b) const {
871 if constexpr (std::is_same_v<ValueType, storm::Interval>) {
872 STORM_LOG_THROW(false, storm::exceptions::NotImplementedException, "ViToPi does not handle interval-based models");
873 return false;
874 }
875 // First create an (inprecise) vi solver to get a good initial strategy for the (potentially precise) policy iteration solver.
876 std::vector<storm::storage::sparse::state_type> initialSched;
877 {
878 Environment viEnv = env;
879 viEnv.solver().minMax().setMethod(MinMaxMethod::ValueIteration);
880 viEnv.solver().setForceExact(false);
881 viEnv.solver().setForceSoundness(false);
882 auto impreciseSolver = GeneralMinMaxLinearEquationSolverFactory<double>().create(viEnv, this->A->template toValueType<double>());
883 impreciseSolver->setHasUniqueSolution(this->hasUniqueSolution());
884 impreciseSolver->setTrackScheduler(true);
885 if (this->hasInitialScheduler()) {
886 auto initSched = this->getInitialScheduler();
887 impreciseSolver->setInitialScheduler(std::move(initSched));
888 }
889 auto impreciseSolverReq = impreciseSolver->getRequirements(viEnv, dir);
890 STORM_LOG_THROW(!impreciseSolverReq.hasEnabledCriticalRequirement(), storm::exceptions::UnmetRequirementException,
891 "The value-iteration based solver has an unmet requirement: " << impreciseSolverReq.getEnabledRequirementsAsString());
892 impreciseSolver->setRequirementsChecked(true);
893 auto xVi = storm::utility::vector::convertNumericVector<double>(x);
894 auto bVi = storm::utility::vector::convertNumericVector<double>(b);
895 impreciseSolver->solveEquations(viEnv, dir, xVi, bVi);
896 initialSched = impreciseSolver->getSchedulerChoices();
897 }
898 STORM_LOG_INFO("Found initial policy using Value Iteration. Starting Policy iteration now.");
899 return performPolicyIteration(env, dir, x, b, std::move(initialSched));
900}
901
902template<typename ValueType, typename SolutionType>
903bool IterativeMinMaxLinearEquationSolver<ValueType, SolutionType>::solveEquationsRationalSearch(Environment const& env, OptimizationDirection dir,
904 std::vector<SolutionType>& x,
905 std::vector<ValueType> const& b) const {
906 if constexpr (std::is_same_v<ValueType, storm::Interval>) {
907 STORM_LOG_THROW(false, storm::exceptions::NotImplementedException, "Rational search does not handle interval-based models");
908 return false;
909 } else {
910 // Set up two value iteration operators. One for exact and one for imprecise computations
911 setUpViOperator();
912 std::shared_ptr<helper::ValueIterationOperator<storm::RationalNumber, false>> exactOp;
913 std::shared_ptr<helper::ValueIterationOperator<double, false>> impreciseOp;
914 std::function<bool(uint64_t, uint64_t)> fixedChoicesCallback;
915 if (this->choiceFixedForRowGroup) {
916 // Ignore those rows that are not selected
917 assert(this->initialScheduler);
918 fixedChoicesCallback = [&](uint64_t groupIndex, uint64_t localRowIndex) {
919 return this->choiceFixedForRowGroup->get(groupIndex) && this->initialScheduler->at(groupIndex) != localRowIndex;
920 };
921 }
922
923 if constexpr (std::is_same_v<ValueType, storm::RationalNumber>) {
924 exactOp = viOperatorNontriv;
925 impreciseOp = std::make_shared<helper::ValueIterationOperator<double, false>>();
926 impreciseOp->setMatrixBackwards(this->A->template toValueType<double>(), &this->A->getRowGroupIndices());
927 if (this->choiceFixedForRowGroup) {
928 impreciseOp->setIgnoredRows(true, fixedChoicesCallback);
929 }
930 } else if constexpr (std::is_same_v<ValueType, double>) {
931 impreciseOp = viOperatorNontriv;
932 exactOp = std::make_shared<helper::ValueIterationOperator<storm::RationalNumber, false>>();
933 exactOp->setMatrixBackwards(this->A->template toValueType<storm::RationalNumber>(), &this->A->getRowGroupIndices());
934 if (this->choiceFixedForRowGroup) {
935 exactOp->setIgnoredRows(true, fixedChoicesCallback);
936 }
937 }
938
940 uint64_t numIterations{0};
941 auto rsCallback = [&](SolverStatus const& current) {
942 this->showProgressIterative(numIterations);
943 return this->updateStatus(current, x, SolverGuarantee::None, numIterations, env.solver().minMax().getMaximalNumberOfIterations());
944 };
945 this->startMeasureProgress();
946 auto status = rsHelper.RS(x, b, numIterations, storm::utility::convertNumber<ValueType>(env.solver().minMax().getPrecision()), dir, rsCallback);
947
948 this->reportStatus(status, numIterations);
949
950 // If requested, we store the scheduler for retrieval.
951 if (this->isTrackSchedulerSet()) {
952 this->extractScheduler(x, b, dir, this->isUncertaintyRobust());
953 }
954
955 if (!this->isCachingEnabled()) {
956 clearCache();
957 }
958
959 return status == SolverStatus::Converged || status == SolverStatus::TerminatedEarly;
960 }
961}
962
963template<typename ValueType, typename SolutionType>
965 auxiliaryRowGroupVector.reset();
966 if (viOperatorTriv) {
967 viOperatorTriv.reset();
968 }
969 if (viOperatorNontriv) {
970 viOperatorNontriv.reset();
971 }
973}
974
978
979} // namespace solver
980} // namespace storm
SolverEnvironment & solver()
uint64_t const & getMaximalNumberOfIterations() const
storm::RationalNumber const & getPrecision() const
storm::solver::MinMaxMethod const & getMethod() const
bool const & getRelativeTerminationCriterion() const
std::optional< storm::RationalNumber > const & getUpperBoundGuessingFactor() const
OviSolverEnvironment const & ovi() const
MinMaxSolverEnvironment & minMax()
virtual void clearCache() const override
Clears the currently cached data that has been stored during previous calls of the solver.
IterativeMinMaxLinearEquationSolver(std::unique_ptr< LinearEquationSolverFactory< SolutionType > > &&linearEquationSolverFactory)
virtual bool internalSolveEquations(Environment const &env, OptimizationDirection dir, std::vector< SolutionType > &x, std::vector< ValueType > const &b) const override
virtual MinMaxLinearEquationSolverRequirements getRequirements(Environment const &env, boost::optional< storm::solver::OptimizationDirection > const &direction=boost::none, bool const &hasInitialScheduler=false) const override
Retrieves the requirements of this solver for solving equations with the current settings.
virtual void clearCache() const
Clears the currently cached data that has been stored during previous calls of the solver.
MinMaxLinearEquationSolverRequirements & requireBounds(bool critical=true)
MinMaxLinearEquationSolverRequirements & requireUniqueSolution(bool critical=true)
MinMaxLinearEquationSolverRequirements & requireLowerBounds(bool critical=true)
MinMaxLinearEquationSolverRequirements & requireValidInitialScheduler(bool critical=true)
MinMaxLinearEquationSolverRequirements & requireUpperBounds(bool critical=true)
Helper class to extract optimal scheduler choices from a MinMax equation system solution.
A bit vector that is internally represented as a vector of 64-bit values.
Definition BitVector.h:18
A class that can be used to build a sparse matrix by adding value by value.
A class that holds a possibly non-square matrix in the compressed row storage format.
void convertToEquationSystem()
Transforms the matrix into an equation system.
SparseMatrix selectRowsFromRowGroups(std::vector< index_type > const &rowGroupToRowIndexMapping, bool insertDiagonalEntries=true) const
Selects exactly one row from each row group of this matrix and returns the resulting matrix.
#define STORM_LOG_INFO(message)
Definition logging.h:29
#define STORM_LOG_WARN(message)
Definition logging.h:30
#define STORM_LOG_ERROR(message)
Definition logging.h:31
#define STORM_LOG_ASSERT(cond, message)
Definition macros.h:11
#define STORM_LOG_WARN_COND(cond, message)
Definition macros.h:38
#define STORM_LOG_THROW(cond, exception, message)
Definition macros.h:30
SFTBDDChecker::ValueType ValueType
bool constexpr maximize(OptimizationDirection d)
void preserveOldRelevantValues(std::vector< ValueType > const &allValues, storm::storage::BitVector const &relevantValues, std::vector< ValueType > &oldValues)
std::string toString(GurobiSolverMethod const &method)
Yields a string representation of the GurobiSolverMethod.
ValueType computeMaxAbsDiff(std::vector< ValueType > const &allValues, storm::storage::BitVector const &relevantValues, std::vector< ValueType > const &oldValues)
void selectVectorValues(std::vector< T > &vector, storm::storage::BitVector const &positions, std::vector< T > const &values)
Selects the elements from a vector at the specified positions and writes them consecutively into anot...
Definition vector.h:188
bool hasNonZeroEntry(std::vector< T > const &v)
Definition vector.h:1103
ValueType min(ValueType const &first, ValueType const &second)
bool isZero(ValueType const &a)
Definition constants.cpp:41
LabParser.cpp.
Definition cli.cpp:18