added feature to linear equation solver factories to take posession of the matrix to forward it to the solvers

Former-commit-id: ed183f1820

src/modelchecker/csl/HybridCtmcCslModelChecker.cpp

@@ -15,7 +15,7 @@
 namespace storm {
     namespace modelchecker {
         template<storm::dd::DdType DdType, class ValueType>
-        HybridCtmcCslModelChecker<DdType, ValueType>::HybridCtmcCslModelChecker(storm::models::symbolic::Ctmc<DdType, ValueType> const& model) : SymbolicPropositionalModelChecker<DdType, ValueType>(model), linearEquationSolverFactory(new storm::utility::solver::LinearEquationSolverFactory<ValueType>()) {
+        HybridCtmcCslModelChecker<DdType, ValueType>::HybridCtmcCslModelChecker(storm::models::symbolic::Ctmc<DdType, ValueType> const& model) : SymbolicPropositionalModelChecker<DdType, ValueType>(model), linearEquationSolverFactory(std::make_unique<storm::utility::solver::GeneralLinearEquationSolverFactory<ValueType>>()) {
             // Intentionally left empty.
         }
         

src/modelchecker/csl/SparseCtmcCslModelChecker.cpp

@@ -24,7 +24,7 @@
 namespace storm {
     namespace modelchecker {
         template <typename SparseCtmcModelType>
-        SparseCtmcCslModelChecker<SparseCtmcModelType>::SparseCtmcCslModelChecker(SparseCtmcModelType const& model) : SparsePropositionalModelChecker<SparseCtmcModelType>(model), linearEquationSolverFactory(new storm::utility::solver::LinearEquationSolverFactory<ValueType>()) {
+        SparseCtmcCslModelChecker<SparseCtmcModelType>::SparseCtmcCslModelChecker(SparseCtmcModelType const& model) : SparsePropositionalModelChecker<SparseCtmcModelType>(model), linearEquationSolverFactory(std::make_unique<storm::utility::solver::GeneralLinearEquationSolverFactory<ValueType>>()) {
             // Intentionally left empty.
         }
         

src/modelchecker/prctl/HybridDtmcPrctlModelChecker.cpp

@@ -31,7 +31,7 @@ namespace storm {
         }
         
         template<storm::dd::DdType DdType, typename ValueType>
-        HybridDtmcPrctlModelChecker<DdType, ValueType>::HybridDtmcPrctlModelChecker(storm::models::symbolic::Dtmc<DdType, ValueType> const& model) : SymbolicPropositionalModelChecker<DdType, ValueType>(model), linearEquationSolverFactory(new storm::utility::solver::LinearEquationSolverFactory<ValueType>()) {
+        HybridDtmcPrctlModelChecker<DdType, ValueType>::HybridDtmcPrctlModelChecker(storm::models::symbolic::Dtmc<DdType, ValueType> const& model) : SymbolicPropositionalModelChecker<DdType, ValueType>(model), linearEquationSolverFactory(std::make_unique<storm::utility::solver::GeneralLinearEquationSolverFactory<ValueType>>()) {
             // Intentionally left empty.
         }
         

src/modelchecker/prctl/SparseDtmcPrctlModelChecker.cpp

@@ -29,7 +29,7 @@ namespace storm {
         }
         
         template<typename SparseDtmcModelType>
-        SparseDtmcPrctlModelChecker<SparseDtmcModelType>::SparseDtmcPrctlModelChecker(SparseDtmcModelType const& model) : SparsePropositionalModelChecker<SparseDtmcModelType>(model), linearEquationSolverFactory(new storm::utility::solver::LinearEquationSolverFactory<ValueType>()) {
+        SparseDtmcPrctlModelChecker<SparseDtmcModelType>::SparseDtmcPrctlModelChecker(SparseDtmcModelType const& model) : SparsePropositionalModelChecker<SparseDtmcModelType>(model), linearEquationSolverFactory(std::make_unique<storm::utility::solver::GeneralLinearEquationSolverFactory<ValueType>>()) {
             // Intentionally left empty.
         }
         

src/solver/EliminationLinearEquationSolver.cpp

@@ -8,6 +8,7 @@
 #include "src/solver/stateelimination/StatePriorityQueue.h"
 #include "src/solver/stateelimination/PrioritizedStateEliminator.h"
 
+#include "src/utility/graph.h"
 #include "src/utility/macros.h"
 #include "src/utility/vector.h"
 #include "src/utility/stateelimination.h"
@@ -27,12 +28,18 @@ namespace storm {
         void EliminationLinearEquationSolver<ValueType>::solveEquationSystem(std::vector<ValueType>& x, std::vector<ValueType> const& b, std::vector<ValueType>* multiplyResult) const {
             STORM_LOG_WARN_COND(multiplyResult == nullptr, "Providing scratch memory will not improve the performance of this solver.");
             
+            // FIXME: This solver will not work for all input systems. More concretely, the current implementation will
+            // not work for systems that have a 0 on the diagonal. This is not a restriction of this technique in general
+            // but arbitrary matrices require pivoting, which is not currently implemented.
+            
             STORM_LOG_DEBUG("Solving equation system using elimination.");
             
             // We need to revert the transformation into an equation system matrix, because the elimination procedure
             // and the distance computation is based on the probability matrix instead.
             storm::storage::SparseMatrix<ValueType> transitionMatrix(A);
+            
             transitionMatrix.convertToEquationSystem();
+            
             storm::storage::SparseMatrix<ValueType> backwardTransitions = transitionMatrix.transpose();
 
             // Initialize the solution to the right-hand side of the equation system.
@@ -45,11 +52,11 @@ namespace storm {
             storm::settings::modules::EliminationSettings::EliminationOrder order = storm::settings::getModule<storm::settings::modules::EliminationSettings>().getEliminationOrder();
             boost::optional<std::vector<uint_fast64_t>> distanceBasedPriorities;
             if (eliminationOrderNeedsDistances(order)) {
-                // Compute the distance from the first state.
-                // FIXME: This is not exactly the initial states.
-                storm::storage::BitVector initialStates = storm::storage::BitVector(A.getRowCount());
-                initialStates.set(0);
-                distanceBasedPriorities = getDistanceBasedPriorities(transitionMatrix, backwardTransitions, initialStates, b, eliminationOrderNeedsForwardDistances(order), eliminationOrderNeedsReversedDistances(order));
+                // Since we have no initial states at this point, we determine a representative of every BSCC regarding
+                // the backward transitions, because this means that every row is reachable from this set of rows, which
+                // we require to make sure we cover every row.
+                storm::storage::BitVector initialRows = storm::utility::graph::getBsccCover(backwardTransitions);
+                distanceBasedPriorities = getDistanceBasedPriorities(transitionMatrix, backwardTransitions, initialRows, b, eliminationOrderNeedsForwardDistances(order), eliminationOrderNeedsReversedDistances(order));
             }
             
             std::shared_ptr<StatePriorityQueue> priorityQueue = createStatePriorityQueue<ValueType>(distanceBasedPriorities, flexibleMatrix, flexibleBackwardTransitions, b, storm::storage::BitVector(x.size(), true));

src/solver/NativeLinearEquationSolver.h

@@ -16,7 +16,7 @@ namespace storm {
         std::ostream& operator<<(std::ostream& out, NativeLinearEquationSolverSolutionMethod const& method);
         
         /*!
-         * A class that uses StoRM's native matrix operations to implement the LinearEquationSolver interface.
+         * A class that uses Storm's native matrix operations to implement the LinearEquationSolver interface.
          */
         template<typename ValueType>
         class NativeLinearEquationSolver : public LinearEquationSolver<ValueType> {

src/utility/graph.cpp

@@ -5,18 +5,22 @@
 #include "src/adapters/CarlAdapter.h"
 
 #include "src/storage/sparse/StateType.h"
+#include "src/storage/dd/Bdd.h"
+#include "src/storage/dd/Add.h"
+#include "src/storage/dd/DdManager.h"
+
+#include "src/storage/StronglyConnectedComponentDecomposition.h"
+
 #include "src/models/symbolic/DeterministicModel.h"
 #include "src/models/symbolic/NondeterministicModel.h"
 #include "src/models/symbolic/StandardRewardModel.h"
 #include "src/models/sparse/DeterministicModel.h"
 #include "src/models/sparse/NondeterministicModel.h"
 #include "src/models/sparse/StandardRewardModel.h"
+
 #include "src/utility/constants.h"
-#include "src/exceptions/InvalidArgumentException.h"
-#include "src/storage/dd/Bdd.h"
-#include "src/storage/dd/Add.h"
-#include "src/storage/dd/DdManager.h"
 #include "src/utility/macros.h"
+#include "src/exceptions/InvalidArgumentException.h"
 
 namespace storm {
     namespace utility {
@@ -82,6 +86,61 @@ namespace storm {
                 return reachableStates;
             }
             
+            template<typename T>
+            storm::storage::BitVector getBsccCover(storm::storage::SparseMatrix<T> const& transitionMatrix) {
+                storm::storage::BitVector result(transitionMatrix.getRowGroupCount());
+                storm::storage::StronglyConnectedComponentDecomposition<T> decomposition(transitionMatrix, false, true);
+                
+                // Take the first state out of each BSCC.
+                for (auto const& scc : decomposition) {
+                    result.set(*scc.begin());
+                }
+                
+                return result;
+            }
+            
+            template<typename T>
+            storm::storage::BitVector getTerminalStateCover(storm::storage::SparseMatrix<T> const& transitionMatrix, storm::storage::BitVector const& initialStates) {
+                storm::storage::BitVector terminalStateCandidates(transitionMatrix.getRowGroupCount());
+                storm::storage::BitVector terminalStatesWithoutSuccessors(transitionMatrix.getRowCount());
+                
+                std::queue<storm::storage::sparse::state_type> stateQueue;
+                storm::storage::BitVector statesInQueue(transitionMatrix.getRowGroupCount());
+                
+                for (auto const& initialState : initialStates) {
+                    stateQueue.emplace(initialState);
+                    statesInQueue.set(initialState);
+                }
+                
+                // Perform a BFS.
+                while (!stateQueue.empty()) {
+                    storm::storage::sparse::state_type currentState = stateQueue.front();
+                    stateQueue.pop();
+                    
+                    auto row = transitionMatrix.getRow(currentState);
+                    if (row.empty()) {
+                        terminalStatesWithoutSuccessors.set(currentState);
+                    } else {
+                        bool hasUnvisitedSuccessor = false;
+                        for (auto const& successorEntry : row) {
+                            if (!statesInQueue.get(successorEntry.getColumn())) {
+                                hasUnvisitedSuccessor = true;
+                                stateQueue.emplace(successorEntry.getColumn());
+                                statesInQueue.set(successorEntry.getColumn());
+                            }
+                        }
+                        if (!hasUnvisitedSuccessor) {
+                            terminalStateCandidates.set(currentState);
+                        }
+                    }
+                }
+                
+                // Now that we have an overapproximation of the set of states we want to compute, we check whether we
+                // need to include some of the states that are terminal state candidates or whether the terminal states
+                // without successors are sufficient.
+                
+            }
+            
             template<typename T>
             std::vector<uint_fast64_t> getDistances(storm::storage::SparseMatrix<T> const& transitionMatrix, storm::storage::BitVector const& initialStates, boost::optional<storm::storage::BitVector> const& subsystem) {
                 std::vector<uint_fast64_t> distances(transitionMatrix.getRowGroupCount());
@@ -993,6 +1052,8 @@ namespace storm {
             
             template storm::storage::BitVector getReachableStates(storm::storage::SparseMatrix<double> const& transitionMatrix, storm::storage::BitVector const& initialStates, storm::storage::BitVector const& constraintStates, storm::storage::BitVector const& targetStates, bool useStepBound, uint_fast64_t maximalSteps);
             
+            template storm::storage::BitVector getBsccCover(storm::storage::SparseMatrix<double> const& transitionMatrix);
+            
             template std::vector<uint_fast64_t> getDistances(storm::storage::SparseMatrix<double> const& transitionMatrix, storm::storage::BitVector const& initialStates, boost::optional<storm::storage::BitVector> const& subsystem);
             
             
@@ -1065,6 +1126,8 @@ namespace storm {
             
             template storm::storage::BitVector getReachableStates(storm::storage::SparseMatrix<float> const& transitionMatrix, storm::storage::BitVector const& initialStates, storm::storage::BitVector const& constraintStates, storm::storage::BitVector const& targetStates, bool useStepBound, uint_fast64_t maximalSteps);
             
+            template storm::storage::BitVector getBsccCover(storm::storage::SparseMatrix<float> const& transitionMatrix);
+            
             template std::vector<uint_fast64_t> getDistances(storm::storage::SparseMatrix<float> const& transitionMatrix, storm::storage::BitVector const& initialStates, boost::optional<storm::storage::BitVector> const& subsystem);
             
             
@@ -1119,6 +1182,8 @@ namespace storm {
             // Instantiations for storm::RationalNumber.
             template storm::storage::BitVector getReachableStates(storm::storage::SparseMatrix<storm::RationalNumber> const& transitionMatrix, storm::storage::BitVector const& initialStates, storm::storage::BitVector const& constraintStates, storm::storage::BitVector const& targetStates, bool useStepBound, uint_fast64_t maximalSteps);
             
+            template storm::storage::BitVector getBsccCover(storm::storage::SparseMatrix<storm::RationalNumber> const& transitionMatrix);
+            
             template std::vector<uint_fast64_t> getDistances(storm::storage::SparseMatrix<storm::RationalNumber> const& transitionMatrix, storm::storage::BitVector const& initialStates, boost::optional<storm::storage::BitVector> const& subsystem);
             
             
@@ -1174,6 +1239,8 @@ namespace storm {
 #ifdef STORM_HAVE_CARL
             template storm::storage::BitVector getReachableStates(storm::storage::SparseMatrix<storm::RationalFunction> const& transitionMatrix, storm::storage::BitVector const& initialStates, storm::storage::BitVector const& constraintStates, storm::storage::BitVector const& targetStates, bool useStepBound, uint_fast64_t maximalSteps);
             
+            template storm::storage::BitVector getBsccCover(storm::storage::SparseMatrix<storm::RationalFunction> const& transitionMatrix);
+            
             template std::vector<uint_fast64_t> getDistances(storm::storage::SparseMatrix<storm::RationalFunction> const& transitionMatrix, storm::storage::BitVector const& initialStates, boost::optional<storm::storage::BitVector> const& subsystem);
             
             

src/utility/graph.h

@@ -61,6 +61,15 @@ namespace storm {
             template<typename T>
             storm::storage::BitVector getReachableStates(storm::storage::SparseMatrix<T> const& transitionMatrix, storm::storage::BitVector const& initialStates, storm::storage::BitVector const& constraintStates, storm::storage::BitVector const& targetStates, bool useStepBound = false, uint_fast64_t maximalSteps = 0);
 
+            /*!
+             * Retrieves a set of states that covers als BSCCs of the system in the sense that for every BSCC exactly
+             * one state is included in the cover.
+             *
+             * @param transitionMatrix The transition relation of the graph structure.
+             */
+            template<typename T>
+            storm::storage::BitVector getBsccCover(storm::storage::SparseMatrix<T> const& transitionMatrix);
+            
             /*!
              * Performs a breadth-first search through the underlying graph structure to compute the distance from all
              * states to the starting states of the search.

src/utility/solver.cpp

@@ -47,18 +47,43 @@ namespace storm {
             }
             
             template<typename ValueType>
-            std::unique_ptr<storm::solver::LinearEquationSolver<ValueType>> LinearEquationSolverFactory<ValueType>::create(storm::storage::SparseMatrix<ValueType> const& matrix) const {
+            std::unique_ptr<storm::solver::LinearEquationSolver<ValueType>> LinearEquationSolverFactory<ValueType>::create(storm::storage::SparseMatrix<ValueType>&& matrix) const {
+                return create(matrix);
+            }
+            
+            template<typename ValueType>
+            std::unique_ptr<storm::solver::LinearEquationSolver<ValueType>> GeneralLinearEquationSolverFactory<ValueType>::create(storm::storage::SparseMatrix<ValueType> const& matrix) const {
+                return selectSolver(matrix);
+            }
+            
+            template<typename ValueType>
+            std::unique_ptr<storm::solver::LinearEquationSolver<ValueType>> GeneralLinearEquationSolverFactory<ValueType>::create(storm::storage::SparseMatrix<ValueType>&& matrix) const {
+                return selectSolver(std::move(matrix));
+            }
+            
+            template<typename ValueType>
+            template<typename MatrixType>
+            std::unique_ptr<storm::solver::LinearEquationSolver<ValueType>> GeneralLinearEquationSolverFactory<ValueType>::selectSolver(MatrixType&& matrix) const {
                 storm::solver::EquationSolverType equationSolver = storm::settings::getModule<storm::settings::modules::MarkovChainSettings>().getEquationSolver();
                 switch (equationSolver) {
-                    case storm::solver::EquationSolverType::Gmmxx: return std::make_unique<storm::solver::GmmxxLinearEquationSolver<ValueType>>(matrix);
-                    case storm::solver::EquationSolverType::Native: return std::make_unique<storm::solver::NativeLinearEquationSolver<ValueType>>(matrix);
-                    case storm::solver::EquationSolverType::Eigen: return std::make_unique<storm::solver::EigenLinearEquationSolver<ValueType>>(matrix);
-                    case storm::solver::EquationSolverType::Elimination: return std::make_unique<storm::solver::EliminationLinearEquationSolver<ValueType>>(matrix);
-                    default: return std::make_unique<storm::solver::GmmxxLinearEquationSolver<ValueType>>(matrix);
+                    case storm::solver::EquationSolverType::Gmmxx: return std::make_unique<storm::solver::GmmxxLinearEquationSolver<ValueType>>(std::forward<MatrixType>(matrix));
+                    case storm::solver::EquationSolverType::Native: return std::make_unique<storm::solver::NativeLinearEquationSolver<ValueType>>(std::forward<MatrixType>(matrix));
+                    case storm::solver::EquationSolverType::Eigen: return std::make_unique<storm::solver::EigenLinearEquationSolver<ValueType>>(std::forward<MatrixType>(matrix));
+                    case storm::solver::EquationSolverType::Elimination: return std::make_unique<storm::solver::EliminationLinearEquationSolver<ValueType>>(std::forward<MatrixType>(matrix));
+                    default: return std::make_unique<storm::solver::GmmxxLinearEquationSolver<ValueType>>(std::forward<MatrixType>(matrix));
                 }
             }
 
-            std::unique_ptr<storm::solver::LinearEquationSolver<storm::RationalNumber>> LinearEquationSolverFactory<storm::RationalNumber>::create(storm::storage::SparseMatrix<storm::RationalNumber> const& matrix) const {
+            std::unique_ptr<storm::solver::LinearEquationSolver<storm::RationalNumber>> GeneralLinearEquationSolverFactory<storm::RationalNumber>::create(storm::storage::SparseMatrix<storm::RationalNumber> const& matrix) const {
+                return selectSolver(matrix);
+            }
+
+            std::unique_ptr<storm::solver::LinearEquationSolver<storm::RationalNumber>> GeneralLinearEquationSolverFactory<storm::RationalNumber>::create(storm::storage::SparseMatrix<storm::RationalNumber>&& matrix) const {
+                return selectSolver(std::move(matrix));
+            }
+            
+            template<typename MatrixType>
+            std::unique_ptr<storm::solver::LinearEquationSolver<storm::RationalNumber>> GeneralLinearEquationSolverFactory<storm::RationalNumber>::selectSolver(MatrixType&& matrix) const {
                 storm::solver::EquationSolverType equationSolver = storm::settings::getModule<storm::settings::modules::MarkovChainSettings>().getEquationSolver();
                 switch (equationSolver) {
                     case storm::solver::EquationSolverType::Elimination: return std::make_unique<storm::solver::EliminationLinearEquationSolver<storm::RationalNumber>>(matrix);
@@ -66,7 +91,16 @@ namespace storm {
                 }
             }
             
-            std::unique_ptr<storm::solver::LinearEquationSolver<storm::RationalFunction>> LinearEquationSolverFactory<storm::RationalFunction>::create(storm::storage::SparseMatrix<storm::RationalFunction> const& matrix) const {
+            std::unique_ptr<storm::solver::LinearEquationSolver<storm::RationalFunction>> GeneralLinearEquationSolverFactory<storm::RationalFunction>::create(storm::storage::SparseMatrix<storm::RationalFunction> const& matrix) const {
+                return selectSolver(matrix);
+            }
+
+            std::unique_ptr<storm::solver::LinearEquationSolver<storm::RationalFunction>> GeneralLinearEquationSolverFactory<storm::RationalFunction>::create(storm::storage::SparseMatrix<storm::RationalFunction>&& matrix) const {
+                return selectSolver(std::move(matrix));
+            }
+
+            template<typename MatrixType>
+            std::unique_ptr<storm::solver::LinearEquationSolver<storm::RationalFunction>> GeneralLinearEquationSolverFactory<storm::RationalFunction>::selectSolver(MatrixType&& matrix) const {
                 storm::solver::EquationSolverType equationSolver = storm::settings::getModule<storm::settings::modules::MarkovChainSettings>().getEquationSolver();
                 switch (equationSolver) {
                     case storm::solver::EquationSolverType::Elimination: return std::make_unique<storm::solver::EliminationLinearEquationSolver<storm::RationalFunction>>(matrix);
@@ -222,6 +256,7 @@ namespace storm {
             template class SymbolicGameSolverFactory<storm::dd::DdType::CUDD, double>;
             template class SymbolicGameSolverFactory<storm::dd::DdType::Sylvan, double>;
             template class LinearEquationSolverFactory<double>;
+            template class GeneralLinearEquationSolverFactory<double>;
             template class GmmxxLinearEquationSolverFactory<double>;
             template class EigenLinearEquationSolverFactory<double>;
             template class NativeLinearEquationSolverFactory<double>;
@@ -229,9 +264,11 @@ namespace storm {
             template class GameSolverFactory<double>;
 
             template class LinearEquationSolverFactory<storm::RationalNumber>;
+            template class GeneralLinearEquationSolverFactory<storm::RationalNumber>;
             template class EigenLinearEquationSolverFactory<storm::RationalNumber>;
 
             template class LinearEquationSolverFactory<storm::RationalFunction>;
+            template class GeneralLinearEquationSolverFactory<storm::RationalFunction>;
             template class EigenLinearEquationSolverFactory<storm::RationalFunction>;
         }
     }

src/utility/solver.h

@@ -87,31 +87,49 @@ namespace storm {
                  * @param matrix The matrix that defines the equation system.
                  * @return A pointer to the newly created solver.
                  */
-                virtual std::unique_ptr<storm::solver::LinearEquationSolver<ValueType>> create(storm::storage::SparseMatrix<ValueType> const& matrix) const;
-            };
+                virtual std::unique_ptr<storm::solver::LinearEquationSolver<ValueType>> create(storm::storage::SparseMatrix<ValueType> const& matrix) const = 0;
                 
-            template<>
-            class LinearEquationSolverFactory<storm::RationalNumber> {
-            public:
                 /*!
-                 * Creates a new linear equation solver instance with the given matrix.
+                 * Creates a new linear equation solver instance with the given matrix. The caller gives up posession of the
+                 * matrix by calling this function.
                  *
                  * @param matrix The matrix that defines the equation system.
                  * @return A pointer to the newly created solver.
                  */
-                virtual std::unique_ptr<storm::solver::LinearEquationSolver<storm::RationalNumber>> create(storm::storage::SparseMatrix<storm::RationalNumber> const& matrix) const;
+                virtual std::unique_ptr<storm::solver::LinearEquationSolver<ValueType>> create(storm::storage::SparseMatrix<ValueType>&& matrix) const;
+            };
+            
+            template<typename ValueType>
+            class GeneralLinearEquationSolverFactory : public LinearEquationSolverFactory<ValueType> {
+            public:
+                virtual std::unique_ptr<storm::solver::LinearEquationSolver<ValueType>> create(storm::storage::SparseMatrix<ValueType> const& matrix) const override;
+                virtual std::unique_ptr<storm::solver::LinearEquationSolver<ValueType>> create(storm::storage::SparseMatrix<ValueType>&& matrix) const override;
+                
+            private:
+                template<typename MatrixType>
+                std::unique_ptr<storm::solver::LinearEquationSolver<ValueType>> selectSolver(MatrixType&& matrix) const;
             };
 
             template<>
-            class LinearEquationSolverFactory<storm::RationalFunction> {
+            class GeneralLinearEquationSolverFactory<storm::RationalNumber> : public LinearEquationSolverFactory<storm::RationalNumber> {
             public:
-                /*!
-                 * Creates a new linear equation solver instance with the given matrix.
-                 *
-                 * @param matrix The matrix that defines the equation system.
-                 * @return A pointer to the newly created solver.
-                 */
-                virtual std::unique_ptr<storm::solver::LinearEquationSolver<storm::RationalFunction>> create(storm::storage::SparseMatrix<storm::RationalFunction> const& matrix) const;
+                virtual std::unique_ptr<storm::solver::LinearEquationSolver<storm::RationalNumber>> create(storm::storage::SparseMatrix<storm::RationalNumber> const& matrix) const override;
+                virtual std::unique_ptr<storm::solver::LinearEquationSolver<storm::RationalNumber>> create(storm::storage::SparseMatrix<storm::RationalNumber>&& matrix) const override;
+                
+            private:
+                template<typename MatrixType>
+                std::unique_ptr<storm::solver::LinearEquationSolver<storm::RationalNumber>> selectSolver(MatrixType&& matrix) const;
+            };
+
+            template<>
+            class GeneralLinearEquationSolverFactory<storm::RationalFunction> : public LinearEquationSolverFactory<storm::RationalFunction> {
+            public:
+                virtual std::unique_ptr<storm::solver::LinearEquationSolver<storm::RationalFunction>> create(storm::storage::SparseMatrix<storm::RationalFunction> const& matrix) const override;
+                virtual std::unique_ptr<storm::solver::LinearEquationSolver<storm::RationalFunction>> create(storm::storage::SparseMatrix<storm::RationalFunction>&& matrix) const override;
+
+            private:
+                template<typename MatrixType>
+                std::unique_ptr<storm::solver::LinearEquationSolver<storm::RationalFunction>> selectSolver(MatrixType&& matrix) const;
             };
             
             template<typename ValueType>