|
|
@ -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<typename ValueType> |
|
|
|
class DsMpi { |
|
|
|
public: |
|
|
|
DsMpi(storm::storage::SparseMatrix<ValueType> const& transitionMatrix, std::vector<ValueType> const& rewards, std::vector<ValueType> const& oneStepTargetProbabilities) : transitionMatrix(transitionMatrix), originalRewards(rewards), backwardTransitions(transitionMatrix.transpose()), p(transitionMatrix.getRowCount()), w(transitionMatrix.getRowCount()), rewards(rewards), targetProbabilities(oneStepTargetProbabilities) { |
|
|
|
// Intentionally left empty.
|
|
|
|
} |
|
|
|
|
|
|
|
std::vector<ValueType> computeUpperBounds() { |
|
|
|
sweep(); |
|
|
|
ValueType lambda = computeLambda(); |
|
|
|
|
|
|
|
// Finally compute the upper bounds for the states.
|
|
|
|
std::vector<ValueType> result(transitionMatrix.getRowCount()); |
|
|
|
auto one = storm::utility::one<ValueType>(); |
|
|
|
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<ValueType>(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<ValueType>(); |
|
|
|
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<storm::storage::sparse::state_type, decltype(cmp)> 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<ValueType> const& transitionMatrix; |
|
|
|
std::vector<ValueType> const& originalRewards; |
|
|
|
|
|
|
|
// Derived from input data.
|
|
|
|
storm::storage::SparseMatrix<ValueType> backwardTransitions; |
|
|
|
|
|
|
|
// Data that the algorithm uses internally.
|
|
|
|
std::vector<ValueType> p; |
|
|
|
std::vector<ValueType> w; |
|
|
|
std::vector<ValueType> rewards; |
|
|
|
std::vector<ValueType> targetProbabilities; |
|
|
|
}; |
|
|
|
|
|
|
|
// This function computes an upper bound on the reachability rewards (see Baier et al, CAV'17).
|
|
|
|
template<typename ValueType> |
|
|
|
std::vector<ValueType> computeUpperRewardBounds(storm::storage::SparseMatrix<ValueType> const& transitionMatrix, std::vector<ValueType> const& rewards, std::vector<ValueType> const& oneStepTargetProbabilities) { |
|
|
|
DsMpi<ValueType> dsmpi(transitionMatrix, rewards, oneStepTargetProbabilities); |
|
|
|
return dsmpi.computeUpperBounds(); |
|
|
|
} |
|
|
|
|
|
|
|
template<typename ValueType, typename RewardModelType> |
|
|
|
std::vector<ValueType> SparseDtmcPrctlHelper<ValueType, RewardModelType>::computeReachabilityRewards(storm::storage::SparseMatrix<ValueType> const& transitionMatrix, storm::storage::SparseMatrix<ValueType> const& backwardTransitions, std::function<std::vector<ValueType>(uint_fast64_t, storm::storage::SparseMatrix<ValueType> const&, storm::storage::BitVector const&)> const& totalStateRewardVectorGetter, storm::storage::BitVector const& targetStates, bool qualitative, storm::solver::LinearEquationSolverFactory<ValueType> const& linearEquationSolverFactory, ModelCheckerHint const& hint) { |
|
|
|
|
|
|
@ -264,10 +391,24 @@ namespace storm { |
|
|
|
|
|
|
|
// Prepare the right-hand side of the equation system.
|
|
|
|
std::vector<ValueType> b = totalStateRewardVectorGetter(submatrix.getRowCount(), transitionMatrix, maybeStates); |
|
|
|
|
|
|
|
storm::solver::LinearEquationSolverRequirements requirements = linearEquationSolverFactory.getRequirements(); |
|
|
|
boost::optional<std::vector<ValueType>> 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<storm::solver::LinearEquationSolver<ValueType>> solver = linearEquationSolverFactory.create(std::move(submatrix)); |
|
|
|
solver->setLowerBound(storm::utility::zero<ValueType>()); |
|
|
|
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.
|
|
|
|