From 8ab7cae9741253618453d72b0cbbbfbeb79b2df4 Mon Sep 17 00:00:00 2001 From: TimQu Date: Mon, 21 Dec 2015 19:16:16 +0100 Subject: [PATCH] early termination Former-commit-id: a1752a6fb765c352dbd7288a317937b4c15a46cb --- .../region/ApproximationModel.cpp | 30 ++++++++++-- src/modelchecker/region/ApproximationModel.h | 4 +- src/modelchecker/region/SamplingModel.cpp | 30 +++++++++--- src/modelchecker/region/SamplingModel.h | 4 +- src/solver/AllowEarlyTerminationCondition.cpp | 25 +++++++++- src/solver/AllowEarlyTerminationCondition.h | 12 +++++ src/solver/GameSolver.cpp | 13 +++-- src/solver/GameSolver.h | 6 +++ .../GmmxxMinMaxLinearEquationSolver.cpp | 6 +++ .../NativeMinMaxLinearEquationSolver.cpp | 6 +++ src/utility/policyguessing.cpp | 48 +++---------------- 11 files changed, 123 insertions(+), 61 deletions(-) diff --git a/src/modelchecker/region/ApproximationModel.cpp b/src/modelchecker/region/ApproximationModel.cpp index 5bf1ff8a2..887ad08e2 100644 --- a/src/modelchecker/region/ApproximationModel.cpp +++ b/src/modelchecker/region/ApproximationModel.cpp @@ -12,6 +12,7 @@ #include "src/models/sparse/Mdp.h" #include "src/models/ModelType.h" #include "src/models/sparse/StandardRewardModel.h" +#include "src/solver/AllowEarlyTerminationCondition.h" #include "src/solver/MinMaxLinearEquationSolver.h" #include "src/solver/GameSolver.h" #include "src/utility/macros.h" @@ -42,7 +43,6 @@ namespace storm { } else { STORM_LOG_THROW(false, storm::exceptions::InvalidArgumentException, "Invalid formula: " << formula << ". Approximation model only supports eventually or reachability reward formulae."); } - this->solverData.player1Goal = storm::solver::SolveGoal(storm::logic::isLowerBound(formula->getComparisonType())); STORM_LOG_THROW(parametricModel.hasLabel("target"), storm::exceptions::InvalidArgumentException, "The given Model has no \"target\"-statelabel."); this->targetStates = parametricModel.getStateLabeling().getStates("target"); STORM_LOG_THROW(parametricModel.hasLabel("sink"), storm::exceptions::InvalidArgumentException, "The given Model has no \"sink\"-statelabel."); @@ -76,6 +76,15 @@ namespace storm { this->solverData.initialStateIndex = newIndices[initialState]; this->solverData.lastMinimizingPolicy = Policy(this->matrixData.matrix.getRowGroupCount(), 0); this->solverData.lastMaximizingPolicy = Policy(this->matrixData.matrix.getRowGroupCount(), 0); + //this->solverData.player1Goal = storm::solver::SolveGoal(storm::logic::isLowerBound(formula->getComparisonType())); + storm::storage::BitVector filter(this->solverData.result.size(), false); + filter.set(this->solverData.initialStateIndex, true); + this->solverData.player1Goal = std::make_unique>( + storm::logic::isLowerBound(formula->getComparisonType()) ? storm::solver::OptimizationDirection::Minimize : storm::solver::OptimizationDirection::Maximize, + formula->getComparisonType(), + formula->getBound(), + filter + ); } template @@ -398,11 +407,18 @@ namespace storm { template - void ApproximationModel::invokeSolver(bool computeLowerBounds, Policy& policy){ + void ApproximationModel::invokeSolver(bool computeLowerBounds, Policy& policy, bool allowEarlyTermination){ storm::solver::SolveGoal player2Goal(computeLowerBounds); if(this->typeOfParametricModel == storm::models::ModelType::Dtmc){ //Invoke mdp model checking std::unique_ptr> solver = storm::solver::configureMinMaxLinearEquationSolver(player2Goal, storm::utility::solver::MinMaxLinearEquationSolverFactory(), this->matrixData.matrix); + if(allowEarlyTermination){ + solver->setEarlyTerminationCriterion(std::make_unique>( + this->solverData.player1Goal->relevantColumns(), + this->solverData.player1Goal->thresholdValue(), + computeLowerBounds + )); + } storm::utility::policyguessing::solveMinMaxLinearEquationSystem(*solver, this->solverData.result, this->vectorData.vector, player2Goal.direction(), @@ -412,9 +428,17 @@ namespace storm { } else if(this->typeOfParametricModel == storm::models::ModelType::Mdp){ //Invoke stochastic two player game model checking std::unique_ptr> solver = storm::utility::solver::GameSolverFactory().create(this->solverData.player1Matrix, this->matrixData.matrix); + if(allowEarlyTermination && this->solverData.player1Goal->minimize() == computeLowerBounds){ + //We play MIN-Min or Max-Max + solver->setEarlyTerminationCriterion(std::make_unique>( + this->solverData.player1Goal->relevantColumnVector, + this->solverData.player1Goal->thresholdValue(), + computeLowerBounds + )); + } storm::utility::policyguessing::solveGame(*solver, this->solverData.result, this->vectorData.vector, - this->solverData.player1Goal.direction(), player2Goal.direction(), + this->solverData.player1Goal->direction(), player2Goal.direction(), this->solverData.lastPlayer1Policy, policy, this->matrixData.targetChoices, (this->computeRewards ? storm::utility::infinity() : storm::utility::zero()) ); diff --git a/src/modelchecker/region/ApproximationModel.h b/src/modelchecker/region/ApproximationModel.h index 492a84fc0..21322f026 100644 --- a/src/modelchecker/region/ApproximationModel.h +++ b/src/modelchecker/region/ApproximationModel.h @@ -76,7 +76,7 @@ namespace storm { void initializeRewards(ParametricSparseModelType const& parametricModel, std::vector const& newIndices); void initializePlayer1Matrix(ParametricSparseModelType const& parametricModel); void instantiate(ParameterRegion const& region, bool computeLowerBounds); - void invokeSolver(bool computeLowerBounds, Policy& policy); + void invokeSolver(bool computeLowerBounds, Policy& policy, bool allowEarlyTermination=true); //A flag that denotes whether we compute probabilities or rewards bool computeRewards; @@ -92,7 +92,7 @@ namespace storm { Policy lastMinimizingPolicy, lastMaximizingPolicy; std::size_t initialStateIndex; //The index which represents the result for the initial state in the result vector //Player 1 represents the nondeterminism of the given mdp (so, this is irrelevant if we approximate values of a DTMC) - storm::solver::SolveGoal player1Goal = storm::solver::SolveGoal(true); //No default cunstructor for solve goal... + std::unique_ptr> player1Goal; storm::storage::SparseMatrix player1Matrix; Policy lastPlayer1Policy; } solverData; diff --git a/src/modelchecker/region/SamplingModel.cpp b/src/modelchecker/region/SamplingModel.cpp index 2c6716609..dc361129e 100644 --- a/src/modelchecker/region/SamplingModel.cpp +++ b/src/modelchecker/region/SamplingModel.cpp @@ -20,7 +20,6 @@ #include "src/utility/policyguessing.h" #include "src/exceptions/UnexpectedException.h" #include "src/exceptions/InvalidArgumentException.h" -#include "storage/dd/CuddBdd.h" namespace storm { namespace modelchecker { @@ -41,7 +40,6 @@ namespace storm { } else { STORM_LOG_THROW(false, storm::exceptions::InvalidArgumentException, "Invalid formula: " << formula << ". Sampling model only supports eventually or reachability reward formulae."); } - this->solverData.solveGoal = storm::solver::SolveGoal(storm::logic::isLowerBound(formula->getComparisonType())); STORM_LOG_THROW(parametricModel.hasLabel("target"), storm::exceptions::InvalidArgumentException, "The given Model has no \"target\"-statelabel."); this->targetStates = parametricModel.getStateLabeling().getStates("target"); STORM_LOG_THROW(parametricModel.hasLabel("sink"), storm::exceptions::InvalidArgumentException, "The given Model has no \"sink\"-statelabel."); @@ -70,6 +68,16 @@ namespace storm { this->solverData.result = std::vector(maybeStates.getNumberOfSetBits(), this->computeRewards ? storm::utility::one() : ConstantType(0.5)); this->solverData.initialStateIndex = newIndices[initialState]; this->solverData.lastPolicy = Policy(this->matrixData.matrix.getRowGroupCount(), 0); + + //this->solverData.solveGoal = storm::solver::SolveGoal(storm::logic::isLowerBound(formula->getComparisonType())); + storm::storage::BitVector filter(this->solverData.result.size(), false); + filter.set(this->solverData.initialStateIndex, true); + this->solverData.solveGoal = std::make_unique>( + storm::logic::isLowerBound(formula->getComparisonType()) ? storm::solver::OptimizationDirection::Minimize : storm::solver::OptimizationDirection::Maximize, + formula->getComparisonType(), + formula->getBound(), + filter + ); } template @@ -253,17 +261,25 @@ namespace storm { } template - void SamplingModel::invokeSolver(){ + void SamplingModel::invokeSolver(bool allowEarlyTermination){ if(this->typeOfParametricModel == storm::models::ModelType::Dtmc){ - std::unique_ptr> solver = storm::utility::solver::LinearEquationSolverFactory().create(this->matrixData.matrix); + std::unique_ptr> solver = storm::utility::solver::LinearEquationSolverFactory().create(this->matrixData.matrix); solver->solveEquationSystem(this->solverData.result, this->vectorData.vector); } else if(this->typeOfParametricModel == storm::models::ModelType::Mdp){ - std::unique_ptr> solver = storm::solver::configureMinMaxLinearEquationSolver(this->solverData.solveGoal, storm::utility::solver::MinMaxLinearEquationSolverFactory(), this->matrixData.matrix); + std::unique_ptr> solver = storm::utility::solver::MinMaxLinearEquationSolverFactory().create(this->matrixData.matrix); + solver->setOptimizationDirection(this->solverData.solveGoal->direction()); + if(allowEarlyTermination){ + solver->setEarlyTerminationCriterion(std::make_unique>( + this->solverData.solveGoal->relevantColumns(), + this->solverData.solveGoal->thresholdValue(), + this->solverData.solveGoal->minimize() + )); + } storm::utility::policyguessing::solveMinMaxLinearEquationSystem(*solver, this->solverData.result, this->vectorData.vector, - this->solverData.solveGoal.direction(), + this->solverData.solveGoal->direction(), this->solverData.lastPolicy, - this->matrixData.targetChoices, (this->computeRewards ? storm::utility::infinity() : storm::utility::zero())); + this->matrixData.targetChoices, (this->computeRewards ? storm::utility::infinity() : storm::utility::zero())); } else { STORM_LOG_THROW(false, storm::exceptions::UnexpectedException, "Unexpected Type of model"); } diff --git a/src/modelchecker/region/SamplingModel.h b/src/modelchecker/region/SamplingModel.h index cb2b15e8c..e5c65ff49 100644 --- a/src/modelchecker/region/SamplingModel.h +++ b/src/modelchecker/region/SamplingModel.h @@ -57,7 +57,7 @@ namespace storm { void initializeProbabilities(ParametricSparseModelType const& parametricModel, std::vector const& newIndices); void initializeRewards(ParametricSparseModelType const& parametricModel, std::vector const& newIndices); void instantiate(std::mapconst& point); - void invokeSolver(); + void invokeSolver(bool allowEarlyTermination=true); //A flag that denotes whether we compute probabilities or rewards bool computeRewards; @@ -72,7 +72,7 @@ namespace storm { std::vector result; //Note: result.size==maybeStates.numberOfSetBits std::size_t initialStateIndex; //The index which represents the result for the initial state in the result vector //The following is only relevant if we consider mdps: - storm::solver::SolveGoal solveGoal = storm::solver::SolveGoal(true); //No default cunstructor for solve goal... + std::unique_ptr> solveGoal; Policy lastPolicy; //best policy from the previous instantiation. Serves as first guess for the next call. } solverData; diff --git a/src/solver/AllowEarlyTerminationCondition.cpp b/src/solver/AllowEarlyTerminationCondition.cpp index 2236c7c78..25c9533b2 100644 --- a/src/solver/AllowEarlyTerminationCondition.cpp +++ b/src/solver/AllowEarlyTerminationCondition.cpp @@ -35,13 +35,34 @@ namespace storm { ValueType initVal = terminationThreshold - 1; ValueType currentThreshold = useMinimumAsExtremum ? storm::utility::vector::max_if(currentValues, filter, initVal) : storm::utility::vector::max_if(currentValues, filter, initVal); - return currentThreshold >= this->terminationThreshold; - + return currentThreshold >= this->terminationThreshold; + } + + template + TerminateAfterFilteredExtremumBelowOrAboveThreshold::TerminateAfterFilteredExtremumBelowOrAboveThreshold(storm::storage::BitVector const& filter, ValueType threshold, bool useMinimum) : + terminationThreshold(threshold), filter(filter), useMinimumAsExtremum(useMinimum) + { + // Intentionally left empty. + } + + template + bool TerminateAfterFilteredExtremumBelowOrAboveThreshold::terminateNow(const std::vector& currentValues) const { + assert(currentValues.size() >= filter.size()); + if(useMinimumAsExtremum){ + ValueType initVal = terminationThreshold + 1; + ValueType currentThreshold = storm::utility::vector::min_if(currentValues, filter, initVal); + return currentThreshold < this->terminationThreshold; + } else { + ValueType initVal = terminationThreshold - 1; + ValueType currentThreshold = storm::utility::vector::max_if(currentValues, filter, initVal); + return currentThreshold > this->terminationThreshold; + } } template class TerminateAfterFilteredExtremumPassesThresholdValue; + template class TerminateAfterFilteredExtremumBelowOrAboveThreshold; template class TerminateAfterFilteredSumPassesThresholdValue; } diff --git a/src/solver/AllowEarlyTerminationCondition.h b/src/solver/AllowEarlyTerminationCondition.h index be9084c70..d789a5e41 100644 --- a/src/solver/AllowEarlyTerminationCondition.h +++ b/src/solver/AllowEarlyTerminationCondition.h @@ -41,6 +41,18 @@ namespace storm { storm::storage::BitVector filter; bool useMinimumAsExtremum; }; + + template + class TerminateAfterFilteredExtremumBelowOrAboveThreshold : public AllowEarlyTerminationCondition{ + public: + TerminateAfterFilteredExtremumBelowOrAboveThreshold(storm::storage::BitVector const& filter, ValueType threshold, bool useMinimum); + bool terminateNow(std::vector const& currentValue) const; + + protected: + ValueType terminationThreshold; + storm::storage::BitVector filter; + bool useMinimumAsExtremum; + }; } } diff --git a/src/solver/GameSolver.cpp b/src/solver/GameSolver.cpp index 9d94149d0..dd81fe52f 100644 --- a/src/solver/GameSolver.cpp +++ b/src/solver/GameSolver.cpp @@ -9,12 +9,12 @@ namespace storm { namespace solver { template - GameSolver::GameSolver(storm::storage::SparseMatrix const& player1Matrix, storm::storage::SparseMatrix const& player2Matrix) : AbstractGameSolver(), player1Matrix(player1Matrix), player2Matrix(player2Matrix) { + GameSolver::GameSolver(storm::storage::SparseMatrix const& player1Matrix, storm::storage::SparseMatrix const& player2Matrix) : AbstractGameSolver(), player1Matrix(player1Matrix), player2Matrix(player2Matrix), earlyTermination(new NoEarlyTerminationCondition()) { // Intentionally left empty. } template - GameSolver::GameSolver(storm::storage::SparseMatrix const& player1Matrix, storm::storage::SparseMatrix const& player2Matrix, double precision, uint_fast64_t maximalNumberOfIterations, bool relative) : AbstractGameSolver(precision, maximalNumberOfIterations, relative), player1Matrix(player1Matrix), player2Matrix(player2Matrix) { + GameSolver::GameSolver(storm::storage::SparseMatrix const& player1Matrix, storm::storage::SparseMatrix const& player2Matrix, double precision, uint_fast64_t maximalNumberOfIterations, bool relative) : AbstractGameSolver(precision, maximalNumberOfIterations, relative), player1Matrix(player1Matrix), player2Matrix(player2Matrix) , earlyTermination(new NoEarlyTerminationCondition()) { // Intentionally left empty. } @@ -69,9 +69,9 @@ namespace storm { std::swap(x, tmpResult); ++iterations; - } while (!converged && iterations < maximalNumberOfIterations); + } while (!converged && iterations < maximalNumberOfIterations && !this->earlyTermination->terminateNow(x)); - STORM_LOG_WARN_COND(converged, "Iterative solver for stochastic two player games did not converge after " << iterations << " iterations."); + STORM_LOG_WARN_COND(converged || this->earlyTermination->terminateNow(x), "Iterative solver for stochastic two player games did not converge after " << iterations << " iterations."); if(this->trackPolicies){ this->player2Policy = std::vector(player2Matrix.getRowGroupCount()); @@ -109,6 +109,11 @@ namespace storm { } } } + + template + void GameSolver::setEarlyTerminationCriterion(std::unique_ptr> v) { + earlyTermination = std::move(v); + } template storm::storage::SparseMatrix const& GameSolver::getPlayer1Matrix() const { diff --git a/src/solver/GameSolver.h b/src/solver/GameSolver.h index 699533ee3..f6ad07089 100644 --- a/src/solver/GameSolver.h +++ b/src/solver/GameSolver.h @@ -4,6 +4,7 @@ #include #include "src/solver/AbstractGameSolver.h" +#include "src/solver/AllowEarlyTerminationCondition.h" #include "src/solver/OptimizationDirection.h" #include "src/storage/sparse/StateType.h" @@ -53,6 +54,8 @@ namespace storm { storm::storage::SparseMatrix const& getPlayer2Matrix() const; storm::storage::SparseMatrix const& getPlayer1Matrix() const; + void setEarlyTerminationCriterion(std::unique_ptr> v); + private: // The matrix defining the choices of player 1. storm::storage::SparseMatrix const& player1Matrix; @@ -60,6 +63,9 @@ namespace storm { // The matrix defining the choices of player 2. storm::storage::SparseMatrix const& player2Matrix; + + std::unique_ptr> earlyTermination; + }; } } diff --git a/src/solver/GmmxxMinMaxLinearEquationSolver.cpp b/src/solver/GmmxxMinMaxLinearEquationSolver.cpp index 61e2f9f1b..ad47bdc85 100644 --- a/src/solver/GmmxxMinMaxLinearEquationSolver.cpp +++ b/src/solver/GmmxxMinMaxLinearEquationSolver.cpp @@ -74,6 +74,8 @@ namespace storm { // Check if the solver converged and issue a warning otherwise. if (converged) { LOG4CPLUS_INFO(logger, "Iterative solver converged after " << iterations << " iterations."); + } else if (this->earlyTermination->terminateNow(*currentX)) { + LOG4CPLUS_INFO(logger, "Iterative solver stopped due to early termination. Performed " << iterations << " iterations."); } else { LOG4CPLUS_WARN(logger, "Iterative solver did not converge after " << iterations << " iterations."); } @@ -85,6 +87,10 @@ namespace storm { } if(this->trackPolicy){ + if(iterations==0){ //may happen due to early termination. Then we need to compute x'= A*x+b + gmm::mult(*gmmxxMatrix, x, *multiplyResult); + gmm::add(b, *multiplyResult); + } this->policy = std::vector(x.size()); storm::utility::vector::reduceVectorMinOrMax(dir, *multiplyResult, x, rowGroupIndices, &(this->policy)); } diff --git a/src/solver/NativeMinMaxLinearEquationSolver.cpp b/src/solver/NativeMinMaxLinearEquationSolver.cpp index e2eaa267b..b3e01ec76 100644 --- a/src/solver/NativeMinMaxLinearEquationSolver.cpp +++ b/src/solver/NativeMinMaxLinearEquationSolver.cpp @@ -74,6 +74,8 @@ namespace storm { // Check if the solver converged and issue a warning otherwise. if (converged) { LOG4CPLUS_INFO(logger, "Iterative solver converged after " << iterations << " iterations."); + } else if (this->earlyTermination->terminateNow(*currentX)) { + LOG4CPLUS_INFO(logger, "Iterative solver stopped due to early termination. Performed " << iterations << " iterations."); } else { LOG4CPLUS_WARN(logger, "Iterative solver did not converge after " << iterations << " iterations."); } @@ -85,6 +87,10 @@ namespace storm { } if(this->trackPolicy){ + if(iterations==0){ //may happen due to early termination. Then we need to compute x'= A*x+b + this->A.multiplyWithVector(x, *multiplyResult); + storm::utility::vector::addVectors(*multiplyResult, b, *multiplyResult); + } this->policy = std::vector(x.size()); storm::utility::vector::reduceVectorMinOrMax(dir, *multiplyResult, x, this->A.getRowGroupIndices(), &(this->policy)); } diff --git a/src/utility/policyguessing.cpp b/src/utility/policyguessing.cpp index 9c5b0c459..1981bc42d 100644 --- a/src/utility/policyguessing.cpp +++ b/src/utility/policyguessing.cpp @@ -55,43 +55,7 @@ namespace storm { solveLinearEquationSystem(inducedA, x, inducedB, probGreater0States, prob0Value, solver.getPrecision(), solver.getRelative()); //x = std::vector(x.size(), storm::utility::zero()); } - } - /* - std::size_t p2Precount=0; - std::size_t p2Postcount=0; - std::size_t p1diff =0; - std::size_t p2diff =0; - std::size_t p2RelevantCount=0; - storm::storage::BitVector relevantP2Groups(pl2Policy.size(),false); - for(std::size_t i = 0; igetColumn()); - } - for (std::size_t i : relevantP2Groups){ - if(pl2Policy[i] != player2Policy[i]){ - ++p2RelevantCount; - } - } - for(std::size_t i = 0; i @@ -209,20 +173,22 @@ namespace storm { storm::utility::vector::selectVectorValues(subB, probGreater0States, b); std::unique_ptr> linEqSysSolver = storm::utility::solver::LinearEquationSolverFactory().create(eqSysA); linEqSysSolver->setRelative(relative); - + linEqSysSolver->setIterations(500); std::size_t iterations = 0; std::vector copyX(subX.size()); ValueType newPrecision = precision; do { linEqSysSolver->setPrecision(newPrecision); if(!linEqSysSolver->solveEquationSystem(subX, subB)){ - break; //Solver did not converge.. so we have to go on with the current solution. + // break; //Solver did not converge.. so we have to go on with the current solution. } subA.multiplyWithVector(subX,copyX); storm::utility::vector::addVectors(copyX, subB, copyX); // = Ax + b ++iterations; - newPrecision *= 0.1; - } while(!storm::utility::vector::equalModuloPrecision(subX, copyX, precision, relative) && iterations<30); + newPrecision *= 0.5; + } while(!storm::utility::vector::equalModuloPrecision(subX, copyX, precision*0.5, relative) && iterations<50); + + STORM_LOG_WARN_COND(iterations<50, "Solving linear equation system did not yield a precise result"); STORM_LOG_DEBUG("Required to increase the precision " << iterations << " times in order to obtain a precise result"); //fill in the result