diff --git a/src/modelchecker/csl/SparseMarkovAutomatonCslModelChecker.h b/src/modelchecker/csl/SparseMarkovAutomatonCslModelChecker.h
index dfdb520fd..68491f85b 100644
--- a/src/modelchecker/csl/SparseMarkovAutomatonCslModelChecker.h
+++ b/src/modelchecker/csl/SparseMarkovAutomatonCslModelChecker.h
@@ -8,7 +8,7 @@
 #include "src/models/MarkovAutomaton.h"
 #include "src/storage/BitVector.h"
 #include "src/storage/MaximalEndComponentDecomposition.h"
-#include "src/solver/AbstractNondeterministicLinearEquationSolver.h"
+#include "src/solver/NondeterministicLinearEquationSolver.h"
 #include "src/solver/LpSolver.h"
 #include "src/utility/solver.h"
 #include "src/utility/graph.h"
@@ -21,7 +21,7 @@ namespace storm {
             template<typename ValueType>
             class SparseMarkovAutomatonCslModelChecker : public AbstractModelChecker<ValueType> {
             public:
-                explicit SparseMarkovAutomatonCslModelChecker(storm::models::MarkovAutomaton<ValueType> const& model, std::shared_ptr<storm::solver::AbstractNondeterministicLinearEquationSolver<ValueType>> nondeterministicLinearEquationSolver = storm::utility::solver::getNondeterministicLinearEquationSolver<ValueType>()) : AbstractModelChecker<ValueType>(model), minimumOperatorStack(), nondeterministicLinearEquationSolver(nondeterministicLinearEquationSolver) {
+                explicit SparseMarkovAutomatonCslModelChecker(storm::models::MarkovAutomaton<ValueType> const& model, std::shared_ptr<storm::solver::NondeterministicLinearEquationSolver<ValueType>> nondeterministicLinearEquationSolver = storm::utility::solver::getNondeterministicLinearEquationSolver<ValueType>()) : AbstractModelChecker<ValueType>(model), minimumOperatorStack(), nondeterministicLinearEquationSolver(nondeterministicLinearEquationSolver) {
                     // Intentionally left empty.
                 }
                 
@@ -132,7 +132,7 @@ namespace storm {
                         }
                     }
                     
-                    std::shared_ptr<storm::solver::AbstractNondeterministicLinearEquationSolver<ValueType>> nondeterministiclinearEquationSolver = storm::utility::solver::getNondeterministicLinearEquationSolver<ValueType>();
+                    std::shared_ptr<storm::solver::NondeterministicLinearEquationSolver<ValueType>> nondeterministiclinearEquationSolver = storm::utility::solver::getNondeterministicLinearEquationSolver<ValueType>();
                     
                     // Perform the actual value iteration
                     // * loop until the step bound has been reached
@@ -565,7 +565,7 @@ namespace storm {
                     storm::utility::vector::selectVectorValuesRepeatedly(b, maybeStates, this->getModel().getNondeterministicChoiceIndices(), rewardValues);
                     
                     // Solve the corresponding system of equations.
-                    std::shared_ptr<storm::solver::AbstractNondeterministicLinearEquationSolver<ValueType>> nondeterministiclinearEquationSolver = storm::utility::solver::getNondeterministicLinearEquationSolver<ValueType>();
+                    std::shared_ptr<storm::solver::NondeterministicLinearEquationSolver<ValueType>> nondeterministiclinearEquationSolver = storm::utility::solver::getNondeterministicLinearEquationSolver<ValueType>();
                     nondeterministiclinearEquationSolver->solveEquationSystem(min, submatrix, x, b, subNondeterministicChoiceIndices);
                     
                     // Create resulting vector.
@@ -588,7 +588,7 @@ namespace storm {
                 /*!
                  * A solver that is used for solving systems of linear equations that are the result of nondeterministic choices.
                  */
-                std::shared_ptr<storm::solver::AbstractNondeterministicLinearEquationSolver<ValueType>> nondeterministicLinearEquationSolver;
+                std::shared_ptr<storm::solver::NondeterministicLinearEquationSolver<ValueType>> nondeterministicLinearEquationSolver;
             };
         }
     }
diff --git a/src/modelchecker/prctl/SparseDtmcPrctlModelChecker.h b/src/modelchecker/prctl/SparseDtmcPrctlModelChecker.h
index 0ac15408b..4208501d0 100644
--- a/src/modelchecker/prctl/SparseDtmcPrctlModelChecker.h
+++ b/src/modelchecker/prctl/SparseDtmcPrctlModelChecker.h
@@ -10,7 +10,7 @@
 
 #include "src/modelchecker/prctl/AbstractModelChecker.h"
 #include "src/models/Dtmc.h"
-#include "src/solver/AbstractLinearEquationSolver.h"
+#include "src/solver/LinearEquationSolver.h"
 #include "src/utility/vector.h"
 #include "src/utility/graph.h"
 
@@ -33,7 +33,7 @@ public:
 	 *
 	 * @param model The DTMC to be checked.
 	 */
-	explicit SparseDtmcPrctlModelChecker(storm::models::Dtmc<Type> const& model, storm::solver::AbstractLinearEquationSolver<Type>* linearEquationSolver) : AbstractModelChecker<Type>(model), linearEquationSolver(linearEquationSolver) {
+	explicit SparseDtmcPrctlModelChecker(storm::models::Dtmc<Type> const& model, storm::solver::LinearEquationSolver<Type>* linearEquationSolver) : AbstractModelChecker<Type>(model), linearEquationSolver(linearEquationSolver) {
 		// Intentionally left empty.
 	}
     
@@ -499,7 +499,7 @@ public:
 
 private:
     // An object that is used for solving linear equations and performing matrix-vector multiplication.
-    std::unique_ptr<storm::solver::AbstractLinearEquationSolver<Type>> linearEquationSolver;
+    std::unique_ptr<storm::solver::LinearEquationSolver<Type>> linearEquationSolver;
 };
 
 } // namespace prctl
diff --git a/src/modelchecker/prctl/SparseMdpPrctlModelChecker.h b/src/modelchecker/prctl/SparseMdpPrctlModelChecker.h
index 1b21804b8..f7fe02731 100644
--- a/src/modelchecker/prctl/SparseMdpPrctlModelChecker.h
+++ b/src/modelchecker/prctl/SparseMdpPrctlModelChecker.h
@@ -13,8 +13,8 @@
 #include <fstream>
 
 #include "src/modelchecker/prctl/AbstractModelChecker.h"
-#include "src/solver/AbstractNondeterministicLinearEquationSolver.h"
-#include "src/solver/AbstractLinearEquationSolver.h"
+#include "src/solver/NondeterministicLinearEquationSolver.h"
+#include "src/solver/LinearEquationSolver.h"
 #include "src/models/Mdp.h"
 #include "src/utility/vector.h"
 #include "src/utility/graph.h"
@@ -42,7 +42,7 @@ namespace storm {
                     // Intentionally left empty.
                 }
                 
-				explicit SparseMdpPrctlModelChecker(storm::models::Mdp<Type> const& model, std::shared_ptr<storm::solver::AbstractNondeterministicLinearEquationSolver<Type>> nondeterministicLinearEquationSolver) : AbstractModelChecker<Type>(model), minimumOperatorStack(), nondeterministicLinearEquationSolver(nondeterministicLinearEquationSolver) {
+				explicit SparseMdpPrctlModelChecker(storm::models::Mdp<Type> const& model, std::shared_ptr<storm::solver::NondeterministicLinearEquationSolver<Type>> nondeterministicLinearEquationSolver) : AbstractModelChecker<Type>(model), minimumOperatorStack(), nondeterministicLinearEquationSolver(nondeterministicLinearEquationSolver) {
 					// Intentionally left empty.
 				}
 
@@ -282,7 +282,7 @@ namespace storm {
                  * @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.
                  */
-                static std::pair<std::vector<Type>, storm::storage::TotalScheduler> computeUnboundedUntilProbabilities(bool minimize, storm::storage::SparseMatrix<Type> const& transitionMatrix, std::vector<uint_fast64_t> nondeterministicChoiceIndices, storm::storage::SparseMatrix<Type> const& backwardTransitions, storm::storage::BitVector const& initialStates, storm::storage::BitVector const& phiStates, storm::storage::BitVector const& psiStates, std::shared_ptr<storm::solver::AbstractNondeterministicLinearEquationSolver<Type>> nondeterministicLinearEquationSolver, bool qualitative) {
+                static std::pair<std::vector<Type>, storm::storage::TotalScheduler> computeUnboundedUntilProbabilities(bool minimize, storm::storage::SparseMatrix<Type> const& transitionMatrix, std::vector<uint_fast64_t> nondeterministicChoiceIndices, storm::storage::SparseMatrix<Type> const& backwardTransitions, storm::storage::BitVector const& initialStates, storm::storage::BitVector const& phiStates, storm::storage::BitVector const& psiStates, std::shared_ptr<storm::solver::NondeterministicLinearEquationSolver<Type>> nondeterministicLinearEquationSolver, bool qualitative) {
                     size_t numberOfStates = phiStates.size();
                     
                     // We need to identify the states which have to be taken out of the matrix, i.e.
@@ -590,7 +590,7 @@ namespace storm {
                 /*!
                  * A solver that is used for solving systems of linear equations that are the result of nondeterministic choices.
                  */
-                std::shared_ptr<storm::solver::AbstractNondeterministicLinearEquationSolver<Type>> nondeterministicLinearEquationSolver;
+                std::shared_ptr<storm::solver::NondeterministicLinearEquationSolver<Type>> nondeterministicLinearEquationSolver;
             };
         } // namespace prctl
     } // namespace modelchecker
diff --git a/src/solver/GlpkLpSolver.cpp b/src/solver/GlpkLpSolver.cpp
index 0d51fb590..11aeeda6e 100644
--- a/src/solver/GlpkLpSolver.cpp
+++ b/src/solver/GlpkLpSolver.cpp
@@ -12,15 +12,23 @@
 
 extern log4cplus::Logger logger;
 
+bool GlpkLpSolverOptionsRegistered = storm::settings::Settings::registerNewModule([] (storm::settings::Settings* instance) -> bool {
+	instance->addOption(storm::settings::OptionBuilder("GlpkLpSolver", "glpkoutput", "", "If set, the glpk output will be printed to the command line.").build());
+    
+    return true;
+});
+
 namespace storm {
     namespace solver {
-        GlpkLpSolver::GlpkLpSolver(std::string const& name, ModelSense const& modelSense) : LpSolver(modelSense), lp(nullptr), nextVariableIndex(1), nextConstraintIndex(1), isOptimal(false), rowIndices(), columnIndices(), coefficientValues() {
+        GlpkLpSolver::GlpkLpSolver(std::string const& name, ModelSense const& modelSense) : LpSolver(modelSense), lp(nullptr), nextVariableIndex(1), nextConstraintIndex(1), rowIndices(), columnIndices(), coefficientValues() {
             // Create the LP problem for glpk.
             lp = glp_create_prob();
             
             // Set its name and model sense.
             glp_set_prob_name(lp, name.c_str());
-            this->setModelSense(modelSense);
+            
+            // Set whether the glpk output shall be printed to the command line.
+            glp_term_out(storm::settings::Settings::getInstance()->isSet("debug") || storm::settings::Settings::getInstance()->isSet("glpkoutput") ? GLP_ON : GLP_OFF);
             
             // Because glpk uses 1-based indexing (wtf!?), we need to put dummy elements into the matrix vectors.
             rowIndices.push_back(0);
@@ -33,6 +41,7 @@ namespace storm {
         }
         
         GlpkLpSolver::~GlpkLpSolver() {
+            // Dispose of all objects allocated dynamically by glpk.
             glp_delete_prob(this->lp);
             glp_free_env();
         }
@@ -58,8 +67,7 @@ namespace storm {
             glp_set_obj_coef(this->lp, nextVariableIndex, objectiveFunctionCoefficient);
             ++nextVariableIndex;
             
-            this->isOptimal = false;
-            
+            this->currentModelHasBeenOptimized = false;
             return nextVariableIndex - 1;
         }
         
@@ -112,15 +120,13 @@ namespace storm {
             coefficientValues.insert(coefficientValues.end(), coefficients.begin(), coefficients.end());
             
             ++nextConstraintIndex;
-            this->isOptimal = false;
-        }
-        
-        void GlpkLpSolver::setModelSense(ModelSense const& newModelSense) {
-            glp_set_obj_dir(this->lp, newModelSense == MINIMIZE ? GLP_MIN : GLP_MAX);
-            this->isOptimal = false;
+            this->currentModelHasBeenOptimized = false;
         }
         
         void GlpkLpSolver::optimize() const {
+            // Start by setting the model sense.
+            glp_set_obj_dir(this->lp, this->getModelSense() == MINIMIZE ? GLP_MIN : GLP_MAX);
+            
             glp_load_matrix(this->lp, rowIndices.size() - 1, rowIndices.data(), columnIndices.data(), coefficientValues.data());
             int error = glp_simplex(this->lp, nullptr);
             
@@ -129,7 +135,25 @@ namespace storm {
                 throw storm::exceptions::InvalidStateException() << "Unable to optimize glpk model (" << error << ").";
             }
             
-            this->isOptimal = true;
+            this->currentModelHasBeenOptimized = true;
+        }
+        
+        bool GlpkLpSolver::isInfeasible() const {
+            int status = glp_get_status(this->lp);
+            return status == GLP_INFEAS;
+        }
+        
+        bool GlpkLpSolver::isUnbounded() const {
+            int status = glp_get_status(this->lp);
+            return status == GLP_UNBND;
+        }
+        
+        bool GlpkLpSolver::isOptimal() const {
+            if (!this->currentModelHasBeenOptimized) {
+                return false;
+            }
+            int status = glp_get_status(this->lp);
+            return status == GLP_OPT;
         }
         
         int_fast64_t GlpkLpSolver::getIntegerValue(uint_fast64_t variableIndex) const {
diff --git a/src/solver/GlpkLpSolver.h b/src/solver/GlpkLpSolver.h
index 3741c70e3..9607e59d3 100644
--- a/src/solver/GlpkLpSolver.h
+++ b/src/solver/GlpkLpSolver.h
@@ -14,10 +14,32 @@
 namespace storm {
     namespace solver {
 #ifdef STORM_HAVE_GLPK
+        
+        /*!
+         * A class that implements the LpSolver interface using glpk as the background solver.
+         */
         class GlpkLpSolver : public LpSolver {
         public:
+            /*!
+             * Constructs a solver with the given name and model sense.
+             *
+             * @param name The name of the LP problem.
+             * @param modelSense A value indicating whether the value of the objective function is to be minimized or
+             * maximized.
+             */
             GlpkLpSolver(std::string const& name, ModelSense const& modelSense);
+            
+            /*!
+             * Constructs a solver with the given name. By default the objective function is assumed to be minimized,
+             * but this may be altered later using a call to setModelSense.
+             *
+             * @param name The name of the LP problem.
+             */
             GlpkLpSolver(std::string const& name);
+            
+            /*!
+             * Destructs a solver by freeing the pointers to glpk's structures.
+             */
             virtual ~GlpkLpSolver();
             
             virtual uint_fast64_t createContinuousVariable(std::string const& name, VariableType const& variableType, double lowerBound, double upperBound, double objectiveFunctionCoefficient) override;
@@ -26,9 +48,11 @@ namespace storm {
             
             virtual void addConstraint(std::string const& name, std::vector<uint_fast64_t> const& variables, std::vector<double> const& coefficients, BoundType const& boundType, double rightHandSideValue) override;
             
-            virtual void setModelSense(ModelSense const& newModelSense) override;
             virtual void optimize() const override;
-            
+            virtual bool isInfeasible() const override;
+            virtual bool isUnbounded() const override;
+            virtual bool isOptimal() const override;
+
             virtual int_fast64_t getIntegerValue(uint_fast64_t variableIndex) const override;
             virtual bool getBinaryValue(uint_fast64_t variableIndex) const override;
             virtual double getContinuousValue(uint_fast64_t variableIndex) const override;
@@ -45,9 +69,6 @@ namespace storm {
             // A counter that keeps track of the next free constraint index.
             uint_fast64_t nextConstraintIndex;
             
-            // A flag that stores whether the model was optimized properly.
-            mutable bool isOptimal;
-            
             // The arrays that store the coefficient matrix of the problem.
             std::vector<int> rowIndices;
             std::vector<int> columnIndices;
@@ -57,6 +78,9 @@ namespace storm {
         // If glpk is not available, we provide a stub implementation that emits an error if any of its methods is called.
         class GlpkLpSolver : public LpSolver {
         public:
+            GlpkLpSolver(std::string const& name, ModelSense const& modelSense) {
+                throw storm::exceptions::NotImplementedException() << "This version of StoRM was compiled without support for glpk. Yet, a method was called that requires this support. Please choose a version of support with glpk support.";
+            }
             
             GlpkLpSolver(std::string const& name) : LpSolver(MINIMIZE) {
                 throw storm::exceptions::NotImplementedException() << "This version of StoRM was compiled without support for glpk. Yet, a method was called that requires this support. Please choose a version of support with glpk support.";
@@ -81,12 +105,20 @@ namespace storm {
             virtual void addConstraint(std::string const& name, std::vector<uint_fast64_t> const& variables, std::vector<double> const& coefficients, BoundType const& boundType, double rightHandSideValue) override {
                 throw storm::exceptions::NotImplementedException() << "This version of StoRM was compiled without support for glpk. Yet, a method was called that requires this support. Please choose a version of support with glpk support.";
             }
+                        
+            virtual void optimize() const override {
+                throw storm::exceptions::NotImplementedException() << "This version of StoRM was compiled without support for glpk. Yet, a method was called that requires this support. Please choose a version of support with glpk support.";
+            }
+
+            virtual bool isInfeasible() const override {
+                throw storm::exceptions::NotImplementedException() << "This version of StoRM was compiled without support for glpk. Yet, a method was called that requires this support. Please choose a version of support with glpk support.";
+            }
             
-            virtual void setModelSense(ModelSense const& newModelSense) {
+            virtual bool isUnbounded() const override {
                 throw storm::exceptions::NotImplementedException() << "This version of StoRM was compiled without support for glpk. Yet, a method was called that requires this support. Please choose a version of support with glpk support.";
             }
             
-            virtual void optimize() const override {
+            virtual bool isOptimal() const override {
                 throw storm::exceptions::NotImplementedException() << "This version of StoRM was compiled without support for glpk. Yet, a method was called that requires this support. Please choose a version of support with glpk support.";
             }
             
diff --git a/src/solver/GmmxxLinearEquationSolver.cpp b/src/solver/GmmxxLinearEquationSolver.cpp
index 402a33d59..dd6845659 100644
--- a/src/solver/GmmxxLinearEquationSolver.cpp
+++ b/src/solver/GmmxxLinearEquationSolver.cpp
@@ -85,7 +85,7 @@ namespace storm {
         }
         
         template<typename ValueType>
-        AbstractLinearEquationSolver<ValueType>* GmmxxLinearEquationSolver<ValueType>::clone() const {
+        LinearEquationSolver<ValueType>* GmmxxLinearEquationSolver<ValueType>::clone() const {
             return new GmmxxLinearEquationSolver<ValueType>(*this);
         }
         
diff --git a/src/solver/GmmxxLinearEquationSolver.h b/src/solver/GmmxxLinearEquationSolver.h
index aaab4b31f..a67de0764 100644
--- a/src/solver/GmmxxLinearEquationSolver.h
+++ b/src/solver/GmmxxLinearEquationSolver.h
@@ -1,16 +1,16 @@
 #ifndef STORM_SOLVER_GMMXXLINEAREQUATIONSOLVER_H_
 #define STORM_SOLVER_GMMXXLINEAREQUATIONSOLVER_H_
 
-#include "AbstractLinearEquationSolver.h"
+#include "LinearEquationSolver.h"
 
 namespace storm {
     namespace solver {
 
         /*!
-         * A class that uses the gmm++ library to implement the AbstractLinearEquationSolver interface.
+         * A class that uses the gmm++ library to implement the LinearEquationSolver interface.
          */
         template<typename ValueType>
-        class GmmxxLinearEquationSolver : public AbstractLinearEquationSolver<ValueType> {
+        class GmmxxLinearEquationSolver : public LinearEquationSolver<ValueType> {
         public:
             // An enumeration specifying the available preconditioners.
             enum Preconditioner {
@@ -41,7 +41,7 @@ namespace storm {
              */
             GmmxxLinearEquationSolver();
             
-            virtual AbstractLinearEquationSolver<ValueType>* clone() const override;
+            virtual LinearEquationSolver<ValueType>* clone() const override;
             
             virtual void solveEquationSystem(storm::storage::SparseMatrix<ValueType> const& A, std::vector<ValueType>& x, std::vector<ValueType> const& b, std::vector<ValueType>* multiplyResult = nullptr) const override;
             
diff --git a/src/solver/GmmxxNondeterministicLinearEquationSolver.cpp b/src/solver/GmmxxNondeterministicLinearEquationSolver.cpp
index f7c13526c..f66c520fe 100644
--- a/src/solver/GmmxxNondeterministicLinearEquationSolver.cpp
+++ b/src/solver/GmmxxNondeterministicLinearEquationSolver.cpp
@@ -34,7 +34,7 @@ namespace storm {
 
         
         template<typename ValueType>
-        AbstractNondeterministicLinearEquationSolver<ValueType>* GmmxxNondeterministicLinearEquationSolver<ValueType>::clone() const {
+        NondeterministicLinearEquationSolver<ValueType>* GmmxxNondeterministicLinearEquationSolver<ValueType>::clone() const {
             return new GmmxxNondeterministicLinearEquationSolver<ValueType>(*this);
         }
         
diff --git a/src/solver/GmmxxNondeterministicLinearEquationSolver.h b/src/solver/GmmxxNondeterministicLinearEquationSolver.h
index 54a48efb8..ac5e4d3fc 100644
--- a/src/solver/GmmxxNondeterministicLinearEquationSolver.h
+++ b/src/solver/GmmxxNondeterministicLinearEquationSolver.h
@@ -1,13 +1,16 @@
 #ifndef STORM_SOLVER_GMMXXNONDETERMINISTICLINEAREQUATIONSOLVER_H_
 #define STORM_SOLVER_GMMXXNONDETERMINISTICLINEAREQUATIONSOLVER_H_
 
-#include "src/solver/AbstractNondeterministicLinearEquationSolver.h"
+#include "src/solver/NondeterministicLinearEquationSolver.h"
 
 namespace storm {
     namespace solver {
         
+        /*!
+         * A class that uses the gmm++ library to implement the NondeterministicLinearEquationSolver interface.
+         */
         template<class Type>
-        class GmmxxNondeterministicLinearEquationSolver : public AbstractNondeterministicLinearEquationSolver<Type> {
+        class GmmxxNondeterministicLinearEquationSolver : public NondeterministicLinearEquationSolver<Type> {
         public:
             /*!
              * Constructs a nondeterministic linear equation solver with parameters being set according to the settings
@@ -25,7 +28,7 @@ namespace storm {
              */
             GmmxxNondeterministicLinearEquationSolver(double precision, uint_fast64_t maximalNumberOfIterations, bool relative = true);
 
-            virtual AbstractNondeterministicLinearEquationSolver<Type>* clone() const;
+            virtual NondeterministicLinearEquationSolver<Type>* clone() const;
             
             virtual void performMatrixVectorMultiplication(bool minimize, storm::storage::SparseMatrix<Type> const& A, std::vector<Type>& x, std::vector<uint_fast64_t> const& nondeterministicChoiceIndices, std::vector<Type>* b = nullptr, uint_fast64_t n = 1, std::vector<Type>* multiplyResult = nullptr) const override;
             
diff --git a/src/solver/GurobiLpSolver.cpp b/src/solver/GurobiLpSolver.cpp
index 5b143fe9c..e062e9935 100644
--- a/src/solver/GurobiLpSolver.cpp
+++ b/src/solver/GurobiLpSolver.cpp
@@ -10,13 +10,23 @@
 
 extern log4cplus::Logger logger;
 
+bool GurobiLpSolverOptionsRegistered = storm::settings::Settings::registerNewModule([] (storm::settings::Settings* instance) -> bool {
+    instance->addOption(storm::settings::OptionBuilder("GurobiLpSolver", "gurobithreads", "", "The number of threads that may be used by Gurobi.").addArgument(storm::settings::ArgumentBuilder::createUnsignedIntegerArgument("count", "The number of threads.").setDefaultValueUnsignedInteger(1).build()).build());
+    
+	instance->addOption(storm::settings::OptionBuilder("GurobiLpSolver", "gurobioutput", "", "If set, the Gurobi output will be printed to the command line.").build());
+    
+	instance->addOption(storm::settings::OptionBuilder("GurobiLpSolver", "gurobiinttol", "", "Sets Gurobi's precision for integer variables.").addArgument(storm::settings::ArgumentBuilder::createDoubleArgument("value", "The precision to achieve.").setDefaultValueDouble(1e-06).addValidationFunctionDouble(storm::settings::ArgumentValidators::doubleRangeValidatorExcluding(0.0, 1.0)).build()).build());
+    
+    return true;
+});
+
 namespace storm {
     namespace solver {
         
-        GurobiLpSolver::GurobiLpSolver(std::string const& name, ModelSense const& modelSense) : LpSolver(modelSense), env(nullptr), model(nullptr), nextVariableIndex(0), isOptimal(false) {
+        GurobiLpSolver::GurobiLpSolver(std::string const& name, ModelSense const& modelSense) : LpSolver(modelSense), env(nullptr), model(nullptr), nextVariableIndex(0) {
             // Create the environment.
             int error = GRBloadenv(&env, "");
-            if (error || env == NULL) {
+            if (error || env == nullptr) {
                 LOG4CPLUS_ERROR(logger, "Could not initialize Gurobi (" << GRBgeterrormsg(env) << ").");
                 throw storm::exceptions::InvalidStateException() << "Could not initialize Gurobi environment (" << GRBgeterrormsg(env) << ").";
             }
@@ -30,27 +40,28 @@ namespace storm {
                 LOG4CPLUS_ERROR(logger, "Could not initialize Gurobi model (" << GRBgeterrormsg(env) << ").");
                 throw storm::exceptions::InvalidStateException() << "Could not initialize Gurobi model (" << GRBgeterrormsg(env) << ").";
             }
-            
-            // Set the required sense of the model.
-            this->setModelSense(modelSense);
         }
         
         GurobiLpSolver::GurobiLpSolver(std::string const& name) : GurobiLpSolver(name, MINIMIZE) {
             // Intentionally left empty.
         }
         
+        GurobiLpSolver::GurobiLpSolver() : GurobiLpSolver("", MINIMIZE) {
+            // Intentionally left empty.
+        }
+        
         GurobiLpSolver::~GurobiLpSolver() {
+            // Dispose of the objects allocated inside Gurobi.
             GRBfreemodel(model);
             GRBfreeenv(env);
         }
         
         void GurobiLpSolver::setGurobiEnvironmentProperties() const {
             // Enable the following line to only print the output of Gurobi if the debug flag is set.
-            // int error = error = GRBsetintparam(env, "OutputFlag", storm::settings::Settings::getInstance()->isSet("debug") ? 1 : 0);
-            int error = error = GRBsetintparam(env, "OutputFlag", 1);
+            int error = error = GRBsetintparam(env, "OutputFlag", storm::settings::Settings::getInstance()->isSet("debug") || storm::settings::Settings::getInstance()->isSet("gurobioutput") ? 1 : 0);
             
             // Enable the following line to restrict Gurobi to one thread only.
-            error = GRBsetintparam(env, "Threads", 1);
+            error = GRBsetintparam(env, "Threads", storm::settings::Settings::getInstance()->getOptionByLongName("gurobithreads").getArgument(0).getValueAsUnsignedInteger());
             
             // Enable the following line to force Gurobi to be as precise about the binary variables as required by the given precision option.
             error = GRBsetdblparam(env, "IntFeasTol", storm::settings::Settings::getInstance()->getOptionByLongName("precision").getArgument(0).getValueAsDouble());
@@ -64,7 +75,7 @@ namespace storm {
             }
             
             // Since the model changed, we erase the optimality flag.
-            this->isOptimal = false;
+            this->currentModelHasBeenOptimized = false;
         }
         
         uint_fast64_t GurobiLpSolver::createContinuousVariable(std::string const& name, VariableType const& variableType, double lowerBound, double upperBound, double objectiveFunctionCoefficient) {
@@ -165,43 +176,78 @@ namespace storm {
             this->updateModel();
         }
         
-        void GurobiLpSolver::setModelSense(ModelSense const& newModelSense) {
-            LpSolver::setModelSense(newModelSense);
-            int error = GRBsetintattr(model, "ModelSense", newModelSense == MINIMIZE ? 1 : -1);
+        void GurobiLpSolver::optimize() const {
+            // First incorporate all recent changes.
+            this->updateModel();
+         
+            // Set the most recently set model sense.
+            int error = GRBsetintattr(model, "ModelSense", this->getModelSense() == MINIMIZE ? 1 : -1);
             if (error) {
                 LOG4CPLUS_ERROR(logger, "Unable to set Gurobi model sense (" << GRBgeterrormsg(env) << ").");
                 throw storm::exceptions::InvalidStateException() << "Unable to set Gurobi model sense (" << GRBgeterrormsg(env) << ").";
             }
-        }
-        
-        void GurobiLpSolver::optimize() const {
-            // First incorporate all recent changes.
-            this->updateModel();
             
             // Then we actually optimize the model.
-            int error = GRBoptimize(model);
+            error = GRBoptimize(model);
             if (error) {
                 LOG4CPLUS_ERROR(logger, "Unable to optimize Gurobi model (" << GRBgeterrormsg(env) << ").");
                 throw storm::exceptions::InvalidStateException() << "Unable to optimize Gurobi model (" << GRBgeterrormsg(env) << ").";
             }
-
-            // Finally, we check whether the solution is legal (i.e. bounded, so we can extract values).
+            
+            this->currentModelHasBeenOptimized = true;
+        }
+        
+        bool GurobiLpSolver::isInfeasible() const {
             int optimalityStatus = 0;
-            error = GRBgetintattr(model, GRB_INT_ATTR_STATUS, &optimalityStatus);
-            if (optimalityStatus == GRB_INF_OR_UNBD) {
-                LOG4CPLUS_ERROR(logger, "Unable to get Gurobi solution from infeasible or unbounded model (" << GRBgeterrormsg(env) << ").");
-                throw storm::exceptions::InvalidStateException() << "Unable to get Gurobi solution from infeasible or unbounded MILP (" << GRBgeterrormsg(env) << ").";
-            } else if (optimalityStatus != GRB_OPTIMAL) {
-                LOG4CPLUS_ERROR(logger, "Unable to get Gurobi solution from unoptimized model (" << GRBgeterrormsg(env) << ").");
-                throw storm::exceptions::InvalidStateException() << "Unable to get Gurobi solution from unoptimized model (" << GRBgeterrormsg(env) << ").";
+            
+            int error = GRBgetintattr(model, GRB_INT_ATTR_STATUS, &optimalityStatus);
+            if (error) {
+                LOG4CPLUS_ERROR(logger, "Unable to retrieve optimization status of Gurobi model (" << GRBgeterrormsg(env) << ").");
+                throw storm::exceptions::InvalidStateException() << "Unable to retrieve optimization status of Gurobi model (" << GRBgeterrormsg(env) << ").";
             }
             
-            this->isOptimal = true;
+            return optimalityStatus == GRB_INFEASIBLE;
+        }
+        
+        bool GurobiLpSolver::isUnbounded() const {
+            int optimalityStatus = 0;
+            
+            int error = GRBgetintattr(model, GRB_INT_ATTR_STATUS, &optimalityStatus);
+            if (error) {
+                LOG4CPLUS_ERROR(logger, "Unable to retrieve optimization status of Gurobi model (" << GRBgeterrormsg(env) << ").");
+                throw storm::exceptions::InvalidStateException() << "Unable to retrieve optimization status of Gurobi model (" << GRBgeterrormsg(env) << ").";
+            }
+            
+            return optimalityStatus == GRB_UNBOUNDED;
+        }
+        
+        bool GurobiLpSolver::isOptimal() const {
+            if (!this->currentModelHasBeenOptimized) {
+                return false;
+            }
+            int optimalityStatus = 0;
+            
+            int error = GRBgetintattr(model, GRB_INT_ATTR_STATUS, &optimalityStatus);
+            if (error) {
+                LOG4CPLUS_ERROR(logger, "Unable to retrieve optimization status of Gurobi model (" << GRBgeterrormsg(env) << ").");
+                throw storm::exceptions::InvalidStateException() << "Unable to retrieve optimization status of Gurobi model (" << GRBgeterrormsg(env) << ").";
+            }
+            
+            return optimalityStatus == GRB_OPTIMAL;
         }
         
         int_fast64_t GurobiLpSolver::getIntegerValue(uint_fast64_t variableIndex) const {
-            if (!isOptimal) {
-                throw storm::exceptions::InvalidStateException() << "Unable to get Gurobi solution from unoptimized model (" << GRBgeterrormsg(env) << ").";
+            if (!this->isOptimal()) {
+                if (this->isInfeasible()) {
+                    LOG4CPLUS_ERROR(logger, "Unable to get Gurobi solution from infeasible model (" << GRBgeterrormsg(env) << ").");
+                    throw storm::exceptions::InvalidStateException() << "Unable to get Gurobi solution from infeasible model (" << GRBgeterrormsg(env) << ").";
+                } else if (this->isUnbounded()) {
+                    LOG4CPLUS_ERROR(logger, "Unable to get Gurobi solution from unbounded model (" << GRBgeterrormsg(env) << ").");
+                    throw storm::exceptions::InvalidStateException() << "Unable to get Gurobi solution from unbounded model (" << GRBgeterrormsg(env) << ").";
+                } else {
+                    LOG4CPLUS_ERROR(logger, "Unable to get Gurobi solution from unoptimized model (" << GRBgeterrormsg(env) << ").");
+                    throw storm::exceptions::InvalidStateException() << "Unable to get Gurobi solution from unoptimized model (" << GRBgeterrormsg(env) << ").";
+                }
             }
             
             double value = 0;
@@ -222,8 +268,17 @@ namespace storm {
         }
         
         bool GurobiLpSolver::getBinaryValue(uint_fast64_t variableIndex) const {
-            if (!isOptimal) {
-                throw storm::exceptions::InvalidStateException() << "Unable to get Gurobi solution from unoptimized model (" << GRBgeterrormsg(env) << ").";
+            if (!this->isOptimal()) {
+                if (this->isInfeasible()) {
+                    LOG4CPLUS_ERROR(logger, "Unable to get Gurobi solution from infeasible model (" << GRBgeterrormsg(env) << ").");
+                    throw storm::exceptions::InvalidStateException() << "Unable to get Gurobi solution from infeasible model (" << GRBgeterrormsg(env) << ").";
+                } else if (this->isUnbounded()) {
+                    LOG4CPLUS_ERROR(logger, "Unable to get Gurobi solution from unbounded model (" << GRBgeterrormsg(env) << ").");
+                    throw storm::exceptions::InvalidStateException() << "Unable to get Gurobi solution from unbounded model (" << GRBgeterrormsg(env) << ").";
+                } else {
+                    LOG4CPLUS_ERROR(logger, "Unable to get Gurobi solution from unoptimized model (" << GRBgeterrormsg(env) << ").");
+                    throw storm::exceptions::InvalidStateException() << "Unable to get Gurobi solution from unoptimized model (" << GRBgeterrormsg(env) << ").";
+                }
             }
 
             double value = 0;
@@ -244,8 +299,17 @@ namespace storm {
         }
         
         double GurobiLpSolver::getContinuousValue(uint_fast64_t variableIndex) const {
-            if (!isOptimal) {
-                throw storm::exceptions::InvalidStateException() << "Unable to get Gurobi solution from unoptimized model (" << GRBgeterrormsg(env) << ").";
+            if (!this->isOptimal()) {
+                if (this->isInfeasible()) {
+                    LOG4CPLUS_ERROR(logger, "Unable to get Gurobi solution from infeasible model (" << GRBgeterrormsg(env) << ").");
+                    throw storm::exceptions::InvalidStateException() << "Unable to get Gurobi solution from infeasible model (" << GRBgeterrormsg(env) << ").";
+                } else if (this->isUnbounded()) {
+                    LOG4CPLUS_ERROR(logger, "Unable to get Gurobi solution from unbounded model (" << GRBgeterrormsg(env) << ").");
+                    throw storm::exceptions::InvalidStateException() << "Unable to get Gurobi solution from unbounded model (" << GRBgeterrormsg(env) << ").";
+                } else {
+                    LOG4CPLUS_ERROR(logger, "Unable to get Gurobi solution from unoptimized model (" << GRBgeterrormsg(env) << ").");
+                    throw storm::exceptions::InvalidStateException() << "Unable to get Gurobi solution from unoptimized model (" << GRBgeterrormsg(env) << ").";
+                }
             }
 
             double value = 0;
diff --git a/src/solver/GurobiLpSolver.h b/src/solver/GurobiLpSolver.h
index 0bd5cf0c9..3e361843e 100644
--- a/src/solver/GurobiLpSolver.h
+++ b/src/solver/GurobiLpSolver.h
@@ -19,10 +19,37 @@ namespace storm {
     namespace solver {
 
 #ifdef STORM_HAVE_GUROBI
+        /*!
+         * A class that implements the LpSolver interface using Gurobi.
+         */
         class GurobiLpSolver : public LpSolver {
         public:
+            /*!
+             * Constructs a solver with the given name and model sense.
+             *
+             * @param name The name of the LP problem.
+             * @param modelSense A value indicating whether the value of the objective function is to be minimized or
+             * maximized.
+             */
             GurobiLpSolver(std::string const& name, ModelSense const& modelSense);
+            
+            /*!
+             * Constructs a solver with the given name. By default the objective function is assumed to be minimized,
+             * but this may be altered later using a call to setModelSense.
+             *
+             * @param name The name of the LP problem.
+             */
             GurobiLpSolver(std::string const& name);
+            
+            /*!
+             * Constructs a solver without a name. By default the objective function is assumed to be minimized,
+             * but this may be altered later using a call to setModelSense.
+             */
+            GurobiLpSolver();
+            
+            /*!
+             * Destructs a solver by freeing the pointers to Gurobi's structures.
+             */
             virtual ~GurobiLpSolver();
             
             virtual uint_fast64_t createContinuousVariable(std::string const& name, VariableType const& variableType, double lowerBound, double upperBound, double objectiveFunctionCoefficient) override;
@@ -31,9 +58,11 @@ namespace storm {
             
             virtual void addConstraint(std::string const& name, std::vector<uint_fast64_t> const& variables, std::vector<double> const& coefficients, BoundType const& boundType, double rightHandSideValue) override;
             
-            virtual void setModelSense(ModelSense const& newModelSense) override;
             virtual void optimize() const override;
-            
+            virtual bool isInfeasible() const override;
+            virtual bool isUnbounded() const override;
+            virtual bool isOptimal() const override;
+
             virtual int_fast64_t getIntegerValue(uint_fast64_t variableIndex) const override;
             virtual bool getBinaryValue(uint_fast64_t variableIndex) const override;
             virtual double getContinuousValue(uint_fast64_t variableIndex) const override;
@@ -41,7 +70,14 @@ namespace storm {
             virtual void writeModelToFile(std::string const& filename) const override;
             
         private:
+            /*!
+             * Sets some properties of the Gurobi environment according to parameters given by the options.
+             */
             void setGurobiEnvironmentProperties() const;
+            
+            /*!
+             * Calls Gurobi to incorporate the latest changes to the model.
+             */
             void updateModel() const;
             
             // The Gurobi environment.
@@ -52,23 +88,23 @@ namespace storm {
             
             // A counter that keeps track of the next free variable index.
             uint_fast64_t nextVariableIndex;
-            
-            // A flag that stores whether the model was optimized properly.
-            mutable bool isOptimal;
         };
 #else
         // If Gurobi is not available, we provide a stub implementation that emits an error if any of its methods is called.
         class GurobiLpSolver : public LpSolver {
         public:
-
-        	GurobiLpSolver(std::string const& name) : LpSolver(MINIMIZE) {
-        		throw storm::exceptions::NotImplementedException() << "This version of StoRM was compiled without support for Gurobi. Yet, a method was called that requires this support. Please choose a version of support with Gurobi support.";
-        	}
-
-        	virtual ~GurobiLpSolver() {
-        		throw storm::exceptions::NotImplementedException() << "This version of StoRM was compiled without support for Gurobi. Yet, a method was called that requires this support. Please choose a version of support with Gurobi support.";
-        	}
-
+            GurobiLpSolver(std::string const& name, ModelSense const& modelSense) {
+                throw storm::exceptions::NotImplementedException() << "This version of StoRM was compiled without support for Gurobi. Yet, a method was called that requires this support. Please choose a version of support with Gurobi support.";
+            }
+            
+            GurobiLpSolver(std::string const& name) {
+                throw storm::exceptions::NotImplementedException() << "This version of StoRM was compiled without support for Gurobi. Yet, a method was called that requires this support. Please choose a version of support with Gurobi support.";
+            }
+            
+            GurobiLpSolver() {
+                throw storm::exceptions::NotImplementedException() << "This version of StoRM was compiled without support for Gurobi. Yet, a method was called that requires this support. Please choose a version of support with Gurobi support.";
+            }
+            
             virtual uint_fast64_t createContinuousVariable(std::string const& name, VariableType const& variableType, double lowerBound, double upperBound, double objectiveFunctionCoefficient) override {
                 throw storm::exceptions::NotImplementedException() << "This version of StoRM was compiled without support for Gurobi. Yet, a method was called that requires this support. Please choose a version of support with Gurobi support.";
             }
@@ -85,7 +121,7 @@ namespace storm {
                 throw storm::exceptions::NotImplementedException() << "This version of StoRM was compiled without support for Gurobi. Yet, a method was called that requires this support. Please choose a version of support with Gurobi support.";
             }
             
-            virtual void setModelSense(ModelSense const& newModelSense) {
+            virtual void setModelSense(ModelSense const& modelSense) {
                 throw storm::exceptions::NotImplementedException() << "This version of StoRM was compiled without support for Gurobi. Yet, a method was called that requires this support. Please choose a version of support with Gurobi support.";
             }
             
@@ -93,6 +129,18 @@ namespace storm {
                 throw storm::exceptions::NotImplementedException() << "This version of StoRM was compiled without support for Gurobi. Yet, a method was called that requires this support. Please choose a version of support with Gurobi support.";
             }
             
+            virtual bool isInfeasible() const override {
+                throw storm::exceptions::NotImplementedException() << "This version of StoRM was compiled without support for Gurobi. Yet, a method was called that requires this support. Please choose a version of support with Gurobi support.";
+            }
+            
+            virtual bool isUnbounded() const override {
+                throw storm::exceptions::NotImplementedException() << "This version of StoRM was compiled without support for Gurobi. Yet, a method was called that requires this support. Please choose a version of support with Gurobi support.";
+            }
+            
+            virtual bool isOptimal() const override {
+                throw storm::exceptions::NotImplementedException() << "This version of StoRM was compiled without support for Gurobi. Yet, a method was called that requires this support. Please choose a version of support with Gurobi support.";
+            }
+            
             virtual int_fast64_t getIntegerValue(uint_fast64_t variableIndex) const override {
                 throw storm::exceptions::NotImplementedException() << "This version of StoRM was compiled without support for Gurobi. Yet, a method was called that requires this support. Please choose a version of support with Gurobi support.";
             }
diff --git a/src/solver/AbstractLinearEquationSolver.h b/src/solver/LinearEquationSolver.h
similarity index 86%
rename from src/solver/AbstractLinearEquationSolver.h
rename to src/solver/LinearEquationSolver.h
index 39c5cab37..2a0f6579e 100644
--- a/src/solver/AbstractLinearEquationSolver.h
+++ b/src/solver/LinearEquationSolver.h
@@ -1,5 +1,5 @@
-#ifndef STORM_SOLVER_ABSTRACTLINEAREQUATIONSOLVER_H_
-#define STORM_SOLVER_ABSTRACTLINEAREQUATIONSOLVER_H_
+#ifndef STORM_SOLVER_LINEAREQUATIONSOLVER_H_
+#define STORM_SOLVER_LINEAREQUATIONSOLVER_H_
 
 #include <vector>
 
@@ -9,18 +9,18 @@ namespace storm {
     namespace solver {
         
         /*!
-         * A class that represents an abstract linear equation solver. In addition to solving a system of linear
+         * An interface that represents an abstract linear equation solver. In addition to solving a system of linear
          * equations, the functionality to repeatedly multiply a matrix with a given vector is provided.
          */
         template<class Type>
-        class AbstractLinearEquationSolver {
+        class LinearEquationSolver {
         public:
             /*!
              * Makes a copy of the linear equation solver.
              *
              * @return A pointer to a copy of the linear equation solver.
              */
-            virtual AbstractLinearEquationSolver<Type>* clone() const = 0;
+            virtual LinearEquationSolver<Type>* clone() const = 0;
             
             /*!
              * Solves the equation system A*x = b. The matrix A is required to be square and have a unique solution.
@@ -52,4 +52,4 @@ namespace storm {
     } // namespace solver
 } // namespace storm
 
-#endif /* STORM_SOLVER_ABSTRACTLINEAREQUATIONSOLVER_H_ */
+#endif /* STORM_SOLVER_LINEAREQUATIONSOLVER_H_ */
diff --git a/src/solver/LpSolver.h b/src/solver/LpSolver.h
index e286e22ab..96e8e469c 100644
--- a/src/solver/LpSolver.h
+++ b/src/solver/LpSolver.h
@@ -1,15 +1,19 @@
 #ifndef STORM_SOLVER_LPSOLVER
 #define STORM_SOLVER_LPSOLVER
 
-#include <cstdint>
 #include <string>
 #include <vector>
 
 namespace storm {
     namespace solver {
-        
+
+        /*!
+         * An interface that captures the functionality of an LP solver.
+         */
         class LpSolver {
         public:
+            // An enumeration to represent the possible types of variables. Variables may be either unbounded, have only
+            // a lower or an upper bound, respectively, or be bounded from below and from above.
             enum VariableType {
                 UNBOUNDED,
                 LOWER_BOUND,
@@ -17,6 +21,7 @@ namespace storm {
                 BOUNDED,
             };
             
+            // An enumeration to represent the possible types of constraints.
             enum BoundType {
                 LESS,
                 LESS_EQUAL,
@@ -25,38 +30,183 @@ namespace storm {
                 EQUAL
             };
             
+            // An enumeration to represent whether the objective function is to be minimized or maximized.
             enum ModelSense {
                 MINIMIZE,
                 MAXIMIZE
             };
             
-            LpSolver(ModelSense const& modelSense) : modelSense(modelSense) {
+            /*!
+             * Creates an empty LP solver. By default the objective function is assumed to be minimized.
+             */
+            LpSolver() : currentModelHasBeenOptimized(false), modelSense(MINIMIZE) {
                 // Intentionally left empty.
             }
             
+            /*!
+             * Creates an empty LP solver with the given model sense.
+             *
+             * @param modelSense A value indicating whether the objective function of this model is to be minimized or
+             * maximized.
+             */
+            LpSolver(ModelSense const& modelSense) : currentModelHasBeenOptimized(false), modelSense(modelSense) {
+                // Intentionally left empty.
+            }
+
+            /*!
+             * Creates a continuous variable, i.e. a variable that may take all real values within its bounds.
+             *
+             * @param name The name of the variable.
+             * @param variableType The type of the variable.
+             * @param lowerBound The lower bound of the variable. Note that this parameter is ignored if the variable 
+             * is not bounded from below.
+             * @param upperBound The upper bound of the variable. Note that this parameter is ignored if the variable
+             * is not bounded from above.
+             * @param objectiveFunctionCoefficient The coefficient with which the variable appears in the objective
+             * function. If set to zero, the variable is irrelevant for the objective value.
+             * @return The index of the newly created variable. This index can be used to retrieve the variables value
+             * in a solution after the model has been optimized.
+             */
             virtual uint_fast64_t createContinuousVariable(std::string const& name, VariableType const& variableType, double lowerBound, double upperBound, double objectiveFunctionCoefficient) = 0;
+            
+            /*!
+             * Creates an integer variable.
+             *
+             * @param name The name of the variable.
+             * @param variableType The type of the variable.
+             * @param lowerBound The lower bound of the variable. Note that this parameter is ignored if the variable
+             * is not bounded from below.
+             * @param upperBound The upper bound of the variable. Note that this parameter is ignored if the variable
+             * is not bounded from above.
+             * @param objectiveFunctionCoefficient The coefficient with which the variable appears in the objective
+             * function. If set to zero, the variable is irrelevant for the objective value.
+             * @return The index of the newly created variable. This index can be used to retrieve the variables value
+             * in a solution after the model has been optimized.
+             */
             virtual uint_fast64_t createIntegerVariable(std::string const& name, VariableType const& variableType, double lowerBound, double upperBound, double objectiveFunctionCoefficient) = 0;
+            
+            /*!
+             * Creates a binary variable, i.e. a variable that may be either zero or one.
+             *
+             * @param name The name of the variable.
+             * @param variableType The type of the variable.
+             * @param lowerBound The lower bound of the variable. Note that this parameter is ignored if the variable
+             * is not bounded from below.
+             * @param upperBound The upper bound of the variable. Note that this parameter is ignored if the variable
+             * is not bounded from above.
+             * @param objectiveFunctionCoefficient The coefficient with which the variable appears in the objective
+             * function. If set to zero, the variable is irrelevant for the objective value.
+             * @return The index of the newly created variable. This index can be used to retrieve the variables value
+             * in a solution after the model has been optimized.
+             */
             virtual uint_fast64_t createBinaryVariable(std::string const& name, double objectiveFunctionCoefficient) = 0;
             
-            virtual void addConstraint(std::string const& name, std::vector<uint_fast64_t> const& variables, std::vector<double> const& coefficients, BoundType const& boundType, double rightHandSideValue) = 0;
+            /*!
+             * Adds a constraint of the form
+             *      a_1*x_1 + ... + a_n*x_n  op  c
+             * to the LP problem.
+             *
+             * @param name The name of the constraint. If empty, the constraint has no name.
+             * @param variables A vector of variable indices that define the appearing variables x_1, ..., x_n.
+             * @param coefficients A vector of coefficients that define the a_1, ..., a_n. The i-ith entry specifies the
+             * coefficient of the variable whose index appears at the i-th position in the vector of variable indices.
+             * @param boundType The bound type specifies the operator op used in the constraint.
+             * @param rhs This defines the value of the constant appearing on the right-hand side of the constraint.
+             */
+            virtual void addConstraint(std::string const& name, std::vector<uint_fast64_t> const& variables, std::vector<double> const& coefficients, BoundType const& boundType, double rhs) = 0;
             
+            /*!
+             * Optimizes the LP problem previously constructed. Afterwards, the methods isInfeasible, isUnbounded and
+             * isOptimal can be used to query the optimality status.
+             */
             virtual void optimize() const = 0;
             
+            /*!
+             * Retrieves whether the model was found to be infeasible. This can only be called after the model has been
+             * optimized and not modified afterwards. Otherwise, an exception is thrown.
+             *
+             * @return True if the model was optimized and found to be infeasible.
+             */
+            virtual bool isInfeasible() const = 0;
+            
+            /*!
+             * Retrieves whether the model was found to be infeasible. This can only be called after the model has been
+             * optimized and not modified afterwards. Otherwise, an exception is thrown.
+             *
+             * @return True if the model was optimized and found to be unbounded.
+             */
+            virtual bool isUnbounded() const = 0;
+            
+            /*!
+             * Retrieves whether the model was found to be optimal, i.e. neither infeasible nor unbounded. This can only
+             * be called after the model has been optimized and not modified afterwards. Otherwise, an exception is
+             * thrown.
+             *
+             * @return True if the model was optimized and found to be neither infeasible nor unbounded.
+             */
+            virtual bool isOptimal() const = 0;
+            
+            /*!
+             * Retrieves the value of the integer variable with the given index. Note that this may only be called, if
+             * the model was found to be optimal, i.e. iff isOptimal() returns true.
+             *
+             * @param variableIndex The index of the integer variable whose value to query. If this index does not
+             * belong to a previously declared integer variable, the behaviour is undefined.
+             */
             virtual int_fast64_t getIntegerValue(uint_fast64_t variableIndex) const = 0;
+            
+            /*!
+             * Retrieves the value of the binary variable with the given index. Note that this may only be called, if
+             * the model was found to be optimal, i.e. iff isOptimal() returns true.
+             *
+             * @param variableIndex The index of the binary variable whose value to query. If this index does not
+             * belong to a previously declared binary variable, the behaviour is undefined.
+             */
             virtual bool getBinaryValue(uint_fast64_t variableIndex) const = 0;
+            
+            /*!
+             * Retrieves the value of the continuous variable with the given index. Note that this may only be called,
+             * if the model was found to be optimal, i.e. iff isOptimal() returns true.
+             *
+             * @param variableIndex The index of the continuous variable whose value to query. If this index does not
+             * belong to a previously declared continuous variable, the behaviour is undefined.
+             */
             virtual double getContinuousValue(uint_fast64_t variableIndex) const = 0;
 
+            /*!
+             * Writes the current LP problem to the given file.
+             *
+             * @param filename The file to which to write the string representation of the LP.
+             */
             virtual void writeModelToFile(std::string const& filename) const = 0;
             
-            virtual void setModelSense(ModelSense const& newModelSense) {
-                this->modelSense = newModelSense;
+            /*!
+             * Sets whether the objective function of this model is to be minimized or maximized.
+             *
+             * @param modelSense The model sense to use.
+             */
+            void setModelSense(ModelSense const& modelSense) {
+                if (modelSense != this->modelSense) {
+                    currentModelHasBeenOptimized = false;
+                }
+                this->modelSense = modelSense;
             }
             
+            /*!
+             * Retrieves whether the objective function of this model is to be minimized or maximized.
+             *
+             * @return A value indicating whether the objective function of this model is to be minimized or maximized.
+             */
             ModelSense getModelSense() const {
                 return modelSense;
             }
             
+        protected:
+            // A flag indicating whether the current model has been optimized and not changed afterwards.
+            mutable bool currentModelHasBeenOptimized;
+            
         private:
+            // A flag that indicates the model sense.
             ModelSense modelSense;
         };
     }
diff --git a/src/solver/NativeLinearEquationSolver.cpp b/src/solver/NativeLinearEquationSolver.cpp
index 0176bf24f..21f09203c 100644
--- a/src/solver/NativeLinearEquationSolver.cpp
+++ b/src/solver/NativeLinearEquationSolver.cpp
@@ -44,7 +44,7 @@ namespace storm {
         }
         
         template<typename ValueType>
-        AbstractLinearEquationSolver<ValueType>* NativeLinearEquationSolver<ValueType>::clone() const {
+        LinearEquationSolver<ValueType>* NativeLinearEquationSolver<ValueType>::clone() const {
             return new NativeLinearEquationSolver<ValueType>(*this);
         }
         
diff --git a/src/solver/NativeLinearEquationSolver.h b/src/solver/NativeLinearEquationSolver.h
index ab3ce2d47..af94cdb43 100644
--- a/src/solver/NativeLinearEquationSolver.h
+++ b/src/solver/NativeLinearEquationSolver.h
@@ -1,16 +1,16 @@
 #ifndef STORM_SOLVER_NATIVELINEAREQUATIONSOLVER_H_
 #define STORM_SOLVER_NATIVELINEAREQUATIONSOLVER_H_
 
-#include "AbstractLinearEquationSolver.h"
+#include "LinearEquationSolver.h"
 
 namespace storm {
     namespace solver {
         
         /*!
-         * A class that uses StoRM's native matrix operations to implement the AbstractLinearEquationSolver interface.
+         * A class that uses StoRM's native matrix operations to implement the LinearEquationSolver interface.
          */
         template<typename ValueType>
-        class NativeLinearEquationSolver : public AbstractLinearEquationSolver<ValueType> {
+        class NativeLinearEquationSolver : public LinearEquationSolver<ValueType> {
         public:
             // An enumeration specifying the available solution methods.
             enum SolutionMethod {
@@ -33,7 +33,7 @@ namespace storm {
              */
             NativeLinearEquationSolver(SolutionMethod method, double precision, uint_fast64_t maximalNumberOfIterations, bool relative = true);
             
-            virtual AbstractLinearEquationSolver<ValueType>* clone() const override;
+            virtual LinearEquationSolver<ValueType>* clone() const override;
             
             virtual void solveEquationSystem(storm::storage::SparseMatrix<ValueType> const& A, std::vector<ValueType>& x, std::vector<ValueType> const& b, std::vector<ValueType>* multiplyResult = nullptr) const override;
             
diff --git a/src/solver/NativeNondeterministicLinearEquationSolver.cpp b/src/solver/NativeNondeterministicLinearEquationSolver.cpp
index df7edd137..6f7459ef7 100644
--- a/src/solver/NativeNondeterministicLinearEquationSolver.cpp
+++ b/src/solver/NativeNondeterministicLinearEquationSolver.cpp
@@ -32,7 +32,7 @@ namespace storm {
         }
         
         template<typename ValueType>
-        AbstractNondeterministicLinearEquationSolver<ValueType>* NativeNondeterministicLinearEquationSolver<ValueType>::clone() const {
+        NondeterministicLinearEquationSolver<ValueType>* NativeNondeterministicLinearEquationSolver<ValueType>::clone() const {
             return new NativeNondeterministicLinearEquationSolver<ValueType>(*this);
         }
         
diff --git a/src/solver/NativeNondeterministicLinearEquationSolver.h b/src/solver/NativeNondeterministicLinearEquationSolver.h
index 4e5126254..51c8ddd97 100644
--- a/src/solver/NativeNondeterministicLinearEquationSolver.h
+++ b/src/solver/NativeNondeterministicLinearEquationSolver.h
@@ -1,16 +1,16 @@
 #ifndef STORM_SOLVER_NATIVENONDETERMINISTICLINEAREQUATIONSOLVER_H_
 #define STORM_SOLVER_NATIVENONDETERMINISTICLINEAREQUATIONSOLVER_H_
 
-#include "src/solver/AbstractNondeterministicLinearEquationSolver.h"
+#include "src/solver/NondeterministicLinearEquationSolver.h"
 
 namespace storm {
     namespace solver {
         
         /*!
-         * A class that uses the gmm++ library to implement the AbstractNondeterminsticLinearEquationSolver interface.
+         * A class that uses the gmm++ library to implement the NondeterminsticLinearEquationSolver interface.
          */
         template<class ValueType>
-        class NativeNondeterministicLinearEquationSolver : public AbstractNondeterministicLinearEquationSolver<ValueType> {
+        class NativeNondeterministicLinearEquationSolver : public NondeterministicLinearEquationSolver<ValueType> {
         public:
             /*!
              * Constructs a nondeterministic linear equation solver with parameters being set according to the settings
@@ -28,7 +28,7 @@ namespace storm {
              */
             NativeNondeterministicLinearEquationSolver(double precision, uint_fast64_t maximalNumberOfIterations, bool relative = true);
             
-            virtual AbstractNondeterministicLinearEquationSolver<ValueType>* clone() const override;
+            virtual NondeterministicLinearEquationSolver<ValueType>* clone() const override;
             
             virtual void performMatrixVectorMultiplication(bool minimize, storm::storage::SparseMatrix<ValueType> const& A, std::vector<ValueType>& x, std::vector<uint_fast64_t> const& nondeterministicChoiceIndices, std::vector<ValueType>* b = nullptr, uint_fast64_t n = 1, std::vector<ValueType>* newX = nullptr) const override;
             
diff --git a/src/solver/AbstractNondeterministicLinearEquationSolver.h b/src/solver/NondeterministicLinearEquationSolver.h
similarity index 86%
rename from src/solver/AbstractNondeterministicLinearEquationSolver.h
rename to src/solver/NondeterministicLinearEquationSolver.h
index 08158742b..8c5bc2334 100644
--- a/src/solver/AbstractNondeterministicLinearEquationSolver.h
+++ b/src/solver/NondeterministicLinearEquationSolver.h
@@ -1,5 +1,5 @@
-#ifndef STORM_SOLVER_ABSTRACTNONDETERMINISTICLINEAREQUATIONSOLVER_H_
-#define STORM_SOLVER_ABSTRACTNONDETERMINISTICLINEAREQUATIONSOLVER_H_
+#ifndef STORM_SOLVER_NONDETERMINISTICLINEAREQUATIONSOLVER_H_
+#define STORM_SOLVER_NONDETERMINISTICLINEAREQUATIONSOLVER_H_
 
 #include <vector>
 
@@ -9,18 +9,19 @@ namespace storm {
     namespace solver {
         
         /*!
-         * A class that represents an abstract nondeterministic linear equation solver. In addition to solving linear
-         * equation systems involving min/max operators, .
+         * A interface that represents an abstract nondeterministic linear equation solver. In addition to solving
+         * linear equation systems involving min/max operators, repeated matrix-vector multiplication functionality is
+         * provided.
          */
         template<class ValueType>
-        class AbstractNondeterministicLinearEquationSolver {
+        class NondeterministicLinearEquationSolver {
         public:
             /*!
              * Makes a copy of the nondeterministic linear equation solver.
              *
              * @return A pointer to a copy of the nondeterministic linear equation solver.
              */
-            virtual AbstractNondeterministicLinearEquationSolver<ValueType>* clone() const = 0;
+            virtual NondeterministicLinearEquationSolver<ValueType>* clone() const = 0;
             
             /*!
              * Solves the equation system x = min/max(A*x + b) given by the parameters.
@@ -64,4 +65,4 @@ namespace storm {
     } // namespace solver
 } // namespace storm
 
-#endif /* STORM_SOLVER_ABSTRACTNONDETERMINISTICLINEAREQUATIONSOLVER_H_ */
+#endif /* STORM_SOLVER_NONDETERMINISTICLINEAREQUATIONSOLVER_H_ */
diff --git a/src/storm.cpp b/src/storm.cpp
index 7600f3372..88ecafad4 100644
--- a/src/storm.cpp
+++ b/src/storm.cpp
@@ -466,10 +466,10 @@ int main(const int argc, const char* argv[]) {
                     storm::modelchecker::csl::SparseMarkovAutomatonCslModelChecker<double> mc(*markovAutomaton);
 //                    std::cout << mc.checkExpectedTime(true, markovAutomaton->getLabeledStates("goal")) << std::endl;
 //                    std::cout << mc.checkExpectedTime(false, markovAutomaton->getLabeledStates("goal")) << std::endl;
-//                    std::cout << mc.checkLongRunAverage(true, markovAutomaton->getLabeledStates("goal")) << std::endl;
-//                    std::cout << mc.checkLongRunAverage(false, markovAutomaton->getLabeledStates("goal")) << std::endl;
-                    std::cout << mc.checkTimeBoundedEventually(true, markovAutomaton->getLabeledStates("goal"), 0, 1) << std::endl;
-                    std::cout << mc.checkTimeBoundedEventually(true, markovAutomaton->getLabeledStates("goal"), 1, 2) << std::endl;
+                    std::cout << mc.checkLongRunAverage(true, markovAutomaton->getLabeledStates("goal")) << std::endl;
+                    std::cout << mc.checkLongRunAverage(false, markovAutomaton->getLabeledStates("goal")) << std::endl;
+//                    std::cout << mc.checkTimeBoundedEventually(true, markovAutomaton->getLabeledStates("goal"), 0, 1) << std::endl;
+//                    std::cout << mc.checkTimeBoundedEventually(true, markovAutomaton->getLabeledStates("goal"), 1, 2) << std::endl;
                     break;
                 }
 				case storm::models::Unknown:
diff --git a/src/utility/StormOptions.cpp b/src/utility/StormOptions.cpp
index 5d8eea24c..447057093 100644
--- a/src/utility/StormOptions.cpp
+++ b/src/utility/StormOptions.cpp
@@ -42,7 +42,7 @@ bool storm::utility::StormOptions::optionsRegistered = storm::settings::Settings
     std::vector<std::string> nondeterministicLinearEquationSolver;
 	nondeterministicLinearEquationSolver.push_back("gmm++");
 	nondeterministicLinearEquationSolver.push_back("native");
-	settings->addOption(storm::settings::OptionBuilder("StoRM Main", "linsolver", "", "Sets which solver is preferred for solving systems of linear equations arising from nondeterministic systems.").addArgument(storm::settings::ArgumentBuilder::createStringArgument("name", "The name of the solver to prefer. Available are: gmm++ and native.").addValidationFunctionString(storm::settings::ArgumentValidators::stringInListValidator(nondeterministicLinearEquationSolver)).setDefaultValueString("native").build()).build());
+	settings->addOption(storm::settings::OptionBuilder("StoRM Main", "ndsolver", "", "Sets which solver is preferred for solving systems of linear equations arising from nondeterministic systems.").addArgument(storm::settings::ArgumentBuilder::createStringArgument("name", "The name of the solver to prefer. Available are: gmm++ and native.").addValidationFunctionString(storm::settings::ArgumentValidators::stringInListValidator(nondeterministicLinearEquationSolver)).setDefaultValueString("native").build()).build());
 
     std::vector<std::string> lpSolvers;
 	lpSolvers.push_back("gurobi");
diff --git a/src/utility/solver.cpp b/src/utility/solver.cpp
index 93b811389..5ab8666a4 100644
--- a/src/utility/solver.cpp
+++ b/src/utility/solver.cpp
@@ -26,32 +26,32 @@ namespace storm {
             }
             
             template<typename ValueType>
-            std::shared_ptr<storm::solver::AbstractLinearEquationSolver<ValueType>> getLinearEquationSolver() {
+            std::shared_ptr<storm::solver::LinearEquationSolver<ValueType>> getLinearEquationSolver() {
                 std::string const& linearEquationSolver = storm::settings::Settings::getInstance()->getOptionByLongName("linsolver").getArgument(0).getValueAsString();
                 if (linearEquationSolver == "gmm++") {
-                    return std::shared_ptr<storm::solver::AbstractLinearEquationSolver<ValueType>>(new storm::solver::GmmxxLinearEquationSolver<ValueType>());
+                    return std::shared_ptr<storm::solver::LinearEquationSolver<ValueType>>(new storm::solver::GmmxxLinearEquationSolver<ValueType>());
                 } else if (linearEquationSolver == "native") {
-                    return std::shared_ptr<storm::solver::AbstractLinearEquationSolver<ValueType>>(new storm::solver::NativeLinearEquationSolver<ValueType>());
+                    return std::shared_ptr<storm::solver::LinearEquationSolver<ValueType>>(new storm::solver::NativeLinearEquationSolver<ValueType>());
                 }
                 
                 throw storm::exceptions::InvalidSettingsException() << "No suitable linear equation solver selected.";
             }
             
             template<typename ValueType>
-            std::shared_ptr<storm::solver::AbstractNondeterministicLinearEquationSolver<ValueType>> getNondeterministicLinearEquationSolver() {
+            std::shared_ptr<storm::solver::NondeterministicLinearEquationSolver<ValueType>> getNondeterministicLinearEquationSolver() {
                 std::string const& nondeterministicLinearEquationSolver = storm::settings::Settings::getInstance()->getOptionByLongName("ndsolver").getArgument(0).getValueAsString();
                 if (nondeterministicLinearEquationSolver == "gmm++") {
-                    return std::shared_ptr<storm::solver::AbstractNondeterministicLinearEquationSolver<ValueType>>(new storm::solver::GmmxxNondeterministicLinearEquationSolver<ValueType>());
+                    return std::shared_ptr<storm::solver::NondeterministicLinearEquationSolver<ValueType>>(new storm::solver::GmmxxNondeterministicLinearEquationSolver<ValueType>());
                 } else if (nondeterministicLinearEquationSolver == "native") {
-                    return std::shared_ptr<storm::solver::AbstractNondeterministicLinearEquationSolver<ValueType>>(new storm::solver::NativeNondeterministicLinearEquationSolver<ValueType>());
+                    return std::shared_ptr<storm::solver::NondeterministicLinearEquationSolver<ValueType>>(new storm::solver::NativeNondeterministicLinearEquationSolver<ValueType>());
                 }
                 
                 throw storm::exceptions::InvalidSettingsException() << "No suitable nondeterministic linear equation solver selected.";
             }
             
-            template std::shared_ptr<storm::solver::AbstractLinearEquationSolver<double>> getLinearEquationSolver();
+            template std::shared_ptr<storm::solver::LinearEquationSolver<double>> getLinearEquationSolver();
 
-            template std::shared_ptr<storm::solver::AbstractNondeterministicLinearEquationSolver<double>> getNondeterministicLinearEquationSolver();
+            template std::shared_ptr<storm::solver::NondeterministicLinearEquationSolver<double>> getNondeterministicLinearEquationSolver();
         }
     }
 }
\ No newline at end of file
diff --git a/src/utility/solver.h b/src/utility/solver.h
index 2475456d9..e5e7d456f 100644
--- a/src/utility/solver.h
+++ b/src/utility/solver.h
@@ -1,8 +1,8 @@
 #ifndef STORM_UTILITY_SOLVER_H_
 #define STORM_UTILITY_SOLVER_H_
 
-#include "src/solver/AbstractLinearEquationSolver.h"
-#include "src/solver/AbstractNondeterministicLinearEquationSolver.h"
+#include "src/solver/LinearEquationSolver.h"
+#include "src/solver/NondeterministicLinearEquationSolver.h"
 #include "src/solver/LpSolver.h"
 
 #include "src/exceptions/InvalidSettingsException.h"
@@ -14,10 +14,10 @@ namespace storm {
             std::shared_ptr<storm::solver::LpSolver> getLpSolver(std::string const& name);
             
             template<typename ValueType>
-            std::shared_ptr<storm::solver::AbstractLinearEquationSolver<ValueType>> getLinearEquationSolver();
+            std::shared_ptr<storm::solver::LinearEquationSolver<ValueType>> getLinearEquationSolver();
 
             template<typename ValueType>
-            std::shared_ptr<storm::solver::AbstractNondeterministicLinearEquationSolver<ValueType>> getNondeterministicLinearEquationSolver();
+            std::shared_ptr<storm::solver::NondeterministicLinearEquationSolver<ValueType>> getNondeterministicLinearEquationSolver();
 
         }
     }
diff --git a/src/utility/vector.h b/src/utility/vector.h
index f779f2f4a..dd3bb2a1f 100644
--- a/src/utility/vector.h
+++ b/src/utility/vector.h
@@ -151,7 +151,7 @@ namespace storm {
                                       std::transform(target.begin() + range.begin(), target.begin() + range.end(), secondOperand.begin() + range.begin(), target.begin() + range.begin(), function);
                                   });
 #else
-                std::transform(target.begin(), target.end(), secondOperand.begin(), secondOperand.begin(), function);
+                std::transform(target.begin(), target.end(), secondOperand.begin(), target.begin(), function);
 #endif
             }
             
diff --git a/test/performance/modelchecker/SparseMdpPrctlModelCheckerTest.cpp b/test/performance/modelchecker/SparseMdpPrctlModelCheckerTest.cpp
index a95aba8d7..a4c6064df 100644
--- a/test/performance/modelchecker/SparseMdpPrctlModelCheckerTest.cpp
+++ b/test/performance/modelchecker/SparseMdpPrctlModelCheckerTest.cpp
@@ -166,7 +166,7 @@ TEST(SparseMdpPrctlModelCheckerTest, Consensus) {
     
 	result = mc.checkNoBoundOperator(*rewardFormula);
     
-	ASSERT_LT(std::abs(result[31168] - 2183.1424220082612919213715), s->getOptionByLongName("precision").getArgument(0).getValueAsDouble());
+	ASSERT_LT(std::abs(result[31168] - 2186.97175), s->getOptionByLongName("precision").getArgument(0).getValueAsDouble());
 	delete rewardFormula;
 
 }
\ No newline at end of file