From 15a4d4757f2abd5d4af6692cb23db5964a5b09c4 Mon Sep 17 00:00:00 2001 From: dehnert Date: Thu, 23 Jun 2016 13:11:13 +0200 Subject: [PATCH] added feature to linear equation solver factories to take posession of the matrix to forward it to the solvers Former-commit-id: ed183f182099b97173f1be729f43fb0576028f03 --- .../csl/HybridCtmcCslModelChecker.cpp | 2 +- .../csl/SparseCtmcCslModelChecker.cpp | 2 +- .../prctl/HybridDtmcPrctlModelChecker.cpp | 2 +- .../prctl/SparseDtmcPrctlModelChecker.cpp | 2 +- .../EliminationLinearEquationSolver.cpp | 17 +++-- src/solver/NativeLinearEquationSolver.h | 2 +- src/utility/graph.cpp | 75 ++++++++++++++++++- src/utility/graph.h | 9 +++ src/utility/solver.cpp | 53 +++++++++++-- src/utility/solver.h | 52 ++++++++----- 10 files changed, 177 insertions(+), 39 deletions(-) diff --git a/src/modelchecker/csl/HybridCtmcCslModelChecker.cpp b/src/modelchecker/csl/HybridCtmcCslModelChecker.cpp index 81bbbd00e..b4c59bb0f 100644 --- a/src/modelchecker/csl/HybridCtmcCslModelChecker.cpp +++ b/src/modelchecker/csl/HybridCtmcCslModelChecker.cpp @@ -15,7 +15,7 @@ namespace storm { namespace modelchecker { template - HybridCtmcCslModelChecker::HybridCtmcCslModelChecker(storm::models::symbolic::Ctmc const& model) : SymbolicPropositionalModelChecker(model), linearEquationSolverFactory(new storm::utility::solver::LinearEquationSolverFactory()) { + HybridCtmcCslModelChecker::HybridCtmcCslModelChecker(storm::models::symbolic::Ctmc const& model) : SymbolicPropositionalModelChecker(model), linearEquationSolverFactory(std::make_unique>()) { // Intentionally left empty. } diff --git a/src/modelchecker/csl/SparseCtmcCslModelChecker.cpp b/src/modelchecker/csl/SparseCtmcCslModelChecker.cpp index bf70ca586..abf21fe01 100644 --- a/src/modelchecker/csl/SparseCtmcCslModelChecker.cpp +++ b/src/modelchecker/csl/SparseCtmcCslModelChecker.cpp @@ -24,7 +24,7 @@ namespace storm { namespace modelchecker { template - SparseCtmcCslModelChecker::SparseCtmcCslModelChecker(SparseCtmcModelType const& model) : SparsePropositionalModelChecker(model), linearEquationSolverFactory(new storm::utility::solver::LinearEquationSolverFactory()) { + SparseCtmcCslModelChecker::SparseCtmcCslModelChecker(SparseCtmcModelType const& model) : SparsePropositionalModelChecker(model), linearEquationSolverFactory(std::make_unique>()) { // Intentionally left empty. } diff --git a/src/modelchecker/prctl/HybridDtmcPrctlModelChecker.cpp b/src/modelchecker/prctl/HybridDtmcPrctlModelChecker.cpp index 4f8c38285..75f651e8a 100644 --- a/src/modelchecker/prctl/HybridDtmcPrctlModelChecker.cpp +++ b/src/modelchecker/prctl/HybridDtmcPrctlModelChecker.cpp @@ -31,7 +31,7 @@ namespace storm { } template - HybridDtmcPrctlModelChecker::HybridDtmcPrctlModelChecker(storm::models::symbolic::Dtmc const& model) : SymbolicPropositionalModelChecker(model), linearEquationSolverFactory(new storm::utility::solver::LinearEquationSolverFactory()) { + HybridDtmcPrctlModelChecker::HybridDtmcPrctlModelChecker(storm::models::symbolic::Dtmc const& model) : SymbolicPropositionalModelChecker(model), linearEquationSolverFactory(std::make_unique>()) { // Intentionally left empty. } diff --git a/src/modelchecker/prctl/SparseDtmcPrctlModelChecker.cpp b/src/modelchecker/prctl/SparseDtmcPrctlModelChecker.cpp index ee9195b6e..75fdaca64 100644 --- a/src/modelchecker/prctl/SparseDtmcPrctlModelChecker.cpp +++ b/src/modelchecker/prctl/SparseDtmcPrctlModelChecker.cpp @@ -29,7 +29,7 @@ namespace storm { } template - SparseDtmcPrctlModelChecker::SparseDtmcPrctlModelChecker(SparseDtmcModelType const& model) : SparsePropositionalModelChecker(model), linearEquationSolverFactory(new storm::utility::solver::LinearEquationSolverFactory()) { + SparseDtmcPrctlModelChecker::SparseDtmcPrctlModelChecker(SparseDtmcModelType const& model) : SparsePropositionalModelChecker(model), linearEquationSolverFactory(std::make_unique>()) { // Intentionally left empty. } diff --git a/src/solver/EliminationLinearEquationSolver.cpp b/src/solver/EliminationLinearEquationSolver.cpp index 4a4321ef9..52a355fbd 100644 --- a/src/solver/EliminationLinearEquationSolver.cpp +++ b/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::solveEquationSystem(std::vector& x, std::vector const& b, std::vector* 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 transitionMatrix(A); + transitionMatrix.convertToEquationSystem(); + storm::storage::SparseMatrix 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().getEliminationOrder(); boost::optional> 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 priorityQueue = createStatePriorityQueue(distanceBasedPriorities, flexibleMatrix, flexibleBackwardTransitions, b, storm::storage::BitVector(x.size(), true)); diff --git a/src/solver/NativeLinearEquationSolver.h b/src/solver/NativeLinearEquationSolver.h index 187bc8b5c..9dabd4e65 100644 --- a/src/solver/NativeLinearEquationSolver.h +++ b/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 class NativeLinearEquationSolver : public LinearEquationSolver { diff --git a/src/utility/graph.cpp b/src/utility/graph.cpp index 5705c10b9..58ff2efd8 100644 --- a/src/utility/graph.cpp +++ b/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 + storm::storage::BitVector getBsccCover(storm::storage::SparseMatrix const& transitionMatrix) { + storm::storage::BitVector result(transitionMatrix.getRowGroupCount()); + storm::storage::StronglyConnectedComponentDecomposition decomposition(transitionMatrix, false, true); + + // Take the first state out of each BSCC. + for (auto const& scc : decomposition) { + result.set(*scc.begin()); + } + + return result; + } + + template + storm::storage::BitVector getTerminalStateCover(storm::storage::SparseMatrix const& transitionMatrix, storm::storage::BitVector const& initialStates) { + storm::storage::BitVector terminalStateCandidates(transitionMatrix.getRowGroupCount()); + storm::storage::BitVector terminalStatesWithoutSuccessors(transitionMatrix.getRowCount()); + + std::queue 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 std::vector getDistances(storm::storage::SparseMatrix const& transitionMatrix, storm::storage::BitVector const& initialStates, boost::optional const& subsystem) { std::vector distances(transitionMatrix.getRowGroupCount()); @@ -993,6 +1052,8 @@ namespace storm { template storm::storage::BitVector getReachableStates(storm::storage::SparseMatrix 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 const& transitionMatrix); + template std::vector getDistances(storm::storage::SparseMatrix const& transitionMatrix, storm::storage::BitVector const& initialStates, boost::optional const& subsystem); @@ -1065,6 +1126,8 @@ namespace storm { template storm::storage::BitVector getReachableStates(storm::storage::SparseMatrix 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 const& transitionMatrix); + template std::vector getDistances(storm::storage::SparseMatrix const& transitionMatrix, storm::storage::BitVector const& initialStates, boost::optional const& subsystem); @@ -1119,6 +1182,8 @@ namespace storm { // Instantiations for storm::RationalNumber. template storm::storage::BitVector getReachableStates(storm::storage::SparseMatrix 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 const& transitionMatrix); + template std::vector getDistances(storm::storage::SparseMatrix const& transitionMatrix, storm::storage::BitVector const& initialStates, boost::optional const& subsystem); @@ -1174,6 +1239,8 @@ namespace storm { #ifdef STORM_HAVE_CARL template storm::storage::BitVector getReachableStates(storm::storage::SparseMatrix 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 const& transitionMatrix); + template std::vector getDistances(storm::storage::SparseMatrix const& transitionMatrix, storm::storage::BitVector const& initialStates, boost::optional const& subsystem); diff --git a/src/utility/graph.h b/src/utility/graph.h index 26f73902d..10dcd3e97 100644 --- a/src/utility/graph.h +++ b/src/utility/graph.h @@ -60,6 +60,15 @@ namespace storm { */ template storm::storage::BitVector getReachableStates(storm::storage::SparseMatrix 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 + storm::storage::BitVector getBsccCover(storm::storage::SparseMatrix const& transitionMatrix); /*! * Performs a breadth-first search through the underlying graph structure to compute the distance from all diff --git a/src/utility/solver.cpp b/src/utility/solver.cpp index 237d22b48..cae5bdb98 100644 --- a/src/utility/solver.cpp +++ b/src/utility/solver.cpp @@ -47,18 +47,43 @@ namespace storm { } template - std::unique_ptr> LinearEquationSolverFactory::create(storm::storage::SparseMatrix const& matrix) const { + std::unique_ptr> LinearEquationSolverFactory::create(storm::storage::SparseMatrix&& matrix) const { + return create(matrix); + } + + template + std::unique_ptr> GeneralLinearEquationSolverFactory::create(storm::storage::SparseMatrix const& matrix) const { + return selectSolver(matrix); + } + + template + std::unique_ptr> GeneralLinearEquationSolverFactory::create(storm::storage::SparseMatrix&& matrix) const { + return selectSolver(std::move(matrix)); + } + + template + template + std::unique_ptr> GeneralLinearEquationSolverFactory::selectSolver(MatrixType&& matrix) const { storm::solver::EquationSolverType equationSolver = storm::settings::getModule().getEquationSolver(); switch (equationSolver) { - case storm::solver::EquationSolverType::Gmmxx: return std::make_unique>(matrix); - case storm::solver::EquationSolverType::Native: return std::make_unique>(matrix); - case storm::solver::EquationSolverType::Eigen: return std::make_unique>(matrix); - case storm::solver::EquationSolverType::Elimination: return std::make_unique>(matrix); - default: return std::make_unique>(matrix); + case storm::solver::EquationSolverType::Gmmxx: return std::make_unique>(std::forward(matrix)); + case storm::solver::EquationSolverType::Native: return std::make_unique>(std::forward(matrix)); + case storm::solver::EquationSolverType::Eigen: return std::make_unique>(std::forward(matrix)); + case storm::solver::EquationSolverType::Elimination: return std::make_unique>(std::forward(matrix)); + default: return std::make_unique>(std::forward(matrix)); } } - std::unique_ptr> LinearEquationSolverFactory::create(storm::storage::SparseMatrix const& matrix) const { + std::unique_ptr> GeneralLinearEquationSolverFactory::create(storm::storage::SparseMatrix const& matrix) const { + return selectSolver(matrix); + } + + std::unique_ptr> GeneralLinearEquationSolverFactory::create(storm::storage::SparseMatrix&& matrix) const { + return selectSolver(std::move(matrix)); + } + + template + std::unique_ptr> GeneralLinearEquationSolverFactory::selectSolver(MatrixType&& matrix) const { storm::solver::EquationSolverType equationSolver = storm::settings::getModule().getEquationSolver(); switch (equationSolver) { case storm::solver::EquationSolverType::Elimination: return std::make_unique>(matrix); @@ -66,7 +91,16 @@ namespace storm { } } - std::unique_ptr> LinearEquationSolverFactory::create(storm::storage::SparseMatrix const& matrix) const { + std::unique_ptr> GeneralLinearEquationSolverFactory::create(storm::storage::SparseMatrix const& matrix) const { + return selectSolver(matrix); + } + + std::unique_ptr> GeneralLinearEquationSolverFactory::create(storm::storage::SparseMatrix&& matrix) const { + return selectSolver(std::move(matrix)); + } + + template + std::unique_ptr> GeneralLinearEquationSolverFactory::selectSolver(MatrixType&& matrix) const { storm::solver::EquationSolverType equationSolver = storm::settings::getModule().getEquationSolver(); switch (equationSolver) { case storm::solver::EquationSolverType::Elimination: return std::make_unique>(matrix); @@ -222,6 +256,7 @@ namespace storm { template class SymbolicGameSolverFactory; template class SymbolicGameSolverFactory; template class LinearEquationSolverFactory; + template class GeneralLinearEquationSolverFactory; template class GmmxxLinearEquationSolverFactory; template class EigenLinearEquationSolverFactory; template class NativeLinearEquationSolverFactory; @@ -229,9 +264,11 @@ namespace storm { template class GameSolverFactory; template class LinearEquationSolverFactory; + template class GeneralLinearEquationSolverFactory; template class EigenLinearEquationSolverFactory; template class LinearEquationSolverFactory; + template class GeneralLinearEquationSolverFactory; template class EigenLinearEquationSolverFactory; } } diff --git a/src/utility/solver.h b/src/utility/solver.h index d61c99442..bacc61e5e 100644 --- a/src/utility/solver.h +++ b/src/utility/solver.h @@ -77,7 +77,7 @@ namespace storm { public: virtual std::unique_ptr> create(storm::dd::Add const& A, storm::dd::Bdd const& allRows, std::set const& rowMetaVariables, std::set const& columnMetaVariables, std::vector> const& rowColumnMetaVariablePairs, std::set const& player1Variables, std::set const& player2Variables) const; }; - + template class LinearEquationSolverFactory { public: @@ -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> create(storm::storage::SparseMatrix const& matrix) const; - }; - - template<> - class LinearEquationSolverFactory { - public: + virtual std::unique_ptr> create(storm::storage::SparseMatrix const& matrix) const = 0; + /*! - * 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> create(storm::storage::SparseMatrix const& matrix) const; + virtual std::unique_ptr> create(storm::storage::SparseMatrix&& matrix) const; + }; + + template + class GeneralLinearEquationSolverFactory : public LinearEquationSolverFactory { + public: + virtual std::unique_ptr> create(storm::storage::SparseMatrix const& matrix) const override; + virtual std::unique_ptr> create(storm::storage::SparseMatrix&& matrix) const override; + + private: + template + std::unique_ptr> selectSolver(MatrixType&& matrix) const; }; template<> - class LinearEquationSolverFactory { + class GeneralLinearEquationSolverFactory : public LinearEquationSolverFactory { 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> create(storm::storage::SparseMatrix const& matrix) const; + virtual std::unique_ptr> create(storm::storage::SparseMatrix const& matrix) const override; + virtual std::unique_ptr> create(storm::storage::SparseMatrix&& matrix) const override; + + private: + template + std::unique_ptr> selectSolver(MatrixType&& matrix) const; + }; + + template<> + class GeneralLinearEquationSolverFactory : public LinearEquationSolverFactory { + public: + virtual std::unique_ptr> create(storm::storage::SparseMatrix const& matrix) const override; + virtual std::unique_ptr> create(storm::storage::SparseMatrix&& matrix) const override; + + private: + template + std::unique_ptr> selectSolver(MatrixType&& matrix) const; }; template