From 98b0bcf1872b378e689b2b035a1230953752e368 Mon Sep 17 00:00:00 2001 From: PBerger Date: Mon, 24 Feb 2014 04:43:07 +0100 Subject: [PATCH] Reimplemented the TopologicalValueIterationNondeterministicLinearEquationSolver with splitting into submatrices. Added a dtmc example for tests with the StronglyConnectedComponentDecomposition. Former-commit-id: 0c33793fe6ca844ac775cd29d050846626526a4d --- examples/mdp/scc/scc.pctl | 4 +- ...onNondeterministicLinearEquationSolver.cpp | 124 ++++++++++++++---- ...tronglyConnectedComponentDecomposition.cpp | 2 +- src/storm.cpp | 18 ++- src/utility/OsDetection.h | 5 + ...ValueIterationMdpPrctlModelCheckerTest.cpp | 32 ++++- ...glyConnectedComponentDecompositionTest.cpp | 28 ++++ 7 files changed, 179 insertions(+), 34 deletions(-) diff --git a/examples/mdp/scc/scc.pctl b/examples/mdp/scc/scc.pctl index 8a853a969..501655389 100644 --- a/examples/mdp/scc/scc.pctl +++ b/examples/mdp/scc/scc.pctl @@ -1,2 +1,2 @@ -Pmin=? [ F end ] -Pmax=? [ F end ] \ No newline at end of file +Pmin=? [ (!statetwo) U end ] +Pmax=? [ (!statetwo) U end ] \ No newline at end of file diff --git a/src/solver/TopologicalValueIterationNondeterministicLinearEquationSolver.cpp b/src/solver/TopologicalValueIterationNondeterministicLinearEquationSolver.cpp index 1af63b729..c82a40af5 100644 --- a/src/solver/TopologicalValueIterationNondeterministicLinearEquationSolver.cpp +++ b/src/solver/TopologicalValueIterationNondeterministicLinearEquationSolver.cpp @@ -7,6 +7,11 @@ #include "src/utility/graph.h" #include "src/models/PseudoModel.h" #include "src/storage/StronglyConnectedComponentDecomposition.h" +#include "src/exceptions/IllegalArgumentException.h" + +#include "log4cplus/logger.h" +#include "log4cplus/loggingmacros.h" +extern log4cplus::Logger logger; namespace storm { namespace solver { @@ -39,11 +44,29 @@ namespace storm { //std::vector> stronglyConnectedComponents = storm::utility::graph::performSccDecomposition(this->getModel(), stronglyConnectedComponents, stronglyConnectedComponentsDependencyGraph); //storm::storage::SparseMatrix stronglyConnectedComponentsDependencyGraph = this->getModel().extractSccDependencyGraph(stronglyConnectedComponents); + std::cout << "TopoSolver Input Matrix: " << A.getRowCount() << " x " << A.getColumnCount() << " with " << A.getEntryCount() << " Entries:" << std::endl; + + uint_fast64_t const rowCount = A.getRowCount(); + for (uint_fast64_t row = 0; row < rowCount; ++row) { + std::cout << "Row " << row << ": "; + auto const& rowElement = A.getRow(row); + for (auto rowIt = rowElement.begin(); rowIt != rowElement.end(); ++rowIt) { + std::cout << rowIt->first << " [" << rowIt->second << "], "; + } + std::cout << std::endl; + } + + storm::models::NonDeterministicMatrixBasedPseudoModel pseudoModel(A, nondeterministicChoiceIndices); - storm::storage::StronglyConnectedComponentDecomposition sccDecomposition(*static_cast*>(&pseudoModel), false, false); - storm::storage::SparseMatrix stronglyConnectedComponentsDependencyGraph = pseudoModel.extractPartitionDependencyGraph(sccDecomposition); - + //storm::storage::StronglyConnectedComponentDecomposition sccDecomposition(*static_cast*>(&pseudoModel), false, false); + storm::storage::StronglyConnectedComponentDecomposition sccDecomposition(pseudoModel, false, false); + + if (sccDecomposition.size() == 0) { + LOG4CPLUS_ERROR(logger, "Can not solve given Equation System as the SCC Decomposition returned no SCCs."); + throw storm::exceptions::IllegalArgumentException() << "Can not solve given Equation System as the SCC Decomposition returned no SCCs."; + } + storm::storage::SparseMatrix stronglyConnectedComponentsDependencyGraph = pseudoModel.extractPartitionDependencyGraph(sccDecomposition); std::vector topologicalSort = storm::utility::graph::getTopologicalSort(stronglyConnectedComponentsDependencyGraph); // Set up the environment for the power method. @@ -52,12 +75,12 @@ namespace storm { multiplyResult = new std::vector(A.getRowCount()); multiplyResultMemoryProvided = false; } - std::vector* currentX = &x; - bool xMemoryProvided = true; - if (newX == nullptr) { - newX = new std::vector(x.size()); - xMemoryProvided = false; - } + std::vector* currentX = nullptr; + //bool xMemoryProvided = true; + //if (newX == nullptr) { + // newX = new std::vector(x.size()); + // xMemoryProvided = false; + //} std::vector* swap = nullptr; uint_fast64_t currentMaxLocalIterations = 0; uint_fast64_t localIterations = 0; @@ -68,20 +91,62 @@ namespace storm { // solved after all SCCs it depends on have been solved. int counter = 0; std::cout << "Solving Equation System using the TopologicalValueIterationNon..." << std::endl; + std::cout << "Found " << sccDecomposition.size() << " SCCs." << std::endl; + for (auto sccIndexIt = topologicalSort.begin(); sccIndexIt != topologicalSort.end() && converged; ++sccIndexIt) { storm::storage::StateBlock const& scc = sccDecomposition[*sccIndexIt]; - std::cout << "SCC " << counter << " contains:" << std::endl; - for (auto sccIt = scc.cbegin(); sccIt != scc.cend(); ++sccIt) { - std::cout << *sccIt << ", "; + std::cout << "SCC " << counter << " from Index " << *sccIndexIt << " contains:" << std::endl; + ++counter; + for (auto state : scc) { + std::cout << state << ", "; } std::cout << std::endl; + // Generate a submatrix + storm::storage::BitVector subMatrixIndices(rowCount, scc.cbegin(), scc.cend()); + storm::storage::SparseMatrix sccSubmatrix = A.getSubmatrix(subMatrixIndices, nondeterministicChoiceIndices, true); + std::vector sccSubB(sccSubmatrix.getRowCount()); + storm::utility::vector::selectVectorValues(sccSubB, subMatrixIndices, nondeterministicChoiceIndices, b); + std::vector sccSubX(sccSubmatrix.getColumnCount()); + std::vector sccSubXSwap(sccSubmatrix.getColumnCount()); + + // Prepare the pointers for swapping in the calculation + currentX = &sccSubX; + swap = &sccSubXSwap; + + storm::utility::vector::selectVectorValues(sccSubX, subMatrixIndices, x); // x is getCols() large, where as b and multiplyResult are getRows() (nondet. choices times states) + std::vector sccSubNondeterministicChoiceIndices(sccSubmatrix.getColumnCount() + 1); + sccSubNondeterministicChoiceIndices.at(0) = 0; + + // Preprocess all dependant states + // Remove outgoing transitions and create the ChoiceIndices + uint_fast64_t innerIndex = 0; + for (uint_fast64_t state: scc) { + // Choice Indices + sccSubNondeterministicChoiceIndices.at(innerIndex + 1) = sccSubNondeterministicChoiceIndices.at(innerIndex) + (nondeterministicChoiceIndices[state + 1] - nondeterministicChoiceIndices[state]); + + for (auto rowGroupIt = nondeterministicChoiceIndices[state]; rowGroupIt != nondeterministicChoiceIndices[state + 1]; ++rowGroupIt) { + storm::storage::SparseMatrix::const_rows row = A.getRow(rowGroupIt); + for (auto rowIt = row.begin(); rowIt != row.end(); ++rowIt) { + if (!subMatrixIndices.get(rowIt->first)) { + // This is an outgoing transition of a state in the SCC to a state not included in the SCC + // Subtracting Pr(tau) * x_other from b fixes that + sccSubB.at(innerIndex) = sccSubB.at(innerIndex) - (rowIt->second * x.at(rowIt->first)); + } + } + } + ++innerIndex; + } + // For the current SCC, we need to perform value iteration until convergence. localIterations = 0; converged = false; while (!converged && localIterations < this->maximalNumberOfIterations) { // Compute x' = A*x + b. + sccSubmatrix.multiplyWithVector(*currentX, *multiplyResult); + storm::utility::vector::addVectorsInPlace(*multiplyResult, sccSubB); + //A.multiplyWithVector(scc, nondeterministicChoiceIndices, *currentX, multiplyResult); //storm::utility::addVectors(scc, nondeterministicChoiceIndices, multiplyResult, b); @@ -93,43 +158,46 @@ namespace storm { // Reduce the vector x' by applying min/max for all non-deterministic choices. if (minimize) { - //storm::utility::reduceVectorMin(*multiplyResult, *newX, scc, nondeterministicChoiceIndices); + storm::utility::vector::reduceVectorMin(*multiplyResult, *swap, sccSubNondeterministicChoiceIndices); } else { - //storm::utility::reduceVectorMax(*multiplyResult, *newX, scc, nondeterministicChoiceIndices); + storm::utility::vector::reduceVectorMax(*multiplyResult, *swap, sccSubNondeterministicChoiceIndices); } // Determine whether the method converged. // TODO: It seems that the equalModuloPrecision call that compares all values should have a higher // running time. In fact, it is faster. This has to be investigated. // converged = storm::utility::equalModuloPrecision(*currentX, *newX, scc, precision, relative); - converged = storm::utility::vector::equalModuloPrecision(*currentX, *newX, this->precision, this->relative); + converged = storm::utility::vector::equalModuloPrecision(*currentX, *swap, this->precision, this->relative); // Update environment variables. - swap = currentX; - currentX = newX; - newX = swap; + std::swap(currentX, swap); + ++localIterations; ++globalIterations; } + // The Result of this SCC has to be taken back into the main result vector + innerIndex = 0; + for (uint_fast64_t state: scc) { + x.at(state) = currentX->at(innerIndex); + ++innerIndex; + } + + // Since the pointers for swapping in the calculation point to temps they should not be valide anymore + currentX = nullptr; + swap = nullptr; + // As the "number of iterations" of the full method is the maximum of the local iterations, we need to keep // track of the maximum. if (localIterations > currentMaxLocalIterations) { currentMaxLocalIterations = localIterations; } } - - // If we performed an odd number of global iterations, we need to swap the x and currentX, because the newest - // result is currently stored in currentX, but x is the output vector. - // TODO: Check whether this is correct or should be put into the for-loop over SCCs. - if (globalIterations % 2 == 1) { - std::swap(x, *currentX); - } - if (!xMemoryProvided) { - delete newX; - } + //if (!xMemoryProvided) { + // delete newX; + //} if (!multiplyResultMemoryProvided) { delete multiplyResult; diff --git a/src/storage/StronglyConnectedComponentDecomposition.cpp b/src/storage/StronglyConnectedComponentDecomposition.cpp index b46a5706f..b1100f369 100644 --- a/src/storage/StronglyConnectedComponentDecomposition.cpp +++ b/src/storage/StronglyConnectedComponentDecomposition.cpp @@ -35,7 +35,7 @@ namespace storm { template StronglyConnectedComponentDecomposition::StronglyConnectedComponentDecomposition(storm::models::AbstractModel const& model, storm::storage::BitVector const& subsystem, bool dropNaiveSccs, bool onlyBottomSccs) { storm::models::ModelBasedPseudoModel encapsulatedModel(model); - performSccDecomposition(*static_cast*>(&encapsulatedModel), subsystem, dropNaiveSccs, onlyBottomSccs); + performSccDecomposition(encapsulatedModel, subsystem, dropNaiveSccs, onlyBottomSccs); } template diff --git a/src/storm.cpp b/src/storm.cpp index 9149b2c99..fadcb3db1 100644 --- a/src/storm.cpp +++ b/src/storm.cpp @@ -17,6 +17,7 @@ #include #include #include +#include #include #include @@ -29,6 +30,7 @@ #include "src/models/AtomicPropositionsLabeling.h" #include "src/modelchecker/prctl/SparseDtmcPrctlModelChecker.h" #include "src/modelchecker/prctl/SparseMdpPrctlModelChecker.h" +#include "src/modelchecker/prctl/TopologicalValueIterationMdpPrctlModelChecker.h" #include "src/solver/GmmxxLinearEquationSolver.h" #include "src/solver/GmmxxNondeterministicLinearEquationSolver.h" #include "src/solver/GurobiLpSolver.h" @@ -132,6 +134,16 @@ void setUpFileLogging() { logger.addAppender(fileLogAppender); } +/*! +* Gives the current working directory +* +* @return std::string The path of the current working directory +*/ +std::string getCurrentWorkingDirectory() { + char temp[512]; + return (getcwd(temp, 512 - 1) ? std::string(temp) : std::string("")); +} + /*! * Prints the header. */ @@ -146,7 +158,8 @@ void printHeader(const int argc, const char* argv[]) { for (int i = 0; i < argc; ++i) { commandStream << argv[i] << " "; } - std::cout << "Command line: " << commandStream.str() << std::endl << std::endl; + std::cout << "Command line: " << commandStream.str() << std::endl; + std::cout << "Current working directory: " << getCurrentWorkingDirectory() << std::endl << std::endl; } /*! @@ -234,7 +247,8 @@ storm::modelchecker::prctl::AbstractModelChecker* createPrctlModelChecke */ storm::modelchecker::prctl::AbstractModelChecker* createPrctlModelChecker(storm::models::Mdp& mdp) { // Create the appropriate model checker. - return new storm::modelchecker::prctl::SparseMdpPrctlModelChecker(mdp); + //return new storm::modelchecker::prctl::SparseMdpPrctlModelChecker(mdp); + return new storm::modelchecker::prctl::TopologicalValueIterationMdpPrctlModelChecker(mdp); } /*! diff --git a/src/utility/OsDetection.h b/src/utility/OsDetection.h index 1969a6de7..37be3f0e7 100644 --- a/src/utility/OsDetection.h +++ b/src/utility/OsDetection.h @@ -4,10 +4,12 @@ #if defined __linux__ || defined __linux # define LINUX # include +# include #include // Required by ErrorHandling.h #include // Required by ErrorHandling.h #include // Required by storm.cpp, Memory Usage #include // Required by storm.cpp, Memory Usage +# define GetCurrentDir getcwd #elif defined TARGET_OS_MAC || defined __apple__ || defined __APPLE__ # define MACOSX # define _DARWIN_USE_64_BIT_INODE @@ -17,6 +19,7 @@ # include // Required by ErrorHandling.h # include // Required by storm.cpp, Memory Usage # include // Required by storm.cpp, Memory Usage +# define GetCurrentDir getcwd #elif defined _WIN32 || defined _WIN64 # define WINDOWS # ifndef NOMINMAX @@ -28,8 +31,10 @@ # include # include # include +# include # define strncpy strncpy_s # define sscanf sscanf_s +# define GetCurrentDir _getcwd // This disables Warning C4250 - Diamond Inheritance Dominance #pragma warning(disable:4250) diff --git a/test/functional/modelchecker/TopologicalValueIterationMdpPrctlModelCheckerTest.cpp b/test/functional/modelchecker/TopologicalValueIterationMdpPrctlModelCheckerTest.cpp index 874170a4f..c2956f779 100644 --- a/test/functional/modelchecker/TopologicalValueIterationMdpPrctlModelCheckerTest.cpp +++ b/test/functional/modelchecker/TopologicalValueIterationMdpPrctlModelCheckerTest.cpp @@ -7,6 +7,36 @@ #include "src/modelchecker/prctl/TopologicalValueIterationMdpPrctlModelChecker.h" #include "src/parser/AutoParser.h" +TEST(TopologicalValueIterationMdpPrctlModelCheckerTest, SmallLinEqSystem) { + storm::storage::SparseMatrixBuilder matrixBuilder(3, 3); + ASSERT_NO_THROW(matrixBuilder.addNextValue(0, 2, 1.0)); + ASSERT_NO_THROW(matrixBuilder.addNextValue(1, 0, 4.0)); + ASSERT_NO_THROW(matrixBuilder.addNextValue(1, 1, 7.0)); + ASSERT_NO_THROW(matrixBuilder.addNextValue(2, 1, 7.0)); + ASSERT_NO_THROW(matrixBuilder.addNextValue(2, 2, -1.0)); + + storm::storage::SparseMatrix matrix; + ASSERT_NO_THROW(matrix = matrixBuilder.build()); + + ASSERT_EQ(3, matrix.getRowCount()); + ASSERT_EQ(3, matrix.getColumnCount()); + ASSERT_EQ(5, matrix.getEntryCount()); + + // Solve the Linear Equation System + storm::solver::TopologicalValueIterationNondeterministicLinearEquationSolver topoSolver; + + std::vector x(3); + std::vector b = { 5, 8, 2 }; + std::vector choices = { 0, 1, 2, 3 }; + + ASSERT_NO_THROW(topoSolver.solveEquationSystem(true, matrix, x, b, choices)); + + storm::settings::Settings* s = storm::settings::Settings::getInstance(); + ASSERT_LT(std::abs(x.at(0) - 0.25), s->getOptionByLongName("precision").getArgument(0).getValueAsDouble()); + ASSERT_LT(std::abs(x.at(1) - 1.0), s->getOptionByLongName("precision").getArgument(0).getValueAsDouble()); + ASSERT_LT(std::abs(x.at(2) - 5.0), s->getOptionByLongName("precision").getArgument(0).getValueAsDouble()); +} + TEST(TopologicalValueIterationMdpPrctlModelCheckerTest, Dice) { storm::settings::Settings* s = storm::settings::Settings::getInstance(); //storm::parser::AutoParser parser(STORM_CPP_BASE_PATH "/examples/mdp/two_dice/two_dice.tra", STORM_CPP_BASE_PATH "/examples/mdp/two_dice/two_dice.lab", "", STORM_CPP_BASE_PATH "/examples/mdp/two_dice/two_dice.flip.trans.rew"); @@ -17,7 +47,7 @@ TEST(TopologicalValueIterationMdpPrctlModelCheckerTest, Dice) { std::shared_ptr> mdp = parser.getModel>(); ASSERT_EQ(mdp->getNumberOfStates(), 11ull); - ASSERT_EQ(mdp->getNumberOfTransitions(), 17ull); + ASSERT_EQ(mdp->getNumberOfTransitions(), 18ull); storm::modelchecker::prctl::TopologicalValueIterationMdpPrctlModelChecker mc(*mdp); diff --git a/test/functional/storage/StronglyConnectedComponentDecompositionTest.cpp b/test/functional/storage/StronglyConnectedComponentDecompositionTest.cpp index 22f04d0d1..0279b72fb 100644 --- a/test/functional/storage/StronglyConnectedComponentDecompositionTest.cpp +++ b/test/functional/storage/StronglyConnectedComponentDecompositionTest.cpp @@ -44,3 +44,31 @@ TEST(StronglyConnectedComponentDecomposition, FullSystem2) { markovAutomaton = nullptr; } + +TEST(StronglyConnectedComponentDecomposition, MatrixBasedSystem) { + storm::parser::AutoParser parser(STORM_CPP_BASE_PATH "/examples/dtmc/scc/scc.tra", STORM_CPP_BASE_PATH "/examples/dtmc/scc/scc.lab", "", ""); + std::shared_ptr> dtmc = parser.getModel>(); + + storm::storage::StronglyConnectedComponentDecomposition sccDecomposition; + ASSERT_NO_THROW(sccDecomposition = storm::storage::StronglyConnectedComponentDecomposition(*dtmc, true, false)); + + ASSERT_EQ(sccDecomposition.size(), 3); + + // Now, because there is no ordering we have to check the contents of the MECs in a symmetrical way. + storm::storage::StateBlock const& scc1 = sccDecomposition[0]; + storm::storage::StateBlock const& scc2 = sccDecomposition[1]; + storm::storage::StateBlock const& scc3 = sccDecomposition[2]; + + std::vector correctScc1 = { 1, 2, 3, 4 }; + std::vector correctScc2 = { 5, 6, 7, 8 }; + std::vector correctScc3 = { 0 }; + + ASSERT_TRUE(scc1 == storm::storage::StateBlock(correctScc1.begin(), correctScc1.end()) || scc1 == storm::storage::StateBlock(correctScc2.begin(), correctScc2.end()) || scc1 == storm::storage::StateBlock(correctScc3.begin(), correctScc3.end())); + ASSERT_TRUE(scc2 == storm::storage::StateBlock(correctScc1.begin(), correctScc1.end()) || scc2 == storm::storage::StateBlock(correctScc2.begin(), correctScc2.end()) || scc2 == storm::storage::StateBlock(correctScc3.begin(), correctScc3.end())); + ASSERT_TRUE(scc3 == storm::storage::StateBlock(correctScc1.begin(), correctScc1.end()) || scc3 == storm::storage::StateBlock(correctScc2.begin(), correctScc2.end()) || scc3 == storm::storage::StateBlock(correctScc3.begin(), correctScc3.end())); + + ASSERT_NO_THROW(sccDecomposition = storm::storage::StronglyConnectedComponentDecomposition(*dtmc, true, true)); + ASSERT_EQ(2, sccDecomposition.size()); + + dtmc = nullptr; +} \ No newline at end of file