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

@@ -56,7 +56,7 @@ public:
 	BoundedUntil(PctlStateFormula<T>* left, PctlStateFormula<T>* right,
 					 uint_fast64_t bound) {
 		this->left = left;
-		this->right = right;;
+		this->right = right;
 		this->bound = bound;
 	}
 

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 <string>
+
+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 T>
+class CumulativeReward : public PctlPathFormula<T> {
+
+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<T>* 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<T> *check(const storm::modelChecker::DtmcPrctlModelChecker<T>& modelChecker) const {
+	  return modelChecker.checkCumulativeReward(*this);
+	}
+
+private:
+	uint_fast64_t bound;
+};
+
+} //namespace formula
+
+} //namespace storm
+
+#endif /* STORM_FORMULA_INSTANTANEOUSREWARD_H_ */

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 T>
+class Eventually : public PctlPathFormula<T> {
+
+public:
+	/*!
+	 * Empty constructor
+	 */
+	Eventually() {
+		this->child = nullptr;
+	}
+
+	/*!
+	 * Constructor
+	 *
+	 * @param child The child node
+	 */
+	Eventually(PctlStateFormula<T>* 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<T>& getChild() const {
+		return *child;
+	}
+
+	/*!
+	 * Sets the subtree
+	 * @param child the new child node
+	 */
+	void setChild(PctlStateFormula<T>* 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<T>* clone() const {
+		Eventually<T>* result = new Eventually<T>();
+		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<T> *check(const storm::modelChecker::DtmcPrctlModelChecker<T>& modelChecker) const {
+	  return modelChecker.checkEventually(*this);
+	}
+
+private:
+	PctlStateFormula<T>* child;
+};
+
+} //namespace formula
+
+} //namespace storm
+
+#endif /* STORM_FORMULA_EVENTUALLY_H_ */

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_ */

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 T>
+class Globally : public PctlPathFormula<T> {
+
+public:
+	/*!
+	 * Empty constructor
+	 */
+	Globally() {
+		this->child = nullptr;
+	}
+
+	/*!
+	 * Constructor
+	 *
+	 * @param child The child node
+	 */
+	Globally(PctlStateFormula<T>* 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<T>& getChild() const {
+		return *child;
+	}
+
+	/*!
+	 * Sets the subtree
+	 * @param child the new child node
+	 */
+	void setChild(PctlStateFormula<T>* 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<T>* clone() const {
+		Next<T>* result = new Next<T>();
+		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<T> *check(const storm::modelChecker::DtmcPrctlModelChecker<T>& modelChecker) const {
+	  return modelChecker.checkGlobally(*this);
+	}
+
+private:
+	PctlStateFormula<T>* child;
+};
+
+} //namespace formula
+
+} //namespace storm
+
+#endif /* STORM_FORMULA_GLOBALLY_H_ */

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 <string>
+
+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 T>
+class InstantaneousReward : public PctlPathFormula<T> {
+
+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<T>* 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<T> *check(const storm::modelChecker::DtmcPrctlModelChecker<T>& modelChecker) const {
+	  return modelChecker.checkInstantaneousReward(*this);
+	}
+
+private:
+	uint_fast64_t bound;
+};
+
+} //namespace formula
+
+} //namespace storm
+
+#endif /* STORM_FORMULA_INSTANTANEOUSREWARD_H_ */

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

src/formula/ProbabilisticOperator.h

@@ -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 T>
-class ProbabilisticOperator : public storm::formula::PctlStateFormula<T> {
-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<T>* 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<T>& 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<T>* 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<T>* clone() const {
-		ProbabilisticOperator<T>* result = new ProbabilisticOperator<T>();
-		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<T>& 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<T>* pathFormula;
-};
-
-} /* namespace formula */
-} /* namespace storm */
-
-#endif /* STORM_FORMULA_PROBABILISTICOPERATOR_H_ */

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 T>
+class ReachabilityReward : public PctlPathFormula<T> {
+
+public:
+	/*!
+	 * Empty constructor
+	 */
+	ReachabilityReward() {
+		this->child = nullptr;
+	}
+
+	/*!
+	 * Constructor
+	 *
+	 * @param child The child node
+	 */
+	ReachabilityReward(PctlStateFormula<T>* 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<T>& getChild() const {
+		return *child;
+	}
+
+	/*!
+	 * Sets the subtree
+	 * @param child the new child node
+	 */
+	void setChild(PctlStateFormula<T>* 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<T>* clone() const {
+		ReachabilityReward<T>* result = new ReachabilityReward<T>();
+		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<T> *check(const storm::modelChecker::DtmcPrctlModelChecker<T>& modelChecker) const {
+	  return modelChecker.checkReachabilityReward(*this);
+	}
+
+private:
+	PctlStateFormula<T>* child;
+};
+
+} //namespace formula
+
+} //namespace storm
+
+#endif /* STORM_FORMULA_REACHABILITYREWARD_H_ */

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 T>
+class RewardIntervalOperator : public PctlStateFormula<T> {
+
+public:
+	/*!
+	 * Empty constructor
+	 */
+	RewardIntervalOperator() {
+		upper = storm::utility::constGetZero<T>();
+		lower = storm::utility::constGetZero<T>();
+		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<T>& 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<T>& 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<T>* 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<T>* clone() const {
+		RewardIntervalOperator<T>* result = new RewardIntervalOperator<T>();
+		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<T>& modelChecker) const {
+	  return modelChecker.checkRewardIntervalOperator(*this);
+	}
+
+private:
+	T lower;
+	T upper;
+	PctlPathFormula<T>* pathFormula;
+};
+
+} //namespace formula
+
+} //namespace storm
+
+#endif /* STORM_FORMULA_REWARDINTERVALOPERATOR_H_ */

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 T>
+class RewardNoBoundsOperator: public storm::formula::PctlFormula<T> {
+public:
+	/*!
+	 * Empty constructor
+	 */
+	RewardNoBoundsOperator() {
+		this->pathFormula = nullptr;
+	}
+
+	/*!
+	 * Constructor
+	 *
+	 * @param pathFormula The child node.
+	 */
+	RewardNoBoundsOperator(PctlPathFormula<T>* pathFormula) {
+		this->pathFormula = pathFormula;
+	}
+
+	/*!
+	 * Destructor
+	 */
+	virtual ~RewardNoBoundsOperator() {
+		if (pathFormula != nullptr) {
+			delete pathFormula;
+		}
+	}
+
+	/*!
+	 * @returns the child node (representation of a PCTL path formula)
+	 */
+	const PctlPathFormula<T>& getPathFormula () const {
+		return *pathFormula;
+	}
+
+	/*!
+	 * Sets the child node
+	 *
+	 * @param pathFormula the path formula that becomes the new child node
+	 */
+	void setPathFormula(PctlPathFormula<T>* 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<T>* pathFormula;
+};
+
+} /* namespace formula */
+
+} /* namespace storm */
+
+#endif /* STORM_FORMULA_REWARDNOBOUNDSOPERATOR_H_ */

src/modelChecker/DtmcPrctlModelChecker.h

@@ -98,13 +98,19 @@ public:
 	 * states.
 	 * @param stateFormula The formula to be checked.
 	 */
-	void check(const storm::formula::PctlStateFormula<Type>& stateFormula) {
-		LOG4CPLUS_INFO(logger, "Model checking formula " << stateFormula.toString());
+	void check(const storm::formula::PctlStateFormula<Type>& 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<Type>& probabilisticNoBoundsFormula) {
-		LOG4CPLUS_INFO(logger, "Model checking formula " << probabilisticNoBoundsFormula.toString());
-		std::cout << "Model checking formula: " << probabilisticNoBoundsFormula.toString() << std::endl;
+	void check(const storm::formula::ProbabilisticNoBoundsOperator<Type>& 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<Type>* 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<Type>& 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<Type>* 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<Type>& formula) const {
+	storm::storage::BitVector* checkProbabilisticIntervalOperator(
+			const storm::formula::ProbabilisticIntervalOperator<Type>& formula) const {
+		// First, we need to compute the probability for satisfying the path formula for each state.
 		std::vector<Type>* 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<Type>& formula) const {
+	storm::storage::BitVector* checkRewardIntervalOperator(
+			const storm::formula::RewardIntervalOperator<Type>& formula) const {
 		// First, we need to compute the probability for satisfying the path formula for each state.
-		std::vector<Type>* probabilisticResult = this->checkPathFormula(formula.getPathFormula());
+		std::vector<Type>* 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<Type>* checkRewardNoBoundsOperator(
+			const storm::formula::RewardNoBoundsOperator<Type>& 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<Type>* checkNext(const storm::formula::Next<Type>& 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<Type>* checkEventually(const storm::formula::Eventually<Type>& formula) const {
+		// Create equivalent temporary until formula and check it.
+		storm::formula::Until<Type> temporaryUntilFormula(new storm::formula::Ap<Type>("true"), formula.getChild().clone());
+		std::vector<Type>* 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<Type>* checkGlobally(const storm::formula::Globally<Type>& formula) const {
+		// Create "equivalent" temporary eventually formula and check it.
+		storm::formula::Eventually<Type> temporaryEventuallyFormula(new storm::formula::Not<Type>(formula.getChild().clone()));
+		std::vector<Type>* 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<Type>* checkUntil(const storm::formula::Until<Type>& 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<Type>* checkInstantaneousReward(const storm::formula::InstantaneousReward<Type>& 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<Type>* checkCumulativeReward(const storm::formula::CumulativeReward<Type>& 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<Type>* checkReachabilityReward(const storm::formula::ReachabilityReward<Type>& formula) const = 0;
+
 private:
 	storm::models::Dtmc<Type>& model;
 };

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<Type>* result = new std::vector<Type>(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<Type> 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<double>("precision"), 0, s->get<unsigned>("lemaxiter"));
+			// Set values of resulting vector according to result.
+			storm::utility::setVectorValues<Type>(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<gmm::csr_matrix<Type>>(*gmmxxMatrix), iter);
-				} else if (precond == "diagonal") {
-					gmm::bicgstab(*gmmxxMatrix, x, b, gmm::diagonal_precond<gmm::csr_matrix<Type>>(*gmmxxMatrix), iter);
-				} else if (precond == "ildlt") {
-					gmm::bicgstab(*gmmxxMatrix, x, b, gmm::ildlt_precond<gmm::csr_matrix<Type>>(*gmmxxMatrix), iter);
-				} else if (precond == "none") {
-					gmm::bicgstab(*gmmxxMatrix, x, b, gmm::identity_matrix(), iter);
+		// Set values of resulting vector that are known exactly.
+		storm::utility::setVectorValues<Type>(result, notExistsPhiUntilPsiStates, storm::utility::constGetZero<Type>());
+		storm::utility::setVectorValues<Type>(result, alwaysPhiUntilPsiStates, storm::utility::constGetOne<Type>());
+
+		return result;
 	}
-			// 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<gmm::csr_matrix<Type>>(*gmmxxMatrix), s->get<unsigned>("restart"), iter);
-				} else if (precond == "diagonal") {
-					gmm::gmres(*gmmxxMatrix, x, b, gmm::diagonal_precond<gmm::csr_matrix<Type>>(*gmmxxMatrix), s->get<unsigned>("restart"), iter);
-				} else if (precond == "ildlt") {
-					gmm::gmres(*gmmxxMatrix, x, b, gmm::ildlt_precond<gmm::csr_matrix<Type>>(*gmmxxMatrix), s->get<unsigned>("restart"), iter);
-				} else if (precond == "none") {
-					gmm::gmres(*gmmxxMatrix, x, b, gmm::identity_matrix(), s->get<unsigned>("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<gmm::csr_matrix<Type>>(*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<gmm::csr_matrix<Type>>(*gmmxxMatrix), iter);
-				} else if (precond == "ildlt") {
-					gmm::qmr(*gmmxxMatrix, x, b, gmm::ildlt_precond<gmm::csr_matrix<Type>>(*gmmxxMatrix), iter);
-				} else if (precond == "none") {
-					gmm::qmr(*gmmxxMatrix, x, b, gmm::identity_matrix(), iter);
+
+	virtual std::vector<Type>* checkInstantaneousReward(const storm::formula::InstantaneousReward<Type>& 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<Type>* gmmxxMatrix = storm::adapters::GmmxxAdapter::toGmmxxSparseMatrix<Type>(*this->getModel().getTransitionProbabilityMatrix());
+
+		// Initialize result to state rewards of the model.
+		std::vector<Type>* result = new std::vector<Type>(*this->getModel().getStateRewards());
+
+		// Now perform matrix-vector multiplication as long as we meet the bound of the formula.
+		std::vector<Type>* swap = nullptr;
+		std::vector<Type>* tmpResult = new std::vector<Type>(this->getModel().getNumberOfStates());
+		for (uint_fast64_t i = 0; i < formula.getBound(); ++i) {
+			gmm::mult(*gmmxxMatrix, *result, *tmpResult);
+			swap = tmpResult;
+			tmpResult = result;
+			result = swap;
 		}
 
-			// Check if the solver converged and issue a warning otherwise.
-			if (iter.converged()) {
-				LOG4CPLUS_INFO(logger, "Iterative solver converged after " << iter.get_iteration() << " iterations.");
+		// Delete temporary variables and return result.
+		delete tmpResult;
+		delete gmmxxMatrix;
+		return result;
+	}
+
+	virtual std::vector<Type>* checkCumulativeReward(const storm::formula::CumulativeReward<Type>& 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<Type>* gmmxxMatrix = storm::adapters::GmmxxAdapter::toGmmxxSparseMatrix<Type>(*this->getModel().getTransitionProbabilityMatrix());
+
+		// Compute the reward vector to add in each step based on the available reward models.
+		std::vector<Type>* 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 {
-				LOG4CPLUS_WARN(logger, "Iterative solver did not converge.");
+			totalRewardVector = new std::vector<Type>(*this->getModel().getStateRewards());
 		}
 
+		std::vector<Type>* result = new std::vector<Type>(this->getModel().getNumberOfStates());
+
+		// Now perform matrix-vector multiplication as long as we meet the bound of the formula.
+		std::vector<Type>* swap = nullptr;
+		std::vector<Type>* tmpResult = new std::vector<Type>(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<Type>* checkReachabilityReward(const storm::formula::ReachabilityReward<Type>& 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<Type>* result = new std::vector<Type>(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<Type>* 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<Type>* gmmxxMatrix = storm::adapters::GmmxxAdapter::toGmmxxSparseMatrix<Type>(*submatrix);
+			delete submatrix;
+
+			// Initialize the x vector with 1 for each element. This is the initial guess for
+			// the iterative solvers.
+			std::vector<Type> x(maybeStates.getNumberOfSetBits(), storm::utility::constGetOne<Type>());
+
+			// Prepare the right-hand side of the equation system.
+			std::vector<Type>* b = new std::vector<Type>(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<Type>* 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<Type>* subStateRewards = new std::vector<Type>(maybeStates.getNumberOfSetBits());
+					storm::utility::setVectorValues(subStateRewards, maybeStates, *this->getModel().getStateRewards());
+					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().getStateRewards());
+			}
+
+			// Solve the corresponding system of linear equations.
+			this->solveLinearEquationSystem(*gmmxxMatrix, x, *b);
+
 			// Set values of resulting vector according to result.
 			storm::utility::setVectorValues<Type>(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<Type>(result, notExistsPhiUntilPsiStates, storm::utility::constGetZero<Type>());
-		storm::utility::setVectorValues<Type>(result, alwaysPhiUntilPsiStates, storm::utility::constGetOne<Type>());
+		storm::utility::setVectorValues(result, *targetStates, storm::utility::constGetZero<Type>());
+		storm::utility::setVectorValues(result, infinityStates, storm::utility::constGetInfinity<Type>());
 
+		// 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<Type> const& A, std::vector<Type>& x, std::vector<Type> 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<double>("precision"), 0, s->get<unsigned>("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<gmm::csr_matrix<Type>>(A), iter);
+			} else if (precond == "diagonal") {
+				gmm::bicgstab(A, x, b, gmm::diagonal_precond<gmm::csr_matrix<Type>>(A), iter);
+			} else if (precond == "ildlt") {
+				gmm::bicgstab(A, x, b, gmm::ildlt_precond<gmm::csr_matrix<Type>>(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<gmm::csr_matrix<Type>>(A), s->get<unsigned>("restart"), iter);
+			} else if (precond == "diagonal") {
+				gmm::gmres(A, x, b, gmm::diagonal_precond<gmm::csr_matrix<Type>>(A), s->get<unsigned>("restart"), iter);
+			} else if (precond == "ildlt") {
+				gmm::gmres(A, x, b, gmm::ildlt_precond<gmm::csr_matrix<Type>>(A), s->get<unsigned>("restart"), iter);
+			} else if (precond == "none") {
+				gmm::gmres(A, x, b, gmm::identity_matrix(), s->get<unsigned>("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<gmm::csr_matrix<Type>>(A), iter);
+			} else if (precond == "diagonal") {
+				gmm::qmr(A, x, b, gmm::diagonal_precond<gmm::csr_matrix<Type>>(A), iter);
+			} else if (precond == "ildlt") {
+				gmm::qmr(A, x, b, gmm::ildlt_precond<gmm::csr_matrix<Type>>(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

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:

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 <class T>
+	static void getAlwaysPhiUntilPsiStates(storm::models::Dtmc<T>& 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

src/storage/SquareSparseMatrix.h

@@ -827,6 +827,38 @@ public:
 		return new storm::storage::JacobiDecomposition<T>(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<T>* getPointwiseProductRowSumVector(storm::storage::SquareSparseMatrix<T> const& otherMatrix) {
+		// Prepare result.
+		std::vector<T>* result = new std::vector<T>(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.

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<storm::models::Dtmc<double>> dtmc = dtmcParser.getDtmc();
 
 	dtmc->printModelInformationToStream(std::cout);
 
-	storm::formula::Ap<double>* trueFormula = new storm::formula::Ap<double>("true");
-	storm::formula::Ap<double>* observe0Greater1Formula = new storm::formula::Ap<double>("observe0Greater1");
-	storm::formula::Until<double>* untilFormula = new storm::formula::Until<double>(trueFormula, observe0Greater1Formula);
-	storm::formula::ProbabilisticNoBoundsOperator<double>* probFormula = new storm::formula::ProbabilisticNoBoundsOperator<double>(untilFormula);
+	storm::formula::Ap<double>* observe0Greater1Formula = new storm::formula::Ap<double>("one");
+	storm::formula::Eventually<double>* eventuallyFormula = new storm::formula::Eventually<double>(observe0Greater1Formula);
+	storm::formula::ProbabilisticNoBoundsOperator<double>* probFormula = new storm::formula::ProbabilisticNoBoundsOperator<double>(eventuallyFormula);
+
+	storm::formula::Ap<double>* done = new storm::formula::Ap<double>("done");
+	storm::formula::ReachabilityReward<double>* reachabilityRewardFormula = new storm::formula::ReachabilityReward<double>(done);
+	storm::formula::RewardNoBoundsOperator<double>* rewardFormula = new storm::formula::RewardNoBoundsOperator<double>(reachabilityRewardFormula);
+
 	storm::modelChecker::GmmxxDtmcPrctlModelChecker<double>* mc = new storm::modelChecker::GmmxxDtmcPrctlModelChecker<double>(*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;

src/utility/CommandLine.cpp

@@ -0,0 +1,20 @@
+/*
+ * CommandLine.cpp
+ *
+ *  Created on: 26.12.2012
+ *      Author: Christian Dehnert
+ */
+
+#include <ostream>
+
+namespace storm {
+
+namespace utility {
+
+void printSeparationLine(std::ostream& out) {
+	out << "------------------------------------------------------" << std::endl;
+}
+
+} // namespace utility
+
+} // namespace storm

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_ */

src/utility/ConstTemplates.h

@@ -8,6 +8,8 @@
 #ifndef STORM_UTILITY_CONSTTEMPLATES_H_
 #define STORM_UTILITY_CONSTTEMPLATES_H_
 
+#include <limits>
+
 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<int>();)
+ *
+ * @return Value Infinity, fit to the return type
+ */
+template<typename _Scalar>
+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<int_fast32_t>::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<int_fast64_t>::max();
+}
+
+/*!
+ * Template specification for double
+ * @return Value Infinity, fit to the type double
+ */
+template<>
+inline double constGetInfinity() {
+   return std::numeric_limits<double>::infinity();
+}
+
+/*! @endcond */
+
 } //namespace utility
 
 } //namespace storm

src/utility/Vector.h

@@ -10,12 +10,14 @@
 
 #include "Eigen/src/Core/Matrix.h"
 
+#include <iostream>
+
 namespace storm {
 
 namespace utility {
 
 template<class T>
-void setVectorValues(std::vector<T>* vector, const storm::storage::BitVector& positions, const std::vector<T> values) {
+void setVectorValues(std::vector<T>* vector, const storm::storage::BitVector& positions, std::vector<T> const& values) {
 	uint_fast64_t oldPosition = 0;
 	for (auto position : positions) {
 		(*vector)[position] = values[oldPosition++];
@@ -36,6 +38,21 @@ void setVectorValues(Eigen::Matrix<T, -1, 1, 0, -1, 1>* eigenVector, const storm
 	}
 }
 
+template<class T>
+void selectVectorValues(std::vector<T>* vector, const storm::storage::BitVector& positions, std::vector<T> const& values) {
+	uint_fast64_t oldPosition = 0;
+	for (auto position : positions) {
+		(*vector)[oldPosition++] = values[position];
+	}
+}
+
+template<class T>
+void subtractFromConstantOneVector(std::vector<T>* vector) {
+	for (auto it = vector->begin(); it != vector->end(); ++it) {
+		*it = storm::utility::constGetOne<T>() - *it;
+	}
+}
+
 } //namespace utility
 
 } //namespace storm