From cb849a9ab832db14fc8f22a6764b7ebe06fce05b Mon Sep 17 00:00:00 2001 From: dehnert Date: Sun, 17 Sep 2017 22:49:21 +0200 Subject: [PATCH] started on computing upper bounds for rewards for interval value iteration --- .../prctl/helper/SparseDtmcPrctlHelper.cpp | 143 +++++++++++++++++- .../prctl/helper/SparseMdpPrctlHelper.cpp | 4 +- .../solver/EigenLinearEquationSolver.cpp | 2 +- .../IterativeMinMaxLinearEquationSolver.cpp | 4 +- src/storm/solver/LinearEquationSolver.cpp | 86 +++++++++-- src/storm/solver/LinearEquationSolver.h | 58 ++++++- .../LinearEquationSolverRequirements.cpp | 32 ++-- .../solver/LinearEquationSolverRequirements.h | 20 +-- ...MinMaxLinearEquationSolverRequirements.cpp | 37 +++-- .../MinMaxLinearEquationSolverRequirements.h | 21 +-- .../solver/NativeLinearEquationSolver.cpp | 8 +- src/storm/utility/vector.h | 30 ++++ .../EliminationLinearEquationSolverTest.cpp | 1 + 13 files changed, 366 insertions(+), 80 deletions(-) diff --git a/src/storm/modelchecker/prctl/helper/SparseDtmcPrctlHelper.cpp b/src/storm/modelchecker/prctl/helper/SparseDtmcPrctlHelper.cpp index fe1f120de..90ed634db 100644 --- a/src/storm/modelchecker/prctl/helper/SparseDtmcPrctlHelper.cpp +++ b/src/storm/modelchecker/prctl/helper/SparseDtmcPrctlHelper.cpp @@ -6,6 +6,8 @@ #include "storm/utility/vector.h" #include "storm/utility/graph.h" +#include "storm/storage/StronglyConnectedComponentDecomposition.h" + #include "storm/solver/LinearEquationSolver.h" #include "storm/modelchecker/results/ExplicitQuantitativeCheckResult.h" @@ -14,6 +16,7 @@ #include "storm/exceptions/InvalidStateException.h" #include "storm/exceptions/InvalidPropertyException.h" #include "storm/exceptions/IllegalArgumentException.h" +#include "storm/exceptions/UncheckedRequirementException.h" namespace storm { namespace modelchecker { @@ -211,6 +214,130 @@ namespace storm { targetStates, qualitative, linearEquationSolverFactory, hint); } + template + class DsMpi { + public: + DsMpi(storm::storage::SparseMatrix const& transitionMatrix, std::vector const& rewards, std::vector const& oneStepTargetProbabilities) : transitionMatrix(transitionMatrix), originalRewards(rewards), backwardTransitions(transitionMatrix.transpose()), p(transitionMatrix.getRowCount()), w(transitionMatrix.getRowCount()), rewards(rewards), targetProbabilities(oneStepTargetProbabilities) { + // Intentionally left empty. + } + + std::vector computeUpperBounds() { + sweep(); + ValueType lambda = computeLambda(); + + // Finally compute the upper bounds for the states. + std::vector result(transitionMatrix.getRowCount()); + auto one = storm::utility::one(); + for (storm::storage::sparse::state_type state = 0; state < result.size(); ++state) { + result[state] = w[state] + (one - p[state]) * lambda; + } + return result; + } + + private: + ValueType computeLambda() { + ValueType lambda = storm::utility::convertNumber(0.0); + for (storm::storage::sparse::state_type state = 0; state < targetProbabilities.size(); ++state) { + // Check whether condition (I) or (II) applies. + ValueType sum = storm::utility::zero(); + for (auto const& e : transitionMatrix.getRow(state)) { + sum += e.getValue() * p[e.getColumn()]; + } + if (p[state] < sum) { + // Condition (I) applies. + ValueType localLambda = sum - p[state]; + ValueType nominator = originalRewards[state]; + for (auto const& e : transitionMatrix.getRow(state)) { + nominator += e.getValue() * w[e.getColumn()]; + } + nominator -= w[state]; + localLambda = nominator / localLambda; + lambda = std::max(lambda, localLambda); + } else { + // Here, condition (II) automatically applies and as the resulting local lambda is 0, we + // don't need to consider it. + } + } + return lambda; + } + + void sweep() { + // Create a priority queue that allows for easy retrieval of the currently best state. + auto cmp = [this](storm::storage::sparse::state_type const& a, storm::storage::sparse::state_type const& b) { + ValueType pa = p[a]; + ValueType pb = p[b]; + if (pa > pb) { + return true; + } else if (pa == pb) { + return w[a] < w[b]; + } + return false; + }; + std::set queue; + + + storm::storage::BitVector visited(p.size()); + + for (storm::storage::sparse::state_type state = 0; state < targetProbabilities.size(); ++state) { + if (!storm::utility::isZero(targetProbabilities[state])) { + queue.insert(state); + } + } + + while (!queue.empty()) { + // Get first entry in queue. + storm::storage::sparse::state_type currentState = *queue.begin(); + queue.erase(queue.begin()); + + // Mark state as visited. + visited.set(currentState); + + // Set weight and probability for the state. + w[currentState] = rewards[currentState]; + p[currentState] = targetProbabilities[currentState]; + + for (auto const& e : backwardTransitions.getRow(currentState)) { + if (visited.get(e.getColumn())) { + continue; + } + + // Delete element from the priority queue if it was in there. + auto it = queue.find(e.getColumn()); + if (it != queue.end()) { + queue.erase(it); + } + + // Update reward/probability values. + rewards[e.getColumn()] += e.getValue() * w[currentState]; + targetProbabilities[e.getColumn()] += e.getValue() * p[currentState]; + + // (Re-)insert the state with the new rewards/target probabilities. + queue.insert(e.getColumn()); + } + } + } + + // References to input data. + storm::storage::SparseMatrix const& transitionMatrix; + std::vector const& originalRewards; + + // Derived from input data. + storm::storage::SparseMatrix backwardTransitions; + + // Data that the algorithm uses internally. + std::vector p; + std::vector w; + std::vector rewards; + std::vector targetProbabilities; + }; + + // This function computes an upper bound on the reachability rewards (see Baier et al, CAV'17). + template + std::vector computeUpperRewardBounds(storm::storage::SparseMatrix const& transitionMatrix, std::vector const& rewards, std::vector const& oneStepTargetProbabilities) { + DsMpi dsmpi(transitionMatrix, rewards, oneStepTargetProbabilities); + return dsmpi.computeUpperBounds(); + } + template std::vector SparseDtmcPrctlHelper::computeReachabilityRewards(storm::storage::SparseMatrix const& transitionMatrix, storm::storage::SparseMatrix const& backwardTransitions, std::function(uint_fast64_t, storm::storage::SparseMatrix const&, storm::storage::BitVector const&)> const& totalStateRewardVectorGetter, storm::storage::BitVector const& targetStates, bool qualitative, storm::solver::LinearEquationSolverFactory const& linearEquationSolverFactory, ModelCheckerHint const& hint) { @@ -264,10 +391,24 @@ namespace storm { // Prepare the right-hand side of the equation system. std::vector b = totalStateRewardVectorGetter(submatrix.getRowCount(), transitionMatrix, maybeStates); + + storm::solver::LinearEquationSolverRequirements requirements = linearEquationSolverFactory.getRequirements(); + boost::optional> upperRewardBounds; + requirements.clearLowerBounds(); + if (requirements.requiresUpperBounds()) { + upperRewardBounds = computeUpperRewardBounds(submatrix, b, transitionMatrix.getConstrainedRowSumVector(maybeStates, targetStates)); + requirements.clearUpperBounds(); + } + STORM_LOG_THROW(requirements.empty(), storm::exceptions::UncheckedRequirementException, "There are unchecked requirements of the solver."); - // Now solve the resulting equation system. + // Create the solvers and provide std::unique_ptr> solver = linearEquationSolverFactory.create(std::move(submatrix)); solver->setLowerBound(storm::utility::zero()); + if (upperRewardBounds) { + solver->setUpperBounds(std::move(upperRewardBounds.get())); + } + + // Now solve the resulting equation system. solver->solveEquations(x, b); // Set values of resulting vector according to result. diff --git a/src/storm/modelchecker/prctl/helper/SparseMdpPrctlHelper.cpp b/src/storm/modelchecker/prctl/helper/SparseMdpPrctlHelper.cpp index 14da77296..5ed84eb25 100644 --- a/src/storm/modelchecker/prctl/helper/SparseMdpPrctlHelper.cpp +++ b/src/storm/modelchecker/prctl/helper/SparseMdpPrctlHelper.cpp @@ -210,7 +210,9 @@ namespace storm { requirements.clearValidInitialScheduler(); } if (type == storm::solver::EquationSystemType::UntilProbabilities) { - requirements.clearGlobalUpperBound(); + requirements.clearBounds(); + } else if (type == storm::solver::EquationSystemType::ReachabilityRewards) { + requirements.clearLowerBounds(); } STORM_LOG_THROW(requirements.empty(), storm::exceptions::UncheckedRequirementException, "There are unchecked requirements of the solver."); } diff --git a/src/storm/solver/EigenLinearEquationSolver.cpp b/src/storm/solver/EigenLinearEquationSolver.cpp index 002352d88..f2df524b3 100644 --- a/src/storm/solver/EigenLinearEquationSolver.cpp +++ b/src/storm/solver/EigenLinearEquationSolver.cpp @@ -276,7 +276,7 @@ namespace storm { } } - // Make sure that all results conform to the bounds. + // Make sure that all results conform to the (global) bounds. storm::utility::vector::clip(x, this->lowerBound, this->upperBound); // Check if the solver converged and issue a warning otherwise. diff --git a/src/storm/solver/IterativeMinMaxLinearEquationSolver.cpp b/src/storm/solver/IterativeMinMaxLinearEquationSolver.cpp index 7acc811de..e9913c3a1 100644 --- a/src/storm/solver/IterativeMinMaxLinearEquationSolver.cpp +++ b/src/storm/solver/IterativeMinMaxLinearEquationSolver.cpp @@ -256,7 +256,7 @@ namespace storm { if (!direction || direction.get() == OptimizationDirection::Maximize) { requirements.requireNoEndComponents(); } - requirements.requireGlobalUpperBound(); + requirements.requireUpperBounds(); } // Then add our requirements on top of that. @@ -273,7 +273,7 @@ namespace storm { } if (this->getSettings().getSolutionMethod() == IterativeMinMaxLinearEquationSolverSettings::SolutionMethod::ValueIteration && this->getSettings().getForceSoundness()) { - requirements.requireGlobalUpperBound(); + requirements.requireUpperBounds(); } return requirements; diff --git a/src/storm/solver/LinearEquationSolver.cpp b/src/storm/solver/LinearEquationSolver.cpp index 41de60629..88e06293c 100644 --- a/src/storm/solver/LinearEquationSolver.cpp +++ b/src/storm/solver/LinearEquationSolver.cpp @@ -131,11 +131,35 @@ namespace storm { cachedRowVector.reset(); } + template + bool LinearEquationSolver::hasLowerBound(BoundType const& type) const { + if (type == BoundType::Any) { + return static_cast(lowerBound) || static_cast(lowerBounds); + } else if (type == BoundType::Global) { + return static_cast(lowerBound); + } else if (type == BoundType::Local) { + return static_cast(lowerBounds); + } + return false; + } + + template + bool LinearEquationSolver::hasUpperBound(BoundType const& type) const { + if (type == BoundType::Any) { + return static_cast(upperBound) || static_cast(upperBounds); + } else if (type == BoundType::Global) { + return static_cast(upperBound); + } else if (type == BoundType::Local) { + return static_cast(upperBounds); + } + return false; + } + template void LinearEquationSolver::setLowerBound(ValueType const& value) { lowerBound = value; } - + template void LinearEquationSolver::setUpperBound(ValueType const& value) { upperBound = value; @@ -148,23 +172,67 @@ namespace storm { } template - bool LinearEquationSolver::hasLowerBound() const { - return static_cast(lowerBound); + ValueType const& LinearEquationSolver::getLowerBound() const { + return lowerBound.get(); } template - bool LinearEquationSolver::hasUpperBound() const { - return static_cast(upperBound); + ValueType const& LinearEquationSolver::getUpperBound() const { + return upperBound.get(); } template - ValueType const& LinearEquationSolver::getLowerBound() const { - return lowerBound.get(); + std::vector const& LinearEquationSolver::getLowerBounds() const { + return lowerBounds.get(); } template - ValueType const& LinearEquationSolver::getUpperBound() const { - return upperBound.get(); + std::vector const& LinearEquationSolver::getUpperBounds() const { + return upperBounds.get(); + } + + template + void LinearEquationSolver::setLowerBounds(std::vector const& values) { + lowerBounds = values; + } + + template + void LinearEquationSolver::setUpperBounds(std::vector const& values) { + upperBounds = values; + } + + template + void LinearEquationSolver::setUpperBounds(std::vector&& values) { + upperBounds = std::move(values); + } + + template + void LinearEquationSolver::setBounds(std::vector const& lower, std::vector const& upper) { + setLowerBounds(lower); + setUpperBounds(upper); + } + + template + void LinearEquationSolver::createUpperBoundsVector(std::unique_ptr>& upperBoundsVector) const { + if (!upperBoundsVector) { + if (this->hasUpperBound(BoundType::Local)) { + upperBoundsVector = std::make_unique>(this->getUpperBounds()); + } else { + upperBoundsVector = std::make_unique>(getMatrixRowCount(), this->getUpperBound()); + } + } else { + if (this->hasUpperBound(BoundType::Local)) { + for (auto& e : *upperBoundsVector) { + e = this->getUpperBound(); + } + } else { + auto upperBoundsIt = this->getUpperBounds().begin(); + for (auto& e : *upperBoundsVector) { + e = *upperBoundsIt; + ++upperBoundsIt; + } + } + } } template diff --git a/src/storm/solver/LinearEquationSolver.h b/src/storm/solver/LinearEquationSolver.h index f17c9d865..dc80bd688 100644 --- a/src/storm/solver/LinearEquationSolver.h +++ b/src/storm/solver/LinearEquationSolver.h @@ -150,6 +150,22 @@ namespace storm { */ virtual void clearCache() const; + enum class BoundType { + Global, + Local, + Any + }; + + /*! + * Retrieves whether this solver has a lower bound. + */ + bool hasLowerBound(BoundType const& type = BoundType::Any) const; + + /*! + * Retrieves whether this solver has an upper bound. + */ + bool hasUpperBound(BoundType const& type = BoundType::Any) const; + /*! * Sets a lower bound for the solution that can potentially be used by the solver. */ @@ -160,34 +176,56 @@ namespace storm { */ void setUpperBound(ValueType const& value); + /*! + * Sets bounds for the solution that can potentially be used by the solver. + */ + void setBounds(ValueType const& lower, ValueType const& upper); + /*! * Retrieves the lower bound (if there is any). */ ValueType const& getLowerBound() const; /*! - * Retrieves the lower bound (if there is any). + * Retrieves the upper bound (if there is any). */ ValueType const& getUpperBound() const; /*! - * Retrieves whether this solver has a lower bound. + * Retrieves a vector containing the lower bounds (if there are any). */ - bool hasLowerBound() const; + std::vector const& getLowerBounds() const; + + /*! + * Retrieves a vector containing the upper bounds (if there are any). + */ + std::vector const& getUpperBounds() const; + + /*! + * Sets lower bounds for the solution that can potentially be used by the solver. + */ + void setLowerBounds(std::vector const& values); + + /*! + * Sets upper bounds for the solution that can potentially be used by the solver. + */ + void setUpperBounds(std::vector const& values); /*! - * Retrieves whether this solver has an upper bound. + * Sets upper bounds for the solution that can potentially be used by the solver. */ - bool hasUpperBound() const; + void setUpperBounds(std::vector&& values); /*! * Sets bounds for the solution that can potentially be used by the solver. */ - void setBounds(ValueType const& lower, ValueType const& upper); - + void setBounds(std::vector const& lower, std::vector const& upper); + protected: virtual bool internalSolveEquations(std::vector& x, std::vector const& b) const = 0; + void createUpperBoundsVector(std::unique_ptr>& upperBoundsVector) const; + // auxiliary storage. If set, this vector has getMatrixRowCount() entries. mutable std::unique_ptr> cachedRowVector; @@ -196,6 +234,12 @@ namespace storm { // An upper bound if one was set. boost::optional upperBound; + + // Lower bounds if they were set. + boost::optional> lowerBounds; + + // Lower bounds if they were set. + boost::optional> upperBounds; private: /*! diff --git a/src/storm/solver/LinearEquationSolverRequirements.cpp b/src/storm/solver/LinearEquationSolverRequirements.cpp index 3d1b34ee9..0e4ed44fc 100644 --- a/src/storm/solver/LinearEquationSolverRequirements.cpp +++ b/src/storm/solver/LinearEquationSolverRequirements.cpp @@ -3,45 +3,45 @@ namespace storm { namespace solver { - LinearEquationSolverRequirements::LinearEquationSolverRequirements() : globalLowerBound(false), globalUpperBound(false) { + LinearEquationSolverRequirements::LinearEquationSolverRequirements() : lowerBounds(false), upperBounds(false) { // Intentionally left empty. } - LinearEquationSolverRequirements& LinearEquationSolverRequirements::requireGlobalLowerBound() { - globalLowerBound = true; + LinearEquationSolverRequirements& LinearEquationSolverRequirements::requireLowerBounds() { + lowerBounds = true; return *this; } - LinearEquationSolverRequirements& LinearEquationSolverRequirements::requireGlobalUpperBound() { - globalUpperBound = true; + LinearEquationSolverRequirements& LinearEquationSolverRequirements::requireUpperBounds() { + upperBounds = true; return *this; } - bool LinearEquationSolverRequirements::requiresGlobalLowerBound() const { - return globalLowerBound; + bool LinearEquationSolverRequirements::requiresLowerBounds() const { + return lowerBounds; } - bool LinearEquationSolverRequirements::requiresGlobalUpperBound() const { - return globalUpperBound; + bool LinearEquationSolverRequirements::requiresUpperBounds() const { + return upperBounds; } bool LinearEquationSolverRequirements::requires(Element const& element) const { switch (element) { - case Element::GlobalLowerBound: return globalLowerBound; break; - case Element::GlobalUpperBound: return globalUpperBound; break; + case Element::LowerBounds: return lowerBounds; break; + case Element::UpperBounds: return upperBounds; break; } } - void LinearEquationSolverRequirements::clearGlobalLowerBound() { - globalLowerBound = false; + void LinearEquationSolverRequirements::clearLowerBounds() { + lowerBounds = false; } - void LinearEquationSolverRequirements::clearGlobalUpperBound() { - globalUpperBound = false; + void LinearEquationSolverRequirements::clearUpperBounds() { + upperBounds = false; } bool LinearEquationSolverRequirements::empty() const { - return !globalLowerBound && !globalUpperBound; + return !lowerBounds && !upperBounds; } } diff --git a/src/storm/solver/LinearEquationSolverRequirements.h b/src/storm/solver/LinearEquationSolverRequirements.h index d07c6a764..168af8005 100644 --- a/src/storm/solver/LinearEquationSolverRequirements.h +++ b/src/storm/solver/LinearEquationSolverRequirements.h @@ -8,27 +8,27 @@ namespace storm { // The different requirements a solver can have. enum class Element { // Requirements that are related to bounds for the actual solution. - GlobalLowerBound, - GlobalUpperBound + LowerBounds, + UpperBounds }; LinearEquationSolverRequirements(); - LinearEquationSolverRequirements& requireGlobalLowerBound(); - LinearEquationSolverRequirements& requireGlobalUpperBound(); + LinearEquationSolverRequirements& requireLowerBounds(); + LinearEquationSolverRequirements& requireUpperBounds(); - bool requiresGlobalLowerBound() const; - bool requiresGlobalUpperBound() const; + bool requiresLowerBounds() const; + bool requiresUpperBounds() const; bool requires(Element const& element) const; - void clearGlobalLowerBound(); - void clearGlobalUpperBound(); + void clearLowerBounds(); + void clearUpperBounds(); bool empty() const; private: - bool globalLowerBound; - bool globalUpperBound; + bool lowerBounds; + bool upperBounds; }; } diff --git a/src/storm/solver/MinMaxLinearEquationSolverRequirements.cpp b/src/storm/solver/MinMaxLinearEquationSolverRequirements.cpp index a90907ff8..41c02731c 100644 --- a/src/storm/solver/MinMaxLinearEquationSolverRequirements.cpp +++ b/src/storm/solver/MinMaxLinearEquationSolverRequirements.cpp @@ -3,7 +3,7 @@ namespace storm { namespace solver { - MinMaxLinearEquationSolverRequirements::MinMaxLinearEquationSolverRequirements(LinearEquationSolverRequirements const& linearEquationSolverRequirements) : noEndComponents(false), validInitialScheduler(false), globalLowerBound(linearEquationSolverRequirements.requiresGlobalLowerBound()), globalUpperBound(linearEquationSolverRequirements.requiresGlobalUpperBound()) { + MinMaxLinearEquationSolverRequirements::MinMaxLinearEquationSolverRequirements(LinearEquationSolverRequirements const& linearEquationSolverRequirements) : noEndComponents(false), validInitialScheduler(false), lowerBounds(linearEquationSolverRequirements.requiresLowerBounds()), upperBounds(linearEquationSolverRequirements.requiresUpperBounds()) { // Intentionally left empty. } @@ -17,13 +17,13 @@ namespace storm { return *this; } - MinMaxLinearEquationSolverRequirements& MinMaxLinearEquationSolverRequirements::requireGlobalLowerBound() { - globalLowerBound = true; + MinMaxLinearEquationSolverRequirements& MinMaxLinearEquationSolverRequirements::requireLowerBounds() { + lowerBounds = true; return *this; } - MinMaxLinearEquationSolverRequirements& MinMaxLinearEquationSolverRequirements::requireGlobalUpperBound() { - globalUpperBound = true; + MinMaxLinearEquationSolverRequirements& MinMaxLinearEquationSolverRequirements::requireUpperBounds() { + upperBounds = true; return *this; } @@ -35,20 +35,20 @@ namespace storm { return validInitialScheduler; } - bool MinMaxLinearEquationSolverRequirements::requiresGlobalLowerBound() const { - return globalLowerBound; + bool MinMaxLinearEquationSolverRequirements::requiresLowerBounds() const { + return lowerBounds; } - bool MinMaxLinearEquationSolverRequirements::requiresGlobalUpperBound() const { - return globalUpperBound; + bool MinMaxLinearEquationSolverRequirements::requiresUpperBounds() const { + return upperBounds; } bool MinMaxLinearEquationSolverRequirements::requires(Element const& element) const { switch (element) { case Element::NoEndComponents: return noEndComponents; break; case Element::ValidInitialScheduler: return validInitialScheduler; break; - case Element::GlobalLowerBound: return globalLowerBound; break; - case Element::GlobalUpperBound: return globalUpperBound; break; + case Element::LowerBounds: return lowerBounds; break; + case Element::UpperBounds: return upperBounds; break; } } @@ -61,16 +61,21 @@ namespace storm { validInitialScheduler = false; } - void MinMaxLinearEquationSolverRequirements::clearGlobalLowerBound() { - globalLowerBound = false; + void MinMaxLinearEquationSolverRequirements::clearLowerBounds() { + lowerBounds = false; } - void MinMaxLinearEquationSolverRequirements::clearGlobalUpperBound() { - globalUpperBound = false; + void MinMaxLinearEquationSolverRequirements::clearUpperBounds() { + upperBounds = false; + } + + void MinMaxLinearEquationSolverRequirements::clearBounds() { + clearLowerBounds(); + clearUpperBounds(); } bool MinMaxLinearEquationSolverRequirements::empty() const { - return !noEndComponents && !validInitialScheduler && !globalLowerBound && !globalUpperBound; + return !noEndComponents && !validInitialScheduler && !lowerBounds && !upperBounds; } } diff --git a/src/storm/solver/MinMaxLinearEquationSolverRequirements.h b/src/storm/solver/MinMaxLinearEquationSolverRequirements.h index 4dfb9f2aa..5b11c2041 100644 --- a/src/storm/solver/MinMaxLinearEquationSolverRequirements.h +++ b/src/storm/solver/MinMaxLinearEquationSolverRequirements.h @@ -16,35 +16,36 @@ namespace storm { ValidInitialScheduler, // Requirements that are related to bounds for the actual solution. - GlobalLowerBound, - GlobalUpperBound + LowerBounds, + UpperBounds }; MinMaxLinearEquationSolverRequirements(LinearEquationSolverRequirements const& linearEquationSolverRequirements = LinearEquationSolverRequirements()); MinMaxLinearEquationSolverRequirements& requireNoEndComponents(); MinMaxLinearEquationSolverRequirements& requireValidInitialScheduler(); - MinMaxLinearEquationSolverRequirements& requireGlobalLowerBound(); - MinMaxLinearEquationSolverRequirements& requireGlobalUpperBound(); + MinMaxLinearEquationSolverRequirements& requireLowerBounds(); + MinMaxLinearEquationSolverRequirements& requireUpperBounds(); bool requiresNoEndComponents() const; bool requiresValidIntialScheduler() const; - bool requiresGlobalLowerBound() const; - bool requiresGlobalUpperBound() const; + bool requiresLowerBounds() const; + bool requiresUpperBounds() const; bool requires(Element const& element) const; void clearNoEndComponents(); void clearValidInitialScheduler(); - void clearGlobalLowerBound(); - void clearGlobalUpperBound(); + void clearLowerBounds(); + void clearUpperBounds(); + void clearBounds(); bool empty() const; private: bool noEndComponents; bool validInitialScheduler; - bool globalLowerBound; - bool globalUpperBound; + bool lowerBounds; + bool upperBounds; }; } diff --git a/src/storm/solver/NativeLinearEquationSolver.cpp b/src/storm/solver/NativeLinearEquationSolver.cpp index 6402f6046..c4cdab06c 100644 --- a/src/storm/solver/NativeLinearEquationSolver.cpp +++ b/src/storm/solver/NativeLinearEquationSolver.cpp @@ -434,13 +434,7 @@ namespace storm { STORM_LOG_INFO("Solving linear equation system (" << x.size() << " rows) with NativeLinearEquationSolver (SoundPower)"); std::vector* lowerX = &x; - if (!this->cachedRowVector) { - this->cachedRowVector = std::make_unique>(getMatrixRowCount(), this->getUpperBound()); - } else { - for (auto& e : *this->cachedRowVector) { - e = this->getUpperBound(); - } - } + this->createUpperBoundsVector(this->cachedRowVector); std::vector* upperX = this->cachedRowVector.get(); bool useGaussSeidelMultiplication = this->getSettings().getPowerMethodMultiplicationStyle() == storm::solver::MultiplicationStyle::GaussSeidel; diff --git a/src/storm/utility/vector.h b/src/storm/utility/vector.h index 4aa684557..069c11759 100644 --- a/src/storm/utility/vector.h +++ b/src/storm/utility/vector.h @@ -883,6 +883,36 @@ namespace storm { } } + /*! + * Takes the input vector and ensures that all entries conform to the bounds. + */ + template + void clip(std::vector& x, ValueType const& bound, bool boundFromBelow) { + for (auto& entry : x) { + if (boundFromBelow && entry < bound) { + entry = bound; + } else if (!boundFromBelow && entry > bound) { + entry = bound; + } + } + } + + /*! + * Takes the input vector and ensures that all entries conform to the bounds. + */ + template + void clip(std::vector& x, std::vector const& bounds, bool boundFromBelow) { + auto boundsIt = bounds.begin(); + for (auto& entry : x) { + if (boundFromBelow && entry < *boundsIt) { + entry = *boundsIt; + } else if (!boundFromBelow && entry > *boundsIt) { + entry = *boundsIt; + } + ++boundsIt; + } + } + /*! * Takes the given offset vector and applies the given contraint. That is, it produces another offset vector that contains * the relative offsets of the entries given by the constraint. diff --git a/src/test/storm/solver/EliminationLinearEquationSolverTest.cpp b/src/test/storm/solver/EliminationLinearEquationSolverTest.cpp index e09d37325..0094b1bcb 100644 --- a/src/test/storm/solver/EliminationLinearEquationSolverTest.cpp +++ b/src/test/storm/solver/EliminationLinearEquationSolverTest.cpp @@ -21,6 +21,7 @@ TEST(EliminationLinearEquationSolver, Solve) { storm::storage::SparseMatrix A; ASSERT_NO_THROW(A = builder.build()); + ASSERT_NO_THROW(A.convertToEquationSystem()); std::vector x(3); std::vector b = {16, -4, -7};