From c1917ce6d9ee1f262a3672806cab1bf5f984cdf6 Mon Sep 17 00:00:00 2001 From: dehnert Date: Tue, 7 Apr 2015 17:28:51 +0200 Subject: [PATCH] Finalized hybrid DTMC model checker. It now passes its tests. Former-commit-id: 99d79e1bc6250695fe057692128d4ab51bf9ad2f --- src/builder/DdPrismModelBuilder.cpp | 26 ++- src/logic/BoundedUntilFormula.cpp | 14 +- src/logic/BoundedUntilFormula.h | 5 +- .../csl/SparseCtmcCslModelChecker.cpp | 4 +- .../prctl/HybridDtmcPrctlModelChecker.cpp | 212 +++++++++++++++++- .../prctl/HybridDtmcPrctlModelChecker.h | 14 +- .../prctl/SparseDtmcPrctlModelChecker.cpp | 38 ++-- .../prctl/SparseMdpPrctlModelChecker.cpp | 31 ++- .../results/HybridQuantitativeCheckResult.cpp | 26 +++ .../results/HybridQuantitativeCheckResult.h | 4 + .../SymbolicQuantitativeCheckResult.cpp | 18 +- .../results/SymbolicQuantitativeCheckResult.h | 7 + src/solver/NativeLinearEquationSolver.cpp | 2 - src/storage/dd/CuddAdd.cpp | 2 +- src/storage/dd/CuddAdd.h | 2 +- src/storage/dd/CuddBdd.cpp | 2 +- src/storage/dd/CuddBdd.h | 2 +- src/utility/graph.h | 44 +++- test/functional/builder/die.pm | 1 + .../GmmxxHybridDtmcPrctlModelCheckerTest.cpp | 165 ++++++++++++++ .../SparseDtmcPrctlModelCheckerTest.cpp | 129 +++++++++++ 21 files changed, 665 insertions(+), 83 deletions(-) create mode 100644 test/functional/modelchecker/GmmxxHybridDtmcPrctlModelCheckerTest.cpp create mode 100644 test/functional/modelchecker/SparseDtmcPrctlModelCheckerTest.cpp diff --git a/src/builder/DdPrismModelBuilder.cpp b/src/builder/DdPrismModelBuilder.cpp index 329e1dfe3..6378a3c56 100644 --- a/src/builder/DdPrismModelBuilder.cpp +++ b/src/builder/DdPrismModelBuilder.cpp @@ -546,21 +546,26 @@ namespace storm { storm::dd::Add states = generationInfo.rowExpressionAdapter->translateExpression(transitionReward.getStatePredicateExpression()); storm::dd::Add rewards = generationInfo.rowExpressionAdapter->translateExpression(transitionReward.getRewardValueExpression()); - storm::dd::Add synchronization; + storm::dd::Add synchronization = generationInfo.manager->getAddOne(); + storm::dd::Add transitions; if (transitionReward.isLabeled()) { - synchronization = getSynchronizationDecisionDiagram(generationInfo, transitionReward.getActionIndex()); + if (generationInfo.program.getModelType() == storm::prism::Program::ModelType::MDP) { + synchronization = getSynchronizationDecisionDiagram(generationInfo, transitionReward.getActionIndex()); + } transitions = globalModule.synchronizingActionToDecisionDiagramMap.at(transitionReward.getActionIndex()).transitionsDd; } else { - synchronization = getSynchronizationDecisionDiagram(generationInfo); + if (generationInfo.program.getModelType() == storm::prism::Program::ModelType::MDP) { + synchronization = getSynchronizationDecisionDiagram(generationInfo); + } transitions = globalModule.independentAction.transitionsDd; } storm::dd::Add transitionRewardDd = synchronization * states * rewards; if (generationInfo.program.getModelType() == storm::prism::Program::ModelType::MDP) { - transitionRewardDd += transitions.notZero().toAdd() * transitionRewardDd; + transitionRewardDd = transitions.notZero().toAdd() * transitionRewardDd; } else { - transitionRewardDd += transitions * transitionRewardDd; + transitionRewardDd = transitions * transitionRewardDd; } // Perform some sanity checks. @@ -575,7 +580,7 @@ namespace storm { if (generationInfo.program.getModelType() == storm::prism::Program::ModelType::DTMC) { transitionRewards /= fullTransitionMatrix; } - + return std::make_pair(stateRewards, transitionRewards); } @@ -636,14 +641,19 @@ namespace storm { stateAndTransitionRewards = createRewardDecisionDiagrams(generationInfo, rewardModel, globalModule, transitionMatrix); } - // Cut the transition matrix to the reachable fragment of the state space. + // Cut the transitions and rewards to the reachable fragment of the state space. storm::dd::Bdd initialStates = createInitialStatesDecisionDiagram(generationInfo); storm::dd::Bdd transitionMatrixBdd = transitionMatrix.notZero(); if (program.getModelType() == storm::prism::Program::ModelType::MDP) { transitionMatrixBdd = transitionMatrixBdd.existsAbstract(generationInfo.allNondeterminismVariables); } storm::dd::Bdd reachableStates = computeReachableStates(generationInfo, initialStates, transitionMatrixBdd); - transitionMatrix *= reachableStates.toAdd(); + storm::dd::Add reachableStatesAdd = reachableStates.toAdd(); + transitionMatrix *= reachableStatesAdd; + if (stateAndTransitionRewards) { + stateAndTransitionRewards.get().first *= reachableStatesAdd; + stateAndTransitionRewards.get().second *= reachableStatesAdd; + } // Detect deadlocks and 1) fix them if requested 2) throw an error otherwise. storm::dd::Bdd statesWithTransition = transitionMatrixBdd.existsAbstract(generationInfo.columnMetaVariables); diff --git a/src/logic/BoundedUntilFormula.cpp b/src/logic/BoundedUntilFormula.cpp index 5c148bc1f..be9ea4f9b 100644 --- a/src/logic/BoundedUntilFormula.cpp +++ b/src/logic/BoundedUntilFormula.cpp @@ -22,16 +22,12 @@ namespace storm { return true; } - bool BoundedUntilFormula::isIntervalBounded() const { - return bounds.which() == 1; - } - - bool BoundedUntilFormula::isIntegerUpperBounded() const { + bool BoundedUntilFormula::hasDiscreteTimeBound() const { return bounds.which() == 0; } bool BoundedUntilFormula::isPctlPathFormula() const { - return this->isIntegerUpperBounded() && this->getLeftSubformula().isPctlStateFormula() && this->getRightSubformula().isPctlStateFormula(); + return this->hasDiscreteTimeBound() && this->getLeftSubformula().isPctlStateFormula() && this->getRightSubformula().isPctlStateFormula(); } bool BoundedUntilFormula::isCslPathFormula() const { @@ -42,7 +38,7 @@ namespace storm { return boost::get>(bounds); } - uint_fast64_t BoundedUntilFormula::getUpperBound() const { + uint_fast64_t BoundedUntilFormula::getDiscreteTimeBound() const { return boost::get(bounds); } @@ -50,11 +46,11 @@ namespace storm { this->getLeftSubformula().writeToStream(out); out << " U"; - if (this->isIntervalBounded()) { + if (!this->hasDiscreteTimeBound()) { std::pair const& intervalBounds = getIntervalBounds(); out << "[" << intervalBounds.first << "," << intervalBounds.second << "] "; } else { - out << "<=" << getUpperBound() << " "; + out << "<=" << getDiscreteTimeBound() << " "; } this->getRightSubformula().writeToStream(out); diff --git a/src/logic/BoundedUntilFormula.h b/src/logic/BoundedUntilFormula.h index 0ed97524a..ffb0f0d32 100644 --- a/src/logic/BoundedUntilFormula.h +++ b/src/logic/BoundedUntilFormula.h @@ -16,11 +16,10 @@ namespace storm { virtual bool containsBoundedUntilFormula() const override; - bool isIntervalBounded() const; - bool isIntegerUpperBounded() const; + bool hasDiscreteTimeBound() const; std::pair const& getIntervalBounds() const; - uint_fast64_t getUpperBound() const; + uint_fast64_t getDiscreteTimeBound() const; virtual bool isPctlPathFormula() const override; virtual bool isCslPathFormula() const override; diff --git a/src/modelchecker/csl/SparseCtmcCslModelChecker.cpp b/src/modelchecker/csl/SparseCtmcCslModelChecker.cpp index 4f24cde22..74540327b 100644 --- a/src/modelchecker/csl/SparseCtmcCslModelChecker.cpp +++ b/src/modelchecker/csl/SparseCtmcCslModelChecker.cpp @@ -41,12 +41,12 @@ namespace storm { ExplicitQualitativeCheckResult const& rightResult = rightResultPointer->asExplicitQualitativeCheckResult(); double lowerBound = 0; double upperBound = 0; - if (pathFormula.isIntervalBounded()) { + if (!pathFormula.hasDiscreteTimeBound()) { std::pair const& intervalBounds = pathFormula.getIntervalBounds(); lowerBound = intervalBounds.first; upperBound = intervalBounds.second; } else { - upperBound = pathFormula.getUpperBound(); + upperBound = pathFormula.getDiscreteTimeBound(); } std::unique_ptr result = std::unique_ptr(new ExplicitQuantitativeCheckResult(this->computeBoundedUntilProbabilitiesHelper(leftResult.getTruthValuesVector(), rightResult.getTruthValuesVector(), this->getModel().getExitRateVector(), qualitative, lowerBound, upperBound))); diff --git a/src/modelchecker/prctl/HybridDtmcPrctlModelChecker.cpp b/src/modelchecker/prctl/HybridDtmcPrctlModelChecker.cpp index 9c959be09..4ddda3339 100644 --- a/src/modelchecker/prctl/HybridDtmcPrctlModelChecker.cpp +++ b/src/modelchecker/prctl/HybridDtmcPrctlModelChecker.cpp @@ -36,9 +36,6 @@ namespace storm { std::pair, storm::dd::Bdd> statesWithProbability01 = storm::utility::graph::performProb01(model, transitionMatrix, phiStates, psiStates); storm::dd::Bdd maybeStates = !statesWithProbability01.first && !statesWithProbability01.second && model.getReachableStates(); - // Create the ODD for the translation between symbolic and explicit storage. - storm::dd::Odd odd(maybeStates); - // Perform some logging. STORM_LOG_INFO("Found " << statesWithProbability01.first.getNonZeroCount() << " 'no' states."); STORM_LOG_INFO("Found " << statesWithProbability01.second.getNonZeroCount() << " 'yes' states."); @@ -51,6 +48,9 @@ namespace storm { } 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(); @@ -70,23 +70,17 @@ namespace storm { submatrix *= maybeStatesAdd.swapVariables(model.getRowColumnMetaVariablePairs()); submatrix = (model.getRowColumnIdentity() * maybeStatesAdd) - submatrix; - submatrix.exportToDot("submatrix.dot"); - // 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. - STORM_LOG_DEBUG("Converting the symbolic matrix to a sparse matrix."); storm::storage::SparseMatrix explicitSubmatrix = submatrix.toMatrix(odd, odd); - - STORM_LOG_DEBUG("Converting the symbolic vector to an explicit vector."); std::vector b = subvector.template toVector(odd); - STORM_LOG_DEBUG("Solving explicit linear equation system."); std::unique_ptr> solver = linearEquationSolverFactory.create(explicitSubmatrix); solver->solveEquationSystem(x, b); - // Now that we have the explicit solution of the system, we need to transform it to a symbolic representation. + // 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())); @@ -103,6 +97,204 @@ namespace storm { return this->computeUntilProbabilitiesHelper(this->getModel(), this->getModel().getTransitionMatrix(), leftResult.getTruthValuesVector(), rightResult.getTruthValuesVector(), qualitative, *this->linearEquationSolverFactory); } + template + std::unique_ptr HybridDtmcPrctlModelChecker::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 HybridDtmcPrctlModelChecker::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()); + } + + template + std::unique_ptr HybridDtmcPrctlModelChecker::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()); + 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); + } + + template + std::unique_ptr HybridDtmcPrctlModelChecker::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::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 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()); + + // Create the solution vector. + std::vector x(maybeStates.getNonZeroCount(), storm::utility::zero()); + + // Translate the symbolic matrix/vector to their explicit representations. + storm::storage::SparseMatrix explicitSubmatrix = submatrix.toMatrix(odd, odd); + std::vector b = subvector.template toVector(odd); + + std::unique_ptr> solver = linearEquationSolverFactory.create(explicitSubmatrix); + solver->performMatrixVectorMultiplication(x, &b, 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 HybridDtmcPrctlModelChecker::computeCumulativeRewards(storm::logic::CumulativeRewardFormula const& rewardPathFormula, 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 HybridDtmcPrctlModelChecker::computeCumulativeRewardsHelper(storm::models::symbolic::Model const& model, storm::dd::Add const& transitionMatrix, uint_fast64_t stepBound, storm::utility::solver::LinearEquationSolverFactory 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(odd, odd); + std::vector b = totalRewardVector.template toVector(odd); + + // Perform the matrix-vector multiplication. + std::unique_ptr> solver = linearEquationSolverFactory.create(explicitMatrix); + solver->performMatrixVectorMultiplication(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 HybridDtmcPrctlModelChecker::computeInstantaneousRewards(storm::logic::InstantaneousRewardFormula const& rewardPathFormula, 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); + } + + template + std::unique_ptr HybridDtmcPrctlModelChecker::computeInstantaneousRewardsHelper(storm::models::symbolic::Model const& model, storm::dd::Add const& transitionMatrix, uint_fast64_t stepBound, storm::utility::solver::LinearEquationSolverFactory 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()); + + // Create the solution vector (and initialize it to the state rewards of the model). + std::vector x = model.getStateRewardVector().template toVector(odd); + + // Translate the symbolic matrix to its explicit representations. + storm::storage::SparseMatrix explicitMatrix = transitionMatrix.toMatrix(odd, odd); + + // Perform the matrix-vector multiplication. + std::unique_ptr> solver = linearEquationSolverFactory.create(explicitMatrix); + solver->performMatrixVectorMultiplication(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 HybridDtmcPrctlModelChecker::computeReachabilityRewards(storm::logic::ReachabilityRewardFormula const& rewardPathFormula, 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(), subResult.getTruthValuesVector(), *this->linearEquationSolverFactory, qualitative); + } + + template + std::unique_ptr HybridDtmcPrctlModelChecker::computeReachabilityRewardsHelper(storm::models::symbolic::Model const& model, storm::dd::Add const& transitionMatrix, storm::dd::Bdd const& targetStates, storm::utility::solver::LinearEquationSolverFactory const& linearEquationSolverFactory, bool qualitative) { + // Only compute the result if the model has at least one reward model. + STORM_LOG_THROW(model.hasStateRewards() || model.hasTransitionRewards(), 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 = model.hasStateRewards() ? maybeStatesAdd * model.getStateRewardVector() : model.getManager().getAddZero(); + if (model.hasTransitionRewards()) { + subvector += (submatrix * model.getTransitionRewardMatrix()).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; + + // Create the solution vector. + std::vector x(maybeStates.getNonZeroCount(), ValueType(0.5)); + + // Translate the symbolic matrix/vector to their explicit representations. + storm::storage::SparseMatrix explicitSubmatrix = submatrix.toMatrix(odd, odd); + std::vector b = subvector.template toVector(odd); + + // Now solve the resulting equation system. + std::unique_ptr> solver = linearEquationSolverFactory.create(explicitSubmatrix); + solver->solveEquationSystem(x, b); + + // 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::Dtmc const& HybridDtmcPrctlModelChecker::getModel() const { return this->template getModelAs>(); diff --git a/src/modelchecker/prctl/HybridDtmcPrctlModelChecker.h b/src/modelchecker/prctl/HybridDtmcPrctlModelChecker.h index bc722f57f..bb7a263a9 100644 --- a/src/modelchecker/prctl/HybridDtmcPrctlModelChecker.h +++ b/src/modelchecker/prctl/HybridDtmcPrctlModelChecker.h @@ -15,15 +15,25 @@ namespace storm { // The implemented methods of the AbstractModelChecker interface. virtual bool canHandle(storm::logic::Formula const& formula) const override; + 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; + 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::LinearEquationSolverFactory 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::LinearEquationSolverFactory 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::LinearEquationSolverFactory 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::LinearEquationSolverFactory const& linearEquationSolverFactory); + static std::unique_ptr computeReachabilityRewardsHelper(storm::models::symbolic::Model const& model, storm::dd::Add const& transitionMatrix, storm::dd::Bdd const& targetStates, storm::utility::solver::LinearEquationSolverFactory 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 46ba68183..a7775139a 100644 --- a/src/modelchecker/prctl/SparseDtmcPrctlModelChecker.cpp +++ b/src/modelchecker/prctl/SparseDtmcPrctlModelChecker.cpp @@ -35,44 +35,41 @@ namespace storm { 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 statesWithProbabilityGreater0 = storm::utility::graph::performProbGreater0(this->getModel().getBackwardTransitions(), phiStates, psiStates, true, stepBound); - STORM_LOG_INFO("Found " << statesWithProbabilityGreater0.getNumberOfSetBits() << " 'maybe' states."); + 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 (!statesWithProbabilityGreater0.empty()) { + 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, statesWithProbabilityGreater0, statesWithProbabilityGreater0, true); - - // Compute the new set of target states in the reduced system. - storm::storage::BitVector rightStatesInReducedSystem = psiStates % statesWithProbabilityGreater0; + storm::storage::SparseMatrix submatrix = this->getModel().getTransitionMatrix().getSubmatrix(true, maybeStates, maybeStates, true); - // Make all rows absorbing that satisfy the second sub-formula. - submatrix.makeRowsAbsorbing(rightStatesInReducedSystem); + // 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(statesWithProbabilityGreater0.getNumberOfSetBits()); - storm::utility::vector::setVectorValues(subresult, rightStatesInReducedSystem, storm::utility::one()); + 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, nullptr, stepBound); + solver->performMatrixVectorMultiplication(subresult, &b, stepBound); // Set the values of the resulting vector accordingly. - storm::utility::vector::setVectorValues(result, statesWithProbabilityGreater0, subresult); - storm::utility::vector::setVectorValues(result, ~statesWithProbabilityGreater0, storm::utility::zero()); + 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& rightResult = rightResultPointer->asExplicitQualitativeCheckResult(); - std::unique_ptr result = std::unique_ptr(new ExplicitQuantitativeCheckResult(this->computeBoundedUntilProbabilitiesHelper(leftResult.getTruthValuesVector(), rightResult.getTruthValuesVector(), pathFormula.getUpperBound()))); - + std::unique_ptr result = std::unique_ptr(new ExplicitQuantitativeCheckResult(this->computeBoundedUntilProbabilitiesHelper(leftResult.getTruthValuesVector(), rightResult.getTruthValuesVector(), pathFormula.getDiscreteTimeBound()))); return result; } @@ -163,7 +160,7 @@ namespace storm { template std::vector SparseDtmcPrctlModelChecker::computeCumulativeRewardsHelper(uint_fast64_t stepBound) const { - // Only compute the result if the model has at least one reward this->getModel(). + // Only compute the result if the model has at least one reward model. STORM_LOG_THROW(this->getModel().hasStateRewards() || this->getModel().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. @@ -204,7 +201,7 @@ namespace storm { // 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."); - // Initialize result to state rewards of the this->getModel(). + // Initialize result to state rewards of the model. std::vector result(this->getModel().getStateRewardVector()); // Perform the matrix vector multiplication as often as required by the formula bound. @@ -223,7 +220,7 @@ namespace storm { template std::vector SparseDtmcPrctlModelChecker::computeReachabilityRewardsHelper(storm::storage::SparseMatrix const& transitionMatrix, boost::optional> const& stateRewardVector, boost::optional> const& transitionRewardMatrix, storm::storage::SparseMatrix const& backwardTransitions, storm::storage::BitVector const& targetStates, storm::utility::solver::LinearEquationSolverFactory const& linearEquationSolverFactory, bool qualitative) { - // Only compute the result if the model has at least one reward this->getModel(). + // Only compute the result if the model has 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. @@ -236,7 +233,7 @@ namespace storm { STORM_LOG_INFO("Found " << maybeStates.getNumberOfSetBits() << " 'maybe' states."); // Create resulting vector. - std::vector result(transitionMatrix.getRowCount()); + std::vector result(transitionMatrix.getRowCount(), storm::utility::zero()); // Check whether we need to compute exact rewards for some states. if (qualitative) { @@ -291,7 +288,6 @@ namespace storm { } // 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; diff --git a/src/modelchecker/prctl/SparseMdpPrctlModelChecker.cpp b/src/modelchecker/prctl/SparseMdpPrctlModelChecker.cpp index 664a15fdd..bc4bb8a8e 100644 --- a/src/modelchecker/prctl/SparseMdpPrctlModelChecker.cpp +++ b/src/modelchecker/prctl/SparseMdpPrctlModelChecker.cpp @@ -35,36 +35,31 @@ namespace storm { std::vector result(this->getModel().getNumberOfStates(), storm::utility::zero()); // Determine the states that have 0 probability of reaching the target states. - storm::storage::BitVector statesWithProbabilityGreater0; + storm::storage::BitVector maybeStates; if (minimize) { - statesWithProbabilityGreater0 = storm::utility::graph::performProbGreater0A(this->getModel().getTransitionMatrix(), this->getModel().getTransitionMatrix().getRowGroupIndices(), this->getModel().getBackwardTransitions(), phiStates, psiStates, true, stepBound); + maybeStates = storm::utility::graph::performProbGreater0A(this->getModel().getTransitionMatrix(), this->getModel().getTransitionMatrix().getRowGroupIndices(), this->getModel().getBackwardTransitions(), phiStates, psiStates, true, stepBound); } else { - statesWithProbabilityGreater0 = storm::utility::graph::performProbGreater0E(this->getModel().getTransitionMatrix(), this->getModel().getTransitionMatrix().getRowGroupIndices(), this->getModel().getBackwardTransitions(), phiStates, psiStates, true, stepBound); + maybeStates = storm::utility::graph::performProbGreater0E(this->getModel().getTransitionMatrix(), this->getModel().getTransitionMatrix().getRowGroupIndices(), this->getModel().getBackwardTransitions(), phiStates, psiStates, true, stepBound); } - STORM_LOG_INFO("Found " << statesWithProbabilityGreater0.getNumberOfSetBits() << " 'maybe' states."); + maybeStates &= ~psiStates; + STORM_LOG_INFO("Found " << maybeStates.getNumberOfSetBits() << " 'maybe' states."); - if (!statesWithProbabilityGreater0.empty()) { + 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, statesWithProbabilityGreater0, statesWithProbabilityGreater0, false); - - // Compute the new set of target states in the reduced system. - storm::storage::BitVector rightStatesInReducedSystem = psiStates % statesWithProbabilityGreater0; - - // Make all rows absorbing that satisfy the second sub-formula. - submatrix.makeRowGroupsAbsorbing(rightStatesInReducedSystem); + 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(statesWithProbabilityGreater0.getNumberOfSetBits()); - storm::utility::vector::setVectorValues(subresult, rightStatesInReducedSystem, storm::utility::one()); + 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, nullptr, stepBound); + solver->performMatrixVectorMultiplication(minimize, subresult, &b, stepBound); // Set the values of the resulting vector accordingly. - storm::utility::vector::setVectorValues(result, statesWithProbabilityGreater0, subresult); - storm::utility::vector::setVectorValues(result, ~statesWithProbabilityGreater0, storm::utility::zero()); + storm::utility::vector::setVectorValues(result, maybeStates, subresult); } + storm::utility::vector::setVectorValues(result, psiStates, storm::utility::one()); return result; } @@ -76,7 +71,7 @@ namespace storm { std::unique_ptr rightResultPointer = this->check(pathFormula.getRightSubformula()); ExplicitQualitativeCheckResult const& leftResult = leftResultPointer->asExplicitQualitativeCheckResult(); ExplicitQualitativeCheckResult const& rightResult = rightResultPointer->asExplicitQualitativeCheckResult(); - std::unique_ptr result = std::unique_ptr(new ExplicitQuantitativeCheckResult(this->computeBoundedUntilProbabilitiesHelper(optimalityType.get() == storm::logic::OptimalityType::Minimize, leftResult.getTruthValuesVector(), rightResult.getTruthValuesVector(), pathFormula.getUpperBound()))); + std::unique_ptr result = std::unique_ptr(new ExplicitQuantitativeCheckResult(this->computeBoundedUntilProbabilitiesHelper(optimalityType.get() == storm::logic::OptimalityType::Minimize, leftResult.getTruthValuesVector(), rightResult.getTruthValuesVector(), pathFormula.getDiscreteTimeBound()))); return result; } diff --git a/src/modelchecker/results/HybridQuantitativeCheckResult.cpp b/src/modelchecker/results/HybridQuantitativeCheckResult.cpp index f6bb4136d..8132b75b6 100644 --- a/src/modelchecker/results/HybridQuantitativeCheckResult.cpp +++ b/src/modelchecker/results/HybridQuantitativeCheckResult.cpp @@ -1,5 +1,6 @@ #include "src/modelchecker/results/HybridQuantitativeCheckResult.h" #include "src/modelchecker/results/SymbolicQualitativeCheckResult.h" +#include "src/storage/dd/CuddDdManager.h" #include "src/exceptions/InvalidOperationException.h" @@ -122,6 +123,31 @@ namespace storm { this->odd = newOdd; } + template + double HybridQuantitativeCheckResult::getMin() const { + // In order to not get false zeros, we need to set the values of all states whose values is not stored + // symbolically to infinity. + storm::dd::Add tmp = symbolicStates.toAdd().ite(this->symbolicValues, reachableStates.getDdManager()->getConstant(storm::utility::infinity())); + double min = tmp.getMin(); + if (!explicitStates.isZero()) { + for (auto const& element : explicitValues) { + min = std::min(min, element); + } + } + return min; + } + + template + double HybridQuantitativeCheckResult::getMax() const { + double max = this->symbolicValues.getMin(); + if (!explicitStates.isZero()) { + for (auto const& element : explicitValues) { + max = std::max(max, element); + } + } + return max; + } + // Explicitly instantiate the class. template class HybridQuantitativeCheckResult; } diff --git a/src/modelchecker/results/HybridQuantitativeCheckResult.h b/src/modelchecker/results/HybridQuantitativeCheckResult.h index d8239f8c9..160bbbb6d 100644 --- a/src/modelchecker/results/HybridQuantitativeCheckResult.h +++ b/src/modelchecker/results/HybridQuantitativeCheckResult.h @@ -44,6 +44,10 @@ namespace storm { virtual void filter(QualitativeCheckResult const& filter) override; + virtual double getMin() const; + + virtual double getMax() const; + private: // The set of all reachable states. storm::dd::Bdd reachableStates; diff --git a/src/modelchecker/results/SymbolicQuantitativeCheckResult.cpp b/src/modelchecker/results/SymbolicQuantitativeCheckResult.cpp index 088b15a63..c9630a4c3 100644 --- a/src/modelchecker/results/SymbolicQuantitativeCheckResult.cpp +++ b/src/modelchecker/results/SymbolicQuantitativeCheckResult.cpp @@ -2,12 +2,13 @@ #include "src/modelchecker/results/SymbolicQualitativeCheckResult.h" #include "src/storage/dd/CuddDd.h" +#include "src/storage/dd/CuddDdManager.h" #include "src/exceptions/InvalidOperationException.h" namespace storm { namespace modelchecker { template - SymbolicQuantitativeCheckResult::SymbolicQuantitativeCheckResult(storm::dd::Bdd const& reachableStates, storm::dd::Add const& values) : reachableStates(reachableStates), values(values) { + SymbolicQuantitativeCheckResult::SymbolicQuantitativeCheckResult(storm::dd::Bdd const& reachableStates, storm::dd::Add const& values) : reachableStates(reachableStates), states(reachableStates), values(values) { // Intentionally left empty. } @@ -33,7 +34,7 @@ namespace storm { template bool SymbolicQuantitativeCheckResult::isResultForAllStates() const { - return true; + return states == reachableStates; } template @@ -69,9 +70,22 @@ namespace storm { template void SymbolicQuantitativeCheckResult::filter(QualitativeCheckResult const& filter) { STORM_LOG_THROW(filter.isSymbolicQualitativeCheckResult(), storm::exceptions::InvalidOperationException, "Cannot filter symbolic check result with non-symbolic filter."); + this->states &= filter.asSymbolicQualitativeCheckResult().getTruthValuesVector(); this->values *= filter.asSymbolicQualitativeCheckResult().getTruthValuesVector().toAdd(); } + template + double SymbolicQuantitativeCheckResult::getMin() const { + // In order to not get false zeros, we need to set the values of all states whose values is not stored + // symbolically to infinity. + return states.toAdd().ite(this->values, states.getDdManager()->getConstant(storm::utility::infinity())).getMin(); + } + + template + double SymbolicQuantitativeCheckResult::getMax() const { + return this->values.getMax(); + } + // Explicitly instantiate the class. template class SymbolicQuantitativeCheckResult; } diff --git a/src/modelchecker/results/SymbolicQuantitativeCheckResult.h b/src/modelchecker/results/SymbolicQuantitativeCheckResult.h index 5f0ba0600..21a8f3db3 100644 --- a/src/modelchecker/results/SymbolicQuantitativeCheckResult.h +++ b/src/modelchecker/results/SymbolicQuantitativeCheckResult.h @@ -34,10 +34,17 @@ namespace storm { virtual void filter(QualitativeCheckResult const& filter) override; + virtual double getMin() const; + + virtual double getMax() const; + private: // The set of all reachable states. storm::dd::Bdd reachableStates; + // The set of states for which this check result contains values. + storm::dd::Bdd states; + // The values of the quantitative check result. storm::dd::Add values; }; diff --git a/src/solver/NativeLinearEquationSolver.cpp b/src/solver/NativeLinearEquationSolver.cpp index 9206ae301..ab4a0355b 100644 --- a/src/solver/NativeLinearEquationSolver.cpp +++ b/src/solver/NativeLinearEquationSolver.cpp @@ -36,8 +36,6 @@ namespace storm { // Get a Jacobi decomposition of the matrix A. std::pair, std::vector> jacobiDecomposition = A.getJacobiDecomposition(); - std::cout << "LU has " << jacobiDecomposition.first.getNonzeroEntryCount() << " nonzeros." << std::endl; - // To avoid copying the contents of the vector in the loop, we create a temporary x to swap with. bool multiplyResultProvided = true; std::vector* nextX = multiplyResult; diff --git a/src/storage/dd/CuddAdd.cpp b/src/storage/dd/CuddAdd.cpp index 71ec05758..8152083cd 100644 --- a/src/storage/dd/CuddAdd.cpp +++ b/src/storage/dd/CuddAdd.cpp @@ -245,7 +245,7 @@ namespace storm { } } - Add Add::swapVariables(std::vector> const& metaVariablePairs) { + Add Add::swapVariables(std::vector> const& metaVariablePairs) const { std::set newContainedMetaVariables; std::vector from; std::vector to; diff --git a/src/storage/dd/CuddAdd.h b/src/storage/dd/CuddAdd.h index 65de2dd10..8d01839d6 100644 --- a/src/storage/dd/CuddAdd.h +++ b/src/storage/dd/CuddAdd.h @@ -316,7 +316,7 @@ namespace storm { * @param metaVariablePairs A vector of meta variable pairs that are to be swapped for one another. * @return The resulting ADD. */ - Add swapVariables(std::vector> const& metaVariablePairs); + Add swapVariables(std::vector> const& metaVariablePairs) const; /*! * Multiplies the current ADD (representing a matrix) with the given matrix by summing over the given meta diff --git a/src/storage/dd/CuddBdd.cpp b/src/storage/dd/CuddBdd.cpp index 6d920f05d..ca93c1742 100644 --- a/src/storage/dd/CuddBdd.cpp +++ b/src/storage/dd/CuddBdd.cpp @@ -147,7 +147,7 @@ namespace storm { return Bdd(this->getDdManager(), this->getCuddBdd().UnivAbstract(cubeBdd.getCuddBdd()), newMetaVariables); } - Bdd Bdd::swapVariables(std::vector> const& metaVariablePairs) { + Bdd Bdd::swapVariables(std::vector> const& metaVariablePairs) const { std::set newContainedMetaVariables; std::vector from; std::vector to; diff --git a/src/storage/dd/CuddBdd.h b/src/storage/dd/CuddBdd.h index 1868edf2e..326072eb9 100644 --- a/src/storage/dd/CuddBdd.h +++ b/src/storage/dd/CuddBdd.h @@ -166,7 +166,7 @@ namespace storm { * @param metaVariablePairs A vector of meta variable pairs that are to be swapped for one another. * @return The resulting BDD. */ - Bdd swapVariables(std::vector> const& metaVariablePairs); + Bdd swapVariables(std::vector> const& metaVariablePairs) const; /*! * Computes the logical and of the current and the given BDD and existentially abstracts from the given set diff --git a/src/utility/graph.h b/src/utility/graph.h index a49d9a2a7..6e950bcd2 100644 --- a/src/utility/graph.h +++ b/src/utility/graph.h @@ -258,13 +258,14 @@ namespace storm { * through phi states before. * * @param model The (symbolic) model for which to compute the set of states. - * @param transitionMatrixBdd The transition matrix of the model as a BDD. + * @param transitionMatrix The transition matrix of the model as a BDD. * @param phiStates The BDD containing all phi states of the model. * @param psiStates The BDD containing all psi states of the model. + * @param stepBound If given, this number indicates the maximal amount of steps allowed. * @return All states with positive probability. */ template - storm::dd::Bdd performProbGreater0(storm::models::symbolic::Model const& model, storm::dd::Bdd const& transitionMatrix, storm::dd::Bdd const& phiStates, storm::dd::Bdd const& psiStates) { + storm::dd::Bdd performProbGreater0(storm::models::symbolic::Model const& model, storm::dd::Bdd const& transitionMatrix, storm::dd::Bdd const& phiStates, storm::dd::Bdd const& psiStates, boost::optional const& stepBound = boost::optional()) { // Initialize environment for backward search. storm::dd::DdManager const& manager = model.getManager(); storm::dd::Bdd lastIterationStates = manager.getBddZero(); @@ -272,6 +273,10 @@ namespace storm { uint_fast64_t iterations = 0; while (lastIterationStates != statesWithProbabilityGreater0) { + if (stepBound && iterations >= stepBound.get()) { + break; + } + lastIterationStates = statesWithProbabilityGreater0; statesWithProbabilityGreater0 = statesWithProbabilityGreater0.swapVariables(model.getRowColumnMetaVariablePairs()); statesWithProbabilityGreater0 &= transitionMatrix; @@ -284,6 +289,41 @@ namespace storm { return statesWithProbabilityGreater0; } + /*! + * Computes the set of states that have a probability of one of reaching psi states after only passing + * through phi states before. + * + * @param model The (symbolic) model for which to compute the set of states. + * @param transitionMatrix The transition matrix of the model as a BDD. + * @param phiStates The BDD containing all phi states of the model. + * @param psiStates The BDD containing all psi states of the model. + * @param statesWithProbabilityGreater0 The set of states with a positive probability of satisfying phi + * until psi as a BDD. + * @return All states with probability 1. + */ + template + storm::dd::Bdd performProb1(storm::models::symbolic::Model const& model, storm::dd::Bdd const& transitionMatrix, storm::dd::Bdd const& phiStates, storm::dd::Bdd const& psiStates, storm::dd::Bdd const& statesWithProbabilityGreater0) { + storm::dd::Bdd statesWithProbability1 = performProbGreater0(model, transitionMatrix, !psiStates && model.getReachableStates(), !statesWithProbabilityGreater0 && model.getReachableStates()); + statesWithProbability1 = !statesWithProbability1 && model.getReachableStates(); + return statesWithProbability1; + } + + /*! + * Computes the set of states that have a probability of one of reaching psi states after only passing + * through phi states before. + * + * @param model The (symbolic) model for which to compute the set of states. + * @param transitionMatrix The transition matrix of the model as a BDD. + * @param phiStates The BDD containing all phi states of the model. + * @param psiStates The BDD containing all psi states of the model. + * @return All states with probability 1. + */ + template + storm::dd::Bdd performProb1(storm::models::symbolic::Model const& model, storm::dd::Bdd const& transitionMatrix, storm::dd::Bdd const& phiStates, storm::dd::Bdd const& psiStates) { + storm::dd::Bdd statesWithProbabilityGreater0 = performProbGreater0(model, transitionMatrix, phiStates, psiStates); + return performProb1(model, transitionMatrix, phiStates, psiStates, statesWithProbabilityGreater0); + } + /*! * Computes the sets of states that have probability 0 or 1, respectively, of satisfying phi until psi in a * deterministic model. diff --git a/test/functional/builder/die.pm b/test/functional/builder/die.pm index dacde8a79..af0797cff 100644 --- a/test/functional/builder/die.pm +++ b/test/functional/builder/die.pm @@ -29,3 +29,4 @@ label "three" = s=7&d=3; label "four" = s=7&d=4; label "five" = s=7&d=5; label "six" = s=7&d=6; +label "done" = s=7; diff --git a/test/functional/modelchecker/GmmxxHybridDtmcPrctlModelCheckerTest.cpp b/test/functional/modelchecker/GmmxxHybridDtmcPrctlModelCheckerTest.cpp new file mode 100644 index 000000000..cd1a014a9 --- /dev/null +++ b/test/functional/modelchecker/GmmxxHybridDtmcPrctlModelCheckerTest.cpp @@ -0,0 +1,165 @@ +#include "gtest/gtest.h" +#include "storm-config.h" + +#include "src/logic/Formulas.h" +#include "src/utility/solver.h" +#include "src/modelchecker/prctl/HybridDtmcPrctlModelChecker.h" +#include "src/modelchecker/results/HybridQuantitativeCheckResult.h" +#include "src/modelchecker/results/SymbolicQualitativeCheckResult.h" +#include "src/modelchecker/results/SymbolicQuantitativeCheckResult.h" +#include "src/parser/PrismParser.h" +#include "src/builder/DdPrismModelBuilder.h" +#include "src/models/symbolic/Dtmc.h" +#include "src/settings/SettingsManager.h" + +TEST(GmmxxHybridDtmcPrctlModelCheckerTest, Die) { + storm::prism::Program program = storm::parser::PrismParser::parse(STORM_CPP_TESTS_BASE_PATH "/functional/builder/die.pm"); + + // Build the die model with its reward model. + typename storm::builder::DdPrismModelBuilder::Options options; + options.buildRewards = true; + options.rewardModelName = "coin_flips"; + std::shared_ptr> model = storm::builder::DdPrismModelBuilder::translateProgram(program, options); + EXPECT_EQ(13, model->getNumberOfStates()); + EXPECT_EQ(20, model->getNumberOfTransitions()); + + ASSERT_EQ(model->getType(), storm::models::ModelType::Dtmc); + + std::shared_ptr> dtmc = model->as>(); + + storm::modelchecker::HybridDtmcPrctlModelChecker checker(*dtmc, std::unique_ptr>(new storm::utility::solver::GmmxxLinearEquationSolverFactory())); + + auto labelFormula = std::make_shared("one"); + auto eventuallyFormula = std::make_shared(labelFormula); + + std::unique_ptr result = checker.check(*eventuallyFormula); + result->filter(storm::modelchecker::SymbolicQualitativeCheckResult(model->getReachableStates(), model->getInitialStates())); + storm::modelchecker::HybridQuantitativeCheckResult& quantitativeResult1 = result->asHybridQuantitativeCheckResult(); + + EXPECT_NEAR(1.0/6.0, quantitativeResult1.getMin(), storm::settings::gmmxxEquationSolverSettings().getPrecision()); + EXPECT_NEAR(1.0/6.0, quantitativeResult1.getMax(), storm::settings::gmmxxEquationSolverSettings().getPrecision()); + + labelFormula = std::make_shared("two"); + eventuallyFormula = std::make_shared(labelFormula); + + result = checker.check(*eventuallyFormula); + result->filter(storm::modelchecker::SymbolicQualitativeCheckResult(model->getReachableStates(), model->getInitialStates())); + storm::modelchecker::HybridQuantitativeCheckResult& quantitativeResult2 = result->asHybridQuantitativeCheckResult(); + + EXPECT_NEAR(1.0/6.0, quantitativeResult2.getMin(), storm::settings::gmmxxEquationSolverSettings().getPrecision()); + EXPECT_NEAR(1.0/6.0, quantitativeResult2.getMax(), storm::settings::gmmxxEquationSolverSettings().getPrecision()); + + labelFormula = std::make_shared("three"); + eventuallyFormula = std::make_shared(labelFormula); + + result = checker.check(*eventuallyFormula); + result->filter(storm::modelchecker::SymbolicQualitativeCheckResult(model->getReachableStates(), model->getInitialStates())); + storm::modelchecker::HybridQuantitativeCheckResult& quantitativeResult3 = result->asHybridQuantitativeCheckResult(); + + EXPECT_NEAR(1.0/6.0, quantitativeResult3.getMin(), storm::settings::gmmxxEquationSolverSettings().getPrecision()); + EXPECT_NEAR(1.0/6.0, quantitativeResult3.getMax(), storm::settings::gmmxxEquationSolverSettings().getPrecision()); + + auto done = std::make_shared("done"); + auto reachabilityRewardFormula = std::make_shared(done); + + result = checker.check(*reachabilityRewardFormula); + result->filter(storm::modelchecker::SymbolicQualitativeCheckResult(model->getReachableStates(), model->getInitialStates())); + storm::modelchecker::HybridQuantitativeCheckResult& quantitativeResult4 = result->asHybridQuantitativeCheckResult(); + + EXPECT_NEAR(11.0/3.0, quantitativeResult4.getMin(), storm::settings::gmmxxEquationSolverSettings().getPrecision()); + EXPECT_NEAR(11.0/3.0, quantitativeResult4.getMax(), storm::settings::gmmxxEquationSolverSettings().getPrecision()); +} + +TEST(GmmxxHybridDtmcPrctlModelCheckerTest, Crowds) { + storm::prism::Program program = storm::parser::PrismParser::parse(STORM_CPP_TESTS_BASE_PATH "/functional/builder/crowds-5-5.pm"); + + std::shared_ptr> model = storm::builder::DdPrismModelBuilder::translateProgram(program); + EXPECT_EQ(8607, model->getNumberOfStates()); + EXPECT_EQ(15113, model->getNumberOfTransitions()); + + ASSERT_EQ(model->getType(), storm::models::ModelType::Dtmc); + + std::shared_ptr> dtmc = model->as>(); + + storm::modelchecker::HybridDtmcPrctlModelChecker checker(*dtmc, std::unique_ptr>(new storm::utility::solver::GmmxxLinearEquationSolverFactory())); + + auto labelFormula = std::make_shared("observe0Greater1"); + auto eventuallyFormula = std::make_shared(labelFormula); + + std::unique_ptr result = checker.check(*eventuallyFormula); + result->filter(storm::modelchecker::SymbolicQualitativeCheckResult(model->getReachableStates(), model->getInitialStates())); + storm::modelchecker::HybridQuantitativeCheckResult& quantitativeResult1 = result->asHybridQuantitativeCheckResult(); + + EXPECT_NEAR(0.3328800375801578281, quantitativeResult1.getMin(), storm::settings::gmmxxEquationSolverSettings().getPrecision()); + EXPECT_NEAR(0.3328800375801578281, quantitativeResult1.getMax(), storm::settings::gmmxxEquationSolverSettings().getPrecision()); + + labelFormula = std::make_shared("observeIGreater1"); + eventuallyFormula = std::make_shared(labelFormula); + + result = checker.check(*eventuallyFormula); + result->filter(storm::modelchecker::SymbolicQualitativeCheckResult(model->getReachableStates(), model->getInitialStates())); + storm::modelchecker::HybridQuantitativeCheckResult& quantitativeResult2 = result->asHybridQuantitativeCheckResult(); + + EXPECT_NEAR(0.1522194965, quantitativeResult2.getMin(), storm::settings::gmmxxEquationSolverSettings().getPrecision()); + EXPECT_NEAR(0.1522194965, quantitativeResult2.getMax(), storm::settings::gmmxxEquationSolverSettings().getPrecision()); + + labelFormula = std::make_shared("observeOnlyTrueSender"); + eventuallyFormula = std::make_shared(labelFormula); + + result = checker.check(*eventuallyFormula); + result->filter(storm::modelchecker::SymbolicQualitativeCheckResult(model->getReachableStates(), model->getInitialStates())); + storm::modelchecker::HybridQuantitativeCheckResult& quantitativeResult3 = result->asHybridQuantitativeCheckResult(); + + EXPECT_NEAR(0.32153724292835045, quantitativeResult3.getMin(), storm::settings::gmmxxEquationSolverSettings().getPrecision()); + EXPECT_NEAR(0.32153724292835045, quantitativeResult3.getMax(), storm::settings::gmmxxEquationSolverSettings().getPrecision()); +} + +TEST(GmmxxHybridDtmcPrctlModelCheckerTest, SynchronousLeader) { + storm::prism::Program program = storm::parser::PrismParser::parse(STORM_CPP_TESTS_BASE_PATH "/functional/builder/leader-3-5.pm"); + + // Build the die model with its reward model. + typename storm::builder::DdPrismModelBuilder::Options options; + options.buildRewards = true; + options.rewardModelName = "num_rounds"; + std::shared_ptr> model = storm::builder::DdPrismModelBuilder::translateProgram(program, options); + EXPECT_EQ(273, model->getNumberOfStates()); + EXPECT_EQ(397, model->getNumberOfTransitions()); + + ASSERT_EQ(model->getType(), storm::models::ModelType::Dtmc); + + std::shared_ptr> dtmc = model->as>(); + + storm::modelchecker::HybridDtmcPrctlModelChecker checker(*dtmc, std::unique_ptr>(new storm::utility::solver::GmmxxLinearEquationSolverFactory())); + + auto labelFormula = std::make_shared("elected"); + auto eventuallyFormula = std::make_shared(labelFormula); + + std::unique_ptr result = checker.check(*eventuallyFormula); + result->filter(storm::modelchecker::SymbolicQualitativeCheckResult(model->getReachableStates(), model->getInitialStates())); + storm::modelchecker::SymbolicQuantitativeCheckResult& quantitativeResult1 = result->asSymbolicQuantitativeCheckResult(); + + EXPECT_NEAR(1.0, quantitativeResult1.getMin(), storm::settings::gmmxxEquationSolverSettings().getPrecision()); + EXPECT_NEAR(1.0, quantitativeResult1.getMax(), storm::settings::gmmxxEquationSolverSettings().getPrecision()); + + labelFormula = std::make_shared("elected"); + auto trueFormula = std::make_shared(true); + auto boundedUntilFormula = std::make_shared(trueFormula, labelFormula, 20); + + result = checker.check(*boundedUntilFormula); + result->filter(storm::modelchecker::SymbolicQualitativeCheckResult(model->getReachableStates(), model->getInitialStates())); + storm::modelchecker::HybridQuantitativeCheckResult& quantitativeResult2 = result->asHybridQuantitativeCheckResult(); + + EXPECT_NEAR(0.99999989760000074, quantitativeResult2.getMin(), storm::settings::gmmxxEquationSolverSettings().getPrecision()); + EXPECT_NEAR(0.99999989760000074, quantitativeResult2.getMax(), storm::settings::gmmxxEquationSolverSettings().getPrecision()); + + labelFormula = std::make_shared("elected"); + auto reachabilityRewardFormula = std::make_shared(labelFormula); + + result = checker.check(*reachabilityRewardFormula); + result->filter(storm::modelchecker::SymbolicQualitativeCheckResult(model->getReachableStates(), model->getInitialStates())); + storm::modelchecker::HybridQuantitativeCheckResult& quantitativeResult3 = result->asHybridQuantitativeCheckResult(); + + EXPECT_NEAR(1.0416666666666643, quantitativeResult3.getMin(), storm::settings::gmmxxEquationSolverSettings().getPrecision()); + EXPECT_NEAR(1.0416666666666643, quantitativeResult3.getMax(), storm::settings::gmmxxEquationSolverSettings().getPrecision()); +} + diff --git a/test/functional/modelchecker/SparseDtmcPrctlModelCheckerTest.cpp b/test/functional/modelchecker/SparseDtmcPrctlModelCheckerTest.cpp new file mode 100644 index 000000000..c204d36dd --- /dev/null +++ b/test/functional/modelchecker/SparseDtmcPrctlModelCheckerTest.cpp @@ -0,0 +1,129 @@ +#include "gtest/gtest.h" +#include "storm-config.h" + +#include "src/logic/Formulas.h" +#include "src/utility/solver.h" +#include "src/modelchecker/prctl/SparseDtmcPrctlModelChecker.h" +#include "src/modelchecker/results/ExplicitQuantitativeCheckResult.h" +#include "src/settings/SettingsManager.h" +#include "src/settings/SettingMemento.h" +#include "src/parser/AutoParser.h" + +TEST(SparseDtmcPrctlModelCheckerTest, Die) { + std::shared_ptr> abstractModel = storm::parser::AutoParser::parseModel(STORM_CPP_BASE_PATH "/examples/dtmc/die/die.tra", STORM_CPP_BASE_PATH "/examples/dtmc/die/die.lab", "", STORM_CPP_BASE_PATH "/examples/dtmc/die/die.coin_flips.trans.rew"); + + ASSERT_EQ(abstractModel->getType(), storm::models::ModelType::Dtmc); + + std::shared_ptr> dtmc = abstractModel->as>(); + + ASSERT_EQ(dtmc->getNumberOfStates(), 13ull); + ASSERT_EQ(dtmc->getNumberOfTransitions(), 20ull); + + storm::modelchecker::SparseDtmcPrctlModelChecker checker(*dtmc, std::unique_ptr>(new storm::utility::solver::NativeLinearEquationSolverFactory())); + + auto labelFormula = std::make_shared("one"); + auto eventuallyFormula = std::make_shared(labelFormula); + + std::unique_ptr result = checker.check(*eventuallyFormula); + storm::modelchecker::ExplicitQuantitativeCheckResult& quantitativeResult1 = result->asExplicitQuantitativeCheckResult(); + + EXPECT_NEAR(1.0/6.0, quantitativeResult1[0], storm::settings::nativeEquationSolverSettings().getPrecision()); + + labelFormula = std::make_shared("two"); + eventuallyFormula = std::make_shared(labelFormula); + + result = checker.check(*eventuallyFormula); + storm::modelchecker::ExplicitQuantitativeCheckResult& quantitativeResult2 = result->asExplicitQuantitativeCheckResult(); + + EXPECT_NEAR(1.0/6.0, quantitativeResult2[0], storm::settings::nativeEquationSolverSettings().getPrecision()); + + labelFormula = std::make_shared("three"); + eventuallyFormula = std::make_shared(labelFormula); + + result = checker.check(*eventuallyFormula); + storm::modelchecker::ExplicitQuantitativeCheckResult& quantitativeResult3 = result->asExplicitQuantitativeCheckResult(); + + EXPECT_NEAR(1.0/6.0, quantitativeResult3[0], storm::settings::nativeEquationSolverSettings().getPrecision()); + + auto done = std::make_shared("done"); + auto reachabilityRewardFormula = std::make_shared(done); + + result = checker.check(*reachabilityRewardFormula); + storm::modelchecker::ExplicitQuantitativeCheckResult& quantitativeResult4 = result->asExplicitQuantitativeCheckResult(); + + EXPECT_NEAR(3.6666650772094727, quantitativeResult4[0], storm::settings::nativeEquationSolverSettings().getPrecision()); +} + +TEST(SparseDtmcPrctlModelCheckerTest, Crowds) { + std::shared_ptr> abstractModel = storm::parser::AutoParser::parseModel(STORM_CPP_BASE_PATH "/examples/dtmc/crowds/crowds5_5.tra", STORM_CPP_BASE_PATH "/examples/dtmc/crowds/crowds5_5.lab", "", ""); + + ASSERT_EQ(abstractModel->getType(), storm::models::ModelType::Dtmc); + + std::shared_ptr> dtmc = abstractModel->as>(); + + ASSERT_EQ(8607ull, dtmc->getNumberOfStates()); + ASSERT_EQ(15113ull, dtmc->getNumberOfTransitions()); + + storm::modelchecker::SparseDtmcPrctlModelChecker checker(*dtmc, std::unique_ptr>(new storm::utility::solver::NativeLinearEquationSolverFactory())); + + auto labelFormula = std::make_shared("observe0Greater1"); + auto eventuallyFormula = std::make_shared(labelFormula); + + std::unique_ptr result = checker.check(*eventuallyFormula); + storm::modelchecker::ExplicitQuantitativeCheckResult& quantitativeResult1 = result->asExplicitQuantitativeCheckResult(); + + EXPECT_NEAR(0.33288205191646525, quantitativeResult1[0], storm::settings::nativeEquationSolverSettings().getPrecision()); + + labelFormula = std::make_shared("observeIGreater1"); + eventuallyFormula = std::make_shared(labelFormula); + + result = checker.check(*eventuallyFormula); + storm::modelchecker::ExplicitQuantitativeCheckResult& quantitativeResult2 = result->asExplicitQuantitativeCheckResult(); + + EXPECT_NEAR(0.15222066094730619, quantitativeResult2[0], storm::settings::nativeEquationSolverSettings().getPrecision()); + + labelFormula = std::make_shared("observeOnlyTrueSender"); + eventuallyFormula = std::make_shared(labelFormula); + + result = checker.check(*eventuallyFormula); + storm::modelchecker::ExplicitQuantitativeCheckResult& quantitativeResult3 = result->asExplicitQuantitativeCheckResult(); + + EXPECT_NEAR(0.32153900158185761, quantitativeResult3[0], storm::settings::nativeEquationSolverSettings().getPrecision()); +} + +TEST(SparseDtmcPrctlModelCheckerTest, SynchronousLeader) { + std::shared_ptr> abstractModel = storm::parser::AutoParser::parseModel(STORM_CPP_BASE_PATH "/examples/dtmc/synchronous_leader/leader4_8.tra", STORM_CPP_BASE_PATH "/examples/dtmc/synchronous_leader/leader4_8.lab", "", STORM_CPP_BASE_PATH "/examples/dtmc/synchronous_leader/leader4_8.pick.trans.rew"); + + ASSERT_EQ(abstractModel->getType(), storm::models::ModelType::Dtmc); + std::shared_ptr> dtmc = abstractModel->as>(); + + ASSERT_EQ(12400ull, dtmc->getNumberOfStates()); + ASSERT_EQ(16495ull, dtmc->getNumberOfTransitions()); + + storm::modelchecker::SparseDtmcPrctlModelChecker checker(*dtmc, std::unique_ptr>(new storm::utility::solver::NativeLinearEquationSolverFactory())); + + auto labelFormula = std::make_shared("elected"); + auto eventuallyFormula = std::make_shared(labelFormula); + + std::unique_ptr result = checker.check(*eventuallyFormula); + storm::modelchecker::ExplicitQuantitativeCheckResult& quantitativeResult1 = result->asExplicitQuantitativeCheckResult(); + + EXPECT_NEAR(1.0, quantitativeResult1[0], storm::settings::nativeEquationSolverSettings().getPrecision()); + + labelFormula = std::make_shared("elected"); + auto trueFormula = std::make_shared(true); + auto boundedUntilFormula = std::make_shared(trueFormula, labelFormula, 20); + + result = checker.check(*boundedUntilFormula); + storm::modelchecker::ExplicitQuantitativeCheckResult& quantitativeResult2 = result->asExplicitQuantitativeCheckResult(); + + EXPECT_NEAR(0.9999965911265462636, quantitativeResult2[0], storm::settings::nativeEquationSolverSettings().getPrecision()); + + labelFormula = std::make_shared("elected"); + auto reachabilityRewardFormula = std::make_shared(labelFormula); + + result = checker.check(*reachabilityRewardFormula); + storm::modelchecker::ExplicitQuantitativeCheckResult& quantitativeResult3 = result->asExplicitQuantitativeCheckResult(); + + EXPECT_NEAR(1.0448979589010925, quantitativeResult3[0], storm::settings::nativeEquationSolverSettings().getPrecision()); +}