diff --git a/src/modelchecker/csl/SparseMarkovAutomatonCslModelChecker.cpp b/src/modelchecker/csl/SparseMarkovAutomatonCslModelChecker.cpp
index 10236bdcd..ab3e8f32f 100644
--- a/src/modelchecker/csl/SparseMarkovAutomatonCslModelChecker.cpp
+++ b/src/modelchecker/csl/SparseMarkovAutomatonCslModelChecker.cpp
@@ -54,7 +54,8 @@ namespace storm {
             } else {
                 upperBound = pathFormula.getDiscreteTimeBound();
             }
-            std::vector<ValueType> result = storm::modelchecker::helper::SparseMarkovAutomatonCslHelper<ValueType>::computeBoundedUntilProbabilities(checkTask.getOptimizationDirection(), this->getModel().getTransitionMatrix(), this->getModel().getExitRates(), this->getModel().getMarkovianStates(), rightResult.getTruthValuesVector(), lowerBound, upperBound, *minMaxLinearEquationSolverFactory);
+            
+            std::vector<ValueType> result = storm::modelchecker::helper::SparseMarkovAutomatonCslHelper<ValueType>::computeBoundedUntilProbabilities(checkTask.getOptimizationDirection(), this->getModel().getTransitionMatrix(), this->getModel().getExitRates(), this->getModel().getMarkovianStates(), rightResult.getTruthValuesVector(), std::make_pair(lowerBound, upperBound), *minMaxLinearEquationSolverFactory);
             return std::unique_ptr<CheckResult>(new ExplicitQuantitativeCheckResult<ValueType>(std::move(result)));
         }
                 
diff --git a/src/modelchecker/csl/helper/SparseMarkovAutomatonCslHelper.cpp b/src/modelchecker/csl/helper/SparseMarkovAutomatonCslHelper.cpp
index 0139cdfcd..0c55b6556 100644
--- a/src/modelchecker/csl/helper/SparseMarkovAutomatonCslHelper.cpp
+++ b/src/modelchecker/csl/helper/SparseMarkovAutomatonCslHelper.cpp
@@ -32,6 +32,7 @@ namespace storm {
 
             template<typename ValueType>
             void SparseMarkovAutomatonCslHelper<ValueType>::computeBoundedReachabilityProbabilities(OptimizationDirection dir, storm::storage::SparseMatrix<ValueType> const& transitionMatrix, std::vector<ValueType> const& exitRates, storm::storage::BitVector const& markovianStates, storm::storage::BitVector const& goalStates, storm::storage::BitVector const& markovianNonGoalStates, storm::storage::BitVector const& probabilisticNonGoalStates, std::vector<ValueType>& markovianNonGoalValues, std::vector<ValueType>& probabilisticNonGoalValues, ValueType delta, uint_fast64_t numberOfSteps, storm::utility::solver::MinMaxLinearEquationSolverFactory<ValueType> const& minMaxLinearEquationSolverFactory) {
+                
                 // Start by computing four sparse matrices:
                 // * a matrix aMarkovian with all (discretized) transitions from Markovian non-goal states to all Markovian non-goal states.
                 // * a matrix aMarkovianToProbabilistic with all (discretized) transitions from Markovian non-goal states to all probabilistic non-goal states.
@@ -116,12 +117,16 @@ namespace storm {
                 storm::utility::vector::addVectors(bProbabilistic, bProbabilisticFixed, bProbabilistic);
                 solver->solveEquationSystem(dir, probabilisticNonGoalValues, bProbabilistic, &multiplicationResultScratchMemory, &aProbabilisticScratchMemory);
             }
-            
 
             template<typename ValueType>
-            std::vector<ValueType> SparseMarkovAutomatonCslHelper<ValueType>::computeBoundedUntilProbabilities(OptimizationDirection dir, storm::storage::SparseMatrix<ValueType> const& transitionMatrix, std::vector<ValueType> const& exitRateVector, storm::storage::BitVector const& markovianStates, storm::storage::BitVector const& psiStates, double lowerBound, double upperBound, storm::utility::solver::MinMaxLinearEquationSolverFactory<ValueType> const& minMaxLinearEquationSolverFactory) {
+            std::vector<ValueType> SparseMarkovAutomatonCslHelper<ValueType>::computeBoundedUntilProbabilities(OptimizationDirection dir, storm::storage::SparseMatrix<ValueType> const& transitionMatrix, std::vector<ValueType> const& exitRateVector, storm::storage::BitVector const& markovianStates, storm::storage::BitVector const& psiStates, std::pair<double, double> const& boundsPair, storm::utility::solver::MinMaxLinearEquationSolverFactory<ValueType> const& minMaxLinearEquationSolverFactory) {
+
                 uint_fast64_t numberOfStates = transitionMatrix.getRowGroupCount();
-                
+               
+                // 'Unpack' the bounds to make them more easily accessible.
+                double lowerBound = boundsPair.first;
+                double upperBound = boundsPair.second;
+
                 // (1) Compute the accuracy we need to achieve the required error bound.
                 ValueType maxExitRate = 0;
                 for (auto value : exitRateVector) {
@@ -505,4 +510,4 @@ namespace storm {
             
         }
     }
-}
\ No newline at end of file
+}
diff --git a/src/modelchecker/csl/helper/SparseMarkovAutomatonCslHelper.h b/src/modelchecker/csl/helper/SparseMarkovAutomatonCslHelper.h
index 96429d727..534f1d272 100644
--- a/src/modelchecker/csl/helper/SparseMarkovAutomatonCslHelper.h
+++ b/src/modelchecker/csl/helper/SparseMarkovAutomatonCslHelper.h
@@ -14,7 +14,7 @@ namespace storm {
             template <typename ValueType>
             class SparseMarkovAutomatonCslHelper {
             public:
-                static std::vector<ValueType> computeBoundedUntilProbabilities(OptimizationDirection dir, storm::storage::SparseMatrix<ValueType> const& transitionMatrix, std::vector<ValueType> const& exitRateVector, storm::storage::BitVector const& markovianStates, storm::storage::BitVector const& psiStates, double lowerBound, double upperBound, storm::utility::solver::MinMaxLinearEquationSolverFactory<ValueType> const& minMaxLinearEquationSolverFactory);
+                static std::vector<ValueType> computeBoundedUntilProbabilities(OptimizationDirection dir, storm::storage::SparseMatrix<ValueType> const& transitionMatrix, std::vector<ValueType> const& exitRateVector, storm::storage::BitVector const& markovianStates, storm::storage::BitVector const& psiStates, std::pair<double, double> const& boundsPair, storm::utility::solver::MinMaxLinearEquationSolverFactory<ValueType> const& minMaxLinearEquationSolverFactory);
                 
                 static std::vector<ValueType> computeUntilProbabilities(OptimizationDirection dir, 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::utility::solver::MinMaxLinearEquationSolverFactory<ValueType> const& minMaxLinearEquationSolverFactory);
                 
diff --git a/src/modelchecker/prctl/helper/SparseDtmcPrctlHelper.cpp b/src/modelchecker/prctl/helper/SparseDtmcPrctlHelper.cpp
index 9e1b6ef14..98cbcf2ea 100644
--- a/src/modelchecker/prctl/helper/SparseDtmcPrctlHelper.cpp
+++ b/src/modelchecker/prctl/helper/SparseDtmcPrctlHelper.cpp
@@ -367,6 +367,7 @@ namespace storm {
                         // Now compute reachability probabilities in the transformed model.
                         storm::storage::SparseMatrix<ValueType> const& newTransitionMatrix = transformedModel.transitionMatrix.get();
                         std::vector<ValueType> conditionalProbabilities = computeUntilProbabilities(newTransitionMatrix, newTransitionMatrix.transpose(), storm::storage::BitVector(newTransitionMatrix.getRowCount(), true), transformedModel.targetStates.get(), qualitative, linearEquationSolverFactory);
+                        
                         storm::utility::vector::setVectorValues(result, transformedModel.beforeStates, conditionalProbabilities);
                     }
                 }
diff --git a/src/modelchecker/prctl/helper/SparseDtmcPrctlHelper.h b/src/modelchecker/prctl/helper/SparseDtmcPrctlHelper.h
index 637c30cb5..0be762298 100644
--- a/src/modelchecker/prctl/helper/SparseDtmcPrctlHelper.h
+++ b/src/modelchecker/prctl/helper/SparseDtmcPrctlHelper.h
@@ -45,6 +45,10 @@ namespace storm {
                 static std::vector<ValueType> computeReachabilityRewards(storm::storage::SparseMatrix<ValueType> const& transitionMatrix, storm::storage::SparseMatrix<ValueType> const& backwardTransitions, std::function<std::vector<ValueType>(uint_fast64_t, storm::storage::SparseMatrix<ValueType> const&, storm::storage::BitVector const&)> const& totalStateRewardVectorGetter, storm::storage::BitVector const& targetStates, bool qualitative, storm::utility::solver::LinearEquationSolverFactory<ValueType> const& linearEquationSolverFactory);
                 
                 struct BaierTransformedModel {
+                    BaierTransformedModel() : noTargetStates(false) {
+                        // Intentionally left empty.
+                    }
+                    
                     storm::storage::BitVector beforeStates;
                     boost::optional<storm::storage::SparseMatrix<ValueType>> transitionMatrix;
                     boost::optional<storm::storage::BitVector> targetStates;
diff --git a/src/models/sparse/MarkovAutomaton.cpp b/src/models/sparse/MarkovAutomaton.cpp
index da55901d3..b3b416cfc 100644
--- a/src/models/sparse/MarkovAutomaton.cpp
+++ b/src/models/sparse/MarkovAutomaton.cpp
@@ -19,9 +19,8 @@ namespace storm {
                                                         std::vector<ValueType> const& exitRates,
                                                         std::unordered_map<std::string, RewardModelType> const& rewardModels,
                                                         boost::optional<std::vector<LabelSet>> const& optionalChoiceLabeling)
-            : NondeterministicModel<ValueType, RewardModelType>(storm::models::ModelType::MarkovAutomaton, transitionMatrix, stateLabeling, rewardModels, optionalChoiceLabeling), markovianStates(markovianStates), exitRates(exitRates) {
+            : NondeterministicModel<ValueType, RewardModelType>(storm::models::ModelType::MarkovAutomaton, transitionMatrix, stateLabeling, rewardModels, optionalChoiceLabeling), markovianStates(markovianStates), exitRates(exitRates), closed(this->checkIsClosed()) {
                 this->turnRatesToProbabilities();
-                closed = !hasHybridState();
             }
             
             template <typename ValueType, typename RewardModelType>
@@ -31,9 +30,8 @@ namespace storm {
                                                         std::vector<ValueType> const& exitRates,
                                                         std::unordered_map<std::string, RewardModelType>&& rewardModels,
                                                         boost::optional<std::vector<LabelSet>>&& optionalChoiceLabeling)
-            : NondeterministicModel<ValueType, RewardModelType>(storm::models::ModelType::MarkovAutomaton, std::move(transitionMatrix), std::move(stateLabeling), std::move(rewardModels), std::move(optionalChoiceLabeling)), markovianStates(markovianStates), exitRates(std::move(exitRates)) {
+            : NondeterministicModel<ValueType, RewardModelType>(storm::models::ModelType::MarkovAutomaton, std::move(transitionMatrix), std::move(stateLabeling), std::move(rewardModels), std::move(optionalChoiceLabeling)), markovianStates(markovianStates), exitRates(std::move(exitRates)), closed(this->checkIsClosed()) {
                 this->turnRatesToProbabilities();
-                closed = !hasHybridState();
             }
             
             template <typename ValueType, typename RewardModelType>
@@ -44,10 +42,9 @@ namespace storm {
                                                                          bool probabilities,
                                                                          std::unordered_map<std::string, RewardModelType>&& rewardModels,
                                                                          boost::optional<std::vector<LabelSet>>&& optionalChoiceLabeling)
-            : NondeterministicModel<ValueType, RewardModelType>(storm::models::ModelType::MarkovAutomaton, std::move(transitionMatrix), std::move(stateLabeling), std::move(rewardModels), std::move(optionalChoiceLabeling)), markovianStates(markovianStates), exitRates(std::move(exitRates)) {
+            : NondeterministicModel<ValueType, RewardModelType>(storm::models::ModelType::MarkovAutomaton, std::move(transitionMatrix), std::move(stateLabeling), std::move(rewardModels), std::move(optionalChoiceLabeling)), markovianStates(markovianStates), exitRates(std::move(exitRates)), closed(this->checkIsClosed()) {
                 assert(probabilities);
                 assert(this->getTransitionMatrix().isProbabilistic());
-                closed = !hasHybridState();
             }
             
             template <typename ValueType, typename RewardModelType>
@@ -70,16 +67,6 @@ namespace storm {
                 return !this->markovianStates.get(state);
             }
             
-            template <typename ValueType, typename RewardModelType>
-            bool MarkovAutomaton<ValueType, RewardModelType>::hasHybridState() const {
-                for (uint_fast64_t state = 0; state < this->getNumberOfStates(); ++state) {
-                    if (this->isHybridState(state)) {
-                        return true;
-                    }
-                }
-                return false;
-            }
-            
             template <typename ValueType, typename RewardModelType>
             std::vector<ValueType> const& MarkovAutomaton<ValueType, RewardModelType>::getExitRates() const {
                 return this->exitRates;
@@ -127,21 +114,19 @@ namespace storm {
                     // Now copy over all choices that need to be kept.
                     uint_fast64_t currentChoice = 0;
                     for (uint_fast64_t state = 0; state < this->getNumberOfStates(); ++state) {
-                        // If the state is a hybrid state, closing it will make it a probabilistic state, so we remove the Markovian marking.
-                        if (this->isHybridState(state)) {
-                            this->markovianStates.set(state, false);
-                        }
-                        
                         // Record the new beginning of choices of this state.
                         newTransitionMatrixBuilder.newRowGroup(currentChoice);
-                        
-                        // If we are currently treating a hybrid state, we need to skip its first choice.
+
+                        // If the state is a hybrid state, closing it will make it a probabilistic state, so we remove the Markovian marking.
+                        // Additionally, we need to remember whether we need to skip the first choice of the state when
+                        // we assemble the new transition matrix.
+                        uint_fast64_t offset = 0;
                         if (this->isHybridState(state)) {
-                            // Remove the Markovian state marking.
                             this->markovianStates.set(state, false);
+                            offset = 1;
                         }
                         
-                        for (uint_fast64_t row = this->getTransitionMatrix().getRowGroupIndices()[state] + (this->isHybridState(state) ? 1 : 0); row < this->getTransitionMatrix().getRowGroupIndices()[state + 1]; ++row) {
+                        for (uint_fast64_t row = this->getTransitionMatrix().getRowGroupIndices()[state] + offset; row < this->getTransitionMatrix().getRowGroupIndices()[state + 1]; ++row) {
                             for (auto const& entry : this->getTransitionMatrix().getRow(row)) {
                                 newTransitionMatrixBuilder.addNextValue(currentChoice, entry.getColumn(), entry.getValue());
                             }
@@ -246,7 +231,7 @@ namespace storm {
             template <typename ValueType, typename RewardModelType>
             void MarkovAutomaton<ValueType, RewardModelType>::turnRatesToProbabilities() {
                 for (auto state : this->markovianStates) {
-                    for (auto& transition : this->getTransitionMatrix().getRowGroup(state)) {
+                    for (auto& transition : this->getTransitionMatrix().getRow(this->getTransitionMatrix().getRowGroupIndices()[state])) {
                         transition.setValue(transition.getValue() / this->exitRates[state]);
                     }
                 }
@@ -269,6 +254,16 @@ namespace storm {
                 return true;
             }
             
+            template <typename ValueType, typename RewardModelType>
+            bool MarkovAutomaton<ValueType, RewardModelType>::checkIsClosed() const {
+                for (auto state : markovianStates) {
+                    if (this->getTransitionMatrix().getRowGroupSize(state) > 1) {
+                        return false;
+                    }
+                }
+                return true;
+            }
+            
             template <typename ValueType, typename RewardModelType>
             std::shared_ptr<storm::models::sparse::Ctmc<ValueType, RewardModelType>> MarkovAutomaton<ValueType, RewardModelType>::convertToCTMC() {
                 STORM_LOG_TRACE("MA matrix:" << std::endl << this->getTransitionMatrix());
@@ -336,4 +331,4 @@ namespace storm {
 
         } // namespace sparse
     } // namespace models
-} // namespace storm
\ No newline at end of file
+} // namespace storm
diff --git a/src/models/sparse/MarkovAutomaton.h b/src/models/sparse/MarkovAutomaton.h
index 4251640fc..9864f4f85 100644
--- a/src/models/sparse/MarkovAutomaton.h
+++ b/src/models/sparse/MarkovAutomaton.h
@@ -165,11 +165,9 @@ namespace storm {
                 void turnRatesToProbabilities();
                 
                 /*!
-                 * Check if at least one hybrid state exists.
-                 *
-                 * @return True, if at least one hybrid state exists, false if none exists.
+                 * Checks whether the automaton is closed by actually looking at the transition information.
                  */
-                bool hasHybridState() const;
+                bool checkIsClosed() const;
                 
                 // A bit vector representing the set of Markovian states.
                 storm::storage::BitVector markovianStates;
diff --git a/src/parser/MarkovAutomatonSparseTransitionParser.cpp b/src/parser/MarkovAutomatonSparseTransitionParser.cpp
index 5a9fcfb84..43f8aed64 100644
--- a/src/parser/MarkovAutomatonSparseTransitionParser.cpp
+++ b/src/parser/MarkovAutomatonSparseTransitionParser.cpp
@@ -6,6 +6,7 @@
 #include "src/exceptions/FileIoException.h"
 #include "src/parser/MappedFile.h"
 #include "src/utility/cstring.h"
+#include "src/utility/constants.h"
 #include "src/utility/macros.h"
 
 namespace storm {
@@ -152,6 +153,16 @@ namespace storm {
                     }
                 } while (!encounteredEOF && !encounteredNewDistribution);
             }
+            
+            // If there are some states with indices that are behind the last source for which no transition was specified,
+            // we need to reserve some space for introducing self-loops later.
+            if (!dontFixDeadlocks) {
+                result.numberOfNonzeroEntries += result.highestStateIndex - lastsource;
+                result.numberOfChoices += result.highestStateIndex - lastsource;
+            } else {
+                STORM_LOG_ERROR("Found deadlock states (e.g. " << lastsource + 1 << ") during parsing. Please fix them or set the appropriate flag.");
+                throw storm::exceptions::WrongFormatException() << "Found deadlock states (e.g. " << lastsource + 1 << ") during parsing. Please fix them or set the appropriate flag.";
+            }
 
             return result;
         }
@@ -259,6 +270,21 @@ namespace storm {
 
                 ++currentChoice;
             }
+            
+            // If there are some states with indices that are behind the last source for which no transition was specified,
+            // we need to insert the self-loops now. Note that we assume all these states to be Markovian.
+            if (!dontFixDeadlocks) {
+                for (uint_fast64_t index = lastsource + 1; index <= firstPassResult.highestStateIndex; ++index) {
+                    result.markovianStates.set(index, true);
+                    result.exitRates[index] = storm::utility::one<ValueType>();
+                    result.transitionMatrixBuilder.newRowGroup(currentChoice);
+                    result.transitionMatrixBuilder.addNextValue(currentChoice, index, storm::utility::one<ValueType>());
+                    ++currentChoice;
+                }
+            } else {
+                STORM_LOG_ERROR("Found deadlock states (e.g. " << lastsource + 1 << ") during parsing. Please fix them or set the appropriate flag.");
+                throw storm::exceptions::WrongFormatException() << "Found deadlock states (e.g. " << lastsource + 1 << ") during parsing. Please fix them or set the appropriate flag.";
+            }
 
             return result;
         }
diff --git a/src/storage/bisimulation/BisimulationDecomposition.cpp b/src/storage/bisimulation/BisimulationDecomposition.cpp
index fb6a6c759..1da648981 100644
--- a/src/storage/bisimulation/BisimulationDecomposition.cpp
+++ b/src/storage/bisimulation/BisimulationDecomposition.cpp
@@ -240,6 +240,9 @@ namespace storm {
                 ++iterations;
 
                 // Get and prepare the next splitter.
+                // Sort the splitters according to their sizes to prefer small splitters. That is just a heuristic, but
+                // tends to work well.
+                std::sort(splitterQueue.begin(), splitterQueue.end(), [] (Block<BlockDataType> const* b1, Block<BlockDataType> const* b2) { return b1->getNumberOfStates() < b2->getNumberOfStates(); } );
                 Block<BlockDataType>* splitter = splitterQueue.front();
                 splitterQueue.pop_front();
                 splitter->data().setSplitter(false);
diff --git a/src/utility/graph.cpp b/src/utility/graph.cpp
index de920ad21..bd10cfd03 100644
--- a/src/utility/graph.cpp
+++ b/src/utility/graph.cpp
@@ -310,8 +310,46 @@ namespace storm {
             }
             
             template <typename T>
-            storm::storage::PartialScheduler computeSchedulerProbGreater0E(storm::storage::BitVector const& probGreater0EStates, storm::storage::SparseMatrix<T> const& transitionMatrix) {
-                return computeSchedulerWithOneSuccessorInStates(probGreater0EStates, transitionMatrix);
+            storm::storage::PartialScheduler computeSchedulerProbGreater0E(storm::storage::SparseMatrix<T> const& transitionMatrix, storm::storage::SparseMatrix<T> const& backwardTransitions, storm::storage::BitVector const& phiStates, storm::storage::BitVector const& psiStates, boost::optional<storm::storage::BitVector> const& rowFilter) {
+                //Perform backwards DFS from psiStates and find a valid choice for each visited state.
+                
+                storm::storage::PartialScheduler result;
+                std::vector<uint_fast64_t> stack;
+                storm::storage::BitVector currentStates(psiStates); //the states that are either psiStates or for which we have found a valid choice.
+                stack.insert(stack.end(), currentStates.begin(), currentStates.end());
+                uint_fast64_t currentState = 0;
+                
+                while (!stack.empty()) {
+                    currentState = stack.back();
+                    stack.pop_back();
+                    
+                    for (typename storm::storage::SparseMatrix<T>::const_iterator predecessorEntryIt = backwardTransitions.begin(currentState), predecessorEntryIte = backwardTransitions.end(currentState); predecessorEntryIt != predecessorEntryIte; ++predecessorEntryIt) {
+                        uint_fast64_t const& predecessor = predecessorEntryIt->getColumn();
+                        if (phiStates.get(predecessor) && !currentStates.get(predecessor)) {
+                            //The predecessor is a probGreater0E state that has not been considered yet. Let's find the right choice that leads to a state in currentStates.
+                            bool foundValidChoice = false;
+                            for (uint_fast64_t row = transitionMatrix.getRowGroupIndices()[predecessor]; row < transitionMatrix.getRowGroupIndices()[predecessor + 1]; ++row) {
+                                if(rowFilter && !rowFilter->get(row)){
+                                    continue;
+                                }
+                                for (typename storm::storage::SparseMatrix<T>::const_iterator successorEntryIt = transitionMatrix.begin(row), successorEntryIte = transitionMatrix.end(row); successorEntryIt != successorEntryIte; ++successorEntryIt) {
+                                    if(currentStates.get(successorEntryIt->getColumn())){
+                                        foundValidChoice = true;
+                                        break;
+                                    }
+                                }
+                                if(foundValidChoice){
+                                    result.setChoice(predecessor, row - transitionMatrix.getRowGroupIndices()[predecessor]);
+                                    currentStates.set(predecessor, true);
+                                    stack.push_back(predecessor);
+                                    break;
+                                }
+                            }
+                            STORM_LOG_INFO_COND(foundValidChoice, "Could not find a valid choice for ProbGreater0E state " << predecessor << ".");
+                        }
+                    }
+                }
+                return result;
             }
             
             template <typename T>
@@ -972,7 +1010,7 @@ namespace storm {
             
             
             
-            template storm::storage::PartialScheduler computeSchedulerProbGreater0E(storm::storage::BitVector const& probGreater0EStates, storm::storage::SparseMatrix<double> const& transitionMatrix);
+            template storm::storage::PartialScheduler computeSchedulerProbGreater0E(storm::storage::SparseMatrix<double> const& transitionMatrix, storm::storage::SparseMatrix<double> const& backwardTransitions, storm::storage::BitVector const& phiStates, storm::storage::BitVector const& psiStates, boost::optional<storm::storage::BitVector> const& rowFilter);
             
             template storm::storage::PartialScheduler computeSchedulerProb0E(storm::storage::BitVector const& prob0EStates, storm::storage::SparseMatrix<double> const& transitionMatrix);
             
diff --git a/src/utility/graph.h b/src/utility/graph.h
index 913a8aaab..54e541ac9 100644
--- a/src/utility/graph.h
+++ b/src/utility/graph.h
@@ -234,13 +234,20 @@ namespace storm {
             storm::storage::PartialScheduler computeSchedulerWithOneSuccessorInStates(storm::storage::BitVector const& states, storm::storage::SparseMatrix<T> const& transitionMatrix);
             
             /*!
-             * Computes a scheduler for the given states that have a scheduler that has a probability greater 0.
+             * Computes a scheduler for the ProbGreater0E-States such that in the induced system the given psiStates are reachable via phiStates
              *
-             * @param probGreater0EStates The states that have a scheduler achieving a probablity greater 0.
              * @param transitionMatrix The transition matrix of the system.
+             * @param backwardTransitions The reversed transition relation.
+             * @param phiStates The set of states satisfying phi.
+             * @param psiStates The set of states satisfying psi.
+             * @param rowFilter If given, the returned scheduler will only pick choices such that rowFilter is true for the corresponding matrixrow.
+             * @return A Scheduler for the ProbGreater0E-States
+             *
+             * @note No choice is defined for ProbGreater0E-States if all the probGreater0-choices violate the row filter.
+             *       This also holds for states that only reach psi via such states.
              */
             template <typename T>
-            storm::storage::PartialScheduler computeSchedulerProbGreater0E(storm::storage::BitVector const& probGreater0EStates, storm::storage::SparseMatrix<T> const& transitionMatrix);
+            storm::storage::PartialScheduler computeSchedulerProbGreater0E(storm::storage::SparseMatrix<T> const& transitionMatrix, storm::storage::SparseMatrix<T> const& backwardTransitions, storm::storage::BitVector const& phiStates, storm::storage::BitVector const& psiStates, boost::optional<storm::storage::BitVector> const& rowFilter = boost::none);
 
             /*!
              * Computes a scheduler for the given states that have a scheduler that has a probability 0.
@@ -252,10 +259,14 @@ namespace storm {
             storm::storage::PartialScheduler computeSchedulerProb0E(storm::storage::BitVector const& prob0EStates, storm::storage::SparseMatrix<T> const& transitionMatrix);
 
             /*!
-             * Computes a scheduler for the given states that have a scheduler that has a probability 0.
+             * Computes a scheduler for the given prob1EStates such that in the induced system the given psiStates are reached with probability 1.
              *
              * @param prob1EStates The states that have a scheduler achieving probablity 1.
              * @param transitionMatrix The transition matrix of the system.
+             * @param backwardTransitions The reversed transition relation.
+             * @param phiStates The set of states satisfying phi.
+             * @param psiStates The set of states satisfying psi.
+             * @return A scheduler for the Prob1E-States
              */
             template <typename T>
             storm::storage::PartialScheduler computeSchedulerProb1E(storm::storage::BitVector const& prob1EStates, storm::storage::SparseMatrix<T> const& transitionMatrix, storm::storage::SparseMatrix<T> const& backwardTransitions, storm::storage::BitVector const& phiStates, storm::storage::BitVector const& psiStates);
diff --git a/src/utility/storm.h b/src/utility/storm.h
index df84e79d6..297afcd56 100644
--- a/src/utility/storm.h
+++ b/src/utility/storm.h
@@ -64,6 +64,7 @@
 #include "src/modelchecker/csl/helper/SparseCtmcCslHelper.h"
 #include "src/modelchecker/csl/SparseMarkovAutomatonCslModelChecker.h"
 #include "src/modelchecker/csl/HybridCtmcCslModelChecker.h"
+#include "src/modelchecker/csl/SparseMarkovAutomatonCslModelChecker.h"
 #include "src/modelchecker/results/ExplicitQualitativeCheckResult.h"
 #include "src/modelchecker/results/SymbolicQualitativeCheckResult.h"
 
@@ -312,6 +313,10 @@ namespace storm {
             result = modelchecker.check(task);
         } else if (model->getType() == storm::models::ModelType::MarkovAutomaton) {
             std::shared_ptr<storm::models::sparse::MarkovAutomaton<ValueType>> ma = model->template as<storm::models::sparse::MarkovAutomaton<ValueType>>();
+            // Close the MA, if it is not already closed.
+            if (!ma->isClosed()) {
+                ma->close();
+            }
             storm::modelchecker::SparseMarkovAutomatonCslModelChecker<storm::models::sparse::MarkovAutomaton<ValueType>> modelchecker(*ma);
             result = modelchecker.check(task);
         } else {
diff --git a/test/functional/modelchecker/GmmxxMdpPrctlModelCheckerTest.cpp b/test/functional/modelchecker/GmmxxMdpPrctlModelCheckerTest.cpp
index 319f03066..6220276ba 100644
--- a/test/functional/modelchecker/GmmxxMdpPrctlModelCheckerTest.cpp
+++ b/test/functional/modelchecker/GmmxxMdpPrctlModelCheckerTest.cpp
@@ -240,3 +240,36 @@ TEST(GmmxxMdpPrctlModelCheckerTest, SchedulerGeneration) {
     EXPECT_EQ(0, scheduler2.getChoice(2));
     EXPECT_EQ(0, scheduler2.getChoice(3));
 }
+
+TEST(GmmxxMdpPrctlModelCheckerTest, tiny_rewards) {
+    storm::prism::Program program = storm::parser::PrismParser::parse(STORM_CPP_TESTS_BASE_PATH "/functional/modelchecker/tiny_rewards.nm");
+
+    // A parser that we use for conveniently constructing the formulas.
+    storm::parser::FormulaParser formulaParser;
+    
+    std::shared_ptr<storm::models::sparse::Model<double>> model = storm::builder::ExplicitPrismModelBuilder<double>(program).translate();
+    EXPECT_EQ(3ul, model->getNumberOfStates());
+    EXPECT_EQ(4ul, model->getNumberOfTransitions());
+    
+    ASSERT_EQ(model->getType(), storm::models::ModelType::Mdp);
+
+    std::shared_ptr<storm::models::sparse::Mdp<double>> mdp = model->as<storm::models::sparse::Mdp<double>>();
+
+    EXPECT_EQ(4ul, mdp->getNumberOfChoices());
+
+    auto solverFactory = std::make_unique<storm::utility::solver::MinMaxLinearEquationSolverFactory<double>>(storm::solver::EquationSolverTypeSelection::Gmmxx);
+    solverFactory->setPreferredTechnique(storm::solver::MinMaxTechniqueSelection::ValueIteration);
+    storm::modelchecker::SparseMdpPrctlModelChecker<storm::models::sparse::Mdp<double>> checker(*mdp, std::move(solverFactory));
+    
+    std::shared_ptr<const storm::logic::Formula> formula = formulaParser.parseSingleFormulaFromString("Rmin=? [F \"target\"]");
+    
+    storm::modelchecker::CheckTask<storm::logic::Formula> checkTask(*formula);
+    
+    std::unique_ptr<storm::modelchecker::CheckResult> result = checker.check(checkTask);
+    
+    ASSERT_TRUE(result->isExplicitQuantitativeCheckResult());
+    EXPECT_NEAR(1,result->asExplicitQuantitativeCheckResult<double>().getValueVector()[0], storm::settings::nativeEquationSolverSettings().getPrecision());
+    EXPECT_NEAR(1,result->asExplicitQuantitativeCheckResult<double>().getValueVector()[1], storm::settings::nativeEquationSolverSettings().getPrecision());
+    EXPECT_NEAR(0,result->asExplicitQuantitativeCheckResult<double>().getValueVector()[2], storm::settings::nativeEquationSolverSettings().getPrecision());
+    
+}
diff --git a/test/functional/modelchecker/tiny_rewards.nm b/test/functional/modelchecker/tiny_rewards.nm
new file mode 100644
index 000000000..5119f6306
--- /dev/null
+++ b/test/functional/modelchecker/tiny_rewards.nm
@@ -0,0 +1,17 @@
+mdp
+
+module mod1
+
+s : [0..2] init 0;
+[] s=0 -> true;
+[] s=0 -> (s'=1);
+[] s=1 -> (s'=2);
+[] s=2 -> (s'=2);
+
+endmodule
+
+rewards
+ [] s=1 : 1;
+endrewards
+
+label "target" = s=2;