From b56766e9936ef0c8fa01286df5e39d9fab974bbc Mon Sep 17 00:00:00 2001 From: dehnert Date: Fri, 14 Aug 2015 16:53:56 +0200 Subject: [PATCH] more work on reward model that turned out to be refactoring in disguise Former-commit-id: 31a7fa480158a47a60f7458559ef988cc1a6f66b --- CMakeLists.txt | 8 +- .../csl/HybridCtmcCslModelChecker.cpp | 18 +- .../csl/SparseCtmcCslModelChecker.cpp | 39 +- .../csl/SparseCtmcCslModelChecker.h | 15 +- .../csl/helper/SparseCtmcCslHelper.cpp | 0 .../csl/helper/SparseCtmcCslHelper.h | 25 ++ .../csl/helper/SymbolicCtmcCslHelper.cpp | 0 .../csl/helper/SymbolicCtmcCslHelper.h | 0 .../prctl/HybridDtmcPrctlModelChecker.h | 6 +- .../prctl/HybridMdpPrctlModelChecker.cpp | 254 +---------- .../prctl/HybridMdpPrctlModelChecker.h | 14 +- .../prctl/SparseDtmcPrctlModelChecker.cpp | 226 +--------- .../prctl/SparseDtmcPrctlModelChecker.h | 19 +- .../prctl/SparseMdpPrctlModelChecker.cpp | 399 ----------------- .../prctl/SparseMdpPrctlModelChecker.h | 22 +- .../prctl/SymbolicDtmcPrctlModelChecker.cpp | 211 +-------- .../prctl/SymbolicDtmcPrctlModelChecker.h | 19 +- .../prctl/helper/HybridDtmcPrctlHelper.cpp | 0 .../prctl/helper/HybridDtmcPrctlHelper.h | 0 .../prctl/helper/HybridMdpPrctlHelper.cpp | 248 +++++++++++ .../prctl/helper/HybridMdpPrctlHelper.h | 34 ++ .../prctl/helper/SparseDtmcPrctlHelper.cpp | 196 +++++++++ .../prctl/helper/SparseDtmcPrctlHelper.h | 39 ++ .../prctl/helper/SparseMdpPrctlHelper.cpp | 407 ++++++++++++++++++ .../prctl/helper/SparseMdpPrctlHelper.h | 45 ++ .../prctl/helper/SymbolicDtmcPrctlHelper.cpp | 195 +++++++++ .../prctl/helper/SymbolicDtmcPrctlHelper.h | 35 ++ .../prctl/helper/SymbolicMdpPrctlHelper.cpp | 0 .../prctl/helper/SymbolicMdpPrctlHelper.h | 0 .../SparsePropositionalModelChecker.h | 3 + src/models/sparse/Ctmc.h | 3 - src/settings/ArgumentBuilder.h | 2 +- .../MaximalEndComponentDecomposition.cpp | 29 +- .../MaximalEndComponentDecomposition.h | 13 +- ...tronglyConnectedComponentDecomposition.cpp | 12 + .../StronglyConnectedComponentDecomposition.h | 24 ++ 36 files changed, 1387 insertions(+), 1173 deletions(-) create mode 100644 src/modelchecker/csl/helper/SparseCtmcCslHelper.cpp create mode 100644 src/modelchecker/csl/helper/SparseCtmcCslHelper.h create mode 100644 src/modelchecker/csl/helper/SymbolicCtmcCslHelper.cpp create mode 100644 src/modelchecker/csl/helper/SymbolicCtmcCslHelper.h create mode 100644 src/modelchecker/prctl/helper/HybridDtmcPrctlHelper.cpp create mode 100644 src/modelchecker/prctl/helper/HybridDtmcPrctlHelper.h create mode 100644 src/modelchecker/prctl/helper/HybridMdpPrctlHelper.cpp create mode 100644 src/modelchecker/prctl/helper/HybridMdpPrctlHelper.h create mode 100644 src/modelchecker/prctl/helper/SparseDtmcPrctlHelper.cpp create mode 100644 src/modelchecker/prctl/helper/SparseDtmcPrctlHelper.h create mode 100644 src/modelchecker/prctl/helper/SparseMdpPrctlHelper.cpp create mode 100644 src/modelchecker/prctl/helper/SparseMdpPrctlHelper.h create mode 100644 src/modelchecker/prctl/helper/SymbolicDtmcPrctlHelper.cpp create mode 100644 src/modelchecker/prctl/helper/SymbolicDtmcPrctlHelper.h create mode 100644 src/modelchecker/prctl/helper/SymbolicMdpPrctlHelper.cpp create mode 100644 src/modelchecker/prctl/helper/SymbolicMdpPrctlHelper.h diff --git a/CMakeLists.txt b/CMakeLists.txt index 6bdaa29a2..96b174c64 100644 --- a/CMakeLists.txt +++ b/CMakeLists.txt @@ -377,8 +377,10 @@ file(GLOB_RECURSE STORM_BUILDER_FILES ${PROJECT_SOURCE_DIR}/src/builder/*.h ${PR file(GLOB_RECURSE STORM_EXCEPTIONS_FILES ${PROJECT_SOURCE_DIR}/src/exceptions/*.h ${PROJECT_SOURCE_DIR}/src/exceptions/*.cpp) file(GLOB_RECURSE STORM_LOGIC_FILES ${PROJECT_SOURCE_DIR}/src/logic/*.h ${PROJECT_SOURCE_DIR}/src/logic/*.cpp) file(GLOB STORM_MODELCHECKER_FILES ${PROJECT_SOURCE_DIR}/src/modelchecker/*.h ${PROJECT_SOURCE_DIR}/src/modelchecker/*.cpp) -file(GLOB_RECURSE STORM_MODELCHECKER_PRCTL_FILES ${PROJECT_SOURCE_DIR}/src/modelchecker/prctl/*.h ${PROJECT_SOURCE_DIR}/src/modelchecker/prctl/*.cpp) -file(GLOB_RECURSE STORM_MODELCHECKER_CSL_FILES ${PROJECT_SOURCE_DIR}/src/modelchecker/csl/*.h ${PROJECT_SOURCE_DIR}/src/modelchecker/csl/*.cpp) +file(GLOB STORM_MODELCHECKER_PRCTL_FILES ${PROJECT_SOURCE_DIR}/src/modelchecker/prctl/*.h ${PROJECT_SOURCE_DIR}/src/modelchecker/prctl/*.cpp) +file(GLOB_RECURSE STORM_MODELCHECKER_PRCTL_HELPER_FILES ${PROJECT_SOURCE_DIR}/src/modelchecker/prctl/helper/*.h ${PROJECT_SOURCE_DIR}/src/modelchecker/prctl/helper/*.cpp) +file(GLOB STORM_MODELCHECKER_CSL_FILES ${PROJECT_SOURCE_DIR}/src/modelchecker/csl/*.h ${PROJECT_SOURCE_DIR}/src/modelchecker/csl/*.cpp) +file(GLOB_RECURSE STORM_MODELCHECKER_CSL_HELPER_FILES ${PROJECT_SOURCE_DIR}/src/modelchecker/csl/helper/*.h ${PROJECT_SOURCE_DIR}/src/modelchecker/csl/helper/*.cpp) file(GLOB_RECURSE STORM_MODELCHECKER_REACHABILITY_FILES ${PROJECT_SOURCE_DIR}/src/modelchecker/reachability/*.h ${PROJECT_SOURCE_DIR}/src/modelchecker/reachability/*.cpp) file(GLOB_RECURSE STORM_MODELCHECKER_PROPOSITIONAL_FILES ${PROJECT_SOURCE_DIR}/src/modelchecker/propositional/*.h ${PROJECT_SOURCE_DIR}/src/modelchecker/propositional/*.cpp) file(GLOB_RECURSE STORM_MODELCHECKER_RESULTS_FILES ${PROJECT_SOURCE_DIR}/src/modelchecker/results/*.h ${PROJECT_SOURCE_DIR}/src/modelchecker/results/*.cpp) @@ -416,7 +418,9 @@ source_group(logic FILES ${STORM_LOGIC_FILES}) source_group(generated FILES ${STORM_BUILD_HEADERS} ${STORM_BUILD_SOURCES}) source_group(modelchecker FILES ${STORM_MODELCHECKER_FILES}) source_group(modelchecker\\prctl FILES ${STORM_MODELCHECKER_PRCTL_FILES}) +source_group(modelchecker\\prctl\\helper FILES ${STORM_MODELCHECKER_PRCTL_HELPER_FILES}) source_group(modelchecker\\csl FILES ${STORM_MODELCHECKER_CSL_FILES}) +source_group(modelchecker\\csl\\helper FILES ${STORM_MODELCHECKER_CSL_HELPER_FILES}) source_group(modelchecker\\reachability FILES ${STORM_MODELCHECKER_REACHABILITY_FILES}) source_group(modelchecker\\propositional FILES ${STORM_MODELCHECKER_PROPOSITIONAL_FILES}) source_group(modelchecker\\results FILES ${STORM_MODELCHECKER_RESULTS_FILES}) diff --git a/src/modelchecker/csl/HybridCtmcCslModelChecker.cpp b/src/modelchecker/csl/HybridCtmcCslModelChecker.cpp index 498a2fc02..428584385 100644 --- a/src/modelchecker/csl/HybridCtmcCslModelChecker.cpp +++ b/src/modelchecker/csl/HybridCtmcCslModelChecker.cpp @@ -1,5 +1,5 @@ #include "src/modelchecker/csl/HybridCtmcCslModelChecker.h" -#include "src/modelchecker/csl/SparseCtmcCslModelChecker.h" +#include "src/modelchecker/csl/helper/SparseCtmcCslHelper.h" #include "src/modelchecker/prctl/HybridDtmcPrctlModelChecker.h" #include "src/storage/dd/CuddOdd.h" @@ -155,7 +155,7 @@ namespace storm { // Finally compute the transient probabilities. std::vector values(statesWithProbabilityGreater0NonPsi.getNonZeroCount(), storm::utility::zero()); - std::vector subresult = SparseCtmcCslModelChecker::computeTransientProbabilities(explicitUniformizedMatrix, &explicitB, upperBound, uniformizationRate, values, *this->linearEquationSolverFactory); + std::vector subresult = SparseCtmcCslHelper::computeTransientProbabilities(explicitUniformizedMatrix, &explicitB, upperBound, uniformizationRate, values, *this->linearEquationSolverFactory); return std::unique_ptr(new HybridQuantitativeCheckResult(this->getModel().getReachableStates(), (psiStates || !statesWithProbabilityGreater0) && this->getModel().getReachableStates(), @@ -195,7 +195,7 @@ namespace storm { storm::storage::SparseMatrix explicitUniformizedMatrix = uniformizedMatrix.toMatrix(odd, odd); // Compute the transient probabilities. - result = SparseCtmcCslModelChecker::computeTransientProbabilities(explicitUniformizedMatrix, nullptr, lowerBound, uniformizationRate, result, *this->linearEquationSolverFactory); + result = SparseCtmcCslHelper::computeTransientProbabilities(explicitUniformizedMatrix, nullptr, lowerBound, uniformizationRate, result, *this->linearEquationSolverFactory); return std::unique_ptr(new HybridQuantitativeCheckResult(this->getModel().getReachableStates(), !relevantStates && this->getModel().getReachableStates(), this->getModel().getManager().getAddZero(), relevantStates, odd, result)); } else { @@ -221,7 +221,7 @@ namespace storm { // Compute the transient probabilities. std::vector values(statesWithProbabilityGreater0NonPsi.getNonZeroCount(), storm::utility::zero()); - std::vector subResult = SparseCtmcCslModelChecker::computeTransientProbabilities(explicitUniformizedMatrix, &explicitB, upperBound - lowerBound, uniformizationRate, values, *this->linearEquationSolverFactory); + std::vector subResult = SparseCtmcCslHelper::computeTransientProbabilities(explicitUniformizedMatrix, &explicitB, upperBound - lowerBound, uniformizationRate, values, *this->linearEquationSolverFactory); // Transform the explicit result to a hybrid check result, so we can easily convert it to // a symbolic qualitative format. @@ -256,7 +256,7 @@ namespace storm { uniformizedMatrix = this->computeUniformizedMatrix(this->getModel(), this->getModel().getTransitionMatrix(), exitRates, relevantStates, uniformizationRate); explicitUniformizedMatrix = uniformizedMatrix.toMatrix(odd, odd); - newSubresult = SparseCtmcCslModelChecker::computeTransientProbabilities(explicitUniformizedMatrix, nullptr, lowerBound, uniformizationRate, newSubresult, *this->linearEquationSolverFactory); + newSubresult = SparseCtmcCslHelper::computeTransientProbabilities(explicitUniformizedMatrix, nullptr, lowerBound, uniformizationRate, newSubresult, *this->linearEquationSolverFactory); return std::unique_ptr(new HybridQuantitativeCheckResult(this->getModel().getReachableStates(), !relevantStates && this->getModel().getReachableStates(), this->getModel().getManager().getAddZero(), relevantStates, odd, newSubresult)); } else { @@ -276,7 +276,7 @@ namespace storm { storm::dd::Add uniformizedMatrix = this->computeUniformizedMatrix(this->getModel(), this->getModel().getTransitionMatrix(), exitRates, statesWithProbabilityGreater0, uniformizationRate); storm::storage::SparseMatrix explicitUniformizedMatrix = uniformizedMatrix.toMatrix(odd, odd); - newSubresult = SparseCtmcCslModelChecker::computeTransientProbabilities(explicitUniformizedMatrix, nullptr, lowerBound, uniformizationRate, newSubresult, *this->linearEquationSolverFactory); + newSubresult = SparseCtmcCslHelper::computeTransientProbabilities(explicitUniformizedMatrix, nullptr, lowerBound, uniformizationRate, newSubresult, *this->linearEquationSolverFactory); return std::unique_ptr(new HybridQuantitativeCheckResult(this->getModel().getReachableStates(), !statesWithProbabilityGreater0 && this->getModel().getReachableStates(), this->getModel().getManager().getAddZero(), statesWithProbabilityGreater0, odd, newSubresult)); } @@ -311,7 +311,7 @@ namespace storm { storm::dd::Add uniformizedMatrix = this->computeUniformizedMatrix(this->getModel(), this->getModel().getTransitionMatrix(), this->getModel().getExitRateVector(), this->getModel().getReachableStates(), uniformizationRate); storm::storage::SparseMatrix explicitUniformizedMatrix = uniformizedMatrix.toMatrix(odd, odd); - result = SparseCtmcCslModelChecker::computeTransientProbabilities(explicitUniformizedMatrix, nullptr, timeBound, uniformizationRate, result, *this->linearEquationSolverFactory); + result = SparseCtmcCslHelper::computeTransientProbabilities(explicitUniformizedMatrix, nullptr, timeBound, uniformizationRate, result, *this->linearEquationSolverFactory); } return std::unique_ptr(new HybridQuantitativeCheckResult(this->getModel().getReachableStates(), this->getModel().getManager().getBddZero(), this->getModel().getManager().getAddZero(), this->getModel().getReachableStates(), odd, result)); @@ -353,7 +353,7 @@ namespace storm { std::vector explicitTotalRewardVector = totalRewardVector.template toVector(odd); // Finally, compute the transient probabilities. - std::vector result = SparseCtmcCslModelChecker::template computeTransientProbabilities(explicitUniformizedMatrix, nullptr, timeBound, uniformizationRate, explicitTotalRewardVector, *this->linearEquationSolverFactory); + std::vector result = SparseCtmcCslHelper::template computeTransientProbabilities(explicitUniformizedMatrix, nullptr, timeBound, uniformizationRate, explicitTotalRewardVector, *this->linearEquationSolverFactory); return std::unique_ptr(new HybridQuantitativeCheckResult(this->getModel().getReachableStates(), this->getModel().getManager().getBddZero(), this->getModel().getManager().getAddZero(), this->getModel().getReachableStates(), std::move(odd), std::move(result))); } @@ -370,7 +370,7 @@ namespace storm { storm::storage::SparseMatrix explicitProbabilityMatrix = probabilityMatrix.toMatrix(odd, odd); std::vector explicitExitRateVector = this->getModel().getExitRateVector().template toVector(odd); - std::vector result = SparseCtmcCslModelChecker::computeLongRunAverageHelper(explicitProbabilityMatrix, subResult.getTruthValuesVector().toVector(odd), &explicitExitRateVector, qualitative, *this->linearEquationSolverFactory); + std::vector result = SparseCtmcCslHelper::computeLongRunAverage(explicitProbabilityMatrix, subResult.getTruthValuesVector().toVector(odd), &explicitExitRateVector, qualitative, *this->linearEquationSolverFactory); return std::unique_ptr(new HybridQuantitativeCheckResult(this->getModel().getReachableStates(), this->getModel().getManager().getBddZero(), this->getModel().getManager().getAddZero(), this->getModel().getReachableStates(), std::move(odd), std::move(result))); } diff --git a/src/modelchecker/csl/SparseCtmcCslModelChecker.cpp b/src/modelchecker/csl/SparseCtmcCslModelChecker.cpp index ccffd7b3d..1dba53698 100644 --- a/src/modelchecker/csl/SparseCtmcCslModelChecker.cpp +++ b/src/modelchecker/csl/SparseCtmcCslModelChecker.cpp @@ -21,12 +21,12 @@ namespace storm { namespace modelchecker { template - SparseCtmcCslModelChecker::SparseCtmcCslModelChecker(SparseCtmcModelType const& model) : SparsePropositionalModelChecker(model), linearEquationSolverFactory(new storm::utility::solver::LinearEquationSolverFactory()) { + SparseCtmcCslModelChecker::SparseCtmcCslModelChecker(SparseCtmcModelType const& model) : SparsePropositionalModelChecker(model), linearEquationSolverFactory(new storm::utility::solver::LinearEquationSolverFactory()) { // Intentionally left empty. } template - SparseCtmcCslModelChecker::SparseCtmcCslModelChecker(SparseCtmcModelType const& model, std::unique_ptr>&& linearEquationSolverFactory) : SparsePropositionalModelChecker(model), linearEquationSolverFactory(std::move(linearEquationSolverFactory)) { + SparseCtmcCslModelChecker::SparseCtmcCslModelChecker(SparseCtmcModelType const& model, std::unique_ptr>&& linearEquationSolverFactory) : SparsePropositionalModelChecker(model), linearEquationSolverFactory(std::move(linearEquationSolverFactory)) { // Intentionally left empty. } @@ -58,7 +58,7 @@ namespace storm { std::unique_ptr SparseCtmcCslModelChecker::computeNextProbabilities(storm::logic::NextFormula const& pathFormula, bool qualitative, boost::optional const& optimalityType) { std::unique_ptr subResultPointer = this->check(pathFormula.getSubformula()); ExplicitQualitativeCheckResult const& subResult = subResultPointer->asExplicitQualitativeCheckResult(); - std::vector result = SparseDtmcPrctlModelChecker::computeNextProbabilitiesHelper(this->computeProbabilityMatrix(this->getModel().getTransitionMatrix(), this->getModel().getExitRateVector()), subResult.getTruthValuesVector(), *this->linearEquationSolverFactory); + std::vector result = SparseDtmcPrctlModelChecker::computeNextProbabilitiesHelper(this->computeProbabilityMatrix(this->getModel().getTransitionMatrix(), this->getModel().getExitRateVector()), subResult.getTruthValuesVector(), *this->linearEquationSolverFactory); return std::unique_ptr(new ExplicitQuantitativeCheckResult(std::move(result))); } @@ -72,12 +72,7 @@ namespace storm { } template - storm::models::sparse::Ctmc const& SparseCtmcCslModelChecker::getModel() const { - return this->template getModelAs>(); - } - - template - std::vector SparseCtmcCslModelChecker::computeBoundedUntilProbabilitiesHelper(storm::storage::BitVector const& phiStates, storm::storage::BitVector const& psiStates, std::vector const& exitRates, bool qualitative, double lowerBound, double upperBound) const { + std::vector SparseCtmcCslModelChecker::computeBoundedUntilProbabilitiesHelper(storm::storage::BitVector const& phiStates, storm::storage::BitVector const& psiStates, std::vector const& exitRates, bool qualitative, double lowerBound, double upperBound) const { // If the time bounds are [0, inf], we rather call untimed reachability. storm::utility::ConstantsComparator comparator; if (comparator.isZero(lowerBound) && comparator.isInfinity(upperBound)) { @@ -241,7 +236,7 @@ namespace storm { } template - storm::storage::SparseMatrix SparseCtmcCslModelChecker::computeUniformizedMatrix(storm::storage::SparseMatrix const& transitionMatrix, storm::storage::BitVector const& maybeStates, ValueType uniformizationRate, std::vector const& exitRates) { + storm::storage::SparseMatrix SparseCtmcCslModelChecker::computeUniformizedMatrix(storm::storage::SparseMatrix const& transitionMatrix, storm::storage::BitVector const& maybeStates, ValueType uniformizationRate, std::vector const& exitRates) { STORM_LOG_DEBUG("Computing uniformized matrix using uniformization rate " << uniformizationRate << "."); STORM_LOG_DEBUG("Keeping " << maybeStates.getNumberOfSetBits() << " rows."); @@ -269,7 +264,7 @@ namespace storm { template template - std::vector SparseCtmcCslModelChecker::computeTransientProbabilities(storm::storage::SparseMatrix const& uniformizedMatrix, std::vector const* addVector, ValueType timeBound, ValueType uniformizationRate, std::vector values, storm::utility::solver::LinearEquationSolverFactory const& linearEquationSolverFactory) { + std::vector SparseCtmcCslModelChecker::computeTransientProbabilities(storm::storage::SparseMatrix const& uniformizedMatrix, std::vector const* addVector, ValueType timeBound, ValueType uniformizationRate, std::vector values, storm::utility::solver::LinearEquationSolverFactory const& linearEquationSolverFactory) { ValueType lambda = timeBound * uniformizationRate; @@ -351,17 +346,17 @@ namespace storm { } template - std::vector SparseCtmcCslModelChecker::computeUntilProbabilitiesHelper(storm::storage::SparseMatrix const& transitionMatrix, storm::storage::SparseMatrix const& backwardTransitions, storm::storage::BitVector const& phiStates, storm::storage::BitVector const& psiStates, bool qualitative, storm::utility::solver::LinearEquationSolverFactory const& linearEquationSolverFactory) { - return SparseDtmcPrctlModelChecker::computeUntilProbabilitiesHelper(transitionMatrix, backwardTransitions, phiStates, psiStates, qualitative, linearEquationSolverFactory); + std::vector SparseCtmcCslModelChecker::computeUntilProbabilitiesHelper(storm::storage::SparseMatrix const& transitionMatrix, storm::storage::SparseMatrix const& backwardTransitions, storm::storage::BitVector const& phiStates, storm::storage::BitVector const& psiStates, bool qualitative, storm::utility::solver::LinearEquationSolverFactory const& linearEquationSolverFactory) { + return SparseDtmcPrctlModelChecker::computeUntilProbabilitiesHelper(transitionMatrix, backwardTransitions, phiStates, psiStates, qualitative, linearEquationSolverFactory); } template - std::unique_ptr SparseCtmcCslModelChecker::computeInstantaneousRewards(storm::logic::InstantaneousRewardFormula const& rewardPathFormula, bool qualitative, boost::optional const& optimalityType) { + std::unique_ptr SparseCtmcCslModelChecker::computeInstantaneousRewards(storm::logic::InstantaneousRewardFormula const& rewardPathFormula, boost::optional const& rewardModelName, bool qualitative, boost::optional const& optimalityType) { return std::unique_ptr(new ExplicitQuantitativeCheckResult(this->computeInstantaneousRewardsHelper(rewardPathFormula.getContinuousTimeBound()))); } template - std::vector SparseCtmcCslModelChecker::computeInstantaneousRewardsHelper(double timeBound) const { + std::vector SparseCtmcCslModelChecker::computeInstantaneousRewardsHelper(double timeBound) const { // Only compute the result if the model has a state-based reward this->getModel(). STORM_LOG_THROW(this->getModel().hasStateRewards(), storm::exceptions::InvalidPropertyException, "Missing reward model for formula. Skipping formula."); @@ -385,12 +380,12 @@ namespace storm { } template - std::unique_ptr SparseCtmcCslModelChecker::computeCumulativeRewards(storm::logic::CumulativeRewardFormula const& rewardPathFormula, bool qualitative, boost::optional const& optimalityType) { + std::unique_ptr SparseCtmcCslModelChecker::computeCumulativeRewards(storm::logic::CumulativeRewardFormula const& rewardPathFormula, boost::optional const& rewardModelName, bool qualitative, boost::optional const& optimalityType) { return std::unique_ptr(new ExplicitQuantitativeCheckResult(this->computeCumulativeRewardsHelper(rewardPathFormula.getContinuousTimeBound()))); } template - std::vector SparseCtmcCslModelChecker::computeCumulativeRewardsHelper(double timeBound) const { + std::vector SparseCtmcCslModelChecker::computeCumulativeRewardsHelper(double timeBound) const { // Only compute the result if the model has a state-based reward this->getModel(). STORM_LOG_THROW(this->getModel().hasStateRewards() || this->getModel().hasTransitionRewards(), storm::exceptions::InvalidPropertyException, "Missing reward model for formula. Skipping formula."); @@ -427,7 +422,7 @@ namespace storm { } template - storm::storage::SparseMatrix SparseCtmcCslModelChecker::computeProbabilityMatrix(storm::storage::SparseMatrix const& rateMatrix, std::vector const& exitRates) { + storm::storage::SparseMatrix SparseCtmcCslModelChecker::computeProbabilityMatrix(storm::storage::SparseMatrix const& rateMatrix, std::vector const& exitRates) { // Turn the rates into probabilities by scaling each row with the exit rate of the state. storm::storage::SparseMatrix result(rateMatrix); for (uint_fast64_t row = 0; row < result.getRowCount(); ++row) { @@ -439,7 +434,7 @@ namespace storm { } template - storm::storage::SparseMatrix SparseCtmcCslModelChecker::computeGeneratorMatrix(storm::storage::SparseMatrix const& rateMatrix, std::vector const& exitRates) { + storm::storage::SparseMatrix SparseCtmcCslModelChecker::computeGeneratorMatrix(storm::storage::SparseMatrix const& rateMatrix, std::vector const& exitRates) { storm::storage::SparseMatrix generatorMatrix(rateMatrix, true); // Place the negative exit rate on the diagonal. @@ -455,7 +450,7 @@ namespace storm { } template - std::unique_ptr SparseCtmcCslModelChecker::computeReachabilityRewards(storm::logic::ReachabilityRewardFormula const& rewardPathFormula, bool qualitative, boost::optional const& optimalityType) { + std::unique_ptr SparseCtmcCslModelChecker::computeReachabilityRewards(storm::logic::ReachabilityRewardFormula const& rewardPathFormula, boost::optional const& rewardModelName, bool qualitative, boost::optional const& optimalityType) { std::unique_ptr subResultPointer = this->check(rewardPathFormula.getSubformula()); ExplicitQualitativeCheckResult const& subResult = subResultPointer->asExplicitQualitativeCheckResult(); storm::storage::SparseMatrix probabilityMatrix = computeProbabilityMatrix(this->getModel().getTransitionMatrix(), this->getModel().getExitRateVector()); @@ -482,8 +477,8 @@ namespace storm { return std::unique_ptr(new ExplicitQuantitativeCheckResult(computeLongRunAverageHelper(probabilityMatrix, subResult.getTruthValuesVector(), &this->getModel().getExitRateVector(), qualitative, *linearEquationSolverFactory))); } - template - std::vector SparseCtmcCslModelChecker::computeLongRunAverageHelper(storm::storage::SparseMatrix const& transitionMatrix, storm::storage::BitVector const& psiStates, std::vector const* exitRateVector, bool qualitative, storm::utility::solver::LinearEquationSolverFactory const& linearEquationSolverFactory) { + template + std::vector SparseCtmcCslModelChecker::computeLongRunAverageHelper(storm::storage::SparseMatrix const& transitionMatrix, storm::storage::BitVector const& psiStates, std::vector const* exitRateVector, bool qualitative, storm::utility::solver::LinearEquationSolverFactory const& linearEquationSolverFactory) { // If there are no goal states, we avoid the computation and directly return zero. uint_fast64_t numOfStates = transitionMatrix.getRowCount(); if (psiStates.empty()) { diff --git a/src/modelchecker/csl/SparseCtmcCslModelChecker.h b/src/modelchecker/csl/SparseCtmcCslModelChecker.h index 498eea05b..d5fadf772 100644 --- a/src/modelchecker/csl/SparseCtmcCslModelChecker.h +++ b/src/modelchecker/csl/SparseCtmcCslModelChecker.h @@ -20,9 +20,7 @@ namespace storm { class SparseCtmcCslModelChecker : public SparsePropositionalModelChecker { public: typedef typename SparseCtmcModelType::ValueType ValueType; - - friend class HybridCtmcCslModelChecker; - friend class SparseDtmcPrctlModelChecker; + typedef typename SparseCtmcModelType::RewardModelType RewardModelType; explicit SparseCtmcCslModelChecker(SparseCtmcModelType const& model); explicit SparseCtmcCslModelChecker(SparseCtmcModelType const& model, std::unique_ptr>&& linearEquationSolverFactory); @@ -32,17 +30,16 @@ namespace storm { virtual std::unique_ptr computeBoundedUntilProbabilities(storm::logic::BoundedUntilFormula const& pathFormula, bool qualitative = false, boost::optional const& optimalityType = boost::optional()) override; virtual std::unique_ptr computeNextProbabilities(storm::logic::NextFormula const& pathFormula, bool qualitative = false, boost::optional const& optimalityType = boost::optional()) override; virtual std::unique_ptr computeUntilProbabilities(storm::logic::UntilFormula const& pathFormula, bool qualitative = false, boost::optional const& optimalityType = boost::optional()) override; - virtual std::unique_ptr computeInstantaneousRewards(storm::logic::InstantaneousRewardFormula const& rewardPathFormula, bool qualitative = false, boost::optional const& optimalityType = boost::optional()) override; - virtual std::unique_ptr computeCumulativeRewards(storm::logic::CumulativeRewardFormula const& rewardPathFormula, bool qualitative = false, boost::optional const& optimalityType = boost::optional()) override; - virtual std::unique_ptr computeReachabilityRewards(storm::logic::ReachabilityRewardFormula const& rewardPathFormula, bool qualitative = false, boost::optional const& optimalityType = boost::optional()) override; + virtual std::unique_ptr computeCumulativeRewards(storm::logic::CumulativeRewardFormula const& rewardPathFormula, boost::optional const& rewardModelName = boost::optional(), bool qualitative = false, boost::optional const& optimalityType = boost::optional()) override; + virtual std::unique_ptr computeInstantaneousRewards(storm::logic::InstantaneousRewardFormula const& rewardPathFormula, boost::optional const& rewardModelName = boost::optional(), bool qualitative = false, boost::optional const& optimalityType = boost::optional()) override; + virtual std::unique_ptr computeReachabilityRewards(storm::logic::ReachabilityRewardFormula const& rewardPathFormula, boost::optional const& rewardModelName = boost::optional(), bool qualitative = false, boost::optional const& optimalityType = boost::optional()) override; virtual std::unique_ptr computeLongRunAverage(storm::logic::StateFormula const& stateFormula, bool qualitative = false, boost::optional const& optimalityType = boost::optional()) override; private: - // The methods that perform the actual checking. std::vector computeBoundedUntilProbabilitiesHelper(storm::storage::BitVector const& phiStates, storm::storage::BitVector const& psiStates, std::vector const& exitRates, bool qualitative, double lowerBound, double upperBound) const; static std::vector computeUntilProbabilitiesHelper(storm::storage::SparseMatrix const& transitionMatrix, storm::storage::SparseMatrix const& backwardTransitions, storm::storage::BitVector const& phiStates, storm::storage::BitVector const& psiStates, bool qualitative, storm::utility::solver::LinearEquationSolverFactory const& linearEquationSolverFactory); - std::vector computeInstantaneousRewardsHelper(double timeBound) const; - std::vector computeCumulativeRewardsHelper(double timeBound) const; + std::vector computeInstantaneousRewardsHelper(RewardModelType const& rewardModel, double timeBound) const; + std::vector computeCumulativeRewardsHelper(RewardModelType const& rewardModel, double timeBound) const; static std::vector computeLongRunAverageHelper(storm::storage::SparseMatrix const& transitionMatrix, storm::storage::BitVector const& psiStates, std::vector const* exitRateVector, bool qualitative, storm::utility::solver::LinearEquationSolverFactory const& linearEquationSolverFactory); /*! * Computes the matrix representing the transitions of the uniformized CTMC. diff --git a/src/modelchecker/csl/helper/SparseCtmcCslHelper.cpp b/src/modelchecker/csl/helper/SparseCtmcCslHelper.cpp new file mode 100644 index 000000000..e69de29bb diff --git a/src/modelchecker/csl/helper/SparseCtmcCslHelper.h b/src/modelchecker/csl/helper/SparseCtmcCslHelper.h new file mode 100644 index 000000000..d2edcbd18 --- /dev/null +++ b/src/modelchecker/csl/helper/SparseCtmcCslHelper.h @@ -0,0 +1,25 @@ +#ifndef STORM_MODELCHECKER_SPARSE_CTMC_CSL_MODELCHECKER_HELPER_H_ +#define STORM_MODELCHECKER_SPARSE_CTMC_CSL_MODELCHECKER_HELPER_H_ + +#include + +#include "src/models/sparse/StandardRewardModel.h" + +#include "src/storage/SparseMatrix.h" +#include "src/storage/BitVector.h" + +#include "src/utility/solver.h" + +namespace storm { + namespace modelchecker { + namespace helper { + template > + class SparseCtmcCslHelper { + public: + + } + } + } +} + +#endif /* STORM_MODELCHECKER_SPARSE_CTMC_CSL_MODELCHECKER_HELPER_H_ */ \ No newline at end of file diff --git a/src/modelchecker/csl/helper/SymbolicCtmcCslHelper.cpp b/src/modelchecker/csl/helper/SymbolicCtmcCslHelper.cpp new file mode 100644 index 000000000..e69de29bb diff --git a/src/modelchecker/csl/helper/SymbolicCtmcCslHelper.h b/src/modelchecker/csl/helper/SymbolicCtmcCslHelper.h new file mode 100644 index 000000000..e69de29bb diff --git a/src/modelchecker/prctl/HybridDtmcPrctlModelChecker.h b/src/modelchecker/prctl/HybridDtmcPrctlModelChecker.h index 8964d29d7..016e2898a 100644 --- a/src/modelchecker/prctl/HybridDtmcPrctlModelChecker.h +++ b/src/modelchecker/prctl/HybridDtmcPrctlModelChecker.h @@ -23,9 +23,9 @@ namespace storm { virtual std::unique_ptr computeBoundedUntilProbabilities(storm::logic::BoundedUntilFormula const& pathFormula, bool qualitative = false, boost::optional const& optimalityType = boost::optional()) override; virtual std::unique_ptr computeNextProbabilities(storm::logic::NextFormula const& pathFormula, bool qualitative = false, boost::optional const& optimalityType = boost::optional()) override; virtual std::unique_ptr computeUntilProbabilities(storm::logic::UntilFormula const& pathFormula, bool qualitative = false, boost::optional const& optimalityType = boost::optional()) override; - virtual std::unique_ptr computeCumulativeRewards(storm::logic::CumulativeRewardFormula const& rewardPathFormula, bool qualitative = false, boost::optional const& optimalityType = boost::optional()) override; - virtual std::unique_ptr computeInstantaneousRewards(storm::logic::InstantaneousRewardFormula const& rewardPathFormula, bool qualitative = false, boost::optional const& optimalityType = boost::optional()) override; - virtual std::unique_ptr computeReachabilityRewards(storm::logic::ReachabilityRewardFormula const& rewardPathFormula, bool qualitative = false, boost::optional const& optimalityType = boost::optional()) override; + virtual std::unique_ptr computeCumulativeRewards(storm::logic::CumulativeRewardFormula const& rewardPathFormula, boost::optional const& rewardModelName = boost::optional(), bool qualitative = false, boost::optional const& optimalityType = boost::optional()) override; + virtual std::unique_ptr computeInstantaneousRewards(storm::logic::InstantaneousRewardFormula const& rewardPathFormula, boost::optional const& rewardModelName = boost::optional(), bool qualitative = false, boost::optional const& optimalityType = boost::optional()) override; + virtual std::unique_ptr computeReachabilityRewards(storm::logic::ReachabilityRewardFormula const& rewardPathFormula, boost::optional const& rewardModelName = boost::optional(), bool qualitative = false, boost::optional const& optimalityType = boost::optional()) override; virtual std::unique_ptr computeLongRunAverage(storm::logic::StateFormula const& stateFormula, bool qualitative = false, boost::optional const& optimalityType = boost::optional()) override; protected: diff --git a/src/modelchecker/prctl/HybridMdpPrctlModelChecker.cpp b/src/modelchecker/prctl/HybridMdpPrctlModelChecker.cpp index f54bed203..a8af6665c 100644 --- a/src/modelchecker/prctl/HybridMdpPrctlModelChecker.cpp +++ b/src/modelchecker/prctl/HybridMdpPrctlModelChecker.cpp @@ -28,71 +28,7 @@ namespace storm { bool HybridMdpPrctlModelChecker::canHandle(storm::logic::Formula const& formula) const { return formula.isPctlStateFormula() || formula.isPctlPathFormula() || formula.isRewardPathFormula(); } - - template - std::unique_ptr HybridMdpPrctlModelChecker::computeUntilProbabilitiesHelper(bool minimize, storm::models::symbolic::NondeterministicModel const& model, storm::dd::Add const& transitionMatrix, storm::dd::Bdd const& phiStates, storm::dd::Bdd const& psiStates, bool qualitative, storm::utility::solver::MinMaxLinearEquationSolverFactory const& linearEquationSolverFactory) { - // We need to identify the states which have to be taken out of the matrix, i.e. all states that have - // probability 0 and 1 of satisfying the until-formula. - std::pair, storm::dd::Bdd> statesWithProbability01; - if (minimize) { - statesWithProbability01 = storm::utility::graph::performProb01Min(model, phiStates, psiStates); - } else { - statesWithProbability01 = storm::utility::graph::performProb01Max(model, phiStates, psiStates); - } - storm::dd::Bdd maybeStates = !statesWithProbability01.first && !statesWithProbability01.second && model.getReachableStates(); - - // Perform some logging. - STORM_LOG_INFO("Found " << statesWithProbability01.first.getNonZeroCount() << " 'no' states."); - STORM_LOG_INFO("Found " << statesWithProbability01.second.getNonZeroCount() << " 'yes' states."); - STORM_LOG_INFO("Found " << maybeStates.getNonZeroCount() << " 'maybe' states."); - - // Check whether we need to compute exact probabilities for some states. - if (qualitative) { - // Set the values for all maybe-states to 0.5 to indicate that their probability values are neither 0 nor 1. - return std::unique_ptr(new storm::modelchecker::SymbolicQuantitativeCheckResult(model.getReachableStates(), statesWithProbability01.second.toAdd() + maybeStates.toAdd() * model.getManager().getConstant(0.5))); - } else { - // If there are maybe states, we need to solve an equation system. - if (!maybeStates.isZero()) { - // Create the ODD for the translation between symbolic and explicit storage. - storm::dd::Odd odd(maybeStates); - - // Create the matrix and the vector for the equation system. - storm::dd::Add maybeStatesAdd = maybeStates.toAdd(); - - // Start by cutting away all rows that do not belong to maybe states. Note that this leaves columns targeting - // non-maybe states in the matrix. - storm::dd::Add submatrix = transitionMatrix * maybeStatesAdd; - - // Then compute the vector that contains the one-step probabilities to a state with probability 1 for all - // maybe states. - storm::dd::Add prob1StatesAsColumn = statesWithProbability01.second.toAdd(); - prob1StatesAsColumn = prob1StatesAsColumn.swapVariables(model.getRowColumnMetaVariablePairs()); - storm::dd::Add subvector = submatrix * prob1StatesAsColumn; - subvector = subvector.sumAbstract(model.getColumnVariables()); - - // Before cutting the non-maybe columns, we need to compute the sizes of the row groups. - std::vector rowGroupSizes = submatrix.notZero().existsAbstract(model.getColumnVariables()).toAdd().sumAbstract(model.getNondeterminismVariables()).template toVector(odd); - - // Finally cut away all columns targeting non-maybe states. - submatrix *= maybeStatesAdd.swapVariables(model.getRowColumnMetaVariablePairs()); - - // Create the solution vector. - std::vector x(maybeStates.getNonZeroCount(), ValueType(0.5)); - - // Translate the symbolic matrix/vector to their explicit representations and solve the equation system. - std::pair, std::vector> explicitRepresentation = submatrix.toMatrixVector(subvector, std::move(rowGroupSizes), model.getNondeterminismVariables(), odd, odd); - - std::unique_ptr> solver = linearEquationSolverFactory.create(explicitRepresentation.first); - solver->solveEquationSystem(minimize, x, explicitRepresentation.second); - - // Return a hybrid check result that stores the numerical values explicitly. - return std::unique_ptr(new storm::modelchecker::HybridQuantitativeCheckResult(model.getReachableStates(), model.getReachableStates() && !maybeStates, statesWithProbability01.second.toAdd(), maybeStates, odd, x)); - } else { - return std::unique_ptr(new storm::modelchecker::SymbolicQuantitativeCheckResult(model.getReachableStates(), statesWithProbability01.second.toAdd())); - } - } - } - + template std::unique_ptr HybridMdpPrctlModelChecker::computeUntilProbabilities(storm::logic::UntilFormula const& pathFormula, bool qualitative, boost::optional const& optimalityType) { STORM_LOG_THROW(optimalityType, storm::exceptions::InvalidArgumentException, "Formula needs to specify whether minimal or maximal values are to be computed on nondeterministic model."); @@ -110,13 +46,7 @@ namespace storm { SymbolicQualitativeCheckResult const& subResult = subResultPointer->asSymbolicQualitativeCheckResult(); return std::unique_ptr(new SymbolicQuantitativeCheckResult(this->getModel().getReachableStates(), this->computeNextProbabilitiesHelper(optimalityType.get() == storm::logic::OptimalityType::Minimize, this->getModel(), this->getModel().getTransitionMatrix(), subResult.getTruthValuesVector()))); } - - template - storm::dd::Add HybridMdpPrctlModelChecker::computeNextProbabilitiesHelper(bool minimize, storm::models::symbolic::NondeterministicModel const& model, storm::dd::Add const& transitionMatrix, storm::dd::Bdd const& nextStates) { - storm::dd::Add result = transitionMatrix * nextStates.swapVariables(model.getRowColumnMetaVariablePairs()).toAdd(); - return result.sumAbstract(model.getColumnVariables()); - } - + template std::unique_ptr HybridMdpPrctlModelChecker::computeBoundedUntilProbabilities(storm::logic::BoundedUntilFormula const& pathFormula, bool qualitative, boost::optional const& optimalityType) { STORM_LOG_THROW(optimalityType, storm::exceptions::InvalidArgumentException, "Formula needs to specify whether minimal or maximal values are to be computed on nondeterministic model."); @@ -129,197 +59,27 @@ namespace storm { } template - std::unique_ptr HybridMdpPrctlModelChecker::computeBoundedUntilProbabilitiesHelper(bool minimize, storm::models::symbolic::NondeterministicModel const& model, storm::dd::Add const& transitionMatrix, storm::dd::Bdd const& phiStates, storm::dd::Bdd const& psiStates, uint_fast64_t stepBound, storm::utility::solver::MinMaxLinearEquationSolverFactory const& linearEquationSolverFactory) { - // We need to identify the states which have to be taken out of the matrix, i.e. all states that have - // probability 0 or 1 of satisfying the until-formula. - storm::dd::Bdd statesWithProbabilityGreater0; - if (minimize) { - statesWithProbabilityGreater0 = storm::utility::graph::performProbGreater0A(model, transitionMatrix.notZero(), phiStates, psiStates); - } else { - statesWithProbabilityGreater0 = storm::utility::graph::performProbGreater0E(model, transitionMatrix.notZero(), phiStates, psiStates); - } - storm::dd::Bdd maybeStates = statesWithProbabilityGreater0 && !psiStates && model.getReachableStates(); - - // If there are maybe states, we need to perform matrix-vector multiplications. - if (!maybeStates.isZero()) { - // Create the ODD for the translation between symbolic and explicit storage. - storm::dd::Odd odd(maybeStates); - - // Create the matrix and the vector for the equation system. - storm::dd::Add maybeStatesAdd = maybeStates.toAdd(); - - // Start by cutting away all rows that do not belong to maybe states. Note that this leaves columns targeting - // non-maybe states in the matrix. - storm::dd::Add submatrix = transitionMatrix * maybeStatesAdd; - - // Then compute the vector that contains the one-step probabilities to a state with probability 1 for all - // maybe states. - storm::dd::Add prob1StatesAsColumn = psiStates.toAdd().swapVariables(model.getRowColumnMetaVariablePairs()); - storm::dd::Add subvector = (submatrix * prob1StatesAsColumn).sumAbstract(model.getColumnVariables()); - - // Before cutting the non-maybe columns, we need to compute the sizes of the row groups. - std::vector rowGroupSizes = submatrix.notZero().existsAbstract(model.getColumnVariables()).toAdd().sumAbstract(model.getNondeterminismVariables()).template toVector(odd); - - // Finally cut away all columns targeting non-maybe states. - submatrix *= maybeStatesAdd.swapVariables(model.getRowColumnMetaVariablePairs()); - - // Create the solution vector. - std::vector x(maybeStates.getNonZeroCount(), storm::utility::zero()); - - // Translate the symbolic matrix/vector to their explicit representations. - std::pair, std::vector> explicitRepresentation = submatrix.toMatrixVector(subvector, std::move(rowGroupSizes), model.getNondeterminismVariables(), odd, odd); - - std::unique_ptr> solver = linearEquationSolverFactory.create(explicitRepresentation.first); - solver->performMatrixVectorMultiplication(minimize, x, &explicitRepresentation.second, stepBound); - - // Return a hybrid check result that stores the numerical values explicitly. - return std::unique_ptr(new storm::modelchecker::HybridQuantitativeCheckResult(model.getReachableStates(), model.getReachableStates() && !maybeStates, psiStates.toAdd(), maybeStates, odd, x)); - } else { - return std::unique_ptr(new storm::modelchecker::SymbolicQuantitativeCheckResult(model.getReachableStates(), psiStates.toAdd())); - } - } - - template - std::unique_ptr HybridMdpPrctlModelChecker::computeCumulativeRewards(storm::logic::CumulativeRewardFormula const& rewardPathFormula, bool qualitative, boost::optional const& optimalityType) { + std::unique_ptr HybridMdpPrctlModelChecker::computeCumulativeRewards(storm::logic::CumulativeRewardFormula const& rewardPathFormula, boost::optional const& rewardModelName, bool qualitative, boost::optional const& optimalityType) { STORM_LOG_THROW(optimalityType, storm::exceptions::InvalidArgumentException, "Formula needs to specify whether minimal or maximal values are to be computed on nondeterministic model."); STORM_LOG_THROW(rewardPathFormula.hasDiscreteTimeBound(), storm::exceptions::InvalidArgumentException, "Formula needs to have a discrete time bound."); return this->computeCumulativeRewardsHelper(optimalityType.get() == storm::logic::OptimalityType::Minimize, this->getModel(), this->getModel().getTransitionMatrix(), rewardPathFormula.getDiscreteTimeBound(), *this->linearEquationSolverFactory); } - - template - std::unique_ptr HybridMdpPrctlModelChecker::computeCumulativeRewardsHelper(bool minimize, storm::models::symbolic::NondeterministicModel const& model, storm::dd::Add const& transitionMatrix, uint_fast64_t stepBound, storm::utility::solver::MinMaxLinearEquationSolverFactory const& linearEquationSolverFactory) { - // Only compute the result if the model has at least one reward this->getModel(). - STORM_LOG_THROW(model.hasStateRewards() || model.hasTransitionRewards(), storm::exceptions::InvalidPropertyException, "Missing reward model for formula. Skipping formula."); - - // Compute the reward vector to add in each step based on the available reward models. - storm::dd::Add totalRewardVector = model.hasStateRewards() ? model.getStateRewardVector() : model.getManager().getAddZero(); - if (model.hasTransitionRewards()) { - totalRewardVector += (transitionMatrix * model.getTransitionRewardMatrix()).sumAbstract(model.getColumnVariables()); - } - - // Create the ODD for the translation between symbolic and explicit storage. - storm::dd::Odd odd(model.getReachableStates()); - - // Create the solution vector. - std::vector x(model.getNumberOfStates(), storm::utility::zero()); - - // Translate the symbolic matrix/vector to their explicit representations. - storm::storage::SparseMatrix explicitMatrix = transitionMatrix.toMatrix(model.getNondeterminismVariables(), odd, odd); - std::vector b = totalRewardVector.template toVector(model.getNondeterminismVariables(), odd, explicitMatrix.getRowGroupIndices()); - - // Perform the matrix-vector multiplication. - std::unique_ptr> solver = linearEquationSolverFactory.create(explicitMatrix); - solver->performMatrixVectorMultiplication(minimize, x, &b, stepBound); - - // Return a hybrid check result that stores the numerical values explicitly. - return std::unique_ptr(new HybridQuantitativeCheckResult(model.getReachableStates(), model.getManager().getBddZero(), model.getManager().getAddZero(), model.getReachableStates(), odd, x)); - } - + template - std::unique_ptr HybridMdpPrctlModelChecker::computeInstantaneousRewards(storm::logic::InstantaneousRewardFormula const& rewardPathFormula, bool qualitative, boost::optional const& optimalityType) { + std::unique_ptr HybridMdpPrctlModelChecker::computeInstantaneousRewards(storm::logic::InstantaneousRewardFormula const& rewardPathFormula, boost::optional const& rewardModelName, bool qualitative, boost::optional const& optimalityType) { STORM_LOG_THROW(optimalityType, storm::exceptions::InvalidArgumentException, "Formula needs to specify whether minimal or maximal values are to be computed on nondeterministic model."); STORM_LOG_THROW(rewardPathFormula.hasDiscreteTimeBound(), storm::exceptions::InvalidArgumentException, "Formula needs to have a discrete time bound."); return this->computeInstantaneousRewardsHelper(optimalityType.get() == storm::logic::OptimalityType::Minimize, this->getModel(), this->getModel().getTransitionMatrix(), rewardPathFormula.getDiscreteTimeBound(), *this->linearEquationSolverFactory); } - - template - std::unique_ptr HybridMdpPrctlModelChecker::computeInstantaneousRewardsHelper(bool minimize, storm::models::symbolic::NondeterministicModel const& model, storm::dd::Add const& transitionMatrix, uint_fast64_t stepBound, storm::utility::solver::MinMaxLinearEquationSolverFactory const& linearEquationSolverFactory) { - // Only compute the result if the model has at least one reward this->getModel(). - STORM_LOG_THROW(model.hasStateRewards(), storm::exceptions::InvalidPropertyException, "Missing reward model for formula. Skipping formula."); - - // Create the ODD for the translation between symbolic and explicit storage. - storm::dd::Odd odd(model.getReachableStates()); - - // Translate the symbolic matrix to its explicit representations. - storm::storage::SparseMatrix explicitMatrix = transitionMatrix.toMatrix(model.getNondeterminismVariables(), odd, odd); - - // Create the solution vector (and initialize it to the state rewards of the model). - std::vector x = model.getStateRewardVector().template toVector(model.getNondeterminismVariables(), odd, explicitMatrix.getRowGroupIndices()); - - // Perform the matrix-vector multiplication. - std::unique_ptr> solver = linearEquationSolverFactory.create(explicitMatrix); - solver->performMatrixVectorMultiplication(minimize, x, nullptr, stepBound); - - // Return a hybrid check result that stores the numerical values explicitly. - return std::unique_ptr(new HybridQuantitativeCheckResult(model.getReachableStates(), model.getManager().getBddZero(), model.getManager().getAddZero(), model.getReachableStates(), odd, x)); - } - + template - std::unique_ptr HybridMdpPrctlModelChecker::computeReachabilityRewards(storm::logic::ReachabilityRewardFormula const& rewardPathFormula, bool qualitative, boost::optional const& optimalityType) { + std::unique_ptr HybridMdpPrctlModelChecker::computeReachabilityRewards(storm::logic::ReachabilityRewardFormula const& rewardPathFormula, boost::optional const& rewardModelName, bool qualitative, boost::optional const& optimalityType) { STORM_LOG_THROW(optimalityType, storm::exceptions::InvalidArgumentException, "Formula needs to specify whether minimal or maximal values are to be computed on nondeterministic model."); std::unique_ptr subResultPointer = this->check(rewardPathFormula.getSubformula()); SymbolicQualitativeCheckResult const& subResult = subResultPointer->asSymbolicQualitativeCheckResult(); return this->computeReachabilityRewardsHelper(optimalityType.get() == storm::logic::OptimalityType::Minimize, this->getModel(), this->getModel().getTransitionMatrix(), this->getModel().getOptionalStateRewardVector(), this->getModel().getOptionalTransitionRewardMatrix(), subResult.getTruthValuesVector(), *this->linearEquationSolverFactory, qualitative); } - template - std::unique_ptr HybridMdpPrctlModelChecker::computeReachabilityRewardsHelper(bool minimize, storm::models::symbolic::NondeterministicModel const& model, storm::dd::Add const& transitionMatrix, boost::optional> const& stateRewardVector, boost::optional> const& transitionRewardMatrix, storm::dd::Bdd const& targetStates, storm::utility::solver::MinMaxLinearEquationSolverFactory const& linearEquationSolverFactory, bool qualitative) { - - // Only compute the result if there is at least one reward model. - STORM_LOG_THROW(stateRewardVector || transitionRewardMatrix, storm::exceptions::InvalidPropertyException, "Missing reward model for formula. Skipping formula."); - - // Determine which states have a reward of infinity by definition. - storm::dd::Bdd infinityStates; - storm::dd::Bdd transitionMatrixBdd = transitionMatrix.notZero(); - if (minimize) { - infinityStates = storm::utility::graph::performProb1A(model, transitionMatrixBdd, model.getReachableStates(), targetStates, storm::utility::graph::performProbGreater0A(model, transitionMatrixBdd, model.getReachableStates(), targetStates)); - } else { - infinityStates = storm::utility::graph::performProb1E(model, transitionMatrixBdd, model.getReachableStates(), targetStates, storm::utility::graph::performProbGreater0E(model, transitionMatrixBdd, model.getReachableStates(), targetStates)); - } - infinityStates = !infinityStates && model.getReachableStates(); - storm::dd::Bdd maybeStates = (!targetStates && !infinityStates) && model.getReachableStates(); - STORM_LOG_INFO("Found " << infinityStates.getNonZeroCount() << " 'infinity' states."); - STORM_LOG_INFO("Found " << targetStates.getNonZeroCount() << " 'target' states."); - STORM_LOG_INFO("Found " << maybeStates.getNonZeroCount() << " 'maybe' states."); - - // Check whether we need to compute exact rewards for some states. - if (qualitative) { - // Set the values for all maybe-states to 1 to indicate that their reward values - // are neither 0 nor infinity. - return std::unique_ptr(new SymbolicQuantitativeCheckResult(model.getReachableStates(), infinityStates.toAdd() * model.getManager().getConstant(storm::utility::infinity()) + maybeStates.toAdd() * model.getManager().getConstant(storm::utility::one()))); - } else { - // If there are maybe states, we need to solve an equation system. - if (!maybeStates.isZero()) { - // Create the ODD for the translation between symbolic and explicit storage. - storm::dd::Odd odd(maybeStates); - - // Create the matrix and the vector for the equation system. - storm::dd::Add maybeStatesAdd = maybeStates.toAdd(); - - // Start by cutting away all rows that do not belong to maybe states. Note that this leaves columns targeting - // non-maybe states in the matrix. - storm::dd::Add submatrix = transitionMatrix * maybeStatesAdd; - - // Then compute the state reward vector to use in the computation. - storm::dd::Add subvector = stateRewardVector ? maybeStatesAdd * stateRewardVector.get() : model.getManager().getAddZero(); - if (transitionRewardMatrix) { - subvector += (submatrix * transitionRewardMatrix.get()).sumAbstract(model.getColumnVariables()); - } - - // Before cutting the non-maybe columns, we need to compute the sizes of the row groups. - std::vector rowGroupSizes = submatrix.notZero().existsAbstract(model.getColumnVariables()).toAdd().sumAbstract(model.getNondeterminismVariables()).template toVector(odd); - - // Finally cut away all columns targeting non-maybe states. - submatrix *= maybeStatesAdd.swapVariables(model.getRowColumnMetaVariablePairs()); - - // Create the solution vector. - std::vector x(maybeStates.getNonZeroCount(), ValueType(0.5)); - - // Translate the symbolic matrix/vector to their explicit representations. - std::pair, std::vector> explicitRepresentation = submatrix.toMatrixVector(subvector, std::move(rowGroupSizes), model.getNondeterminismVariables(), odd, odd); - - // Now solve the resulting equation system. - std::unique_ptr> solver = linearEquationSolverFactory.create(explicitRepresentation.first); - solver->solveEquationSystem(minimize, x, explicitRepresentation.second); - - // Return a hybrid check result that stores the numerical values explicitly. - return std::unique_ptr(new storm::modelchecker::HybridQuantitativeCheckResult(model.getReachableStates(), model.getReachableStates() && !maybeStates, infinityStates.toAdd() * model.getManager().getConstant(storm::utility::infinity()), maybeStates, odd, x)); - } else { - return std::unique_ptr(new storm::modelchecker::SymbolicQuantitativeCheckResult(model.getReachableStates(), infinityStates.toAdd() * model.getManager().getConstant(storm::utility::infinity()))); - } - } - } - template storm::models::symbolic::Mdp const& HybridMdpPrctlModelChecker::getModel() const { return this->template getModelAs>(); diff --git a/src/modelchecker/prctl/HybridMdpPrctlModelChecker.h b/src/modelchecker/prctl/HybridMdpPrctlModelChecker.h index 24a33ee79..ee22ff766 100644 --- a/src/modelchecker/prctl/HybridMdpPrctlModelChecker.h +++ b/src/modelchecker/prctl/HybridMdpPrctlModelChecker.h @@ -18,22 +18,14 @@ namespace storm { virtual std::unique_ptr computeBoundedUntilProbabilities(storm::logic::BoundedUntilFormula const& pathFormula, bool qualitative = false, boost::optional const& optimalityType = boost::optional()) override; virtual std::unique_ptr computeNextProbabilities(storm::logic::NextFormula const& pathFormula, bool qualitative = false, boost::optional const& optimalityType = boost::optional()) override; virtual std::unique_ptr computeUntilProbabilities(storm::logic::UntilFormula const& pathFormula, bool qualitative = false, boost::optional const& optimalityType = boost::optional()) override; - virtual std::unique_ptr computeCumulativeRewards(storm::logic::CumulativeRewardFormula const& rewardPathFormula, bool qualitative = false, boost::optional const& optimalityType = boost::optional()) override; - virtual std::unique_ptr computeInstantaneousRewards(storm::logic::InstantaneousRewardFormula const& rewardPathFormula, bool qualitative = false, boost::optional const& optimalityType = boost::optional()) override; - virtual std::unique_ptr computeReachabilityRewards(storm::logic::ReachabilityRewardFormula const& rewardPathFormula, bool qualitative = false, boost::optional const& optimalityType = boost::optional()) override; + virtual std::unique_ptr computeCumulativeRewards(storm::logic::CumulativeRewardFormula const& rewardPathFormula, boost::optional const& rewardModelName = boost::optional(), bool qualitative = false, boost::optional const& optimalityType = boost::optional()) override; + virtual std::unique_ptr computeInstantaneousRewards(storm::logic::InstantaneousRewardFormula const& rewardPathFormula, boost::optional const& rewardModelName = boost::optional(), bool qualitative = false, boost::optional const& optimalityType = boost::optional()) override; + virtual std::unique_ptr computeReachabilityRewards(storm::logic::ReachabilityRewardFormula const& rewardPathFormula, boost::optional const& rewardModelName = boost::optional(), bool qualitative = false, boost::optional const& optimalityType = boost::optional()) override; protected: storm::models::symbolic::Mdp const& getModel() const override; private: - // The methods that perform the actual checking. - static std::unique_ptr computeBoundedUntilProbabilitiesHelper(bool minimize, storm::models::symbolic::NondeterministicModel const& model, storm::dd::Add const& transitionMatrix, storm::dd::Bdd const& phiStates, storm::dd::Bdd const& psiStates, uint_fast64_t stepBound, storm::utility::solver::MinMaxLinearEquationSolverFactory const& linearEquationSolverFactory); - static storm::dd::Add computeNextProbabilitiesHelper(bool minimize, storm::models::symbolic::NondeterministicModel const& model, storm::dd::Add const& transitionMatrix, storm::dd::Bdd const& nextStates); - static std::unique_ptr computeUntilProbabilitiesHelper(bool minimize, storm::models::symbolic::NondeterministicModel const& model, storm::dd::Add const& transitionMatrix, storm::dd::Bdd const& phiStates, storm::dd::Bdd const& psiStates, bool qualitative, storm::utility::solver::MinMaxLinearEquationSolverFactory const& linearEquationSolverFactory); - static std::unique_ptr computeCumulativeRewardsHelper(bool minimize, storm::models::symbolic::NondeterministicModel const& model, storm::dd::Add const& transitionMatrix, uint_fast64_t stepBound, storm::utility::solver::MinMaxLinearEquationSolverFactory const& linearEquationSolverFactory); - static std::unique_ptr computeInstantaneousRewardsHelper(bool minimize, storm::models::symbolic::NondeterministicModel const& model, storm::dd::Add const& transitionMatrix, uint_fast64_t stepBound, storm::utility::solver::MinMaxLinearEquationSolverFactory const& linearEquationSolverFactory); - static std::unique_ptr computeReachabilityRewardsHelper(bool minimize, storm::models::symbolic::NondeterministicModel const& model, storm::dd::Add const& transitionMatrix, boost::optional> const& stateRewardVector, boost::optional> const& transitionRewardMatrix, storm::dd::Bdd const& targetStates, storm::utility::solver::MinMaxLinearEquationSolverFactory const& linearEquationSolverFactory, bool qualitative); - // An object that is used for retrieving linear equation solvers. std::unique_ptr> linearEquationSolverFactory; }; diff --git a/src/modelchecker/prctl/SparseDtmcPrctlModelChecker.cpp b/src/modelchecker/prctl/SparseDtmcPrctlModelChecker.cpp index d8d2f1378..ec9743ffe 100644 --- a/src/modelchecker/prctl/SparseDtmcPrctlModelChecker.cpp +++ b/src/modelchecker/prctl/SparseDtmcPrctlModelChecker.cpp @@ -4,17 +4,11 @@ #include #include "src/utility/macros.h" -#include "src/utility/vector.h" -#include "src/utility/graph.h" -#include "src/utility/solver.h" -#include "src/modelchecker/csl/SparseCtmcCslModelChecker.h" #include "src/modelchecker/results/ExplicitQualitativeCheckResult.h" #include "src/modelchecker/results/ExplicitQuantitativeCheckResult.h" -#include "src/exceptions/InvalidStateException.h" - -#include "src/exceptions/InvalidPropertyException.h" +#include "src/modelchecker/prctl/helper/SparseDtmcPrctlHelper.h" namespace storm { namespace modelchecker { @@ -33,58 +27,15 @@ namespace storm { return formula.isPctlStateFormula() || formula.isPctlPathFormula() || formula.isRewardPathFormula(); } - template - std::vector::ValueType> SparseDtmcPrctlModelChecker::computeBoundedUntilProbabilitiesHelper(storm::storage::BitVector const& phiStates, storm::storage::BitVector const& psiStates, uint_fast64_t stepBound) const { - std::vector result(this->getModel().getNumberOfStates(), storm::utility::zero()); - - // If we identify the states that have probability 0 of reaching the target states, we can exclude them in the further analysis. - storm::storage::BitVector maybeStates = storm::utility::graph::performProbGreater0(this->getModel().getBackwardTransitions(), phiStates, psiStates, true, stepBound); - maybeStates &= ~psiStates; - STORM_LOG_INFO("Found " << maybeStates.getNumberOfSetBits() << " 'maybe' states."); - - if (!maybeStates.empty()) { - // We can eliminate the rows and columns from the original transition probability matrix that have probability 0. - storm::storage::SparseMatrix submatrix = this->getModel().getTransitionMatrix().getSubmatrix(true, maybeStates, maybeStates, true); - - // Create the vector of one-step probabilities to go to target states. - std::vector b = this->getModel().getTransitionMatrix().getConstrainedRowSumVector(maybeStates, psiStates); - - // Create the vector with which to multiply. - std::vector subresult(maybeStates.getNumberOfSetBits()); - - // Perform the matrix vector multiplication as often as required by the formula bound. - STORM_LOG_THROW(linearEquationSolverFactory != nullptr, storm::exceptions::InvalidStateException, "No valid linear equation solver available."); - std::unique_ptr> solver = linearEquationSolverFactory->create(submatrix); - solver->performMatrixVectorMultiplication(subresult, &b, stepBound); - - // Set the values of the resulting vector accordingly. - storm::utility::vector::setVectorValues(result, maybeStates, subresult); - } - storm::utility::vector::setVectorValues(result, psiStates, storm::utility::one()); - - return result; - } - template std::unique_ptr SparseDtmcPrctlModelChecker::computeBoundedUntilProbabilities(storm::logic::BoundedUntilFormula const& pathFormula, bool qualitative, boost::optional const& optimalityType) { STORM_LOG_THROW(pathFormula.hasDiscreteTimeBound(), storm::exceptions::InvalidArgumentException, "Formula needs to have a discrete time bound."); std::unique_ptr leftResultPointer = this->check(pathFormula.getLeftSubformula()); std::unique_ptr rightResultPointer = this->check(pathFormula.getRightSubformula()); - ExplicitQualitativeCheckResult const& leftResult = leftResultPointer->asExplicitQualitativeCheckResult();; + ExplicitQualitativeCheckResult const& leftResult = leftResultPointer->asExplicitQualitativeCheckResult(); ExplicitQualitativeCheckResult const& rightResult = rightResultPointer->asExplicitQualitativeCheckResult(); - std::unique_ptr result = std::unique_ptr(new ExplicitQuantitativeCheckResult(this->computeBoundedUntilProbabilitiesHelper(leftResult.getTruthValuesVector(), rightResult.getTruthValuesVector(), pathFormula.getDiscreteTimeBound()))); - return result; - } - - template - std::vector::ValueType> SparseDtmcPrctlModelChecker::computeNextProbabilitiesHelper(storm::storage::SparseMatrix const& transitionMatrix, storm::storage::BitVector const& nextStates, storm::utility::solver::LinearEquationSolverFactory const& linearEquationSolverFactory) { - // Create the vector with which to multiply and initialize it correctly. - std::vector result(transitionMatrix.getRowCount()); - storm::utility::vector::setVectorValues(result, nextStates, storm::utility::one()); - - // Perform one single matrix-vector multiplication. - std::unique_ptr> solver = linearEquationSolverFactory.create(transitionMatrix); - solver->performMatrixVectorMultiplication(result); + std::vector numericResult = storm::modelchecker::helper::SparseDtmcPrctlHelper::computeBoundedUntilProbabilities(this->getModel().getTransitionMatrix(), this->getModel().getBackwardTransitions(), leftResult.getTruthValuesVector(), rightResult.getTruthValuesVector(), pathFormula.getDiscreteTimeBound(), *linearEquationSolverFactory); + std::unique_ptr result = std::unique_ptr(new ExplicitQuantitativeCheckResult(std::move(numericResult))); return result; } @@ -92,189 +43,50 @@ namespace storm { std::unique_ptr SparseDtmcPrctlModelChecker::computeNextProbabilities(storm::logic::NextFormula const& pathFormula, bool qualitative, boost::optional const& optimalityType) { std::unique_ptr subResultPointer = this->check(pathFormula.getSubformula()); ExplicitQualitativeCheckResult const& subResult = subResultPointer->asExplicitQualitativeCheckResult(); - return std::unique_ptr(new ExplicitQuantitativeCheckResult(this->computeNextProbabilitiesHelper(this->getModel().getTransitionMatrix(), subResult.getTruthValuesVector(), *this->linearEquationSolverFactory))); - } - - template - std::vector::ValueType> SparseDtmcPrctlModelChecker::computeUntilProbabilitiesHelper(storm::storage::SparseMatrix const& transitionMatrix, storm::storage::SparseMatrix const& backwardTransitions, storm::storage::BitVector const& phiStates, storm::storage::BitVector const& psiStates, bool qualitative, storm::utility::solver::LinearEquationSolverFactory const& linearEquationSolverFactory) { - // We need to identify the states which have to be taken out of the matrix, i.e. - // all states that have probability 0 and 1 of satisfying the until-formula. - std::pair statesWithProbability01 = storm::utility::graph::performProb01(backwardTransitions, phiStates, psiStates); - storm::storage::BitVector statesWithProbability0 = std::move(statesWithProbability01.first); - storm::storage::BitVector statesWithProbability1 = std::move(statesWithProbability01.second); - - // Perform some logging. - storm::storage::BitVector maybeStates = ~(statesWithProbability0 | statesWithProbability1); - STORM_LOG_INFO("Found " << statesWithProbability0.getNumberOfSetBits() << " 'no' states."); - STORM_LOG_INFO("Found " << statesWithProbability1.getNumberOfSetBits() << " 'yes' states."); - STORM_LOG_INFO("Found " << maybeStates.getNumberOfSetBits() << " 'maybe' states."); - - // Create resulting vector. - std::vector result(transitionMatrix.getRowCount()); - - // Check whether we need to compute exact probabilities for some states. - if (qualitative) { - // Set the values for all maybe-states to 0.5 to indicate that their probability values are neither 0 nor 1. - storm::utility::vector::setVectorValues(result, maybeStates, ValueType(0.5)); - } else { - if (!maybeStates.empty()) { - // In this case we have have to compute the probabilities. - - // We can eliminate the rows and columns from the original transition probability matrix. - storm::storage::SparseMatrix submatrix = transitionMatrix.getSubmatrix(true, maybeStates, maybeStates, true); - - // Converting the matrix from the fixpoint notation to the form needed for the equation - // system. That is, we go from x = A*x + b to (I-A)x = b. - submatrix.convertToEquationSystem(); - - // Initialize the x vector with 0.5 for each element. This is the initial guess for - // the iterative solvers. It should be safe as for all 'maybe' states we know that the - // probability is strictly larger than 0. - std::vector x(maybeStates.getNumberOfSetBits(), ValueType(0.5)); - - // Prepare the right-hand side of the equation system. For entry i this corresponds to - // the accumulated probability of going from state i to some 'yes' state. - std::vector b = transitionMatrix.getConstrainedRowSumVector(maybeStates, statesWithProbability1); - - // Now solve the created system of linear equations. - std::unique_ptr> solver = linearEquationSolverFactory.create(submatrix); - solver->solveEquationSystem(x, b); - - // Set values of resulting vector according to result. - storm::utility::vector::setVectorValues(result, maybeStates, x); - } - } - - // Set values of resulting vector that are known exactly. - storm::utility::vector::setVectorValues(result, statesWithProbability0, storm::utility::zero()); - storm::utility::vector::setVectorValues(result, statesWithProbability1, storm::utility::one()); - - return result; + std::vector numericResult = storm::modelchecker::helper::SparseDtmcPrctlHelper::computeNextProbabilities(this->getModel().getTransitionMatrix(), subResult.getTruthValuesVector(), *linearEquationSolverFactory); + return std::unique_ptr(new ExplicitQuantitativeCheckResult(std::move(numericResult))); } template std::unique_ptr SparseDtmcPrctlModelChecker::computeUntilProbabilities(storm::logic::UntilFormula const& pathFormula, bool qualitative, boost::optional const& optimalityType) { std::unique_ptr leftResultPointer = this->check(pathFormula.getLeftSubformula()); std::unique_ptr rightResultPointer = this->check(pathFormula.getRightSubformula()); - ExplicitQualitativeCheckResult const& leftResult = leftResultPointer->asExplicitQualitativeCheckResult();; - ExplicitQualitativeCheckResult const& rightResult = rightResultPointer->asExplicitQualitativeCheckResult();; - return std::unique_ptr(new ExplicitQuantitativeCheckResult(this->computeUntilProbabilitiesHelper(this->getModel().getTransitionMatrix(), this->getModel().getBackwardTransitions(), leftResult.getTruthValuesVector(), rightResult.getTruthValuesVector(), qualitative, *this->linearEquationSolverFactory))); - } - - template - std::vector::ValueType> SparseDtmcPrctlModelChecker::computeCumulativeRewardsHelper(RewardModelType const& rewardModel, uint_fast64_t stepBound) const { - // Compute the reward vector to add in each step based on the available reward models. - std::vector totalRewardVector = rewardModel.getTotalRewardVector(this->getModel().getTransitionMatrix()); - - // Initialize result to either the state rewards of the model or the null vector. - std::vector result = rewardModel.getTotalStateActionRewardVector(this->getModel().getNumberOfStates(), this->getModel().getTransitionMatrix().getRowGroupIndices()); - - // Perform the matrix vector multiplication as often as required by the formula bound. - STORM_LOG_THROW(linearEquationSolverFactory != nullptr, storm::exceptions::InvalidStateException, "No valid linear equation solver available."); - std::unique_ptr> solver = linearEquationSolverFactory->create(this->getModel().getTransitionMatrix()); - solver->performMatrixVectorMultiplication(result, &totalRewardVector, stepBound); - - return result; + ExplicitQualitativeCheckResult const& leftResult = leftResultPointer->asExplicitQualitativeCheckResult(); + ExplicitQualitativeCheckResult const& rightResult = rightResultPointer->asExplicitQualitativeCheckResult(); + std::vector numericResult = storm::modelchecker::helper::SparseDtmcPrctlHelper::computeUntilProbabilities(this->getModel().getTransitionMatrix(), this->getModel().getBackwardTransitions(), leftResult.getTruthValuesVector(), rightResult.getTruthValuesVector(), qualitative, *linearEquationSolverFactory); + return std::unique_ptr(new ExplicitQuantitativeCheckResult(std::move(numericResult))); } template std::unique_ptr SparseDtmcPrctlModelChecker::computeCumulativeRewards(storm::logic::CumulativeRewardFormula const& rewardPathFormula, boost::optional const& rewardModelName, bool qualitative, boost::optional const& optimalityType) { STORM_LOG_THROW(rewardPathFormula.hasDiscreteTimeBound(), storm::exceptions::InvalidArgumentException, "Formula needs to have a discrete time bound."); - return std::unique_ptr(new ExplicitQuantitativeCheckResult(this->computeCumulativeRewardsHelper(rewardModelName ? this->getModel().getRewardModel(rewardModelName.get()) : this->getModel().getRewardModel(""), rewardPathFormula.getDiscreteTimeBound()))); - } - - template - std::vector::ValueType> SparseDtmcPrctlModelChecker::computeInstantaneousRewardsHelper(RewardModelType const& rewardModel, uint_fast64_t stepCount) const { - // Only compute the result if the model has a state-based reward this->getModel(). - STORM_LOG_THROW(rewardModel.hasStateRewards() || rewardModel.hasStateActionRewards(), storm::exceptions::InvalidPropertyException, "Missing reward model for formula. Skipping formula."); - - // Initialize result to state rewards of the model. - std::vector result = rewardModel.getTotalStateActionRewardVector(this->getModel().getNumberOfStates(), this->getModel().getTransitionMatrix().getRowGroupIndices()); - - // Perform the matrix vector multiplication as often as required by the formula bound. - STORM_LOG_THROW(linearEquationSolverFactory != nullptr, storm::exceptions::InvalidStateException, "No valid linear equation solver available."); - std::unique_ptr> solver = linearEquationSolverFactory->create(this->getModel().getTransitionMatrix()); - solver->performMatrixVectorMultiplication(result, nullptr, stepCount); - - return result; + std::vector numericResult = storm::modelchecker::helper::SparseDtmcPrctlHelper::computeCumulativeRewards(this->getModel().getTransitionMatrix(), rewardModelName ? this->getModel().getRewardModel(rewardModelName.get()) : this->getModel().getRewardModel(""), rewardPathFormula.getDiscreteTimeBound(), *linearEquationSolverFactory); + return std::unique_ptr(new ExplicitQuantitativeCheckResult(std::move(numericResult))); } template std::unique_ptr SparseDtmcPrctlModelChecker::computeInstantaneousRewards(storm::logic::InstantaneousRewardFormula const& rewardPathFormula, boost::optional const& rewardModelName, bool qualitative, boost::optional const& optimalityType) { STORM_LOG_THROW(rewardPathFormula.hasDiscreteTimeBound(), storm::exceptions::InvalidArgumentException, "Formula needs to have a discrete time bound."); - return std::unique_ptr(new ExplicitQuantitativeCheckResult(this->computeInstantaneousRewardsHelper(rewardModelName ? this->getModel().getRewardModel(rewardModelName.get()) : this->getModel().getRewardModel(""), rewardPathFormula.getDiscreteTimeBound()))); - } - - template - std::vector::ValueType> SparseDtmcPrctlModelChecker::computeReachabilityRewardsHelper(RewardModelType const& rewardModel, storm::storage::SparseMatrix const& transitionMatrix, storm::storage::SparseMatrix const& backwardTransitions, storm::storage::BitVector const& targetStates, storm::utility::solver::LinearEquationSolverFactory const& linearEquationSolverFactory, bool qualitative) { - - // Determine which states have a reward of infinity by definition. - storm::storage::BitVector trueStates(transitionMatrix.getRowCount(), true); - storm::storage::BitVector infinityStates = storm::utility::graph::performProb1(backwardTransitions, trueStates, targetStates); - infinityStates.complement(); - storm::storage::BitVector maybeStates = ~targetStates & ~infinityStates; - STORM_LOG_INFO("Found " << infinityStates.getNumberOfSetBits() << " 'infinity' states."); - STORM_LOG_INFO("Found " << targetStates.getNumberOfSetBits() << " 'target' states."); - STORM_LOG_INFO("Found " << maybeStates.getNumberOfSetBits() << " 'maybe' states."); - - // Create resulting vector. - std::vector result(transitionMatrix.getRowCount(), storm::utility::zero()); - - // Check whether we need to compute exact rewards for some states. - if (qualitative) { - // Set the values for all maybe-states to 1 to indicate that their reward values - // are neither 0 nor infinity. - storm::utility::vector::setVectorValues(result, maybeStates, storm::utility::one()); - } else { - // In this case we have to compute the reward values for the remaining states. - // We can eliminate the rows and columns from the original transition probability matrix. - storm::storage::SparseMatrix submatrix = transitionMatrix.getSubmatrix(true, maybeStates, maybeStates, true); - - // Converting the matrix from the fixpoint notation to the form needed for the equation - // system. That is, we go from x = A*x + b to (I-A)x = b. - submatrix.convertToEquationSystem(); - - // Initialize the x vector with 1 for each element. This is the initial guess for - // the iterative solvers. - std::vector x(submatrix.getColumnCount(), storm::utility::one()); - - // Prepare the right-hand side of the equation system. - std::vector b = rewardModel.getTotalRewardVector(submatrix.getRowCount(), transitionMatrix, maybeStates); - - // Now solve the resulting equation system. - std::unique_ptr> solver = linearEquationSolverFactory.create(submatrix); - solver->solveEquationSystem(x, b); - - // Set values of resulting vector according to result. - storm::utility::vector::setVectorValues(result, maybeStates, x); - } - - // Set values of resulting vector that are known exactly. - storm::utility::vector::setVectorValues(result, infinityStates, storm::utility::infinity()); - - return result; + std::vector numericResult = storm::modelchecker::helper::SparseDtmcPrctlHelper::computeInstantaneousRewards(this->getModel().getTransitionMatrix(), rewardModelName ? this->getModel().getRewardModel(rewardModelName.get()) : this->getModel().getRewardModel(""), rewardPathFormula.getDiscreteTimeBound(), *linearEquationSolverFactory); + return std::unique_ptr(new ExplicitQuantitativeCheckResult(std::move(numericResult))); } template std::unique_ptr SparseDtmcPrctlModelChecker::computeReachabilityRewards(storm::logic::ReachabilityRewardFormula const& rewardPathFormula, boost::optional const& rewardModelName, bool qualitative, boost::optional const& optimalityType) { std::unique_ptr subResultPointer = this->check(rewardPathFormula.getSubformula()); ExplicitQualitativeCheckResult const& subResult = subResultPointer->asExplicitQualitativeCheckResult(); - return std::unique_ptr(new ExplicitQuantitativeCheckResult(this->computeReachabilityRewardsHelper(rewardModelName ? this->getModel().getRewardModel(rewardModelName.get()) : this->getModel().getRewardModel(""), this->getModel().getTransitionMatrix(), this->getModel().getBackwardTransitions(), subResult.getTruthValuesVector(), *this->linearEquationSolverFactory, qualitative))); + std::vector numericResult = storm::modelchecker::helper::SparseDtmcPrctlHelper::computeReachabilityRewards(this->getModel().getTransitionMatrix(), this->getModel().getBackwardTransitions(), rewardModelName ? this->getModel().getRewardModel(rewardModelName.get()) : this->getModel().getRewardModel(""), subResult.getTruthValuesVector(), qualitative, *linearEquationSolverFactory); + return std::unique_ptr(new ExplicitQuantitativeCheckResult(std::move(numericResult))); } template std::unique_ptr SparseDtmcPrctlModelChecker::computeLongRunAverage(storm::logic::StateFormula const& stateFormula, bool qualitative, boost::optional const& optimalityType) { std::unique_ptr subResultPointer = this->check(stateFormula); ExplicitQualitativeCheckResult const& subResult = subResultPointer->asExplicitQualitativeCheckResult(); - - return std::unique_ptr(new ExplicitQuantitativeCheckResult(computeLongRunAverageHelper(this->getModel().getTransitionMatrix(), subResult.getTruthValuesVector(), qualitative, *linearEquationSolverFactory))); + std::vector numericResult = storm::modelchecker::helper::SparseCtmcCslHelper::computeLongRunAverage(this->getModel().getTransitionMatrix(), this->getModel().getTransitionMatrix(), subResult.getTruthValuesVector(), qualitative, *linearEquationSolverFactory); + return std::unique_ptr(new ExplicitQuantitativeCheckResult(std::move(numericResult))); } - - template - std::vector::ValueType> SparseDtmcPrctlModelChecker::computeLongRunAverageHelper(storm::storage::SparseMatrix const& transitionMatrix, storm::storage::BitVector const& psiStates, bool qualitative, storm::utility::solver::LinearEquationSolverFactory const& linearEquationSolverFactory) { - return SparseCtmcCslModelChecker::computeLongRunAverageHelper(transitionMatrix, psiStates, nullptr, qualitative, linearEquationSolverFactory); - } - + template class SparseDtmcPrctlModelChecker>; } } \ No newline at end of file diff --git a/src/modelchecker/prctl/SparseDtmcPrctlModelChecker.h b/src/modelchecker/prctl/SparseDtmcPrctlModelChecker.h index f1d9744fd..cc574d289 100644 --- a/src/modelchecker/prctl/SparseDtmcPrctlModelChecker.h +++ b/src/modelchecker/prctl/SparseDtmcPrctlModelChecker.h @@ -12,9 +12,6 @@ namespace storm { template class HybridDtmcPrctlModelChecker; - // Forward-declare CTMC model checker so we can make it a friend. - template class SparseCtmcCslModelChecker; - template class SparseDtmcPrctlModelChecker : public SparsePropositionalModelChecker { public: @@ -22,7 +19,6 @@ namespace storm { typedef typename SparseDtmcModelType::RewardModelType RewardModelType; friend class HybridDtmcPrctlModelChecker; - friend class SparseCtmcCslModelChecker; explicit SparseDtmcPrctlModelChecker(SparseDtmcModelType const& model); explicit SparseDtmcPrctlModelChecker(SparseDtmcModelType const& model, std::unique_ptr>&& linearEquationSolverFactory); @@ -32,21 +28,12 @@ namespace storm { virtual std::unique_ptr computeBoundedUntilProbabilities(storm::logic::BoundedUntilFormula const& pathFormula, bool qualitative = false, boost::optional const& optimalityType = boost::optional()) override; virtual std::unique_ptr computeNextProbabilities(storm::logic::NextFormula const& pathFormula, bool qualitative = false, boost::optional const& optimalityType = boost::optional()) override; virtual std::unique_ptr computeUntilProbabilities(storm::logic::UntilFormula const& pathFormula, bool qualitative = false, boost::optional const& optimalityType = boost::optional()) override; - virtual std::unique_ptr computeCumulativeRewards(storm::logic::CumulativeRewardFormula const& rewardPathFormula, boost::optional const& rewardModelName = boost::optional(), bool qualitative = false, boost::optional const& optimalityType = boost::optional()); - virtual std::unique_ptr computeInstantaneousRewards(storm::logic::InstantaneousRewardFormula const& rewardPathFormula, boost::optional const& rewardModelName = boost::optional(), bool qualitative = false, boost::optional const& optimalityType = boost::optional()); - virtual std::unique_ptr computeReachabilityRewards(storm::logic::ReachabilityRewardFormula const& rewardPathFormula, boost::optional const& rewardModelName = boost::optional(), bool qualitative = false, boost::optional const& optimalityType = boost::optional()); + virtual std::unique_ptr computeCumulativeRewards(storm::logic::CumulativeRewardFormula const& rewardPathFormula, boost::optional const& rewardModelName = boost::optional(), bool qualitative = false, boost::optional const& optimalityType = boost::optional()) override; + virtual std::unique_ptr computeInstantaneousRewards(storm::logic::InstantaneousRewardFormula const& rewardPathFormula, boost::optional const& rewardModelName = boost::optional(), bool qualitative = false, boost::optional const& optimalityType = boost::optional()) override; + virtual std::unique_ptr computeReachabilityRewards(storm::logic::ReachabilityRewardFormula const& rewardPathFormula, boost::optional const& rewardModelName = boost::optional(), bool qualitative = false, boost::optional const& optimalityType = boost::optional()) override; virtual std::unique_ptr computeLongRunAverage(storm::logic::StateFormula const& stateFormula, bool qualitative = false, boost::optional const& optimalityType = boost::optional()) override; private: - // The methods that perform the actual checking. - std::vector computeBoundedUntilProbabilitiesHelper(storm::storage::BitVector const& phiStates, storm::storage::BitVector const& psiStates, uint_fast64_t stepBound) const; - static std::vector computeNextProbabilitiesHelper(storm::storage::SparseMatrix const& transitionMatrix, storm::storage::BitVector const& nextStates, storm::utility::solver::LinearEquationSolverFactory const& linearEquationSolverFactory); - static std::vector computeUntilProbabilitiesHelper(storm::storage::SparseMatrix const& transitionMatrix, storm::storage::SparseMatrix const& backwardTransitions, storm::storage::BitVector const& phiStates, storm::storage::BitVector const& psiStates, bool qualitative, storm::utility::solver::LinearEquationSolverFactory const& linearEquationSolverFactory); - std::vector computeInstantaneousRewardsHelper(RewardModelType const& rewardModel, uint_fast64_t stepCount) const; - std::vector computeCumulativeRewardsHelper(RewardModelType const& rewardModel, uint_fast64_t stepBound) const; - static std::vector computeReachabilityRewardsHelper(RewardModelType const& rewardModel, storm::storage::SparseMatrix const& transitionMatrix, storm::storage::SparseMatrix const& backwardTransitions, storm::storage::BitVector const& targetStates, storm::utility::solver::LinearEquationSolverFactory const& linearEquationSolverFactory, bool qualitative); - static std::vector computeLongRunAverageHelper(storm::storage::SparseMatrix const& transitionMatrix, storm::storage::BitVector const& psiStates, bool qualitative, storm::utility::solver::LinearEquationSolverFactory const& linearEquationSolverFactory); - // An object that is used for retrieving linear equation solvers. std::unique_ptr> linearEquationSolverFactory; }; diff --git a/src/modelchecker/prctl/SparseMdpPrctlModelChecker.cpp b/src/modelchecker/prctl/SparseMdpPrctlModelChecker.cpp index b079ff63c..adceed6f8 100644 --- a/src/modelchecker/prctl/SparseMdpPrctlModelChecker.cpp +++ b/src/modelchecker/prctl/SparseMdpPrctlModelChecker.cpp @@ -36,40 +36,6 @@ namespace storm { return formula.isPctlStateFormula() || formula.isPctlPathFormula() || formula.isRewardPathFormula(); } - template - std::vector SparseMdpPrctlModelChecker::computeBoundedUntilProbabilitiesHelper(bool minimize, storm::storage::BitVector const& phiStates, storm::storage::BitVector const& psiStates, uint_fast64_t stepBound) const { - std::vector result(this->getModel().getNumberOfStates(), storm::utility::zero()); - - // Determine the states that have 0 probability of reaching the target states. - storm::storage::BitVector maybeStates; - if (minimize) { - maybeStates = storm::utility::graph::performProbGreater0A(this->getModel().getTransitionMatrix(), this->getModel().getTransitionMatrix().getRowGroupIndices(), this->getModel().getBackwardTransitions(), phiStates, psiStates, true, stepBound); - } else { - maybeStates = storm::utility::graph::performProbGreater0E(this->getModel().getTransitionMatrix(), this->getModel().getTransitionMatrix().getRowGroupIndices(), this->getModel().getBackwardTransitions(), phiStates, psiStates, true, stepBound); - } - maybeStates &= ~psiStates; - STORM_LOG_INFO("Found " << maybeStates.getNumberOfSetBits() << " 'maybe' states."); - - if (!maybeStates.empty()) { - // We can eliminate the rows and columns from the original transition probability matrix that have probability 0. - storm::storage::SparseMatrix submatrix = this->getModel().getTransitionMatrix().getSubmatrix(true, maybeStates, maybeStates, false); - std::vector b = this->getModel().getTransitionMatrix().getConstrainedRowGroupSumVector(maybeStates, psiStates); - - // Create the vector with which to multiply. - std::vector subresult(maybeStates.getNumberOfSetBits()); - - STORM_LOG_THROW(MinMaxLinearEquationSolverFactory != nullptr, storm::exceptions::InvalidStateException, "No valid equation solver available."); - std::unique_ptr> solver = MinMaxLinearEquationSolverFactory->create(submatrix); - solver->performMatrixVectorMultiplication(minimize, subresult, &b, stepBound); - - // Set the values of the resulting vector accordingly. - storm::utility::vector::setVectorValues(result, maybeStates, subresult); - } - storm::utility::vector::setVectorValues(result, psiStates, storm::utility::one()); - - return result; - } - template std::unique_ptr SparseMdpPrctlModelChecker::computeBoundedUntilProbabilities(storm::logic::BoundedUntilFormula const& pathFormula, bool qualitative, boost::optional const& optimalityType) { STORM_LOG_THROW(optimalityType, storm::exceptions::InvalidArgumentException, "Formula needs to specify whether minimal or maximal values are to be computed on nondeterministic model."); @@ -81,19 +47,6 @@ namespace storm { return result; } - template - std::vector SparseMdpPrctlModelChecker::computeNextProbabilitiesHelper(bool minimize, storm::storage::BitVector const& nextStates) { - // Create the vector with which to multiply and initialize it correctly. - std::vector result(this->getModel().getNumberOfStates()); - storm::utility::vector::setVectorValues(result, nextStates, storm::utility::one()); - - STORM_LOG_THROW(MinMaxLinearEquationSolverFactory != nullptr, storm::exceptions::InvalidStateException, "No valid equation solver available."); - std::unique_ptr> solver = MinMaxLinearEquationSolverFactory->create(this->getModel().getTransitionMatrix()); - solver->performMatrixVectorMultiplication(minimize, result); - - return result; - } - template std::unique_ptr SparseMdpPrctlModelChecker::computeNextProbabilities(storm::logic::NextFormula const& pathFormula, bool qualitative, boost::optional const& optimalityType) { STORM_LOG_THROW(optimalityType, storm::exceptions::InvalidArgumentException, "Formula needs to specify whether minimal or maximal values are to be computed on nondeterministic model."); @@ -102,68 +55,6 @@ namespace storm { return std::unique_ptr(new ExplicitQuantitativeCheckResult(this->computeNextProbabilitiesHelper(optimalityType.get() == storm::logic::OptimalityType::Minimize, subResult.getTruthValuesVector()))); } - template - std::vector SparseMdpPrctlModelChecker::computeUntilProbabilitiesHelper(bool minimize, storm::storage::BitVector const& phiStates, storm::storage::BitVector const& psiStates, bool qualitative) const { - return computeUntilProbabilitiesHelper(minimize, this->getModel().getTransitionMatrix(), this->getModel().getBackwardTransitions(), phiStates, psiStates, *MinMaxLinearEquationSolverFactory, qualitative); - } - - template - std::vector SparseMdpPrctlModelChecker::computeUntilProbabilitiesHelper(bool minimize, storm::storage::SparseMatrix const& transitionMatrix, storm::storage::SparseMatrix const& backwardTransitions, storm::storage::BitVector const& phiStates, storm::storage::BitVector const& psiStates, storm::utility::solver::MinMaxLinearEquationSolverFactory const& MinMaxLinearEquationSolverFactory, bool qualitative) { - size_t numberOfStates = phiStates.size(); - - // We need to identify the states which have to be taken out of the matrix, i.e. - // all states that have probability 0 and 1 of satisfying the until-formula. - std::pair statesWithProbability01; - if (minimize) { - statesWithProbability01 = storm::utility::graph::performProb01Min(transitionMatrix, transitionMatrix.getRowGroupIndices(), backwardTransitions, phiStates, psiStates); - } else { - statesWithProbability01 = storm::utility::graph::performProb01Max(transitionMatrix, transitionMatrix.getRowGroupIndices(), backwardTransitions, phiStates, psiStates); - } - storm::storage::BitVector statesWithProbability0 = std::move(statesWithProbability01.first); - storm::storage::BitVector statesWithProbability1 = std::move(statesWithProbability01.second); - storm::storage::BitVector maybeStates = ~(statesWithProbability0 | statesWithProbability1); - LOG4CPLUS_INFO(logger, "Found " << statesWithProbability0.getNumberOfSetBits() << " 'no' states."); - LOG4CPLUS_INFO(logger, "Found " << statesWithProbability1.getNumberOfSetBits() << " 'yes' states."); - LOG4CPLUS_INFO(logger, "Found " << maybeStates.getNumberOfSetBits() << " 'maybe' states."); - - // Create resulting vector. - std::vector result(numberOfStates); - - // Check whether we need to compute exact probabilities for some states. - if (qualitative) { - // Set the values for all maybe-states to 0.5 to indicate that their probability values are neither 0 nor 1. - storm::utility::vector::setVectorValues(result, maybeStates, ValueType(0.5)); - } else { - if (!maybeStates.empty()) { - // In this case we have have to compute the probabilities. - - // First, we can eliminate the rows and columns from the original transition probability matrix for states - // whose probabilities are already known. - storm::storage::SparseMatrix submatrix = transitionMatrix.getSubmatrix(true, maybeStates, maybeStates, false); - - // Prepare the right-hand side of the equation system. For entry i this corresponds to - // the accumulated probability of going from state i to some 'yes' state. - std::vector b = transitionMatrix.getConstrainedRowGroupSumVector(maybeStates, statesWithProbability1); - - // Create vector for results for maybe states. - std::vector x(maybeStates.getNumberOfSetBits()); - - // Solve the corresponding system of equations. - std::unique_ptr> solver = MinMaxLinearEquationSolverFactory.create(submatrix); - solver->solveEquationSystem(minimize, x, b); - - // Set values of resulting vector according to result. - storm::utility::vector::setVectorValues(result, maybeStates, x); - } - } - - // Set values of resulting vector that are known exactly. - storm::utility::vector::setVectorValues(result, statesWithProbability0, storm::utility::zero()); - storm::utility::vector::setVectorValues(result, statesWithProbability1, storm::utility::one()); - - return result; - } - template std::unique_ptr SparseMdpPrctlModelChecker::computeUntilProbabilities(storm::logic::UntilFormula const& pathFormula, bool qualitative, boost::optional const& optimalityType) { STORM_LOG_THROW(optimalityType, storm::exceptions::InvalidArgumentException, "Formula needs to specify whether minimal or maximal values are to be computed on nondeterministic model."); @@ -174,28 +65,6 @@ namespace storm { return std::unique_ptr(new ExplicitQuantitativeCheckResult(SparseMdpPrctlModelChecker::computeUntilProbabilitiesHelper(optimalityType.get() == storm::logic::OptimalityType::Minimize, this->getModel().getTransitionMatrix(), this->getModel().getBackwardTransitions(), leftResult.getTruthValuesVector(), rightResult.getTruthValuesVector(), *MinMaxLinearEquationSolverFactory, qualitative))); } - template - std::vector SparseMdpPrctlModelChecker::computeCumulativeRewardsHelper(RewardModelType const& rewardModel, bool minimize, uint_fast64_t stepBound) const { - // Only compute the result if the model has at least one reward this->getModel(). - STORM_LOG_THROW(!rewardModel.empty(), storm::exceptions::InvalidPropertyException, "Missing reward model for formula. Skipping formula."); - - // Compute the reward vector to add in each step based on the available reward models. - std::vector totalRewardVector = rewardModel.getTotalRewardVector(this->getModel().getTransitionMatrix()); - - // Initialize result to either the state rewards of the model or the null vector. - std::vector result; - if (rewardModel.hasStateRewards()) { - result = rewardModel.getStateRewardVector(); - } else { - result.resize(this->getModel().getNumberOfStates()); - } - - std::unique_ptr> solver = MinMaxLinearEquationSolverFactory->create(this->getModel().getTransitionMatrix()); - solver->performMatrixVectorMultiplication(minimize, result, &totalRewardVector, stepBound); - - return result; - } - template std::unique_ptr SparseMdpPrctlModelChecker::computeCumulativeRewards(storm::logic::CumulativeRewardFormula const& rewardPathFormula, boost::optional const& rewardModelName, bool qualitative, boost::optional const& optimalityType) { STORM_LOG_THROW(optimalityType, storm::exceptions::InvalidArgumentException, "Formula needs to specify whether minimal or maximal values are to be computed on nondeterministic model."); @@ -203,85 +72,13 @@ namespace storm { return std::unique_ptr(new ExplicitQuantitativeCheckResult(this->computeCumulativeRewardsHelper(rewardModelName ? this->getModel().getRewardModel(rewardModelName.get()) : this->getModel().getRewardModel(""), optimalityType.get() == storm::logic::OptimalityType::Minimize, rewardPathFormula.getDiscreteTimeBound()))); } - template - std::vector SparseMdpPrctlModelChecker::computeInstantaneousRewardsHelper(RewardModelType const& rewardModel, bool minimize, uint_fast64_t stepCount) const { - // Only compute the result if the model has a state-based reward this->getModel(). - STORM_LOG_THROW(rewardModel.hasStateRewards(), storm::exceptions::InvalidPropertyException, "Missing reward model for formula. Skipping formula."); - - // Initialize result to state rewards of the this->getModel(). - std::vector result(rewardModel.getStateRewardVector()); - - STORM_LOG_THROW(MinMaxLinearEquationSolverFactory != nullptr, storm::exceptions::InvalidStateException, "No valid linear equation solver available."); - std::unique_ptr> solver = MinMaxLinearEquationSolverFactory->create(this->getModel().getTransitionMatrix()); - solver->performMatrixVectorMultiplication(minimize, result, nullptr, stepCount); - - return result; - } - template std::unique_ptr SparseMdpPrctlModelChecker::computeInstantaneousRewards(storm::logic::InstantaneousRewardFormula const& rewardPathFormula, boost::optional const& rewardModelName, bool qualitative, boost::optional const& optimalityType) { STORM_LOG_THROW(optimalityType, storm::exceptions::InvalidArgumentException, "Formula needs to specify whether minimal or maximal values are to be computed on nondeterministic model."); STORM_LOG_THROW(rewardPathFormula.hasDiscreteTimeBound(), storm::exceptions::InvalidArgumentException, "Formula needs to have a discrete time bound."); return std::unique_ptr(new ExplicitQuantitativeCheckResult(this->computeInstantaneousRewardsHelper(rewardModelName ? this->getModel().getRewardModel(rewardModelName.get()) : this->getModel().getRewardModel(""), optimalityType.get() == storm::logic::OptimalityType::Minimize, rewardPathFormula.getDiscreteTimeBound()))); } - - template - std::vector SparseMdpPrctlModelChecker::computeReachabilityRewardsHelper(RewardModelType const& rewardModel, bool minimize, storm::storage::SparseMatrix const& transitionMatrix, storm::storage::SparseMatrix const& backwardTransitions, storm::storage::BitVector const& targetStates, storm::utility::solver::MinMaxLinearEquationSolverFactory const& MinMaxLinearEquationSolverFactory, bool qualitative) const { - // Only compute the result if the model has at least one reward this->getModel(). - STORM_LOG_THROW(!rewardModel.empty(), storm::exceptions::InvalidPropertyException, "Missing reward model for formula. Skipping formula."); - - // Determine which states have a reward of infinity by definition. - storm::storage::BitVector infinityStates; - storm::storage::BitVector trueStates(transitionMatrix.getRowCount(), true); - if (minimize) { - infinityStates = std::move(storm::utility::graph::performProb1A(transitionMatrix, transitionMatrix.getRowGroupIndices(), backwardTransitions, trueStates, targetStates)); - } else { - infinityStates = std::move(storm::utility::graph::performProb1E(transitionMatrix, transitionMatrix.getRowGroupIndices(), backwardTransitions, trueStates, targetStates)); - } - infinityStates.complement(); - storm::storage::BitVector maybeStates = ~targetStates & ~infinityStates; - LOG4CPLUS_INFO(logger, "Found " << infinityStates.getNumberOfSetBits() << " 'infinity' states."); - LOG4CPLUS_INFO(logger, "Found " << targetStates.getNumberOfSetBits() << " 'target' states."); - LOG4CPLUS_INFO(logger, "Found " << maybeStates.getNumberOfSetBits() << " 'maybe' states."); - - // Create resulting vector. - std::vector result(transitionMatrix.getRowCount()); - - // Check whether we need to compute exact rewards for some states. - if (qualitative) { - LOG4CPLUS_INFO(logger, "The rewards for the initial states were determined in a preprocessing step. No exact rewards were computed."); - // Set the values for all maybe-states to 1 to indicate that their reward values - // are neither 0 nor infinity. - storm::utility::vector::setVectorValues(result, maybeStates, storm::utility::one()); - } else { - // In this case we have to compute the reward values for the remaining states. - - // We can eliminate the rows and columns from the original transition probability matrix for states - // whose reward values are already known. - storm::storage::SparseMatrix submatrix = transitionMatrix.getSubmatrix(true, maybeStates, maybeStates, false); - // Prepare the right-hand side of the equation system. - std::vector b = rewardModel.getTotalRewardVector(submatrix.getRowCount(), transitionMatrix, maybeStates); - - // Create vector for results for maybe states. - std::vector x(maybeStates.getNumberOfSetBits()); - - // Solve the corresponding system of equations. - std::unique_ptr> solver = MinMaxLinearEquationSolverFactory.create(submatrix); - solver->solveEquationSystem(minimize, x, b); - - - // Set values of resulting vector according to result. - storm::utility::vector::setVectorValues(result, maybeStates, x); - } - - // Set values of resulting vector that are known exactly. - storm::utility::vector::setVectorValues(result, targetStates, storm::utility::zero()); - storm::utility::vector::setVectorValues(result, infinityStates, storm::utility::infinity()); - - return result; - } - template std::unique_ptr SparseMdpPrctlModelChecker::computeReachabilityRewards(storm::logic::ReachabilityRewardFormula const& rewardPathFormula, boost::optional const& rewardModelName, bool qualitative, boost::optional const& optimalityType) { STORM_LOG_THROW(optimalityType, storm::exceptions::InvalidArgumentException, "Formula needs to specify whether minimal or maximal values are to be computed on nondeterministic model."); @@ -289,158 +86,6 @@ namespace storm { ExplicitQualitativeCheckResult const& subResult = subResultPointer->asExplicitQualitativeCheckResult(); return std::unique_ptr(new ExplicitQuantitativeCheckResult(this->computeReachabilityRewardsHelper(rewardModelName ? this->getModel().getRewardModel(rewardModelName.get()) : this->getModel().getRewardModel(""), optimalityType.get() == storm::logic::OptimalityType::Minimize, this->getModel().getTransitionMatrix(), this->getModel().getBackwardTransitions(), subResult.getTruthValuesVector(), *this->MinMaxLinearEquationSolverFactory, qualitative))); } - - - template - std::vector SparseMdpPrctlModelChecker::computeLongRunAverageHelper(bool minimize, storm::storage::BitVector const& psiStates, bool qualitative) const { - // If there are no goal states, we avoid the computation and directly return zero. - auto numOfStates = this->getModel().getNumberOfStates(); - if (psiStates.empty()) { - return std::vector(numOfStates, storm::utility::zero()); - } - - // Likewise, if all bits are set, we can avoid the computation and set. - if ((~psiStates).empty()) { - return std::vector(numOfStates, storm::utility::one()); - } - - // Start by decomposing the MDP into its MECs. - storm::storage::MaximalEndComponentDecomposition mecDecomposition(this->getModel()); - - // Get some data members for convenience. - typename storm::storage::SparseMatrix const& transitionMatrix = this->getModel().getTransitionMatrix(); - std::vector const& nondeterministicChoiceIndices = this->getModel().getNondeterministicChoiceIndices(); - ValueType zero = storm::utility::zero(); - - //first calculate LRA for the Maximal End Components. - storm::storage::BitVector statesInMecs(numOfStates); - std::vector stateToMecIndexMap(transitionMatrix.getColumnCount()); - std::vector lraValuesForEndComponents(mecDecomposition.size(), zero); - - for (uint_fast64_t currentMecIndex = 0; currentMecIndex < mecDecomposition.size(); ++currentMecIndex) { - storm::storage::MaximalEndComponent const& mec = mecDecomposition[currentMecIndex]; - - lraValuesForEndComponents[currentMecIndex] = computeLraForMaximalEndComponent(minimize, transitionMatrix, psiStates, mec); - - // Gather information for later use. - for (auto const& stateChoicesPair : mec) { - statesInMecs.set(stateChoicesPair.first); - stateToMecIndexMap[stateChoicesPair.first] = currentMecIndex; - } - } - - // For fast transition rewriting, we build some auxiliary data structures. - storm::storage::BitVector statesNotContainedInAnyMec = ~statesInMecs; - uint_fast64_t firstAuxiliaryStateIndex = statesNotContainedInAnyMec.getNumberOfSetBits(); - uint_fast64_t lastStateNotInMecs = 0; - uint_fast64_t numberOfStatesNotInMecs = 0; - std::vector statesNotInMecsBeforeIndex; - statesNotInMecsBeforeIndex.reserve(this->getModel().getNumberOfStates()); - for (auto state : statesNotContainedInAnyMec) { - while (lastStateNotInMecs <= state) { - statesNotInMecsBeforeIndex.push_back(numberOfStatesNotInMecs); - ++lastStateNotInMecs; - } - ++numberOfStatesNotInMecs; - } - - // Finally, we are ready to create the SSP matrix and right-hand side of the SSP. - std::vector b; - typename storm::storage::SparseMatrixBuilder sspMatrixBuilder(0, 0, 0, false, true, numberOfStatesNotInMecs + mecDecomposition.size()); - - // If the source state is not contained in any MEC, we copy its choices (and perform the necessary modifications). - uint_fast64_t currentChoice = 0; - for (auto state : statesNotContainedInAnyMec) { - sspMatrixBuilder.newRowGroup(currentChoice); - - for (uint_fast64_t choice = nondeterministicChoiceIndices[state]; choice < nondeterministicChoiceIndices[state + 1]; ++choice, ++currentChoice) { - std::vector auxiliaryStateToProbabilityMap(mecDecomposition.size()); - b.push_back(storm::utility::zero()); - - for (auto element : transitionMatrix.getRow(choice)) { - if (statesNotContainedInAnyMec.get(element.getColumn())) { - // If the target state is not contained in an MEC, we can copy over the entry. - sspMatrixBuilder.addNextValue(currentChoice, statesNotInMecsBeforeIndex[element.getColumn()], element.getValue()); - } else { - // If the target state is contained in MEC i, we need to add the probability to the corresponding field in the vector - // so that we are able to write the cumulative probability to the MEC into the matrix. - auxiliaryStateToProbabilityMap[stateToMecIndexMap[element.getColumn()]] += element.getValue(); - } - } - - // Now insert all (cumulative) probability values that target an MEC. - for (uint_fast64_t mecIndex = 0; mecIndex < auxiliaryStateToProbabilityMap.size(); ++mecIndex) { - if (auxiliaryStateToProbabilityMap[mecIndex] != 0) { - sspMatrixBuilder.addNextValue(currentChoice, firstAuxiliaryStateIndex + mecIndex, auxiliaryStateToProbabilityMap[mecIndex]); - } - } - } - } - - // Now we are ready to construct the choices for the auxiliary states. - for (uint_fast64_t mecIndex = 0; mecIndex < mecDecomposition.size(); ++mecIndex) { - storm::storage::MaximalEndComponent const& mec = mecDecomposition[mecIndex]; - sspMatrixBuilder.newRowGroup(currentChoice); - - for (auto const& stateChoicesPair : mec) { - uint_fast64_t state = stateChoicesPair.first; - boost::container::flat_set const& choicesInMec = stateChoicesPair.second; - - for (uint_fast64_t choice = nondeterministicChoiceIndices[state]; choice < nondeterministicChoiceIndices[state + 1]; ++choice) { - std::vector auxiliaryStateToProbabilityMap(mecDecomposition.size()); - - // If the choice is not contained in the MEC itself, we have to add a similar distribution to the auxiliary state. - if (choicesInMec.find(choice) == choicesInMec.end()) { - b.push_back(storm::utility::zero()); - - for (auto element : transitionMatrix.getRow(choice)) { - if (statesNotContainedInAnyMec.get(element.getColumn())) { - // If the target state is not contained in an MEC, we can copy over the entry. - sspMatrixBuilder.addNextValue(currentChoice, statesNotInMecsBeforeIndex[element.getColumn()], element.getValue()); - } else { - // If the target state is contained in MEC i, we need to add the probability to the corresponding field in the vector - // so that we are able to write the cumulative probability to the MEC into the matrix. - auxiliaryStateToProbabilityMap[stateToMecIndexMap[element.getColumn()]] += element.getValue(); - } - } - - // Now insert all (cumulative) probability values that target an MEC. - for (uint_fast64_t targetMecIndex = 0; targetMecIndex < auxiliaryStateToProbabilityMap.size(); ++targetMecIndex) { - if (auxiliaryStateToProbabilityMap[targetMecIndex] != 0) { - sspMatrixBuilder.addNextValue(currentChoice, firstAuxiliaryStateIndex + targetMecIndex, auxiliaryStateToProbabilityMap[targetMecIndex]); - } - } - - ++currentChoice; - } - } - } - - // For each auxiliary state, there is the option to achieve the reward value of the LRA associated with the MEC. - ++currentChoice; - b.push_back(lraValuesForEndComponents[mecIndex]); - } - - // Finalize the matrix and solve the corresponding system of equations. - storm::storage::SparseMatrix sspMatrix = sspMatrixBuilder.build(currentChoice); - - std::vector sspResult(numberOfStatesNotInMecs + mecDecomposition.size()); - std::unique_ptr> solver = MinMaxLinearEquationSolverFactory->create(sspMatrix); - solver->solveEquationSystem(minimize, sspResult, b); - - // Prepare result vector. - std::vector result(this->getModel().getNumberOfStates()); - - // Set the values for states not contained in MECs. - storm::utility::vector::setVectorValues(result, statesNotContainedInAnyMec, sspResult); - - // Set the values for all states in MECs. - for (auto state : statesInMecs) { - result[state] = sspResult[firstAuxiliaryStateIndex + stateToMecIndexMap[state]]; - } - - return result; - } template std::unique_ptr SparseMdpPrctlModelChecker::computeLongRunAverage(storm::logic::StateFormula const& stateFormula, bool qualitative, boost::optional const& optimalityType) { @@ -451,50 +96,6 @@ namespace storm { return std::unique_ptr(new ExplicitQuantitativeCheckResult(this->computeLongRunAverageHelper(optimalityType.get() == storm::logic::OptimalityType::Minimize, subResult.getTruthValuesVector(), qualitative))); } - - template - typename SparseMdpModelType::ValueType SparseMdpPrctlModelChecker::computeLraForMaximalEndComponent(bool minimize, storm::storage::SparseMatrix const& transitionMatrix, storm::storage::BitVector const& psiStates, storm::storage::MaximalEndComponent const& mec) { - std::shared_ptr solver = storm::utility::solver::getLpSolver("LRA for MEC"); - solver->setModelSense(minimize ? storm::solver::LpSolver::ModelSense::Maximize : storm::solver::LpSolver::ModelSense::Minimize); - - //// First, we need to create the variables for the problem. - std::map stateToVariableMap; - for (auto const& stateChoicesPair : mec) { - std::string variableName = "h" + std::to_string(stateChoicesPair.first); - stateToVariableMap[stateChoicesPair.first] = solver->addUnboundedContinuousVariable(variableName); - } - storm::expressions::Variable lambda = solver->addUnboundedContinuousVariable("L", 1); - solver->update(); - - // Now we encode the problem as constraints. - for (auto const& stateChoicesPair : mec) { - uint_fast64_t state = stateChoicesPair.first; - - // Now, based on the type of the state, create a suitable constraint. - for (auto choice : stateChoicesPair.second) { - storm::expressions::Expression constraint = -lambda; - ValueType r = 0; - - for (auto element : transitionMatrix.getRow(choice)) { - constraint = constraint + stateToVariableMap.at(element.getColumn()) * solver->getConstant(element.getValue()); - if (psiStates.get(element.getColumn())) { - r += element.getValue(); - } - } - constraint = solver->getConstant(r) + constraint; - - if (minimize) { - constraint = stateToVariableMap.at(state) <= constraint; - } else { - constraint = stateToVariableMap.at(state) >= constraint; - } - solver->addConstraint("state" + std::to_string(state) + "," + std::to_string(choice), constraint); - } - } - - solver->optimize(); - return solver->getContinuousValue(lambda); - } template class SparseMdpPrctlModelChecker>; } diff --git a/src/modelchecker/prctl/SparseMdpPrctlModelChecker.h b/src/modelchecker/prctl/SparseMdpPrctlModelChecker.h index 4223c7e62..2bd9c091a 100644 --- a/src/modelchecker/prctl/SparseMdpPrctlModelChecker.h +++ b/src/modelchecker/prctl/SparseMdpPrctlModelChecker.h @@ -5,8 +5,6 @@ #include "src/models/sparse/Mdp.h" #include "src/utility/solver.h" #include "src/solver/MinMaxLinearEquationSolver.h" -#include "src/storage/MaximalEndComponent.h" - namespace storm { namespace counterexamples { @@ -18,18 +16,12 @@ namespace storm { } namespace modelchecker { - - // Forward-declare other model checkers to make them friend classes. - template - class SparseMarkovAutomatonCslModelChecker; - template class SparseMdpPrctlModelChecker : public SparsePropositionalModelChecker { public: typedef typename SparseMdpModelType::ValueType ValueType; typedef typename SparseMdpModelType::RewardModelType RewardModelType; - friend class SparseMarkovAutomatonCslModelChecker; friend class storm::counterexamples::SMTMinimalCommandSetGenerator; friend class storm::counterexamples::MILPMinimalLabelSetGenerator; @@ -47,20 +39,8 @@ namespace storm { virtual std::unique_ptr computeLongRunAverage(storm::logic::StateFormula const& stateFormula, bool qualitative = false, boost::optional const& optimalityType = boost::optional()) override; private: - // The methods that perform the actual checking. - std::vector computeBoundedUntilProbabilitiesHelper(bool minimize, storm::storage::BitVector const& phiStates, storm::storage::BitVector const& psiStates, uint_fast64_t stepBound) const; - std::vector computeNextProbabilitiesHelper(bool minimize, storm::storage::BitVector const& nextStates); - std::vector computeUntilProbabilitiesHelper(bool minimize, storm::storage::BitVector const& phiStates, storm::storage::BitVector const& psiStates, bool qualitative) const; - static std::vector computeUntilProbabilitiesHelper(bool minimize, storm::storage::SparseMatrix const& transitionMatrix, storm::storage::SparseMatrix const& backwardTransitions, storm::storage::BitVector const& phiStates, storm::storage::BitVector const& psiStates, storm::utility::solver::MinMaxLinearEquationSolverFactory const& MinMaxLinearEquationSolverFactory, bool qualitative); - std::vector computeInstantaneousRewardsHelper(RewardModelType const& rewardModel, bool minimize, uint_fast64_t stepCount) const; - std::vector computeCumulativeRewardsHelper(RewardModelType const& rewardModel, bool minimize, uint_fast64_t stepBound) const; - std::vector computeReachabilityRewardsHelper(RewardModelType const& rewardModel, bool minimize, storm::storage::SparseMatrix const& transitionMatrix, storm::storage::SparseMatrix const& backwardTransitions, storm::storage::BitVector const& targetStates, storm::utility::solver::MinMaxLinearEquationSolverFactory const& MinMaxLinearEquationSolverFactory, bool qualitative) const; - std::vector computeLongRunAverageHelper(bool minimize, storm::storage::BitVector const& psiStates, bool qualitative) const; - - static ValueType computeLraForMaximalEndComponent(bool minimize, storm::storage::SparseMatrix const& transitionMatrix, storm::storage::BitVector const& goalStates, storm::storage::MaximalEndComponent const& mec); - // An object that is used for retrieving solvers for systems of linear equations that are the result of nondeterministic choices. - std::unique_ptr> MinMaxLinearEquationSolverFactory; + std::unique_ptr> minMaxLinearEquationSolverFactory; }; } // namespace modelchecker } // namespace storm diff --git a/src/modelchecker/prctl/SymbolicDtmcPrctlModelChecker.cpp b/src/modelchecker/prctl/SymbolicDtmcPrctlModelChecker.cpp index 022eca0af..891061edc 100644 --- a/src/modelchecker/prctl/SymbolicDtmcPrctlModelChecker.cpp +++ b/src/modelchecker/prctl/SymbolicDtmcPrctlModelChecker.cpp @@ -1,13 +1,13 @@ #include "src/modelchecker/prctl/SymbolicDtmcPrctlModelChecker.h" -#include "src/storage/dd/CuddOdd.h" +#include "src/modelchecker/prctl/helper/SymbolicDtmcPrctlHelper.h" + +#include "src/storage/dd/Add.h" #include "src/utility/macros.h" -#include "src/utility/graph.h" #include "src/modelchecker/results/SymbolicQualitativeCheckResult.h" #include "src/modelchecker/results/SymbolicQuantitativeCheckResult.h" -#include "src/modelchecker/results/HybridQuantitativeCheckResult.h" #include "src/exceptions/InvalidStateException.h" #include "src/exceptions/InvalidPropertyException.h" @@ -29,80 +29,24 @@ namespace storm { return formula.isPctlStateFormula() || formula.isPctlPathFormula() || formula.isRewardPathFormula(); } - template - std::unique_ptr SymbolicDtmcPrctlModelChecker::computeUntilProbabilitiesHelper(storm::models::symbolic::Model const& model, storm::dd::Add const& transitionMatrix, storm::dd::Bdd const& phiStates, storm::dd::Bdd const& psiStates, bool qualitative, storm::utility::solver::SymbolicLinearEquationSolverFactory const& linearEquationSolverFactory) { - // We need to identify the states which have to be taken out of the matrix, i.e. all states that have - // probability 0 and 1 of satisfying the until-formula. - std::pair, storm::dd::Bdd> statesWithProbability01 = storm::utility::graph::performProb01(model, transitionMatrix, phiStates, psiStates); - storm::dd::Bdd maybeStates = !statesWithProbability01.first && !statesWithProbability01.second && model.getReachableStates(); - - // Perform some logging. - STORM_LOG_INFO("Found " << statesWithProbability01.first.getNonZeroCount() << " 'no' states."); - STORM_LOG_INFO("Found " << statesWithProbability01.second.getNonZeroCount() << " 'yes' states."); - STORM_LOG_INFO("Found " << maybeStates.getNonZeroCount() << " 'maybe' states."); - - // Check whether we need to compute exact probabilities for some states. - if (qualitative) { - // Set the values for all maybe-states to 0.5 to indicate that their probability values are neither 0 nor 1. - return std::unique_ptr(new storm::modelchecker::SymbolicQuantitativeCheckResult(model.getReachableStates(), statesWithProbability01.second.toAdd() + maybeStates.toAdd() * model.getManager().getConstant(0.5))); - } else { - // If there are maybe states, we need to solve an equation system. - if (!maybeStates.isZero()) { - // Create the ODD for the translation between symbolic and explicit storage. - storm::dd::Odd odd(maybeStates); - - // Create the matrix and the vector for the equation system. - storm::dd::Add maybeStatesAdd = maybeStates.toAdd(); - - // Start by cutting away all rows that do not belong to maybe states. Note that this leaves columns targeting - // non-maybe states in the matrix. - storm::dd::Add submatrix = transitionMatrix * maybeStatesAdd; - - // Then compute the vector that contains the one-step probabilities to a state with probability 1 for all - // maybe states. - storm::dd::Add prob1StatesAsColumn = statesWithProbability01.second.toAdd(); - prob1StatesAsColumn = prob1StatesAsColumn.swapVariables(model.getRowColumnMetaVariablePairs()); - storm::dd::Add subvector = submatrix * prob1StatesAsColumn; - subvector = subvector.sumAbstract(model.getColumnVariables()); - - // Finally cut away all columns targeting non-maybe states and convert the matrix into the matrix needed - // for solving the equation system (i.e. compute (I-A)). - submatrix *= maybeStatesAdd.swapVariables(model.getRowColumnMetaVariablePairs()); - submatrix = (model.getRowColumnIdentity() * maybeStatesAdd) - submatrix; - - // Solve the equation system. - std::unique_ptr> solver = linearEquationSolverFactory.create(submatrix, maybeStates, model.getRowVariables(), model.getColumnVariables(), model.getRowColumnMetaVariablePairs()); - storm::dd::Add result = solver->solveEquationSystem(model.getManager().getConstant(0.5) * maybeStatesAdd, subvector); - - return std::unique_ptr(new storm::modelchecker::SymbolicQuantitativeCheckResult(model.getReachableStates(), statesWithProbability01.second.toAdd() + result)); - } else { - return std::unique_ptr(new storm::modelchecker::SymbolicQuantitativeCheckResult(model.getReachableStates(), statesWithProbability01.second.toAdd())); - } - } - } - template std::unique_ptr SymbolicDtmcPrctlModelChecker::computeUntilProbabilities(storm::logic::UntilFormula const& pathFormula, bool qualitative, boost::optional const& optimalityType) { std::unique_ptr leftResultPointer = this->check(pathFormula.getLeftSubformula()); std::unique_ptr rightResultPointer = this->check(pathFormula.getRightSubformula()); SymbolicQualitativeCheckResult const& leftResult = leftResultPointer->asSymbolicQualitativeCheckResult(); SymbolicQualitativeCheckResult const& rightResult = rightResultPointer->asSymbolicQualitativeCheckResult(); - return this->computeUntilProbabilitiesHelper(this->getModel(), this->getModel().getTransitionMatrix(), leftResult.getTruthValuesVector(), rightResult.getTruthValuesVector(), qualitative, *this->linearEquationSolverFactory); + storm::dd::Add numericResult = storm::modelchecker::helper::SymbolicDtmcPrctlHelper::computeUntilProbabilities(this->getModel(), this->getModel().getTransitionMatrix(), leftResult.getTruthValuesVector(), rightResult.getTruthValuesVector(), qualitative, *this->linearEquationSolverFactory); + return std::unique_ptr>(new SymbolicQuantitativeCheckResult(this->getModel().getReachableStates(), numericResult)); } template std::unique_ptr SymbolicDtmcPrctlModelChecker::computeNextProbabilities(storm::logic::NextFormula const& pathFormula, bool qualitative, boost::optional const& optimalityType) { std::unique_ptr subResultPointer = this->check(pathFormula.getSubformula()); SymbolicQualitativeCheckResult const& subResult = subResultPointer->asSymbolicQualitativeCheckResult(); - return std::unique_ptr(new SymbolicQuantitativeCheckResult(this->getModel().getReachableStates(), this->computeNextProbabilitiesHelper(this->getModel(), this->getModel().getTransitionMatrix(), subResult.getTruthValuesVector()))); - } - - template - storm::dd::Add SymbolicDtmcPrctlModelChecker::computeNextProbabilitiesHelper(storm::models::symbolic::Model const& model, storm::dd::Add const& transitionMatrix, storm::dd::Bdd const& nextStates) { - storm::dd::Add result = transitionMatrix * nextStates.swapVariables(model.getRowColumnMetaVariablePairs()).toAdd(); - return result.sumAbstract(model.getColumnVariables()); + storm::dd::Add numericResult = storm::modelchecker::helper::SymbolicDtmcPrctlHelper::computeNextProbabilities(this->getModel(), this->getModel().getTransitionMatrix(), subResult.getTruthValuesVector()); + return std::unique_ptr>(new SymbolicQuantitativeCheckResult(this->getModel().getReachableStates(), numericResult)); } - + template std::unique_ptr SymbolicDtmcPrctlModelChecker::computeBoundedUntilProbabilities(storm::logic::BoundedUntilFormula const& pathFormula, bool qualitative, boost::optional const& optimalityType) { STORM_LOG_THROW(pathFormula.hasDiscreteTimeBound(), storm::exceptions::InvalidArgumentException, "Formula needs to have a discrete time bound."); @@ -110,147 +54,30 @@ namespace storm { std::unique_ptr rightResultPointer = this->check(pathFormula.getRightSubformula()); SymbolicQualitativeCheckResult const& leftResult = leftResultPointer->asSymbolicQualitativeCheckResult(); SymbolicQualitativeCheckResult const& rightResult = rightResultPointer->asSymbolicQualitativeCheckResult(); - return this->computeBoundedUntilProbabilitiesHelper(this->getModel(), this->getModel().getTransitionMatrix(), leftResult.getTruthValuesVector(), rightResult.getTruthValuesVector(), pathFormula.getDiscreteTimeBound(), *this->linearEquationSolverFactory); + storm::dd::Add numericResult = storm::modelchecker::helper::SymbolicDtmcPrctlHelper::computeBoundedUntilProbabilities(this->getModel(), this->getModel().getTransitionMatrix(), leftResult.getTruthValuesVector(), rightResult.getTruthValuesVector(), pathFormula.getDiscreteTimeBound(), *this->linearEquationSolverFactory); + return std::unique_ptr>(new SymbolicQuantitativeCheckResult(this->getModel().getReachableStates(), numericResult)); } template - std::unique_ptr SymbolicDtmcPrctlModelChecker::computeBoundedUntilProbabilitiesHelper(storm::models::symbolic::Model const& model, storm::dd::Add const& transitionMatrix, storm::dd::Bdd const& phiStates, storm::dd::Bdd const& psiStates, uint_fast64_t stepBound, storm::utility::solver::SymbolicLinearEquationSolverFactory const& linearEquationSolverFactory) { - // We need to identify the states which have to be taken out of the matrix, i.e. all states that have - // probability 0 or 1 of satisfying the until-formula. - storm::dd::Bdd statesWithProbabilityGreater0 = storm::utility::graph::performProbGreater0(model, transitionMatrix.notZero(), phiStates, psiStates, stepBound); - storm::dd::Bdd maybeStates = statesWithProbabilityGreater0 && !psiStates && model.getReachableStates(); - - // If there are maybe states, we need to perform matrix-vector multiplications. - if (!maybeStates.isZero()) { - // Create the ODD for the translation between symbolic and explicit storage. - storm::dd::Odd odd(maybeStates); - - // Create the matrix and the vector for the equation system. - storm::dd::Add maybeStatesAdd = maybeStates.toAdd(); - - // Start by cutting away all rows that do not belong to maybe states. Note that this leaves columns targeting - // non-maybe states in the matrix. - storm::dd::Add submatrix = transitionMatrix * maybeStatesAdd; - - // Then compute the vector that contains the one-step probabilities to a state with probability 1 for all - // maybe states. - storm::dd::Add prob1StatesAsColumn = psiStates.toAdd().swapVariables(model.getRowColumnMetaVariablePairs()); - storm::dd::Add subvector = (submatrix * prob1StatesAsColumn).sumAbstract(model.getColumnVariables()); - - // Finally cut away all columns targeting non-maybe states. - submatrix *= maybeStatesAdd.swapVariables(model.getRowColumnMetaVariablePairs()); - - // Perform the matrix-vector multiplication. - std::unique_ptr> solver = linearEquationSolverFactory.create(submatrix, maybeStates, model.getRowVariables(), model.getColumnVariables(), model.getRowColumnMetaVariablePairs()); - storm::dd::Add result = solver->performMatrixVectorMultiplication(model.getManager().getAddZero(), &subvector, stepBound); - - return std::unique_ptr(new storm::modelchecker::SymbolicQuantitativeCheckResult(model.getReachableStates(), psiStates.toAdd() + result)); - } else { - return std::unique_ptr(new storm::modelchecker::SymbolicQuantitativeCheckResult(model.getReachableStates(), psiStates.toAdd())); - } - } - - template - std::unique_ptr SymbolicDtmcPrctlModelChecker::computeCumulativeRewards(storm::logic::CumulativeRewardFormula const& rewardPathFormula, bool qualitative, boost::optional const& optimalityType) { + std::unique_ptr SymbolicDtmcPrctlModelChecker::computeCumulativeRewards(storm::logic::CumulativeRewardFormula const& rewardPathFormula, boost::optional const& rewardModelName, bool qualitative, boost::optional const& optimalityType) { STORM_LOG_THROW(rewardPathFormula.hasDiscreteTimeBound(), storm::exceptions::InvalidArgumentException, "Formula needs to have a discrete time bound."); - return this->computeCumulativeRewardsHelper(this->getModel(), this->getModel().getTransitionMatrix(), rewardPathFormula.getDiscreteTimeBound(), *this->linearEquationSolverFactory); - } - - template - std::unique_ptr SymbolicDtmcPrctlModelChecker::computeCumulativeRewardsHelper(storm::models::symbolic::Model const& model, storm::dd::Add const& transitionMatrix, uint_fast64_t stepBound, storm::utility::solver::SymbolicLinearEquationSolverFactory const& linearEquationSolverFactory) { - // Only compute the result if the model has at least one reward this->getModel(). - STORM_LOG_THROW(model.hasStateRewards() || model.hasTransitionRewards(), storm::exceptions::InvalidPropertyException, "Missing reward model for formula. Skipping formula."); - - // Compute the reward vector to add in each step based on the available reward models. - storm::dd::Add totalRewardVector = model.hasStateRewards() ? model.getStateRewardVector() : model.getManager().getAddZero(); - if (model.hasTransitionRewards()) { - totalRewardVector += (transitionMatrix * model.getTransitionRewardMatrix()).sumAbstract(model.getColumnVariables()); - } - - // Perform the matrix-vector multiplication. - std::unique_ptr> solver = linearEquationSolverFactory.create(transitionMatrix, model.getReachableStates(), model.getRowVariables(), model.getColumnVariables(), model.getRowColumnMetaVariablePairs()); - storm::dd::Add result = solver->performMatrixVectorMultiplication(model.getManager().getAddZero(), &totalRewardVector, stepBound); - - return std::unique_ptr(new storm::modelchecker::SymbolicQuantitativeCheckResult(model.getReachableStates(), result)); + storm::dd::Add numericResult = storm::modelchecker::helper::SymbolicDtmcPrctlHelper::computeCumulativeRewards(this->getModel(), this->getModel().getTransitionMatrix(), rewardPathFormula.getDiscreteTimeBound(), *this->linearEquationSolverFactory); + return std::unique_ptr>(new SymbolicQuantitativeCheckResult(this->getModel().getReachableStates(), numericResult)); } template - std::unique_ptr SymbolicDtmcPrctlModelChecker::computeInstantaneousRewards(storm::logic::InstantaneousRewardFormula const& rewardPathFormula, bool qualitative, boost::optional const& optimalityType) { + std::unique_ptr SymbolicDtmcPrctlModelChecker::computeInstantaneousRewards(storm::logic::InstantaneousRewardFormula const& rewardPathFormula, boost::optional const& rewardModelName, bool qualitative, boost::optional const& optimalityType) { STORM_LOG_THROW(rewardPathFormula.hasDiscreteTimeBound(), storm::exceptions::InvalidArgumentException, "Formula needs to have a discrete time bound."); - return this->computeInstantaneousRewardsHelper(this->getModel(), this->getModel().getTransitionMatrix(), rewardPathFormula.getDiscreteTimeBound(), *this->linearEquationSolverFactory); + storm::dd::Add numericResult = storm::modelchecker::helper::SymbolicDtmcPrctlHelper::computeInstantaneousRewards(this->getModel(), this->getModel().getTransitionMatrix(), rewardPathFormula.getDiscreteTimeBound(), *this->linearEquationSolverFactory); + return std::unique_ptr>(new SymbolicQuantitativeCheckResult(this->getModel().getReachableStates(), numericResult)); } template - std::unique_ptr SymbolicDtmcPrctlModelChecker::computeInstantaneousRewardsHelper(storm::models::symbolic::Model const& model, storm::dd::Add const& transitionMatrix, uint_fast64_t stepBound, storm::utility::solver::SymbolicLinearEquationSolverFactory const& linearEquationSolverFactory) { - // Only compute the result if the model has at least one reward this->getModel(). - STORM_LOG_THROW(model.hasStateRewards(), storm::exceptions::InvalidPropertyException, "Missing reward model for formula. Skipping formula."); - - // Perform the matrix-vector multiplication. - std::unique_ptr> solver = linearEquationSolverFactory.create(transitionMatrix, model.getReachableStates(), model.getRowVariables(), model.getColumnVariables(), model.getRowColumnMetaVariablePairs()); - storm::dd::Add result = solver->performMatrixVectorMultiplication(model.getStateRewardVector(), nullptr, stepBound); - - return std::unique_ptr(new storm::modelchecker::SymbolicQuantitativeCheckResult(model.getReachableStates(), result)); - } - - template - std::unique_ptr SymbolicDtmcPrctlModelChecker::computeReachabilityRewards(storm::logic::ReachabilityRewardFormula const& rewardPathFormula, bool qualitative, boost::optional const& optimalityType) { + std::unique_ptr SymbolicDtmcPrctlModelChecker::computeReachabilityRewards(storm::logic::ReachabilityRewardFormula const& rewardPathFormula, boost::optional const& rewardModelName, bool qualitative, boost::optional const& optimalityType) { std::unique_ptr subResultPointer = this->check(rewardPathFormula.getSubformula()); SymbolicQualitativeCheckResult const& subResult = subResultPointer->asSymbolicQualitativeCheckResult(); - return this->computeReachabilityRewardsHelper(this->getModel(), this->getModel().getTransitionMatrix(), this->getModel().getOptionalStateRewardVector(), this->getModel().getOptionalTransitionRewardMatrix(), subResult.getTruthValuesVector(), *this->linearEquationSolverFactory, qualitative); - } - - template - std::unique_ptr SymbolicDtmcPrctlModelChecker::computeReachabilityRewardsHelper(storm::models::symbolic::Model const& model, storm::dd::Add const& transitionMatrix, boost::optional> const& stateRewardVector, boost::optional> const& transitionRewardMatrix, storm::dd::Bdd const& targetStates, storm::utility::solver::SymbolicLinearEquationSolverFactory const& linearEquationSolverFactory, bool qualitative) { - - // Only compute the result if there is at least one reward model. - STORM_LOG_THROW(stateRewardVector || transitionRewardMatrix, storm::exceptions::InvalidPropertyException, "Missing reward model for formula. Skipping formula."); - - // Determine which states have a reward of infinity by definition. - storm::dd::Bdd infinityStates = storm::utility::graph::performProb1(model, transitionMatrix.notZero(), model.getReachableStates(), targetStates); - infinityStates = !infinityStates && model.getReachableStates(); - storm::dd::Bdd maybeStates = (!targetStates && !infinityStates) && model.getReachableStates(); - STORM_LOG_INFO("Found " << infinityStates.getNonZeroCount() << " 'infinity' states."); - STORM_LOG_INFO("Found " << targetStates.getNonZeroCount() << " 'target' states."); - STORM_LOG_INFO("Found " << maybeStates.getNonZeroCount() << " 'maybe' states."); - - // Check whether we need to compute exact rewards for some states. - if (qualitative) { - // Set the values for all maybe-states to 1 to indicate that their reward values - // are neither 0 nor infinity. - return std::unique_ptr(new SymbolicQuantitativeCheckResult(model.getReachableStates(), infinityStates.toAdd() * model.getManager().getConstant(storm::utility::infinity()) + maybeStates.toAdd() * model.getManager().getConstant(storm::utility::one()))); - } else { - // If there are maybe states, we need to solve an equation system. - if (!maybeStates.isZero()) { - // Create the ODD for the translation between symbolic and explicit storage. - storm::dd::Odd odd(maybeStates); - - // Create the matrix and the vector for the equation system. - storm::dd::Add maybeStatesAdd = maybeStates.toAdd(); - - // Start by cutting away all rows that do not belong to maybe states. Note that this leaves columns targeting - // non-maybe states in the matrix. - storm::dd::Add submatrix = transitionMatrix * maybeStatesAdd; - - // Then compute the state reward vector to use in the computation. - storm::dd::Add subvector = stateRewardVector ? maybeStatesAdd * stateRewardVector.get() : model.getManager().getAddZero(); - if (transitionRewardMatrix) { - subvector += (submatrix * transitionRewardMatrix.get()).sumAbstract(model.getColumnVariables()); - } - - // Finally cut away all columns targeting non-maybe states and convert the matrix into the matrix needed - // for solving the equation system (i.e. compute (I-A)). - submatrix *= maybeStatesAdd.swapVariables(model.getRowColumnMetaVariablePairs()); - submatrix = (model.getRowColumnIdentity() * maybeStatesAdd) - submatrix; - - // Solve the equation system. - std::unique_ptr> solver = linearEquationSolverFactory.create(submatrix, maybeStates, model.getRowVariables(), model.getColumnVariables(), model.getRowColumnMetaVariablePairs()); - storm::dd::Add result = solver->solveEquationSystem(model.getManager().getConstant(0.5) * maybeStatesAdd, subvector); - - return std::unique_ptr(new storm::modelchecker::SymbolicQuantitativeCheckResult(model.getReachableStates(), infinityStates.toAdd() * model.getManager().getConstant(storm::utility::infinity()) + result)); - } else { - return std::unique_ptr(new storm::modelchecker::SymbolicQuantitativeCheckResult(model.getReachableStates(), infinityStates.toAdd() * model.getManager().getConstant(storm::utility::infinity()))); - } - } + storm::dd::Add numericResult = storm::modelchecker::helper::SymbolicDtmcPrctlHelper::computeReachabilityRewards(this->getModel(), this->getModel().getTransitionMatrix(), this->getModel().getOptionalStateRewardVector(), this->getModel().getOptionalTransitionRewardMatrix(), subResult.getTruthValuesVector(), qualitative, *this->linearEquationSolverFactory); + return std::unique_ptr>(new SymbolicQuantitativeCheckResult(this->getModel().getReachableStates(), numericResult)); } template diff --git a/src/modelchecker/prctl/SymbolicDtmcPrctlModelChecker.h b/src/modelchecker/prctl/SymbolicDtmcPrctlModelChecker.h index 9dca0335b..761427f5e 100644 --- a/src/modelchecker/prctl/SymbolicDtmcPrctlModelChecker.h +++ b/src/modelchecker/prctl/SymbolicDtmcPrctlModelChecker.h @@ -7,14 +7,9 @@ namespace storm { namespace modelchecker { - template - class SymbolicCtmcCslModelChecker; - template class SymbolicDtmcPrctlModelChecker : public SymbolicPropositionalModelChecker { public: - friend class SymbolicCtmcCslModelChecker; - explicit SymbolicDtmcPrctlModelChecker(storm::models::symbolic::Dtmc const& model); explicit SymbolicDtmcPrctlModelChecker(storm::models::symbolic::Dtmc const& model, std::unique_ptr>&& linearEquationSolverFactory); @@ -23,22 +18,14 @@ namespace storm { virtual std::unique_ptr computeBoundedUntilProbabilities(storm::logic::BoundedUntilFormula const& pathFormula, bool qualitative = false, boost::optional const& optimalityType = boost::optional()) override; virtual std::unique_ptr computeNextProbabilities(storm::logic::NextFormula const& pathFormula, bool qualitative = false, boost::optional const& optimalityType = boost::optional()) override; virtual std::unique_ptr computeUntilProbabilities(storm::logic::UntilFormula const& pathFormula, bool qualitative = false, boost::optional const& optimalityType = boost::optional()) override; - virtual std::unique_ptr computeCumulativeRewards(storm::logic::CumulativeRewardFormula const& rewardPathFormula, bool qualitative = false, boost::optional const& optimalityType = boost::optional()) override; - virtual std::unique_ptr computeInstantaneousRewards(storm::logic::InstantaneousRewardFormula const& rewardPathFormula, bool qualitative = false, boost::optional const& optimalityType = boost::optional()) override; - virtual std::unique_ptr computeReachabilityRewards(storm::logic::ReachabilityRewardFormula const& rewardPathFormula, bool qualitative = false, boost::optional const& optimalityType = boost::optional()) override; + virtual std::unique_ptr computeCumulativeRewards(storm::logic::CumulativeRewardFormula const& rewardPathFormula, boost::optional const& rewardModelName = boost::optional(), bool qualitative = false, boost::optional const& optimalityType = boost::optional()) override; + virtual std::unique_ptr computeInstantaneousRewards(storm::logic::InstantaneousRewardFormula const& rewardPathFormula, boost::optional const& rewardModelName = boost::optional(), bool qualitative = false, boost::optional const& optimalityType = boost::optional()) override; + virtual std::unique_ptr computeReachabilityRewards(storm::logic::ReachabilityRewardFormula const& rewardPathFormula, boost::optional const& rewardModelName = boost::optional(), bool qualitative = false, boost::optional const& optimalityType = boost::optional()) override; protected: storm::models::symbolic::Dtmc const& getModel() const override; private: - // The methods that perform the actual checking. - static std::unique_ptr computeBoundedUntilProbabilitiesHelper(storm::models::symbolic::Model const& model, storm::dd::Add const& transitionMatrix, storm::dd::Bdd const& phiStates, storm::dd::Bdd const& psiStates, uint_fast64_t stepBound, storm::utility::solver::SymbolicLinearEquationSolverFactory const& linearEquationSolverFactory); - static storm::dd::Add computeNextProbabilitiesHelper(storm::models::symbolic::Model const& model, storm::dd::Add const& transitionMatrix, storm::dd::Bdd const& nextStates); - static std::unique_ptr computeUntilProbabilitiesHelper(storm::models::symbolic::Model const& model, storm::dd::Add const& transitionMatrix, storm::dd::Bdd const& phiStates, storm::dd::Bdd const& psiStates, bool qualitative, storm::utility::solver::SymbolicLinearEquationSolverFactory const& linearEquationSolverFactory); - static std::unique_ptr computeCumulativeRewardsHelper(storm::models::symbolic::Model const& model, storm::dd::Add const& transitionMatrix, uint_fast64_t stepBound, storm::utility::solver::SymbolicLinearEquationSolverFactory const& linearEquationSolverFactory); - static std::unique_ptr computeInstantaneousRewardsHelper(storm::models::symbolic::Model const& model, storm::dd::Add const& transitionMatrix, uint_fast64_t stepBound, storm::utility::solver::SymbolicLinearEquationSolverFactory const& linearEquationSolverFactory); - static std::unique_ptr computeReachabilityRewardsHelper(storm::models::symbolic::Model const& model, storm::dd::Add const& transitionMatrix, boost::optional> const& stateRewardVector, boost::optional> const& transitionRewardMatrix, storm::dd::Bdd const& targetStates, storm::utility::solver::SymbolicLinearEquationSolverFactory const& linearEquationSolverFactory, bool qualitative); - // An object that is used for retrieving linear equation solvers. std::unique_ptr> linearEquationSolverFactory; }; diff --git a/src/modelchecker/prctl/helper/HybridDtmcPrctlHelper.cpp b/src/modelchecker/prctl/helper/HybridDtmcPrctlHelper.cpp new file mode 100644 index 000000000..e69de29bb diff --git a/src/modelchecker/prctl/helper/HybridDtmcPrctlHelper.h b/src/modelchecker/prctl/helper/HybridDtmcPrctlHelper.h new file mode 100644 index 000000000..e69de29bb diff --git a/src/modelchecker/prctl/helper/HybridMdpPrctlHelper.cpp b/src/modelchecker/prctl/helper/HybridMdpPrctlHelper.cpp new file mode 100644 index 000000000..23cebc26a --- /dev/null +++ b/src/modelchecker/prctl/helper/HybridMdpPrctlHelper.cpp @@ -0,0 +1,248 @@ +#include "src/modelchecker/prctl/helper/HybridMdpPrctlHelper.h" + +namespace storm { + namespace modelchecker { + namespace helper { + + template + std::unique_ptr HybridMdpPrctlModelChecker::computeUntilProbabilities(bool minimize, storm::models::symbolic::NondeterministicModel const& model, storm::dd::Add const& transitionMatrix, storm::dd::Bdd const& phiStates, storm::dd::Bdd const& psiStates, bool qualitative, storm::utility::solver::MinMaxLinearEquationSolverFactory const& linearEquationSolverFactory) { + // We need to identify the states which have to be taken out of the matrix, i.e. all states that have + // probability 0 and 1 of satisfying the until-formula. + std::pair, storm::dd::Bdd> statesWithProbability01; + if (minimize) { + statesWithProbability01 = storm::utility::graph::performProb01Min(model, phiStates, psiStates); + } else { + statesWithProbability01 = storm::utility::graph::performProb01Max(model, phiStates, psiStates); + } + storm::dd::Bdd maybeStates = !statesWithProbability01.first && !statesWithProbability01.second && model.getReachableStates(); + + // Perform some logging. + STORM_LOG_INFO("Found " << statesWithProbability01.first.getNonZeroCount() << " 'no' states."); + STORM_LOG_INFO("Found " << statesWithProbability01.second.getNonZeroCount() << " 'yes' states."); + STORM_LOG_INFO("Found " << maybeStates.getNonZeroCount() << " 'maybe' states."); + + // Check whether we need to compute exact probabilities for some states. + if (qualitative) { + // Set the values for all maybe-states to 0.5 to indicate that their probability values are neither 0 nor 1. + return std::unique_ptr(new storm::modelchecker::SymbolicQuantitativeCheckResult(model.getReachableStates(), statesWithProbability01.second.toAdd() + maybeStates.toAdd() * model.getManager().getConstant(0.5))); + } else { + // If there are maybe states, we need to solve an equation system. + if (!maybeStates.isZero()) { + // Create the ODD for the translation between symbolic and explicit storage. + storm::dd::Odd odd(maybeStates); + + // Create the matrix and the vector for the equation system. + storm::dd::Add maybeStatesAdd = maybeStates.toAdd(); + + // Start by cutting away all rows that do not belong to maybe states. Note that this leaves columns targeting + // non-maybe states in the matrix. + storm::dd::Add submatrix = transitionMatrix * maybeStatesAdd; + + // Then compute the vector that contains the one-step probabilities to a state with probability 1 for all + // maybe states. + storm::dd::Add prob1StatesAsColumn = statesWithProbability01.second.toAdd(); + prob1StatesAsColumn = prob1StatesAsColumn.swapVariables(model.getRowColumnMetaVariablePairs()); + storm::dd::Add subvector = submatrix * prob1StatesAsColumn; + subvector = subvector.sumAbstract(model.getColumnVariables()); + + // Before cutting the non-maybe columns, we need to compute the sizes of the row groups. + std::vector rowGroupSizes = submatrix.notZero().existsAbstract(model.getColumnVariables()).toAdd().sumAbstract(model.getNondeterminismVariables()).template toVector(odd); + + // Finally cut away all columns targeting non-maybe states. + submatrix *= maybeStatesAdd.swapVariables(model.getRowColumnMetaVariablePairs()); + + // Create the solution vector. + std::vector x(maybeStates.getNonZeroCount(), ValueType(0.5)); + + // Translate the symbolic matrix/vector to their explicit representations and solve the equation system. + std::pair, std::vector> explicitRepresentation = submatrix.toMatrixVector(subvector, std::move(rowGroupSizes), model.getNondeterminismVariables(), odd, odd); + + std::unique_ptr> solver = linearEquationSolverFactory.create(explicitRepresentation.first); + solver->solveEquationSystem(minimize, x, explicitRepresentation.second); + + // Return a hybrid check result that stores the numerical values explicitly. + return std::unique_ptr(new storm::modelchecker::HybridQuantitativeCheckResult(model.getReachableStates(), model.getReachableStates() && !maybeStates, statesWithProbability01.second.toAdd(), maybeStates, odd, x)); + } else { + return std::unique_ptr(new storm::modelchecker::SymbolicQuantitativeCheckResult(model.getReachableStates(), statesWithProbability01.second.toAdd())); + } + } + } + + template + storm::dd::Add HybridMdpPrctlModelChecker::computeNextProbabilitiesHelper(bool minimize, storm::models::symbolic::NondeterministicModel const& model, storm::dd::Add const& transitionMatrix, storm::dd::Bdd const& nextStates) { + storm::dd::Add result = transitionMatrix * nextStates.swapVariables(model.getRowColumnMetaVariablePairs()).toAdd(); + return result.sumAbstract(model.getColumnVariables()); + } + + template + std::unique_ptr HybridMdpPrctlModelChecker::computeBoundedUntilProbabilitiesHelper(bool minimize, storm::models::symbolic::NondeterministicModel const& model, storm::dd::Add const& transitionMatrix, storm::dd::Bdd const& phiStates, storm::dd::Bdd const& psiStates, uint_fast64_t stepBound, storm::utility::solver::MinMaxLinearEquationSolverFactory const& linearEquationSolverFactory) { + // We need to identify the states which have to be taken out of the matrix, i.e. all states that have + // probability 0 or 1 of satisfying the until-formula. + storm::dd::Bdd statesWithProbabilityGreater0; + if (minimize) { + statesWithProbabilityGreater0 = storm::utility::graph::performProbGreater0A(model, transitionMatrix.notZero(), phiStates, psiStates); + } else { + statesWithProbabilityGreater0 = storm::utility::graph::performProbGreater0E(model, transitionMatrix.notZero(), phiStates, psiStates); + } + storm::dd::Bdd maybeStates = statesWithProbabilityGreater0 && !psiStates && model.getReachableStates(); + + // If there are maybe states, we need to perform matrix-vector multiplications. + if (!maybeStates.isZero()) { + // Create the ODD for the translation between symbolic and explicit storage. + storm::dd::Odd odd(maybeStates); + + // Create the matrix and the vector for the equation system. + storm::dd::Add maybeStatesAdd = maybeStates.toAdd(); + + // Start by cutting away all rows that do not belong to maybe states. Note that this leaves columns targeting + // non-maybe states in the matrix. + storm::dd::Add submatrix = transitionMatrix * maybeStatesAdd; + + // Then compute the vector that contains the one-step probabilities to a state with probability 1 for all + // maybe states. + storm::dd::Add prob1StatesAsColumn = psiStates.toAdd().swapVariables(model.getRowColumnMetaVariablePairs()); + storm::dd::Add subvector = (submatrix * prob1StatesAsColumn).sumAbstract(model.getColumnVariables()); + + // Before cutting the non-maybe columns, we need to compute the sizes of the row groups. + std::vector rowGroupSizes = submatrix.notZero().existsAbstract(model.getColumnVariables()).toAdd().sumAbstract(model.getNondeterminismVariables()).template toVector(odd); + + // Finally cut away all columns targeting non-maybe states. + submatrix *= maybeStatesAdd.swapVariables(model.getRowColumnMetaVariablePairs()); + + // Create the solution vector. + std::vector x(maybeStates.getNonZeroCount(), storm::utility::zero()); + + // Translate the symbolic matrix/vector to their explicit representations. + std::pair, std::vector> explicitRepresentation = submatrix.toMatrixVector(subvector, std::move(rowGroupSizes), model.getNondeterminismVariables(), odd, odd); + + std::unique_ptr> solver = linearEquationSolverFactory.create(explicitRepresentation.first); + solver->performMatrixVectorMultiplication(minimize, x, &explicitRepresentation.second, stepBound); + + // Return a hybrid check result that stores the numerical values explicitly. + return std::unique_ptr(new storm::modelchecker::HybridQuantitativeCheckResult(model.getReachableStates(), model.getReachableStates() && !maybeStates, psiStates.toAdd(), maybeStates, odd, x)); + } else { + return std::unique_ptr(new storm::modelchecker::SymbolicQuantitativeCheckResult(model.getReachableStates(), psiStates.toAdd())); + } + } + + template + std::unique_ptr HybridMdpPrctlModelChecker::computeInstantaneousRewardsHelper(bool minimize, storm::models::symbolic::NondeterministicModel const& model, storm::dd::Add const& transitionMatrix, uint_fast64_t stepBound, storm::utility::solver::MinMaxLinearEquationSolverFactory const& linearEquationSolverFactory) { + // Only compute the result if the model has at least one reward this->getModel(). + STORM_LOG_THROW(model.hasStateRewards(), storm::exceptions::InvalidPropertyException, "Missing reward model for formula. Skipping formula."); + + // Create the ODD for the translation between symbolic and explicit storage. + storm::dd::Odd odd(model.getReachableStates()); + + // Translate the symbolic matrix to its explicit representations. + storm::storage::SparseMatrix explicitMatrix = transitionMatrix.toMatrix(model.getNondeterminismVariables(), odd, odd); + + // Create the solution vector (and initialize it to the state rewards of the model). + std::vector x = model.getStateRewardVector().template toVector(model.getNondeterminismVariables(), odd, explicitMatrix.getRowGroupIndices()); + + // Perform the matrix-vector multiplication. + std::unique_ptr> solver = linearEquationSolverFactory.create(explicitMatrix); + solver->performMatrixVectorMultiplication(minimize, x, nullptr, stepBound); + + // Return a hybrid check result that stores the numerical values explicitly. + return std::unique_ptr(new HybridQuantitativeCheckResult(model.getReachableStates(), model.getManager().getBddZero(), model.getManager().getAddZero(), model.getReachableStates(), odd, x)); + } + + template + std::unique_ptr HybridMdpPrctlModelChecker::computeCumulativeRewardsHelper(bool minimize, storm::models::symbolic::NondeterministicModel const& model, storm::dd::Add const& transitionMatrix, uint_fast64_t stepBound, storm::utility::solver::MinMaxLinearEquationSolverFactory const& linearEquationSolverFactory) { + // Only compute the result if the model has at least one reward this->getModel(). + STORM_LOG_THROW(model.hasStateRewards() || model.hasTransitionRewards(), storm::exceptions::InvalidPropertyException, "Missing reward model for formula. Skipping formula."); + + // Compute the reward vector to add in each step based on the available reward models. + storm::dd::Add totalRewardVector = model.hasStateRewards() ? model.getStateRewardVector() : model.getManager().getAddZero(); + if (model.hasTransitionRewards()) { + totalRewardVector += (transitionMatrix * model.getTransitionRewardMatrix()).sumAbstract(model.getColumnVariables()); + } + + // Create the ODD for the translation between symbolic and explicit storage. + storm::dd::Odd odd(model.getReachableStates()); + + // Create the solution vector. + std::vector x(model.getNumberOfStates(), storm::utility::zero()); + + // Translate the symbolic matrix/vector to their explicit representations. + storm::storage::SparseMatrix explicitMatrix = transitionMatrix.toMatrix(model.getNondeterminismVariables(), odd, odd); + std::vector b = totalRewardVector.template toVector(model.getNondeterminismVariables(), odd, explicitMatrix.getRowGroupIndices()); + + // Perform the matrix-vector multiplication. + std::unique_ptr> solver = linearEquationSolverFactory.create(explicitMatrix); + solver->performMatrixVectorMultiplication(minimize, x, &b, stepBound); + + // Return a hybrid check result that stores the numerical values explicitly. + return std::unique_ptr(new HybridQuantitativeCheckResult(model.getReachableStates(), model.getManager().getBddZero(), model.getManager().getAddZero(), model.getReachableStates(), odd, x)); + } + + template + std::unique_ptr HybridMdpPrctlModelChecker::computeReachabilityRewardsHelper(bool minimize, storm::models::symbolic::NondeterministicModel const& model, storm::dd::Add const& transitionMatrix, boost::optional> const& stateRewardVector, boost::optional> const& transitionRewardMatrix, storm::dd::Bdd const& targetStates, storm::utility::solver::MinMaxLinearEquationSolverFactory const& linearEquationSolverFactory, bool qualitative) { + + // Only compute the result if there is at least one reward model. + STORM_LOG_THROW(stateRewardVector || transitionRewardMatrix, storm::exceptions::InvalidPropertyException, "Missing reward model for formula. Skipping formula."); + + // Determine which states have a reward of infinity by definition. + storm::dd::Bdd infinityStates; + storm::dd::Bdd transitionMatrixBdd = transitionMatrix.notZero(); + if (minimize) { + infinityStates = storm::utility::graph::performProb1A(model, transitionMatrixBdd, model.getReachableStates(), targetStates, storm::utility::graph::performProbGreater0A(model, transitionMatrixBdd, model.getReachableStates(), targetStates)); + } else { + infinityStates = storm::utility::graph::performProb1E(model, transitionMatrixBdd, model.getReachableStates(), targetStates, storm::utility::graph::performProbGreater0E(model, transitionMatrixBdd, model.getReachableStates(), targetStates)); + } + infinityStates = !infinityStates && model.getReachableStates(); + storm::dd::Bdd maybeStates = (!targetStates && !infinityStates) && model.getReachableStates(); + STORM_LOG_INFO("Found " << infinityStates.getNonZeroCount() << " 'infinity' states."); + STORM_LOG_INFO("Found " << targetStates.getNonZeroCount() << " 'target' states."); + STORM_LOG_INFO("Found " << maybeStates.getNonZeroCount() << " 'maybe' states."); + + // Check whether we need to compute exact rewards for some states. + if (qualitative) { + // Set the values for all maybe-states to 1 to indicate that their reward values + // are neither 0 nor infinity. + return std::unique_ptr(new SymbolicQuantitativeCheckResult(model.getReachableStates(), infinityStates.toAdd() * model.getManager().getConstant(storm::utility::infinity()) + maybeStates.toAdd() * model.getManager().getConstant(storm::utility::one()))); + } else { + // If there are maybe states, we need to solve an equation system. + if (!maybeStates.isZero()) { + // Create the ODD for the translation between symbolic and explicit storage. + storm::dd::Odd odd(maybeStates); + + // Create the matrix and the vector for the equation system. + storm::dd::Add maybeStatesAdd = maybeStates.toAdd(); + + // Start by cutting away all rows that do not belong to maybe states. Note that this leaves columns targeting + // non-maybe states in the matrix. + storm::dd::Add submatrix = transitionMatrix * maybeStatesAdd; + + // Then compute the state reward vector to use in the computation. + storm::dd::Add subvector = stateRewardVector ? maybeStatesAdd * stateRewardVector.get() : model.getManager().getAddZero(); + if (transitionRewardMatrix) { + subvector += (submatrix * transitionRewardMatrix.get()).sumAbstract(model.getColumnVariables()); + } + + // Before cutting the non-maybe columns, we need to compute the sizes of the row groups. + std::vector rowGroupSizes = submatrix.notZero().existsAbstract(model.getColumnVariables()).toAdd().sumAbstract(model.getNondeterminismVariables()).template toVector(odd); + + // Finally cut away all columns targeting non-maybe states. + submatrix *= maybeStatesAdd.swapVariables(model.getRowColumnMetaVariablePairs()); + + // Create the solution vector. + std::vector x(maybeStates.getNonZeroCount(), ValueType(0.5)); + + // Translate the symbolic matrix/vector to their explicit representations. + std::pair, std::vector> explicitRepresentation = submatrix.toMatrixVector(subvector, std::move(rowGroupSizes), model.getNondeterminismVariables(), odd, odd); + + // Now solve the resulting equation system. + std::unique_ptr> solver = linearEquationSolverFactory.create(explicitRepresentation.first); + solver->solveEquationSystem(minimize, x, explicitRepresentation.second); + + // Return a hybrid check result that stores the numerical values explicitly. + return std::unique_ptr(new storm::modelchecker::HybridQuantitativeCheckResult(model.getReachableStates(), model.getReachableStates() && !maybeStates, infinityStates.toAdd() * model.getManager().getConstant(storm::utility::infinity()), maybeStates, odd, x)); + } else { + return std::unique_ptr(new storm::modelchecker::SymbolicQuantitativeCheckResult(model.getReachableStates(), infinityStates.toAdd() * model.getManager().getConstant(storm::utility::infinity()))); + } + } + } + } + } +} \ No newline at end of file diff --git a/src/modelchecker/prctl/helper/HybridMdpPrctlHelper.h b/src/modelchecker/prctl/helper/HybridMdpPrctlHelper.h new file mode 100644 index 000000000..e55b5762b --- /dev/null +++ b/src/modelchecker/prctl/helper/HybridMdpPrctlHelper.h @@ -0,0 +1,34 @@ +#ifndef STORM_MODELCHECKER_SPARSE_MDP_PRCTL_MODELCHECKER_HELPER_H_ +#define STORM_MODELCHECKER_SPARSE_MDP_PRCTL_MODELCHECKER_HELPER_H_ + +#include "src/models/symbolic/Model.h" + +#include "src/storage/dd/Add.h" +#include "src/storage/dd/Bdd.h" + +#include "src/utility/solver.h" + +namespace storm { + namespace modelchecker { + namespace helper { + + template + class HybridMdpPrctlHelper { + public: + static std::unique_ptr computeBoundedUntilProbabilities(bool minimize, storm::models::symbolic::NondeterministicModel const& model, storm::dd::Add const& transitionMatrix, storm::dd::Bdd const& phiStates, storm::dd::Bdd const& psiStates, uint_fast64_t stepBound, storm::utility::solver::MinMaxLinearEquationSolverFactory const& linearEquationSolverFactory); + + static storm::dd::Add computeNextProbabilities(bool minimize, storm::models::symbolic::NondeterministicModel const& model, storm::dd::Add const& transitionMatrix, storm::dd::Bdd const& nextStates); + + static std::unique_ptr computeUntilProbabilities(bool minimize, storm::models::symbolic::NondeterministicModel const& model, storm::dd::Add const& transitionMatrix, storm::dd::Bdd const& phiStates, storm::dd::Bdd const& psiStates, bool qualitative, storm::utility::solver::MinMaxLinearEquationSolverFactory const& linearEquationSolverFactory); + + static std::unique_ptr computeCumulativeRewards(bool minimize, storm::models::symbolic::NondeterministicModel const& model, storm::dd::Add const& transitionMatrix, uint_fast64_t stepBound, storm::utility::solver::MinMaxLinearEquationSolverFactory const& linearEquationSolverFactory); + + static std::unique_ptr computeInstantaneousRewards(bool minimize, storm::models::symbolic::NondeterministicModel const& model, storm::dd::Add const& transitionMatrix, uint_fast64_t stepBound, storm::utility::solver::MinMaxLinearEquationSolverFactory const& linearEquationSolverFactory); + + static std::unique_ptr computeReachabilityRewards(bool minimize, storm::models::symbolic::NondeterministicModel const& model, storm::dd::Add const& transitionMatrix, boost::optional> const& stateRewardVector, boost::optional> const& transitionRewardMatrix, storm::dd::Bdd const& targetStates, bool qualitative, storm::utility::solver::MinMaxLinearEquationSolverFactory const& linearEquationSolverFactory); + } + } + } +} + +#endif /* STORM_MODELCHECKER_SPARSE_MDP_PRCTL_MODELCHECKER_HELPER_H_ */ diff --git a/src/modelchecker/prctl/helper/SparseDtmcPrctlHelper.cpp b/src/modelchecker/prctl/helper/SparseDtmcPrctlHelper.cpp new file mode 100644 index 000000000..0cb95666f --- /dev/null +++ b/src/modelchecker/prctl/helper/SparseDtmcPrctlHelper.cpp @@ -0,0 +1,196 @@ +#include "src/modelchecker/prctl/helper/SparseDtmcPrctlHelper.h" + +#include "src/utility/macros.h" +#include "src/utility/vector.h" +#include "src/utility/graph.h" + +#include "src/exceptions/InvalidStateException.h" +#include "src/exceptions/InvalidPropertyException.h" + +namespace storm { + namespace modelchecker { + namespace helper { + template + std::vector SparseDtmcPrctlHelper::computeBoundedUntilProbabilities(storm::storage::SparseMatrix const& transitionMatrix, storm::storage::SparseMatrix const& backwardTransitions, storm::storage::BitVector const& phiStates, storm::storage::BitVector const& psiStates, uint_fast64_t stepBound, storm::utility::solver::LinearEquationSolverFactory const& linearEquationSolverFactory) { + std::vector result(transitionMatrix.getRowCount(), storm::utility::zero()); + + // If we identify the states that have probability 0 of reaching the target states, we can exclude them in the further analysis. + storm::storage::BitVector maybeStates = storm::utility::graph::performProbGreater0(backwardTransitions, phiStates, psiStates, true, stepBound); + maybeStates &= ~psiStates; + STORM_LOG_INFO("Found " << maybeStates.getNumberOfSetBits() << " 'maybe' states."); + + if (!maybeStates.empty()) { + // We can eliminate the rows and columns from the original transition probability matrix that have probability 0. + storm::storage::SparseMatrix submatrix = transitionMatrix.getSubmatrix(true, maybeStates, maybeStates, true); + + // Create the vector of one-step probabilities to go to target states. + std::vector b = transitionMatrix.getConstrainedRowSumVector(maybeStates, psiStates); + + // Create the vector with which to multiply. + std::vector subresult(maybeStates.getNumberOfSetBits()); + + // Perform the matrix vector multiplication as often as required by the formula bound. + std::unique_ptr> solver = linearEquationSolverFactory.create(submatrix); + solver->performMatrixVectorMultiplication(subresult, &b, stepBound); + + // Set the values of the resulting vector accordingly. + storm::utility::vector::setVectorValues(result, maybeStates, subresult); + } + storm::utility::vector::setVectorValues(result, psiStates, storm::utility::one()); + + return result; + } + + template + std::vector SparseDtmcPrctlHelper::computeUntilProbabilities(storm::storage::SparseMatrix const& transitionMatrix, storm::storage::SparseMatrix const& backwardTransitions, storm::storage::BitVector const& phiStates, storm::storage::BitVector const& psiStates, bool qualitative, storm::utility::solver::LinearEquationSolverFactory const& linearEquationSolverFactory) { + // We need to identify the states which have to be taken out of the matrix, i.e. + // all states that have probability 0 and 1 of satisfying the until-formula. + std::pair statesWithProbability01 = storm::utility::graph::performProb01(backwardTransitions, phiStates, psiStates); + storm::storage::BitVector statesWithProbability0 = std::move(statesWithProbability01.first); + storm::storage::BitVector statesWithProbability1 = std::move(statesWithProbability01.second); + + // Perform some logging. + storm::storage::BitVector maybeStates = ~(statesWithProbability0 | statesWithProbability1); + STORM_LOG_INFO("Found " << statesWithProbability0.getNumberOfSetBits() << " 'no' states."); + STORM_LOG_INFO("Found " << statesWithProbability1.getNumberOfSetBits() << " 'yes' states."); + STORM_LOG_INFO("Found " << maybeStates.getNumberOfSetBits() << " 'maybe' states."); + + // Create resulting vector. + std::vector result(transitionMatrix.getRowCount()); + + // Check whether we need to compute exact probabilities for some states. + if (qualitative) { + // Set the values for all maybe-states to 0.5 to indicate that their probability values are neither 0 nor 1. + storm::utility::vector::setVectorValues(result, maybeStates, ValueType(0.5)); + } else { + if (!maybeStates.empty()) { + // In this case we have have to compute the probabilities. + + // We can eliminate the rows and columns from the original transition probability matrix. + storm::storage::SparseMatrix submatrix = transitionMatrix.getSubmatrix(true, maybeStates, maybeStates, true); + + // Converting the matrix from the fixpoint notation to the form needed for the equation + // system. That is, we go from x = A*x + b to (I-A)x = b. + submatrix.convertToEquationSystem(); + + // Initialize the x vector with 0.5 for each element. This is the initial guess for + // the iterative solvers. It should be safe as for all 'maybe' states we know that the + // probability is strictly larger than 0. + std::vector x(maybeStates.getNumberOfSetBits(), ValueType(0.5)); + + // Prepare the right-hand side of the equation system. For entry i this corresponds to + // the accumulated probability of going from state i to some 'yes' state. + std::vector b = transitionMatrix.getConstrainedRowSumVector(maybeStates, statesWithProbability1); + + // Now solve the created system of linear equations. + std::unique_ptr> solver = linearEquationSolverFactory.create(submatrix); + solver->solveEquationSystem(x, b); + + // Set values of resulting vector according to result. + storm::utility::vector::setVectorValues(result, maybeStates, x); + } + } + + // Set values of resulting vector that are known exactly. + storm::utility::vector::setVectorValues(result, statesWithProbability0, storm::utility::zero()); + storm::utility::vector::setVectorValues(result, statesWithProbability1, storm::utility::one()); + + return result; + } + + template + std::vector SparseDtmcPrctlHelper::computeNextProbabilities(storm::storage::SparseMatrix const& transitionMatrix, storm::storage::BitVector const& nextStates, storm::utility::solver::LinearEquationSolverFactory const& linearEquationSolverFactory) { + // Create the vector with which to multiply and initialize it correctly. + std::vector result(transitionMatrix.getRowCount()); + storm::utility::vector::setVectorValues(result, nextStates, storm::utility::one()); + + // Perform one single matrix-vector multiplication. + std::unique_ptr> solver = linearEquationSolverFactory.create(transitionMatrix); + solver->performMatrixVectorMultiplication(result); + return result; + } + + template + std::vector SparseDtmcPrctlHelper::computeCumulativeRewards(storm::storage::SparseMatrix const& transitionMatrix, RewardModelType const& rewardModel, uint_fast64_t stepBound, storm::utility::solver::LinearEquationSolverFactory const& linearEquationSolverFactory) { + // Compute the reward vector to add in each step based on the available reward models. + std::vector totalRewardVector = rewardModel.getTotalRewardVector(transitionMatrix); + + // Initialize result to either the state rewards of the model or the null vector. + std::vector result = rewardModel.getTotalStateActionRewardVector(transitionMatrix.getRowCount(), transitionMatrix.getRowGroupIndices()); + + // Perform the matrix vector multiplication as often as required by the formula bound. + std::unique_ptr> solver = linearEquationSolverFactory.create(transitionMatrix); + solver->performMatrixVectorMultiplication(result, &totalRewardVector, stepBound); + + return result; + } + + template + std::vector SparseDtmcPrctlHelper::computeInstantaneousRewards(storm::storage::SparseMatrix const& transitionMatrix, RewardModelType const& rewardModel, uint_fast64_t stepCount, storm::utility::solver::LinearEquationSolverFactory const& linearEquationSolverFactory) { + // Only compute the result if the model has a state-based reward this->getModel(). + STORM_LOG_THROW(rewardModel.hasStateRewards() || rewardModel.hasStateActionRewards(), storm::exceptions::InvalidPropertyException, "Missing reward model for formula. Skipping formula."); + + // Initialize result to state rewards of the model. + std::vector result = rewardModel.getTotalStateActionRewardVector(transitionMatrix.getRowCount(), transitionMatrix.getRowGroupIndices()); + + // Perform the matrix vector multiplication as often as required by the formula bound. + std::unique_ptr> solver = linearEquationSolverFactory.create(transitionMatrix); + solver->performMatrixVectorMultiplication(result, nullptr, stepCount); + + return result; + } + + template + std::vector SparseDtmcPrctlHelper::computeReachabilityRewards(storm::storage::SparseMatrix const& transitionMatrix, storm::storage::SparseMatrix const& backwardTransitions, RewardModelType const& rewardModel, storm::storage::BitVector const& targetStates, bool qualitative, storm::utility::solver::LinearEquationSolverFactory const& linearEquationSolverFactory) { + + // Determine which states have a reward of infinity by definition. + storm::storage::BitVector trueStates(transitionMatrix.getRowCount(), true); + storm::storage::BitVector infinityStates = storm::utility::graph::performProb1(backwardTransitions, trueStates, targetStates); + infinityStates.complement(); + storm::storage::BitVector maybeStates = ~targetStates & ~infinityStates; + STORM_LOG_INFO("Found " << infinityStates.getNumberOfSetBits() << " 'infinity' states."); + STORM_LOG_INFO("Found " << targetStates.getNumberOfSetBits() << " 'target' states."); + STORM_LOG_INFO("Found " << maybeStates.getNumberOfSetBits() << " 'maybe' states."); + + // Create resulting vector. + std::vector result(transitionMatrix.getRowCount(), storm::utility::zero()); + + // Check whether we need to compute exact rewards for some states. + if (qualitative) { + // Set the values for all maybe-states to 1 to indicate that their reward values + // are neither 0 nor infinity. + storm::utility::vector::setVectorValues(result, maybeStates, storm::utility::one()); + } else { + // In this case we have to compute the reward values for the remaining states. + // We can eliminate the rows and columns from the original transition probability matrix. + storm::storage::SparseMatrix submatrix = transitionMatrix.getSubmatrix(true, maybeStates, maybeStates, true); + + // Converting the matrix from the fixpoint notation to the form needed for the equation + // system. That is, we go from x = A*x + b to (I-A)x = b. + submatrix.convertToEquationSystem(); + + // Initialize the x vector with 1 for each element. This is the initial guess for + // the iterative solvers. + std::vector x(submatrix.getColumnCount(), storm::utility::one()); + + // Prepare the right-hand side of the equation system. + std::vector b = rewardModel.getTotalRewardVector(submatrix.getRowCount(), transitionMatrix, maybeStates); + + // Now solve the resulting equation system. + std::unique_ptr> solver = linearEquationSolverFactory.create(submatrix); + solver->solveEquationSystem(x, b); + + // Set values of resulting vector according to result. + storm::utility::vector::setVectorValues(result, maybeStates, x); + } + + // Set values of resulting vector that are known exactly. + storm::utility::vector::setVectorValues(result, infinityStates, storm::utility::infinity()); + + return result; + } + + template class SparseDtmcPrctlHelper; + } + } +} diff --git a/src/modelchecker/prctl/helper/SparseDtmcPrctlHelper.h b/src/modelchecker/prctl/helper/SparseDtmcPrctlHelper.h new file mode 100644 index 000000000..57d2a9de7 --- /dev/null +++ b/src/modelchecker/prctl/helper/SparseDtmcPrctlHelper.h @@ -0,0 +1,39 @@ +#ifndef STORM_MODELCHECKER_SPARSE_DTMC_PRCTL_MODELCHECKER_HELPER_H_ +#define STORM_MODELCHECKER_SPARSE_DTMC_PRCTL_MODELCHECKER_HELPER_H_ + +#include + +#include "src/models/sparse/StandardRewardModel.h" + +#include "src/storage/SparseMatrix.h" +#include "src/storage/BitVector.h" + +#include "src/utility/solver.h" + +namespace storm { + namespace modelchecker { + namespace helper { + + template > + class SparseDtmcPrctlHelper { + public: + static std::vector computeBoundedUntilProbabilities(storm::storage::SparseMatrix const& transitionMatrix, storm::storage::SparseMatrix const& backwardTransitions, storm::storage::BitVector const& phiStates, storm::storage::BitVector const& psiStates, uint_fast64_t stepBound, storm::utility::solver::LinearEquationSolverFactory const& linearEquationSolverFactory); + + static std::vector computeNextProbabilities(storm::storage::SparseMatrix const& transitionMatrix, storm::storage::BitVector const& nextStates, storm::utility::solver::LinearEquationSolverFactory const& linearEquationSolverFactory); + + static std::vector computeUntilProbabilities(storm::storage::SparseMatrix const& transitionMatrix, storm::storage::SparseMatrix const& backwardTransitions, storm::storage::BitVector const& phiStates, storm::storage::BitVector const& psiStates, bool qualitative, storm::utility::solver::LinearEquationSolverFactory const& linearEquationSolverFactory); + + static std::vector computeCumulativeRewards(storm::storage::SparseMatrix const& transitionMatrix, RewardModelType const& rewardModel, uint_fast64_t stepBound, storm::utility::solver::LinearEquationSolverFactory const& linearEquationSolverFactory); + + static std::vector computeInstantaneousRewards(storm::storage::SparseMatrix const& transitionMatrix, RewardModelType const& rewardModel, uint_fast64_t stepCount, storm::utility::solver::LinearEquationSolverFactory const& linearEquationSolverFactory); + + static std::vector computeReachabilityRewards(storm::storage::SparseMatrix const& transitionMatrix, storm::storage::SparseMatrix const& backwardTransitions, RewardModelType const& rewardModel, storm::storage::BitVector const& targetStates, bool qualitative, storm::utility::solver::LinearEquationSolverFactory const& linearEquationSolverFactory); + + static std::vector computeLongRunAverage(storm::storage::SparseMatrix const& transitionMatrix, storm::storage::BitVector const& psiStates, bool qualitative, storm::utility::solver::LinearEquationSolverFactory const& linearEquationSolverFactory); + }; + + } + } +} + +#endif /* STORM_MODELCHECKER_SPARSE_DTMC_PRCTL_MODELCHECKER_HELPER_H_ */ \ No newline at end of file diff --git a/src/modelchecker/prctl/helper/SparseMdpPrctlHelper.cpp b/src/modelchecker/prctl/helper/SparseMdpPrctlHelper.cpp new file mode 100644 index 000000000..ff278e20e --- /dev/null +++ b/src/modelchecker/prctl/helper/SparseMdpPrctlHelper.cpp @@ -0,0 +1,407 @@ +#include "src/modelchecker/prctl/helper/SparseMdpPrctlHelper.h" + +#include "src/storage/MaximalEndComponentDecomposition.h" + +#include "src/utility/macros.h" +#include "src/utility/vector.h" +#include "src/utility/graph.h" + +#include "src/exceptions/InvalidStateException.h" +#include "src/exceptions/InvalidPropertyException.h" + +namespace storm { + namespace modelchecker { + namespace helper { + template + std::vector SparseMdpPrctlHelper::computeBoundedUntilProbabilities(bool minimize, storm::storage::SparseMatrix const& transitionMatrix, storm::storage::SparseMatrix const& backwardTransitions, storm::storage::BitVector const& phiStates, storm::storage::BitVector const& psiStates, uint_fast64_t stepBound, storm::utility::solver::MinMaxLinearEquationSolverFactory const& minMaxLinearEquationSolverFactory) { + std::vector result(transitionMatrix.getRowCount(), storm::utility::zero()); + + // Determine the states that have 0 probability of reaching the target states. + storm::storage::BitVector maybeStates; + if (minimize) { + maybeStates = storm::utility::graph::performProbGreater0A(transitionMatrix, transitionMatrix.getRowGroupIndices(), backwardTransitions, phiStates, psiStates, true, stepBound); + } else { + maybeStates = storm::utility::graph::performProbGreater0E(transitionMatrix, transitionMatrix.getRowGroupIndices(), backwardTransitions, phiStates, psiStates, true, stepBound); + } + maybeStates &= ~psiStates; + STORM_LOG_INFO("Found " << maybeStates.getNumberOfSetBits() << " 'maybe' states."); + + if (!maybeStates.empty()) { + // We can eliminate the rows and columns from the original transition probability matrix that have probability 0. + storm::storage::SparseMatrix submatrix = transitionMatrix.getSubmatrix(true, maybeStates, maybeStates, false); + std::vector b = transitionMatrix.getConstrainedRowGroupSumVector(maybeStates, psiStates); + + // Create the vector with which to multiply. + std::vector subresult(maybeStates.getNumberOfSetBits()); + + std::unique_ptr> solver = minMaxLinearEquationSolverFactory.create(submatrix); + solver->performMatrixVectorMultiplication(minimize, subresult, &b, stepBound); + + // Set the values of the resulting vector accordingly. + storm::utility::vector::setVectorValues(result, maybeStates, subresult); + } + storm::utility::vector::setVectorValues(result, psiStates, storm::utility::one()); + + return result; + } + + template + std::vector SparseMdpPrctlHelper::computeNextProbabilities(bool minimize, storm::storage::SparseMatrix const& transitionMatrix, storm::storage::BitVector const& nextStates, storm::utility::solver::MinMaxLinearEquationSolverFactory const& minMaxLinearEquationSolverFactory) { + // Create the vector with which to multiply and initialize it correctly. + std::vector result(transitionMatrix.getRowCount()); + storm::utility::vector::setVectorValues(result, nextStates, storm::utility::one()); + + std::unique_ptr> solver = minMaxLinearEquationSolverFactory.create(transitionMatrix); + solver->performMatrixVectorMultiplication(minimize, result); + + return result; + } + + template + std::vector SparseMdpPrctlHelper::computeUntilProbabilities(bool minimize, storm::storage::SparseMatrix const& transitionMatrix, storm::storage::SparseMatrix const& backwardTransitions, storm::storage::BitVector const& phiStates, storm::storage::BitVector const& psiStates, bool qualitative, storm::utility::solver::MinMaxLinearEquationSolverFactory const& minMaxLinearEquationSolverFactory) { + uint_fast64_t numberOfStates = transitionMatrix.getRowCount(); + + // We need to identify the states which have to be taken out of the matrix, i.e. + // all states that have probability 0 and 1 of satisfying the until-formula. + std::pair statesWithProbability01; + if (minimize) { + statesWithProbability01 = storm::utility::graph::performProb01Min(transitionMatrix, transitionMatrix.getRowGroupIndices(), backwardTransitions, phiStates, psiStates); + } else { + statesWithProbability01 = storm::utility::graph::performProb01Max(transitionMatrix, transitionMatrix.getRowGroupIndices(), backwardTransitions, phiStates, psiStates); + } + storm::storage::BitVector statesWithProbability0 = std::move(statesWithProbability01.first); + storm::storage::BitVector statesWithProbability1 = std::move(statesWithProbability01.second); + storm::storage::BitVector maybeStates = ~(statesWithProbability0 | statesWithProbability1); + LOG4CPLUS_INFO(logger, "Found " << statesWithProbability0.getNumberOfSetBits() << " 'no' states."); + LOG4CPLUS_INFO(logger, "Found " << statesWithProbability1.getNumberOfSetBits() << " 'yes' states."); + LOG4CPLUS_INFO(logger, "Found " << maybeStates.getNumberOfSetBits() << " 'maybe' states."); + + // Create resulting vector. + std::vector result(numberOfStates); + + // Check whether we need to compute exact probabilities for some states. + if (qualitative) { + // Set the values for all maybe-states to 0.5 to indicate that their probability values are neither 0 nor 1. + storm::utility::vector::setVectorValues(result, maybeStates, ValueType(0.5)); + } else { + if (!maybeStates.empty()) { + // In this case we have have to compute the probabilities. + + // First, we can eliminate the rows and columns from the original transition probability matrix for states + // whose probabilities are already known. + storm::storage::SparseMatrix submatrix = transitionMatrix.getSubmatrix(true, maybeStates, maybeStates, false); + + // Prepare the right-hand side of the equation system. For entry i this corresponds to + // the accumulated probability of going from state i to some 'yes' state. + std::vector b = transitionMatrix.getConstrainedRowGroupSumVector(maybeStates, statesWithProbability1); + + // Create vector for results for maybe states. + std::vector x(maybeStates.getNumberOfSetBits()); + + // Solve the corresponding system of equations. + std::unique_ptr> solver = minMaxLinearEquationSolverFactory.create(submatrix); + solver->solveEquationSystem(minimize, x, b); + + // Set values of resulting vector according to result. + storm::utility::vector::setVectorValues(result, maybeStates, x); + } + } + + // Set values of resulting vector that are known exactly. + storm::utility::vector::setVectorValues(result, statesWithProbability0, storm::utility::zero()); + storm::utility::vector::setVectorValues(result, statesWithProbability1, storm::utility::one()); + + return result; + } + + template + std::vector SparseMdpPrctlHelper::computeInstantaneousRewards(storm::storage::SparseMatrix const& transitionMatrix, RewardModelType const& rewardModel, bool minimize, uint_fast64_t stepCount, storm::utility::solver::MinMaxLinearEquationSolverFactory const& minMaxLinearEquationSolverFactory) { + // Only compute the result if the model has a state-based reward this->getModel(). + STORM_LOG_THROW(rewardModel.hasStateRewards(), storm::exceptions::InvalidPropertyException, "Missing reward model for formula. Skipping formula."); + + // Initialize result to state rewards of the this->getModel(). + std::vector result(rewardModel.getStateRewardVector()); + + std::unique_ptr> solver = minMaxLinearEquationSolverFactory.create(transitionMatrix); + solver->performMatrixVectorMultiplication(minimize, result, nullptr, stepCount); + + return result; + } + + template + std::vector SparseMdpPrctlHelper::computeCumulativeRewards(storm::storage::SparseMatrix const& transitionMatrix, RewardModelType const& rewardModel, bool minimize, uint_fast64_t stepBound, storm::utility::solver::MinMaxLinearEquationSolverFactory const& minMaxLinearEquationSolverFactory) { + // Only compute the result if the model has at least one reward this->getModel(). + STORM_LOG_THROW(!rewardModel.empty(), storm::exceptions::InvalidPropertyException, "Missing reward model for formula. Skipping formula."); + + // Compute the reward vector to add in each step based on the available reward models. + std::vector totalRewardVector = rewardModel.getTotalRewardVector(transitionMatrix); + + // Initialize result to either the state rewards of the model or the null vector. + std::vector result; + if (rewardModel.hasStateRewards()) { + result = rewardModel.getStateRewardVector(); + } else { + result.resize(transitionMatrix.getRowCount()); + } + + std::unique_ptr> solver = minMaxLinearEquationSolverFactory.create(transitionMatrix); + solver->performMatrixVectorMultiplication(minimize, result, &totalRewardVector, stepBound); + + return result; + } + + template + std::vector SparseMdpPrctlHelper::computeReachabilityRewards(bool minimize, storm::storage::SparseMatrix const& transitionMatrix, storm::storage::SparseMatrix const& backwardTransitions, RewardModelType const& rewardModel, storm::storage::BitVector const& targetStates, bool qualitative, storm::utility::solver::MinMaxLinearEquationSolverFactory const& minMaxLinearEquationSolverFactory) { + // Only compute the result if the model has at least one reward this->getModel(). + STORM_LOG_THROW(!rewardModel.empty(), storm::exceptions::InvalidPropertyException, "Missing reward model for formula. Skipping formula."); + + // Determine which states have a reward of infinity by definition. + storm::storage::BitVector infinityStates; + storm::storage::BitVector trueStates(transitionMatrix.getRowCount(), true); + if (minimize) { + infinityStates = std::move(storm::utility::graph::performProb1A(transitionMatrix, transitionMatrix.getRowGroupIndices(), backwardTransitions, trueStates, targetStates)); + } else { + infinityStates = std::move(storm::utility::graph::performProb1E(transitionMatrix, transitionMatrix.getRowGroupIndices(), backwardTransitions, trueStates, targetStates)); + } + infinityStates.complement(); + storm::storage::BitVector maybeStates = ~targetStates & ~infinityStates; + LOG4CPLUS_INFO(logger, "Found " << infinityStates.getNumberOfSetBits() << " 'infinity' states."); + LOG4CPLUS_INFO(logger, "Found " << targetStates.getNumberOfSetBits() << " 'target' states."); + LOG4CPLUS_INFO(logger, "Found " << maybeStates.getNumberOfSetBits() << " 'maybe' states."); + + // Create resulting vector. + std::vector result(transitionMatrix.getRowCount()); + + // Check whether we need to compute exact rewards for some states. + if (qualitative) { + LOG4CPLUS_INFO(logger, "The rewards for the initial states were determined in a preprocessing step. No exact rewards were computed."); + // Set the values for all maybe-states to 1 to indicate that their reward values + // are neither 0 nor infinity. + storm::utility::vector::setVectorValues(result, maybeStates, storm::utility::one()); + } else { + // In this case we have to compute the reward values for the remaining states. + + // We can eliminate the rows and columns from the original transition probability matrix for states + // whose reward values are already known. + storm::storage::SparseMatrix submatrix = transitionMatrix.getSubmatrix(true, maybeStates, maybeStates, false); + + // Prepare the right-hand side of the equation system. + std::vector b = rewardModel.getTotalRewardVector(submatrix.getRowCount(), transitionMatrix, maybeStates); + + // Create vector for results for maybe states. + std::vector x(maybeStates.getNumberOfSetBits()); + + // Solve the corresponding system of equations. + std::unique_ptr> solver = minMaxLinearEquationSolverFactory.create(submatrix); + solver->solveEquationSystem(minimize, x, b); + + + // Set values of resulting vector according to result. + storm::utility::vector::setVectorValues(result, maybeStates, x); + } + + // Set values of resulting vector that are known exactly. + storm::utility::vector::setVectorValues(result, targetStates, storm::utility::zero()); + storm::utility::vector::setVectorValues(result, infinityStates, storm::utility::infinity()); + + return result; + } + + template + std::vector SparseMdpPrctlHelper::computeLongRunAverage(bool minimize, storm::storage::SparseMatrix const& transitionMatrix, storm::storage::SparseMatrix const& backwardTransitions, storm::storage::BitVector const& psiStates, bool qualitative, storm::utility::solver::MinMaxLinearEquationSolverFactory const& minMaxLinearEquationSolverFactory) { + // If there are no goal states, we avoid the computation and directly return zero. + uint_fast64_t numberOfStates = transitionMatrix.getRowCount(); + if (psiStates.empty()) { + return std::vector(numberOfStates, storm::utility::zero()); + } + + // Likewise, if all bits are set, we can avoid the computation and set. + if ((~psiStates).empty()) { + return std::vector(numberOfStates, storm::utility::one()); + } + + // Start by decomposing the MDP into its MECs. + storm::storage::MaximalEndComponentDecomposition mecDecomposition(transitionMatrix, backwardTransitions); + + // Get some data members for convenience. + std::vector const& nondeterministicChoiceIndices = transitionMatrix.getRowGroupIndices(); + ValueType zero = storm::utility::zero(); + + //first calculate LRA for the Maximal End Components. + storm::storage::BitVector statesInMecs(numberOfStates); + std::vector stateToMecIndexMap(transitionMatrix.getColumnCount()); + std::vector lraValuesForEndComponents(mecDecomposition.size(), zero); + + for (uint_fast64_t currentMecIndex = 0; currentMecIndex < mecDecomposition.size(); ++currentMecIndex) { + storm::storage::MaximalEndComponent const& mec = mecDecomposition[currentMecIndex]; + + lraValuesForEndComponents[currentMecIndex] = computeLraForMaximalEndComponent(minimize, transitionMatrix, psiStates, mec); + + // Gather information for later use. + for (auto const& stateChoicesPair : mec) { + statesInMecs.set(stateChoicesPair.first); + stateToMecIndexMap[stateChoicesPair.first] = currentMecIndex; + } + } + + // For fast transition rewriting, we build some auxiliary data structures. + storm::storage::BitVector statesNotContainedInAnyMec = ~statesInMecs; + uint_fast64_t firstAuxiliaryStateIndex = statesNotContainedInAnyMec.getNumberOfSetBits(); + uint_fast64_t lastStateNotInMecs = 0; + uint_fast64_t numberOfStatesNotInMecs = 0; + std::vector statesNotInMecsBeforeIndex; + statesNotInMecsBeforeIndex.reserve(numberOfStates); + for (auto state : statesNotContainedInAnyMec) { + while (lastStateNotInMecs <= state) { + statesNotInMecsBeforeIndex.push_back(numberOfStatesNotInMecs); + ++lastStateNotInMecs; + } + ++numberOfStatesNotInMecs; + } + + // Finally, we are ready to create the SSP matrix and right-hand side of the SSP. + std::vector b; + typename storm::storage::SparseMatrixBuilder sspMatrixBuilder(0, 0, 0, false, true, numberOfStatesNotInMecs + mecDecomposition.size()); + + // If the source state is not contained in any MEC, we copy its choices (and perform the necessary modifications). + uint_fast64_t currentChoice = 0; + for (auto state : statesNotContainedInAnyMec) { + sspMatrixBuilder.newRowGroup(currentChoice); + + for (uint_fast64_t choice = nondeterministicChoiceIndices[state]; choice < nondeterministicChoiceIndices[state + 1]; ++choice, ++currentChoice) { + std::vector auxiliaryStateToProbabilityMap(mecDecomposition.size()); + b.push_back(storm::utility::zero()); + + for (auto element : transitionMatrix.getRow(choice)) { + if (statesNotContainedInAnyMec.get(element.getColumn())) { + // If the target state is not contained in an MEC, we can copy over the entry. + sspMatrixBuilder.addNextValue(currentChoice, statesNotInMecsBeforeIndex[element.getColumn()], element.getValue()); + } else { + // If the target state is contained in MEC i, we need to add the probability to the corresponding field in the vector + // so that we are able to write the cumulative probability to the MEC into the matrix. + auxiliaryStateToProbabilityMap[stateToMecIndexMap[element.getColumn()]] += element.getValue(); + } + } + + // Now insert all (cumulative) probability values that target an MEC. + for (uint_fast64_t mecIndex = 0; mecIndex < auxiliaryStateToProbabilityMap.size(); ++mecIndex) { + if (auxiliaryStateToProbabilityMap[mecIndex] != 0) { + sspMatrixBuilder.addNextValue(currentChoice, firstAuxiliaryStateIndex + mecIndex, auxiliaryStateToProbabilityMap[mecIndex]); + } + } + } + } + + // Now we are ready to construct the choices for the auxiliary states. + for (uint_fast64_t mecIndex = 0; mecIndex < mecDecomposition.size(); ++mecIndex) { + storm::storage::MaximalEndComponent const& mec = mecDecomposition[mecIndex]; + sspMatrixBuilder.newRowGroup(currentChoice); + + for (auto const& stateChoicesPair : mec) { + uint_fast64_t state = stateChoicesPair.first; + boost::container::flat_set const& choicesInMec = stateChoicesPair.second; + + for (uint_fast64_t choice = nondeterministicChoiceIndices[state]; choice < nondeterministicChoiceIndices[state + 1]; ++choice) { + std::vector auxiliaryStateToProbabilityMap(mecDecomposition.size()); + + // If the choice is not contained in the MEC itself, we have to add a similar distribution to the auxiliary state. + if (choicesInMec.find(choice) == choicesInMec.end()) { + b.push_back(storm::utility::zero()); + + for (auto element : transitionMatrix.getRow(choice)) { + if (statesNotContainedInAnyMec.get(element.getColumn())) { + // If the target state is not contained in an MEC, we can copy over the entry. + sspMatrixBuilder.addNextValue(currentChoice, statesNotInMecsBeforeIndex[element.getColumn()], element.getValue()); + } else { + // If the target state is contained in MEC i, we need to add the probability to the corresponding field in the vector + // so that we are able to write the cumulative probability to the MEC into the matrix. + auxiliaryStateToProbabilityMap[stateToMecIndexMap[element.getColumn()]] += element.getValue(); + } + } + + // Now insert all (cumulative) probability values that target an MEC. + for (uint_fast64_t targetMecIndex = 0; targetMecIndex < auxiliaryStateToProbabilityMap.size(); ++targetMecIndex) { + if (auxiliaryStateToProbabilityMap[targetMecIndex] != 0) { + sspMatrixBuilder.addNextValue(currentChoice, firstAuxiliaryStateIndex + targetMecIndex, auxiliaryStateToProbabilityMap[targetMecIndex]); + } + } + + ++currentChoice; + } + } + } + + // For each auxiliary state, there is the option to achieve the reward value of the LRA associated with the MEC. + ++currentChoice; + b.push_back(lraValuesForEndComponents[mecIndex]); + } + + // Finalize the matrix and solve the corresponding system of equations. + storm::storage::SparseMatrix sspMatrix = sspMatrixBuilder.build(currentChoice); + + std::vector sspResult(numberOfStatesNotInMecs + mecDecomposition.size()); + std::unique_ptr> solver = minMaxLinearEquationSolverFactory.create(sspMatrix); + solver->solveEquationSystem(minimize, sspResult, b); + + // Prepare result vector. + std::vector result(numberOfStates, zero); + + // Set the values for states not contained in MECs. + storm::utility::vector::setVectorValues(result, statesNotContainedInAnyMec, sspResult); + + // Set the values for all states in MECs. + for (auto state : statesInMecs) { + result[state] = sspResult[firstAuxiliaryStateIndex + stateToMecIndexMap[state]]; + } + + return result; + } + + template + ValueType SparseMdpPrctlHelper::computeLraForMaximalEndComponent(bool minimize, storm::storage::SparseMatrix const& transitionMatrix, storm::storage::BitVector const& psiStates, storm::storage::MaximalEndComponent const& mec) { + std::shared_ptr solver = storm::utility::solver::getLpSolver("LRA for MEC"); + solver->setModelSense(minimize ? storm::solver::LpSolver::ModelSense::Maximize : storm::solver::LpSolver::ModelSense::Minimize); + + //// First, we need to create the variables for the problem. + std::map stateToVariableMap; + for (auto const& stateChoicesPair : mec) { + std::string variableName = "h" + std::to_string(stateChoicesPair.first); + stateToVariableMap[stateChoicesPair.first] = solver->addUnboundedContinuousVariable(variableName); + } + storm::expressions::Variable lambda = solver->addUnboundedContinuousVariable("L", 1); + solver->update(); + + // Now we encode the problem as constraints. + for (auto const& stateChoicesPair : mec) { + uint_fast64_t state = stateChoicesPair.first; + + // Now, based on the type of the state, create a suitable constraint. + for (auto choice : stateChoicesPair.second) { + storm::expressions::Expression constraint = -lambda; + ValueType r = 0; + + for (auto element : transitionMatrix.getRow(choice)) { + constraint = constraint + stateToVariableMap.at(element.getColumn()) * solver->getConstant(element.getValue()); + if (psiStates.get(element.getColumn())) { + r += element.getValue(); + } + } + constraint = solver->getConstant(r) + constraint; + + if (minimize) { + constraint = stateToVariableMap.at(state) <= constraint; + } else { + constraint = stateToVariableMap.at(state) >= constraint; + } + solver->addConstraint("state" + std::to_string(state) + "," + std::to_string(choice), constraint); + } + } + + solver->optimize(); + return solver->getContinuousValue(lambda); + } + + template class SparseMdpPrctlHelper; + } + } +} diff --git a/src/modelchecker/prctl/helper/SparseMdpPrctlHelper.h b/src/modelchecker/prctl/helper/SparseMdpPrctlHelper.h new file mode 100644 index 000000000..2c3aa99e4 --- /dev/null +++ b/src/modelchecker/prctl/helper/SparseMdpPrctlHelper.h @@ -0,0 +1,45 @@ +#ifndef STORM_MODELCHECKER_SPARSE_MDP_PRCTL_MODELCHECKER_HELPER_H_ +#define STORM_MODELCHECKER_SPARSE_MDP_PRCTL_MODELCHECKER_HELPER_H_ + +#include + +#include "src/models/sparse/StandardRewardModel.h" + +#include "src/storage/SparseMatrix.h" +#include "src/storage/BitVector.h" +#include "src/storage/MaximalEndComponent.h" + +#include "src/utility/solver.h" + +namespace storm { + namespace modelchecker { + namespace helper { + + template > + class SparseMdpPrctlHelper { + public: + static std::vector computeBoundedUntilProbabilities(bool minimize, storm::storage::SparseMatrix const& transitionMatrix, storm::storage::SparseMatrix const& backwardTransitions, storm::storage::BitVector const& phiStates, storm::storage::BitVector const& psiStates, uint_fast64_t stepBound, storm::utility::solver::MinMaxLinearEquationSolverFactory const& minMaxLinearEquationSolverFactory); + + static std::vector computeNextProbabilities(bool minimize, storm::storage::SparseMatrix const& transitionMatrix, storm::storage::BitVector const& nextStates, storm::utility::solver::MinMaxLinearEquationSolverFactory const& minMaxLinearEquationSolverFactory); + + static std::vector computeUntilProbabilities(bool minimize, storm::storage::SparseMatrix const& transitionMatrix, storm::storage::SparseMatrix const& backwardTransitions, storm::storage::BitVector const& phiStates, storm::storage::BitVector const& psiStates, bool qualitative, storm::utility::solver::MinMaxLinearEquationSolverFactory const& minMaxLinearEquationSolverFactory); + + static std::vector computeUntilProbabilities(bool minimize, storm::storage::SparseMatrix const& transitionMatrix, storm::storage::SparseMatrix const& backwardTransitions, storm::storage::BitVector const& phiStates, storm::storage::BitVector const& psiStates, storm::utility::solver::MinMaxLinearEquationSolverFactory const& MinMaxLinearEquationSolverFactory, bool qualitative); + + static std::vector computeInstantaneousRewards(storm::storage::SparseMatrix const& transitionMatrix, RewardModelType const& rewardModel, bool minimize, uint_fast64_t stepCount, storm::utility::solver::MinMaxLinearEquationSolverFactory const& minMaxLinearEquationSolverFactory); + + static std::vector computeCumulativeRewards(storm::storage::SparseMatrix const& transitionMatrix, RewardModelType const& rewardModel, bool minimize, uint_fast64_t stepBound, storm::utility::solver::MinMaxLinearEquationSolverFactory const& minMaxLinearEquationSolverFactory); + + static std::vector computeReachabilityRewards(bool minimize, storm::storage::SparseMatrix const& transitionMatrix, storm::storage::SparseMatrix const& backwardTransitions, RewardModelType const& rewardModel, storm::storage::BitVector const& targetStates, bool qualitative, storm::utility::solver::MinMaxLinearEquationSolverFactory const& minMaxLinearEquationSolverFactory); + + static std::vector computeLongRunAverage(bool minimize, storm::storage::SparseMatrix const& transitionMatrix, storm::storage::SparseMatrix const& backwardTransitions, storm::storage::BitVector const& psiStates, bool qualitative, storm::utility::solver::MinMaxLinearEquationSolverFactory const& minMaxLinearEquationSolverFactory); + + private: + static ValueType computeLraForMaximalEndComponent(bool minimize, storm::storage::SparseMatrix const& transitionMatrix, storm::storage::BitVector const& goalStates, storm::storage::MaximalEndComponent const& mec); + }; + + } + } +} + +#endif /* STORM_MODELCHECKER_SPARSE_MDP_PRCTL_MODELCHECKER_HELPER_H_ */ \ No newline at end of file diff --git a/src/modelchecker/prctl/helper/SymbolicDtmcPrctlHelper.cpp b/src/modelchecker/prctl/helper/SymbolicDtmcPrctlHelper.cpp new file mode 100644 index 000000000..66c523331 --- /dev/null +++ b/src/modelchecker/prctl/helper/SymbolicDtmcPrctlHelper.cpp @@ -0,0 +1,195 @@ +#include "src/modelchecker/prctl/helper/SymbolicDtmcPrctlHelper.h" + +#include "src/storage/dd/DdType.h" +#include "src/storage/dd/CuddAdd.h" +#include "src/storage/dd/CuddBdd.h" +#include "src/storage/dd/CuddOdd.h" + +#include "src/utility/graph.h" + +#include "src/exceptions/InvalidPropertyException.h" + +namespace storm { + namespace modelchecker { + namespace helper { + + template + storm::dd::Add SymbolicDtmcPrctlHelper::computeUntilProbabilities(storm::models::symbolic::Model const& model, storm::dd::Add const& transitionMatrix, storm::dd::Bdd const& phiStates, storm::dd::Bdd const& psiStates, bool qualitative, storm::utility::solver::SymbolicLinearEquationSolverFactory const& linearEquationSolverFactory) { + // We need to identify the states which have to be taken out of the matrix, i.e. all states that have + // probability 0 and 1 of satisfying the until-formula. + std::pair, storm::dd::Bdd> statesWithProbability01 = storm::utility::graph::performProb01(model, transitionMatrix, phiStates, psiStates); + storm::dd::Bdd maybeStates = !statesWithProbability01.first && !statesWithProbability01.second && model.getReachableStates(); + + // Perform some logging. + STORM_LOG_INFO("Found " << statesWithProbability01.first.getNonZeroCount() << " 'no' states."); + STORM_LOG_INFO("Found " << statesWithProbability01.second.getNonZeroCount() << " 'yes' states."); + STORM_LOG_INFO("Found " << maybeStates.getNonZeroCount() << " 'maybe' states."); + + // Check whether we need to compute exact probabilities for some states. + if (qualitative) { + // Set the values for all maybe-states to 0.5 to indicate that their probability values are neither 0 nor 1. + return statesWithProbability01.second.toAdd() + maybeStates.toAdd() * model.getManager().getConstant(0.5); + } else { + // If there are maybe states, we need to solve an equation system. + if (!maybeStates.isZero()) { + // Create the ODD for the translation between symbolic and explicit storage. + storm::dd::Odd odd(maybeStates); + + // Create the matrix and the vector for the equation system. + storm::dd::Add maybeStatesAdd = maybeStates.toAdd(); + + // Start by cutting away all rows that do not belong to maybe states. Note that this leaves columns targeting + // non-maybe states in the matrix. + storm::dd::Add submatrix = transitionMatrix * maybeStatesAdd; + + // Then compute the vector that contains the one-step probabilities to a state with probability 1 for all + // maybe states. + storm::dd::Add prob1StatesAsColumn = statesWithProbability01.second.toAdd(); + prob1StatesAsColumn = prob1StatesAsColumn.swapVariables(model.getRowColumnMetaVariablePairs()); + storm::dd::Add subvector = submatrix * prob1StatesAsColumn; + subvector = subvector.sumAbstract(model.getColumnVariables()); + + // Finally cut away all columns targeting non-maybe states and convert the matrix into the matrix needed + // for solving the equation system (i.e. compute (I-A)). + submatrix *= maybeStatesAdd.swapVariables(model.getRowColumnMetaVariablePairs()); + submatrix = (model.getRowColumnIdentity() * maybeStatesAdd) - submatrix; + + // Solve the equation system. + std::unique_ptr> solver = linearEquationSolverFactory.create(submatrix, maybeStates, model.getRowVariables(), model.getColumnVariables(), model.getRowColumnMetaVariablePairs()); + storm::dd::Add result = solver->solveEquationSystem(model.getManager().getConstant(0.5) * maybeStatesAdd, subvector); + + return statesWithProbability01.second.toAdd() + result; + } else { + return statesWithProbability01.second.toAdd(); + } + } + } + + template + storm::dd::Add SymbolicDtmcPrctlHelper::computeNextProbabilities(storm::models::symbolic::Model const& model, storm::dd::Add const& transitionMatrix, storm::dd::Bdd const& nextStates) { + storm::dd::Add result = transitionMatrix * nextStates.swapVariables(model.getRowColumnMetaVariablePairs()).toAdd(); + return result.sumAbstract(model.getColumnVariables()); + } + + template + storm::dd::Add SymbolicDtmcPrctlHelper::computeBoundedUntilProbabilities(storm::models::symbolic::Model const& model, storm::dd::Add const& transitionMatrix, storm::dd::Bdd const& phiStates, storm::dd::Bdd const& psiStates, uint_fast64_t stepBound, storm::utility::solver::SymbolicLinearEquationSolverFactory const& linearEquationSolverFactory) { + // We need to identify the states which have to be taken out of the matrix, i.e. all states that have + // probability 0 or 1 of satisfying the until-formula. + storm::dd::Bdd statesWithProbabilityGreater0 = storm::utility::graph::performProbGreater0(model, transitionMatrix.notZero(), phiStates, psiStates, stepBound); + storm::dd::Bdd maybeStates = statesWithProbabilityGreater0 && !psiStates && model.getReachableStates(); + + // If there are maybe states, we need to perform matrix-vector multiplications. + if (!maybeStates.isZero()) { + // Create the ODD for the translation between symbolic and explicit storage. + storm::dd::Odd odd(maybeStates); + + // Create the matrix and the vector for the equation system. + storm::dd::Add maybeStatesAdd = maybeStates.toAdd(); + + // Start by cutting away all rows that do not belong to maybe states. Note that this leaves columns targeting + // non-maybe states in the matrix. + storm::dd::Add submatrix = transitionMatrix * maybeStatesAdd; + + // Then compute the vector that contains the one-step probabilities to a state with probability 1 for all + // maybe states. + storm::dd::Add prob1StatesAsColumn = psiStates.toAdd().swapVariables(model.getRowColumnMetaVariablePairs()); + storm::dd::Add subvector = (submatrix * prob1StatesAsColumn).sumAbstract(model.getColumnVariables()); + + // Finally cut away all columns targeting non-maybe states. + submatrix *= maybeStatesAdd.swapVariables(model.getRowColumnMetaVariablePairs()); + + // Perform the matrix-vector multiplication. + std::unique_ptr> solver = linearEquationSolverFactory.create(submatrix, maybeStates, model.getRowVariables(), model.getColumnVariables(), model.getRowColumnMetaVariablePairs()); + storm::dd::Add result = solver->performMatrixVectorMultiplication(model.getManager().getAddZero(), &subvector, stepBound); + + return psiStates.toAdd() + result; + } else { + return psiStates.toAdd(); + } + } + + template + storm::dd::Add SymbolicDtmcPrctlHelper::computeCumulativeRewards(storm::models::symbolic::Model const& model, storm::dd::Add const& transitionMatrix, uint_fast64_t stepBound, storm::utility::solver::SymbolicLinearEquationSolverFactory const& linearEquationSolverFactory) { + // Only compute the result if the model has at least one reward this->getModel(). + STORM_LOG_THROW(model.hasStateRewards() || model.hasTransitionRewards(), storm::exceptions::InvalidPropertyException, "Missing reward model for formula. Skipping formula."); + + // Compute the reward vector to add in each step based on the available reward models. + storm::dd::Add totalRewardVector = model.hasStateRewards() ? model.getStateRewardVector() : model.getManager().getAddZero(); + if (model.hasTransitionRewards()) { + totalRewardVector += (transitionMatrix * model.getTransitionRewardMatrix()).sumAbstract(model.getColumnVariables()); + } + + // Perform the matrix-vector multiplication. + std::unique_ptr> solver = linearEquationSolverFactory.create(transitionMatrix, model.getReachableStates(), model.getRowVariables(), model.getColumnVariables(), model.getRowColumnMetaVariablePairs()); + return solver->performMatrixVectorMultiplication(model.getManager().getAddZero(), &totalRewardVector, stepBound); + } + + template + storm::dd::Add SymbolicDtmcPrctlHelper::computeInstantaneousRewards(storm::models::symbolic::Model const& model, storm::dd::Add const& transitionMatrix, uint_fast64_t stepBound, storm::utility::solver::SymbolicLinearEquationSolverFactory const& linearEquationSolverFactory) { + // Only compute the result if the model has at least one reward this->getModel(). + STORM_LOG_THROW(model.hasStateRewards(), storm::exceptions::InvalidPropertyException, "Missing reward model for formula. Skipping formula."); + + // Perform the matrix-vector multiplication. + std::unique_ptr> solver = linearEquationSolverFactory.create(transitionMatrix, model.getReachableStates(), model.getRowVariables(), model.getColumnVariables(), model.getRowColumnMetaVariablePairs()); + return solver->performMatrixVectorMultiplication(model.getStateRewardVector(), nullptr, stepBound); + } + + template + storm::dd::Add SymbolicDtmcPrctlHelper::computeReachabilityRewards(storm::models::symbolic::Model const& model, storm::dd::Add const& transitionMatrix, boost::optional> const& stateRewardVector, boost::optional> const& transitionRewardMatrix, storm::dd::Bdd const& targetStates, bool qualitative, storm::utility::solver::SymbolicLinearEquationSolverFactory const& linearEquationSolverFactory) { + + // Only compute the result if there is at least one reward model. + STORM_LOG_THROW(stateRewardVector || transitionRewardMatrix, storm::exceptions::InvalidPropertyException, "Missing reward model for formula. Skipping formula."); + + // Determine which states have a reward of infinity by definition. + storm::dd::Bdd infinityStates = storm::utility::graph::performProb1(model, transitionMatrix.notZero(), model.getReachableStates(), targetStates); + infinityStates = !infinityStates && model.getReachableStates(); + storm::dd::Bdd maybeStates = (!targetStates && !infinityStates) && model.getReachableStates(); + STORM_LOG_INFO("Found " << infinityStates.getNonZeroCount() << " 'infinity' states."); + STORM_LOG_INFO("Found " << targetStates.getNonZeroCount() << " 'target' states."); + STORM_LOG_INFO("Found " << maybeStates.getNonZeroCount() << " 'maybe' states."); + + // Check whether we need to compute exact rewards for some states. + if (qualitative) { + // Set the values for all maybe-states to 1 to indicate that their reward values + // are neither 0 nor infinity. + return infinityStates.toAdd() * model.getManager().getConstant(storm::utility::infinity()) + maybeStates.toAdd() * model.getManager().getConstant(storm::utility::one()); + } else { + // If there are maybe states, we need to solve an equation system. + if (!maybeStates.isZero()) { + // Create the ODD for the translation between symbolic and explicit storage. + storm::dd::Odd odd(maybeStates); + + // Create the matrix and the vector for the equation system. + storm::dd::Add maybeStatesAdd = maybeStates.toAdd(); + + // Start by cutting away all rows that do not belong to maybe states. Note that this leaves columns targeting + // non-maybe states in the matrix. + storm::dd::Add submatrix = transitionMatrix * maybeStatesAdd; + + // Then compute the state reward vector to use in the computation. + storm::dd::Add subvector = stateRewardVector ? maybeStatesAdd * stateRewardVector.get() : model.getManager().getAddZero(); + if (transitionRewardMatrix) { + subvector += (submatrix * transitionRewardMatrix.get()).sumAbstract(model.getColumnVariables()); + } + + // Finally cut away all columns targeting non-maybe states and convert the matrix into the matrix needed + // for solving the equation system (i.e. compute (I-A)). + submatrix *= maybeStatesAdd.swapVariables(model.getRowColumnMetaVariablePairs()); + submatrix = (model.getRowColumnIdentity() * maybeStatesAdd) - submatrix; + + // Solve the equation system. + std::unique_ptr> solver = linearEquationSolverFactory.create(submatrix, maybeStates, model.getRowVariables(), model.getColumnVariables(), model.getRowColumnMetaVariablePairs()); + storm::dd::Add result = solver->solveEquationSystem(model.getManager().getConstant(0.5) * maybeStatesAdd, subvector); + + return infinityStates.toAdd() * model.getManager().getConstant(storm::utility::infinity()) + result; + } else { + return infinityStates.toAdd() * model.getManager().getConstant(storm::utility::infinity()); + } + } + } + + template class SymbolicDtmcPrctlHelper; + + } + } +} \ No newline at end of file diff --git a/src/modelchecker/prctl/helper/SymbolicDtmcPrctlHelper.h b/src/modelchecker/prctl/helper/SymbolicDtmcPrctlHelper.h new file mode 100644 index 000000000..5f1519ab0 --- /dev/null +++ b/src/modelchecker/prctl/helper/SymbolicDtmcPrctlHelper.h @@ -0,0 +1,35 @@ +#ifndef STORM_MODELCHECKER_SPARSE_DTMC_PRCTL_MODELCHECKER_HELPER_H_ +#define STORM_MODELCHECKER_SPARSE_DTMC_PRCTL_MODELCHECKER_HELPER_H_ + +#include "src/models/symbolic/Model.h" + +#include "src/storage/dd/Add.h" +#include "src/storage/dd/Bdd.h" + +#include "src/utility/solver.h" + +namespace storm { + namespace modelchecker { + namespace helper { + + template + class SymbolicDtmcPrctlHelper { + public: + static storm::dd::Add computeBoundedUntilProbabilities(storm::models::symbolic::Model const& model, storm::dd::Add const& transitionMatrix, storm::dd::Bdd const& phiStates, storm::dd::Bdd const& psiStates, uint_fast64_t stepBound, storm::utility::solver::SymbolicLinearEquationSolverFactory const& linearEquationSolverFactory); + + static storm::dd::Add computeNextProbabilities(storm::models::symbolic::Model const& model, storm::dd::Add const& transitionMatrix, storm::dd::Bdd const& nextStates); + + static storm::dd::Add computeUntilProbabilities(storm::models::symbolic::Model const& model, storm::dd::Add const& transitionMatrix, storm::dd::Bdd const& phiStates, storm::dd::Bdd const& psiStates, bool qualitative, storm::utility::solver::SymbolicLinearEquationSolverFactory const& linearEquationSolverFactory); + + static storm::dd::Add computeCumulativeRewards(storm::models::symbolic::Model const& model, storm::dd::Add const& transitionMatrix, uint_fast64_t stepBound, storm::utility::solver::SymbolicLinearEquationSolverFactory const& linearEquationSolverFactory); + + static storm::dd::Add computeInstantaneousRewards(storm::models::symbolic::Model const& model, storm::dd::Add const& transitionMatrix, uint_fast64_t stepBound, storm::utility::solver::SymbolicLinearEquationSolverFactory const& linearEquationSolverFactory); + + static storm::dd::Add computeReachabilityRewards(storm::models::symbolic::Model const& model, storm::dd::Add const& transitionMatrix, boost::optional> const& stateRewardVector, boost::optional> const& transitionRewardMatrix, storm::dd::Bdd const& targetStates, bool qualitative, storm::utility::solver::SymbolicLinearEquationSolverFactory const& linearEquationSolverFactory); + }; + + } + } +} + +#endif /* STORM_MODELCHECKER_SPARSE_DTMC_PRCTL_MODELCHECKER_HELPER_H_ */ \ No newline at end of file diff --git a/src/modelchecker/prctl/helper/SymbolicMdpPrctlHelper.cpp b/src/modelchecker/prctl/helper/SymbolicMdpPrctlHelper.cpp new file mode 100644 index 000000000..e69de29bb diff --git a/src/modelchecker/prctl/helper/SymbolicMdpPrctlHelper.h b/src/modelchecker/prctl/helper/SymbolicMdpPrctlHelper.h new file mode 100644 index 000000000..e69de29bb diff --git a/src/modelchecker/propositional/SparsePropositionalModelChecker.h b/src/modelchecker/propositional/SparsePropositionalModelChecker.h index 5270f4b37..f24ece9ab 100644 --- a/src/modelchecker/propositional/SparsePropositionalModelChecker.h +++ b/src/modelchecker/propositional/SparsePropositionalModelChecker.h @@ -11,6 +11,9 @@ namespace storm { template class SparsePropositionalModelChecker : public AbstractModelChecker { public: + typedef typename SparseModelType::ValueType ValueType; + typedef typename SparseModelType::RewardModelType RewardModelType; + explicit SparsePropositionalModelChecker(SparseModelType const& model); // The implemented methods of the AbstractModelChecker interface. diff --git a/src/models/sparse/Ctmc.h b/src/models/sparse/Ctmc.h index 0712cdb75..728c645ae 100644 --- a/src/models/sparse/Ctmc.h +++ b/src/models/sparse/Ctmc.h @@ -1,9 +1,6 @@ #ifndef STORM_MODELS_SPARSE_CTMC_H_ #define STORM_MODELS_SPARSE_CTMC_H_ -#include -#include - #include "src/models/sparse/DeterministicModel.h" #include "src/utility/OsDetection.h" diff --git a/src/settings/ArgumentBuilder.h b/src/settings/ArgumentBuilder.h index 18f03334d..8239286c5 100644 --- a/src/settings/ArgumentBuilder.h +++ b/src/settings/ArgumentBuilder.h @@ -190,7 +190,7 @@ return *this; \ ArgumentBuilder(ArgumentType type, std::string const& name, std::string const& description) : hasBeenBuilt(false), type(type), name(name), description(description), isOptional(false), hasDefaultValue(false), defaultValue_String(), defaultValue_Integer(), defaultValue_UnsignedInteger(), defaultValue_Double(), defaultValue_Boolean(), userValidationFunctions_String(), userValidationFunctions_Integer(), userValidationFunctions_UnsignedInteger(), userValidationFunctions_Double(), userValidationFunctions_Boolean() { // Intentionally left empty. } - + // A flag that stores whether an argument has been built using this builder. bool hasBeenBuilt; diff --git a/src/storage/MaximalEndComponentDecomposition.cpp b/src/storage/MaximalEndComponentDecomposition.cpp index f32db1b00..cc809e03f 100644 --- a/src/storage/MaximalEndComponentDecomposition.cpp +++ b/src/storage/MaximalEndComponentDecomposition.cpp @@ -14,12 +14,17 @@ namespace storm { template MaximalEndComponentDecomposition::MaximalEndComponentDecomposition(storm::models::sparse::NondeterministicModel const& model) { - performMaximalEndComponentDecomposition(model, storm::storage::BitVector(model.getNumberOfStates(), true)); + performMaximalEndComponentDecomposition(model.getTransitionMatrix(), model.getBackwardTransitions(), storm::storage::BitVector(model.getNumberOfStates(), true)); + } + + template + MaximalEndComponentDecomposition::MaximalEndComponentDecomposition(storm::storage::SparseMatrix const& transitionMatrix, storm::storage::SparseMatrix const& backwardTransitions) { + performMaximalEndComponentDecomposition(transitionMatrix, backwardTransitions, storm::storage::BitVector(transitionMatrix.getRowGroupCount(), true)); } template MaximalEndComponentDecomposition::MaximalEndComponentDecomposition(storm::models::sparse::NondeterministicModel const& model, storm::storage::BitVector const& subsystem) { - performMaximalEndComponentDecomposition(model, subsystem); + performMaximalEndComponentDecomposition(model.getTransitionMatrix(), model.getBackwardTransitions(), subsystem); } template @@ -45,30 +50,28 @@ namespace storm { } template - void MaximalEndComponentDecomposition::performMaximalEndComponentDecomposition(storm::models::sparse::NondeterministicModel const& model, storm::storage::BitVector const& subsystem) { - // Get some references for convenient access. - storm::storage::SparseMatrix backwardTransitions = model.getBackwardTransitions(); - std::vector const& nondeterministicChoiceIndices = model.getNondeterministicChoiceIndices(); - storm::storage::SparseMatrix const& transitionMatrix = model.getTransitionMatrix(); + void MaximalEndComponentDecomposition::performMaximalEndComponentDecomposition(storm::storage::SparseMatrix const& transitionMatrix, storm::storage::SparseMatrix backwardTransitions, storm::storage::BitVector const& subsystem) { + // Get some data for convenient access. + uint_fast64_t numberOfStates = transitionMatrix.getRowGroupCount(); + std::vector const& nondeterministicChoiceIndices = transitionMatrix.getRowGroupIndices(); // Initialize the maximal end component list to be the full state space. std::list endComponentStateSets; endComponentStateSets.emplace_back(subsystem.begin(), subsystem.end()); - storm::storage::BitVector statesToCheck(model.getNumberOfStates()); + storm::storage::BitVector statesToCheck(numberOfStates); // The iterator used here should really be a const_iterator. // However, gcc 4.8 (and assorted libraries) does not provide an erase(const_iterator) method for std::list // but only an erase(iterator). This is in compliance with the c++11 draft N3337, which specifies the change // from iterator to const_iterator only for "set, multiset, map [and] multimap". - // FIXME: As soon as gcc provides an erase(const_iterator) method, change this iterator back to a const_iterator. - for (std::list::iterator mecIterator = endComponentStateSets.begin(); mecIterator != endComponentStateSets.end();) { + for (std::list::const_iterator mecIterator = endComponentStateSets.begin(); mecIterator != endComponentStateSets.end();) { StateBlock const& mec = *mecIterator; // Keep track of whether the MEC changed during this iteration. bool mecChanged = false; // Get an SCC decomposition of the current MEC candidate. - StronglyConnectedComponentDecomposition sccs(model, mec, true); + StronglyConnectedComponentDecomposition sccs(transitionMatrix, mec, true); // We need to do another iteration in case we have either more than once SCC or the SCC is smaller than // the MEC canditate itself. @@ -79,7 +82,7 @@ namespace storm { statesToCheck.set(scc.begin(), scc.end()); while (!statesToCheck.empty()) { - storm::storage::BitVector statesToRemove(model.getNumberOfStates()); + storm::storage::BitVector statesToRemove(numberOfStates); for (auto state : statesToCheck) { bool keepStateInMEC = false; @@ -132,7 +135,7 @@ namespace storm { } } - std::list::iterator eraseIterator(mecIterator); + std::list::const_iterator eraseIterator(mecIterator); ++mecIterator; endComponentStateSets.erase(eraseIterator); } else { diff --git a/src/storage/MaximalEndComponentDecomposition.h b/src/storage/MaximalEndComponentDecomposition.h index d2dfa077a..672163e52 100644 --- a/src/storage/MaximalEndComponentDecomposition.h +++ b/src/storage/MaximalEndComponentDecomposition.h @@ -26,6 +26,14 @@ namespace storm { */ MaximalEndComponentDecomposition(storm::models::sparse::NondeterministicModel const& model); + /* + * Creates an MEC decomposition of the given model (represented by a row-grouped matrix). + * + * @param transitionMatrix The transition relation of model to decompose into MECs. + * @param backwardTransition The reversed transition relation. + */ + MaximalEndComponentDecomposition(storm::storage::SparseMatrix const& transitionMatrix, storm::storage::SparseMatrix const& backwardTransitions); + /*! * Creates an MEC decomposition of the given subsystem in the given model. * @@ -67,10 +75,11 @@ namespace storm { * Performs the actual decomposition of the given subsystem in the given model into MECs. As a side-effect * this stores the MECs found in the current decomposition. * - * @param model The model whose subsystem to decompose into MECs. + * @param transitionMatrix The transition matrix representing the system whose subsystem to decompose into MECs. + * @param backwardTransitions The reversed transition relation. * @param subsystem The subsystem to decompose. */ - void performMaximalEndComponentDecomposition(storm::models::sparse::NondeterministicModel const& model, storm::storage::BitVector const& subsystem); + void performMaximalEndComponentDecomposition(storm::storage::SparseMatrix const& transitionMatrix, storm::storage::SparseMatrix backwardTransitions, storm::storage::BitVector const& subsystem); }; } } diff --git a/src/storage/StronglyConnectedComponentDecomposition.cpp b/src/storage/StronglyConnectedComponentDecomposition.cpp index 9d032d709..a03a22b63 100644 --- a/src/storage/StronglyConnectedComponentDecomposition.cpp +++ b/src/storage/StronglyConnectedComponentDecomposition.cpp @@ -25,6 +25,18 @@ namespace storm { performSccDecomposition(model.getTransitionMatrix(), subsystem, dropNaiveSccs, onlyBottomSccs); } + template + StronglyConnectedComponentDecomposition::StronglyConnectedComponentDecomposition(storm::storage::SparseMatrix const& transitionMatrix, StateBlock const& block, bool dropNaiveSccs, bool onlyBottomSccs) { + storm::storage::BitVector subsystem(transitionMatrix.getRowGroupCount(), block.begin(), block.end()); + performSccDecomposition(transitionMatrix, subsystem, dropNaiveSccs, onlyBottomSccs); + } + + + template + StronglyConnectedComponentDecomposition::StronglyConnectedComponentDecomposition(storm::storage::SparseMatrix const& transitionMatrix, bool dropNaiveSccs, bool onlyBottomSccs) { + performSccDecomposition(transitionMatrix, storm::storage::BitVector(transitionMatrix.getRowGroupCount(), true), dropNaiveSccs, onlyBottomSccs); + } + template StronglyConnectedComponentDecomposition::StronglyConnectedComponentDecomposition(storm::storage::SparseMatrix const& transitionMatrix, storm::storage::BitVector const& subsystem, bool dropNaiveSccs, bool onlyBottomSccs) { performSccDecomposition(transitionMatrix, subsystem, dropNaiveSccs, onlyBottomSccs); diff --git a/src/storage/StronglyConnectedComponentDecomposition.h b/src/storage/StronglyConnectedComponentDecomposition.h index 9f1e3d3b1..052c35cbe 100644 --- a/src/storage/StronglyConnectedComponentDecomposition.h +++ b/src/storage/StronglyConnectedComponentDecomposition.h @@ -62,6 +62,30 @@ namespace storm { */ StronglyConnectedComponentDecomposition(storm::models::sparse::Model const& model, storm::storage::BitVector const& subsystem, bool dropNaiveSccs = false, bool onlyBottomSccs = false); + /* + * Creates an SCC decomposition of the given subsystem in the given system (whose transition relation is + * given by a sparse matrix). + * + * @param transitionMatrix The transition matrix of the system to decompose. + * @param block The block to decompose into SCCs. + * @param dropNaiveSccs A flag that indicates whether trivial SCCs (i.e. SCCs consisting of just one state + * without a self-loop) are to be kept in the decomposition. + * @param onlyBottomSccs If set to true, only bottom SCCs, i.e. SCCs in which all states have no way of + * leaving the SCC), are kept. + */ + StronglyConnectedComponentDecomposition(storm::storage::SparseMatrix const& transitionMatrix, StateBlock const& block, bool dropNaiveSccs = false, bool onlyBottomSccs = false); + + /* + * Creates an SCC decomposition of the given system (whose transition relation is given by a sparse matrix). + * + * @param transitionMatrix The transition matrix of the system to decompose. + * @param dropNaiveSccs A flag that indicates whether trivial SCCs (i.e. SCCs consisting of just one state + * without a self-loop) are to be kept in the decomposition. + * @param onlyBottomSccs If set to true, only bottom SCCs, i.e. SCCs in which all states have no way of + * leaving the SCC), are kept. + */ + StronglyConnectedComponentDecomposition(storm::storage::SparseMatrix const& transitionMatrix, bool dropNaiveSccs = false, bool onlyBottomSccs = false); + /* * Creates an SCC decomposition of the given subsystem in the given system (whose transition relation is * given by a sparse matrix).