Merge remote-tracking branch 'upstream/master'

CHANGELOG.md

@@ -7,6 +7,9 @@ The releases of major and minor versions contain an overview of changes since th
 Version 1.2.x
 -------------
 
+### Version 1.2.1 (to be released)
+- Multi-dimensional reward bounded reachability properties for DTMCs.
+
 ### Version 1.2.0 (2017/12)
 - C++ api changes: Building model takes `BuilderOptions` instead of extended list of Booleans, does not depend on settings anymore.
 - `storm-cli-utilities` now contains cli related stuff, instead of `storm-lib`

doc/checklist_new_release.md

@@ -13,6 +13,7 @@ Note that in most case a simultaneous release of [carl](https://github.com/smtra
 
 3. Check that storm builds without errors and all tests are successful
    * [Travis](https://travis-ci.org/moves-rwth/storm) should run successfully
+   * Invoke the script test_build_configurations.py to build and check different cmake configurations
 
 4. Set new storm version:
    * Set new storm version in `version.cmake`

doc/scripts/test_build_configurations.py

@@ -0,0 +1,118 @@
+import subprocess
+import os
+import shutil
+import time
+import math
+import sys
+
+logfileDir = "test_build_configurations_logs"
+pathToConfigFile = ""
+globalCmakeArguments = ""
+globalMakeArguments = ""
+
+if len(sys.argv) == 1 or len(sys.argv) > 5:
+	print "Usage: " + sys.argv[0] + "/path/to/storm /path/to/configurations.txt \"optional make arguments\" \"optional cmake arguments\""
+	print "Example: " + sys.argv[0] + " ~/storm test_build_configurations.txt \"-j 48 -k\" \"-DZ3_ROOT=/path/to/z3/\""
+	sys.exit(0)
+
+pathToStorm = sys.argv[1]
+pathToConfigFile = sys.argv[2]
+
+if len(sys.argv) > 3:
+	globalMakeArguments = sys.argv[3]
+	print globalMakeArguments
+if len(sys.argv) > 4:
+	globalCmakeArguments = sys.argv[4]
+	print globalCmakeArguments
+
+# create directory for log files
+if not os.path.exists(logfileDir):
+    os.makedirs(logfileDir)
+
+unsuccessfulConfigs = ""
+globalStartTime = time.time()
+
+with open(pathToConfigFile) as configfile:
+    localStartTime = time.time()
+    configId=0
+    for localCmakeArguments in configfile:
+        
+        cmakeArguments = globalCmakeArguments + " " + localCmakeArguments.strip('\n')
+        buildDir = os.path.join(pathToStorm, "build{}".format(configId))
+        logfilename = os.path.join(logfileDir, "build{}".format(configId) + "_" + time.strftime("%Y-%m-%d-%H-%M-%S") + ".log")
+        
+        
+        print "Building configuration {} with cmake options ".format(configId) + cmakeArguments
+        print "\tCreating log file " + logfilename + " ..."
+        with open(logfilename, "w") as logfile:
+            success=True
+            logfile.write("Log for test configuration " + cmakeArguments + "\n")
+                          
+            print "\tCreating build directory" + buildDir + " ..."
+            if os.path.exists(buildDir):
+                print "\t\tRemoving existing directory " + buildDir
+                shutil.rmtree(buildDir, ignore_errors=True)
+            os.makedirs(buildDir)
+            logfile.write("Build directory is " + buildDir + "\n")
+            print "\t\tdone"
+            
+            if success:
+                print "\tInvoking cmake ..."
+                cmakeCommand = "cmake .. {}".format(cmakeArguments)
+                logfile.write ("\n\n CALLING CMAKE \n\n" + cmakeCommand + "\n")
+                try:
+                    cmakeOutput = subprocess.check_output("cd " + buildDir + "; " + cmakeCommand , shell=True,stderr=subprocess.STDOUT)
+                    print "\t\tdone"
+                    logfile.write(cmakeOutput)
+                except subprocess.CalledProcessError as e:
+                    success=False
+                    print "\t\tfail"
+                    print e.output
+                    logfile.write(e.output)
+                                
+            if success:
+                print "\tInvoking make ..."
+                makeCommand = "make {}".format(globalMakeArguments)
+                logfile.write ("\n\n CALLING MAKE \n\n" + makeCommand + "\n")
+                try:
+                    makeOutput = subprocess.check_output("cd " + buildDir + "; " + makeCommand , shell=True,stderr=subprocess.STDOUT)
+                    print "\t\tdone"
+                    logfile.write(makeOutput)
+                except subprocess.CalledProcessError as e:
+                    success=False
+                    print "\t\tfail"
+                    print e.output
+                    logfile.write(e.output)
+            
+            if success:
+                print "\tInvoking make check..."
+                makeCheckCommand = "make check"
+                logfile.write ("\n\n CALLING MAKE CHECK \n\n" + makeCheckCommand + "\n")
+                try:
+                    makeCheckOutput = subprocess.check_output("cd " + buildDir + "; " + makeCheckCommand , shell=True,stderr=subprocess.STDOUT)
+                    print "\t\tdone"
+                    logfile.write(makeCheckOutput)
+                except subprocess.CalledProcessError as e:
+                    success=False
+                    print "\t\tfail"
+                    print e.output
+                    logfile.write(e.output)
+                    
+            localEndTime = time.time()
+            if success:
+                print "\tConfiguration build and tested successfully within {} minutes.".format(int((localEndTime - localStartTime)/60))
+            else:
+                print "\tAn error occurred for this configuration."
+                unsuccessfulConfigs += buildDir + " with arguments " + cmakeArguments + "\n"
+
+            configId += 1
+globalEndTime = time.time()
+print "All tests completed after {} minutes.".format(int((globalEndTime - globalStartTime)/60))
+if unsuccessfulConfigs == "":
+    print "All configurations were build and tested successfully."
+else:
+    print "The following configurations failed: \n" + unsuccessfulConfigs
+    
+                
+
+

doc/scripts/test_build_configurations.txt

@@ -0,0 +1,5 @@
+-DSTORM_USE_CLN_EA=OFF -DSTORM_USE_CLN_RF=OFF
+-DSTORM_USE_CLN_EA=OFF -DSTORM_USE_CLN_RF=ON
+-DSTORM_USE_CLN_EA=ON -DSTORM_USE_CLN_RF=OFF
+-DSTORM_USE_CLN_EA=ON -DSTORM_USE_CLN_RF=ON
+-DCMAKE_BUILD_TYPE=DEBUG -DSTORM_DEVELOPER=ON

resources/examples/testfiles/dtmc/crowds_cost_bounded.pm

@@ -0,0 +1,260 @@
+// CROWDS [Reiter,Rubin]
+// Vitaly Shmatikov, 2002
+
+// Note:
+// Change everything marked CWDSIZ when changing the size of the crowd
+// Change everything marked CWDMAX when increasing max size of the crowd
+
+dtmc
+
+// Probability of forwarding
+const double PF = 0.8;
+
+// Probability that a crowd member is bad
+const double  badC = 0.091;
+// const double  badC = 0.167;
+
+const int CrowdSize; // CWDSIZ: actual number of good crowd members
+const int MaxGood=20; // CWDMAX: maximum number of good crowd members
+
+// Process definitions
+module crowds
+
+	// Auxiliary variables
+	launch:   bool init true;       // Start modeling?
+	new:      bool init false;      // Initialize a new protocol instance?
+	start:    bool init false;      // Start the protocol?
+	run:      bool init false;      // Run the protocol?
+	lastSeen: [0..MaxGood] init MaxGood;   // Last crowd member to touch msg
+	good:     bool init false;      // Crowd member is good?
+	bad:      bool init false;      //              ... bad?
+	recordLast: bool init false;    // Record last seen crowd member?
+	badObserve: bool init false;    // Bad members observes who sent msg?
+	deliver:  bool init false;      // Deliver message to destination?
+	done:     bool init false;      // Protocol instance finished?
+	
+	[] launch -> (new'=true) & (launch'=false);
+	// Set up a new protocol instance
+	[newrun] new -> (new'=false) & (start'=true);
+	
+	// SENDER
+	// Start the protocol
+	[] start -> (lastSeen'=0) & (run'=true) & (deliver'=false) & (start'=false);
+	
+	// CROWD MEMBERS
+	// Good or bad crowd member?
+	[] !good & !bad & !deliver & run ->
+	             1-badC : (good'=true) & (recordLast'=true) & (run'=false) +
+	               badC : (bad'=true)  & (badObserve'=true) & (run'=false);
+
+	// GOOD MEMBERS
+	// Forward with probability PF, else deliver
+	[] good & !deliver & run -> PF : (good'=false) + 1-PF : (deliver'=true);
+	// Record the last crowd member who touched the msg;
+	// all good members may appear with equal probability
+	//    Note: This is backward.  In the real protocol, each honest
+	//          forwarder randomly chooses the next forwarder.
+	//          Here, the identity of an honest forwarder is randomly
+	//          chosen *after* it has forwarded the message.
+	[] recordLast & CrowdSize=2 ->
+	        1/2 : (lastSeen'=0) & (recordLast'=false) & (run'=true) +
+	        1/2 : (lastSeen'=1) & (recordLast'=false) & (run'=true);
+	[] recordLast & CrowdSize=4 ->
+	        1/4 : (lastSeen'=0) & (recordLast'=false) & (run'=true) +
+	        1/4 : (lastSeen'=1) & (recordLast'=false) & (run'=true) +
+	        1/4 : (lastSeen'=2) & (recordLast'=false) & (run'=true) +
+	        1/4 : (lastSeen'=3) & (recordLast'=false) & (run'=true);
+	[] recordLast & CrowdSize=5 ->
+	        1/5 : (lastSeen'=0) & (recordLast'=false) & (run'=true) +
+	        1/5 : (lastSeen'=1) & (recordLast'=false) & (run'=true) +
+	        1/5 : (lastSeen'=2) & (recordLast'=false) & (run'=true) +
+	        1/5 : (lastSeen'=3) & (recordLast'=false) & (run'=true) +
+	        1/5 : (lastSeen'=4) & (recordLast'=false) & (run'=true);
+	[] recordLast & CrowdSize=10 ->
+	        1/10 : (lastSeen'=0) & (recordLast'=false) & (run'=true) +
+	        1/10 : (lastSeen'=1) & (recordLast'=false) & (run'=true) +
+	        1/10 : (lastSeen'=2) & (recordLast'=false) & (run'=true) +
+	        1/10 : (lastSeen'=3) & (recordLast'=false) & (run'=true) +
+	        1/10 : (lastSeen'=4) & (recordLast'=false) & (run'=true) +
+	        1/10 : (lastSeen'=5) & (recordLast'=false) & (run'=true) +
+	        1/10 : (lastSeen'=6) & (recordLast'=false) & (run'=true) +
+	        1/10 : (lastSeen'=7) & (recordLast'=false) & (run'=true) +
+	        1/10 : (lastSeen'=8) & (recordLast'=false) & (run'=true) +
+	        1/10 : (lastSeen'=9) & (recordLast'=false) & (run'=true);
+	[] recordLast & CrowdSize=15 ->
+	        1/15 : (lastSeen'=0)  & (recordLast'=false) & (run'=true) +
+	        1/15 : (lastSeen'=1)  & (recordLast'=false) & (run'=true) +
+	        1/15 : (lastSeen'=2)  & (recordLast'=false) & (run'=true) +
+	        1/15 : (lastSeen'=3)  & (recordLast'=false) & (run'=true) +
+	        1/15 : (lastSeen'=4)  & (recordLast'=false) & (run'=true) +
+	        1/15 : (lastSeen'=5)  & (recordLast'=false) & (run'=true) +
+	        1/15 : (lastSeen'=6)  & (recordLast'=false) & (run'=true) +
+	        1/15 : (lastSeen'=7)  & (recordLast'=false) & (run'=true) +
+	        1/15 : (lastSeen'=8)  & (recordLast'=false) & (run'=true) +
+	        1/15 : (lastSeen'=9)  & (recordLast'=false) & (run'=true) +
+	        1/15 : (lastSeen'=10) & (recordLast'=false) & (run'=true) +
+	        1/15 : (lastSeen'=11) & (recordLast'=false) & (run'=true) +
+	        1/15 : (lastSeen'=12) & (recordLast'=false) & (run'=true) +
+	        1/15 : (lastSeen'=13) & (recordLast'=false) & (run'=true) +
+	        1/15 : (lastSeen'=14) & (recordLast'=false) & (run'=true);
+	[] recordLast & CrowdSize=20 ->
+	        1/20 : (lastSeen'=0)  & (recordLast'=false) & (run'=true) +
+	        1/20 : (lastSeen'=1)  & (recordLast'=false) & (run'=true) +
+	        1/20 : (lastSeen'=2)  & (recordLast'=false) & (run'=true) +
+	        1/20 : (lastSeen'=3)  & (recordLast'=false) & (run'=true) +
+	        1/20 : (lastSeen'=4)  & (recordLast'=false) & (run'=true) +
+	        1/20 : (lastSeen'=5)  & (recordLast'=false) & (run'=true) +
+	        1/20 : (lastSeen'=6)  & (recordLast'=false) & (run'=true) +
+	        1/20 : (lastSeen'=7)  & (recordLast'=false) & (run'=true) +
+	        1/20 : (lastSeen'=8)  & (recordLast'=false) & (run'=true) +
+	        1/20 : (lastSeen'=9)  & (recordLast'=false) & (run'=true) +
+	        1/20 : (lastSeen'=10) & (recordLast'=false) & (run'=true) +
+	        1/20 : (lastSeen'=11) & (recordLast'=false) & (run'=true) +
+	        1/20 : (lastSeen'=12) & (recordLast'=false) & (run'=true) +
+	        1/20 : (lastSeen'=13) & (recordLast'=false) & (run'=true) +
+	        1/20 : (lastSeen'=14) & (recordLast'=false) & (run'=true) +
+	        1/20 : (lastSeen'=15) & (recordLast'=false) & (run'=true) +
+	        1/20 : (lastSeen'=16) & (recordLast'=false) & (run'=true) +
+	        1/20 : (lastSeen'=17) & (recordLast'=false) & (run'=true) +
+	        1/20 : (lastSeen'=18) & (recordLast'=false) & (run'=true) +
+	        1/20 : (lastSeen'=19) & (recordLast'=false) & (run'=true);
+	
+	// BAD MEMBERS
+	// Remember from whom the message was received and deliver
+	// CWDMAX: 1 rule per each good crowd member
+	[obs0]  lastSeen=0  & badObserve -> (deliver'=true) & (run'=true) & (badObserve'=false);
+	[obs1]  lastSeen=1  & badObserve -> (deliver'=true) & (run'=true) & (badObserve'=false);
+	[obs2]  lastSeen=2  & badObserve -> (deliver'=true) & (run'=true) & (badObserve'=false);
+	[obs3]  lastSeen=3  & badObserve -> (deliver'=true) & (run'=true) & (badObserve'=false);
+	[obs4]  lastSeen=4  & badObserve -> (deliver'=true) & (run'=true) & (badObserve'=false);
+	[obs5]  lastSeen=5  & badObserve -> (deliver'=true) & (run'=true) & (badObserve'=false);
+	[obs6]  lastSeen=6  & badObserve -> (deliver'=true) & (run'=true) & (badObserve'=false);
+	[obs7]  lastSeen=7  & badObserve -> (deliver'=true) & (run'=true) & (badObserve'=false);
+	[obs8]  lastSeen=8  & badObserve -> (deliver'=true) & (run'=true) & (badObserve'=false);
+	[obs9]  lastSeen=9  & badObserve -> (deliver'=true) & (run'=true) & (badObserve'=false);
+	[obs10] lastSeen=10 & badObserve -> (deliver'=true) & (run'=true) & (badObserve'=false);
+	[obs11] lastSeen=11 & badObserve -> (deliver'=true) & (run'=true) & (badObserve'=false);
+	[obs12] lastSeen=12 & badObserve -> (deliver'=true) & (run'=true) & (badObserve'=false);
+	[obs13] lastSeen=13 & badObserve -> (deliver'=true) & (run'=true) & (badObserve'=false);
+	[obs14] lastSeen=14 & badObserve -> (deliver'=true) & (run'=true) & (badObserve'=false);
+	[obs15] lastSeen=15 & badObserve -> (deliver'=true) & (run'=true) & (badObserve'=false);
+	[obs16] lastSeen=16 & badObserve -> (deliver'=true) & (run'=true) & (badObserve'=false);
+	[obs17] lastSeen=17 & badObserve -> (deliver'=true) & (run'=true) & (badObserve'=false);
+	[obs18] lastSeen=18 & badObserve -> (deliver'=true) & (run'=true) & (badObserve'=false);
+	[obs19] lastSeen=19 & badObserve -> (deliver'=true) & (run'=true) & (badObserve'=false);
+
+	// RECIPIENT
+	// Delivery to destination
+	[] deliver & run -> (done'=true) & (deliver'=false) & (run'=false) & (good'=false) & (bad'=false);
+	// Start a new instance
+	[] done -> (new'=true) & (done'=false) & (run'=false) & (lastSeen'=MaxGood);
+	
+endmodule
+
+rewards "num_runs"
+	[newrun] true : 1;
+endrewards
+
+rewards "observe0"
+	[obs0]  true : 1;
+endrewards
+
+rewards "observe1"
+	[obs1]  true : 1;
+endrewards
+
+rewards "observe2"
+	[obs2]  true : 1;
+endrewards
+
+rewards "observe3"
+	[obs3]  true : 1;
+endrewards
+
+rewards "observe4"
+	[obs4]  true : 1;
+endrewards
+
+rewards "observe5"
+	[obs5]  true : 1;
+endrewards
+
+rewards "observe6"
+	[obs6]  true : 1;
+endrewards
+
+rewards "observe7"
+	[obs7]  true : 1;
+endrewards
+
+rewards "observe8"
+	[obs8]  true : 1;
+endrewards
+
+rewards "observe9"
+	[obs9]  true : 1;
+endrewards
+
+rewards "observe10"
+	[obs10] true : 1;
+endrewards
+
+rewards "observe11"
+	[obs11] true : 1;
+endrewards
+
+rewards "observe12"
+	[obs12] true : 1;
+endrewards
+
+rewards "observe13"
+	[obs13] true : 1;
+endrewards
+
+rewards "observe14"
+	[obs14] true : 1;
+endrewards
+
+rewards "observe15"
+	[obs15] true : 1;
+endrewards
+
+rewards "observe16"
+	[obs16] true : 1;
+endrewards
+
+rewards "observe17"
+	[obs17] true : 1;
+endrewards
+
+rewards "observe18"
+	[obs18] true : 1;
+endrewards
+
+rewards "observe19"
+	[obs19] true : 1;
+endrewards
+
+rewards "observeI"
+	[obs1] true : 1;
+	[obs2] true : 1;
+	[obs3] true : 1;
+	[obs4] true : 1;
+	[obs5] true : 1;
+	[obs6] true : 1;
+	[obs7] true : 1;
+	[obs8] true : 1;
+	[obs9] true : 1;
+	[obs10] true : 1;
+	[obs11] true : 1;
+	[obs12] true : 1;
+	[obs13] true : 1;
+	[obs14] true : 1;
+	[obs15] true : 1;
+	[obs16] true : 1;
+	[obs17] true : 1;
+	[obs18] true : 1;
+	[obs19] true : 1;
+endrewards
+
+

src/storm/modelchecker/prctl/SparseDtmcPrctlModelChecker.cpp

@@ -36,21 +36,32 @@ namespace storm {
         template<typename SparseDtmcModelType>
         bool SparseDtmcPrctlModelChecker<SparseDtmcModelType>::canHandle(CheckTask<storm::logic::Formula, ValueType> const& checkTask) const {
             storm::logic::Formula const& formula = checkTask.getFormula();
-            return formula.isInFragment(storm::logic::prctl().setLongRunAverageRewardFormulasAllowed(true).setLongRunAverageProbabilitiesAllowed(true).setConditionalProbabilityFormulasAllowed(true).setConditionalRewardFormulasAllowed(true).setOnlyEventuallyFormuluasInConditionalFormulasAllowed(true));
+            return formula.isInFragment(storm::logic::prctl().setLongRunAverageRewardFormulasAllowed(true).setLongRunAverageProbabilitiesAllowed(true).setConditionalProbabilityFormulasAllowed(true).setConditionalRewardFormulasAllowed(true).setOnlyEventuallyFormuluasInConditionalFormulasAllowed(true).setRewardBoundedUntilFormulasAllowed(true).setRewardBoundedCumulativeRewardFormulasAllowed(true).setMultiDimensionalBoundedUntilFormulasAllowed(true).setMultiDimensionalCumulativeRewardFormulasAllowed(true));
         }
         
         template<typename SparseDtmcModelType>
         std::unique_ptr<CheckResult> SparseDtmcPrctlModelChecker<SparseDtmcModelType>::computeBoundedUntilProbabilities(Environment const& env, CheckTask<storm::logic::BoundedUntilFormula, ValueType> const& checkTask) {
             storm::logic::BoundedUntilFormula const& pathFormula = checkTask.getFormula();
-            STORM_LOG_THROW(!pathFormula.hasLowerBound() && pathFormula.hasUpperBound(), storm::exceptions::InvalidPropertyException, "Formula needs to have single upper time bound.");
-            STORM_LOG_THROW(pathFormula.hasIntegerUpperBound(), storm::exceptions::InvalidPropertyException, "Formula needs to have discrete upper time bound.");
-            std::unique_ptr<CheckResult> leftResultPointer = this->check(env, pathFormula.getLeftSubformula());
-            std::unique_ptr<CheckResult> rightResultPointer = this->check(env, pathFormula.getRightSubformula());
-            ExplicitQualitativeCheckResult const& leftResult = leftResultPointer->asExplicitQualitativeCheckResult();
-            ExplicitQualitativeCheckResult const& rightResult = rightResultPointer->asExplicitQualitativeCheckResult();
-            std::vector<ValueType> numericResult = storm::modelchecker::helper::SparseDtmcPrctlHelper<ValueType>::computeBoundedUntilProbabilities(env, storm::solver::SolveGoal<ValueType>(this->getModel(), checkTask), this->getModel().getTransitionMatrix(), this->getModel().getBackwardTransitions(), leftResult.getTruthValuesVector(), rightResult.getTruthValuesVector(), pathFormula.getNonStrictUpperBound<uint64_t>(), *linearEquationSolverFactory, checkTask.getHint());
-            std::unique_ptr<CheckResult> result = std::unique_ptr<CheckResult>(new ExplicitQuantitativeCheckResult<ValueType>(std::move(numericResult)));
-            return result;
+            if (pathFormula.isMultiDimensional() || pathFormula.getTimeBoundReference().isRewardBound()) {
+                STORM_LOG_THROW(checkTask.isOnlyInitialStatesRelevantSet(), storm::exceptions::InvalidOperationException, "Checking non-trivial bounded until probabilities can only be computed for the initial states of the model.");
+                storm::logic::OperatorInformation opInfo;
+                if (checkTask.isBoundSet()) {
+                    opInfo.bound = checkTask.getBound();
+                }
+                auto formula = std::make_shared<storm::logic::ProbabilityOperatorFormula>(checkTask.getFormula().asSharedPointer(), opInfo);
+                auto numericResult = storm::modelchecker::helper::SparseDtmcPrctlHelper<ValueType>::computeRewardBoundedValues(env, this->getModel(), formula, *linearEquationSolverFactory);
+                return std::unique_ptr<CheckResult>(new ExplicitQuantitativeCheckResult<ValueType>(std::move(numericResult)));
+            } else {
+                STORM_LOG_THROW(!pathFormula.hasLowerBound() && pathFormula.hasUpperBound(), storm::exceptions::InvalidPropertyException, "Formula needs to have single upper time bound.");
+                STORM_LOG_THROW(pathFormula.hasIntegerUpperBound(), storm::exceptions::InvalidPropertyException, "Formula needs to have discrete upper time bound.");
+                std::unique_ptr<CheckResult> leftResultPointer = this->check(env, pathFormula.getLeftSubformula());
+                std::unique_ptr<CheckResult> rightResultPointer = this->check(env, pathFormula.getRightSubformula());
+                ExplicitQualitativeCheckResult const& leftResult = leftResultPointer->asExplicitQualitativeCheckResult();
+                ExplicitQualitativeCheckResult const& rightResult = rightResultPointer->asExplicitQualitativeCheckResult();
+                std::vector<ValueType> numericResult = storm::modelchecker::helper::SparseDtmcPrctlHelper<ValueType>::computeStepBoundedUntilProbabilities(env, storm::solver::SolveGoal<ValueType>(this->getModel(), checkTask), this->getModel().getTransitionMatrix(), this->getModel().getBackwardTransitions(), leftResult.getTruthValuesVector(), rightResult.getTruthValuesVector(), pathFormula.getNonStrictUpperBound<uint64_t>(), *linearEquationSolverFactory, checkTask.getHint());
+                std::unique_ptr<CheckResult> result = std::unique_ptr<CheckResult>(new ExplicitQuantitativeCheckResult<ValueType>(std::move(numericResult)));
+               return result;
+            }
         }
         
         template<typename SparseDtmcModelType>
@@ -85,9 +96,20 @@ namespace storm {
         template<typename SparseDtmcModelType>
         std::unique_ptr<CheckResult> SparseDtmcPrctlModelChecker<SparseDtmcModelType>::computeCumulativeRewards(Environment const& env, storm::logic::RewardMeasureType, CheckTask<storm::logic::CumulativeRewardFormula, ValueType> const& checkTask) {
             storm::logic::CumulativeRewardFormula const& rewardPathFormula = checkTask.getFormula();
-            STORM_LOG_THROW(rewardPathFormula.hasIntegerBound(), storm::exceptions::InvalidPropertyException, "Formula needs to have a discrete time bound.");
-            std::vector<ValueType> numericResult = storm::modelchecker::helper::SparseDtmcPrctlHelper<ValueType>::computeCumulativeRewards(env, storm::solver::SolveGoal<ValueType>(this->getModel(), checkTask), this->getModel().getTransitionMatrix(), checkTask.isRewardModelSet() ? this->getModel().getRewardModel(checkTask.getRewardModel()) : this->getModel().getRewardModel(""), rewardPathFormula.getNonStrictBound<uint64_t>(), *linearEquationSolverFactory);
-            return std::unique_ptr<CheckResult>(new ExplicitQuantitativeCheckResult<ValueType>(std::move(numericResult)));
+            if (rewardPathFormula.isMultiDimensional() || rewardPathFormula.getTimeBoundReference().isRewardBound()) {
+                STORM_LOG_THROW(checkTask.isOnlyInitialStatesRelevantSet(), storm::exceptions::InvalidOperationException, "Checking non-trivial bounded until probabilities can only be computed for the initial states of the model.");
+                storm::logic::OperatorInformation opInfo;
+                if (checkTask.isBoundSet()) {
+                    opInfo.bound = checkTask.getBound();
+                }
+                auto formula = std::make_shared<storm::logic::RewardOperatorFormula>(checkTask.getFormula().asSharedPointer(), checkTask.getRewardModel(), opInfo);
+                auto numericResult = storm::modelchecker::helper::SparseDtmcPrctlHelper<ValueType>::computeRewardBoundedValues(env, this->getModel(), formula, *linearEquationSolverFactory);
+                return std::unique_ptr<CheckResult>(new ExplicitQuantitativeCheckResult<ValueType>(std::move(numericResult)));
+            } else {
+                STORM_LOG_THROW(rewardPathFormula.hasIntegerBound(), storm::exceptions::InvalidPropertyException, "Formula needs to have a discrete time bound.");
+                std::vector<ValueType> numericResult = storm::modelchecker::helper::SparseDtmcPrctlHelper<ValueType>::computeCumulativeRewards(env, storm::solver::SolveGoal<ValueType>(this->getModel(), checkTask), this->getModel().getTransitionMatrix(), checkTask.isRewardModelSet() ? this->getModel().getRewardModel(checkTask.getRewardModel()) : this->getModel().getRewardModel(""), rewardPathFormula.getNonStrictBound<uint64_t>(), *linearEquationSolverFactory);
+                return std::unique_ptr<CheckResult>(new ExplicitQuantitativeCheckResult<ValueType>(std::move(numericResult)));
+            }
         }
         
         template<typename SparseDtmcModelType>

src/storm/modelchecker/prctl/helper/SparseDtmcPrctlHelper.cpp

@@ -15,6 +15,18 @@
 #include "storm/modelchecker/results/ExplicitQuantitativeCheckResult.h"
 #include "storm/modelchecker/hints/ExplicitModelCheckerHint.h"
 #include "storm/modelchecker/prctl/helper/DsMpiUpperRewardBoundsComputer.h"
+#include "storm/modelchecker/prctl/helper/rewardbounded/MultiDimensionalRewardUnfolding.h"
+
+#include "storm/environment/solver/SolverEnvironment.h"
+
+#include "storm/settings/SettingsManager.h"
+#include "storm/settings/modules/GeneralSettings.h"
+#include "storm/settings/modules/CoreSettings.h"
+#include "storm/settings/modules/IOSettings.h"
+
+#include "storm/utility/Stopwatch.h"
+#include "storm/utility/ProgressMeasurement.h"
+#include "storm/utility/export.h"
 
 #include "storm/utility/macros.h"
 #include "storm/utility/ConstantsComparator.h"
@@ -29,12 +41,12 @@ namespace storm {
     namespace modelchecker {
         namespace helper {
             template<typename ValueType, typename RewardModelType>
-            std::vector<ValueType> SparseDtmcPrctlHelper<ValueType, RewardModelType>::computeBoundedUntilProbabilities(Environment const& env, storm::solver::SolveGoal<ValueType>&& goal, storm::storage::SparseMatrix<ValueType> const& transitionMatrix, storm::storage::SparseMatrix<ValueType> const& backwardTransitions, storm::storage::BitVector const& phiStates, storm::storage::BitVector const& psiStates, uint_fast64_t stepBound, storm::solver::LinearEquationSolverFactory<ValueType> const& linearEquationSolverFactory, ModelCheckerHint const& hint) {
+            std::vector<ValueType> SparseDtmcPrctlHelper<ValueType, RewardModelType>::computeStepBoundedUntilProbabilities(Environment const& env, storm::solver::SolveGoal<ValueType>&& goal, storm::storage::SparseMatrix<ValueType> const& transitionMatrix, storm::storage::SparseMatrix<ValueType> const& backwardTransitions, storm::storage::BitVector const& phiStates, storm::storage::BitVector const& psiStates, uint_fast64_t stepBound, storm::solver::LinearEquationSolverFactory<ValueType> const& linearEquationSolverFactory, ModelCheckerHint const& hint) {
                 std::vector<ValueType> result(transitionMatrix.getRowCount(), storm::utility::zero<ValueType>());
                 
                 // If we identify the states that have probability 0 of reaching the target states, we can exclude them in the further analysis.
                 storm::storage::BitVector maybeStates;
-                  
+                
                 if (hint.isExplicitModelCheckerHint() && hint.template asExplicitModelCheckerHint<ValueType>().getComputeOnlyMaybeStates()) {
                     maybeStates = hint.template asExplicitModelCheckerHint<ValueType>().getMaybeStates();
                 } else {
@@ -68,8 +80,168 @@ namespace storm {
                 return result;
             }
             
+            template<typename ValueType>
+            std::vector<ValueType> analyzeTrivialDtmcEpochModel(typename rewardbounded::MultiDimensionalRewardUnfolding<ValueType, true>::EpochModel& epochModel) {
+                
+                std::vector<ValueType> epochResult;
+                epochResult.reserve(epochModel.epochInStates.getNumberOfSetBits());
+                auto stepSolutionIt = epochModel.stepSolutions.begin();
+                auto stepChoiceIt = epochModel.stepChoices.begin();
+                for (auto const& state : epochModel.epochInStates) {
+                    while (*stepChoiceIt < state) {
+                        ++stepChoiceIt;
+                        ++stepSolutionIt;
+                    }
+                    if (epochModel.objectiveRewardFilter.front().get(state)) {
+                        if (*stepChoiceIt == state) {
+                            epochResult.push_back(epochModel.objectiveRewards.front()[state] + *stepSolutionIt);
+                        } else {
+                            epochResult.push_back(epochModel.objectiveRewards.front()[state]);
+                        }
+                    } else {
+                        if (*stepChoiceIt == state) {
+                            epochResult.push_back(*stepSolutionIt);
+                        } else {
+                            epochResult.push_back(storm::utility::zero<ValueType>());
+                        }
+                    }
+                }
+                return epochResult;
+            }
+            
+            template<typename ValueType>
+            std::vector<ValueType> analyzeNonTrivialDtmcEpochModel(Environment const& env, typename rewardbounded::MultiDimensionalRewardUnfolding<ValueType, true>::EpochModel& epochModel, std::vector<ValueType>& x, std::vector<ValueType>& b, std::unique_ptr<storm::solver::LinearEquationSolver<ValueType>>& linEqSolver, storm::solver::LinearEquationSolverFactory<ValueType> const& linearEquationSolverFactory, boost::optional<ValueType> const& lowerBound, boost::optional<ValueType> const& upperBound) {
+ 
+                // Update some data for the case that the Matrix has changed
+                if (epochModel.epochMatrixChanged) {
+                    x.assign(epochModel.epochMatrix.getRowGroupCount(), storm::utility::zero<ValueType>());
+                    linEqSolver = linearEquationSolverFactory.create(env, epochModel.epochMatrix, storm::solver::LinearEquationSolverTask::SolveEquations);
+                    linEqSolver->setCachingEnabled(true);
+                    auto req = linEqSolver->getRequirements(env, storm::solver::LinearEquationSolverTask::SolveEquations);
+                    if (lowerBound) {
+                        linEqSolver->setLowerBound(lowerBound.get());
+                        req.clearLowerBounds();
+                    }
+                    if (upperBound) {
+                        linEqSolver->setUpperBound(upperBound.get());
+                        req.clearUpperBounds();
+                    }
+                    STORM_LOG_THROW(req.empty(), storm::exceptions::UncheckedRequirementException, "At least one requirement was not checked.");
+                }
+                
+                // Prepare the right hand side of the equation system
+                b.assign(epochModel.epochMatrix.getRowCount(), storm::utility::zero<ValueType>());
+                std::vector<ValueType> const& objectiveValues = epochModel.objectiveRewards.front();
+                for (auto const& choice : epochModel.objectiveRewardFilter.front()) {
+                    b[choice] = objectiveValues[choice];
+                }
+                auto stepSolutionIt = epochModel.stepSolutions.begin();
+                for (auto const& choice : epochModel.stepChoices) {
+                    b[choice] += *stepSolutionIt;
+                    ++stepSolutionIt;
+                }
+                assert(stepSolutionIt == epochModel.stepSolutions.end());
+                
+                // Solve the minMax equation system
+                linEqSolver->solveEquations(env, x, b);
+                
+                return storm::utility::vector::filterVector(x, epochModel.epochInStates);
+            }
+            
+            template<>
+            std::map<storm::storage::sparse::state_type, storm::RationalFunction> SparseDtmcPrctlHelper<storm::RationalFunction>::computeRewardBoundedValues(Environment const& env, storm::models::sparse::Dtmc<storm::RationalFunction> const& model, std::shared_ptr<storm::logic::OperatorFormula const> rewardBoundedFormula, storm::solver::LinearEquationSolverFactory<storm::RationalFunction> const& linearEquationSolverFactory) {
+                STORM_LOG_THROW(false, storm::exceptions::NotSupportedException, "The specified property is not supported by this value type.");
+                return std::map<storm::storage::sparse::state_type, storm::RationalFunction>();
+            }
+            
             template<typename ValueType, typename RewardModelType>
+            std::map<storm::storage::sparse::state_type, ValueType> SparseDtmcPrctlHelper<ValueType, RewardModelType>::computeRewardBoundedValues(Environment const& env, storm::models::sparse::Dtmc<ValueType> const& model, std::shared_ptr<storm::logic::OperatorFormula const> rewardBoundedFormula, storm::solver::LinearEquationSolverFactory<ValueType> const& linearEquationSolverFactory) {
+                storm::utility::Stopwatch swAll(true), swBuild, swCheck;
+                
+                storm::modelchecker::helper::rewardbounded::MultiDimensionalRewardUnfolding<ValueType, true> rewardUnfolding(model, rewardBoundedFormula);
+                
+                // Get lower and upper bounds for the solution.
+                auto lowerBound = rewardUnfolding.getLowerObjectiveBound();
+                auto upperBound = rewardUnfolding.getUpperObjectiveBound();
+                
+                // Initialize epoch models
+                auto initEpoch = rewardUnfolding.getStartEpoch();
+                auto epochOrder = rewardUnfolding.getEpochComputationOrder(initEpoch);
+                
+                // initialize data that will be needed for each epoch
+                std::vector<ValueType> x, b;
+                std::unique_ptr<storm::solver::LinearEquationSolver<ValueType>> linEqSolver;
 
+                Environment preciseEnv = env;
+                ValueType precision = rewardUnfolding.getRequiredEpochModelPrecision(initEpoch, storm::utility::convertNumber<ValueType>(storm::settings::getModule<storm::settings::modules::GeneralSettings>().getPrecision()));
+                preciseEnv.solver().setLinearEquationSolverPrecision(storm::utility::convertNumber<storm::RationalNumber>(precision));
+                
+                // In case of cdf export we store the necessary data.
+                std::vector<std::vector<ValueType>> cdfData;
+
+                // Set the correct equation problem format
+                rewardUnfolding.setEquationSystemFormatForEpochModel(linearEquationSolverFactory.getEquationProblemFormat(preciseEnv));
+                bool convertToEquationSystem = linearEquationSolverFactory.getEquationProblemFormat(preciseEnv) == solver::LinearEquationSolverProblemFormat::EquationSystem;
+                
+                storm::utility::ProgressMeasurement progress("epochs");
+                progress.setMaxCount(epochOrder.size());
+                progress.startNewMeasurement(0);
+                uint64_t numCheckedEpochs = 0;
+                for (auto const& epoch : epochOrder) {
+                    swBuild.start();
+                    auto& epochModel = rewardUnfolding.setCurrentEpoch(epoch);
+                    swBuild.stop(); swCheck.start();
+                    // If the epoch matrix is empty we do not need to solve a linear equation system
+                    if ((convertToEquationSystem && epochModel.epochMatrix.isIdentityMatrix()) || (!convertToEquationSystem && epochModel.epochMatrix.getEntryCount() == 0)) {
+                        rewardUnfolding.setSolutionForCurrentEpoch(analyzeTrivialDtmcEpochModel<ValueType>(epochModel));
+                    } else {
+                        rewardUnfolding.setSolutionForCurrentEpoch(analyzeNonTrivialDtmcEpochModel<ValueType>(preciseEnv, epochModel, x, b, linEqSolver, linearEquationSolverFactory, lowerBound, upperBound));
+                    }
+                    swCheck.stop();
+                    if (storm::settings::getModule<storm::settings::modules::IOSettings>().isExportCdfSet() && !rewardUnfolding.getEpochManager().hasBottomDimension(epoch)) {
+                        std::vector<ValueType> cdfEntry;
+                        for (uint64_t i = 0; i < rewardUnfolding.getEpochManager().getDimensionCount(); ++i) {
+                            uint64_t offset = rewardUnfolding.getDimension(i).isUpperBounded ? 0 : 1;
+                            cdfEntry.push_back(storm::utility::convertNumber<ValueType>(rewardUnfolding.getEpochManager().getDimensionOfEpoch(epoch, i) + offset) * rewardUnfolding.getDimension(i).scalingFactor);
+                        }
+                        cdfEntry.push_back(rewardUnfolding.getInitialStateResult(epoch));
+                        cdfData.push_back(std::move(cdfEntry));
+                    }
+                    ++numCheckedEpochs;
+                    progress.updateProgress(numCheckedEpochs);
+                }
+                
+                std::map<storm::storage::sparse::state_type, ValueType> result;
+                for (auto const& initState : model.getInitialStates()) {
+                    result[initState] = rewardUnfolding.getInitialStateResult(initEpoch, initState);
+                }
+                
+                swAll.stop();
+                
+                if (storm::settings::getModule<storm::settings::modules::IOSettings>().isExportCdfSet()) {
+                    std::vector<std::string> headers;
+                    for (uint64_t i = 0; i < rewardUnfolding.getEpochManager().getDimensionCount(); ++i) {
+                        headers.push_back(rewardUnfolding.getDimension(i).formula->toString());
+                    }
+                    headers.push_back("Result");
+                    storm::utility::exportDataToCSVFile<ValueType, std::string, std::string>(storm::settings::getModule<storm::settings::modules::IOSettings>().getExportCdfDirectory() + "cdf.csv", cdfData, headers);
+                }
+
+                if (storm::settings::getModule<storm::settings::modules::CoreSettings>().isShowStatisticsSet()) {
+                    STORM_PRINT_AND_LOG("---------------------------------" << std::endl);
+                    STORM_PRINT_AND_LOG("Statistics:" << std::endl);
+                    STORM_PRINT_AND_LOG("---------------------------------" << std::endl);
+                    STORM_PRINT_AND_LOG("          #checked epochs: " << epochOrder.size() << "." << std::endl);
+                    STORM_PRINT_AND_LOG("             overall Time: " << swAll << "." << std::endl);
+                    STORM_PRINT_AND_LOG("Epoch Model building Time: " << swBuild << "." << std::endl);
+                    STORM_PRINT_AND_LOG("Epoch Model checking Time: " << swCheck << "." << std::endl);
+                    STORM_PRINT_AND_LOG("---------------------------------" << std::endl);
+                }
+                
+                return  result;
+            }
+            
+            template<typename ValueType, typename RewardModelType>
             std::vector<ValueType> SparseDtmcPrctlHelper<ValueType, RewardModelType>::computeUntilProbabilities(Environment const& env, storm::solver::SolveGoal<ValueType>&& goal, storm::storage::SparseMatrix<ValueType> const& transitionMatrix, storm::storage::SparseMatrix<ValueType> const& backwardTransitions, storm::storage::BitVector const& phiStates, storm::storage::BitVector const& psiStates, bool qualitative, storm::solver::LinearEquationSolverFactory<ValueType> const& linearEquationSolverFactory, ModelCheckerHint const& hint) {
                 
                 std::vector<ValueType> result(transitionMatrix.getRowCount(), storm::utility::zero<ValueType>());

src/storm/modelchecker/prctl/helper/SparseDtmcPrctlHelper.h

@@ -5,6 +5,7 @@
 
 #include <boost/optional.hpp>
 
+#include "storm/models/sparse/Dtmc.h"
 #include "storm/models/sparse/StandardRewardModel.h"
 #include "storm/modelchecker/hints/ModelCheckerHint.h"
 
@@ -26,7 +27,10 @@ namespace storm {
             template <typename ValueType, typename RewardModelType = storm::models::sparse::StandardRewardModel<ValueType>>
             class SparseDtmcPrctlHelper {
             public:
-                static std::vector<ValueType> computeBoundedUntilProbabilities(Environment const& env, storm::solver::SolveGoal<ValueType>&& goal, storm::storage::SparseMatrix<ValueType> const& transitionMatrix, storm::storage::SparseMatrix<ValueType> const& backwardTransitions, storm::storage::BitVector const& phiStates, storm::storage::BitVector const& psiStates, uint_fast64_t stepBound, storm::solver::LinearEquationSolverFactory<ValueType> const& linearEquationSolverFactory, ModelCheckerHint const& hint = ModelCheckerHint());
+
+                static std::vector<ValueType> computeStepBoundedUntilProbabilities(Environment const& env, storm::solver::SolveGoal<ValueType>&& goal, storm::storage::SparseMatrix<ValueType> const& transitionMatrix, storm::storage::SparseMatrix<ValueType> const& backwardTransitions, storm::storage::BitVector const& phiStates, storm::storage::BitVector const& psiStates, uint_fast64_t stepBound, storm::solver::LinearEquationSolverFactory<ValueType> const& linearEquationSolverFactory, ModelCheckerHint const& hint = ModelCheckerHint());
+                
+                static std::map<storm::storage::sparse::state_type, ValueType> computeRewardBoundedValues(Environment const& env, storm::models::sparse::Dtmc<ValueType> const& model, std::shared_ptr<storm::logic::OperatorFormula const> rewardBoundedFormula, storm::solver::LinearEquationSolverFactory<ValueType> const& linearEquationSolverFactory);
                 
                 static std::vector<ValueType> computeNextProbabilities(Environment const& env, storm::storage::SparseMatrix<ValueType> const& transitionMatrix, storm::storage::BitVector const& nextStates, storm::solver::LinearEquationSolverFactory<ValueType> const& linearEquationSolverFactory);
                 

src/storm/modelchecker/prctl/helper/SparseMdpPrctlHelper.cpp

@@ -87,7 +87,7 @@ namespace storm {
             }
 
             template<typename ValueType>
-            std::vector<ValueType> analyzeTrivialEpochModel(OptimizationDirection dir, typename rewardbounded::MultiDimensionalRewardUnfolding<ValueType, true>::EpochModel& epochModel) {
+            std::vector<ValueType> analyzeTrivialMdpEpochModel(OptimizationDirection dir, typename rewardbounded::MultiDimensionalRewardUnfolding<ValueType, true>::EpochModel& epochModel) {
                 // Assert that the epoch model is indeed trivial
                 assert(epochModel.epochMatrix.getEntryCount() == 0);
                 
@@ -137,7 +137,7 @@ namespace storm {
             }
             
             template<typename ValueType>
-            std::vector<ValueType> analyzeNonTrivialEpochModel(Environment const& env, OptimizationDirection dir, typename rewardbounded::MultiDimensionalRewardUnfolding<ValueType, true>::EpochModel& epochModel, std::vector<ValueType>& x, std::vector<ValueType>& b, std::unique_ptr<storm::solver::MinMaxLinearEquationSolver<ValueType>>& minMaxSolver, storm::solver::MinMaxLinearEquationSolverFactory<ValueType> const& minMaxLinearEquationSolverFactory, boost::optional<ValueType> const& lowerBound, boost::optional<ValueType> const& upperBound) {
+            std::vector<ValueType> analyzeNonTrivialMdpEpochModel(Environment const& env, OptimizationDirection dir, typename rewardbounded::MultiDimensionalRewardUnfolding<ValueType, true>::EpochModel& epochModel, std::vector<ValueType>& x, std::vector<ValueType>& b, std::unique_ptr<storm::solver::MinMaxLinearEquationSolver<ValueType>>& minMaxSolver, storm::solver::MinMaxLinearEquationSolverFactory<ValueType> const& minMaxLinearEquationSolverFactory, boost::optional<ValueType> const& lowerBound, boost::optional<ValueType> const& upperBound) {
  
                 // Update some data for the case that the Matrix has changed
                 if (epochModel.epochMatrixChanged) {
@@ -215,9 +215,9 @@ namespace storm {
                     swBuild.stop(); swCheck.start();
                     // If the epoch matrix is empty we do not need to solve a linear equation system
                     if (epochModel.epochMatrix.getEntryCount() == 0) {
-                        rewardUnfolding.setSolutionForCurrentEpoch(analyzeTrivialEpochModel<ValueType>(dir, epochModel));
+                        rewardUnfolding.setSolutionForCurrentEpoch(analyzeTrivialMdpEpochModel<ValueType>(dir, epochModel));
                     } else {
-                        rewardUnfolding.setSolutionForCurrentEpoch(analyzeNonTrivialEpochModel<ValueType>(preciseEnv, dir, epochModel, x, b, minMaxSolver, minMaxLinearEquationSolverFactory, lowerBound, upperBound));
+                        rewardUnfolding.setSolutionForCurrentEpoch(analyzeNonTrivialMdpEpochModel<ValueType>(preciseEnv, dir, epochModel, x, b, minMaxSolver, minMaxLinearEquationSolverFactory, lowerBound, upperBound));
                     }
                     swCheck.stop();
                     if (storm::settings::getModule<storm::settings::modules::IOSettings>().isExportCdfSet() && !rewardUnfolding.getEpochManager().hasBottomDimension(epoch)) {

src/storm/modelchecker/prctl/helper/rewardbounded/MultiDimensionalRewardUnfolding.cpp

@@ -10,6 +10,8 @@
 #include "storm/modelchecker/propositional/SparsePropositionalModelChecker.h"
 #include "storm/modelchecker/results/ExplicitQualitativeCheckResult.h"
 #include "storm/modelchecker/prctl/helper/BaierUpperRewardBoundsComputer.h"
+#include "storm/models/sparse/Mdp.h"
+#include "storm/models/sparse/Dtmc.h"
 
 #include "storm/transformer/EndComponentEliminator.h"
 
@@ -24,14 +26,13 @@ namespace storm {
             namespace rewardbounded {
     
                 template<typename ValueType, bool SingleObjectiveMode>
-                MultiDimensionalRewardUnfolding<ValueType, SingleObjectiveMode>::MultiDimensionalRewardUnfolding(storm::models::sparse::Mdp<ValueType> const& model, std::vector<storm::modelchecker::multiobjective::Objective<ValueType>> const& objectives) : model(model), objectives(objectives) {
+                MultiDimensionalRewardUnfolding<ValueType, SingleObjectiveMode>::MultiDimensionalRewardUnfolding(storm::models::sparse::Model<ValueType> const& model, std::vector<storm::modelchecker::multiobjective::Objective<ValueType>> const& objectives) : model(model), objectives(objectives) {
                     initialize();
                 }
                 
                 template<typename ValueType, bool SingleObjectiveMode>
-                MultiDimensionalRewardUnfolding<ValueType, SingleObjectiveMode>::MultiDimensionalRewardUnfolding(storm::models::sparse::Mdp<ValueType> const& model, std::shared_ptr<storm::logic::OperatorFormula const> objectiveFormula) : model(model) {
+                MultiDimensionalRewardUnfolding<ValueType, SingleObjectiveMode>::MultiDimensionalRewardUnfolding(storm::models::sparse::Model<ValueType> const& model, std::shared_ptr<storm::logic::OperatorFormula const> objectiveFormula) : model(model) {
                     
-                    STORM_LOG_THROW(objectiveFormula->hasOptimalityType(), storm::exceptions::InvalidPropertyException, "Formula needs to specify whether minimal or maximal values are to be computed on nondeterministic model.");
                     if (objectiveFormula->isProbabilityOperatorFormula()) {
                         if (objectiveFormula->getSubformula().isMultiObjectiveFormula()) {
                             for (auto const& subFormula : objectiveFormula->getSubformula().asMultiObjectiveFormula().getSubformulas()) {
@@ -91,7 +92,7 @@ namespace storm {
                                 // lower bounded until formulas with non-trivial left hand side are excluded as this would require some additional effort (in particular the ProductModel::transformMemoryState method).
                                 STORM_LOG_THROW(dimension.isUpperBounded || subformula.getLeftSubformula(dim).isTrueFormula(), storm::exceptions::NotSupportedException, "Lower bounded until formulas are only supported by this method if the left subformula is 'true'. Got " << subformula << " instead.");
                                 if (subformula.getTimeBoundReference(dim).isTimeBound() || subformula.getTimeBoundReference(dim).isStepBound()) {
-                                    dimensionWiseEpochSteps.push_back(std::vector<uint64_t>(model.getNumberOfChoices(), 1));
+                                    dimensionWiseEpochSteps.push_back(std::vector<uint64_t>(model.getTransitionMatrix().getRowCount(), 1));
                                     dimension.scalingFactor = storm::utility::one<ValueType>();
                                 } else {
                                     STORM_LOG_ASSERT(subformula.getTimeBoundReference(dim).isRewardBound(), "Unexpected type of time bound.");
@@ -114,7 +115,7 @@ namespace storm {
                                 dimension.objectiveIndex = objIndex;
                                 dimension.isUpperBounded = true;
                                 if (subformula.getTimeBoundReference(dim).isTimeBound() || subformula.getTimeBoundReference(dim).isStepBound()) {
-                                    dimensionWiseEpochSteps.push_back(std::vector<uint64_t>(model.getNumberOfChoices(), 1));
+                                    dimensionWiseEpochSteps.push_back(std::vector<uint64_t>(model.getTransitionMatrix().getRowCount(), 1));
                                     dimension.scalingFactor = storm::utility::one<ValueType>();
                                 } else {
                                     STORM_LOG_ASSERT(subformula.getTimeBoundReference(dim).isRewardBound(), "Unexpected type of time bound.");
@@ -166,8 +167,8 @@ namespace storm {
                     epochManager = EpochManager(dimensions.size());
                     
                     // Convert the epoch steps to a choice-wise representation
-                    epochSteps.reserve(model.getNumberOfChoices());
-                    for (uint64_t choice = 0; choice < model.getNumberOfChoices(); ++choice) {
+                    epochSteps.reserve(model.getTransitionMatrix().getRowCount());
+                    for (uint64_t choice = 0; choice < model.getTransitionMatrix().getRowCount(); ++choice) {
                         Epoch step;
                         uint64_t dim = 0;
                         for (auto const& dimensionSteps : dimensionWiseEpochSteps) {
@@ -396,7 +397,7 @@ namespace storm {
                     // std::cout << "Setting epoch class for epoch " << epochManager.toString(epoch) << std::endl;
                     auto productObjectiveRewards = productModel->computeObjectiveRewards(epochClass, objectives);
                     
-                    storm::storage::BitVector stepChoices(productModel->getProduct().getNumberOfChoices(), false);
+                    storm::storage::BitVector stepChoices(productModel->getProduct().getTransitionMatrix().getRowCount(), false);
                     uint64_t choice = 0;
                     for (auto const& step : productModel->getSteps()) {
                         if (!epochManager.isZeroEpoch(step) && epochManager.getSuccessorEpoch(epoch, step) != epoch) {
@@ -421,7 +422,7 @@ namespace storm {
                         }
                     }
                     
-                    storm::storage::BitVector zeroObjRewardChoices(productModel->getProduct().getNumberOfChoices(), true);
+                    storm::storage::BitVector zeroObjRewardChoices(productModel->getProduct().getTransitionMatrix().getRowCount(), true);
                     for (uint64_t objIndex = 0; objIndex < objectives.size(); ++objIndex) {
                         if (violatedLowerBoundedDimensions.isDisjointFrom(objectiveDimensions[objIndex])) {
                             zeroObjRewardChoices &= storm::utility::vector::filterZero(productObjectiveRewards[objIndex]);
@@ -433,15 +434,68 @@ namespace storm {
                     storm::storage::BitVector productInStates = productModel->getInStates(epochClass);
                     // The epoch model only needs to consider the states that are reachable from a relevant state
                     storm::storage::BitVector consideredStates = storm::utility::graph::getReachableStates(epochModel.epochMatrix, productInStates, allProductStates, ~allProductStates);
-                    // std::cout << "numInStates = " << productInStates.getNumberOfSetBits() << std::endl;
-                    // std::cout << "numConsideredStates = " << consideredStates.getNumberOfSetBits() << std::endl;
                     
                     // We assume that there is no end component in which objective reward is earned
                     STORM_LOG_ASSERT(!storm::utility::graph::checkIfECWithChoiceExists(epochModel.epochMatrix, epochModel.epochMatrix.transpose(true), allProductStates, ~zeroObjRewardChoices & ~stepChoices), "There is a scheduler that yields infinite reward for one objective. This case should be excluded");
-                    auto ecElimResult = storm::transformer::EndComponentEliminator<ValueType>::transform(epochModel.epochMatrix, consideredStates, zeroObjRewardChoices & ~stepChoices, consideredStates);
-                    epochModel.epochMatrix = std::move(ecElimResult.matrix);
-                    epochModelToProductChoiceMap = std::move(ecElimResult.newToOldRowMapping);
                     
+                    // Create the epoch model matrix
+                    std::vector<uint64_t> productToEpochModelStateMapping;
+                    if (model.isOfType(storm::models::ModelType::Dtmc)) {
+                        assert(zeroObjRewardChoices.size() == productModel->getProduct().getNumberOfStates());
+                        assert(stepChoices.size() == productModel->getProduct().getNumberOfStates());
+                        STORM_LOG_ASSERT(equationSolverProblemFormatForEpochModel.is_initialized(), "Linear equation problem format was not set.");
+                        bool convertToEquationSystem = equationSolverProblemFormatForEpochModel.get() == storm::solver::LinearEquationSolverProblemFormat::EquationSystem;
+                        // For DTMCs we consider the subsystem induced by the considered states.
+                        // The transitions for states with zero reward are filtered out to guarantee a unique solution of the eq-system.
+                        auto backwardTransitions = epochModel.epochMatrix.transpose(true);
+                        storm::storage::BitVector nonZeroRewardStates = storm::utility::graph::performProbGreater0(backwardTransitions, consideredStates, consideredStates & (~zeroObjRewardChoices | stepChoices));
+                        // If there is at least one considered state with reward zero, we have to add a 'zero-reward-state' to the epoch model.
+                        bool requiresZeroRewardState = nonZeroRewardStates != consideredStates;
+                        uint64_t numEpochModelStates = nonZeroRewardStates.getNumberOfSetBits();
+                        uint64_t zeroRewardInState = numEpochModelStates;
+                        if (requiresZeroRewardState) {
+                            ++numEpochModelStates;
+                        }
+                        storm::storage::SparseMatrixBuilder<ValueType> builder;
+                        if (!nonZeroRewardStates.empty()) {
+                            builder = storm::storage::SparseMatrixBuilder<ValueType>(epochModel.epochMatrix.getSubmatrix(true, nonZeroRewardStates, nonZeroRewardStates, convertToEquationSystem));
+                        }
+                        if (requiresZeroRewardState) {
+                            if (convertToEquationSystem) {
+                                // add a diagonal entry
+                                builder.addNextValue(zeroRewardInState, zeroRewardInState, storm::utility::zero<ValueType>());
+                            }
+                            epochModel.epochMatrix = builder.build(numEpochModelStates, numEpochModelStates);
+                        } else {
+                            assert (!nonZeroRewardStates.empty());
+                            epochModel.epochMatrix = builder.build();
+                        }
+                        if (convertToEquationSystem) {
+                            epochModel.epochMatrix.convertToEquationSystem();
+                        }
+                        
+                        epochModelToProductChoiceMap.clear();
+                        epochModelToProductChoiceMap.reserve(numEpochModelStates);
+                        productToEpochModelStateMapping.assign(nonZeroRewardStates.size(), zeroRewardInState);
+                        for (auto const& productState : nonZeroRewardStates) {
+                            productToEpochModelStateMapping[productState] = epochModelToProductChoiceMap.size();
+                            epochModelToProductChoiceMap.push_back(productState);
+                        }
+                        if (requiresZeroRewardState) {
+                            uint64_t zeroRewardProductState = (consideredStates & ~nonZeroRewardStates).getNextSetIndex(0);
+                            assert(zeroRewardProductState < consideredStates.size());
+                            epochModelToProductChoiceMap.push_back(zeroRewardProductState);
+                        }
+                    } else if (model.isOfType(storm::models::ModelType::Mdp)) {
+                        // Eliminate zero-reward end components
+                        auto ecElimResult = storm::transformer::EndComponentEliminator<ValueType>::transform(epochModel.epochMatrix, consideredStates, zeroObjRewardChoices & ~stepChoices, consideredStates);
+                        epochModel.epochMatrix = std::move(ecElimResult.matrix);
+                        epochModelToProductChoiceMap = std::move(ecElimResult.newToOldRowMapping);
+                        productToEpochModelStateMapping = std::move(ecElimResult.oldToNewStateMapping);
+                    } else {
+                        STORM_LOG_THROW(false, storm::exceptions::UnexpectedException, "Unsupported model type.");
+                    }
+                
                     epochModel.stepChoices = storm::storage::BitVector(epochModel.epochMatrix.getRowCount(), false);
                     for (uint64_t choice = 0; choice < epochModel.epochMatrix.getRowCount(); ++choice) {
                         if (stepChoices.get(epochModelToProductChoiceMap[choice])) {
@@ -466,14 +520,14 @@ namespace storm {
                     
                     epochModel.epochInStates = storm::storage::BitVector(epochModel.epochMatrix.getRowGroupCount(), false);
                     for (auto const& productState : productInStates) {
-                        STORM_LOG_ASSERT(ecElimResult.oldToNewStateMapping[productState] < epochModel.epochMatrix.getRowGroupCount(), "Selected product state does not exist in the epoch model.");
-                        epochModel.epochInStates.set(ecElimResult.oldToNewStateMapping[productState], true);
+                        STORM_LOG_ASSERT(productToEpochModelStateMapping[productState] < epochModel.epochMatrix.getRowGroupCount(), "Selected product state does not exist in the epoch model.");
+                        epochModel.epochInStates.set(productToEpochModelStateMapping[productState], true);
                     }
                     
                     std::vector<uint64_t> toEpochModelInStatesMap(productModel->getProduct().getNumberOfStates(), std::numeric_limits<uint64_t>::max());
                     std::vector<uint64_t> epochModelStateToInStateMap = epochModel.epochInStates.getNumberOfSetBitsBeforeIndices();
                     for (auto const& productState : productInStates) {
-                        toEpochModelInStatesMap[productState] = epochModelStateToInStateMap[ecElimResult.oldToNewStateMapping[productState]];
+                        toEpochModelInStatesMap[productState] = epochModelStateToInStateMap[productToEpochModelStateMapping[productState]];
                     }
                     productStateToEpochModelInStateMap = std::make_shared<std::vector<uint64_t> const>(std::move(toEpochModelInStatesMap));
                     
@@ -487,6 +541,13 @@ namespace storm {
                 }
                 
      
+                template<typename ValueType, bool SingleObjectiveMode>
+                void MultiDimensionalRewardUnfolding<ValueType, SingleObjectiveMode>::setEquationSystemFormatForEpochModel(storm::solver::LinearEquationSolverProblemFormat eqSysFormat) {
+                    STORM_LOG_ASSERT(model.isOfType(storm::models::ModelType::Dtmc), "Trying to set the equation problem format although the model is not deterministic.");
+                    equationSolverProblemFormatForEpochModel = eqSysFormat;
+                }
+
+                
                 template<typename ValueType, bool SingleObjectiveMode>
                 template<bool SO, typename std::enable_if<SO, int>::type>
                 typename MultiDimensionalRewardUnfolding<ValueType, SingleObjectiveMode>::SolutionType MultiDimensionalRewardUnfolding<ValueType, SingleObjectiveMode>::getScaledSolution(SolutionType const& solution, ValueType const& scalingFactor) const {

src/storm/modelchecker/prctl/helper/rewardbounded/MultiDimensionalRewardUnfolding.h

@@ -8,7 +8,8 @@
 #include "storm/modelchecker/prctl/helper/rewardbounded/EpochManager.h"
 #include "storm/modelchecker/prctl/helper/rewardbounded/ProductModel.h"
 #include "storm/modelchecker/prctl/helper/rewardbounded/Dimension.h"
-#include "storm/models/sparse/Mdp.h"
+#include "storm/models/sparse/Model.h"
+#include "storm/solver/LinearEquationSolverProblemFormat.h"
 #include "storm/utility/vector.h"
 #include "storm/storage/memorystructure/MemoryStructure.h"
 #include "storm/utility/Stopwatch.h"
@@ -43,8 +44,8 @@ namespace storm {
                      * @param objectives The (preprocessed) objectives
                      *
                      */
-                    MultiDimensionalRewardUnfolding(storm::models::sparse::Mdp<ValueType> const& model, std::vector<storm::modelchecker::multiobjective::Objective<ValueType>> const& objectives);
-                    MultiDimensionalRewardUnfolding(storm::models::sparse::Mdp<ValueType> const& model, std::shared_ptr<storm::logic::OperatorFormula const> objectiveFormula);
+                    MultiDimensionalRewardUnfolding(storm::models::sparse::Model<ValueType> const& model, std::vector<storm::modelchecker::multiobjective::Objective<ValueType>> const& objectives);
+                    MultiDimensionalRewardUnfolding(storm::models::sparse::Model<ValueType> const& model, std::shared_ptr<storm::logic::OperatorFormula const> objectiveFormula);
                     
                     ~MultiDimensionalRewardUnfolding() = default;
                     
@@ -53,6 +54,8 @@ namespace storm {
                     
                     EpochModel& setCurrentEpoch(Epoch const& epoch);
                     
+                    void setEquationSystemFormatForEpochModel(storm::solver::LinearEquationSolverProblemFormat eqSysFormat);
+                    
                     /*!
                      * Returns the precision required for the analyzis of each epoch model in order to achieve the given overall precision
                      */
@@ -106,7 +109,7 @@ namespace storm {
                     EpochSolution const& getEpochSolution(std::map<Epoch, EpochSolution const*> const& solutions, Epoch const& epoch);
                     SolutionType const& getStateSolution(EpochSolution const& epochSolution, uint64_t const& productState);
                     
-                    storm::models::sparse::Mdp<ValueType> const& model;
+                    storm::models::sparse::Model<ValueType> const& model;
                     std::vector<storm::modelchecker::multiobjective::Objective<ValueType>> objectives;
                     
                     std::unique_ptr<ProductModel<ValueType>> productModel;
@@ -118,6 +121,9 @@ namespace storm {
                     EpochModel epochModel;
                     boost::optional<Epoch> currentEpoch;
                     
+                    // In case of DTMCs we have different options for the equation problem format the epoch model will have.
+                    boost::optional<storm::solver::LinearEquationSolverProblemFormat> equationSolverProblemFormatForEpochModel;
+                    
                     EpochManager epochManager;
                     
                     std::vector<Dimension<ValueType>> dimensions;

src/storm/modelchecker/prctl/helper/rewardbounded/ProductModel.cpp

@@ -9,6 +9,8 @@
 
 #include "storm/modelchecker/propositional/SparsePropositionalModelChecker.h"
 #include "storm/modelchecker/results/ExplicitQualitativeCheckResult.h"
+#include "storm/models/sparse/Mdp.h"
+#include "storm/models/sparse/Dtmc.h"
 
 #include "storm/exceptions/UnexpectedException.h"
 #include "storm/exceptions/NotSupportedException.h"
@@ -19,7 +21,7 @@ namespace storm {
             namespace rewardbounded {
                 
                 template<typename ValueType>
-                ProductModel<ValueType>::ProductModel(storm::models::sparse::Mdp<ValueType> const& model, std::vector<storm::modelchecker::multiobjective::Objective<ValueType>> const& objectives, std::vector<Dimension<ValueType>> const& dimensions, std::vector<storm::storage::BitVector> const& objectiveDimensions, EpochManager const& epochManager, std::vector<Epoch> const& originalModelSteps) : dimensions(dimensions), objectiveDimensions(objectiveDimensions), epochManager(epochManager), memoryStateManager(dimensions.size()) {
+                ProductModel<ValueType>::ProductModel(storm::models::sparse::Model<ValueType> const& model, std::vector<storm::modelchecker::multiobjective::Objective<ValueType>> const& objectives, std::vector<Dimension<ValueType>> const& dimensions, std::vector<storm::storage::BitVector> const& objectiveDimensions, EpochManager const& epochManager, std::vector<Epoch> const& originalModelSteps) : dimensions(dimensions), objectiveDimensions(objectiveDimensions), epochManager(epochManager), memoryStateManager(dimensions.size()) {
                     
                     for (uint64_t dim = 0; dim < dimensions.size(); ++dim) {
                         if (!dimensions[dim].memoryLabel) {
@@ -34,7 +36,7 @@ namespace storm {
                     storm::storage::SparseModelMemoryProduct<ValueType> productBuilder(memory.product(model));
                     
                     setReachableProductStates(productBuilder, originalModelSteps, memoryStateMap);
-                    product = productBuilder.build()->template as<storm::models::sparse::Mdp<ValueType>>();
+                    product = productBuilder.build();
                     
                     uint64_t numModelStates = productBuilder.getOriginalModel().getNumberOfStates();
                     MemoryState upperMemStateBound = memoryStateManager.getUpperMemoryStateBound();
@@ -57,7 +59,7 @@ namespace storm {
                     }
                     
                     // Map choice indices of the product to the state where it origins
-                    choiceToStateMap.reserve(getProduct().getNumberOfChoices());
+                    choiceToStateMap.reserve(getProduct().getTransitionMatrix().getRowCount());
                     for (uint64_t productState = 0; productState < numProductStates; ++productState) {
                         uint64_t groupSize = getProduct().getTransitionMatrix().getRowGroupSize(productState);
                         for (uint64_t i = 0; i < groupSize; ++i) {
@@ -66,7 +68,7 @@ namespace storm {
                     }
                     
                     // Compute the epoch steps for the product
-                    steps.resize(getProduct().getNumberOfChoices(), 0);
+                    steps.resize(getProduct().getTransitionMatrix().getRowCount(), 0);
                     for (uint64_t modelState = 0; modelState < numModelStates; ++modelState) {
                         uint64_t numChoices = productBuilder.getOriginalModel().getTransitionMatrix().getRowGroupSize(modelState);
                         uint64_t firstChoice = productBuilder.getOriginalModel().getTransitionMatrix().getRowGroupIndices()[modelState];
@@ -92,9 +94,9 @@ namespace storm {
                 
                 
                 template<typename ValueType>
-                storm::storage::MemoryStructure ProductModel<ValueType>::computeMemoryStructure(storm::models::sparse::Mdp<ValueType> const& model, std::vector<storm::modelchecker::multiobjective::Objective<ValueType>> const& objectives) const {
+                storm::storage::MemoryStructure ProductModel<ValueType>::computeMemoryStructure(storm::models::sparse::Model<ValueType> const& model, std::vector<storm::modelchecker::multiobjective::Objective<ValueType>> const& objectives) const {
                     
-                    storm::modelchecker::SparsePropositionalModelChecker<storm::models::sparse::Mdp<ValueType>> mc(model);
+                    storm::modelchecker::SparsePropositionalModelChecker<storm::models::sparse::Model<ValueType>> mc(model);
                     
                     // Create a memory structure that remembers whether (sub)objectives are satisfied
                     storm::storage::MemoryStructure memory = storm::storage::MemoryStructureBuilder<ValueType>::buildTrivialMemoryStructure(model);
@@ -285,7 +287,7 @@ namespace storm {
                 }
     
                 template<typename ValueType>
-                storm::models::sparse::Mdp<ValueType> const& ProductModel<ValueType>::getProduct() const {
+                storm::models::sparse::Model<ValueType> const& ProductModel<ValueType>::getProduct() const {
                     return *product;
                 }
                 
@@ -341,7 +343,7 @@ namespace storm {
                     for (uint64_t objIndex = 0; objIndex < objectives.size(); ++objIndex) {
                         auto const& formula = *objectives[objIndex].formula;
                         if (formula.isProbabilityOperatorFormula()) {
-                            storm::modelchecker::SparsePropositionalModelChecker<storm::models::sparse::Mdp<ValueType>> mc(getProduct());
+                            storm::modelchecker::SparsePropositionalModelChecker<storm::models::sparse::Model<ValueType>> mc(getProduct());
                             std::vector<uint64_t> dimensionIndexMap;
                             for (auto const& globalDimensionIndex : objectiveDimensions[objIndex]) {
                                 dimensionIndexMap.push_back(globalDimensionIndex);

src/storm/modelchecker/prctl/helper/rewardbounded/ProductModel.h

@@ -7,7 +7,7 @@
 #include "storm/modelchecker/prctl/helper/rewardbounded/EpochManager.h"
 #include "storm/modelchecker/prctl/helper/rewardbounded/MemoryStateManager.h"
 #include "storm/modelchecker/prctl/helper/rewardbounded/Dimension.h"
-#include "storm/models/sparse/Mdp.h"
+#include "storm/models/sparse/Model.h"
 #include "storm/utility/vector.h"
 #include "storm/storage/memorystructure/MemoryStructure.h"
 #include "storm/storage/memorystructure/SparseModelMemoryProduct.h"
@@ -26,9 +26,9 @@ namespace storm {
                     typedef typename MemoryStateManager::MemoryState MemoryState;
     
                     
-                    ProductModel(storm::models::sparse::Mdp<ValueType> const& model, std::vector<storm::modelchecker::multiobjective::Objective<ValueType>> const& objectives, std::vector<Dimension<ValueType>> const& dimensions, std::vector<storm::storage::BitVector> const& objectiveDimensions, EpochManager const& epochManager, std::vector<Epoch> const& originalModelSteps);
+                    ProductModel(storm::models::sparse::Model<ValueType> const& model, std::vector<storm::modelchecker::multiobjective::Objective<ValueType>> const& objectives, std::vector<Dimension<ValueType>> const& dimensions, std::vector<storm::storage::BitVector> const& objectiveDimensions, EpochManager const& epochManager, std::vector<Epoch> const& originalModelSteps);
                     
-                    storm::models::sparse::Mdp<ValueType> const& getProduct() const;
+                    storm::models::sparse::Model<ValueType> const& getProduct() const;
                     std::vector<Epoch> const& getSteps() const;
                     
                     bool productStateExists(uint64_t const& modelState, uint64_t const& memoryState) const;
@@ -50,7 +50,7 @@ namespace storm {
                     
                 private:
                     
-                    storm::storage::MemoryStructure computeMemoryStructure(storm::models::sparse::Mdp<ValueType> const& model, std::vector<storm::modelchecker::multiobjective::Objective<ValueType>> const& objectives) const;
+                    storm::storage::MemoryStructure computeMemoryStructure(storm::models::sparse::Model<ValueType> const& model, std::vector<storm::modelchecker::multiobjective::Objective<ValueType>> const& objectives) const;
                     std::vector<MemoryState> computeMemoryStateMap(storm::storage::MemoryStructure const& memory) const;
     
                     
@@ -66,7 +66,7 @@ namespace storm {
                     EpochManager const& epochManager;
                     MemoryStateManager memoryStateManager;
     
-                    std::shared_ptr<storm::models::sparse::Mdp<ValueType>> product;
+                    std::shared_ptr<storm::models::sparse::Model<ValueType>> product;
                     std::vector<Epoch> steps;
                     std::map<EpochClass, storm::storage::BitVector> reachableStates;
                     std::map<EpochClass, storm::storage::BitVector> inStates;

src/storm/storage/SparseMatrix.cpp

@@ -2015,6 +2015,35 @@ namespace storm {
             return true;
         }
         
+        template<typename ValueType>
+        bool SparseMatrix<ValueType>::isIdentityMatrix() const {
+            if (this->getRowCount() != this->getColumnCount()) {
+                return false;
+            }
+            if (this->getNonzeroEntryCount() != this->getRowCount()) {
+                return false;
+            }
+            for (uint64_t row = 0; row < this->getRowCount(); ++row) {
+                bool rowHasEntry = false;
+                for (auto const& entry : this->getRow(row)) {
+                    if (entry.getColumn() == row) {
+                        if (!storm::utility::isOne(entry.getValue())) {
+                            return false;
+                        }
+                        rowHasEntry = true;
+                    } else {
+                        if (!storm::utility::isZero(entry.getValue())) {
+                            return false;
+                        }
+                    }
+                }
+                if (!rowHasEntry) {
+                    return false;
+                }
+            }
+            return true;
+        }
+
         template<typename ValueType>
         std::ostream& operator<<(std::ostream& out, SparseMatrix<ValueType> const& matrix) {
             // Print column numbers in header.

src/storm/storage/SparseMatrix.h

@@ -929,6 +929,9 @@ namespace storm {
             template<typename OtherValueType>
             bool isSubmatrixOf(SparseMatrix<OtherValueType> const& matrix) const;
             
+            // Returns true if the matrix is the identity matrix
+            bool isIdentityMatrix() const;
+            
             template<typename TPrime>
             friend std::ostream& operator<<(std::ostream& out, SparseMatrix<TPrime> const& matrix);
 

src/storm/storage/sparse/StateStorage.cpp

@@ -5,7 +5,7 @@ namespace storm {
         namespace sparse {
                         
             template <typename StateType>
-            StateStorage<StateType>::StateStorage(uint64_t bitsPerState) : stateToId(bitsPerState, 10000000), initialStateIndices(), deadlockStateIndices(), bitsPerState(bitsPerState) {
+            StateStorage<StateType>::StateStorage(uint64_t bitsPerState) : stateToId(bitsPerState, 100000), initialStateIndices(), deadlockStateIndices(), bitsPerState(bitsPerState) {
                 // Intentionally left empty.
             }
 

src/test/storm/modelchecker/SparseDtmcMultiDimensionalRewardUnfoldingTest.cpp

@@ -0,0 +1,99 @@
+#include "gtest/gtest.h"
+#include "storm-config.h"
+
+#include "storm/modelchecker/results/ExplicitQuantitativeCheckResult.h"
+#include "storm/models/sparse/Dtmc.h"
+#include "storm/settings/modules/GeneralSettings.h"
+#include "storm/settings/SettingsManager.h"
+#include "storm/utility/constants.h"
+#include "storm/api/storm.h"
+#include "storm/environment/Environment.h"
+
+TEST(SparseDtmcMultiDimensionalRewardUnfoldingTest, cost_bounded_die) {
+    storm::Environment env;
+    std::string programFile = STORM_TEST_RESOURCES_DIR "/dtmc/die.pm";
+    std::string formulasAsString = "P=? [ F{\"coin_flips\"}<=2 \"two\" ] ";
+    formulasAsString += "; P=? [ F{\"coin_flips\"}<=3 \"two\" ] ";
+    formulasAsString += "; P=? [ F{\"coin_flips\"}<=8 \"two\" ] ";
+
+    // programm, model,  formula
+    storm::prism::Program program = storm::api::parseProgram(programFile);
+    program = storm::utility::prism::preprocess(program, "");
+    std::vector<std::shared_ptr<storm::logic::Formula const>> formulas = storm::api::extractFormulasFromProperties(storm::api::parsePropertiesForPrismProgram(formulasAsString, program));
+    std::shared_ptr<storm::models::sparse::Dtmc<storm::RationalNumber>> dtmc = storm::api::buildSparseModel<storm::RationalNumber>(program, formulas)->as<storm::models::sparse::Dtmc<storm::RationalNumber>>();
+    uint_fast64_t const initState = *dtmc->getInitialStates().begin();;
+    std::unique_ptr<storm::modelchecker::CheckResult> result;
+    
+    result = storm::api::verifyWithSparseEngine(dtmc, storm::api::createTask<storm::RationalNumber>(formulas[0], true));
+    ASSERT_TRUE(result->isExplicitQuantitativeCheckResult());
+    EXPECT_EQ(storm::utility::convertNumber<storm::RationalNumber>(std::string("0")), result->asExplicitQuantitativeCheckResult<storm::RationalNumber>()[initState]);
+
+    result = storm::api::verifyWithSparseEngine(dtmc, storm::api::createTask<storm::RationalNumber>(formulas[1], true));
+    ASSERT_TRUE(result->isExplicitQuantitativeCheckResult());
+    EXPECT_EQ(storm::utility::convertNumber<storm::RationalNumber>(std::string("1/8")), result->asExplicitQuantitativeCheckResult<storm::RationalNumber>()[initState]);
+
+    result = storm::api::verifyWithSparseEngine(dtmc, storm::api::createTask<storm::RationalNumber>(formulas[2], true));
+    ASSERT_TRUE(result->isExplicitQuantitativeCheckResult());
+    EXPECT_EQ(storm::utility::convertNumber<storm::RationalNumber>(std::string("21/128")), result->asExplicitQuantitativeCheckResult<storm::RationalNumber>()[initState]);
+}
+
+TEST(SparseDtmcMultiDimensionalRewardUnfoldingTest, cost_bounded_leader) {
+    storm::Environment env;
+    std::string programFile = STORM_TEST_RESOURCES_DIR "/dtmc/leader-3-5.pm";
+    std::string formulasAsString = "P=? [ F{\"num_rounds\"}<=1 \"elected\" ] ";
+    formulasAsString += "; P=? [ F{\"num_rounds\"}<=2 \"elected\" ] ";
+    formulasAsString += "; P=? [ F{\"num_rounds\"}>2 \"elected\" ] ";
+    formulasAsString += "; P=? [ F{\"num_rounds\"}>=2,{\"num_rounds\"}<3 \"elected\" ] ";
+
+    // programm, model,  formula
+    storm::prism::Program program = storm::api::parseProgram(programFile);
+    program = storm::utility::prism::preprocess(program, "");
+    std::vector<std::shared_ptr<storm::logic::Formula const>> formulas = storm::api::extractFormulasFromProperties(storm::api::parsePropertiesForPrismProgram(formulasAsString, program));
+    std::shared_ptr<storm::models::sparse::Dtmc<storm::RationalNumber>> dtmc = storm::api::buildSparseModel<storm::RationalNumber>(program, formulas)->as<storm::models::sparse::Dtmc<storm::RationalNumber>>();
+    uint_fast64_t const initState = *dtmc->getInitialStates().begin();;
+    std::unique_ptr<storm::modelchecker::CheckResult> result;
+    
+    result = storm::api::verifyWithSparseEngine(dtmc, storm::api::createTask<storm::RationalNumber>(formulas[0], true));
+    ASSERT_TRUE(result->isExplicitQuantitativeCheckResult());
+    EXPECT_EQ(storm::utility::convertNumber<storm::RationalNumber>(std::string("24/25")), result->asExplicitQuantitativeCheckResult<storm::RationalNumber>()[initState]);
+
+    result = storm::api::verifyWithSparseEngine(dtmc, storm::api::createTask<storm::RationalNumber>(formulas[1], true));
+    ASSERT_TRUE(result->isExplicitQuantitativeCheckResult());
+    EXPECT_EQ(storm::utility::convertNumber<storm::RationalNumber>(std::string("624/625")), result->asExplicitQuantitativeCheckResult<storm::RationalNumber>()[initState]);
+
+    result = storm::api::verifyWithSparseEngine(dtmc, storm::api::createTask<storm::RationalNumber>(formulas[2], true));
+    ASSERT_TRUE(result->isExplicitQuantitativeCheckResult());
+    EXPECT_EQ(storm::utility::convertNumber<storm::RationalNumber>(std::string("1/625")), result->asExplicitQuantitativeCheckResult<storm::RationalNumber>()[initState]);
+
+    result = storm::api::verifyWithSparseEngine(dtmc, storm::api::createTask<storm::RationalNumber>(formulas[3], true));
+    ASSERT_TRUE(result->isExplicitQuantitativeCheckResult());
+    EXPECT_EQ(storm::utility::convertNumber<storm::RationalNumber>(std::string("24/625")), result->asExplicitQuantitativeCheckResult<storm::RationalNumber>()[initState]);
+}
+
+TEST(SparseDtmcMultiDimensionalRewardUnfoldingTest, cost_bounded_crowds) {
+    storm::Environment env;
+    std::string programFile = STORM_TEST_RESOURCES_DIR "/dtmc/crowds_cost_bounded.pm";
+    std::string formulasAsString = "P=? [F{\"num_runs\"}<=3,{\"observe0\"}>1 true]";
+    formulasAsString += "; P=? [F{\"num_runs\"}<=3,{\"observe1\"}>1 true]";
+    formulasAsString += "; R{\"observe0\"}=? [C{\"num_runs\"}<=3]";
+
+    // programm, model,  formula
+    storm::prism::Program program = storm::api::parseProgram(programFile);
+    program = storm::utility::prism::preprocess(program, "CrowdSize=4");
+    std::vector<std::shared_ptr<storm::logic::Formula const>> formulas = storm::api::extractFormulasFromProperties(storm::api::parsePropertiesForPrismProgram(formulasAsString, program));
+    std::shared_ptr<storm::models::sparse::Dtmc<storm::RationalNumber>> dtmc = storm::api::buildSparseModel<storm::RationalNumber>(program, formulas)->as<storm::models::sparse::Dtmc<storm::RationalNumber>>();
+    uint_fast64_t const initState = *dtmc->getInitialStates().begin();;
+    std::unique_ptr<storm::modelchecker::CheckResult> result;
+    
+    result = storm::api::verifyWithSparseEngine(dtmc, storm::api::createTask<storm::RationalNumber>(formulas[0], true));
+    ASSERT_TRUE(result->isExplicitQuantitativeCheckResult());
+    EXPECT_EQ(storm::utility::convertNumber<storm::RationalNumber>(std::string("78686542099694893/1268858272000000000")), result->asExplicitQuantitativeCheckResult<storm::RationalNumber>()[initState]);
+
+    result = storm::api::verifyWithSparseEngine(dtmc, storm::api::createTask<storm::RationalNumber>(formulas[1], true));
+    ASSERT_TRUE(result->isExplicitQuantitativeCheckResult());
+    EXPECT_EQ(storm::utility::convertNumber<storm::RationalNumber>(std::string("13433618626105041/1268858272000000000")), result->asExplicitQuantitativeCheckResult<storm::RationalNumber>()[initState]);
+
+    result = storm::api::verifyWithSparseEngine(dtmc, storm::api::createTask<storm::RationalNumber>(formulas[2], true));
+    ASSERT_TRUE(result->isExplicitQuantitativeCheckResult());
+    EXPECT_EQ(storm::utility::convertNumber<storm::RationalNumber>(std::string("620529/1364000")), result->asExplicitQuantitativeCheckResult<storm::RationalNumber>()[initState]);
+}