Storm
A Modern Probabilistic Model Checker
Loading...
Searching...
No Matches
SymbolicMinMaxLinearEquationSolver.cpp
Go to the documentation of this file.
2
4
7
9
11
15#include "storm/utility/dd.h"
17
18namespace storm {
19namespace solver {
20
21template<storm::dd::DdType DdType, typename ValueType>
23 : SymbolicEquationSolver<DdType, ValueType>(), uniqueSolution(false), requirementsChecked(false) {
24 // Intentionally left empty.
25}
26
27template<storm::dd::DdType DdType, typename ValueType>
29 storm::dd::Add<DdType, ValueType> const& A, storm::dd::Bdd<DdType> const& allRows, storm::dd::Bdd<DdType> const& illegalMask,
30 std::set<storm::expressions::Variable> const& rowMetaVariables, std::set<storm::expressions::Variable> const& columnMetaVariables,
31 std::set<storm::expressions::Variable> const& choiceVariables,
32 std::vector<std::pair<storm::expressions::Variable, storm::expressions::Variable>> const& rowColumnMetaVariablePairs,
33 std::unique_ptr<SymbolicLinearEquationSolverFactory<DdType, ValueType>>&& linearEquationSolverFactory)
34 : SymbolicEquationSolver<DdType, ValueType>(allRows),
35 A(A),
36 illegalMask(illegalMask),
37 illegalMaskAdd(illegalMask.ite(A.getDdManager().getConstant(storm::utility::infinity<ValueType>()), A.getDdManager().template getAddZero<ValueType>())),
38 rowMetaVariables(rowMetaVariables),
39 columnMetaVariables(columnMetaVariables),
40 choiceVariables(choiceVariables),
41 rowColumnMetaVariablePairs(rowColumnMetaVariablePairs),
42 linearEquationSolverFactory(std::move(linearEquationSolverFactory)),
43 uniqueSolution(false),
44 requirementsChecked(false) {
45 // Intentionally left empty.
46}
47
48template<storm::dd::DdType DdType, typename ValueType>
49MinMaxMethod SymbolicMinMaxLinearEquationSolver<DdType, ValueType>::getMethod(Environment const& env, bool isExactMode) const {
50 // Adjust the method if none was specified and we are using rational numbers.
51 auto method = env.solver().minMax().getMethod();
52
53 if (isExactMode && method != MinMaxMethod::RationalSearch) {
54 if (env.solver().minMax().isMethodSetFromDefault()) {
56 "Selecting 'rational search' as the solution technique to guarantee exact results. If you want to override this, please explicitly specify a "
57 "different method.");
58 method = MinMaxMethod::RationalSearch;
59 } else {
60 STORM_LOG_WARN("The selected solution method does not guarantee exact results.");
61 }
62 }
63 if (method != MinMaxMethod::ValueIteration && method != MinMaxMethod::PolicyIteration && method != MinMaxMethod::RationalSearch) {
64 STORM_LOG_WARN("Selected method is not supported for this solver, switching to value iteration.");
65 method = MinMaxMethod::ValueIteration;
66 }
67
68 return method;
69}
70
71template<storm::dd::DdType DdType, typename ValueType>
75 storm::dd::Add<DdType, ValueType> const& b) const {
76 STORM_LOG_WARN_COND_DEBUG(this->isRequirementsCheckedSet(),
77 "The requirements of the solver have not been marked as checked. Please provide the appropriate check or mark the requirements "
78 "as checked (if applicable).");
79
80 switch (getMethod(env, std::is_same<ValueType, storm::RationalNumber>::value)) {
81 case MinMaxMethod::ValueIteration:
82 return solveEquationsValueIteration(env, dir, x, b);
83 break;
84 case MinMaxMethod::PolicyIteration:
85 return solveEquationsPolicyIteration(env, dir, x, b);
86 break;
87 case MinMaxMethod::RationalSearch:
88 return solveEquationsRationalSearch(env, dir, x, b);
89 break;
90 default:
91 STORM_LOG_THROW(false, storm::exceptions::InvalidEnvironmentException, "The selected min max technique is not supported by this solver.");
92 }
94}
95
96template<storm::dd::DdType DdType, typename ValueType>
100 storm::dd::Add<DdType, ValueType> const& b, ValueType const& precision,
101 bool relativeTerminationCriterion, uint64_t maximalIterations) const {
102 // Set up local variables.
104 uint64_t iterations = 0;
105
106 // Value iteration loop.
108 while (status == SolverStatus::InProgress && iterations < maximalIterations) {
109 // Compute tmp = A * x + b
110 storm::dd::Add<DdType, ValueType> localXAsColumn = localX.swapVariables(this->rowColumnMetaVariablePairs);
111 storm::dd::Add<DdType, ValueType> tmp = this->A.multiplyMatrix(localXAsColumn, this->columnMetaVariables);
112 tmp += b;
113
114 if (dir == storm::solver::OptimizationDirection::Minimize) {
115 tmp += illegalMaskAdd;
116 tmp = tmp.minAbstract(this->choiceVariables);
117 } else {
118 tmp = tmp.maxAbstract(this->choiceVariables);
119 }
120
121 // Now check if the process already converged within our precision.
122 if (localX.equalModuloPrecision(tmp, precision, relativeTerminationCriterion)) {
124 }
125
126 // Set up next iteration.
127 localX = tmp;
128 ++iterations;
130 status = SolverStatus::Aborted;
131 }
132 }
133
134 if (status == SolverStatus::InProgress && iterations < maximalIterations) {
136 }
137
138 return SymbolicMinMaxLinearEquationSolver<DdType, ValueType>::ValueIterationResult(status, iterations, localX);
139}
140
141template<storm::dd::DdType DdType, typename ValueType>
143 storm::dd::Add<DdType, ValueType> const& b) const {
144 storm::dd::Add<DdType, ValueType> xAsColumn = x.swapVariables(this->rowColumnMetaVariablePairs);
145 storm::dd::Add<DdType, ValueType> tmp = this->A.multiplyMatrix(xAsColumn, this->columnMetaVariables);
146 tmp += b;
147
148 if (dir == storm::solver::OptimizationDirection::Minimize) {
149 tmp += illegalMaskAdd;
150 tmp = tmp.minAbstract(this->choiceVariables);
151 } else {
152 tmp = tmp.maxAbstract(this->choiceVariables);
153 }
154
155 return x == tmp;
156}
157
158template<storm::dd::DdType DdType, typename ValueType>
159template<typename RationalType, typename ImpreciseType>
161 OptimizationDirection dir, uint64_t precision, SymbolicMinMaxLinearEquationSolver<DdType, RationalType> const& rationalSolver,
162 storm::dd::Add<DdType, ImpreciseType> const& x, storm::dd::Add<DdType, RationalType> const& rationalB, bool& isSolution) {
164
165 for (uint64_t p = 1; p < precision; ++p) {
166 sharpenedX = x.sharpenKwekMehlhorn(p);
167 isSolution = rationalSolver.isSolution(dir, sharpenedX, rationalB);
168
169 if (isSolution) {
170 break;
171 }
172 }
173
174 return sharpenedX;
175}
176
177template<storm::dd::DdType DdType, typename ValueType>
178template<typename RationalType, typename ImpreciseType>
179storm::dd::Add<DdType, RationalType> SymbolicMinMaxLinearEquationSolver<DdType, ValueType>::solveEquationsRationalSearchHelper(
180 Environment const& env, OptimizationDirection dir, SymbolicMinMaxLinearEquationSolver<DdType, RationalType> const& rationalSolver,
181 SymbolicMinMaxLinearEquationSolver<DdType, ImpreciseType> const& impreciseSolver, storm::dd::Add<DdType, RationalType> const& rationalB,
183 // Storage for the rational sharpened vector and the power iteration intermediate vector.
186
187 // The actual rational search.
188 uint64_t overallIterations = 0;
189 uint64_t valueIterationInvocations = 0;
190 ValueType precision = storm::utility::convertNumber<ValueType>(env.solver().minMax().getPrecision());
191 uint64_t maxIter = env.solver().minMax().getMaximalNumberOfIterations();
192 bool relative = env.solver().minMax().getRelativeTerminationCriterion();
194 while (status == SolverStatus::InProgress && overallIterations < maxIter) {
195 typename SymbolicMinMaxLinearEquationSolver<DdType, ImpreciseType>::ValueIterationResult viResult =
196 impreciseSolver.performValueIteration(dir, currentX, b, storm::utility::convertNumber<ImpreciseType, ValueType>(precision), relative, maxIter);
197
198 ++valueIterationInvocations;
199 STORM_LOG_TRACE("Completed " << valueIterationInvocations << " value iteration invocations, the last one with precision " << precision
200 << " completed in " << viResult.iterations << " iterations.");
201
202 // Count the iterations.
203 overallIterations += viResult.iterations;
204
205 // Compute maximal precision until which to sharpen.
206 uint64_t p =
207 storm::utility::convertNumber<uint64_t>(storm::utility::ceil(storm::utility::log10<ValueType>(storm::utility::one<ValueType>() / precision)));
208
209 bool isSolution = false;
210 sharpenedX = sharpen<RationalType, ImpreciseType>(dir, p, rationalSolver, viResult.values, rationalB, isSolution);
211
212 if (isSolution) {
214 } else {
215 currentX = viResult.values;
216 precision /= storm::utility::convertNumber<ValueType, uint64_t>(10);
217 }
219 status = SolverStatus::Aborted;
220 }
221 }
222
223 if (status == SolverStatus::InProgress && overallIterations < maxIter) {
225 }
226
227 if (status == SolverStatus::Converged) {
228 STORM_LOG_INFO("Iterative solver (rational search) converged in " << overallIterations << " iterations.");
229 } else {
230 STORM_LOG_WARN("Iterative solver (rational search) did not converge in " << overallIterations << " iterations.");
231 }
232
233 return sharpenedX;
234}
235
236template<storm::dd::DdType DdType, typename ValueType>
237template<typename ImpreciseType>
238typename std::enable_if<std::is_same<ValueType, ImpreciseType>::value && storm::NumberTraits<ValueType>::IsExact, storm::dd::Add<DdType, ValueType>>::type
239SymbolicMinMaxLinearEquationSolver<DdType, ValueType>::solveEquationsRationalSearchHelper(Environment const& env,
242 storm::dd::Add<DdType, ValueType> const& b) const {
243 return solveEquationsRationalSearchHelper<ValueType, ValueType>(env, dir, *this, *this, b, this->getLowerBoundsVector(), b);
244}
245
246template<storm::dd::DdType DdType, typename ValueType>
247template<typename ImpreciseType>
248typename std::enable_if<std::is_same<ValueType, ImpreciseType>::value && !storm::NumberTraits<ValueType>::IsExact, storm::dd::Add<DdType, ValueType>>::type
249SymbolicMinMaxLinearEquationSolver<DdType, ValueType>::solveEquationsRationalSearchHelper(Environment const& env,
252 storm::dd::Add<DdType, ValueType> const& b) const {
253 storm::dd::Add<DdType, storm::RationalNumber> rationalB = b.template toValueType<storm::RationalNumber>();
254 SymbolicMinMaxLinearEquationSolver<DdType, storm::RationalNumber> rationalSolver(
255 this->A.template toValueType<storm::RationalNumber>(), this->allRows, this->illegalMask, this->rowMetaVariables, this->columnMetaVariables,
256 this->choiceVariables, this->rowColumnMetaVariablePairs, std::make_unique<GeneralSymbolicLinearEquationSolverFactory<DdType, storm::RationalNumber>>());
257
259 solveEquationsRationalSearchHelper<storm::RationalNumber, ImpreciseType>(env, dir, rationalSolver, *this, rationalB, this->getLowerBoundsVector(), b);
260 return rationalResult.template toValueType<ValueType>();
261}
262
263template<storm::dd::DdType DdType, typename ValueType>
264template<typename ImpreciseType>
265typename std::enable_if<!std::is_same<ValueType, ImpreciseType>::value, storm::dd::Add<DdType, ValueType>>::type
266SymbolicMinMaxLinearEquationSolver<DdType, ValueType>::solveEquationsRationalSearchHelper(Environment const& env,
269 storm::dd::Add<DdType, ValueType> const& b) const {
270 // First try to find a solution using the imprecise value type.
273 try {
274 impreciseX = this->getLowerBoundsVector().template toValueType<ImpreciseType>();
275 storm::dd::Add<DdType, ImpreciseType> impreciseB = b.template toValueType<ImpreciseType>();
276 SymbolicMinMaxLinearEquationSolver<DdType, ImpreciseType> impreciseSolver(
277 this->A.template toValueType<ImpreciseType>(), this->allRows, this->illegalMask, this->rowMetaVariables, this->columnMetaVariables,
278 this->choiceVariables, this->rowColumnMetaVariablePairs, std::make_unique<GeneralSymbolicLinearEquationSolverFactory<DdType, ImpreciseType>>());
279
280 rationalResult = solveEquationsRationalSearchHelper<ValueType, ImpreciseType>(env, dir, *this, impreciseSolver, b, impreciseX, impreciseB);
281 } catch (storm::exceptions::PrecisionExceededException const& e) {
282 STORM_LOG_WARN("Precision of value type was exceeded, trying to recover by switching to rational arithmetic.");
283
284 // Fall back to precise value type if the precision of the imprecise value type was exceeded.
285 rationalResult = solveEquationsRationalSearchHelper<ValueType, ValueType>(env, dir, *this, *this, b, impreciseX.template toValueType<ValueType>(), b);
286 }
287 return rationalResult.template toValueType<ValueType>();
288}
289
290template<storm::dd::DdType DdType, typename ValueType>
291storm::dd::Add<DdType, ValueType> SymbolicMinMaxLinearEquationSolver<DdType, ValueType>::solveEquationsRationalSearch(
292 Environment const& env, storm::solver::OptimizationDirection const& dir, storm::dd::Add<DdType, ValueType> const& x,
293 storm::dd::Add<DdType, ValueType> const& b) const {
294 return solveEquationsRationalSearchHelper<double>(env, dir, x, b);
295}
296
297template<storm::dd::DdType DdType, typename ValueType>
298storm::dd::Add<DdType, ValueType> SymbolicMinMaxLinearEquationSolver<DdType, ValueType>::solveEquationsValueIteration(
299 Environment const& env, storm::solver::OptimizationDirection const& dir, storm::dd::Add<DdType, ValueType> const& x,
300 storm::dd::Add<DdType, ValueType> const& b) const {
301 // Set up the environment.
303
304 if (this->hasUniqueSolution()) {
305 localX = x;
306 } else {
307 // If we were given an initial scheduler, we take its solution as the starting point.
308 if (this->hasInitialScheduler()) {
309 // The linear equation solver should be at least as precise as this solver
310 std::unique_ptr<storm::Environment> environmentOfSolverStorage;
311 auto precOfSolver = env.solver().getPrecisionOfLinearEquationSolver(env.solver().getLinearEquationSolverType());
313 bool changePrecision = precOfSolver.first && precOfSolver.first.get() > env.solver().minMax().getPrecision();
314 bool changeRelative = precOfSolver.second && !precOfSolver.second.get() && env.solver().minMax().getRelativeTerminationCriterion();
315 if (changePrecision || changeRelative) {
316 environmentOfSolverStorage = std::make_unique<storm::Environment>(env);
317 boost::optional<storm::RationalNumber> newPrecision;
318 boost::optional<bool> newRelative;
319 if (changePrecision) {
320 newPrecision = env.solver().minMax().getPrecision();
321 }
322 if (changeRelative) {
323 newRelative = true;
324 }
325 environmentOfSolverStorage->solver().setLinearEquationSolverPrecision(newPrecision, newRelative);
326 }
327 }
328 storm::Environment const& environmentOfSolver = environmentOfSolverStorage ? *environmentOfSolverStorage : env;
329 localX = solveEquationsWithScheduler(environmentOfSolver, this->getInitialScheduler(), x, b);
330 } else {
331 localX = this->getLowerBoundsVector();
332 }
333 }
334
335 ValueType precision = storm::utility::convertNumber<ValueType>(env.solver().minMax().getPrecision());
336 ValueIterationResult viResult = performValueIteration(dir, localX, b, precision, env.solver().minMax().getRelativeTerminationCriterion(),
337 env.solver().minMax().getMaximalNumberOfIterations());
338
339 if (viResult.status == SolverStatus::Converged) {
340 STORM_LOG_INFO("Iterative solver (value iteration) converged in " << viResult.iterations << " iterations.");
341 } else {
342 STORM_LOG_WARN("Iterative solver (value iteration) did not converge in " << viResult.iterations << " iterations.");
343 }
344
345 return viResult.values;
346}
347
348template<storm::dd::DdType DdType, typename ValueType>
349storm::dd::Add<DdType, ValueType> SymbolicMinMaxLinearEquationSolver<DdType, ValueType>::solveEquationsWithScheduler(
350 Environment const& env, storm::dd::Bdd<DdType> const& scheduler, storm::dd::Add<DdType, ValueType> const& x,
351 storm::dd::Add<DdType, ValueType> const& b) const {
352 std::unique_ptr<SymbolicLinearEquationSolver<DdType, ValueType>> solver =
353 linearEquationSolverFactory->create(env, this->allRows, this->rowMetaVariables, this->columnMetaVariables, this->rowColumnMetaVariablePairs);
354 this->forwardBounds(*solver);
355
357 (storm::utility::dd::getRowColumnDiagonal<DdType>(x.getDdManager(), this->rowColumnMetaVariablePairs) && this->allRows).template toAdd<ValueType>();
358 return solveEquationsWithScheduler(env, *solver, scheduler, x, b, diagonal);
359}
360
361template<storm::dd::DdType DdType, typename ValueType>
362storm::dd::Add<DdType, ValueType> SymbolicMinMaxLinearEquationSolver<DdType, ValueType>::solveEquationsWithScheduler(
363 Environment const& env, SymbolicLinearEquationSolver<DdType, ValueType>& solver, storm::dd::Bdd<DdType> const& scheduler,
365 // Apply scheduler to the matrix and vector.
367 scheduler.ite(this->A, scheduler.getDdManager().template getAddZero<ValueType>()).sumAbstract(this->choiceVariables);
368 if (solver.getEquationProblemFormat(env) == storm::solver::LinearEquationSolverProblemFormat::EquationSystem) {
369 schedulerA = diagonal - schedulerA;
370 }
371
373 scheduler.ite(b, scheduler.getDdManager().template getAddZero<ValueType>()).sumAbstract(this->choiceVariables);
374
375 // Set the matrix for the solver.
376 solver.setMatrix(schedulerA);
377
378 // Solve for the value of the scheduler.
379 storm::dd::Add<DdType, ValueType> schedulerX = solver.solveEquations(env, x, schedulerB);
380
381 return schedulerX;
382}
383
384template<storm::dd::DdType DdType, typename ValueType>
385storm::dd::Add<DdType, ValueType> SymbolicMinMaxLinearEquationSolver<DdType, ValueType>::solveEquationsPolicyIteration(
386 Environment const& env, storm::solver::OptimizationDirection const& dir, storm::dd::Add<DdType, ValueType> const& x,
387 storm::dd::Add<DdType, ValueType> const& b) const {
388 // Set up the environment.
389 storm::dd::Add<DdType, ValueType> currentSolution = x;
391 (storm::utility::dd::getRowColumnDiagonal<DdType>(x.getDdManager(), this->rowColumnMetaVariablePairs) && this->allRows).template toAdd<ValueType>();
392 uint_fast64_t iterations = 0;
393 bool converged = false;
394
395 // Choose initial scheduler.
396 storm::dd::Bdd<DdType> scheduler =
397 this->hasInitialScheduler()
398 ? this->getInitialScheduler()
399 : (this->A.sumAbstract(this->columnMetaVariables).notZero() || b.notZero()).existsAbstractRepresentative(this->choiceVariables);
400
401 // Initialize linear equation solver.
402 // It should be at least as precise as this solver.
403 std::unique_ptr<storm::Environment> environmentOfSolverStorage;
404 auto precOfSolver = env.solver().getPrecisionOfLinearEquationSolver(env.solver().getLinearEquationSolverType());
406 bool changePrecision = precOfSolver.first && precOfSolver.first.get() > env.solver().minMax().getPrecision();
407 bool changeRelative = precOfSolver.second && !precOfSolver.second.get() && env.solver().minMax().getRelativeTerminationCriterion();
408 if (changePrecision || changeRelative) {
409 environmentOfSolverStorage = std::make_unique<storm::Environment>(env);
410 boost::optional<storm::RationalNumber> newPrecision;
411 boost::optional<bool> newRelative;
412 if (changePrecision) {
413 newPrecision = env.solver().minMax().getPrecision();
414 }
415 if (changeRelative) {
416 newRelative = true;
417 }
418 environmentOfSolverStorage->solver().setLinearEquationSolverPrecision(newPrecision, newRelative);
419 }
420 }
421 storm::Environment const& environmentOfSolver = environmentOfSolverStorage ? *environmentOfSolverStorage : env;
422
423 std::unique_ptr<SymbolicLinearEquationSolver<DdType, ValueType>> linearEquationSolver = linearEquationSolverFactory->create(
424 environmentOfSolver, this->allRows, this->rowMetaVariables, this->columnMetaVariables, this->rowColumnMetaVariablePairs);
425 this->forwardBounds(*linearEquationSolver);
426
427 // Iteratively solve and improve the scheduler.
428 uint64_t maxIter = env.solver().minMax().getMaximalNumberOfIterations();
429 while (!converged && iterations < maxIter) {
431 solveEquationsWithScheduler(environmentOfSolver, *linearEquationSolver, scheduler, currentSolution, b, diagonal);
432
433 // Policy improvement step.
435 this->A.multiplyMatrix(schedulerX.swapVariables(this->rowColumnMetaVariablePairs), this->columnMetaVariables) + b;
436
437 storm::dd::Bdd<DdType> nextScheduler;
438 if (dir == storm::solver::OptimizationDirection::Minimize) {
439 choiceValues += illegalMaskAdd;
440
441 storm::dd::Add<DdType, ValueType> newStateValues = choiceValues.minAbstract(this->choiceVariables);
442 storm::dd::Bdd<DdType> improvedStates = newStateValues.less(schedulerX);
443
444 nextScheduler = improvedStates.ite(choiceValues.minAbstractRepresentative(this->choiceVariables), scheduler);
445 } else {
446 storm::dd::Add<DdType, ValueType> newStateValues = choiceValues.maxAbstract(this->choiceVariables);
447 storm::dd::Bdd<DdType> improvedStates = newStateValues.greater(schedulerX);
448
449 nextScheduler = improvedStates.ite(choiceValues.maxAbstractRepresentative(this->choiceVariables), scheduler);
450 }
451
452 // Check for convergence.
453 converged = nextScheduler == scheduler;
454
455 // Set up next iteration.
456 if (!converged) {
457 scheduler = nextScheduler;
458 }
459
460 currentSolution = schedulerX;
461 ++iterations;
462 }
463
464 if (converged) {
465 STORM_LOG_INFO("Iterative solver (policy iteration) converged in " << iterations << " iterations.");
466 } else {
467 STORM_LOG_WARN("Iterative solver (policy iteration) did not converge in " << iterations << " iterations.");
468 }
469
470 return currentSolution;
471}
472
473template<storm::dd::DdType DdType, typename ValueType>
477 uint_fast64_t n) const {
479
480 // Perform matrix-vector multiplication while the bound is met.
481 for (uint_fast64_t i = 0; i < n; ++i) {
482 xCopy = xCopy.swapVariables(this->rowColumnMetaVariablePairs);
483 xCopy = this->A.multiplyMatrix(xCopy, this->columnMetaVariables);
484 if (b != nullptr) {
485 xCopy += *b;
486 }
487
488 if (dir == storm::solver::OptimizationDirection::Minimize) {
489 // This is a hack and only here because of the lack of a suitable minAbstract/maxAbstract function
490 // that can properly deal with a restriction of the choices.
491 xCopy += illegalMaskAdd;
492 xCopy = xCopy.minAbstract(this->choiceVariables);
493 } else {
494 xCopy = xCopy.maxAbstract(this->choiceVariables);
495 }
496 }
497
498 return xCopy;
499}
500
501template<storm::dd::DdType DdType, typename ValueType>
503 this->initialScheduler = scheduler;
504}
505
506template<storm::dd::DdType DdType, typename ValueType>
510
511template<storm::dd::DdType DdType, typename ValueType>
513 return static_cast<bool>(initialScheduler);
514}
515
516template<storm::dd::DdType DdType, typename ValueType>
518 Environment const& env, boost::optional<storm::solver::OptimizationDirection> const& direction) const {
520
521 auto method = getMethod(env, std::is_same<ValueType, storm::RationalNumber>::value);
522 if (method == MinMaxMethod::PolicyIteration) {
523 if (!this->hasUniqueSolution()) {
524 requirements.requireValidInitialScheduler();
525 }
526 } else if (method == MinMaxMethod::ValueIteration) {
527 if (!this->hasUniqueSolution()) {
528 if (!direction || direction.get() == storm::solver::OptimizationDirection::Maximize) {
529 requirements.requireLowerBounds();
530 }
531 if (!direction || direction.get() == storm::solver::OptimizationDirection::Minimize) {
532 requirements.requireValidInitialScheduler();
533 }
534 }
535 } else if (method == MinMaxMethod::RationalSearch) {
536 requirements.requireLowerBounds();
537 if (!this->hasUniqueSolution() && (!direction || direction.get() == storm::solver::OptimizationDirection::Minimize)) {
538 requirements.requireUniqueSolution();
539 }
540 } else {
541 STORM_LOG_THROW(false, storm::exceptions::InvalidEnvironmentException, "The selected min max technique is not supported by this solver.");
542 }
543
544 return requirements;
545}
546
547template<storm::dd::DdType DdType, typename ValueType>
549 this->uniqueSolution = value;
550}
551
552template<storm::dd::DdType DdType, typename ValueType>
554 return this->uniqueSolution;
555}
556
557template<storm::dd::DdType DdType, typename ValueType>
559 this->requirementsChecked = value;
560}
561
562template<storm::dd::DdType DdType, typename ValueType>
564 return this->requirementsChecked;
565}
566
567template<storm::dd::DdType DdType, typename ValueType>
569 if (this->hasLowerBound()) {
570 solver.setLowerBound(this->getLowerBound());
571 }
572 if (this->hasLowerBounds()) {
573 solver.setLowerBounds(this->getLowerBounds());
574 }
575 if (this->hasUpperBound()) {
576 solver.setUpperBound(this->getUpperBound());
577 }
578 if (this->hasUpperBounds()) {
579 solver.setUpperBounds(this->getUpperBounds());
580 }
581}
582
583template<storm::dd::DdType DdType, typename ValueType>
585 Environment const& env, bool hasUniqueSolution, boost::optional<storm::solver::OptimizationDirection> const& direction) const {
586 std::unique_ptr<storm::solver::SymbolicMinMaxLinearEquationSolver<DdType, ValueType>> solver = this->create();
587 solver->setHasUniqueSolution(hasUniqueSolution);
588 return solver->getRequirements(env, direction);
589}
590
591template<storm::dd::DdType DdType, typename ValueType>
592std::unique_ptr<storm::solver::SymbolicMinMaxLinearEquationSolver<DdType, ValueType>>
594 storm::dd::Add<DdType, ValueType> const& A, storm::dd::Bdd<DdType> const& allRows, storm::dd::Bdd<DdType> const& illegalMask,
595 std::set<storm::expressions::Variable> const& rowMetaVariables, std::set<storm::expressions::Variable> const& columnMetaVariables,
596 std::set<storm::expressions::Variable> const& choiceVariables,
597 std::vector<std::pair<storm::expressions::Variable, storm::expressions::Variable>> const& rowColumnMetaVariablePairs) const {
598 return std::make_unique<SymbolicMinMaxLinearEquationSolver<DdType, ValueType>>(
599 A, allRows, illegalMask, rowMetaVariables, columnMetaVariables, choiceVariables, rowColumnMetaVariablePairs,
600 std::make_unique<GeneralSymbolicLinearEquationSolverFactory<DdType, ValueType>>());
601}
602
603template<storm::dd::DdType DdType, typename ValueType>
604std::unique_ptr<storm::solver::SymbolicMinMaxLinearEquationSolver<DdType, ValueType>>
606 return std::make_unique<SymbolicMinMaxLinearEquationSolver<DdType, ValueType>>();
607}
608
609template class SymbolicMinMaxLinearEquationSolver<storm::dd::DdType::CUDD, double>;
610template class SymbolicMinMaxLinearEquationSolver<storm::dd::DdType::Sylvan, double>;
611template class SymbolicMinMaxLinearEquationSolver<storm::dd::DdType::Sylvan, storm::RationalNumber>;
612
613template class GeneralSymbolicMinMaxLinearEquationSolverFactory<storm::dd::DdType::CUDD, double>;
614template class GeneralSymbolicMinMaxLinearEquationSolverFactory<storm::dd::DdType::Sylvan, double>;
615template class GeneralSymbolicMinMaxLinearEquationSolverFactory<storm::dd::DdType::Sylvan, storm::RationalNumber>;
616
617} // namespace solver
618} // namespace storm
SolverEnvironment & solver()
storm::solver::MinMaxMethod const & getMethod() const
MinMaxSolverEnvironment & minMax()
Add< LibraryType, ValueType > swapVariables(std::vector< std::pair< storm::expressions::Variable, storm::expressions::Variable > > const &metaVariablePairs) const
Swaps the given pairs of meta variables in the ADD.
Definition Add.cpp:292
Add< LibraryType, storm::RationalNumber > sharpenKwekMehlhorn(uint64_t precision) const
Retrieves the function that sharpens all values in the current ADD with the Kwek-Mehlhorn algorithm.
Definition Add.cpp:151
Bdd< LibraryType > greater(Add< LibraryType, ValueType > const &other) const
Retrieves the function that maps all evaluations to one whose function value in the first ADD are gre...
Definition Add.cpp:116
Add< LibraryType, ValueType > maxAbstract(std::set< storm::expressions::Variable > const &metaVariables) const
Max-abstracts from the given meta variables.
Definition Add.cpp:198
Bdd< LibraryType > minAbstractRepresentative(std::set< storm::expressions::Variable > const &metaVariables) const
Similar to minAbstract, but does not abstract from the variables but rather picks a valuation of each...
Definition Add.cpp:192
Add< LibraryType, ValueType > minAbstract(std::set< storm::expressions::Variable > const &metaVariables) const
Min-abstracts from the given meta variables.
Definition Add.cpp:185
Add< LibraryType, ValueType > multiplyMatrix(Add< LibraryType, ValueType > const &otherMatrix, std::set< storm::expressions::Variable > const &summationMetaVariables) const
Multiplies the current ADD (representing a matrix) with the given matrix by summing over the given me...
Definition Add.cpp:372
bool equalModuloPrecision(Add< LibraryType, ValueType > const &other, ValueType const &precision, bool relative=true) const
Checks whether the current and the given ADD represent the same function modulo some given precision.
Definition Add.cpp:211
Bdd< LibraryType > maxAbstractRepresentative(std::set< storm::expressions::Variable > const &metaVariables) const
Similar to maxAbstract, but does not abstract from the variables but rather picks a valuation of each...
Definition Add.cpp:205
Bdd< LibraryType > less(Add< LibraryType, ValueType > const &other) const
Retrieves the function that maps all evaluations to one whose function value in the first ADD are les...
Definition Add.cpp:106
Bdd< LibraryType > ite(Bdd< LibraryType > const &thenBdd, Bdd< LibraryType > const &elseBdd) const
Performs an if-then-else with the given operands, i.e.
Definition Bdd.cpp:113
DdManager< LibraryType > & getDdManager() const
Retrieves the manager that is responsible for this DD.
Definition Dd.cpp:39
MinMaxLinearEquationSolverRequirements & requireUniqueSolution(bool critical=true)
MinMaxLinearEquationSolverRequirements & requireLowerBounds(bool critical=true)
MinMaxLinearEquationSolverRequirements & requireValidInitialScheduler(bool critical=true)
virtual void setUpperBound(ValueType const &lowerBound)
virtual void setLowerBounds(storm::dd::Add< DdType, ValueType > const &lowerBounds)
virtual void setUpperBounds(storm::dd::Add< DdType, ValueType > const &upperBounds)
virtual void setLowerBound(ValueType const &lowerBound)
An interface that represents an abstract symbolic linear equation solver.
MinMaxLinearEquationSolverRequirements getRequirements(Environment const &env, bool hasUniqueSolution=false, boost::optional< storm::solver::OptimizationDirection > const &direction=boost::none) const
Retrieves the requirements of the solver that would be created when calling create() right now.
An interface that represents an abstract symbolic linear equation solver.
bool hasUniqueSolution() const
Retrieves whether the solution to the min max equation system is assumed to be unique.
bool hasInitialScheduler() const
Retrieves whether an initial scheduler was set.
void setHasUniqueSolution(bool value=true)
Sets whether the solution to the min max equation system is known to be unique.
virtual storm::dd::Add< DdType, ValueType > multiply(storm::solver::OptimizationDirection const &dir, storm::dd::Add< DdType, ValueType > const &x, storm::dd::Add< DdType, ValueType > const *b=nullptr, uint_fast64_t n=1) const
Performs repeated matrix-vector multiplication, using x[0] = x and x[i + 1] = A*x[i] + b.
virtual storm::dd::Add< DdType, ValueType > solveEquations(Environment const &env, storm::solver::OptimizationDirection const &dir, storm::dd::Add< DdType, ValueType > const &x, storm::dd::Add< DdType, ValueType > const &b) const
Solves the equation system A*x = b.
bool isSolution(OptimizationDirection dir, storm::dd::Add< DdType, ValueType > const &x, storm::dd::Add< DdType, ValueType > const &b) const
Determines whether the given vector x satisfies x = min/max Ax + b.
void setInitialScheduler(storm::dd::Bdd< DdType > const &scheduler)
Sets an initial scheduler that is required by some solvers (see requirements).
bool isRequirementsCheckedSet() const
Retrieves whether the solver is aware that the requirements were checked.
storm::dd::Bdd< DdType > const & getInitialScheduler() const
Retrieves the initial scheduler (if there is any).
void setRequirementsChecked(bool value=true)
Notifies the solver that the requirements for solving equations have been checked.
virtual MinMaxLinearEquationSolverRequirements getRequirements(Environment const &env, boost::optional< storm::solver::OptimizationDirection > const &direction=boost::none) const
Retrieves the requirements of the solver.
#define STORM_LOG_INFO(message)
Definition logging.h:29
#define STORM_LOG_WARN(message)
Definition logging.h:30
#define STORM_LOG_TRACE(message)
Definition logging.h:17
#define STORM_LOG_THROW(cond, exception, message)
Definition macros.h:30
#define STORM_LOG_WARN_COND_DEBUG(cond, message)
Definition macros.h:18
SFTBDDChecker::ValueType ValueType
bool isTerminate()
Check whether the program should terminate (due to some abort signal).
ValueType ceil(ValueType const &number)
LabParser.cpp.
Definition cli.cpp:18