From b096180de808916f8dc2a7c845a01b8b015875da Mon Sep 17 00:00:00 2001 From: David_Korzeniewski Date: Fri, 20 Mar 2015 17:15:49 +0100 Subject: [PATCH] LRA on DTMCs implemented Former-commit-id: 633d81323d15dc040410cd73b682d4072bbf65fe --- .../prctl/SparseDtmcPrctlModelChecker.cpp | 129 ++++++++++++++++++ .../prctl/SparseDtmcPrctlModelChecker.h | 7 +- 2 files changed, 135 insertions(+), 1 deletion(-) diff --git a/src/modelchecker/prctl/SparseDtmcPrctlModelChecker.cpp b/src/modelchecker/prctl/SparseDtmcPrctlModelChecker.cpp index 9d12d6166..c7aa308e0 100644 --- a/src/modelchecker/prctl/SparseDtmcPrctlModelChecker.cpp +++ b/src/modelchecker/prctl/SparseDtmcPrctlModelChecker.cpp @@ -1,6 +1,7 @@ #include "src/modelchecker/prctl/SparseDtmcPrctlModelChecker.h" #include +#include #include "src/utility/macros.h" #include "src/utility/vector.h" @@ -9,6 +10,8 @@ #include "src/modelchecker/results/ExplicitQualitativeCheckResult.h" #include "src/modelchecker/results/ExplicitQuantitativeCheckResult.h" +#include "src/storage/StronglyConnectedComponentDecomposition.h" + #include "src/exceptions/InvalidPropertyException.h" namespace storm { @@ -297,6 +300,132 @@ namespace storm { return std::unique_ptr(new ExplicitQuantitativeCheckResult(this->computeReachabilityRewardsHelper(subResult.getTruthValuesVector(), qualitative))); } + template + std::vector SparseDtmcPrctlModelChecker::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 DTMC into its BSCCs. + storm::storage::StronglyConnectedComponentDecomposition bsccDecomposition(this->getModelAs>(), false, true); + + // Get some data members for convenience. + typename storm::storage::SparseMatrix const& transitionMatrix = this->getModel().getTransitionMatrix(); + + // Now start with compute the long-run average for all BSCCs in isolation. + std::vector lraValuesForBsccs; + + // While doing so, we already gather some information for the following steps. + std::vector stateToBsccIndexMap(numOfStates); + storm::storage::BitVector statesInBsccs(numOfStates); + + for (uint_fast64_t currentBsccIndex = 0; currentBsccIndex < bsccDecomposition.size(); ++currentBsccIndex) { + storm::storage::StronglyConnectedComponent const& bscc = bsccDecomposition[currentBsccIndex]; + + // Gather information for later use. + for (auto const& state : bscc) { + statesInBsccs.set(state); + stateToBsccIndexMap[state] = currentBsccIndex; + } + + // Compute the LRA value for the current BSCC + lraValuesForBsccs.push_back(this->computeLraForBSCC(transitionMatrix, psiStates, bscc)); + } + + // For fast transition rewriting, we build some auxiliary data structures. + storm::storage::BitVector statesNotContainedInAnyBscc = ~statesInBsccs; + + // Prepare result vector. + std::vector result(numOfStates); + + //Set the values for all states in BSCCs. + for (auto state : statesInBsccs) { + result[state] = lraValuesForBsccs[stateToBsccIndexMap[state]]; + } + + //for all states not in any bscc set the result to the minimal/maximal value of the reachable BSCCs + //there might be a more efficient way to do this... + for (auto state : statesNotContainedInAnyBscc){ + + //calculate what result values the reachable states in BSCCs have + storm::storage::BitVector currentState(numOfStates); + currentState.set(state); + storm::storage::BitVector reachableStates = storm::utility::graph::getReachableStates( + transitionMatrix, currentState, storm::storage::BitVector(numOfStates, true), statesInBsccs + ); + + storm::storage::BitVector reachableBsccStates = statesInBsccs & reachableStates; + std::vector reachableResults(reachableBsccStates.getNumberOfSetBits()); + storm::utility::vector::selectVectorValues(reachableResults, reachableBsccStates, result); + + //now select the minimal/maximal element + if (minimize){ + result[state] = *std::min_element(reachableResults.begin(), reachableResults.end()); + } else { + result[state] = *std::max_element(reachableResults.begin(), reachableResults.end()); + } + } + + return result; + } + + template + std::unique_ptr SparseDtmcPrctlModelChecker::computeLongRunAverage(storm::logic::StateFormula const& stateFormula, 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(stateFormula); + ExplicitQualitativeCheckResult const& subResult = subResultPointer->asExplicitQualitativeCheckResult(); + + return std::unique_ptr(new ExplicitQuantitativeCheckResult(this->computeLongRunAverageHelper(optimalityType.get() == storm::logic::OptimalityType::Minimize, subResult.getTruthValuesVector(), qualitative))); + } + + + template + ValueType SparseDtmcPrctlModelChecker::computeLraForBSCC(storm::storage::SparseMatrix const& transitionMatrix, storm::storage::BitVector const& psiStates, storm::storage::StronglyConnectedComponent const& bscc) { + //if size is one we can immediately derive the result + if (bscc.size() == 1){ + if (psiStates.get(*(bscc.begin()))) { + return storm::utility::one(); + } else{ + return storm::utility::zero(); + } + } + std::unique_ptr> solver = storm::utility::solver::getLinearEquationSolver(); + + storm::storage::BitVector subsystem = storm::storage::BitVector(transitionMatrix.getRowCount()); + subsystem.set(bscc.begin(), bscc.end()); + + storm::storage::SparseMatrix subsystemMatrix = transitionMatrix.getSubmatrix(false, subsystem, subsystem); + std::vector steadyStateDist(subsystemMatrix.getRowCount(), storm::utility::zero()); + steadyStateDist[0] = storm::utility::one(); + + solver->solveEquationSystem(subsystemMatrix, steadyStateDist, steadyStateDist); + + ValueType length = storm::utility::zero(); + for (auto value : steadyStateDist) { + length += value; + } + storm::utility::vector::scaleVectorInPlace(steadyStateDist, storm::utility::one() / length); + + std::vector steadyStateForPsiStates(transitionMatrix.getRowCount(), storm::utility::zero()); + storm::utility::vector::setVectorValues(steadyStateForPsiStates, psiStates & subsystem, steadyStateDist); + + ValueType result = storm::utility::zero(); + + for (auto value : steadyStateForPsiStates) { + result += value; + } + + return result; + } + template storm::models::Dtmc const& SparseDtmcPrctlModelChecker::getModel() const { return this->template getModelAs>(); diff --git a/src/modelchecker/prctl/SparseDtmcPrctlModelChecker.h b/src/modelchecker/prctl/SparseDtmcPrctlModelChecker.h index f7557d63b..d1e880e1a 100644 --- a/src/modelchecker/prctl/SparseDtmcPrctlModelChecker.h +++ b/src/modelchecker/prctl/SparseDtmcPrctlModelChecker.h @@ -5,6 +5,7 @@ #include "src/models/Dtmc.h" #include "src/utility/solver.h" #include "src/solver/LinearEquationSolver.h" +#include "src/storage/StronglyConnectedComponent.h" namespace storm { namespace modelchecker { @@ -23,6 +24,7 @@ namespace storm { 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 computeLongRunAverage(storm::logic::StateFormula const& stateFormula, bool qualitative = false, boost::optional const& optimalityType = boost::optional()) override; protected: storm::models::Dtmc const& getModel() const override; @@ -35,7 +37,10 @@ namespace storm { std::vector computeInstantaneousRewardsHelper(uint_fast64_t stepCount) const; std::vector computeCumulativeRewardsHelper(uint_fast64_t stepBound) const; std::vector computeReachabilityRewardsHelper(storm::storage::BitVector const& targetStates, bool qualitative) const; - + std::vector computeLongRunAverageHelper(bool minimize, storm::storage::BitVector const& psiStates, bool qualitative) const; + + static ValueType computeLraForBSCC(storm::storage::SparseMatrix const& transitionMatrix, storm::storage::BitVector const& goalStates, storm::storage::StronglyConnectedComponent const& bscc); + // An object that is used for solving linear equations and performing matrix-vector multiplication. std::unique_ptr> linearEquationSolver; };