From 58cf8118fed93b7e6d48b6ee863cdf3d193bd6b0 Mon Sep 17 00:00:00 2001 From: dehnert Date: Wed, 26 Dec 2012 19:30:07 +0100 Subject: [PATCH] Initial version of reward model checking for DTMCs. Added two convenience operators to PCTL (Eventually and Globally) and added missing reward formulas. --- src/formula/BoundedUntil.h | 2 +- src/formula/CumulativeReward.h | 115 +++++++ src/formula/Eventually.h | 125 ++++++++ src/formula/Formulas.h | 9 +- src/formula/Globally.h | 125 ++++++++ src/formula/InstantaneousReward.h | 115 +++++++ src/formula/ProbabilisticIntervalOperator.h | 4 - src/formula/ProbabilisticOperator.h | 158 ---------- src/formula/ReachabilityReward.h | 122 ++++++++ src/formula/RewardIntervalOperator.h | 171 ++++++++++ src/formula/RewardNoBoundsOperator.h | 110 +++++++ src/modelChecker/DtmcPrctlModelChecker.h | 145 +++++++-- src/modelChecker/GmmxxDtmcPrctlModelChecker.h | 296 ++++++++++++++---- src/models/Dtmc.h | 25 +- src/solver/GraphAnalyzer.h | 25 ++ src/storage/SquareSparseMatrix.h | 32 ++ src/storm.cpp | 22 +- src/utility/CommandLine.cpp | 20 ++ src/utility/CommandLine.h | 24 ++ src/utility/ConstTemplates.h | 53 ++++ src/utility/Vector.h | 19 +- 21 files changed, 1457 insertions(+), 260 deletions(-) create mode 100644 src/formula/CumulativeReward.h create mode 100644 src/formula/Eventually.h create mode 100644 src/formula/Globally.h create mode 100644 src/formula/InstantaneousReward.h delete mode 100644 src/formula/ProbabilisticOperator.h create mode 100644 src/formula/ReachabilityReward.h create mode 100644 src/formula/RewardIntervalOperator.h create mode 100644 src/formula/RewardNoBoundsOperator.h create mode 100644 src/utility/CommandLine.cpp create mode 100644 src/utility/CommandLine.h diff --git a/src/formula/BoundedUntil.h b/src/formula/BoundedUntil.h index 36d5c7fc1..1dd8f255f 100644 --- a/src/formula/BoundedUntil.h +++ b/src/formula/BoundedUntil.h @@ -56,7 +56,7 @@ public: BoundedUntil(PctlStateFormula* left, PctlStateFormula* right, uint_fast64_t bound) { this->left = left; - this->right = right;; + this->right = right; this->bound = bound; } diff --git a/src/formula/CumulativeReward.h b/src/formula/CumulativeReward.h new file mode 100644 index 000000000..cbbf1e18b --- /dev/null +++ b/src/formula/CumulativeReward.h @@ -0,0 +1,115 @@ +/* + * InstantaneousReward.h + * + * Created on: 26.12.2012 + * Author: Christian Dehnert + */ + +#ifndef STORM_FORMULA_CUMULATIVEREWARD_H_ +#define STORM_FORMULA_CUMULATIVEREWARD_H_ + +#include "PctlPathFormula.h" +#include "PctlStateFormula.h" +#include "boost/integer/integer_mask.hpp" +#include + +namespace storm { + +namespace formula { + +/*! + * @brief + * Class for a PCTL (path) formula tree with a Cumulative Reward node as root. + * + * The subtrees are seen as part of the object and deleted with the object + * (this behavior can be prevented by setting them to NULL before deletion) + * + * @see PctlPathFormula + * @see PctlFormula + */ +template +class CumulativeReward : public PctlPathFormula { + +public: + /*! + * Empty constructor + */ + CumulativeReward() { + bound = 0; + } + + /*! + * Constructor + * + * @param bound The time bound of the reward formula + */ + CumulativeReward(uint_fast64_t bound) { + this->bound = bound; + } + + /*! + * Empty destructor. + */ + virtual ~CumulativeReward() { + // Intentionally left empty. + } + + /*! + * @returns the time instance for the instantaneous reward operator + */ + uint_fast64_t getBound() const { + return bound; + } + + /*! + * Sets the the time instance for the instantaneous reward operator + * + * @param bound the new bound. + */ + void setBound(uint_fast64_t bound) { + this->bound = bound; + } + + /*! + * @returns a string representation of the formula + */ + virtual std::string toString() const { + std::string result = "C<="; + result += std::to_string(bound); + return result; + } + + /*! + * Clones the called object. + * + * Performs a "deep copy", i.e. the subtrees of the new object are clones of the original ones + * + * @returns a new BoundedUntil-object that is identical the called object. + */ + virtual PctlPathFormula* clone() const { + return new CumulativeReward(bound); + } + + + /*! + * Calls the model checker to check this formula. + * Needed to infer the correct type of formula class. + * + * @note This function should only be called in a generic check function of a model checker class. For other uses, + * the methods of the model checker should be used. + * + * @returns A vector indicating the probability that the formula holds for each state. + */ + virtual std::vector *check(const storm::modelChecker::DtmcPrctlModelChecker& modelChecker) const { + return modelChecker.checkCumulativeReward(*this); + } + +private: + uint_fast64_t bound; +}; + +} //namespace formula + +} //namespace storm + +#endif /* STORM_FORMULA_INSTANTANEOUSREWARD_H_ */ diff --git a/src/formula/Eventually.h b/src/formula/Eventually.h new file mode 100644 index 000000000..1a781fbb7 --- /dev/null +++ b/src/formula/Eventually.h @@ -0,0 +1,125 @@ +/* + * Next.h + * + * Created on: 26.12.2012 + * Author: Christian Dehnert + */ + +#ifndef STORM_FORMULA_EVENTUALLY_H_ +#define STORM_FORMULA_EVENTUALLY_H_ + +#include "PctlPathFormula.h" +#include "PctlStateFormula.h" + +namespace storm { + +namespace formula { + +/*! + * @brief + * Class for a PCTL (path) formula tree with an Eventually node as root. + * + * Has one PCTL state formula as sub formula/tree. + * + * @par Semantics + * The formula holds iff eventually \e child holds. + * + * The subtree is seen as part of the object and deleted with the object + * (this behavior can be prevented by setting them to nullptr before deletion) + * + * @see PctlPathFormula + * @see PctlFormula + */ +template +class Eventually : public PctlPathFormula { + +public: + /*! + * Empty constructor + */ + Eventually() { + this->child = nullptr; + } + + /*! + * Constructor + * + * @param child The child node + */ + Eventually(PctlStateFormula* child) { + this->child = child; + } + + /*! + * Constructor. + * + * Also deletes the subtree. + * (this behaviour can be prevented by setting the subtrees to nullptr before deletion) + */ + virtual ~Eventually() { + if (child != nullptr) { + delete child; + } + } + + /*! + * @returns the child node + */ + const PctlStateFormula& getChild() const { + return *child; + } + + /*! + * Sets the subtree + * @param child the new child node + */ + void setChild(PctlStateFormula* child) { + this->child = child; + } + + /*! + * @returns a string representation of the formula + */ + virtual std::string toString() const { + std::string result = "F "; + result += child->toString(); + return result; + } + + /*! + * Clones the called object. + * + * Performs a "deep copy", i.e. the subtrees of the new object are clones of the original ones + * + * @returns a new BoundedUntil-object that is identical the called object. + */ + virtual PctlPathFormula* clone() const { + Eventually* result = new Eventually(); + if (child != nullptr) { + result->setChild(child); + } + return result; + } + + /*! + * Calls the model checker to check this formula. + * Needed to infer the correct type of formula class. + * + * @note This function should only be called in a generic check function of a model checker class. For other uses, + * the methods of the model checker should be used. + * + * @returns A vector indicating the probability that the formula holds for each state. + */ + virtual std::vector *check(const storm::modelChecker::DtmcPrctlModelChecker& modelChecker) const { + return modelChecker.checkEventually(*this); + } + +private: + PctlStateFormula* child; +}; + +} //namespace formula + +} //namespace storm + +#endif /* STORM_FORMULA_EVENTUALLY_H_ */ diff --git a/src/formula/Formulas.h b/src/formula/Formulas.h index 3cd6ce94c..b93952b9b 100644 --- a/src/formula/Formulas.h +++ b/src/formula/Formulas.h @@ -17,9 +17,16 @@ #include "PctlFormula.h" #include "PctlPathFormula.h" #include "PctlStateFormula.h" -#include "ProbabilisticOperator.h" #include "ProbabilisticNoBoundsOperator.h" #include "ProbabilisticIntervalOperator.h" #include "Until.h" +#include "Eventually.h" +#include "Globally.h" + +#include "InstantaneousReward.h" +#include "CumulativeReward.h" +#include "ReachabilityReward.h" +#include "RewardIntervalOperator.h" +#include "RewardNoBoundsOperator.h" #endif /* STORM_FORMULA_FORMULAS_H_ */ diff --git a/src/formula/Globally.h b/src/formula/Globally.h new file mode 100644 index 000000000..02a18b884 --- /dev/null +++ b/src/formula/Globally.h @@ -0,0 +1,125 @@ +/* + * Next.h + * + * Created on: 26.12.2012 + * Author: Christian Dehnert + */ + +#ifndef STORM_FORMULA_GLOBALLY_H_ +#define STORM_FORMULA_GLOBALLY_H_ + +#include "PctlPathFormula.h" +#include "PctlStateFormula.h" + +namespace storm { + +namespace formula { + +/*! + * @brief + * Class for a PCTL (path) formula tree with a Globally node as root. + * + * Has one PCTL state formula as sub formula/tree. + * + * @par Semantics + * The formula holds iff globally \e child holds. + * + * The subtree is seen as part of the object and deleted with the object + * (this behavior can be prevented by setting them to nullptr before deletion) + * + * @see PctlPathFormula + * @see PctlFormula + */ +template +class Globally : public PctlPathFormula { + +public: + /*! + * Empty constructor + */ + Globally() { + this->child = nullptr; + } + + /*! + * Constructor + * + * @param child The child node + */ + Globally(PctlStateFormula* child) { + this->child = child; + } + + /*! + * Constructor. + * + * Also deletes the subtree. + * (this behaviour can be prevented by setting the subtrees to nullptr before deletion) + */ + virtual ~Globally() { + if (child != nullptr) { + delete child; + } + } + + /*! + * @returns the child node + */ + const PctlStateFormula& getChild() const { + return *child; + } + + /*! + * Sets the subtree + * @param child the new child node + */ + void setChild(PctlStateFormula* child) { + this->child = child; + } + + /*! + * @returns a string representation of the formula + */ + virtual std::string toString() const { + std::string result = "G "; + result += child->toString(); + return result; + } + + /*! + * Clones the called object. + * + * Performs a "deep copy", i.e. the subtrees of the new object are clones of the original ones + * + * @returns a new BoundedUntil-object that is identical the called object. + */ + virtual PctlPathFormula* clone() const { + Next* result = new Next(); + if (child != nullptr) { + result->setChild(child); + } + return result; + } + + /*! + * Calls the model checker to check this formula. + * Needed to infer the correct type of formula class. + * + * @note This function should only be called in a generic check function of a model checker class. For other uses, + * the methods of the model checker should be used. + * + * @returns A vector indicating the probability that the formula holds for each state. + */ + virtual std::vector *check(const storm::modelChecker::DtmcPrctlModelChecker& modelChecker) const { + return modelChecker.checkGlobally(*this); + } + +private: + PctlStateFormula* child; +}; + +} //namespace formula + +} //namespace storm + +#endif /* STORM_FORMULA_GLOBALLY_H_ */ diff --git a/src/formula/InstantaneousReward.h b/src/formula/InstantaneousReward.h new file mode 100644 index 000000000..d6067d344 --- /dev/null +++ b/src/formula/InstantaneousReward.h @@ -0,0 +1,115 @@ +/* + * InstantaneousReward.h + * + * Created on: 26.12.2012 + * Author: Christian Dehnert + */ + +#ifndef STORM_FORMULA_INSTANTANEOUSREWARD_H_ +#define STORM_FORMULA_INSTANTANEOUSREWARD_H_ + +#include "PctlPathFormula.h" +#include "PctlStateFormula.h" +#include "boost/integer/integer_mask.hpp" +#include + +namespace storm { + +namespace formula { + +/*! + * @brief + * Class for a PCTL (path) formula tree with a Instantaneous Reward node as root. + * + * The subtrees are seen as part of the object and deleted with the object + * (this behavior can be prevented by setting them to NULL before deletion) + * + * @see PctlPathFormula + * @see PctlFormula + */ +template +class InstantaneousReward : public PctlPathFormula { + +public: + /*! + * Empty constructor + */ + InstantaneousReward() { + bound = 0; + } + + /*! + * Constructor + * + * @param bound The time instance of the reward formula + */ + InstantaneousReward(uint_fast64_t bound) { + this->bound = bound; + } + + /*! + * Empty destructor. + */ + virtual ~InstantaneousReward() { + // Intentionally left empty. + } + + /*! + * @returns the time instance for the instantaneous reward operator + */ + uint_fast64_t getBound() const { + return bound; + } + + /*! + * Sets the the time instance for the instantaneous reward operator + * + * @param bound the new bound. + */ + void setBound(uint_fast64_t bound) { + this->bound = bound; + } + + /*! + * @returns a string representation of the formula + */ + virtual std::string toString() const { + std::string result = "I="; + result += std::to_string(bound); + return result; + } + + /*! + * Clones the called object. + * + * Performs a "deep copy", i.e. the subtrees of the new object are clones of the original ones + * + * @returns a new BoundedUntil-object that is identical the called object. + */ + virtual PctlPathFormula* clone() const { + return new InstantaneousReward(bound); + } + + + /*! + * Calls the model checker to check this formula. + * Needed to infer the correct type of formula class. + * + * @note This function should only be called in a generic check function of a model checker class. For other uses, + * the methods of the model checker should be used. + * + * @returns A vector indicating the probability that the formula holds for each state. + */ + virtual std::vector *check(const storm::modelChecker::DtmcPrctlModelChecker& modelChecker) const { + return modelChecker.checkInstantaneousReward(*this); + } + +private: + uint_fast64_t bound; +}; + +} //namespace formula + +} //namespace storm + +#endif /* STORM_FORMULA_INSTANTANEOUSREWARD_H_ */ diff --git a/src/formula/ProbabilisticIntervalOperator.h b/src/formula/ProbabilisticIntervalOperator.h index 2ca2b9bbb..175167292 100644 --- a/src/formula/ProbabilisticIntervalOperator.h +++ b/src/formula/ProbabilisticIntervalOperator.h @@ -21,10 +21,6 @@ namespace formula { * Class for a PCTL formula tree with a P (probablistic) operator node over a probability interval * as root. * - * If the probability interval consist just of one single value (i.e. it is [x,x] for some - * real number x), the class ProbabilisticOperator should be used instead. - * - * * Has one PCTL path formula as sub formula/tree. * * @par Semantics diff --git a/src/formula/ProbabilisticOperator.h b/src/formula/ProbabilisticOperator.h deleted file mode 100644 index 118bf1cdb..000000000 --- a/src/formula/ProbabilisticOperator.h +++ /dev/null @@ -1,158 +0,0 @@ -/* - * ProbabilisticOperator.h - * - * Created on: 07.12.2012 - * Author: Thomas Heinemann - */ - -#ifndef STORM_FORMULA_PROBABILISTICOPERATOR_H_ -#define STORM_FORMULA_PROBABILISTICOPERATOR_H_ - -#include "PctlStateFormula.h" - -namespace storm { -namespace formula { - -/*! - * @brief - * Class for a PCTL formula tree with a P (probablistic) operator node over a single real valued - * probability as root. - * - * If the probability interval consist just of one single value (i.e. it is [x,x] for some - * real number x), the class ProbabilisticOperator should be used instead. - * - * - * Has one PCTL path formula as sub formula/tree. - * - * @par Semantics - * The formula holds iff the probability that the path formula holds is equal to the probablility - * specified in this operator - * - * The subtree is seen as part of the object and deleted with it - * (this behavior can be prevented by setting them to NULL before deletion) - * - * - * @see PctlStateFormula - * @see PctlPathFormula - * @see ProbabilisticIntervalOperator - * @see ProbabilisticNoBoundsOperator - * @see PctlFormula - */ -template -class ProbabilisticOperator : public storm::formula::PctlStateFormula { -public: - /*! - * Empty constructor - */ - ProbabilisticOperator() { - this->pathFormula = NULL; - } - - /*! - * Constructor - * - * @param bound The expected value for path formulas - * @param pathFormula The child node - */ - ProbabilisticOperator(T bound, PctlPathFormula* pathFormula) { - this->bound = bound; - this->pathFormula = *pathFormula; - } - - /*! - * Destructor - * - * The subtree is deleted with the object - * (this behavior can be prevented by setting them to NULL before deletion) - */ - virtual ~ProbabilisticOperator() { - if (pathFormula != NULL) { - delete pathFormula; - } - } - - /*! - * @returns the child node (representation of a PCTL path formula) - */ - const PctlPathFormula& getPathFormula () const { - return *pathFormula; - } - - /*! - * @returns the bound for the probability - */ - const T& getBound() const { - return bound; - } - - /*! - * Sets the child node - * - * @param pathFormula the path formula that becomes the new child node - */ - void setPathFormula(PctlPathFormula* pathFormula) { - this->pathFormula = pathFormula; - } - - /*! - * Sets the expected probability that the path formula holds. - * - * @param bound The bound for the probability - */ - void setBound(T bound) { - this->bound = bound; - } - - /*! - * Clones the called object. - * - * Performs a "deep copy", i.e. the subtrees of the new object are clones of the original ones - * - * @returns a new ProbabilisticOperator-object that is identical to the called object. - */ - virtual PctlStateFormula* clone() const { - ProbabilisticOperator* result = new ProbabilisticOperator(); - result->setBound(bound); - if (pathFormula != NULL) { - result->setPathFormula(pathFormula->clone()); - } - return result; - } - - /*! - * Calls the model checker to check this formula. - * Needed to infer the correct type of formula class. - * - * @note This function should only be called in a generic check function of a model checker - * class. For other uses, the methods of the model checker should be used. - * - * @returns A bit vector indicating all states that satisfy the formula represented by the - * called object. - */ - virtual storm::storage::BitVector *check( - const storm::modelChecker::DtmcPrctlModelChecker& modelChecker) const { - return modelChecker.checkProbabilisticOperator(*this); - } - - /*! - * Returns a string representation of this PctlStateFormula - * - * @returns a string representation of this PctlStateFormula - */ - virtual std::string toString() const { - std::string result = " P="; - result += std::to_string(bound); - result += " ("; - result += pathFormula->toString(); - result += ")"; - return result; - } -private: - T bound; - PctlPathFormula* pathFormula; -}; - -} /* namespace formula */ -} /* namespace storm */ - -#endif /* STORM_FORMULA_PROBABILISTICOPERATOR_H_ */ diff --git a/src/formula/ReachabilityReward.h b/src/formula/ReachabilityReward.h new file mode 100644 index 000000000..9f5d6ee43 --- /dev/null +++ b/src/formula/ReachabilityReward.h @@ -0,0 +1,122 @@ +/* + * Next.h + * + * Created on: 26.12.2012 + * Author: Christian Dehnert + */ + +#ifndef STORM_FORMULA_REACHABILITYREWARD_H_ +#define STORM_FORMULA_REACHABILITYREWARD_H_ + +#include "PctlPathFormula.h" +#include "PctlStateFormula.h" + +namespace storm { + +namespace formula { + +/*! + * @brief + * Class for a PCTL (path) formula tree with an Reachability Reward node as root. + * + * Has one PCTL state formula as sub formula/tree. + * + * The subtree is seen as part of the object and deleted with the object + * (this behavior can be prevented by setting them to nullptr before deletion) + * + * @see PctlPathFormula + * @see PctlFormula + */ +template +class ReachabilityReward : public PctlPathFormula { + +public: + /*! + * Empty constructor + */ + ReachabilityReward() { + this->child = nullptr; + } + + /*! + * Constructor + * + * @param child The child node + */ + ReachabilityReward(PctlStateFormula* child) { + this->child = child; + } + + /*! + * Constructor. + * + * Also deletes the subtree. + * (this behaviour can be prevented by setting the subtrees to nullptr before deletion) + */ + virtual ~ReachabilityReward() { + if (child != nullptr) { + delete child; + } + } + + /*! + * @returns the child node + */ + const PctlStateFormula& getChild() const { + return *child; + } + + /*! + * Sets the subtree + * @param child the new child node + */ + void setChild(PctlStateFormula* child) { + this->child = child; + } + + /*! + * @returns a string representation of the formula + */ + virtual std::string toString() const { + std::string result = "F "; + result += child->toString(); + return result; + } + + /*! + * Clones the called object. + * + * Performs a "deep copy", i.e. the subtrees of the new object are clones of the original ones + * + * @returns a new BoundedUntil-object that is identical the called object. + */ + virtual PctlPathFormula* clone() const { + ReachabilityReward* result = new ReachabilityReward(); + if (child != nullptr) { + result->setChild(child); + } + return result; + } + + /*! + * Calls the model checker to check this formula. + * Needed to infer the correct type of formula class. + * + * @note This function should only be called in a generic check function of a model checker class. For other uses, + * the methods of the model checker should be used. + * + * @returns A vector indicating the probability that the formula holds for each state. + */ + virtual std::vector *check(const storm::modelChecker::DtmcPrctlModelChecker& modelChecker) const { + return modelChecker.checkReachabilityReward(*this); + } + +private: + PctlStateFormula* child; +}; + +} //namespace formula + +} //namespace storm + +#endif /* STORM_FORMULA_REACHABILITYREWARD_H_ */ diff --git a/src/formula/RewardIntervalOperator.h b/src/formula/RewardIntervalOperator.h new file mode 100644 index 000000000..d18e3b7d9 --- /dev/null +++ b/src/formula/RewardIntervalOperator.h @@ -0,0 +1,171 @@ +/* + * ProbabilisticOperator.h + * + * Created on: 19.10.2012 + * Author: Thomas Heinemann + */ + +#ifndef STORM_FORMULA_REWARDINTERVALOPERATOR_H_ +#define STORM_FORMULA_REWARDINTERVALOPERATOR_H_ + +#include "PctlStateFormula.h" +#include "PctlPathFormula.h" +#include "utility/ConstTemplates.h" + +namespace storm { + +namespace formula { + +/*! + * @brief + * Class for a PCTL formula tree with a R (reward) operator node over a reward interval as root. + * + * Has a reward path formula as sub formula/tree. + * + * @par Semantics + * The formula holds iff the reward of the reward path formula is inside the bounds + * specified in this operator + * + * The subtree is seen as part of the object and deleted with it + * (this behavior can be prevented by setting them to NULL before deletion) + * + * + * @see PctlStateFormula + * @see PctlPathFormula + * @see ProbabilisticOperator + * @see ProbabilisticNoBoundsOperator + * @see PctlFormula + */ +template +class RewardIntervalOperator : public PctlStateFormula { + +public: + /*! + * Empty constructor + */ + RewardIntervalOperator() { + upper = storm::utility::constGetZero(); + lower = storm::utility::constGetZero(); + pathFormula = nullptr; + } + + /*! + * Constructor + * + * @param lowerBound The lower bound for the probability + * @param upperBound The upper bound for the probability + * @param pathFormula The child node + */ + RewardIntervalOperator(T lowerBound, T upperBound, PctlPathFormula& pathFormula) { + this->lower = lowerBound; + this->upper = upperBound; + this->pathFormula = &pathFormula; + } + + /*! + * Destructor + * + * The subtree is deleted with the object + * (this behavior can be prevented by setting them to NULL before deletion) + */ + virtual ~RewardIntervalOperator() { + if (pathFormula != nullptr) { + delete pathFormula; + } + } + + /*! + * @returns the child node (representation of a PCTL path formula) + */ + const PctlPathFormula& getPathFormula () const { + return *pathFormula; + } + + /*! + * @returns the lower bound for the probability + */ + const T& getLowerBound() const { + return lower; + } + + /*! + * @returns the upper bound for the probability + */ + const T& getUpperBound() const { + return upper; + } + + /*! + * Sets the child node + * + * @param pathFormula the path formula that becomes the new child node + */ + void setPathFormula(PctlPathFormula* pathFormula) { + this->pathFormula = pathFormula; + } + + /*! + * Sets the interval in which the probability that the path formula holds may lie in. + * + * @param lowerBound The lower bound for the probability + * @param upperBound The upper bound for the probability + */ + void setInterval(T lowerBound, T upperBound) { + this->lower = lowerBound; + this->upper = upperBound; + } + + /*! + * @returns a string representation of the formula + */ + virtual std::string toString() const { + std::string result = "R ["; + result += std::to_string(lower); + result += ", "; + result += std::to_string(upper); + result += "] ["; + result += pathFormula->toString(); + result += "]"; + return result; + } + + /*! + * Clones the called object. + * + * Performs a "deep copy", i.e. the subtrees of the new object are clones of the original ones + * + * @returns a new AND-object that is identical the called object. + */ + virtual PctlStateFormula* clone() const { + RewardIntervalOperator* result = new RewardIntervalOperator(); + result->setInterval(lower, upper); + if (pathFormula != nullptr) { + result->setPathFormula(pathFormula->clone()); + } + return result; + } + + /*! + * Calls the model checker to check this formula. + * Needed to infer the correct type of formula class. + * + * @note This function should only be called in a generic check function of a model checker class. For other uses, + * the methods of the model checker should be used. + * + * @returns A bit vector indicating all states that satisfy the formula represented by the called object. + */ + virtual storm::storage::BitVector *check(const storm::modelChecker::DtmcPrctlModelChecker& modelChecker) const { + return modelChecker.checkRewardIntervalOperator(*this); + } + +private: + T lower; + T upper; + PctlPathFormula* pathFormula; +}; + +} //namespace formula + +} //namespace storm + +#endif /* STORM_FORMULA_REWARDINTERVALOPERATOR_H_ */ diff --git a/src/formula/RewardNoBoundsOperator.h b/src/formula/RewardNoBoundsOperator.h new file mode 100644 index 000000000..b1acd4acf --- /dev/null +++ b/src/formula/RewardNoBoundsOperator.h @@ -0,0 +1,110 @@ +/* + * RewardNoBoundsOperator.h + * + * Created on: 12.12.2012 + * Author: thomas + */ + +#ifndef STORM_FORMULA_REWARDNOBOUNDSOPERATOR_H_ +#define STORM_FORMULA_REWARDNOBOUNDSOPERATOR_H_ + +#include "PctlFormula.h" +#include "PctlPathFormula.h" + +namespace storm { + +namespace formula { + +/*! + * @brief + * Class for a PCTL formula tree with a R (reward) operator without declaration of reward values + * as root. + * + * Checking a formula with this operator as root returns the reward for the reward path formula for + * each state + * + * Has one PCTL path formula as sub formula/tree. + * + * @note + * This class is a hybrid of a state and path formula, and may only appear as the outermost operator. + * Hence, it is seen as neither a state nor a path formula, but is directly derived from PctlFormula. + * + * @note + * This class does not contain a check() method like the other formula classes. + * The check method should only be called by the model checker to infer the correct check function for sub + * formulas. As this operator can only appear at the root, the method is not useful here. + * Use the checkRewardNoBoundsOperator method from the DtmcPrctlModelChecker class instead. + * + * The subtree is seen as part of the object and deleted with it + * (this behavior can be prevented by setting them to NULL before deletion) + * + * + * @see PctlStateFormula + * @see PctlPathFormula + * @see ProbabilisticOperator + * @see ProbabilisticIntervalOperator + * @see PctlFormula + */ +template +class RewardNoBoundsOperator: public storm::formula::PctlFormula { +public: + /*! + * Empty constructor + */ + RewardNoBoundsOperator() { + this->pathFormula = nullptr; + } + + /*! + * Constructor + * + * @param pathFormula The child node. + */ + RewardNoBoundsOperator(PctlPathFormula* pathFormula) { + this->pathFormula = pathFormula; + } + + /*! + * Destructor + */ + virtual ~RewardNoBoundsOperator() { + if (pathFormula != nullptr) { + delete pathFormula; + } + } + + /*! + * @returns the child node (representation of a PCTL path formula) + */ + const PctlPathFormula& getPathFormula () const { + return *pathFormula; + } + + /*! + * Sets the child node + * + * @param pathFormula the path formula that becomes the new child node + */ + void setPathFormula(PctlPathFormula* pathFormula) { + this->pathFormula = pathFormula; + } + + /*! + * @returns a string representation of the formula + */ + virtual std::string toString() const { + std::string result = " R=? ["; + result += pathFormula->toString(); + result += "]"; + return result; + } + +private: + PctlPathFormula* pathFormula; +}; + +} /* namespace formula */ + +} /* namespace storm */ + +#endif /* STORM_FORMULA_REWARDNOBOUNDSOPERATOR_H_ */ diff --git a/src/modelChecker/DtmcPrctlModelChecker.h b/src/modelChecker/DtmcPrctlModelChecker.h index 5f48a8e8d..3a5c50a41 100644 --- a/src/modelChecker/DtmcPrctlModelChecker.h +++ b/src/modelChecker/DtmcPrctlModelChecker.h @@ -98,13 +98,19 @@ public: * states. * @param stateFormula The formula to be checked. */ - void check(const storm::formula::PctlStateFormula& stateFormula) { - LOG4CPLUS_INFO(logger, "Model checking formula " << stateFormula.toString()); + void check(const storm::formula::PctlStateFormula& 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 = 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; } + std::cout << std::endl; + storm::utility::printSeparationLine(std::cout); delete result; } @@ -113,9 +119,10 @@ public: * (probability) for all initial states. * @param probabilisticNoBoundsFormula The formula to be checked. */ - void check(const storm::formula::ProbabilisticNoBoundsOperator& probabilisticNoBoundsFormula) { - LOG4CPLUS_INFO(logger, "Model checking formula " << probabilisticNoBoundsFormula.toString()); - std::cout << "Model checking formula: " << probabilisticNoBoundsFormula.toString() << std::endl; + void check(const storm::formula::ProbabilisticNoBoundsOperator& probabilisticNoBoundsFormula) const { + std::cout << std::endl; + LOG4CPLUS_INFO(logger, "Model checking formula\t" << probabilisticNoBoundsFormula.toString()); + std::cout << "Model checking formula:\t" << probabilisticNoBoundsFormula.toString() << std::endl; std::vector* result = checkProbabilisticNoBoundsOperator(probabilisticNoBoundsFormula); LOG4CPLUS_INFO(logger, "Result for initial states:"); std::cout << "Result for initial states:" << std::endl; @@ -123,6 +130,29 @@ public: LOG4CPLUS_INFO(logger, "\t" << initialState << ": " << (*result)[initialState]); std::cout << "\t" << initialState << ": " << (*result)[initialState] << std::endl; } + std::cout << std::endl; + storm::utility::printSeparationLine(std::cout); + delete result; + } + + /*! + * Checks the given reward operator (with no bound) on the DTMC and prints the result + * (reward values) for all initial states. + * @param rewardNoBoundsFormula The formula to be checked. + */ + void check(const storm::formula::RewardNoBoundsOperator& rewardNoBoundsFormula) { + std::cout << std::endl; + LOG4CPLUS_INFO(logger, "Model checking formula\t" << rewardNoBoundsFormula.toString()); + std::cout << "Model checking formula:\t" << rewardNoBoundsFormula.toString() << std::endl; + std::vector* result = checkRewardNoBoundsOperator(rewardNoBoundsFormula); + 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; + } + std::cout << std::endl; + storm::utility::printSeparationLine(std::cout); delete result; } @@ -193,38 +223,45 @@ public: } /*! - * The check method for a state formula with a probabilistic operator node as root in its - * formula tree + * The check method for a state formula with a probabilistic interval 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* checkProbabilisticOperator( - const storm::formula::ProbabilisticOperator& formula) const { + storm::storage::BitVector* checkProbabilisticIntervalOperator( + const storm::formula::ProbabilisticIntervalOperator& formula) const { + // First, we need to compute the probability for satisfying the path formula for each state. std::vector* probabilisticResult = 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()); - Type bound = formula.getBound(); + + // Now, we can compute which states meet the bound specified in this operator, i.e. + // lie in the interval that was given along with this operator, and set the corresponding bits + // to true in the resulting vector. + Type lower = formula.getLowerBound(); + Type upper = formula.getUpperBound(); for (uint_fast64_t i = 0; i < this->getModel().getNumberOfStates(); ++i) { - if ((*probabilisticResult)[i] == bound) result->set(i, true); + if ((*probabilisticResult)[i] >= lower && (*probabilisticResult)[i] <= upper) result->set(i, true); } + // Delete the probabilities computed for the states and return result. delete probabilisticResult; - return result; } /*! - * The check method for a state formula with a probabilistic interval operator node as root in + * The check method for a state formula with a reward interval 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* checkProbabilisticIntervalOperator( - const storm::formula::ProbabilisticIntervalOperator& formula) const { + storm::storage::BitVector* checkRewardIntervalOperator( + const storm::formula::RewardIntervalOperator& formula) const { // First, we need to compute the probability for satisfying the path formula for each state. - std::vector* probabilisticResult = this->checkPathFormula(formula.getPathFormula()); + std::vector* rewardResult = 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()); @@ -235,15 +272,14 @@ public: Type lower = formula.getLowerBound(); Type upper = formula.getUpperBound(); for (uint_fast64_t i = 0; i < this->getModel().getNumberOfStates(); ++i) { - if ((*probabilisticResult)[i] >= lower && (*probabilisticResult)[i] <= upper) result->set(i, true); + if ((*rewardResult)[i] >= lower && (*rewardResult)[i] <= upper) result->set(i, true); } - // Delete the probabilities computed for the states and return result. - delete probabilisticResult; + // Delete the reward values computed for the states and return result. + delete rewardResult; return result; } - /*! * The check method for a state formula with a probabilistic operator node without bounds as root * in its formula tree @@ -256,6 +292,18 @@ public: return formula.getPathFormula().check(*this); } + /*! + * The check method for a state formula with a reward 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* checkRewardNoBoundsOperator( + const storm::formula::RewardNoBoundsOperator& formula) const { + return formula.getPathFormula().check(*this); + } + /*! * The check method for a path formula; Will infer the actual type of formula and delegate it * to the specialized method @@ -283,6 +331,36 @@ public: */ virtual std::vector* checkNext(const storm::formula::Next& formula) const = 0; + /*! + * 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()); + std::vector* result = this->checkUntil(temporaryUntilFormula); + return result; + } + + /*! + * 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 * @@ -291,6 +369,33 @@ public: */ 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; + private: storm::models::Dtmc& model; }; diff --git a/src/modelChecker/GmmxxDtmcPrctlModelChecker.h b/src/modelChecker/GmmxxDtmcPrctlModelChecker.h index ab7bce9b9..17bdfd8e4 100644 --- a/src/modelChecker/GmmxxDtmcPrctlModelChecker.h +++ b/src/modelChecker/GmmxxDtmcPrctlModelChecker.h @@ -17,6 +17,7 @@ #include "src/utility/ConstTemplates.h" #include "src/utility/Settings.h" #include "src/adapters/GmmxxAdapter.h" +#include "src/exceptions/InvalidArgumentException.h" #include "gmm/gmm_matrix.h" #include "gmm/gmm_iter_solvers.h" @@ -71,6 +72,7 @@ public: delete tmpResult; // Delete intermediate results and return result. + delete gmmxxMatrix; delete leftStates; delete rightStates; return result; @@ -122,7 +124,7 @@ public: storm::storage::BitVector maybeStates = ~(notExistsPhiUntilPsiStates | alwaysPhiUntilPsiStates); LOG4CPLUS_INFO(logger, "Found " << maybeStates.getNumberOfSetBits() << " 'maybe' states."); - // Create resulting vector and set values accordingly. + // 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. @@ -147,84 +149,176 @@ public: std::vector b(maybeStates.getNumberOfSetBits()); this->getModel().getTransitionProbabilityMatrix()->getConstrainedRowCountVector(maybeStates, alwaysPhiUntilPsiStates, &b); - // Get the settings object to customize linear solving. - storm::settings::Settings* s = storm::settings::instance(); + // Solve the corresponding system of linear equations. + this->solveLinearEquationSystem(*gmmxxMatrix, x, b); - // Prepare an iteration object that determines the accuracy, maximum number of iterations - // and the like. - gmm::iteration iter(s->get("precision"), 0, s->get("lemaxiter")); + // Set values of resulting vector according to result. + storm::utility::setVectorValues(result, maybeStates, x); - // Now do the actual solving. - LOG4CPLUS_INFO(logger, "Starting iterative solver."); - const std::string& precond = s->getString("precond"); - if (precond == "ilu") { - LOG4CPLUS_INFO(logger, "Using ILU preconditioner."); - } else if (precond == "diagonal") { - LOG4CPLUS_INFO(logger, "Using diagonal preconditioner."); - } else if (precond == "ildlt") { - LOG4CPLUS_INFO(logger, "Using ILDLT preconditioner."); - } else if (precond == "none") { - LOG4CPLUS_INFO(logger, "Using no preconditioner."); - } + // Delete temporary matrix. + delete gmmxxMatrix; + } - if (s->getString("lemethod") == "bicgstab") { - LOG4CPLUS_INFO(logger, "Using BiCGStab method."); - if (precond == "ilu") { - gmm::bicgstab(*gmmxxMatrix, x, b, gmm::ilu_precond>(*gmmxxMatrix), iter); - } else if (precond == "diagonal") { - gmm::bicgstab(*gmmxxMatrix, x, b, gmm::diagonal_precond>(*gmmxxMatrix), iter); - } else if (precond == "ildlt") { - gmm::bicgstab(*gmmxxMatrix, x, b, gmm::ildlt_precond>(*gmmxxMatrix), iter); - } else if (precond == "none") { - gmm::bicgstab(*gmmxxMatrix, x, b, gmm::identity_matrix(), iter); - } - // FIXME: gmres has been disabled, because it triggers gmm++ compilation errors - /* } else if (s->getString("lemethod").compare("gmres") == 0) { - LOG4CPLUS_INFO(logger, "Using GMRES method."); - if (precond.compare("ilu")) { - gmm::gmres(*gmmxxMatrix, x, b, gmm::ilu_precond>(*gmmxxMatrix), s->get("restart"), iter); - } else if (precond == "diagonal") { - gmm::gmres(*gmmxxMatrix, x, b, gmm::diagonal_precond>(*gmmxxMatrix), s->get("restart"), iter); - } else if (precond == "ildlt") { - gmm::gmres(*gmmxxMatrix, x, b, gmm::ildlt_precond>(*gmmxxMatrix), s->get("restart"), iter); - } else if (precond == "none") { - gmm::gmres(*gmmxxMatrix, x, b, gmm::identity_matrix(), s->get("restart"), iter); - } */ - } else if (s->getString("lemethod") == "qmr") { - LOG4CPLUS_INFO(logger, "Using QMR method."); - if (precond == "ilu") { - gmm::qmr(*gmmxxMatrix, x, b, gmm::ilu_precond>(*gmmxxMatrix), iter); - } /* FIXME: The following line throws a warning as there should be brackets around such a construction - * TBH, I don't understand it completely (why the comparison with 0?), so I don't know how to fix it - * (Thomas Heinemann, 2012-12-21) - */ - else if (precond == "diagonal") { - gmm::qmr(*gmmxxMatrix, x, b, gmm::diagonal_precond>(*gmmxxMatrix), iter); - } else if (precond == "ildlt") { - gmm::qmr(*gmmxxMatrix, x, b, gmm::ildlt_precond>(*gmmxxMatrix), iter); - } else if (precond == "none") { - gmm::qmr(*gmmxxMatrix, x, b, gmm::identity_matrix(), iter); - } + // Set values of resulting vector that are known exactly. + storm::utility::setVectorValues(result, notExistsPhiUntilPsiStates, storm::utility::constGetZero()); + storm::utility::setVectorValues(result, alwaysPhiUntilPsiStates, 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::InvalidArgumentException() << "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().getTransitionProbabilityMatrix()); + + // Initialize result to state rewards of the model. + std::vector* result = new std::vector(*this->getModel().getStateRewards()); + + // 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::InvalidArgumentException() << "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().getTransitionProbabilityMatrix()); + + // 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().getTransitionProbabilityMatrix()->getPointwiseProductRowSumVector(*this->getModel().getTransitionRewardMatrix()); + if (this->getModel().hasStateRewards()) { + gmm::add(*this->getModel().getStateRewards(), *totalRewardVector); } + } else { + totalRewardVector = new std::vector(*this->getModel().getStateRewards()); + } + + std::vector* result = new std::vector(this->getModel().getNumberOfStates()); - // Check if the solver converged and issue a warning otherwise. - if (iter.converged()) { - LOG4CPLUS_INFO(logger, "Iterative solver converged after " << iter.get_iteration() << " iterations."); + // 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."); + throw storm::exceptions::InvalidArgumentException() << "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::solver::GraphAnalyzer::getAlwaysPhiUntilPsiStates(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; + if (maybeStates.getNumberOfSetBits() > 0) { + // Now we can eliminate the rows and columns from the original transition probability matrix. + storm::storage::SquareSparseMatrix* submatrix = this->getModel().getTransitionProbabilityMatrix()->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(maybeStates.getNumberOfSetBits(), storm::utility::constGetOne()); + + // Prepare the right-hand side of the equation system. + std::vector* b = new std::vector(maybeStates.getNumberOfSetBits()); + 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().getTransitionProbabilityMatrix()->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(maybeStates.getNumberOfSetBits()); + storm::utility::setVectorValues(subStateRewards, maybeStates, *this->getModel().getStateRewards()); + gmm::add(*subStateRewards, *b); + delete subStateRewards; + } } else { - LOG4CPLUS_WARN(logger, "Iterative solver did not converge."); + // 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().getStateRewards()); } + // 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. + // Delete temporary matrix and right-hand side. delete gmmxxMatrix; + delete b; } // Set values of resulting vector that are known exactly. - storm::utility::setVectorValues(result, notExistsPhiUntilPsiStates, storm::utility::constGetZero()); - storm::utility::setVectorValues(result, alwaysPhiUntilPsiStates, storm::utility::constGetOne()); + 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; } @@ -274,6 +368,82 @@ public: throw exceptions::InvalidSettingsException() << "Argument " << preconditioner << " for option 'precond' is invalid."; } } + +private: + /*! + * Solves the linear equation system Ax=b with the given parameters. + * + * @param A The matrix A specifying the coefficients of the linear 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(); + + // Prepare an iteration object that determines the accuracy, maximum number of iterations + // and the like. + gmm::iteration iter(s->get("precision"), 0, s->get("lemaxiter")); + + // Now do the actual solving. + LOG4CPLUS_INFO(logger, "Starting iterative solver."); + const std::string& precond = s->getString("precond"); + if (precond == "ilu") { + LOG4CPLUS_INFO(logger, "Using ILU preconditioner."); + } else if (precond == "diagonal") { + LOG4CPLUS_INFO(logger, "Using diagonal preconditioner."); + } else if (precond == "ildlt") { + LOG4CPLUS_INFO(logger, "Using ILDLT preconditioner."); + } else if (precond == "none") { + LOG4CPLUS_INFO(logger, "Using no preconditioner."); + } + + if (s->getString("lemethod") == "bicgstab") { + LOG4CPLUS_INFO(logger, "Using BiCGStab method."); + if (precond == "ilu") { + gmm::bicgstab(A, x, b, gmm::ilu_precond>(A), iter); + } else if (precond == "diagonal") { + gmm::bicgstab(A, x, b, gmm::diagonal_precond>(A), iter); + } else if (precond == "ildlt") { + gmm::bicgstab(A, x, b, gmm::ildlt_precond>(A), iter); + } else if (precond == "none") { + gmm::bicgstab(A, x, b, gmm::identity_matrix(), iter); + } + // FIXME: gmres has been disabled, because it triggers gmm++ compilation errors + /* } else if (s->getString("lemethod").compare("gmres") == 0) { + LOG4CPLUS_INFO(logger, "Using GMRES method."); + if (precond.compare("ilu")) { + gmm::gmres(A, x, b, gmm::ilu_precond>(A), s->get("restart"), iter); + } else if (precond == "diagonal") { + gmm::gmres(A, x, b, gmm::diagonal_precond>(A), s->get("restart"), iter); + } else if (precond == "ildlt") { + gmm::gmres(A, x, b, gmm::ildlt_precond>(A), s->get("restart"), iter); + } else if (precond == "none") { + gmm::gmres(A, x, b, gmm::identity_matrix(), s->get("restart"), iter); + } */ + } else if (s->getString("lemethod") == "qmr") { + LOG4CPLUS_INFO(logger, "Using QMR method."); + if (precond == "ilu") { + gmm::qmr(A, x, b, gmm::ilu_precond>(A), iter); + } else if (precond == "diagonal") { + gmm::qmr(A, x, b, gmm::diagonal_precond>(A), iter); + } else if (precond == "ildlt") { + gmm::qmr(A, x, b, gmm::ildlt_precond>(A), iter); + } else if (precond == "none") { + gmm::qmr(A, x, b, gmm::identity_matrix(), iter); + } + } + + // 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 diff --git a/src/models/Dtmc.h b/src/models/Dtmc.h index 69831533b..8b6c8ed43 100644 --- a/src/models/Dtmc.h +++ b/src/models/Dtmc.h @@ -17,6 +17,7 @@ #include "GraphTransitions.h" #include "src/storage/SquareSparseMatrix.h" #include "src/exceptions/InvalidArgumentException.h" +#include "src/utility/CommandLine.h" namespace storm { @@ -148,13 +149,29 @@ public: return *this->backwardTransitions; } + /*! + * Retrieves whether this DTMC has a state reward model. + * @return True if this DTMC has a state reward model. + */ + bool hasStateRewards() { + return this->stateRewards != nullptr; + } + + /*! + * Retrieves whether this DTMC has a transition reward model. + * @return True if this DTMC has a transition reward model. + */ + bool hasTransitionRewards() { + return this->transitionRewardMatrix != nullptr; + } + /*! * Prints information about the model to the specified stream. * @param out The stream the information is to be printed to. */ void printModelInformationToStream(std::ostream& out) const { - out << "-------------------------------------------------------------- " - << std::endl; + storm::utility::printSeparationLine(out); + out << std::endl; out << "Model type: \t\tDTMC" << std::endl; out << "States: \t\t" << this->getNumberOfStates() << std::endl; out << "Transitions: \t\t" << this->getNumberOfTransitions() << std::endl; @@ -163,8 +180,8 @@ public: << (this->probabilityMatrix->getSizeInMemory() + this->stateLabeling->getSizeInMemory() + sizeof(*this))/1024 << " kbytes" << std::endl; - out << "-------------------------------------------------------------- " - << std::endl; + out << std::endl; + storm::utility::printSeparationLine(out); } private: diff --git a/src/solver/GraphAnalyzer.h b/src/solver/GraphAnalyzer.h index 3c99b4426..70b6b10f5 100644 --- a/src/solver/GraphAnalyzer.h +++ b/src/solver/GraphAnalyzer.h @@ -92,6 +92,31 @@ public: alwaysPhiUntilPsiStates->complement(); } + /*! + * Computes the set of states of the given model for which all paths lead to + * the given set of target states and only visit states from the filter set + * before. The results are written to the given bit vector. + * @param model The model whose graph structure to search. + * @param phiStates A bit vector of all states satisfying phi. + * @param psiStates A bit vector of all states satisfying psi. + * @param alwaysPhiUntilPsiStates A pointer to the result of the search for states that only + * have paths satisfying phi until psi. + */ + template + static void getAlwaysPhiUntilPsiStates(storm::models::Dtmc& model, const storm::storage::BitVector& phiStates, const storm::storage::BitVector& psiStates, storm::storage::BitVector* alwaysPhiUntilPsiStates) { + // Check for valid parameter. + if (alwaysPhiUntilPsiStates == nullptr) { + LOG4CPLUS_ERROR(logger, "Parameter 'alwaysPhiUntilPhiStates' must not be null."); + throw storm::exceptions::InvalidArgumentException("Parameter 'alwaysPhiUntilPhiStates' must not be null."); + } + + storm::storage::BitVector existsPhiUntilPsiStates(model.getNumberOfStates()); + GraphAnalyzer::getExistsPhiUntilPsiStates(model, phiStates, psiStates, &existsPhiUntilPsiStates); + GraphAnalyzer::getExistsPhiUntilPsiStates(model, ~psiStates, ~existsPhiUntilPsiStates, alwaysPhiUntilPsiStates); + alwaysPhiUntilPsiStates->complement(); + } + + /*! * Computes the set of states of the given model for which all paths lead to * the given set of target states and only visit states from the filter set diff --git a/src/storage/SquareSparseMatrix.h b/src/storage/SquareSparseMatrix.h index fc657a7df..492b771a2 100644 --- a/src/storage/SquareSparseMatrix.h +++ b/src/storage/SquareSparseMatrix.h @@ -827,6 +827,38 @@ public: return new storm::storage::JacobiDecomposition(resultLU, resultDinv); } + /*! + * Performs a pointwise matrix multiplication of the matrix with the given matrix and returns a + * vector containing the sum of the elements in each row of the resulting matrix. + * @param otherMatrix A reference to the matrix with which to perform the pointwise multiplication. + * This matrix must be a submatrix of the current matrix in the sense that it may not have + * non-zero entries at indices where there is a zero in the current matrix. + * @return A vector containing the sum of the elements in each row of the matrix resulting from + * pointwise multiplication of the current matrix with the given matrix. + */ + std::vector* getPointwiseProductRowSumVector(storm::storage::SquareSparseMatrix const& otherMatrix) { + // Prepare result. + std::vector* result = new std::vector(rowCount); + + // Iterate over all elements of the current matrix and either continue with the next element + // in case the given matrix does not have a non-zero element at this column position, or + // multiply the two entries and add the result to the corresponding position in the vector. + uint_fast64_t otherRow = 0; + for (uint_fast64_t row = 0; row < rowCount; ++row) { + (*result)[row] += diagonalStorage[row] * otherMatrix.diagonalStorage[row]; + for (uint_fast64_t element = otherMatrix.rowIndications[row], nextElement = rowIndications[row]; element < otherMatrix.rowIndications[row + 1]; ++element) { + if (otherMatrix.columnIndications[element] < columnIndications[nextElement]) { + continue; + } else { + (*result)[row] += otherMatrix.valueStorage[element] * valueStorage[nextElement]; + ++nextElement; + } + } + } + + return result; + } + /*! * Returns the size of the matrix in memory measured in bytes. * @return The size of the matrix in memory measured in bytes. diff --git a/src/storm.cpp b/src/storm.cpp index f692179c0..084012bb8 100644 --- a/src/storm.cpp +++ b/src/storm.cpp @@ -66,10 +66,10 @@ void setUpFileLogging() { * Prints the header. */ void printHeader(const int argc, const char* argv[]) { - std::cout << "STORM" << std::endl; + std::cout << "StoRM" << std::endl; std::cout << "====" << std::endl << std::endl; - std::cout << "Version: 1.0" << std::endl; + std::cout << "Version: 1.0 Alpha" << std::endl; // "Compute" the command line argument string with which STORM was invoked. std::stringstream commandStream; for (int i = 0; i < argc; ++i) { @@ -143,20 +143,26 @@ void cleanUp() { */ void testChecking() { storm::settings::Settings* s = storm::settings::instance(); - storm::parser::DtmcParser dtmcParser(s->getString("trafile"), s->getString("labfile")); + storm::parser::DtmcParser dtmcParser(s->getString("trafile"), s->getString("labfile"), s->getString("staterew"), s->getString("transrew")); std::shared_ptr> dtmc = dtmcParser.getDtmc(); dtmc->printModelInformationToStream(std::cout); - storm::formula::Ap* trueFormula = new storm::formula::Ap("true"); - storm::formula::Ap* observe0Greater1Formula = new storm::formula::Ap("observe0Greater1"); - storm::formula::Until* untilFormula = new storm::formula::Until(trueFormula, observe0Greater1Formula); - storm::formula::ProbabilisticNoBoundsOperator* probFormula = new storm::formula::ProbabilisticNoBoundsOperator(untilFormula); + storm::formula::Ap* observe0Greater1Formula = new storm::formula::Ap("one"); + storm::formula::Eventually* eventuallyFormula = new storm::formula::Eventually(observe0Greater1Formula); + storm::formula::ProbabilisticNoBoundsOperator* probFormula = new storm::formula::ProbabilisticNoBoundsOperator(eventuallyFormula); + + storm::formula::Ap* done = new storm::formula::Ap("done"); + storm::formula::ReachabilityReward* reachabilityRewardFormula = new storm::formula::ReachabilityReward(done); + storm::formula::RewardNoBoundsOperator* rewardFormula = new storm::formula::RewardNoBoundsOperator(reachabilityRewardFormula); + storm::modelChecker::GmmxxDtmcPrctlModelChecker* mc = new storm::modelChecker::GmmxxDtmcPrctlModelChecker(*dtmc); mc->check(*probFormula); + mc->check(*rewardFormula); delete mc; delete probFormula; + delete rewardFormula; } /*! @@ -169,7 +175,7 @@ int main(const int argc, const char* argv[]) { } printHeader(argc, argv); - testChecking(); + // testChecking(); cleanUp(); return 0; diff --git a/src/utility/CommandLine.cpp b/src/utility/CommandLine.cpp new file mode 100644 index 000000000..c79948f43 --- /dev/null +++ b/src/utility/CommandLine.cpp @@ -0,0 +1,20 @@ +/* + * CommandLine.cpp + * + * Created on: 26.12.2012 + * Author: Christian Dehnert + */ + +#include + +namespace storm { + +namespace utility { + +void printSeparationLine(std::ostream& out) { + out << "------------------------------------------------------" << std::endl; +} + +} // namespace utility + +} // namespace storm diff --git a/src/utility/CommandLine.h b/src/utility/CommandLine.h new file mode 100644 index 000000000..04599348d --- /dev/null +++ b/src/utility/CommandLine.h @@ -0,0 +1,24 @@ +/* + * CommandLine.h + * + * Created on: 26.12.2012 + * Author: Christian Dehnert + */ + +#ifndef STORM_UTILITY_COMMANDLINE_H_ +#define STORM_UTILITY_COMMANDLINE_H_ + +namespace storm { + +namespace utility { + +/*! + * Prints a separation line on the command line. + */ +void printSeparationLine(std::ostream& out); + +} //namespace utility + +} //namespace storm + +#endif /* STORM_UTILITY_COMMANDLINE_H_ */ diff --git a/src/utility/ConstTemplates.h b/src/utility/ConstTemplates.h index 23fe315ad..1681faf68 100644 --- a/src/utility/ConstTemplates.h +++ b/src/utility/ConstTemplates.h @@ -8,6 +8,8 @@ #ifndef STORM_UTILITY_CONSTTEMPLATES_H_ #define STORM_UTILITY_CONSTTEMPLATES_H_ +#include + namespace storm { namespace utility { @@ -110,6 +112,57 @@ inline double constGetOne() { /*! @endcond */ +/*! + * Returns a constant value of infinity that is fit to the type it is being written to. + * As (at least) gcc has problems to use the correct template by the return value + * only, the function gets a pointer as a parameter to infer the return type. + * + * @note + * The template parameter is just inferred by the return type; GCC is not able to infer this + * automatically, hence the type should always be stated explicitly (e.g. @c constGetOne();) + * + * @return Value Infinity, fit to the return type + */ +template +static inline _Scalar constGetInfinity() { + return std::numeric_limits<_Scalar>::infinity(); +} + +/*! @cond TEMPLATE_SPECIALIZATION + * (By default, the template specifications are not included in the documentation) + */ + +/*! + * Template specification for int_fast32_t + * @return Value Infinity, fit to the type uint_fast32_t + */ +template<> +inline int_fast32_t constGetInfinity() { + throw storm::exceptions::InvalidArgumentException() << "Integer has no infinity."; + return std::numeric_limits::max(); +} + +/*! + * Template specification for uint_fast64_t + * @return Value Infinity, fit to the type uint_fast64_t + */ +template<> +inline uint_fast64_t constGetInfinity() { + throw storm::exceptions::InvalidArgumentException() << "Integer has no infinity."; + return std::numeric_limits::max(); +} + +/*! + * Template specification for double + * @return Value Infinity, fit to the type double + */ +template<> +inline double constGetInfinity() { + return std::numeric_limits::infinity(); +} + +/*! @endcond */ + } //namespace utility } //namespace storm diff --git a/src/utility/Vector.h b/src/utility/Vector.h index c315643b6..187d3c5e2 100644 --- a/src/utility/Vector.h +++ b/src/utility/Vector.h @@ -10,12 +10,14 @@ #include "Eigen/src/Core/Matrix.h" +#include + namespace storm { namespace utility { template -void setVectorValues(std::vector* vector, const storm::storage::BitVector& positions, const std::vector values) { +void setVectorValues(std::vector* vector, const storm::storage::BitVector& positions, std::vector const& values) { uint_fast64_t oldPosition = 0; for (auto position : positions) { (*vector)[position] = values[oldPosition++]; @@ -36,6 +38,21 @@ void setVectorValues(Eigen::Matrix* eigenVector, const storm } } +template +void selectVectorValues(std::vector* vector, const storm::storage::BitVector& positions, std::vector const& values) { + uint_fast64_t oldPosition = 0; + for (auto position : positions) { + (*vector)[oldPosition++] = values[position]; + } +} + +template +void subtractFromConstantOneVector(std::vector* vector) { + for (auto it = vector->begin(); it != vector->end(); ++it) { + *it = storm::utility::constGetOne() - *it; + } +} + } //namespace utility } //namespace storm