diff --git a/src/counterexamples/MILPMinimalLabelSetGenerator.h b/src/counterexamples/MILPMinimalLabelSetGenerator.h
index 611da0f7e..a3387b174 100644
--- a/src/counterexamples/MILPMinimalLabelSetGenerator.h
+++ b/src/counterexamples/MILPMinimalLabelSetGenerator.h
@@ -32,7 +32,7 @@ namespace storm {
          * property in terms of used labels.
          */
         template <class T>
-        class MinimalLabelSetGenerator {
+        class MILPMinimalLabelSetGenerator {
 #ifdef STORM_HAVE_GUROBI
         private:
             /*!
diff --git a/src/counterexamples/SMTMinimalCommandSetGenerator.h b/src/counterexamples/SMTMinimalCommandSetGenerator.h
index 4608f4e8e..40373b8fe 100644
--- a/src/counterexamples/SMTMinimalCommandSetGenerator.h
+++ b/src/counterexamples/SMTMinimalCommandSetGenerator.h
@@ -11,6 +11,15 @@
 // To detect whether the usage of Z3 is possible, this include is neccessary.
 #include "storm-config.h"
 
+// If we have Z3 available, we have to include the C++ header.
+#ifdef STORM_HAVE_Z3
+#include "z3++.h"
+#endif
+
+#include "src/ir/Program.h"
+#include "src/modelchecker/prctl/SparseMdpPrctlModelChecker.h"
+#include "src/solver/GmmxxNondeterministicLinearEquationSolver.h"
+
 namespace storm {
     namespace counterexamples {
         
@@ -19,30 +28,233 @@ namespace storm {
          * property in terms of used labels.
          */
         template <class T>
-        class MinimalLabelSetGenerator {
+        class SMTMinimalCommandSetGenerator {
 #ifdef STORM_HAVE_Z3
         private:
+            struct VariableInformation {
+                std::vector<z3::expr> labelVariables;
+                std::vector<z3::expr> auxiliaryVariables;
+                std::map<uint_fast64_t, uint_fast64_t> labelToIndexMap;
+            };
+            
+            /*!
+             * Computes the set of relevant labels in the model. Relevant labels are choice labels such that there exists
+             * a scheduler that satisfies phi until psi with a nonzero probability.
+             *
+             * @param labeledMdp The MDP to search for relevant labels.
+             * @param phiStates A bit vector representing all states that satisfy phi.
+             * @param psiStates A bit vector representing all states that satisfy psi.
+             * @return A set of relevant labels, where relevant is defined as above.
+             */
+            static std::set<uint_fast64_t> getRelevantLabels(storm::models::Mdp<T> const& labeledMdp, storm::storage::BitVector const& phiStates, storm::storage::BitVector const& psiStates) {
+                // Create result.
+                std::set<uint_fast64_t> relevantLabels;
+                
+                // Compute all relevant states, i.e. states for which there exists a scheduler that has a non-zero
+                // probabilitiy of satisfying phi until psi.
+                storm::storage::SparseMatrix<bool> backwardTransitions = labeledMdp.getBackwardTransitions();
+                storm::storage::BitVector relevantStates = storm::utility::graph::performProbGreater0E(labeledMdp, backwardTransitions, phiStates, psiStates);
+                relevantStates &= ~psiStates;
+
+                // Retrieve some references for convenient access.
+                storm::storage::SparseMatrix<T> const& transitionMatrix = labeledMdp.getTransitionMatrix();
+                std::vector<uint_fast64_t> const& nondeterministicChoiceIndices = labeledMdp.getNondeterministicChoiceIndices();
+                std::vector<std::set<uint_fast64_t>> const& choiceLabeling = labeledMdp.getChoiceLabeling();
+
+                // Now traverse all choices of all relevant states and check whether there is a successor target state.
+                // If so, the associated labels become relevant. Also, if a choice of relevant state has at least one
+                // relevant successor, the choice becomes relevant.
+                for (auto state : relevantStates) {
+                    for (uint_fast64_t row = nondeterministicChoiceIndices[state]; row < nondeterministicChoiceIndices[state + 1]; ++row) {
+                        for (typename storm::storage::SparseMatrix<T>::ConstIndexIterator successorIt = transitionMatrix.constColumnIteratorBegin(row); successorIt != transitionMatrix.constColumnIteratorEnd(row); ++successorIt) {
+                            // If there is a relevant successor, we need to add the labels of the current choice.
+                            if (relevantStates.get(*successorIt) || psiStates.get(*successorIt)) {
+                                for (auto const& label : choiceLabeling[row]) {
+                                    relevantLabels.insert(label);
+                                }
+                            }
+                        }
+                    }
+                }
+                
+                LOG4CPLUS_DEBUG(logger, "Found " << relevantLabels.size() << " relevant labels.");
+                return relevantLabels;
+            }
             
+            /*!
+             * Creates all necessary base expressions for the relevant labels.
+             *
+             * @param context The Z3 context in which to create the expressions.
+             * @param relevantCommands A set of relevant labels for which to create the expressions.
+             * @return A mapping from relevant labels to their corresponding expressions.
+             */
+            static VariableInformation createExpressionsForRelevantLabels(z3::context& context, std::set<uint_fast64_t> const& relevantLabels) {
+                VariableInformation variableInformation;
+                
+                // Create stringstream to build expression names.
+                std::stringstream variableName;
+                
+                for (auto label : relevantLabels) {
+                    variableInformation.labelToIndexMap[label] = variableInformation.labelVariables.size();
+                    
+                    // Clear contents of the stream to construct new expression name.
+                    variableName.clear();
+                    variableName.str("");
+                    variableName << "c" << label;
+                    
+                    variableInformation.labelVariables.push_back(context.bool_const(variableName.str().c_str()));
+                    
+                    // Clear contents of the stream to construct new expression name.
+                    variableName.clear();
+                    variableName.str("");
+                    variableName << "h" << label;
+                    
+                    variableInformation.auxiliaryVariables.push_back(context.bool_const(variableName.str().c_str()));
+                }
+                
+                return variableInformation;
+            }
+
+            /*!
+             * Asserts the constraints that are initially known.
+             *
+             * @param program The program for which to build the constraints.
+             * @param context The Z3 context in which to build the expressions.
+             * @param solver The solver in which to assert the constraints.
+             * @param variableInformation A structure with information about the variables for the labels.
+             */
+            static void assertInitialConstraints(storm::ir::Program const& program, z3::context& context, z3::solver& solver, VariableInformation const& variableInformation) {
+                // Assert that at least one of the labels must be taken.
+                z3::expr formula = variableInformation.labelVariables.at(0);
+                for (uint_fast64_t index = 1; index < variableInformation.labelVariables.size(); ++index) {
+                    formula = formula || variableInformation.labelVariables.at(index);
+                }
+                solver.add(formula);
+                
+                for (uint_fast64_t index = 0; index < variableInformation.labelVariables.size(); ++index) {
+                    solver.add(!variableInformation.labelVariables[index] || variableInformation.auxiliaryVariables[index]);
+                }
+            }
+
+            /*!
+             * Performs one Fu-Malik-Maxsat step.
+             *
+             * @param context The Z3 context in which to build the expressions.
+             * @param solver The solver to use for the satisfiability evaluation.
+             * @param variableInformation A structure with information about the variables for the labels.
+             * @return True iff the constraint system was satisfiable.
+             */
+            static bool fuMalikMaxsatStep(z3::context& context, z3::solver& solver, VariableInformation const& variableInformation) {
+                z3::expr_vector assumptions(context);
+                for (auto const& auxVariable : variableInformation.auxiliaryVariables) {
+                    assumptions.push_back(!auxVariable);
+                }
+                
+                std::cout << solver << std::endl;
+                
+                // Check whether the assumptions are satisfiable.
+                z3::check_result result = solver.check(assumptions);
+                
+                if (result == z3::check_result::sat) {
+                    return true;
+                } else {
+                    z3::expr_vector unsatCore = solver.unsat_core();
+                    
+                    std::vector<z3::expr> blockingVariables;
+                    blockingVariables.reserve(unsatCore.size());
+                }
+                
+                return false;
+            }
+
+            
+            /*!
+             * Finds the smallest set of labels such that the constraint system of the solver is still satisfiable.
+             *
+             * @param context The Z3 context in which to build the expressions.
+             * @param solver The solver to use for the satisfiability evaluation.
+             * @param variableInformation A structure with information about the variables for the labels.
+             * @return The smallest set of labels such that the constraint system of the solver is still satisfiable.
+             */
+            static std::set<uint_fast64_t> findSmallestCommandSet(z3::context& context, z3::solver& solver, VariableInformation const& variableInformation) {
+                for (uint_fast64_t i = 0; i < variableInformation.labelToIndexMap.size(); ++i) {
+                    if (fuMalikMaxsatStep(context, solver, variableInformation)) {
+                        break;
+                    }
+                }
+                
+                // Now we are ready to construct the label set from the model of the solver.
+                
+                std::set<uint_fast64_t> result;
+                
+                for (auto const& labelIndexPair : variableInformation.labelToIndexMap) {
+                    result.insert(labelIndexPair.first);
+                }
+                
+                return result;
+            }
 #endif
             
         public:
-            static std::unordered_set<uint_fast64_t> getMinimalCommandSet(storm::ir::Program const& program, storm::models::Mdp<T> const& labeledMdp, storm::storage::BitVector const& phiStates, storm::storage::BitVector const& psiStates, double probabilityThreshold, bool checkThresholdFeasible = false) {
+            static std::set<uint_fast64_t> getMinimalCommandSet(storm::ir::Program const& program, storm::models::Mdp<T> const& labeledMdp, storm::storage::BitVector const& phiStates, storm::storage::BitVector const& psiStates, double probabilityThreshold, bool checkThresholdFeasible = false) {
 #ifdef STORM_HAVE_Z3
                 // (0) Check whether the MDP is indeed labeled.
                 if (!labeledMdp.hasChoiceLabels()) {
                     throw storm::exceptions::InvalidArgumentException() << "Minimal command set generation is impossible for unlabeled model.";
                 }
                 
+                // (1) FIXME: check whether its possible to exceed the threshold if checkThresholdFeasible is set.
+
+                // (2) Identify all commands that are relevant, because only these need to be considered later.
+                std::set<uint_fast64_t> relevantCommands = getRelevantLabels(labeledMdp, phiStates, psiStates);
+                
+                // (3) Create context for solver.
+                z3::context context;
                 
+                // (4) Create the variables for the relevant commands.
+                VariableInformation variableInformation = createExpressionsForRelevantLabels(context, relevantCommands);
+                
+                // (5) After all variables have been created, create a solver for that context.
+                z3::solver solver(context);
+
+                // (5) Build the initial constraint system.
+                assertInitialConstraints(program, context, solver, variableInformation);
+                
+                // (6) Find the smallest set of commands that satisfies all constraints. If the probability of
+                // satisfying phi until psi exceeds the given threshold, the set of labels is minimal and can be returned.
+                // Otherwise, the current solution has to be ruled out and the next smallest solution is retrieved from
+                // the solver.
+                double maximalReachabilityProbability = 0;
+                std::set<uint_fast64_t> commandSet;
+                do {
+                    commandSet = findSmallestCommandSet(context, solver, variableInformation);
+                    
+                    storm::models::Mdp<T> subMdp = labeledMdp.restrictChoiceLabels(commandSet);
+                    storm::modelchecker::prctl::SparseMdpPrctlModelChecker<T> modelchecker(subMdp, new storm::solver::GmmxxNondeterministicLinearEquationSolver<T>());
+                    std::vector<T> result = modelchecker.checkUntil(false, phiStates, psiStates, false, nullptr);
+                    
+                    // Now determine the maximalReachabilityProbability.
+                    for (auto state : labeledMdp.getInitialStates()) {
+                        maximalReachabilityProbability = std::max(maximalReachabilityProbability, result[state]);
+                    }
+                } while (maximalReachabilityProbability < probabilityThreshold);
+                
+                std::cout << "Achieved probability: " << maximalReachabilityProbability << " with " << commandSet.size() << " commands." << std::endl;
+                std::cout << "Taken commands are:" << std::endl;
+                for (auto label : commandSet) {
+                    std::cout << label << ", ";
+                }
+                std::cout << std::endl;
                 
-                return std::unordered_set<uint_fast64_t>();
+                // (7) Return the resulting command set.
+                return commandSet;
                 
 #else
                 throw storm::exceptions::NotImplementedException() << "This functionality is unavailable since StoRM has been compiled without support for Z3.";
 #endif
             }
             
-        }
+        };
         
     } // namespace counterexamples
 } // namespace storm
diff --git a/src/modelchecker/prctl/SparseMdpPrctlModelChecker.h b/src/modelchecker/prctl/SparseMdpPrctlModelChecker.h
index f399c0e68..9f9683a06 100644
--- a/src/modelchecker/prctl/SparseMdpPrctlModelChecker.h
+++ b/src/modelchecker/prctl/SparseMdpPrctlModelChecker.h
@@ -76,28 +76,27 @@ namespace storm {
                 }
                 
                 /*!
-                 * Checks the given formula that is a bounded-until formula.
+                 * Computes the probability to satisfy phi until psi within a limited number of steps for each state.
                  *
-                 * @param formula The formula to check.
-                 * @param qualitative A flag indicating whether the formula only needs to be evaluated qualitatively, i.e. if the
+                 * @param phiStates A bit vector indicating which states satisfy phi.
+                 * @param psiStates A bit vector indicating which states satisfy psi.
+                 * @param stepBound The upper bound for the number of steps.
+                 * @param qualitative A flag indicating whether the check only needs to be done qualitatively, i.e. if the
                  * results are only compared against the bounds 0 and 1. If set to true, this will most likely results that are only
                  * qualitatively correct, i.e. do not represent the correct value, but only the correct relation with respect to the
                  * bounds 0 and 1.
-                 * @returns The probabilities for the given formula to hold on every state of the model associated with this model
-                 * checker. If the qualitative flag is set, exact probabilities might not be computed.
+                 * @return The probabilities for satisfying phi until psi within a limited number of steps for each state.
+                 * If the qualitative flag is set, exact probabilities might not be computed.
                  */
-                virtual std::vector<Type> checkBoundedUntil(const storm::property::prctl::BoundedUntil<Type>& formula, bool qualitative) const {
-                    // First, we need to compute the states that satisfy the sub-formulas of the until-formula.
-                    storm::storage::BitVector leftStates = formula.getLeft().check(*this);
-                    storm::storage::BitVector rightStates = formula.getRight().check(*this);
+                std::vector<Type> checkBoundedUntil(storm::storage::BitVector const& phiStates, storm::storage::BitVector const& psiStates, uint_fast64_t stepBound, bool qualitative) const {
                     std::vector<Type> result(this->getModel().getNumberOfStates());
-                    
+
                     // Determine the states that have 0 probability of reaching the target states.
                     storm::storage::BitVector statesWithProbabilityGreater0;
                     if (this->minimumOperatorStack.top()) {
-                        statesWithProbabilityGreater0 = storm::utility::graph::performProbGreater0A(this->getModel(), this->getModel().getBackwardTransitions(), leftStates, rightStates, true, formula.getBound());
+                        statesWithProbabilityGreater0 = storm::utility::graph::performProbGreater0A(this->getModel(), this->getModel().getBackwardTransitions(), phiStates, psiStates, true, stepBound);
                     } else {
-                        statesWithProbabilityGreater0 = storm::utility::graph::performProbGreater0E(this->getModel(), this->getModel().getBackwardTransitions(), leftStates, rightStates, true, formula.getBound());
+                        statesWithProbabilityGreater0 = storm::utility::graph::performProbGreater0E(this->getModel(), this->getModel().getBackwardTransitions(), phiStates, psiStates, true, stepBound);
                     }
                     
                     // Check if we already know the result (i.e. probability 0) for all initial states and
@@ -118,7 +117,7 @@ namespace storm {
                         std::vector<uint_fast64_t> subNondeterministicChoiceIndices = this->computeNondeterministicChoiceIndicesForConstraint(statesWithProbabilityGreater0);
                         
                         // Compute the new set of target states in the reduced system.
-                        storm::storage::BitVector rightStatesInReducedSystem = statesWithProbabilityGreater0 % rightStates;
+                        storm::storage::BitVector rightStatesInReducedSystem = statesWithProbabilityGreater0 % psiStates;
                         
                         // Make all rows absorbing that satisfy the second sub-formula.
                         submatrix.makeRowsAbsorbing(rightStatesInReducedSystem, subNondeterministicChoiceIndices);
@@ -128,7 +127,7 @@ namespace storm {
                         storm::utility::vector::setVectorValues(subresult, rightStatesInReducedSystem, storm::utility::constGetOne<Type>());
                         
                         if (linearEquationSolver != nullptr) {
-                            this->linearEquationSolver->performMatrixVectorMultiplication(this->minimumOperatorStack.top(), submatrix, subresult, subNondeterministicChoiceIndices, nullptr, formula.getBound());
+                            this->linearEquationSolver->performMatrixVectorMultiplication(this->minimumOperatorStack.top(), submatrix, subresult, subNondeterministicChoiceIndices, nullptr, stepBound);
                         } else {
                             throw storm::exceptions::InvalidStateException() << "No valid linear equation solver available.";
                         }
@@ -142,20 +141,32 @@ namespace storm {
                 }
                 
                 /*!
-                 * Checks the given formula that is a next formula.
+                 * Checks the given formula that is a bounded-until formula.
                  *
                  * @param formula The formula to check.
                  * @param qualitative A flag indicating whether the formula only needs to be evaluated qualitatively, i.e. if the
                  * results are only compared against the bounds 0 and 1. If set to true, this will most likely results that are only
                  * qualitatively correct, i.e. do not represent the correct value, but only the correct relation with respect to the
                  * bounds 0 and 1.
-                 * @returns The probabilities for the given formula to hold on every state of the model associated with this model
+                 * @return The probabilities for the given formula to hold on every state of the model associated with this model
                  * checker. If the qualitative flag is set, exact probabilities might not be computed.
                  */
-                virtual std::vector<Type> checkNext(const storm::property::prctl::Next<Type>& formula, bool qualitative) const {
-                    // First, we need to compute the states that satisfy the sub-formula of the next-formula.
-                    storm::storage::BitVector nextStates = formula.getChild().check(*this);
-                    
+                virtual std::vector<Type> checkBoundedUntil(storm::property::prctl::BoundedUntil<Type> const& formula, bool qualitative) const {
+                    return checkBoundedUntil(formula.getLeft().check(*this), formula.getRight().check(*this), formula.getBound(), qualitative);
+                }
+                
+                /*!
+                 * Computes the probability to reach the given set of states in the next step for each state.
+                 *
+                 * @param nextStates A bit vector defining the states to reach in the next state.
+                 * @param qualitative A flag indicating whether the formula only needs to be evaluated qualitatively, i.e. if the
+                 * results are only compared against the bounds 0 and 1. If set to true, this will most likely results that are only
+                 * qualitatively correct, i.e. do not represent the correct value, but only the correct relation with respect to the
+                 * bounds 0 and 1.
+                 * @return The probabilities to reach the gien set of states in the next step for each state. If the
+                 * qualitative flag is set, exact probabilities might not be computed.
+                 */
+                virtual std::vector<Type> checkNext(storm::storage::BitVector const& nextStates, bool qualitative) const {
                     // Create the vector with which to multiply and initialize it correctly.
                     std::vector<Type> result(this->getModel().getNumberOfStates());
                     storm::utility::vector::setVectorValues(result, nextStates, storm::utility::constGetOne<Type>());
@@ -169,6 +180,21 @@ namespace storm {
                     return result;
                 }
                 
+                /*!
+                 * Checks the given formula that is a next formula.
+                 *
+                 * @param formula The formula to check.
+                 * @param qualitative A flag indicating whether the formula only needs to be evaluated qualitatively, i.e. if the
+                 * results are only compared against the bounds 0 and 1. If set to true, this will most likely results that are only
+                 * qualitatively correct, i.e. do not represent the correct value, but only the correct relation with respect to the
+                 * bounds 0 and 1.
+                 * @return The probabilities for the given formula to hold on every state of the model associated with this model
+                 * checker. If the qualitative flag is set, exact probabilities might not be computed.
+                 */
+                virtual std::vector<Type> checkNext(const storm::property::prctl::Next<Type>& formula, bool qualitative) const {
+                    return checkNext(formula.getChild().check(*this), qualitative);
+                }
+                
                 /*!
                  * Checks the given formula that is a bounded-eventually formula.
                  *
@@ -239,14 +265,15 @@ namespace storm {
                  * checker. If the qualitative flag is set, exact probabilities might not be computed.
                  */
                 virtual std::vector<Type> checkUntil(const storm::property::prctl::Until<Type>& formula, bool qualitative) const {
-                    return this->checkUntil(this->minimumOperatorStack.top(), formula, qualitative, nullptr);
+                    return this->checkUntil(this->minimumOperatorStack.top(), formula.getLeft().check(*this), formula.getRight().check(*this), qualitative, nullptr);
                 }
                 
                 /*!
-                 * Check the given formula that is an until formula.
+                 * Computes the extremal probability to satisfy phi until psi for each state in the model.
                  *
                  * @param minimize If set, the probability is minimized and maximized otherwise.
-                 * @param formula The formula to check.
+                 * @param phiStates A bit vector indicating which states satisfy phi.
+                 * @param psiStates A bit vector indicating which states satisfy psi.
                  * @param qualitative A flag indicating whether the formula only needs to be evaluated qualitatively, i.e. if the
                  * results are only compared against the bounds 0 and 1. If set to true, this will most likely results that are only
                  * qualitatively correct, i.e. do not represent the correct value, but only the correct relation with respect to the
@@ -254,21 +281,17 @@ namespace storm {
                  * @param scheduler If <code>qualitative</code> is false and this vector is non-null and has as many elements as
                  * there are states in the MDP, this vector will represent a scheduler for the model that achieves the probability
                  * returned by model checking. To this end, the vector will hold the nondeterministic choice made for each state.
-                 * @return The probabilities for the given formula to hold on every state of the model associated with this model
-                 * checker. If the qualitative flag is set, exact probabilities might not be computed.
+                 * @return The probabilities for the satisfying phi until psi for each state of the model. If the
+                 * qualitative flag is set, exact probabilities might not be computed.
                  */
-                virtual std::vector<Type> checkUntil(bool minimize, const storm::property::prctl::Until<Type>& formula, bool qualitative, std::vector<uint_fast64_t>* scheduler) const {
-                    // First, we need to compute the states that satisfy the sub-formulas of the until-formula.
-                    storm::storage::BitVector leftStates = formula.getLeft().check(*this);
-                    storm::storage::BitVector rightStates = formula.getRight().check(*this);
-                    
-                    // Then, we need to identify the states which have to be taken out of the matrix, i.e.
+                std::vector<Type> checkUntil(bool minimize, storm::storage::BitVector const& phiStates, storm::storage::BitVector const& psiStates, bool qualitative, std::vector<uint_fast64_t>* scheduler) const {
+                    // We need to identify the states which have to be taken out of the matrix, i.e.
                     // all states that have probability 0 and 1 of satisfying the until-formula.
                     std::pair<storm::storage::BitVector, storm::storage::BitVector> statesWithProbability01;
                     if (minimize) {
-                        statesWithProbability01 = storm::utility::graph::performProb01Min(this->getModel(), leftStates, rightStates);
+                        statesWithProbability01 = storm::utility::graph::performProb01Min(this->getModel(), phiStates, psiStates);
                     } else {
-                        statesWithProbability01 = storm::utility::graph::performProb01Max(this->getModel(), leftStates, rightStates);
+                        statesWithProbability01 = storm::utility::graph::performProb01Max(this->getModel(), phiStates, psiStates);
                     }
                     storm::storage::BitVector statesWithProbability0 = std::move(statesWithProbability01.first);
                     storm::storage::BitVector statesWithProbability1 = std::move(statesWithProbability01.second);
@@ -340,7 +363,7 @@ namespace storm {
                  * results are only compared against the bound 0. If set to true, this will most likely results that are only
                  * qualitatively correct, i.e. do not represent the correct value, but only the correct relation with respect to the
                  * bound 0.
-                 * @returns The reward values for the given formula for every state of the model associated with this model
+                 * @return The reward values for the given formula for every state of the model associated with this model
                  * checker. If the qualitative flag is set, exact values might not be computed.
                  */
                 virtual std::vector<Type> checkInstantaneousReward(const storm::property::prctl::InstantaneousReward<Type>& formula, bool qualitative) const {
@@ -370,7 +393,7 @@ namespace storm {
                  * results are only compared against the bound 0. If set to true, this will most likely results that are only
                  * qualitatively correct, i.e. do not represent the correct value, but only the correct relation with respect to the
                  * bound 0.
-                 * @returns The reward values for the given formula for every state of the model associated with this model
+                 * @return The reward values for the given formula for every state of the model associated with this model
                  * checker. If the qualitative flag is set, exact values might not be computed.
                  */
                 virtual std::vector<Type> checkCumulativeReward(const storm::property::prctl::CumulativeReward<Type>& formula, bool qualitative) const {
@@ -420,14 +443,14 @@ namespace storm {
                  * checker. If the qualitative flag is set, exact values might not be computed.
                  */
                 virtual std::vector<Type> checkReachabilityReward(const storm::property::prctl::ReachabilityReward<Type>& formula, bool qualitative) const {
-                    return this->checkReachabilityReward(this->minimumOperatorStack.top(), formula, qualitative, nullptr);
+                    return this->checkReachabilityReward(this->minimumOperatorStack.top(), formula.getChild().check(*this), qualitative, nullptr);
                 }
                 
                 /*!
-                 * Checks the given formula that is a reachability reward formula.
+                 * Computes the expected reachability reward that is gained before a target state is reached for each state.
                  *
                  * @param minimize If set, the reward is to be minimized and maximized otherwise.
-                 * @param formula The formula to check.
+                 * @param targetStates The target states before which rewards can be gained.
                  * @param qualitative A flag indicating whether the formula only needs to be evaluated qualitatively, i.e. if the
                  * results are only compared against the bound 0. If set to true, this will most likely results that are only
                  * qualitatively correct, i.e. do not represent the correct value, but only the correct relation with respect to the
@@ -435,19 +458,16 @@ namespace storm {
                  * @param scheduler If <code>qualitative</code> is false and this vector is non-null and has as many elements as
                  * there are states in the MDP, this vector will represent a scheduler for the model that achieves the probability
                  * returned by model checking. To this end, the vector will hold the nondeterministic choice made for each state.
-                 * @return The reward values for the given formula for every state of the model associated with this model
-                 * checker. If the qualitative flag is set, exact values might not be computed.
+                 * @return The expected reward values gained before a target state is reached for each state. If the
+                 * qualitative flag is set, exact values might not be computed.
                  */
-                virtual std::vector<Type> checkReachabilityReward(bool minimize, const storm::property::prctl::ReachabilityReward<Type>& formula, bool qualitative, std::vector<uint_fast64_t>* scheduler) const {
+                virtual std::vector<Type> checkReachabilityReward(bool minimize, storm::storage::BitVector const& targetStates, bool qualitative, std::vector<uint_fast64_t>* scheduler) const {
                     // Only compute the result if the model has at least one reward model.
-                    if (!this->getModel().hasStateRewards() && !this->getModel().hasTransitionRewards()) {
-                        LOG4CPLUS_ERROR(logger, "Missing reward model for formula. Skipping formula");
+                    if (!(this->getModel().hasStateRewards() || this->getModel().hasTransitionRewards())) {
+                        LOG4CPLUS_ERROR(logger, "Missing reward model for formula.");
                         throw storm::exceptions::InvalidPropertyException() << "Missing reward model for formula.";
                     }
                     
-                    // Determine the states for which the target predicate holds.
-                    storm::storage::BitVector targetStates = formula.getChild().check(*this);
-                    
                     // Determine which states have a reward of infinity by definition.
                     storm::storage::BitVector infinityStates;
                     storm::storage::BitVector trueStates(this->getModel().getNumberOfStates(), true);
diff --git a/src/models/AbstractModel.h b/src/models/AbstractModel.h
index f32cb3ad4..deffd9999 100644
--- a/src/models/AbstractModel.h
+++ b/src/models/AbstractModel.h
@@ -217,7 +217,7 @@ class AbstractModel: public std::enable_shared_from_this<AbstractModel<T>> {
             // First, we need to count how many backward transitions each state has.
             for (uint_fast64_t i = 0; i < numberOfStates; ++i) {
                 typename storm::storage::SparseMatrix<T>::Rows rows = this->getRows(i);
-                for (auto& transition : rows) {
+                for (auto const& transition : rows) {
                     if (transition.value() > 0) {
                         ++rowIndications[transition.column() + 1];
                     }
diff --git a/src/models/Mdp.h b/src/models/Mdp.h
index 0f97e304a..48a2584d5 100644
--- a/src/models/Mdp.h
+++ b/src/models/Mdp.h
@@ -127,16 +127,16 @@ public:
 	}
 
     /*!
-     * Constructs an MDP by copying the given MDP and restricting the choices of each state to the ones whose label set
+     * Constructs an MDP by copying the current MDP and restricting the choices of each state to the ones whose label set
      * is contained in the given label set.
      *
-     * @param originalModel The model to restrict.
      * @param enabledChoiceLabels A set of labels that determines which choices of the original model can be taken
      * and which ones need to be ignored.
+     * @return A restricted version of the current MDP that only uses choice labels from the given set.
      */
-    Mdp<T> restrictChoiceLabels(Mdp<T> const& originalModel, std::set<uint_fast64_t> const& enabledChoiceLabels) {
+    Mdp<T> restrictChoiceLabels(std::set<uint_fast64_t> const& enabledChoiceLabels) const {
         // Only perform this operation if the given model has choice labels.
-        if (!originalModel.hasChoiceLabels()) {
+        if (!this->hasChoiceLabels()) {
             throw storm::exceptions::InvalidArgumentException() << "Restriction to label set is impossible for unlabeled model.";
         }
         
@@ -144,26 +144,40 @@ public:
         
         storm::storage::SparseMatrix<T> transitionMatrix;
         transitionMatrix.initialize();
+        std::vector<uint_fast64_t> nondeterministicChoiceIndices;
         
         // Check for each choice of each state, whether the choice labels are fully contained in the given label set.
+        uint_fast64_t currentRow = 0;
         for(uint_fast64_t state = 0; state < this->getNumberOfStates(); ++state) {
+            bool stateHasValidChoice = false;
             for (uint_fast64_t choice = this->getNondeterministicChoiceIndices()[state]; choice < this->getNondeterministicChoiceIndices()[state + 1]; ++choice) {
                 bool choiceValid = storm::utility::set::isSubsetOf(choiceLabeling[state], enabledChoiceLabels);
                 
                 // If the choice is valid, copy over all its elements.
                 if (choiceValid) {
+                    if (!stateHasValidChoice) {
+                        nondeterministicChoiceIndices.push_back(currentRow);
+                    }
+                    stateHasValidChoice = true;
                     typename storm::storage::SparseMatrix<T>::Rows row = this->getTransitionMatrix().getRows(choice, choice);
                     for (typename storm::storage::SparseMatrix<T>::ConstIterator rowIt = row.begin(), rowIte = row.end(); rowIt != rowIte; ++rowIt) {
-                        transitionMatrix.insertNextValue(choice, rowIt.column(), rowIt.value(), true);
+                        transitionMatrix.insertNextValue(currentRow, rowIt.column(), rowIt.value(), true);
                     }
-                } else {
-                    // If the choice may not be taken, we insert a self-loop to the state instead.
-                    transitionMatrix.insertNextValue(choice, state, storm::utility::constGetOne<T>(), true);
-                }
+                    ++currentRow;
+                } 
+            }
+            
+            // If no choice of the current state may be taken, we insert a self-loop to the state instead.
+            if (!stateHasValidChoice) {
+                nondeterministicChoiceIndices.push_back(currentRow);
+                transitionMatrix.insertNextValue(currentRow, state, storm::utility::constGetOne<T>(), true);
+                ++currentRow;
             }
         }
-        
-        Mdp<T> restrictedMdp(std::move(transitionMatrix), storm::models::AtomicPropositionsLabeling(this->getStateLabeling()), std::vector<uint_fast64_t>(this->getNondeterministicChoiceIndices()), this->hasStateRewards() ? boost::optional<std::vector<T>>(this->getStateRewardVector()) : boost::optional<std::vector<T>>(), this->hasTransitionRewards() ? boost::optional<storm::storage::SparseMatrix<T>>(this->getTransitionRewardMatrix()) : boost::optional<storm::storage::SparseMatrix<T>>(), boost::optional<std::vector<std::set<uint_fast64_t>>>(this->getChoiceLabeling()));
+        transitionMatrix.finalize(true);
+        nondeterministicChoiceIndices.push_back(currentRow);
+                
+        Mdp<T> restrictedMdp(std::move(transitionMatrix), storm::models::AtomicPropositionsLabeling(this->getStateLabeling()), std::move(nondeterministicChoiceIndices), this->hasStateRewards() ? boost::optional<std::vector<T>>(this->getStateRewardVector()) : boost::optional<std::vector<T>>(), this->hasTransitionRewards() ? boost::optional<storm::storage::SparseMatrix<T>>(this->getTransitionRewardMatrix()) : boost::optional<storm::storage::SparseMatrix<T>>(), boost::optional<std::vector<std::set<uint_fast64_t>>>(this->getChoiceLabeling()));
         return restrictedMdp;
     }
     
diff --git a/src/storage/SparseMatrix.h b/src/storage/SparseMatrix.h
index 0623090ac..84fe0bf94 100644
--- a/src/storage/SparseMatrix.h
+++ b/src/storage/SparseMatrix.h
@@ -159,7 +159,7 @@ public:
 		 *
 		 * @return The value of the current non-zero element to which this iterator points.
          */
-		T const& value() {
+		T const& value() const {
 			return *valuePtr;
 		}
         
diff --git a/src/storm.cpp b/src/storm.cpp
index 39497409a..9543930da 100644
--- a/src/storm.cpp
+++ b/src/storm.cpp
@@ -27,6 +27,7 @@
 #include "src/solver/GmmxxLinearEquationSolver.h"
 #include "src/solver/GmmxxNondeterministicLinearEquationSolver.h"
 #include "src/counterexamples/MILPMinimalLabelSetGenerator.h"
+#include "src/counterexamples/SMTMinimalCommandSetGenerator.h"
 #include "src/parser/AutoParser.h"
 #include "src/parser/PrctlParser.h"
 #include "src/utility/ErrorHandling.h"
@@ -332,17 +333,27 @@ int main(const int argc, const char* argv[]) {
 		} else if (s->isSet("symbolic")) {
 			std::string const& programFile = s->getOptionByLongName("symbolic").getArgument(0).getValueAsString();
 			std::string const& constants = s->getOptionByLongName("constants").getArgument(0).getValueAsString();
-			std::shared_ptr<storm::models::AbstractModel<storm::storage::LabeledValues<double>>> model = storm::adapters::ExplicitModelAdapter<storm::storage::LabeledValues<double>>::translateProgram(storm::parser::PrismParserFromFile(programFile), constants);
+            storm::ir::Program program = storm::parser::PrismParserFromFile(programFile);
+			std::shared_ptr<storm::models::AbstractModel<double>> model = storm::adapters::ExplicitModelAdapter<double>::translateProgram(program, constants);
 			model->printModelInformationToStream(std::cout);
 
             // Enable the following lines to test the MinimalLabelSetGenerator.
-//            if (model->getType() == storm::models::MDP) {
-//                std::shared_ptr<storm::models::Mdp<storm::storage::LabeledValues<double>>> labeledMdp = model->as<storm::models::Mdp<storm::storage::LabeledValues<double>>>();
-//                storm::storage::BitVector const& finishedStates = labeledMdp->getLabeledStates("finished");
-//                storm::storage::BitVector const& allCoinsEqual1States = labeledMdp->getLabeledStates("all_coins_equal_1");
-//                storm::storage::BitVector targetStates = finishedStates & allCoinsEqual1States;
-//                storm::counterexamples::MinimalLabelSetGenerator<storm::storage::LabeledValues<double>>::getMinimalLabelSet(*labeledMdp, storm::storage::BitVector(labeledMdp->getNumberOfStates(), true), targetStates, 0.3, true, true);
-//            }
+            if (model->getType() == storm::models::MDP) {
+                std::shared_ptr<storm::models::Mdp<double>> labeledMdp = model->as<storm::models::Mdp<double>>();
+                storm::storage::BitVector const& finishedStates = labeledMdp->getLabeledStates("finished");
+                storm::storage::BitVector const& allCoinsEqual1States = labeledMdp->getLabeledStates("all_coins_equal_1");
+                storm::storage::BitVector targetStates = finishedStates & allCoinsEqual1States;
+                storm::counterexamples::MILPMinimalLabelSetGenerator<double>::getMinimalLabelSet(*labeledMdp, storm::storage::BitVector(labeledMdp->getNumberOfStates(), true), targetStates, 0.3, true, true);
+            }
+            
+            // Enable the following lines to test the SMTMinimalCommandSetGenerator.
+            if (model->getType() == storm::models::MDP) {
+                std::shared_ptr<storm::models::Mdp<double>> labeledMdp = model->as<storm::models::Mdp<double>>();
+                storm::storage::BitVector const& finishedStates = labeledMdp->getLabeledStates("finished");
+                storm::storage::BitVector const& allCoinsEqual1States = labeledMdp->getLabeledStates("all_coins_equal_1");
+                storm::storage::BitVector targetStates = finishedStates & allCoinsEqual1States;
+                storm::counterexamples::SMTMinimalCommandSetGenerator<double>::getMinimalCommandSet(program, *labeledMdp, storm::storage::BitVector(labeledMdp->getNumberOfStates(), true), targetStates, 0.3, true);
+            }
 		}
 
         // Perform clean-up and terminate.
diff --git a/src/utility/set.h b/src/utility/set.h
index d264034eb..7d5d85051 100644
--- a/src/utility/set.h
+++ b/src/utility/set.h
@@ -43,6 +43,7 @@ namespace storm {
                         
                     // Otherwise, we have found an equivalent element and can continue with the next one.
                 }
+                return true;
             }
             
         } // namespace set
diff --git a/storm-config.h.in b/storm-config.h.in
index 8f7a83312..08bf0ca38 100644
--- a/storm-config.h.in
+++ b/storm-config.h.in
@@ -19,8 +19,11 @@
 // Whether Gurobi is available and to be used (define/undef)
 #@STORM_CPP_GUROBI_DEF@ STORM_HAVE_GUROBI
 
+// Whether Z3 is available and to be used (define/undef)
+#@STORM_CPP_Z3_DEF@ STORM_HAVE_Z3
+
 // Whether Intel Threading Building Blocks are available and to be used (define/undef)
 #@STORM_CPP_INTELTBB_DEF@ STORM_HAVE_INTELTBB
 
 
-#endif // STORM_GENERATED_STORMCONFIG_H_
\ No newline at end of file
+#endif // STORM_GENERATED_STORMCONFIG_H_