diff --git a/src/formula/NoBoundOperator.h b/src/formula/NoBoundOperator.h index 79578fee2..72f035713 100644 --- a/src/formula/NoBoundOperator.h +++ b/src/formula/NoBoundOperator.h @@ -74,8 +74,8 @@ public: /*! * Empty constructor */ - NoBoundOperator() { - this->pathFormula = NULL; + NoBoundOperator() : optimalityOperator(false), minimumOperator(false) { + this->pathFormula = nullptr; } /*! @@ -83,7 +83,19 @@ public: * * @param pathFormula The child node. */ - NoBoundOperator(AbstractPathFormula* pathFormula) { + NoBoundOperator(AbstractPathFormula* pathFormula) : optimalityOperator(false), minimumOperator(false) { + this->pathFormula = pathFormula; + } + + /*! + * Constructor + * + * @param pathFormula The child node. + * @param minimumOperator A flag indicating whether this operator is a minimizing or a + * maximizing operator. + */ + NoBoundOperator(AbstractPathFormula* pathFormula, bool minimumOperator) + : optimalityOperator(true), minimumOperator(minimumOperator) { this->pathFormula = pathFormula; } @@ -142,8 +154,34 @@ public: return checker.conforms(this->pathFormula); } + /*! + * Retrieves whether the operator is to be interpreted as an optimizing (i.e. min/max) operator. + * @returns True if the operator is an optimizing operator. + */ + bool isOptimalityOperator() const { + return optimalityOperator; + } + + /*! + * Retrieves whether the operator is a minimizing operator given that it is an optimality + * operator. + * @returns True if the operator is an optimizing operator and it is a minimizing operator and + * false otherwise, i.e. if it is either not an optimizing operator or not a minimizing operator. + */ + bool isMinimumOperator() const { + return optimalityOperator && minimumOperator; + } + private: AbstractPathFormula* pathFormula; + + // A flag that indicates whether this operator is meant as an optimizing (i.e. min/max) operator + // over a nondeterministic model. + bool optimalityOperator; + + // In the case this operator is an optimizing operator, this flag indicates whether it is + // looking for the minimum or the maximum value. + bool minimumOperator; }; } /* namespace formula */ diff --git a/src/formula/ProbabilisticNoBoundOperator.h b/src/formula/ProbabilisticNoBoundOperator.h index 40c812cf7..cae4df96e 100644 --- a/src/formula/ProbabilisticNoBoundOperator.h +++ b/src/formula/ProbabilisticNoBoundOperator.h @@ -51,7 +51,7 @@ public: /*! * Empty constructor */ - ProbabilisticNoBoundOperator() : NoBoundOperator(nullptr), optimalityOperator(false), minimumOperator(false) { + ProbabilisticNoBoundOperator() : NoBoundOperator(nullptr) { // Intentionally left empty } @@ -60,18 +60,7 @@ public: * * @param pathFormula The child node. */ - ProbabilisticNoBoundOperator(AbstractPathFormula* pathFormula) : NoBoundOperator(pathFormula), - optimalityOperator(false), minimumOperator(false) { - // Intentionally left empty - } - - /*! - * Constructor - * - * @param pathFormula The child node. - */ - ProbabilisticNoBoundOperator(AbstractPathFormula* pathFormula, bool minimumOperator) : NoBoundOperator(pathFormula), - optimalityOperator(true), minimumOperator(minimumOperator) { + ProbabilisticNoBoundOperator(AbstractPathFormula* pathFormula) : NoBoundOperator(pathFormula) { // Intentionally left empty } @@ -80,8 +69,8 @@ public: */ virtual std::string toString() const { std::string result = "P"; - if (optimalityOperator) { - if (minimumOperator) { + if (this->isOptimalityOperator()) { + if (this->isMinimumOperator()) { result += "min"; } else { result += "max"; @@ -92,15 +81,6 @@ public: result += "]"; return result; } - -private: - // A flag that indicates whether this operator is meant as an optimizing (i.e. min/max) operator - // over a nondeterministic model. - bool optimalityOperator; - - // In the case this operator is an optimizing operator, this flag indicates whether it is - // looking for the minimum or the maximum value. - bool minimumOperator; }; } /* namespace formula */ diff --git a/src/modelchecker/DtmcPrctlModelChecker.h b/src/modelchecker/DtmcPrctlModelChecker.h index 48b10d54b..32e1baa21 100644 --- a/src/modelchecker/DtmcPrctlModelChecker.h +++ b/src/modelchecker/DtmcPrctlModelChecker.h @@ -178,6 +178,7 @@ public: } if (!this->getModel().hasAtomicProposition(formula.getAp())) { + LOG4CPLUS_ERROR(logger, "Atomic proposition '" << formula.getAp() << "' is invalid."); throw storm::exceptions::InvalidPropertyException() << "Atomic proposition '" << formula.getAp() << "' is invalid."; return nullptr; } @@ -246,6 +247,10 @@ public: * @returns The set of states satisfying the formula, represented by a bit vector */ std::vector* checkNoBoundOperator(const storm::formula::NoBoundOperator& formula) const { + // Check if the operator was an optimality operator and report a warning in that case. + if (formula.isOptimalityOperator()) { + LOG4CPLUS_WARN(logger, "Formula contains min/max operator which is not meaningful over deterministic models."); + } return formula.getPathFormula().check(*this); } diff --git a/src/modelchecker/GmmxxDtmcPrctlModelChecker.h b/src/modelchecker/GmmxxDtmcPrctlModelChecker.h index 359462c1f..3d3971b35 100644 --- a/src/modelchecker/GmmxxDtmcPrctlModelChecker.h +++ b/src/modelchecker/GmmxxDtmcPrctlModelChecker.h @@ -329,7 +329,7 @@ public: * @return The name of this module. */ static std::string getModuleName() { - return "gmm++"; + return "gmm++det"; } /*! diff --git a/src/modelchecker/GmmxxMdpPrctlModelChecker.h b/src/modelchecker/GmmxxMdpPrctlModelChecker.h new file mode 100644 index 000000000..d0d940eba --- /dev/null +++ b/src/modelchecker/GmmxxMdpPrctlModelChecker.h @@ -0,0 +1,409 @@ +/* + * GmmxxDtmcPrctlModelChecker.h + * + * Created on: 06.12.2012 + * Author: Christian Dehnert + */ + +#ifndef STORM_MODELCHECKER_GMMXXDTMCPRCTLMODELCHECKER_H_ +#define STORM_MODELCHECKER_GMMXXDTMCPRCTLMODELCHECKER_H_ + +#include + +#include "src/models/Mdp.h" +#include "src/modelchecker/MdpPrctlModelChecker.h" +#include "src/utility/GraphAnalyzer.h" +#include "src/utility/Vector.h" +#include "src/utility/ConstTemplates.h" +#include "src/utility/Settings.h" +#include "src/adapters/GmmxxAdapter.h" +#include "src/exceptions/InvalidPropertyException.h" +#include "src/storage/JacobiDecomposition.h" + +#include "gmm/gmm_matrix.h" +#include "gmm/gmm_iter_solvers.h" + +#include "log4cplus/logger.h" +#include "log4cplus/loggingmacros.h" + +extern log4cplus::Logger logger; + +namespace storm { + +namespace modelChecker { + +/* + * A model checking engine that makes use of the gmm++ backend. + */ +template +class GmmxxDtmcPrctlModelChecker : public MdpPrctlModelChecker { + +public: + explicit GmmxxDtmcPrctlModelChecker(storm::models::Mdp& mdp) : MdpPrctlModelChecker(mdp) { } + + virtual ~GmmxxDtmcPrctlModelChecker() { } + + virtual std::vector* checkBoundedUntil(const storm::formula::BoundedUntil& formula) const { + // First, we need to compute the states that satisfy the sub-formulas of the until-formula. + storm::storage::BitVector* leftStates = this->checkStateFormula(formula.getLeft()); + storm::storage::BitVector* rightStates = this->checkStateFormula(formula.getRight()); + + // Copy the matrix before we make any changes. + storm::storage::SparseMatrix tmpMatrix(*this->getModel().getTransitionMatrix()); + + // Get the starting row indices for the non-deterministic choices to reduce the resulting + // vector properly. + std::shared_ptr> nondeterministicChoiceIndices = this->getModel().getNondeterministicChoiceIndices(); + + // Make all rows absorbing that violate both sub-formulas or satisfy the second sub-formula. + tmpMatrix.makeRowsAbsorbing(~(*leftStates | *rightStates) | *rightStates); + + // Transform the transition probability matrix to the gmm++ format to use its arithmetic. + gmm::csr_matrix* gmmxxMatrix = storm::adapters::GmmxxAdapter::toGmmxxSparseMatrix(tmpMatrix); + + // Create the vector with which to multiply. + std::vector* result = new std::vector(this->getModel().getNumberOfStates()); + storm::utility::setVectorValues(result, *rightStates, storm::utility::constGetOne()); + + // Create vector for result of multiplication, which is reduced to the result vector after + // each multiplication. + std::vector* multiplyResult = new std::vector(this->getModel().getTransitionMatrix().getRowCount()); + + // Now perform matrix-vector multiplication as long as we meet the bound of the formula. + for (uint_fast64_t i = 0; i < formula.getBound(); ++i) { + gmm::mult(*gmmxxMatrix, *result, *multiplyResult); + + if (minimumOperatorStack.top()) { + storm::utility::reduceVectorMin(*multiplyResult, result, *nondeterministicChoiceIndices); + } else { + storm::utility::reduceVectorMax(*multiplyResult, result, *nondeterministicChoiceIndices); + } + } + delete multiplyResult; + + // Delete intermediate results and return result. + delete gmmxxMatrix; + delete leftStates; + delete rightStates; + return result; + } + + virtual std::vector* checkNext(const storm::formula::Next& formula) const { + // First, we need to compute the states that satisfy the sub-formula of the next-formula. + storm::storage::BitVector* nextStates = this->checkStateFormula(formula.getChild()); + + // Transform the transition probability matrix to the gmm++ format to use its arithmetic. + gmm::csr_matrix* gmmxxMatrix = storm::adapters::GmmxxAdapter::toGmmxxSparseMatrix(*this->getModel().getTransitionMatrix()); + + // Create the vector with which to multiply and initialize it correctly. + std::vector* result = new std::vector(this->getModel().getNumberOfStates()); + storm::utility::setVectorValues(result, *nextStates, storm::utility::constGetOne()); + + // Delete obsolete sub-result. + delete nextStates; + + // Create resulting vector. + std::vector* temporaryResult = new std::vector(this->getModel().getTransitionMatrix().getRowCount()); + + // Perform the actual computation, namely matrix-vector multiplication. + gmm::mult(*gmmxxMatrix, *result, *temporaryResult); + + // Get the starting row indices for the non-deterministic choices to reduce the resulting + // vector properly. + std::shared_ptr> nondeterministicChoiceIndices = this->getModel().getNondeterministicChoiceIndices(); + + if (minimumOperatorStack.top()) { + storm::utility::reduceVectorMin(*temporaryResult, result, *nondeterministicChoiceIndices); + } else { + storm::utility::reduceVectorMax(*temporaryResult, result, *nondeterministicChoiceIndices); + } + + // Delete temporary matrix plus temporary result and return result. + delete gmmxxMatrix; + delete temporaryResult; + return result; + } + + virtual std::vector* checkUntil(const storm::formula::Until& formula) const { + // First, we need to compute the states that satisfy the sub-formulas of the until-formula. + storm::storage::BitVector* leftStates = this->checkStateFormula(formula.getLeft()); + storm::storage::BitVector* rightStates = this->checkStateFormula(formula.getRight()); + + // Then, 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. + storm::storage::BitVector statesWithProbability0(this->getModel().getNumberOfStates()); + storm::storage::BitVector statesWithProbability1(this->getModel().getNumberOfStates()); + if (minimumOperatorStack.top()) { + storm::utility::GraphAnalyzer::performProb01Min(this->getModel(), *leftStates, *rightStates, &statesWithProbability0, &statesWithProbability1); + } else { + storm::utility::GraphAnalyzer::performProb01Max(this->getModel(), *leftStates, *rightStates, &statesWithProbability0, &statesWithProbability1); + } + + // Delete sub-results that are obsolete now. + delete leftStates; + delete rightStates; + + LOG4CPLUS_INFO(logger, "Found " << statesWithProbability0.getNumberOfSetBits() << " 'no' states."); + LOG4CPLUS_INFO(logger, "Found " << statesWithProbability1.getNumberOfSetBits() << " 'yes' states."); + storm::storage::BitVector maybeStates = ~(statesWithProbability0 | statesWithProbability1); + LOG4CPLUS_INFO(logger, "Found " << maybeStates.getNumberOfSetBits() << " 'maybe' states."); + + // Create resulting vector. + std::vector* result = new std::vector(this->getModel().getNumberOfStates()); + + // Only try to solve system if there are states for which the probability is unknown. + uint_fast64_t mayBeStatesSetBitCount = maybeStates.getNumberOfSetBits(); + if (mayBeStatesSetBitCount > 0) { + // Now we can eliminate the rows and columns from the original transition probability matrix. + storm::storage::SparseMatrix* submatrix = this->getModel().getTransitionMatrix()->getSubmatrix(maybeStates); + + // Transform the submatrix to the gmm++ format to use its solvers. + gmm::csr_matrix* gmmxxMatrix = storm::adapters::GmmxxAdapter::toGmmxxSparseMatrix(*submatrix); + delete submatrix; + + // Initialize the x vector with 0.5 for each element. This is the initial guess for + // iteratively solving the equations. + std::vector x(mayBeStatesSetBitCount, Type(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(mayBeStatesSetBitCount); + this->getModel().getTransitionMatrix()->getConstrainedRowCountVector(maybeStates, statesWithProbability1, &b); + + // Solve the corresponding system of linear equations. + this->solveEquationSystem(*gmmxxMatrix, x, b); + + // Set values of resulting vector according to result. + storm::utility::setVectorValues(result, maybeStates, x); + + // Delete temporary matrix. + delete gmmxxMatrix; + } + + // Set values of resulting vector that are known exactly. + storm::utility::setVectorValues(result, statesWithProbability0, storm::utility::constGetZero()); + storm::utility::setVectorValues(result, statesWithProbability1, storm::utility::constGetOne()); + + return result; + } + + virtual std::vector* checkInstantaneousReward(const storm::formula::InstantaneousReward& formula) const { + // Only compute the result if the model has a state-based reward model. + if (!this->getModel().hasStateRewards()) { + LOG4CPLUS_ERROR(logger, "Missing (state-based) reward model for formula."); + throw storm::exceptions::InvalidPropertyException() << "Missing (state-based) reward model for formula."; + } + + // Transform the transition probability matrix to the gmm++ format to use its arithmetic. + gmm::csr_matrix* gmmxxMatrix = storm::adapters::GmmxxAdapter::toGmmxxSparseMatrix(*this->getModel().getTransitionMatrix()); + + // Initialize result to state rewards of the model. + std::vector* result = new std::vector(*this->getModel().getStateRewardVector()); + + // Now perform matrix-vector multiplication as long as we meet the bound of the formula. + std::vector* swap = nullptr; + std::vector* tmpResult = new std::vector(this->getModel().getNumberOfStates()); + for (uint_fast64_t i = 0; i < formula.getBound(); ++i) { + gmm::mult(*gmmxxMatrix, *result, *tmpResult); + swap = tmpResult; + tmpResult = result; + result = swap; + } + + // Delete temporary variables and return result. + delete tmpResult; + delete gmmxxMatrix; + return result; + } + + virtual std::vector* checkCumulativeReward(const storm::formula::CumulativeReward& formula) const { + // Only compute the result if the model has at least one reward model. + if (!this->getModel().hasStateRewards() && !this->getModel().hasTransitionRewards()) { + LOG4CPLUS_ERROR(logger, "Missing reward model for formula."); + throw storm::exceptions::InvalidPropertyException() << "Missing reward model for formula."; + } + + // Transform the transition probability matrix to the gmm++ format to use its arithmetic. + gmm::csr_matrix* gmmxxMatrix = storm::adapters::GmmxxAdapter::toGmmxxSparseMatrix(*this->getModel().getTransitionMatrix()); + + // Compute the reward vector to add in each step based on the available reward models. + std::vector* totalRewardVector = nullptr; + if (this->getModel().hasTransitionRewards()) { + totalRewardVector = this->getModel().getTransitionMatrix()->getPointwiseProductRowSumVector(*this->getModel().getTransitionRewardMatrix()); + if (this->getModel().hasStateRewards()) { + gmm::add(*this->getModel().getStateRewardVector(), *totalRewardVector); + } + } else { + totalRewardVector = new std::vector(*this->getModel().getStateRewardVector()); + } + + std::vector* result = new std::vector(this->getModel().getNumberOfStates()); + + // Now perform matrix-vector multiplication as long as we meet the bound of the formula. + std::vector* swap = nullptr; + std::vector* tmpResult = new std::vector(this->getModel().getNumberOfStates()); + for (uint_fast64_t i = 0; i < formula.getBound(); ++i) { + gmm::mult(*gmmxxMatrix, *result, *tmpResult); + swap = tmpResult; + tmpResult = result; + result = swap; + + // Add the reward vector to the result. + gmm::add(*totalRewardVector, *result); + } + + // Delete temporary variables and return result. + delete tmpResult; + delete gmmxxMatrix; + delete totalRewardVector; + return result; + } + + virtual std::vector* checkReachabilityReward(const storm::formula::ReachabilityReward& formula) const { + // Only compute the result if the model has at least one reward model. + if (!this->getModel().hasStateRewards() && !this->getModel().hasTransitionRewards()) { + LOG4CPLUS_ERROR(logger, "Missing reward model for formula. Skipping formula"); + throw storm::exceptions::InvalidPropertyException() << "Missing reward model for formula."; + } + + // Determine the states for which the target predicate holds. + storm::storage::BitVector* targetStates = this->checkStateFormula(formula.getChild()); + + // Determine which states have a reward of infinity by definition. + storm::storage::BitVector infinityStates(this->getModel().getNumberOfStates()); + storm::storage::BitVector trueStates(this->getModel().getNumberOfStates(), true); + storm::utility::GraphAnalyzer::performProb1(this->getModel(), trueStates, *targetStates, &infinityStates); + infinityStates.complement(); + + // Create resulting vector. + std::vector* result = new std::vector(this->getModel().getNumberOfStates()); + + // Check whether there are states for which we have to compute the result. + storm::storage::BitVector maybeStates = ~(*targetStates) & ~infinityStates; + const int maybeStatesSetBitCount = maybeStates.getNumberOfSetBits(); + if (maybeStatesSetBitCount > 0) { + // Now we can eliminate the rows and columns from the original transition probability matrix. + storm::storage::SparseMatrix* submatrix = this->getModel().getTransitionMatrix()->getSubmatrix(maybeStates); + // 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(); + + // Transform the submatrix to the gmm++ format to use its solvers. + gmm::csr_matrix* gmmxxMatrix = storm::adapters::GmmxxAdapter::toGmmxxSparseMatrix(*submatrix); + delete submatrix; + + // Initialize the x vector with 1 for each element. This is the initial guess for + // the iterative solvers. + std::vector x(maybeStatesSetBitCount, storm::utility::constGetOne()); + + // Prepare the right-hand side of the equation system. + std::vector* b = new std::vector(maybeStatesSetBitCount); + if (this->getModel().hasTransitionRewards()) { + // If a transition-based reward model is available, we initialize the right-hand + // side to the vector resulting from summing the rows of the pointwise product + // of the transition probability matrix and the transition reward matrix. + std::vector* pointwiseProductRowSumVector = this->getModel().getTransitionMatrix()->getPointwiseProductRowSumVector(*this->getModel().getTransitionRewardMatrix()); + storm::utility::selectVectorValues(b, maybeStates, *pointwiseProductRowSumVector); + delete pointwiseProductRowSumVector; + + if (this->getModel().hasStateRewards()) { + // If a state-based reward model is also available, we need to add this vector + // as well. As the state reward vector contains entries not just for the states + // that we still consider (i.e. maybeStates), we need to extract these values + // first. + std::vector* subStateRewards = new std::vector(maybeStatesSetBitCount); + storm::utility::setVectorValues(subStateRewards, maybeStates, *this->getModel().getStateRewardVector()); + gmm::add(*subStateRewards, *b); + delete subStateRewards; + } + } else { + // If only a state-based reward model is available, we take this vector as the + // right-hand side. As the state reward vector contains entries not just for the + // states that we still consider (i.e. maybeStates), we need to extract these values + // first. + storm::utility::setVectorValues(b, maybeStates, *this->getModel().getStateRewardVector()); + } + + // Solve the corresponding system of linear equations. + this->solveLinearEquationSystem(*gmmxxMatrix, x, *b); + + // Set values of resulting vector according to result. + storm::utility::setVectorValues(result, maybeStates, x); + + // Delete temporary matrix and right-hand side. + delete gmmxxMatrix; + delete b; + } + + // Set values of resulting vector that are known exactly. + storm::utility::setVectorValues(result, *targetStates, storm::utility::constGetZero()); + storm::utility::setVectorValues(result, infinityStates, storm::utility::constGetInfinity()); + + // Delete temporary storages and return result. + delete targetStates; + return result; + } + + /*! + * Returns the name of this module. + * @return The name of this module. + */ + static std::string getModuleName() { + return "gmm++nondet"; + } + + /*! + * Returns a trigger such that if the option "matrixlib" is set to "gmm++", this model checker + * is to be used. + * @return An option trigger for this module. + */ + static std::pair getOptionTrigger() { + return std::pair("matrixlib", "gmm++"); + } + + /*! + * Registers all options associated with the gmm++ matrix library. + */ + static void putOptions(boost::program_options::options_description* desc) { + desc->add_options()("lemaxiter", boost::program_options::value()->default_value(10000), "Sets the maximal number of iterations for iterative equation solving."); + desc->add_options()("precision", boost::program_options::value()->default_value(10e-6), "Sets the precision for iterative linear equation solving."); + } + +private: + /*! + * Solves the given equation system under the given parameters using the power method. + * + * @param A The matrix A specifying the coefficients of the equations. + * @param x The vector x for which to solve the equations. The initial value of the elements of + * this vector are used as the initial guess and might thus influence performance and convergence. + * @param b The vector b specifying the values on the right-hand-sides of the equations. + * @return The solution of the system of linear equations in form of the elements of the vector + * x. + */ + void solveLinearEquationSystem(gmm::csr_matrix const& A, std::vector& x, std::vector const& b) const { + // Get the settings object to customize linear solving. + storm::settings::Settings* s = storm::settings::instance(); + + // Get relevant user-defined settings for solving the equations + double precison = s->get("precision"); + unsigned maxIterations = s->get("lemaxiter"); + + + + + + // Check if the solver converged and issue a warning otherwise. + if (iter.converged()) { + LOG4CPLUS_INFO(logger, "Iterative solver converged after " << iter.get_iteration() << " iterations."); + } else { + LOG4CPLUS_WARN(logger, "Iterative solver did not converge."); + } + + } +}; + +} //namespace modelChecker + +} //namespace storm + +#endif /* STORM_MODELCHECKER_GMMXXDTMCPRCTLMODELCHECKER_H_ */ diff --git a/src/modelchecker/MdpPrctlModelChecker.h b/src/modelchecker/MdpPrctlModelChecker.h new file mode 100644 index 000000000..e3139062f --- /dev/null +++ b/src/modelchecker/MdpPrctlModelChecker.h @@ -0,0 +1,375 @@ +/* + * MdpPrctlModelChecker.h + * + * Created on: 15.02.2013 + * Author: Christian Dehnert + */ + +#ifndef STORM_MODELCHECKER_MDPPRCTLMODELCHECKER_H_ +#define STORM_MODELCHECKER_MDPPRCTLMODELCHECKER_H_ + +#include "src/formula/Formulas.h" +#include "src/utility/Vector.h" +#include "src/storage/SparseMatrix.h" + +#include "src/models/Mdp.h" +#include "src/storage/BitVector.h" +#include "src/exceptions/InvalidPropertyException.h" +#include "src/utility/Vector.h" +#include "src/modelchecker/AbstractModelChecker.h" +#include + +#include "log4cplus/logger.h" +#include "log4cplus/loggingmacros.h" + +extern log4cplus::Logger logger; + +namespace storm { + +namespace modelChecker { + +/*! + * @brief + * Interface for model checker classes. + * + * This class provides basic functions that are the same for all subclasses, but mainly only declares + * abstract methods that are to be implemented in concrete instances. + * + * @attention This class is abstract. + */ +template +class MdpPrctlModelChecker : + public virtual AbstractModelChecker { +public: + /*! + * Constructor + * + * @param model The dtmc model which is checked. + */ + explicit MdpPrctlModelChecker(storm::models::Mdp& model) : model(model), minimumOperatorStack() { + + } + + /*! + * Copy constructor + * + * @param modelChecker The model checker that is copied. + */ + explicit MdpPrctlModelChecker(const storm::modelChecker::MdpPrctlModelChecker* modelChecker) : model(new storm::models::Mdp(modelChecker->getModel())), minimumOperatorStack() { + + } + + /*! + * Destructor + */ + virtual ~MdpPrctlModelChecker() { + // Intentionally left empty. + } + + /*! + * @returns A reference to the dtmc of the model checker. + */ + storm::models::Mdp& getModel() const { + return this->model; + } + + /*! + * Sets the DTMC model which is checked + * @param model + */ + void setModel(storm::models::Mdp& model) { + this->model = &model; + } + + /*! + * Checks the given state formula on the DTMC and prints the result (true/false) for all initial + * states. + * @param stateFormula The formula to be checked. + */ + void check(const storm::formula::AbstractStateFormula& stateFormula) const { + std::cout << std::endl; + LOG4CPLUS_INFO(logger, "Model checking formula\t" << stateFormula.toString()); + std::cout << "Model checking formula:\t" << stateFormula.toString() << std::endl; + storm::storage::BitVector* result = nullptr; + try { + result = stateFormula.check(*this); + LOG4CPLUS_INFO(logger, "Result for initial states:"); + std::cout << "Result for initial states:" << std::endl; + for (auto initialState : *this->getModel().getLabeledStates("init")) { + LOG4CPLUS_INFO(logger, "\t" << initialState << ": " << (result->get(initialState) ? "satisfied" : "not satisfied")); + std::cout << "\t" << initialState << ": " << (*result)[initialState] << std::endl; + } + delete result; + } catch (std::exception& e) { + std::cout << "Error during computation: " << e.what() << "Skipping property." << std::endl; + if (result != nullptr) { + delete result; + } + } + std::cout << std::endl; + storm::utility::printSeparationLine(std::cout); + } + + /*! + * Checks the given operator (with no bound) on the DTMC and prints the result + * (probability/rewards) for all initial states. + * @param noBoundFormula The formula to be checked. + */ + void check(const storm::formula::NoBoundOperator& noBoundFormula) const { + std::cout << std::endl; + LOG4CPLUS_INFO(logger, "Model checking formula\t" << noBoundFormula.toString()); + std::cout << "Model checking formula:\t" << noBoundFormula.toString() << std::endl; + std::vector* result = nullptr; + try { + result = noBoundFormula.check(*this); + LOG4CPLUS_INFO(logger, "Result for initial states:"); + std::cout << "Result for initial states:" << std::endl; + for (auto initialState : *this->getModel().getLabeledStates("init")) { + LOG4CPLUS_INFO(logger, "\t" << initialState << ": " << (*result)[initialState]); + std::cout << "\t" << initialState << ": " << (*result)[initialState] << std::endl; + } + delete result; + } catch (std::exception& e) { + std::cout << "Error during computation: " << e.what() << " Skipping property." << std::endl; + if (result != nullptr) { + delete result; + } + } + std::cout << std::endl; + storm::utility::printSeparationLine(std::cout); + } + + /*! + * The check method for a state formula; Will infer the actual type of formula and delegate it + * to the specialized method + * + * @param formula The state formula to check + * @returns The set of states satisfying the formula, represented by a bit vector + */ + storm::storage::BitVector* checkStateFormula(const storm::formula::AbstractStateFormula& formula) const { + return formula.check(*this); + } + + /*! + * The check method for a state formula with an And node as root in its formula tree + * + * @param formula The And formula to check + * @returns The set of states satisfying the formula, represented by a bit vector + */ + storm::storage::BitVector* checkAnd(const storm::formula::And& formula) const { + storm::storage::BitVector* result = checkStateFormula(formula.getLeft()); + storm::storage::BitVector* right = checkStateFormula(formula.getRight()); + (*result) &= (*right); + delete right; + return result; + } + + /*! + * The check method for a formula with an AP node as root in its formula tree + * + * @param formula The Ap state formula to check + * @returns The set of states satisfying the formula, represented by a bit vector + */ + storm::storage::BitVector* checkAp(const storm::formula::Ap& formula) const { + if (formula.getAp().compare("true") == 0) { + return new storm::storage::BitVector(this->getModel().getNumberOfStates(), true); + } else if (formula.getAp().compare("false") == 0) { + return new storm::storage::BitVector(this->getModel().getNumberOfStates()); + } + + if (!this->getModel().hasAtomicProposition(formula.getAp())) { + LOG4CPLUS_ERROR(logger, "Atomic proposition '" << formula.getAp() << "' is invalid."); + throw storm::exceptions::InvalidPropertyException() << "Atomic proposition '" << formula.getAp() << "' is invalid."; + return nullptr; + } + + return new storm::storage::BitVector(*this->getModel().getLabeledStates(formula.getAp())); + } + + /*! + * The check method for a formula with a Not node as root in its formula tree + * + * @param formula The Not state formula to check + * @returns The set of states satisfying the formula, represented by a bit vector + */ + storm::storage::BitVector* checkNot(const storm::formula::Not& formula) const { + storm::storage::BitVector* result = checkStateFormula(formula.getChild()); + result->complement(); + return result; + } + + /*! + * The check method for a state formula with an Or node as root in its formula tree + * + * @param formula The Or state formula to check + * @returns The set of states satisfying the formula, represented by a bit vector + */ + virtual storm::storage::BitVector* checkOr(const storm::formula::Or& formula) const { + storm::storage::BitVector* result = checkStateFormula(formula.getLeft()); + storm::storage::BitVector* right = checkStateFormula(formula.getRight()); + (*result) |= (*right); + delete right; + return result; + } + + /*! + * The check method for a state formula with a bound operator node as root in + * its formula tree + * + * @param formula The state formula to check + * @returns The set of states satisfying the formula, represented by a bit vector + */ + storm::storage::BitVector* checkPathBoundOperator(const storm::formula::PathBoundOperator& formula) const { + // First, we need to compute the probability for satisfying the path formula for each state. + std::vector* quantitativeResult = this->checkPathFormula(formula.getPathFormula()); + + // Create resulting bit vector, which will hold the yes/no-answer for every state. + storm::storage::BitVector* result = new storm::storage::BitVector(this->getModel().getNumberOfStates()); + + // Now, we can compute which states meet the bound specified in this operator and set the + // corresponding bits to true in the resulting vector. + for (uint_fast64_t i = 0; i < this->getModel().getNumberOfStates(); ++i) { + if (formula.meetsBound((*quantitativeResult)[i])) { + result->set(i, true); + } + } + + // Delete the probabilities computed for the states and return result. + delete quantitativeResult; + return result; + } + + /*! + * The check method for a state formula with a probabilistic operator node without bounds as root + * in its formula tree + * + * @param formula The state formula to check + * @returns The set of states satisfying the formula, represented by a bit vector + */ + std::vector* checkNoBoundOperator(const storm::formula::NoBoundOperator& formula) const { + // Check if the operator was an non-optimality operator and report an error in that case. + if (!formula.isOptimalityOperator()) { + LOG4CPLUS_ERROR(logger, "Formula does not specify neither min nor max optimality, which is not meaningful over nondeterministic models."); + throw storm::exceptions::InvalidArgumentException() << "Formula does not specify neither min nor max optimality, which is not meaningful over nondeterministic models."; + } + minimumOperatorStack.push(formula.isMinimumOperator()); + std::vector* result = formula.getPathFormula().check(*this); + minimumOperatorStack.pop(); + return result; + } + + /*! + * The check method for a path formula; Will infer the actual type of formula and delegate it + * to the specialized method + * + * @param formula The path formula to check + * @returns for each state the probability that the path formula holds. + */ + std::vector* checkPathFormula(const storm::formula::AbstractPathFormula& formula) const { + return formula.check(*this); + } + + /*! + * The check method for a path formula with a Bounded Until operator node as root in its formula tree + * + * @param formula The Bounded Until path formula to check + * @returns for each state the probability that the path formula holds. + */ + virtual std::vector* checkBoundedUntil(const storm::formula::BoundedUntil& formula) const = 0; + + /*! + * The check method for a path formula with a Next operator node as root in its formula tree + * + * @param formula The Next path formula to check + * @returns for each state the probability that the path formula holds. + */ + virtual std::vector* checkNext(const storm::formula::Next& formula) const = 0; + + /*! + * The check method for a path formula with a Bounded Eventually operator node as root in its + * formula tree + * + * @param formula The Bounded Eventually path formula to check + * @returns for each state the probability that the path formula holds + */ + virtual std::vector* checkBoundedEventually(const storm::formula::BoundedEventually& formula) const { + // Create equivalent temporary bounded until formula and check it. + storm::formula::BoundedUntil temporaryBoundedUntilFormula(new storm::formula::Ap("true"), formula.getChild().clone(), formula.getBound()); + return this->checkBoundedUntil(temporaryBoundedUntilFormula); + } + + /*! + * The check method for a path formula with an Eventually operator node as root in its formula tree + * + * @param formula The Eventually path formula to check + * @returns for each state the probability that the path formula holds + */ + virtual std::vector* checkEventually(const storm::formula::Eventually& formula) const { + // Create equivalent temporary until formula and check it. + storm::formula::Until temporaryUntilFormula(new storm::formula::Ap("true"), formula.getChild().clone()); + return this->checkUntil(temporaryUntilFormula); + } + + /*! + * The check method for a path formula with a Globally operator node as root in its formula tree + * + * @param formula The Globally path formula to check + * @returns for each state the probability that the path formula holds + */ + virtual std::vector* checkGlobally(const storm::formula::Globally& formula) const { + // Create "equivalent" temporary eventually formula and check it. + storm::formula::Eventually temporaryEventuallyFormula(new storm::formula::Not(formula.getChild().clone())); + std::vector* result = this->checkEventually(temporaryEventuallyFormula); + + // Now subtract the resulting vector from the constant one vector to obtain final result. + storm::utility::subtractFromConstantOneVector(result); + return result; + } + + /*! + * The check method for a path formula with an Until operator node as root in its formula tree + * + * @param formula The Until path formula to check + * @returns for each state the probability that the path formula holds. + */ + virtual std::vector* checkUntil(const storm::formula::Until& formula) const = 0; + + /*! + * The check method for a path formula with an Instantaneous Reward operator node as root in its + * formula tree + * + * @param formula The Instantaneous Reward formula to check + * @returns for each state the reward that the instantaneous reward yields + */ + virtual std::vector* checkInstantaneousReward(const storm::formula::InstantaneousReward& formula) const = 0; + + /*! + * The check method for a path formula with a Cumulative Reward operator node as root in its + * formula tree + * + * @param formula The Cumulative Reward formula to check + * @returns for each state the reward that the cumulative reward yields + */ + virtual std::vector* checkCumulativeReward(const storm::formula::CumulativeReward& formula) const = 0; + + /*! + * The check method for a path formula with a Reachability Reward operator node as root in its + * formula tree + * + * @param formula The Reachbility Reward formula to check + * @returns for each state the reward that the reachability reward yields + */ + virtual std::vector* checkReachabilityReward(const storm::formula::ReachabilityReward& formula) const = 0; + +protected: + std::stack minimumOperatorStack; + +private: + storm::models::Mdp& model; +}; + +} //namespace modelChecker + +} //namespace storm + +#endif /* STORM_MODELCHECKER_DTMCPRCTLMODELCHECKER_H_ */ diff --git a/src/utility/Vector.h b/src/utility/Vector.h index a69b87e4a..680990e38 100644 --- a/src/utility/Vector.h +++ b/src/utility/Vector.h @@ -53,6 +53,48 @@ void subtractFromConstantOneVector(std::vector* vector) { } } +template +void reduceVectorMin(std::vector const& source, std::vector* target, std::vector const& filter) { + uint_fast64_t currentSourceRow = 0; + uint_fast64_t currentTargetRow = 0; + for (auto it = source->cbegin(); it != source->cend(); ++it, ++currentSourceRow) { + // Check whether we have considered all from rows for the current to row. + if (filter[currentTargetRow + 1] <= currentSourceRow) { + ++currentTargetRow; + (*target)[currentTargetRow] = (*source)[currentSourceRow]; + continue; + } + + + // We have to minimize the value, so only overwrite the current value if the + // value is actually lower. + if (*it < (*target)[currentTargetRow]) { + (*source)[currentTargetRow] = *it; + } + } +} + +template +void reduceVectorMax(std::vector const& source, std::vector* target, std::vector const& filter) { + uint_fast64_t currentSourceRow = 0; + uint_fast64_t currentTargetRow = 0; + for (auto it = source->cbegin(); it != source->cend(); ++it, ++currentSourceRow) { + // Check whether we have considered all from rows for the current to row. + if (filter[currentTargetRow + 1] <= currentSourceRow) { + ++currentTargetRow; + (*target)[currentTargetRow] = (*source)[currentSourceRow]; + continue; + } + + + // We have to maximize the value, so only overwrite the current value if the + // value is actually greater. + if (*it > (*target)[currentTargetRow]) { + (*source)[currentTargetRow] = *it; + } + } +} + } //namespace utility } //namespace storm