From bcd193dd5711445a00648d84931ff0e5669eb663 Mon Sep 17 00:00:00 2001 From: Tim Quatmann Date: Mon, 2 Dec 2019 13:53:23 +0100 Subject: [PATCH 01/24] Implemented Value iteration based LRA computation for CTMCs. --- .../csl/SparseCtmcCslModelChecker.cpp | 6 +- .../csl/helper/HybridCtmcCslHelper.cpp | 9 +- .../csl/helper/SparseCtmcCslHelper.cpp | 493 +++++++++++------- .../csl/helper/SparseCtmcCslHelper.h | 19 +- src/storm/storage/SparseMatrix.cpp | 46 +- src/storm/storage/SparseMatrix.h | 13 + 6 files changed, 366 insertions(+), 220 deletions(-) diff --git a/src/storm/modelchecker/csl/SparseCtmcCslModelChecker.cpp b/src/storm/modelchecker/csl/SparseCtmcCslModelChecker.cpp index 89a2c9e03..05555b125 100644 --- a/src/storm/modelchecker/csl/SparseCtmcCslModelChecker.cpp +++ b/src/storm/modelchecker/csl/SparseCtmcCslModelChecker.cpp @@ -132,16 +132,14 @@ namespace storm { std::unique_ptr subResultPointer = this->check(env, stateFormula); ExplicitQualitativeCheckResult const& subResult = subResultPointer->asExplicitQualitativeCheckResult(); - storm::storage::SparseMatrix probabilityMatrix = storm::modelchecker::helper::SparseCtmcCslHelper::computeProbabilityMatrix(this->getModel().getTransitionMatrix(), this->getModel().getExitRateVector()); - std::vector numericResult = storm::modelchecker::helper::SparseCtmcCslHelper::computeLongRunAverageProbabilities(env, storm::solver::SolveGoal(this->getModel(), checkTask), probabilityMatrix, subResult.getTruthValuesVector(), &this->getModel().getExitRateVector()); + std::vector numericResult = storm::modelchecker::helper::SparseCtmcCslHelper::computeLongRunAverageProbabilities(env, storm::solver::SolveGoal(this->getModel(), checkTask), this->getModel().getTransitionMatrix(), subResult.getTruthValuesVector(), &this->getModel().getExitRateVector()); return std::unique_ptr(new ExplicitQuantitativeCheckResult(std::move(numericResult))); } template std::unique_ptr SparseCtmcCslModelChecker::computeLongRunAverageRewards(Environment const& env, storm::logic::RewardMeasureType rewardMeasureType, CheckTask const& checkTask) { - storm::storage::SparseMatrix probabilityMatrix = storm::modelchecker::helper::SparseCtmcCslHelper::computeProbabilityMatrix(this->getModel().getTransitionMatrix(), this->getModel().getExitRateVector()); auto rewardModel = storm::utility::createFilteredRewardModel(this->getModel(), checkTask); - std::vector numericResult = storm::modelchecker::helper::SparseCtmcCslHelper::computeLongRunAverageRewards(env, storm::solver::SolveGoal(this->getModel(), checkTask), probabilityMatrix, rewardModel.get(), &this->getModel().getExitRateVector()); + std::vector numericResult = storm::modelchecker::helper::SparseCtmcCslHelper::computeLongRunAverageRewards(env, storm::solver::SolveGoal(this->getModel(), checkTask), this->getModel().getTransitionMatrix(), rewardModel.get(), &this->getModel().getExitRateVector()); return std::unique_ptr(new ExplicitQuantitativeCheckResult(std::move(numericResult))); } diff --git a/src/storm/modelchecker/csl/helper/HybridCtmcCslHelper.cpp b/src/storm/modelchecker/csl/helper/HybridCtmcCslHelper.cpp index d76a28884..3555c7216 100644 --- a/src/storm/modelchecker/csl/helper/HybridCtmcCslHelper.cpp +++ b/src/storm/modelchecker/csl/helper/HybridCtmcCslHelper.cpp @@ -339,19 +339,18 @@ namespace storm { template std::unique_ptr HybridCtmcCslHelper::computeLongRunAverageProbabilities(Environment const& env, storm::models::symbolic::Ctmc const& model, storm::dd::Add const& rateMatrix, storm::dd::Add const& exitRateVector, storm::dd::Bdd const& psiStates) { - storm::dd::Add probabilityMatrix = computeProbabilityMatrix(rateMatrix, exitRateVector); storm::utility::Stopwatch conversionWatch(true); // Create ODD for the translation. storm::dd::Odd odd = model.getReachableStates().createOdd(); - storm::storage::SparseMatrix explicitProbabilityMatrix = probabilityMatrix.toMatrix(odd, odd); + storm::storage::SparseMatrix explicitRateMatrix = rateMatrix.toMatrix(odd, odd); std::vector explicitExitRateVector = exitRateVector.toVector(odd); conversionWatch.stop(); STORM_LOG_INFO("Converting symbolic matrix/vector to explicit representation done in " << conversionWatch.getTimeInMilliseconds() << "ms."); - std::vector result = storm::modelchecker::helper::SparseCtmcCslHelper::computeLongRunAverageProbabilities(env, storm::solver::SolveGoal(), explicitProbabilityMatrix, psiStates.toVector(odd), &explicitExitRateVector); + std::vector result = storm::modelchecker::helper::SparseCtmcCslHelper::computeLongRunAverageProbabilities(env, storm::solver::SolveGoal(), explicitRateMatrix, psiStates.toVector(odd), &explicitExitRateVector); return std::unique_ptr(new HybridQuantitativeCheckResult(model.getReachableStates(), model.getManager().getBddZero(), model.getManager().template getAddZero(), model.getReachableStates(), std::move(odd), std::move(result))); } @@ -367,12 +366,12 @@ namespace storm { // Create ODD for the translation. storm::dd::Odd odd = model.getReachableStates().createOdd(); - storm::storage::SparseMatrix explicitProbabilityMatrix = probabilityMatrix.toMatrix(odd, odd); + storm::storage::SparseMatrix explicitRateMatrix = rateMatrix.toMatrix(odd, odd); std::vector explicitExitRateVector = exitRateVector.toVector(odd); conversionWatch.stop(); STORM_LOG_INFO("Converting symbolic matrix/vector to explicit representation done in " << conversionWatch.getTimeInMilliseconds() << "ms."); - std::vector result = storm::modelchecker::helper::SparseCtmcCslHelper::computeLongRunAverageRewards(env, storm::solver::SolveGoal(), explicitProbabilityMatrix, rewardModel.getTotalRewardVector(probabilityMatrix, model.getColumnVariables(), exitRateVector, true).toVector(odd), &explicitExitRateVector); + std::vector result = storm::modelchecker::helper::SparseCtmcCslHelper::computeLongRunAverageRewards(env, storm::solver::SolveGoal(), explicitRateMatrix, rewardModel.getTotalRewardVector(probabilityMatrix, model.getColumnVariables(), exitRateVector, true).toVector(odd), &explicitExitRateVector); return std::unique_ptr(new HybridQuantitativeCheckResult(model.getReachableStates(), model.getManager().getBddZero(), model.getManager().template getAddZero(), model.getReachableStates(), std::move(odd), std::move(result))); } diff --git a/src/storm/modelchecker/csl/helper/SparseCtmcCslHelper.cpp b/src/storm/modelchecker/csl/helper/SparseCtmcCslHelper.cpp index ed87a9e1a..2125ee677 100644 --- a/src/storm/modelchecker/csl/helper/SparseCtmcCslHelper.cpp +++ b/src/storm/modelchecker/csl/helper/SparseCtmcCslHelper.cpp @@ -25,6 +25,7 @@ #include "storm/exceptions/InvalidPropertyException.h" #include "storm/exceptions/FormatUnsupportedBySolverException.h" #include "storm/exceptions/UncheckedRequirementException.h" +#include "storm/exceptions/NotSupportedException.h" namespace storm { namespace modelchecker { @@ -370,10 +371,10 @@ namespace storm { } template - std::vector SparseCtmcCslHelper::computeLongRunAverageProbabilities(Environment const& env, storm::solver::SolveGoal&& goal, storm::storage::SparseMatrix const& probabilityMatrix, storm::storage::BitVector const& psiStates, std::vector const* exitRateVector) { + std::vector SparseCtmcCslHelper::computeLongRunAverageProbabilities(Environment const& env, storm::solver::SolveGoal&& goal, storm::storage::SparseMatrix const& rateMatrix, storm::storage::BitVector const& psiStates, std::vector const* exitRateVector) { // If there are no goal states, we avoid the computation and directly return zero. - uint_fast64_t numberOfStates = probabilityMatrix.getRowCount(); + uint_fast64_t numberOfStates = rateMatrix.getRowCount(); if (psiStates.empty()) { return std::vector(numberOfStates, storm::utility::zero()); } @@ -386,7 +387,7 @@ namespace storm { ValueType zero = storm::utility::zero(); ValueType one = storm::utility::one(); - return computeLongRunAverages(env, std::move(goal), probabilityMatrix, + return computeLongRunAverages(env, std::move(goal), rateMatrix, [&zero, &one, &psiStates] (storm::storage::sparse::state_type const& state) -> ValueType { if (psiStates.get(state)) { return one; @@ -397,16 +398,33 @@ namespace storm { } template - std::vector SparseCtmcCslHelper::computeLongRunAverageRewards(Environment const& env, storm::solver::SolveGoal&& goal, storm::storage::SparseMatrix const& probabilityMatrix, RewardModelType const& rewardModel, std::vector const* exitRateVector) { + std::vector SparseCtmcCslHelper::computeLongRunAverageRewards(Environment const& env, storm::solver::SolveGoal&& goal, storm::storage::SparseMatrix const& rateMatrix, RewardModelType const& rewardModel, std::vector const* exitRateVector) { // Only compute the result if the model has a state-based reward model. STORM_LOG_THROW(!rewardModel.empty(), storm::exceptions::InvalidPropertyException, "Missing reward model for formula. Skipping formula."); - return computeLongRunAverageRewards(env, std::move(goal), probabilityMatrix, rewardModel.getTotalRewardVector(probabilityMatrix, *exitRateVector), exitRateVector); + return computeLongRunAverages(env, std::move(goal), rateMatrix, + [&] (storm::storage::sparse::state_type const& state) -> ValueType { + ValueType result = rewardModel.hasStateRewards() ? rewardModel.getStateReward(state) : storm::utility::zero(); + if (rewardModel.hasStateActionRewards()) { + // State action rewards are multiplied with the exit rate r(s). Then, multiplying the reward with the expected time we stay at s (i.e. 1/r(s)) yields the original state reward + if (exitRateVector) { + result += rewardModel.getStateActionReward(state) * (*exitRateVector)[state]; + } else { + result += rewardModel.getStateActionReward(state); + } + } + if (rewardModel.hasTransitionRewards()) { + // Transition rewards are already multiplied with the rates + result += rateMatrix.getPointwiseProductRowSum(rewardModel.getTransitionRewardMatrix(), state); + } + return result; + }, + exitRateVector); } template - std::vector SparseCtmcCslHelper::computeLongRunAverageRewards(Environment const& env, storm::solver::SolveGoal&& goal, storm::storage::SparseMatrix const& probabilityMatrix, std::vector const& stateRewardVector, std::vector const* exitRateVector) { - return computeLongRunAverages(env, std::move(goal), probabilityMatrix, + std::vector SparseCtmcCslHelper::computeLongRunAverageRewards(Environment const& env, storm::solver::SolveGoal&& goal, storm::storage::SparseMatrix const& rateMatrix, std::vector const& stateRewardVector, std::vector const* exitRateVector) { + return computeLongRunAverages(env, std::move(goal), rateMatrix, [&stateRewardVector] (storm::storage::sparse::state_type const& state) -> ValueType { return stateRewardVector[state]; }, @@ -414,183 +432,49 @@ namespace storm { } template - std::vector SparseCtmcCslHelper::computeLongRunAverages(Environment const& env, storm::solver::SolveGoal&& goal, storm::storage::SparseMatrix const& probabilityMatrix, std::function const& valueGetter, std::vector const* exitRateVector){ - uint_fast64_t numberOfStates = probabilityMatrix.getRowCount(); - + std::vector SparseCtmcCslHelper::computeLongRunAverages(Environment const& env, storm::solver::SolveGoal&& goal, storm::storage::SparseMatrix const& rateMatrix, std::function const& valueGetter, std::vector const* exitRateVector){ + storm::storage::SparseMatrix probabilityMatrix; + if (exitRateVector) { + probabilityMatrix = computeProbabilityMatrix(rateMatrix, *exitRateVector); + } else { + probabilityMatrix = rateMatrix; + } + uint_fast64_t numberOfStates = rateMatrix.getRowCount(); + // Start by decomposing the CTMC into its BSCCs. - storm::storage::StronglyConnectedComponentDecomposition bsccDecomposition(probabilityMatrix, storm::storage::StronglyConnectedComponentDecompositionOptions().onlyBottomSccs()); + storm::storage::StronglyConnectedComponentDecomposition bsccDecomposition(rateMatrix, storm::storage::StronglyConnectedComponentDecompositionOptions().onlyBottomSccs()); STORM_LOG_DEBUG("Found " << bsccDecomposition.size() << " BSCCs."); - - // Get some data members for convenience. - ValueType one = storm::utility::one(); - ValueType zero = storm::utility::zero(); - + // Prepare the vector holding the LRA values for each of the BSCCs. - std::vector bsccLra(bsccDecomposition.size(), zero); + std::vector bsccLra; + bsccLra.reserve(bsccDecomposition.size()); - // First we check which states are in BSCCs. + // Keep track of the maximal and minimal value occuring in one of the BSCCs + ValueType maxValue, minValue; storm::storage::BitVector statesInBsccs(numberOfStates); - storm::storage::BitVector firstStatesInBsccs(numberOfStates); - - for (uint_fast64_t currentBsccIndex = 0; currentBsccIndex < bsccDecomposition.size(); ++currentBsccIndex) { - storm::storage::StronglyConnectedComponent const& bscc = bsccDecomposition[currentBsccIndex]; - - // Gather information for later use. - bool first = true; + auto backwardTransitions = rateMatrix.transpose(); + for (auto const& bscc : bsccDecomposition) { for (auto const& state : bscc) { statesInBsccs.set(state); - if (first) { - firstStatesInBsccs.set(state); - } - first = false; + } + bsccLra.push_back(computeLongRunAveragesForBscc(env, bscc, rateMatrix, backwardTransitions, valueGetter, exitRateVector)); + if (bsccLra.size() == 1) { + maxValue = bsccLra.back(); + minValue = bsccLra.back(); + } else { + maxValue = std::max(bsccLra.back(), maxValue); + minValue = std::min(bsccLra.back(), minValue); } } - storm::storage::BitVector statesNotInBsccs = ~statesInBsccs; + storm::storage::BitVector statesNotInBsccs = ~statesInBsccs; STORM_LOG_DEBUG("Found " << statesInBsccs.getNumberOfSetBits() << " states in BSCCs."); - // Prepare a vector holding the index within all states that are in BSCCs for every state. - std::vector indexInStatesInBsccs; - - // Prepare a vector that maps the index within the set of all states in BSCCs to the index of the containing BSCC. - std::vector stateToBsccIndexMap; - - if (!statesInBsccs.empty()) { - firstStatesInBsccs = firstStatesInBsccs % statesInBsccs; - - // Then we construct an equation system that yields the steady state probabilities for all states in BSCCs. - storm::storage::SparseMatrix bsccEquationSystem = probabilityMatrix.getSubmatrix(false, statesInBsccs, statesInBsccs, true); - - // Since in the fix point equation, we need to multiply the vector from the left, we convert this to a - // multiplication from the right by transposing the system. - bsccEquationSystem = bsccEquationSystem.transpose(false, true); - - // Create an auxiliary structure that makes it easy to look up the indices within the set of BSCC states. - uint_fast64_t lastIndex = 0; - uint_fast64_t currentNumberOfSetBits = 0; - indexInStatesInBsccs.reserve(probabilityMatrix.getRowCount()); - for (auto index : statesInBsccs) { - while (lastIndex <= index) { - indexInStatesInBsccs.push_back(currentNumberOfSetBits); - ++lastIndex; - } - ++currentNumberOfSetBits; - } - - stateToBsccIndexMap.resize(statesInBsccs.getNumberOfSetBits()); - for (uint_fast64_t currentBsccIndex = 0; currentBsccIndex < bsccDecomposition.size(); ++currentBsccIndex) { - storm::storage::StronglyConnectedComponent const& bscc = bsccDecomposition[currentBsccIndex]; - for (auto const& state : bscc) { - stateToBsccIndexMap[indexInStatesInBsccs[state]] = currentBsccIndex; - } - } - - // Build a different system depending on the problem format of the equation solver. - // Check solver requirements. - storm::solver::GeneralLinearEquationSolverFactory linearEquationSolverFactory; - auto requirements = linearEquationSolverFactory.getRequirements(env); - requirements.clearLowerBounds(); - requirements.clearUpperBounds(); - STORM_LOG_THROW(!requirements.hasEnabledCriticalRequirement(), storm::exceptions::UncheckedRequirementException, "Solver requirements " + requirements.getEnabledRequirementsAsString() + " not checked."); - - bool fixedPointSystem = false; - if (linearEquationSolverFactory.getEquationProblemFormat(env) == storm::solver::LinearEquationSolverProblemFormat::FixedPointSystem) { - fixedPointSystem = true; - } - - // Now build the final equation system matrix, the initial guess and the right-hand side in one go. - std::vector bsccEquationSystemRightSide(bsccEquationSystem.getColumnCount(), zero); - storm::storage::SparseMatrixBuilder builder; - for (uint_fast64_t row = 0; row < bsccEquationSystem.getRowCount(); ++row) { - - // If the current row is the first one belonging to a BSCC, we substitute it by the constraint that the - // values for states of this BSCC must sum to one. However, in order to have a non-zero value on the - // diagonal, we add the constraint of the BSCC that produces a 1 on the diagonal. - if (firstStatesInBsccs.get(row)) { - uint_fast64_t requiredBscc = stateToBsccIndexMap[row]; - storm::storage::StronglyConnectedComponent const& bscc = bsccDecomposition[requiredBscc]; - - if (fixedPointSystem) { - for (auto const& state : bscc) { - if (row == indexInStatesInBsccs[state]) { - builder.addNextValue(row, indexInStatesInBsccs[state], zero); - } else { - builder.addNextValue(row, indexInStatesInBsccs[state], -one); - } - } - } else { - for (auto const& state : bscc) { - builder.addNextValue(row, indexInStatesInBsccs[state], one); - } - } - - bsccEquationSystemRightSide[row] = one; - } else { - // Otherwise, we copy the row, and subtract 1 from the diagonal (only for the equation solver format). - for (auto& entry : bsccEquationSystem.getRow(row)) { - if (fixedPointSystem || entry.getColumn() != row) { - builder.addNextValue(row, entry.getColumn(), entry.getValue()); - } else { - builder.addNextValue(row, entry.getColumn(), entry.getValue() - one); - } - } - } - } - - // Create the initial guess for the LRAs. We take a uniform distribution over all states in a BSCC. - std::vector bsccEquationSystemSolution(bsccEquationSystem.getColumnCount(), zero); - for (uint_fast64_t bsccIndex = 0; bsccIndex < bsccDecomposition.size(); ++bsccIndex) { - storm::storage::StronglyConnectedComponent const& bscc = bsccDecomposition[bsccIndex]; - - for (auto const& state : bscc) { - bsccEquationSystemSolution[indexInStatesInBsccs[state]] = one / bscc.size(); - } - } - - bsccEquationSystem = builder.build(); - { - std::unique_ptr> solver = linearEquationSolverFactory.create(env, std::move(bsccEquationSystem)); - solver->setLowerBound(storm::utility::zero()); - solver->setUpperBound(storm::utility::one()); - solver->solveEquations(env, bsccEquationSystemSolution, bsccEquationSystemRightSide); - } - - // If exit rates were given, we need to 'fix' the results to also account for the timing behaviour. - if (exitRateVector != nullptr) { - std::vector bsccTotalValue(bsccDecomposition.size(), zero); - for (auto stateIter = statesInBsccs.begin(); stateIter != statesInBsccs.end(); ++stateIter) { - bsccTotalValue[stateToBsccIndexMap[indexInStatesInBsccs[*stateIter]]] += bsccEquationSystemSolution[indexInStatesInBsccs[*stateIter]] * (one / (*exitRateVector)[*stateIter]); - } - - for (auto stateIter = statesInBsccs.begin(); stateIter != statesInBsccs.end(); ++stateIter) { - bsccEquationSystemSolution[indexInStatesInBsccs[*stateIter]] = (bsccEquationSystemSolution[indexInStatesInBsccs[*stateIter]] * (one / (*exitRateVector)[*stateIter])) / bsccTotalValue[stateToBsccIndexMap[indexInStatesInBsccs[*stateIter]]]; - } - } - - // Calculate LRA Value for each BSCC from steady state distribution in BSCCs. - for (uint_fast64_t bsccIndex = 0; bsccIndex < bsccDecomposition.size(); ++bsccIndex) { - storm::storage::StronglyConnectedComponent const& bscc = bsccDecomposition[bsccIndex]; - - for (auto const& state : bscc) { - bsccLra[stateToBsccIndexMap[indexInStatesInBsccs[state]]] += valueGetter(state) * bsccEquationSystemSolution[indexInStatesInBsccs[state]]; - } - } - - for (uint_fast64_t bsccIndex = 0; bsccIndex < bsccDecomposition.size(); ++bsccIndex) { - STORM_LOG_DEBUG("Found LRA " << bsccLra[bsccIndex] << " for BSCC " << bsccIndex << "."); - } - } else { - for (uint_fast64_t bsccIndex = 0; bsccIndex < bsccDecomposition.size(); ++bsccIndex) { - storm::storage::StronglyConnectedComponent const& bscc = bsccDecomposition[bsccIndex]; - - // At this point, all BSCCs are known to contain exactly one state, which is why we can set all values - // directly (based on whether or not the contained state is a psi state). - bsccLra[bsccIndex] = valueGetter(*bscc.begin()); - } - - for (uint_fast64_t bsccIndex = 0; bsccIndex < bsccDecomposition.size(); ++bsccIndex) { - STORM_LOG_DEBUG("Found LRA " << bsccLra[bsccIndex] << " for BSCC " << bsccIndex << "."); + std::vector stateToBsccMap(statesInBsccs.size(), -1); + for (uint64_t bsccIndex = 0; bsccIndex < bsccDecomposition.size(); ++bsccIndex) { + for (auto const& state : bsccDecomposition[bsccIndex]) { + stateToBsccMap[state] = bsccIndex; } } @@ -605,30 +489,44 @@ namespace storm { rewardRightSide.reserve(statesNotInBsccs.getNumberOfSetBits()); for (auto state : statesNotInBsccs) { - ValueType reward = zero; - for (auto entry : probabilityMatrix.getRow(state)) { + ValueType reward = storm::utility::zero(); + for (auto entry : rateMatrix.getRow(state)) { if (statesInBsccs.get(entry.getColumn())) { - reward += entry.getValue() * bsccLra[stateToBsccIndexMap[indexInStatesInBsccs[entry.getColumn()]]]; + if (exitRateVector) { + reward += (entry.getValue() / (*exitRateVector)[state]) * bsccLra[stateToBsccMap[entry.getColumn()]]; + } else { + reward += entry.getValue() * bsccLra[stateToBsccMap[entry.getColumn()]]; + } } } rewardRightSide.push_back(reward); } - storm::storage::SparseMatrix rewardEquationSystemMatrix = probabilityMatrix.getSubmatrix(false, statesNotInBsccs, statesNotInBsccs, true); - rewardEquationSystemMatrix.convertToEquationSystem(); - - rewardSolution = std::vector(rewardEquationSystemMatrix.getColumnCount(), one); - - { - storm::solver::GeneralLinearEquationSolverFactory linearEquationSolverFactory; - // Check solver requirements - auto requirements = linearEquationSolverFactory.getRequirements(env); - STORM_LOG_THROW(!requirements.hasEnabledCriticalRequirement(), storm::exceptions::UncheckedRequirementException, "Solver requirements " + requirements.getEnabledRequirementsAsString() + " not checked."); - // Check whether we have the right input format for the solver. - STORM_LOG_THROW(linearEquationSolverFactory.getEquationProblemFormat(env) == storm::solver::LinearEquationSolverProblemFormat::EquationSystem, storm::exceptions::FormatUnsupportedBySolverException, "The selected solver does not support the required format."); - std::unique_ptr> solver = linearEquationSolverFactory.create(env, std::move(rewardEquationSystemMatrix)); - solver->solveEquations(env, rewardSolution, rewardRightSide); + // Compute reachability rewards + storm::solver::GeneralLinearEquationSolverFactory linearEquationSolverFactory; + bool isEqSysFormat = linearEquationSolverFactory.getEquationProblemFormat(env) == storm::solver::LinearEquationSolverProblemFormat::EquationSystem; + storm::storage::SparseMatrix rewardEquationSystemMatrix = rateMatrix.getSubmatrix(false, statesNotInBsccs, statesNotInBsccs, isEqSysFormat); + if (exitRateVector) { + uint64_t localRow = 0; + for (auto const& globalRow : statesNotInBsccs) { + for (auto& entry : rewardEquationSystemMatrix.getRow(localRow)) { + entry.setValue(entry.getValue() / (*exitRateVector)[globalRow]); + } + ++localRow; + } + } + if (isEqSysFormat) { + rewardEquationSystemMatrix.convertToEquationSystem(); } + rewardSolution = std::vector(rewardEquationSystemMatrix.getColumnCount(), (maxValue + minValue) / storm::utility::convertNumber(2)); + // std::cout << rewardEquationSystemMatrix << std::endl; + // std::cout << storm::utility::vector::toString(rewardRightSide) << std::endl; + std::unique_ptr> solver = linearEquationSolverFactory.create(env, std::move(rewardEquationSystemMatrix)); + solver->setBounds(minValue, maxValue); + // Check solver requirements + auto requirements = solver->getRequirements(env); + STORM_LOG_THROW(!requirements.hasEnabledCriticalRequirement(), storm::exceptions::UncheckedRequirementException, "Solver requirements " + requirements.getEnabledRequirementsAsString() + " not checked."); + solver->solveEquations(env, rewardSolution, rewardRightSide); } // Fill the result vector. @@ -653,6 +551,225 @@ namespace storm { return result; } + template + ValueType SparseCtmcCslHelper::computeLongRunAveragesForBscc(Environment const& env, storm::storage::StronglyConnectedComponent const& bscc, storm::storage::SparseMatrix const& rateMatrix, storm::storage::SparseMatrix const& backwardTransitions, std::function const& valueGetter, std::vector const* exitRateVector) { + std::cout << std::endl << "========== BSCC =======" << std::endl; + storm::storage::BitVector bsccStates(rateMatrix.getRowCount()); + for (auto const& s : bscc) { bsccStates.set(s); std::cout << s << ": " << valueGetter(s) << "\t";} + auto subm = rateMatrix.getSubmatrix(false,bsccStates,bsccStates); + std::cout << std::endl << subm << std::endl; + + std::cout << std::endl << "========== EQSYS =======" << std::endl; + + // Catch the trivial case for bsccs of size 1 + if (bscc.size() == 1) { + std::cout << std::endl << "========== Result = " << valueGetter(*bscc.begin()) << " =======" << std::endl; + + return valueGetter(*bscc.begin()); + } + + return computeLongRunAveragesForBsccVi(env, bscc, rateMatrix, backwardTransitions, valueGetter, exitRateVector); + return computeLongRunAveragesForBsccEqSys(env, bscc, rateMatrix, backwardTransitions, valueGetter, exitRateVector); + } + + template <> + storm::RationalFunction SparseCtmcCslHelper::computeLongRunAveragesForBsccVi(Environment const&, storm::storage::StronglyConnectedComponent const&, storm::storage::SparseMatrix const&, storm::storage::SparseMatrix const&, std::function const&, std::vector const*) { + STORM_LOG_THROW(false, storm::exceptions::NotSupportedException, "The requested Method for LRA computation is not supported for parametric models."); + } + + template + ValueType SparseCtmcCslHelper::computeLongRunAveragesForBsccVi(Environment const& env, storm::storage::StronglyConnectedComponent const& bscc, storm::storage::SparseMatrix const& rateMatrix, storm::storage::SparseMatrix const& backwardTransitions, std::function const& valueGetter, std::vector const* exitRateVector) { + + // Initialize data about the bscc + storm::storage::BitVector bsccStates(rateMatrix.getRowGroupCount(), false); + for (auto const& state : bscc) { + bsccStates.set(state); + } + + // Get the uniformization rate + ValueType uniformizationRate = storm::utility::one(); + if (exitRateVector) { + uniformizationRate = storm::utility::vector::max_if(*exitRateVector, bsccStates); + } + // To ensure that the model is aperiodic, we need to make sure that every Markovian state gets a self loop. + // Hence, we increase the uniformization rate a little. + uniformizationRate += storm::utility::one(); // Todo: try other values such as *=1.01 + + // Get the transitions of the submodel + typename storm::storage::SparseMatrix bsccMatrix = rateMatrix.getSubmatrix(true, bsccStates, bsccStates, true); + + // Uniformize the transitions + uint64_t subState = 0; + for (auto state : bsccStates) { + for (auto& entry : bsccMatrix.getRow(subState)) { + if (entry.getColumn() == subState) { + if (exitRateVector) { + entry.setValue(storm::utility::one() + (entry.getValue() - (*exitRateVector)[state]) / uniformizationRate); + } else { + entry.setValue(storm::utility::one() + (entry.getValue() - storm::utility::one()) / uniformizationRate); + } + } else { + entry.setValue(entry.getValue() / uniformizationRate); + } + } + ++subState; + } + + // Compute the rewards obtained in a single uniformization step + std::vector markovianRewards; + markovianRewards.reserve(bsccMatrix.getRowCount()); + for (auto const& state : bsccStates) { + ValueType stateRewardScalingFactor = storm::utility::one() / uniformizationRate; + // ValueType actionRewardScalingFactor = exitRateVector[state] / uniformizationRate; // action rewards are currently not supported + markovianRewards.push_back(valueGetter(state) * stateRewardScalingFactor); + } + + // start the iterations + ValueType precision = storm::utility::convertNumber(1e-6) / uniformizationRate; // TODO env.solver().minMax().getPrecision()) / uniformizationRate; + bool relative = true; // TODO env.solver().minMax().getRelativeTerminationCriterion(); + std::vector v(bsccMatrix.getRowCount(), storm::utility::zero()); + std::vector w = v; + std::unique_ptr> solver; + while (true) { + // Compute the values for the markovian states. We also keep track of the maximal and minimal difference between two values (for convergence checking) + auto vIt = v.begin(); + uint64_t row = 0; + ValueType newValue = markovianRewards[row] + bsccMatrix.multiplyRowWithVector(row, w); + ValueType maxDiff = newValue - *vIt; + ValueType minDiff = maxDiff; + *vIt = newValue; + for (++vIt, ++row; row < bsccMatrix.getRowCount(); ++vIt, ++row) { + newValue = markovianRewards[row] + bsccMatrix.multiplyRowWithVector(row, w); + ValueType diff = newValue - *vIt; + maxDiff = std::max(maxDiff, diff); + minDiff = std::min(minDiff, diff); + *vIt = newValue; + } + + // Check for convergence + if ((maxDiff - minDiff) <= (relative ? (precision * (v.front() + minDiff)) : precision)) { + break; + } + + // update the rhs of the equation system + ValueType referenceValue = v.front(); + storm::utility::vector::applyPointwise(v, w, [&referenceValue] (ValueType const& v_i) -> ValueType { return v_i - referenceValue; }); + } + std::cout << std::endl << "========== VI Result = " << v.front() * uniformizationRate << " =======" << std::endl; + return v.front() * uniformizationRate; + } + + template + ValueType SparseCtmcCslHelper::computeLongRunAveragesForBsccEqSys(Environment const& env, storm::storage::StronglyConnectedComponent const& bscc, storm::storage::SparseMatrix const& rateMatrix, storm::storage::SparseMatrix const& backwardTransitions, std::function const& valueGetter, std::vector const* exitRateVector) { + + // Get a mapping from global state indices to local ones. + std::unordered_map toLocalIndexMap; + uint64_t localIndex = 0; + for (auto const& globalIndex : bscc) { + toLocalIndexMap[globalIndex] = localIndex; + ++localIndex; + } + + // Build the equation system matrix + + + // Build the equation system matrix for A[s,s] = R(s,s) - r(s) & A[s,s'] = R(s,s') (where s != s') + // x_0+...+x_n=1 & x*A=0 <=> x_0+...+x_n=1 & A^t*x=0 [ <=> 1-x_1+...+x_n=x_0 & (1-A^t)*x = x ] + storm::solver::GeneralLinearEquationSolverFactory linearEquationSolverFactory; + storm::storage::SparseMatrixBuilder builder(bscc.size(), bscc.size()); + bool isEquationSystemFormat = linearEquationSolverFactory.getEquationProblemFormat(env) == storm::solver::LinearEquationSolverProblemFormat::EquationSystem; + // The first row asserts that the values sum up to one + uint64_t row = 0; + if (isEquationSystemFormat) { + for (uint64_t state = 0; state < bscc.size(); ++state) { + builder.addNextValue(row, state, storm::utility::one()); + } + } else { + for (uint64_t state = 1; state < bscc.size(); ++state) { + builder.addNextValue(row, state, -storm::utility::one()); + } + } + ++row; + // Build the equation system matrix + // We can skip the first state to make the equation system matrix square + auto stateIt = bscc.begin(); + for (++stateIt; stateIt != bscc.end(); ++stateIt) { + ValueType diagonalValue = (exitRateVector == nullptr) ? -storm::utility::one() : -(*exitRateVector)[*stateIt]; + if (!isEquationSystemFormat) { + diagonalValue = storm::utility::one() - diagonalValue; + } + bool insertedDiagonal = storm::utility::isZero(diagonalValue); + for (auto const& backwardsEntry : backwardTransitions.getRow(*stateIt)) { + auto localIndexIt = toLocalIndexMap.find(backwardsEntry.getColumn()); + if (localIndexIt != toLocalIndexMap.end()) { + ValueType val = backwardsEntry.getValue(); + if (!isEquationSystemFormat) { + val = -val; + } + uint64_t localIndex = localIndexIt->second; + if (!insertedDiagonal && localIndex == row) { + builder.addNextValue(row, localIndex, val + diagonalValue); + insertedDiagonal = true; + } else { + if (!insertedDiagonal && localIndex > row) { + builder.addNextValue(row, row, diagonalValue); + insertedDiagonal = true; + } + builder.addNextValue(row, localIndex, val); + } + } + } + if (!insertedDiagonal) { + builder.addNextValue(row, row, diagonalValue); + } + ++row; + } + + // Create a linear equation solver + auto matrix = builder.build(); + std::cout << matrix << std::endl; + auto solver = linearEquationSolverFactory.create(env, std::move(matrix)); + solver->setBounds(storm::utility::zero(), storm::utility::one()); + // Check solver requirements. + auto requirements = solver->getRequirements(env); + STORM_LOG_THROW(!requirements.hasEnabledCriticalRequirement(), storm::exceptions::UncheckedRequirementException, "Solver requirements " + requirements.getEnabledRequirementsAsString() + " not checked."); + + std::vector bsccEquationSystemRightSide(bscc.size(), storm::utility::zero()); + bsccEquationSystemRightSide.front() = storm::utility::one(); + std::vector bsccEquationSystemSolution(bscc.size(), storm::utility::one() / storm::utility::convertNumber(bscc.size())); + std::cout << storm::utility::vector::toString(bsccEquationSystemRightSide) << std::endl; + solver->solveEquations(env, bsccEquationSystemSolution, bsccEquationSystemRightSide); + std::cout << storm::utility::vector::toString(bsccEquationSystemSolution) << std::endl; + // If exit rates were given, we need to 'fix' the results to also account for the timing behaviour. + if (false && exitRateVector != nullptr) { + ValueType totalValue = storm::utility::zero(); + auto solIt = bsccEquationSystemSolution.begin(); + for (auto const& globalState : bscc) { + totalValue += (*solIt) * (storm::utility::one() / (*exitRateVector)[globalState]); + ++solIt; + } + assert(solIt == bsccEquationSystemSolution.end()); + solIt = bsccEquationSystemSolution.begin(); + for (auto const& globalState : bscc) { + *solIt = ((*solIt) * (storm::utility::one() / (*exitRateVector)[globalState])) / totalValue; + ++solIt; + } + assert(solIt == bsccEquationSystemSolution.end()); + } + + // Calculate final LRA Value + ValueType result = storm::utility::zero(); + auto solIt = bsccEquationSystemSolution.begin(); + for (auto const& globalState : bscc) { + result += valueGetter(globalState) * (*solIt); + ++solIt; + } + assert(solIt == bsccEquationSystemSolution.end()); + std::cout << std::endl << "========== Result = " << result << " =======" << std::endl; + + return result; + } + template ::SupportsExponential, int>::type> std::vector SparseCtmcCslHelper::computeAllTransientProbabilities(Environment const& env, storm::storage::SparseMatrix const& rateMatrix, storm::storage::BitVector const& initialStates, storm::storage::BitVector const& phiStates, storm::storage::BitVector const& psiStates, std::vector const& exitRates, double timeBound) { diff --git a/src/storm/modelchecker/csl/helper/SparseCtmcCslHelper.h b/src/storm/modelchecker/csl/helper/SparseCtmcCslHelper.h index 4fb01d5ff..8a45ea75f 100644 --- a/src/storm/modelchecker/csl/helper/SparseCtmcCslHelper.h +++ b/src/storm/modelchecker/csl/helper/SparseCtmcCslHelper.h @@ -14,6 +14,10 @@ namespace storm { class Environment; + namespace storage { + class StronglyConnectedComponent; + } + namespace modelchecker { namespace helper { class SparseCtmcCslHelper { @@ -52,13 +56,13 @@ namespace storm { static std::vector computeTotalRewards(Environment const& env, storm::solver::SolveGoal&& goal, storm::storage::SparseMatrix const& rateMatrix, storm::storage::SparseMatrix const& backwardTransitions, std::vector const& exitRateVector, RewardModelType const& rewardModel, bool qualitative); template - static std::vector computeLongRunAverageProbabilities(Environment const& env, storm::solver::SolveGoal&& goal, storm::storage::SparseMatrix const& probabilityMatrix, storm::storage::BitVector const& psiStates, std::vector const* exitRateVector); + static std::vector computeLongRunAverageProbabilities(Environment const& env, storm::solver::SolveGoal&& goal, storm::storage::SparseMatrix const& rateMatrix, storm::storage::BitVector const& psiStates, std::vector const* exitRateVector); template - static std::vector computeLongRunAverageRewards(Environment const& env, storm::solver::SolveGoal&& goal, storm::storage::SparseMatrix const& probabilityMatrix, RewardModelType const& rewardModel, std::vector const* exitRateVector); + static std::vector computeLongRunAverageRewards(Environment const& env, storm::solver::SolveGoal&& goal, storm::storage::SparseMatrix const& rateMatrix, RewardModelType const& rewardModel, std::vector const* exitRateVector); template - static std::vector computeLongRunAverageRewards(Environment const& env, storm::solver::SolveGoal&& goal, storm::storage::SparseMatrix const& probabilityMatrix, std::vector const& stateRewardVector, std::vector const* exitRateVector); + static std::vector computeLongRunAverageRewards(Environment const& env, storm::solver::SolveGoal&& goal, storm::storage::SparseMatrix const& rateMatrix, std::vector const& stateRewardVector, std::vector const* exitRateVector); template static std::vector computeReachabilityTimes(Environment const& env, storm::solver::SolveGoal&& goal, storm::storage::SparseMatrix const& rateMatrix, storm::storage::SparseMatrix const& backwardTransitions, std::vector const& exitRateVector, storm::storage::BitVector const& targetStates, bool qualitative); @@ -119,7 +123,14 @@ namespace storm { private: template - static std::vector computeLongRunAverages(Environment const& env, storm::solver::SolveGoal&& goal, storm::storage::SparseMatrix const& probabilityMatrix, std::function const& valueGetter, std::vector const* exitRateVector); + static std::vector computeLongRunAverages(Environment const& env, storm::solver::SolveGoal&& goal, storm::storage::SparseMatrix const& rateMatrix, std::function const& valueGetter, std::vector const* exitRateVector); + template + static ValueType computeLongRunAveragesForBscc(Environment const& env, storm::storage::StronglyConnectedComponent const& bscc, storm::storage::SparseMatrix const& rateMatrix, storm::storage::SparseMatrix const& backwardTransitions, std::function const& valueGetter, std::vector const* exitRateVector); + template + static ValueType computeLongRunAveragesForBsccVi(Environment const& env, storm::storage::StronglyConnectedComponent const& bscc, storm::storage::SparseMatrix const& rateMatrix, storm::storage::SparseMatrix const& backwardTransitions, std::function const& valueGetter, std::vector const* exitRateVector); + template + static ValueType computeLongRunAveragesForBsccEqSys(Environment const& env, storm::storage::StronglyConnectedComponent const& bscc, storm::storage::SparseMatrix const& rateMatrix, storm::storage::SparseMatrix const& backwardTransitions, std::function const& valueGetter, std::vector const* exitRateVector); + }; } } diff --git a/src/storm/storage/SparseMatrix.cpp b/src/storm/storage/SparseMatrix.cpp index b3e8cd8c3..857a8d478 100644 --- a/src/storm/storage/SparseMatrix.cpp +++ b/src/storm/storage/SparseMatrix.cpp @@ -1375,30 +1375,38 @@ namespace storm { } #endif + template + template + ResultValueType SparseMatrix::getPointwiseProductRowSum(storm::storage::SparseMatrix const& otherMatrix, index_type const& row) const { + typename storm::storage::SparseMatrix::const_iterator it1 = this->begin(row); + typename storm::storage::SparseMatrix::const_iterator ite1 = this->end(row); + typename storm::storage::SparseMatrix::const_iterator it2 = otherMatrix.begin(row); + typename storm::storage::SparseMatrix::const_iterator ite2 = otherMatrix.end(row); + + ResultValueType result = storm::utility::zero(); + for (;it1 != ite1 && it2 != ite2; ++it1) { + if (it1->getColumn() < it2->getColumn()) { + continue; + } else { + // If the precondition of this method (i.e. that the given matrix is a submatrix + // of the current one) was fulfilled, we know now that the two elements are in + // the same column, so we can multiply and add them to the row sum vector. + STORM_LOG_ASSERT(it1->getColumn() == it2->getColumn(), "The given matrix is not a submatrix of this one."); + result += it2->getValue() * OtherValueType(it1->getValue()); + ++it2; + } + } + return result; + } + template template std::vector SparseMatrix::getPointwiseProductRowSumVector(storm::storage::SparseMatrix const& otherMatrix) const { - std::vector result(rowCount, storm::utility::zero()); - - // Iterate over all elements of the current matrix and either continue with the next element in case the - // given matrix does not have a non-zero element at this column position, or multiply the two entries and - // add the result to the corresponding position in the vector. + std::vector result; + result.reserve(rowCount); for (index_type row = 0; row < rowCount && row < otherMatrix.getRowCount(); ++row) { - typename storm::storage::SparseMatrix::const_iterator it2 = otherMatrix.begin(row); - typename storm::storage::SparseMatrix::const_iterator ite2 = otherMatrix.end(row); - for (const_iterator it1 = this->begin(row), ite1 = this->end(row); it1 != ite1 && it2 != ite2; ++it1) { - if (it1->getColumn() < it2->getColumn()) { - continue; - } else { - // If the precondition of this method (i.e. that the given matrix is a submatrix - // of the current one) was fulfilled, we know now that the two elements are in - // the same column, so we can multiply and add them to the row sum vector. - result[row] += it2->getValue() * OtherValueType(it1->getValue()); - ++it2; - } - } + result.push_back(getPointwiseProductRowSum(otherMatrix, row)); } - return result; } diff --git a/src/storm/storage/SparseMatrix.h b/src/storm/storage/SparseMatrix.h index 125b3bd39..880c1ee36 100644 --- a/src/storm/storage/SparseMatrix.h +++ b/src/storm/storage/SparseMatrix.h @@ -819,6 +819,19 @@ namespace storm { */ std::pair, std::vector> getJacobiDecomposition() const; + /*! + * Performs a pointwise multiplication of the entries in the given row of this matrix and the entries of + * the given row of the other matrix and returns the sum. + * + * @param otherMatrix A reference to the matrix with which to perform the pointwise multiplication. This + * matrix must be a submatrix of the current matrix in the sense that it may not have entries at indices + * where there is no entry in the current matrix. + * @return the sum of the product of the entries in the given row. + */ + template + ResultValueType getPointwiseProductRowSum(storm::storage::SparseMatrix const& otherMatrix, index_type const& row) const; + + /*! * Performs a pointwise matrix multiplication of the matrix with the given matrix and returns a vector * containing the sum of the entries in each row of the resulting matrix. From 7017fc1ab0b829b1f51ad713aadd60f855ef4f6c Mon Sep 17 00:00:00 2001 From: Tim Quatmann Date: Tue, 3 Dec 2019 10:55:14 +0100 Subject: [PATCH 02/24] Added LRA settings. --- src/storm-pars/settings/ParsSettings.cpp | 2 + .../csl/helper/SparseCtmcCslHelper.cpp | 244 ++++++++++-------- src/storm/settings/SettingsManager.cpp | 2 + .../modules/LongRunAverageSolverSettings.cpp | 101 ++++++++ .../modules/LongRunAverageSolverSettings.h | 92 +++++++ .../dtmc/LraDtmcPrctlModelCheckerTest.cpp | 4 +- 6 files changed, 332 insertions(+), 113 deletions(-) create mode 100644 src/storm/settings/modules/LongRunAverageSolverSettings.cpp create mode 100644 src/storm/settings/modules/LongRunAverageSolverSettings.h diff --git a/src/storm-pars/settings/ParsSettings.cpp b/src/storm-pars/settings/ParsSettings.cpp index 8e8592c5e..c0aab397d 100644 --- a/src/storm-pars/settings/ParsSettings.cpp +++ b/src/storm-pars/settings/ParsSettings.cpp @@ -18,6 +18,7 @@ #include "storm/settings/modules/NativeEquationSolverSettings.h" #include "storm/settings/modules/TopologicalEquationSolverSettings.h" #include "storm/settings/modules/EliminationSettings.h" +#include "storm/settings/modules/LongRunAverageSolverSettings.h" #include "storm/settings/modules/MinMaxEquationSolverSettings.h" #include "storm/settings/modules/GameSolverSettings.h" #include "storm/settings/modules/BisimulationSettings.h" @@ -48,6 +49,7 @@ namespace storm { storm::settings::addModule(); storm::settings::addModule(); storm::settings::addModule(); + storm::settings::addModule(); storm::settings::addModule(); storm::settings::addModule(); storm::settings::addModule(); diff --git a/src/storm/modelchecker/csl/helper/SparseCtmcCslHelper.cpp b/src/storm/modelchecker/csl/helper/SparseCtmcCslHelper.cpp index 2125ee677..4d381e885 100644 --- a/src/storm/modelchecker/csl/helper/SparseCtmcCslHelper.cpp +++ b/src/storm/modelchecker/csl/helper/SparseCtmcCslHelper.cpp @@ -14,6 +14,7 @@ #include "storm/storage/StronglyConnectedComponentDecomposition.h" #include "storm/adapters/RationalFunctionAdapter.h" +#include "storm/environment/solver/LongRunAverageSolverEnvironment.h" #include "storm/utility/macros.h" #include "storm/utility/vector.h" @@ -559,17 +560,30 @@ namespace storm { auto subm = rateMatrix.getSubmatrix(false,bsccStates,bsccStates); std::cout << std::endl << subm << std::endl; - std::cout << std::endl << "========== EQSYS =======" << std::endl; - // Catch the trivial case for bsccs of size 1 - if (bscc.size() == 1) { - std::cout << std::endl << "========== Result = " << valueGetter(*bscc.begin()) << " =======" << std::endl; - - return valueGetter(*bscc.begin()); + // Catch the case where all values are the same (this includes the special case where the bscc is of size 1) + auto it = bscc.begin(); + ValueType val = valueGetter(*it); + for (++it; it != bscc.end(); ++it) { + if (valueGetter(*it) != val) { + break; + } + } + if (it == bscc.end()) { + // All entries have the same LRA + std::cout << std::endl << "==== TRIVIAL: RESULT=" << val << "=======" << std::endl; + return val; } - return computeLongRunAveragesForBsccVi(env, bscc, rateMatrix, backwardTransitions, valueGetter, exitRateVector); - return computeLongRunAveragesForBsccEqSys(env, bscc, rateMatrix, backwardTransitions, valueGetter, exitRateVector); + std::cout << std::endl << "========== SOLVING =======" << std::endl; + storm::solver::LraMethod method = env.solver().lra().getDetLraMethod(); + if (method == storm::solver::LraMethod::ValueIteration) { + return computeLongRunAveragesForBsccVi(env, bscc, rateMatrix, backwardTransitions, valueGetter, exitRateVector); + } else if (method == storm::solver::LraMethod::LraDistributionEquations) { + //return computeLongRunAveragesForBsccLraDistrEqSys(env, bscc, rateMatrix, backwardTransitions, valueGetter, exitRateVector); + } + STORM_LOG_WARN_COND(method == storm::solver::LraMethod::GainBiasEquations, "Unsupported lra method selected. Defaulting to gain-bias-equations."); + return computeLongRunAveragesForBsccEqSys(env, bscc, rateMatrix, backwardTransitions, valueGetter, exitRateVector); } template <> @@ -593,7 +607,7 @@ namespace storm { } // To ensure that the model is aperiodic, we need to make sure that every Markovian state gets a self loop. // Hence, we increase the uniformization rate a little. - uniformizationRate += storm::utility::one(); // Todo: try other values such as *=1.01 + uniformizationRate *= (storm::utility::one() + storm::utility::convertNumber(env.solver().lra().getAperiodicFactor())); // Get the transitions of the submodel typename storm::storage::SparseMatrix bsccMatrix = rateMatrix.getSubmatrix(true, bsccStates, bsccStates, true); @@ -619,48 +633,69 @@ namespace storm { std::vector markovianRewards; markovianRewards.reserve(bsccMatrix.getRowCount()); for (auto const& state : bsccStates) { - ValueType stateRewardScalingFactor = storm::utility::one() / uniformizationRate; - // ValueType actionRewardScalingFactor = exitRateVector[state] / uniformizationRate; // action rewards are currently not supported - markovianRewards.push_back(valueGetter(state) * stateRewardScalingFactor); + markovianRewards.push_back(valueGetter(state) / uniformizationRate); } // start the iterations - ValueType precision = storm::utility::convertNumber(1e-6) / uniformizationRate; // TODO env.solver().minMax().getPrecision()) / uniformizationRate; - bool relative = true; // TODO env.solver().minMax().getRelativeTerminationCriterion(); - std::vector v(bsccMatrix.getRowCount(), storm::utility::zero()); - std::vector w = v; - std::unique_ptr> solver; - while (true) { + ValueType precision = storm::utility::convertNumber(env.solver().lra().getPrecision()) / uniformizationRate; + bool relative = env.solver().lra().getRelativeTerminationCriterion(); + if (!relative) { + precision /= uniformizationRate; + } + std::vector x(bsccMatrix.getRowCount(), storm::utility::zero()); + std::vector xPrime(x.size()); + + auto multiplier = storm::solver::MultiplierFactory().create(env, bsccMatrix); + ValueType maxDiff, minDiff; + uint64_t iter = 0; + boost::optional maxIter; + if (env.solver().lra().isMaximalIterationCountSet()) { + maxIter = env.solver().lra().getMaximalIterationCount(); + } + while (!maxIter.is_initialized() || iter < maxIter.get()) { + ++iter; // Compute the values for the markovian states. We also keep track of the maximal and minimal difference between two values (for convergence checking) - auto vIt = v.begin(); - uint64_t row = 0; - ValueType newValue = markovianRewards[row] + bsccMatrix.multiplyRowWithVector(row, w); - ValueType maxDiff = newValue - *vIt; - ValueType minDiff = maxDiff; - *vIt = newValue; - for (++vIt, ++row; row < bsccMatrix.getRowCount(); ++vIt, ++row) { - newValue = markovianRewards[row] + bsccMatrix.multiplyRowWithVector(row, w); - ValueType diff = newValue - *vIt; + multiplier->multiply(env, x, &markovianRewards, xPrime); + + // update xPrime and check for convergence + // to avoid large (and numerically unstable) x-values, we substract a reference value. + auto xIt = x.begin(); + auto xPrimeIt = xPrime.begin(); + ValueType refVal = *xPrimeIt; + maxDiff = *xPrimeIt - *xIt; + minDiff = maxDiff; + *xPrimeIt -= refVal; + *xIt = *xPrimeIt; + for (++xIt, ++xPrimeIt; xIt != x.end(); ++xIt, ++xPrimeIt) { + ValueType diff = *xPrimeIt - *xIt; maxDiff = std::max(maxDiff, diff); minDiff = std::min(minDiff, diff); - *vIt = newValue; + *xPrimeIt -= refVal; + *xIt = *xPrimeIt; } - - // Check for convergence - if ((maxDiff - minDiff) <= (relative ? (precision * (v.front() + minDiff)) : precision)) { + + // Check for convergence. The uniformization rate is already incorporated into the precision parameter + if ((maxDiff - minDiff) <= (relative ? (precision * minDiff) : precision)) { break; } - // update the rhs of the equation system - ValueType referenceValue = v.front(); - storm::utility::vector::applyPointwise(v, w, [&referenceValue] (ValueType const& v_i) -> ValueType { return v_i - referenceValue; }); } - std::cout << std::endl << "========== VI Result = " << v.front() * uniformizationRate << " =======" << std::endl; - return v.front() * uniformizationRate; + return (maxDiff + minDiff) * uniformizationRate / (storm::utility::convertNumber(2.0)); + } + + template <> + storm::RationalFunction SparseCtmcCslHelper::computeLongRunAveragesForBsccEqSys(Environment const&, storm::storage::StronglyConnectedComponent const&, storm::storage::SparseMatrix const&, storm::storage::SparseMatrix const&, std::function const&, std::vector const*) { + STORM_LOG_THROW(false, storm::exceptions::NotSupportedException, "The requested Method for LRA computation is not supported for parametric models."); + // This is not possible due to uniformization } template ValueType SparseCtmcCslHelper::computeLongRunAveragesForBsccEqSys(Environment const& env, storm::storage::StronglyConnectedComponent const& bscc, storm::storage::SparseMatrix const& rateMatrix, storm::storage::SparseMatrix const& backwardTransitions, std::function const& valueGetter, std::vector const* exitRateVector) { + // We build the equation system as in Line 3 of Algorithm 3 from + // Kretinsky, Meggendorfer: Efficient Strategy Iteration for Mean Payoff in Markov Decision Processes (ATVA 2017) + // The first variable corresponds to the gain of the bscc whereas the subsequent variables yield the bias for each state s_1, s_2, .... + // No bias variable for s_0 is needed since it is always set to zero, yielding an nxn equation system matrix + // To properly deal with CMTC, we also need to uniformize the transitions. // Get a mapping from global state indices to local ones. std::unordered_map toLocalIndexMap; @@ -670,104 +705,91 @@ namespace storm { ++localIndex; } - // Build the equation system matrix - + // Get the uniformization rate + ValueType uniformizationRate = storm::utility::one(); + if (exitRateVector) { + uniformizationRate = *std::max_element(exitRateVector->begin(), exitRateVector->end()); + } - // Build the equation system matrix for A[s,s] = R(s,s) - r(s) & A[s,s'] = R(s,s') (where s != s') - // x_0+...+x_n=1 & x*A=0 <=> x_0+...+x_n=1 & A^t*x=0 [ <=> 1-x_1+...+x_n=x_0 & (1-A^t)*x = x ] + // Build the equation system matrix and vector. storm::solver::GeneralLinearEquationSolverFactory linearEquationSolverFactory; - storm::storage::SparseMatrixBuilder builder(bscc.size(), bscc.size()); bool isEquationSystemFormat = linearEquationSolverFactory.getEquationProblemFormat(env) == storm::solver::LinearEquationSolverProblemFormat::EquationSystem; - // The first row asserts that the values sum up to one + storm::storage::SparseMatrixBuilder builder(bscc.size(), bscc.size()); + std::vector eqsysVector; + eqsysVector.reserve(bscc.size()); + // The first row asserts that the weighted bias variables and the reward at s_0 sum up to the gain uint64_t row = 0; - if (isEquationSystemFormat) { - for (uint64_t state = 0; state < bscc.size(); ++state) { - builder.addNextValue(row, state, storm::utility::one()); - } - } else { - for (uint64_t state = 1; state < bscc.size(); ++state) { - builder.addNextValue(row, state, -storm::utility::one()); + ValueType entryValue; + for (auto const& globalState : bscc) { + // Coefficient for the gain variable + if (isEquationSystemFormat) { + // '1-0' in row 0 and -(-1) in other rows + builder.addNextValue(row, 0, storm::utility::one()); + } else if (row > 0) { + // No coeficient in row 0, othwerise substract the gain + builder.addNextValue(row, 0, -storm::utility::one()); } - } - ++row; - // Build the equation system matrix - // We can skip the first state to make the equation system matrix square - auto stateIt = bscc.begin(); - for (++stateIt; stateIt != bscc.end(); ++stateIt) { - ValueType diagonalValue = (exitRateVector == nullptr) ? -storm::utility::one() : -(*exitRateVector)[*stateIt]; - if (!isEquationSystemFormat) { - diagonalValue = storm::utility::one() - diagonalValue; + // Compute weighted sum over successor state. As this is a BSCC, each successor state will again be in the BSCC. + bool needDiagonalEntry = ((exitRateVector && (*exitRateVector)[globalState] != uniformizationRate) || isEquationSystemFormat) && (row > 0); + ValueType diagonalValue; + if (needDiagonalEntry) { + if (isEquationSystemFormat) { + diagonalValue = (*exitRateVector)[globalState] / uniformizationRate; // 1 - (1 - r(s)/unifrate) + } else { + diagonalValue = storm::utility::one() - (*exitRateVector)[globalState] / uniformizationRate; + } } - bool insertedDiagonal = storm::utility::isZero(diagonalValue); - for (auto const& backwardsEntry : backwardTransitions.getRow(*stateIt)) { - auto localIndexIt = toLocalIndexMap.find(backwardsEntry.getColumn()); - if (localIndexIt != toLocalIndexMap.end()) { - ValueType val = backwardsEntry.getValue(); - if (!isEquationSystemFormat) { - val = -val; - } - uint64_t localIndex = localIndexIt->second; - if (!insertedDiagonal && localIndex == row) { - builder.addNextValue(row, localIndex, val + diagonalValue); - insertedDiagonal = true; - } else { - if (!insertedDiagonal && localIndex > row) { - builder.addNextValue(row, row, diagonalValue); - insertedDiagonal = true; - } - builder.addNextValue(row, localIndex, val); + for (auto const& entry : rateMatrix.getRow(globalState)) { + uint64_t col = toLocalIndexMap[entry.getColumn()]; + if (col == 0) { + //Skip transition to state_0. This corresponds to setting the bias of state_0 to zero + continue; + } + entryValue = entry.getValue(); + if (exitRateVector) { + // entryValue /= uniformizationRate; + } + if (isEquationSystemFormat) { + entryValue = -entryValue; + } + if (needDiagonalEntry && col >= row) { + if (col == row) { + entryValue += diagonalValue; + } else { // col > row + builder.addNextValue(row, row, diagonalValue); } + needDiagonalEntry = false; } + builder.addNextValue(row, col, entryValue); } - if (!insertedDiagonal) { + if (needDiagonalEntry) { builder.addNextValue(row, row, diagonalValue); } + if (exitRateVector) { + eqsysVector.push_back(valueGetter(globalState) / uniformizationRate); + } else { + eqsysVector.push_back(valueGetter(globalState)); + } ++row; } // Create a linear equation solver auto matrix = builder.build(); std::cout << matrix << std::endl; + std::cout << "RHS: " << storm::utility::vector::toString(eqsysVector) << std::endl; auto solver = linearEquationSolverFactory.create(env, std::move(matrix)); - solver->setBounds(storm::utility::zero(), storm::utility::one()); // Check solver requirements. auto requirements = solver->getRequirements(env); STORM_LOG_THROW(!requirements.hasEnabledCriticalRequirement(), storm::exceptions::UncheckedRequirementException, "Solver requirements " + requirements.getEnabledRequirementsAsString() + " not checked."); + // Todo: Find bounds on the bias variable. Just inserting the maximal value from the vector probably does not work. - std::vector bsccEquationSystemRightSide(bscc.size(), storm::utility::zero()); - bsccEquationSystemRightSide.front() = storm::utility::one(); - std::vector bsccEquationSystemSolution(bscc.size(), storm::utility::one() / storm::utility::convertNumber(bscc.size())); - std::cout << storm::utility::vector::toString(bsccEquationSystemRightSide) << std::endl; - solver->solveEquations(env, bsccEquationSystemSolution, bsccEquationSystemRightSide); - std::cout << storm::utility::vector::toString(bsccEquationSystemSolution) << std::endl; - // If exit rates were given, we need to 'fix' the results to also account for the timing behaviour. - if (false && exitRateVector != nullptr) { - ValueType totalValue = storm::utility::zero(); - auto solIt = bsccEquationSystemSolution.begin(); - for (auto const& globalState : bscc) { - totalValue += (*solIt) * (storm::utility::one() / (*exitRateVector)[globalState]); - ++solIt; - } - assert(solIt == bsccEquationSystemSolution.end()); - solIt = bsccEquationSystemSolution.begin(); - for (auto const& globalState : bscc) { - *solIt = ((*solIt) * (storm::utility::one() / (*exitRateVector)[globalState])) / totalValue; - ++solIt; - } - assert(solIt == bsccEquationSystemSolution.end()); - } + std::vector eqSysSol(bscc.size(), storm::utility::zero()); + // Take the mean of the rewards as an initial guess for the gain + eqSysSol.front() = std::accumulate(eqsysVector.begin(), eqsysVector.end(), storm::utility::zero()) / storm::utility::convertNumber(bscc.size()); + solver->solveEquations(env, eqSysSol, eqsysVector); - // Calculate final LRA Value - ValueType result = storm::utility::zero(); - auto solIt = bsccEquationSystemSolution.begin(); - for (auto const& globalState : bscc) { - result += valueGetter(globalState) * (*solIt); - ++solIt; - } - assert(solIt == bsccEquationSystemSolution.end()); - std::cout << std::endl << "========== Result = " << result << " =======" << std::endl; - - return result; + std::cout << std::endl << "========== Result = " << eqSysSol.front() << "*" << uniformizationRate << " =======" << std::endl; + return eqSysSol.front();// * uniformizationRate; } template ::SupportsExponential, int>::type> diff --git a/src/storm/settings/SettingsManager.cpp b/src/storm/settings/SettingsManager.cpp index 74bb4d0a7..7088f924e 100644 --- a/src/storm/settings/SettingsManager.cpp +++ b/src/storm/settings/SettingsManager.cpp @@ -24,6 +24,7 @@ #include "storm/settings/modules/GmmxxEquationSolverSettings.h" #include "storm/settings/modules/NativeEquationSolverSettings.h" #include "storm/settings/modules/EliminationSettings.h" +#include "storm/settings/modules/LongRunAverageSolverSettings.h" #include "storm/settings/modules/MinMaxEquationSolverSettings.h" #include "storm/settings/modules/GameSolverSettings.h" #include "storm/settings/modules/BisimulationSettings.h" @@ -654,6 +655,7 @@ namespace storm { storm::settings::addModule(); storm::settings::addModule(); storm::settings::addModule(); + storm::settings::addModule(); storm::settings::addModule(); storm::settings::addModule(); storm::settings::addModule(); diff --git a/src/storm/settings/modules/LongRunAverageSolverSettings.cpp b/src/storm/settings/modules/LongRunAverageSolverSettings.cpp new file mode 100644 index 000000000..9826bdca8 --- /dev/null +++ b/src/storm/settings/modules/LongRunAverageSolverSettings.cpp @@ -0,0 +1,101 @@ +#include "storm/settings/modules/LongRunAverageSolverSettings.h" + +#include "storm/settings/Option.h" +#include "storm/settings/ArgumentBuilder.h" +#include "storm/settings/OptionBuilder.h" + +#include "storm/utility/macros.h" +#include "storm/exceptions/IllegalArgumentValueException.h" + +namespace storm { + namespace settings { + namespace modules { + + const std::string LongRunAverageSolverSettings::moduleName = "lra"; + const std::string LongRunAverageSolverSettings::detLraMethodOptionName = "detmethod"; + const std::string LongRunAverageSolverSettings::nondetLraMethodOptionName = "nondetmethod"; + const std::string LongRunAverageSolverSettings::maximalIterationsOptionName = "maxiter"; + const std::string LongRunAverageSolverSettings::maximalIterationsOptionShortName = "i"; + const std::string LongRunAverageSolverSettings::precisionOptionName = "precision"; + const std::string LongRunAverageSolverSettings::absoluteOptionName = "absolute"; + const std::string LongRunAverageSolverSettings::aperiodicFactorOptionName = "aperiodicfactor"; + + LongRunAverageSolverSettings::LongRunAverageSolverSettings() : ModuleSettings(moduleName) { + + std::vector detLraMethods = {"gb", "gain-bias-equations", "distr", "lra-distribution-equations", "vi", "value-iteration"}; + this->addOption(storm::settings::OptionBuilder(moduleName, detLraMethodOptionName, true, "Sets which method is preferred for computing long run averages on deterministic models.").setIsAdvanced() + .addArgument(storm::settings::ArgumentBuilder::createStringArgument("name", "The name of a long run average computation method.").addValidatorString(ArgumentValidatorFactory::createMultipleChoiceValidator(detLraMethods)).setDefaultValueString("gb").build()).build()); + + std::vector nondetLraMethods = {"vi", "value-iteration", "linear-programming", "lp"}; + this->addOption(storm::settings::OptionBuilder(moduleName, nondetLraMethodOptionName, true, "Sets which method is preferred for computing long run averages on models with nondeterminism.").setIsAdvanced() + .addArgument(storm::settings::ArgumentBuilder::createStringArgument("name", "The name of a long run average computation method.").addValidatorString(ArgumentValidatorFactory::createMultipleChoiceValidator(nondetLraMethods)).setDefaultValueString("vi").build()).build()); + + this->addOption(storm::settings::OptionBuilder(moduleName, maximalIterationsOptionName, false, "The maximal number of iterations to perform before iterative solving is aborted.").setShortName(maximalIterationsOptionShortName).setIsAdvanced().addArgument(storm::settings::ArgumentBuilder::createUnsignedIntegerArgument("count", "The maximal iteration count.").build()).build()); + + this->addOption(storm::settings::OptionBuilder(moduleName, precisionOptionName, false, "The precision used for detecting convergence of iterative methods.").setIsAdvanced().addArgument(storm::settings::ArgumentBuilder::createDoubleArgument("value", "The precision to achieve.").setDefaultValueDouble(1e-06).addValidatorDouble(ArgumentValidatorFactory::createDoubleRangeValidatorExcluding(0.0, 1.0)).build()).build()); + + this->addOption(storm::settings::OptionBuilder(moduleName, absoluteOptionName, false, "Sets whether the relative or the absolute error is considered for detecting convergence.").setIsAdvanced().build()); + + this->addOption(storm::settings::OptionBuilder(moduleName, aperiodicFactorOptionName, true, "If required by the selected method (e.g. vi), this factor controls how the system is made aperiodic").setIsAdvanced().addArgument(storm::settings::ArgumentBuilder::createDoubleArgument("value", "The factor.").setDefaultValueDouble(0.125).addValidatorDouble(ArgumentValidatorFactory::createDoubleRangeValidatorExcluding(0.0, 1.0)).build()).build()); + + } + + storm::solver::LraMethod LongRunAverageSolverSettings::getDetLraMethod() const { + std::string lraMethodString = this->getOption(nondetLraMethodOptionName).getArgumentByName("name").getValueAsString(); + if (lraMethodString == "gain-bias-equations" || lraMethodString == "gb") { + return storm::solver::LraMethod::GainBiasEquations; + } + if (lraMethodString == "lra-distribution-equations" || lraMethodString == "distr") { + return storm::solver::LraMethod::LraDistributionEquations; + } + if (lraMethodString == "value-iteration" || lraMethodString == "vi") { + return storm::solver::LraMethod::ValueIteration; + } + STORM_LOG_THROW(false, storm::exceptions::IllegalArgumentValueException, "Unknown lra solving technique for deterministic models:'" << lraMethodString << "'."); + } + + bool LongRunAverageSolverSettings::isDetLraMethodSetFromDefaultValue() const { + return !this->getOption(detLraMethodOptionName).getArgumentByName("name").getHasBeenSet() || this->getOption(detLraMethodOptionName).getArgumentByName("name").wasSetFromDefaultValue(); + } + + storm::solver::LraMethod LongRunAverageSolverSettings::getNondetLraMethod() const { + std::string lraMethodString = this->getOption(nondetLraMethodOptionName).getArgumentByName("name").getValueAsString(); + if (lraMethodString == "value-iteration" || lraMethodString == "vi") { + return storm::solver::LraMethod::ValueIteration; + } else if (lraMethodString == "linear-programming" || lraMethodString == "lp") { + return storm::solver::LraMethod::LinearProgramming; + } + STORM_LOG_THROW(false, storm::exceptions::IllegalArgumentValueException, "Unknown lra solving technique for nondeterministic models:'" << lraMethodString << "'."); + } + + bool LongRunAverageSolverSettings::isNondetLraMethodSetFromDefaultValue() const { + return !this->getOption(nondetLraMethodOptionName).getArgumentByName("name").getHasBeenSet() || this->getOption(nondetLraMethodOptionName).getArgumentByName("name").wasSetFromDefaultValue(); + } + + bool LongRunAverageSolverSettings::isMaximalIterationCountSet() const { + return this->getOption(maximalIterationsOptionName).getHasOptionBeenSet(); + } + + uint_fast64_t LongRunAverageSolverSettings::getMaximalIterationCount() const { + return this->getOption(maximalIterationsOptionName).getArgumentByName("count").getValueAsUnsignedInteger(); + } + + bool LongRunAverageSolverSettings::isPrecisionSet() const { + return this->getOption(precisionOptionName).getHasOptionBeenSet(); + } + + double LongRunAverageSolverSettings::getPrecision() const { + return this->getOption(precisionOptionName).getArgumentByName("value").getValueAsDouble(); + } + + bool LongRunAverageSolverSettings::isRelativePrecision() const { + return this->getOption(absoluteOptionName).getHasOptionBeenSet(); + } + + double LongRunAverageSolverSettings::getAperiodicFactor() const { + return this->getOption(aperiodicFactorOptionName).getArgumentByName("value").getValueAsDouble(); + } + + } + } +} diff --git a/src/storm/settings/modules/LongRunAverageSolverSettings.h b/src/storm/settings/modules/LongRunAverageSolverSettings.h new file mode 100644 index 000000000..413adf34d --- /dev/null +++ b/src/storm/settings/modules/LongRunAverageSolverSettings.h @@ -0,0 +1,92 @@ +#pragma once + +#include "storm-config.h" +#include "storm/settings/modules/ModuleSettings.h" + +#include "storm/solver/SolverSelectionOptions.h" +#include "storm/solver/MultiplicationStyle.h" + +namespace storm { + namespace settings { + namespace modules { + + /*! + * This class represents the LRA solver settings. + */ + class LongRunAverageSolverSettings : public ModuleSettings { + public: + + LongRunAverageSolverSettings(); + + /*! + * Retrieves the selected long run average method for deterministic models. + */ + storm::solver::LraMethod getDetLraMethod() const; + + /*! + * Retrieves whether the LraMethod for deterministic models was set from a default value. + */ + bool isDetLraMethodSetFromDefaultValue() const; + + /*! + * Retrieves the selected long run average method for nondeterministic models. + */ + storm::solver::LraMethod getNondetLraMethod() const; + + /*! + * Retrieves whether the LraMethod for nondeterministic models was set from a default value. + */ + bool isNondetLraMethodSetFromDefaultValue() const; + + /*! + * Retrieves whether a maximal iteration count for iterative methods was set. + */ + bool isMaximalIterationCountSet() const; + + /*! + * Retrieves the maximal number of iterations to perform until giving up on converging. + * + * @return The maximal iteration count. + */ + uint_fast64_t getMaximalIterationCount() const; + + /*! + * Retrieves whether the precision has been set. + * + * @return True iff the precision has been set. + */ + bool isPrecisionSet() const; + + /*! + * Retrieves the precision that is used for detecting convergence. + * + * @return The precision to use for detecting convergence. + */ + double getPrecision() const; + + /*! + * Retrieves whether the convergence criterion has been set to relative. + */ + bool isRelativePrecision() const; + + /*! + * Retrieves a factor that describes how the system is made aperiodic (if necessary by the method) + */ + double getAperiodicFactor() const; + + // The name of the module. + static const std::string moduleName; + + private: + static const std::string detLraMethodOptionName; + static const std::string nondetLraMethodOptionName; + static const std::string maximalIterationsOptionName; + static const std::string maximalIterationsOptionShortName; + static const std::string precisionOptionName; + static const std::string absoluteOptionName; + static const std::string aperiodicFactorOptionName; + }; + + } + } +} diff --git a/src/test/storm/modelchecker/prctl/dtmc/LraDtmcPrctlModelCheckerTest.cpp b/src/test/storm/modelchecker/prctl/dtmc/LraDtmcPrctlModelCheckerTest.cpp index 7caec1167..6c63a1d21 100755 --- a/src/test/storm/modelchecker/prctl/dtmc/LraDtmcPrctlModelCheckerTest.cpp +++ b/src/test/storm/modelchecker/prctl/dtmc/LraDtmcPrctlModelCheckerTest.cpp @@ -52,7 +52,7 @@ namespace { class EigenRationalLUEnvironment { public: - typedef double ValueType; + typedef storm::RationalNumber ValueType; static const bool isExact = true; static storm::Environment createEnvironment() { storm::Environment env; @@ -70,7 +70,7 @@ namespace { storm::Environment env; env.solver().setLinearEquationSolverType(storm::solver::EquationSolverType::Native); env.solver().native().setMethod(storm::solver::NativeLinearEquationSolverMethod::SOR); - env.solver().native().setSorOmega(storm::utility::convertNumber(0.9)); + env.solver().native().setSorOmega(storm::utility::convertNumber(0.8)); // A test fails if this is set to 0.9... env.solver().native().setPrecision(storm::utility::convertNumber(1e-8)); return env; } From 7a026922b738cd777385b7afcd38ed1f3e019f08 Mon Sep 17 00:00:00 2001 From: Tim Quatmann Date: Tue, 3 Dec 2019 10:55:46 +0100 Subject: [PATCH 03/24] Added LRA Environment --- src/storm/environment/SubEnvironment.cpp | 2 + .../LongRunAverageSolverEnvironment.cpp | 94 +++++++++++++++++++ .../solver/LongRunAverageSolverEnvironment.h | 49 ++++++++++ .../environment/solver/SolverEnvironment.cpp | 9 ++ .../environment/solver/SolverEnvironment.h | 4 + src/storm/solver/SolverSelectionOptions.h | 2 +- 6 files changed, 159 insertions(+), 1 deletion(-) create mode 100644 src/storm/environment/solver/LongRunAverageSolverEnvironment.cpp create mode 100644 src/storm/environment/solver/LongRunAverageSolverEnvironment.h diff --git a/src/storm/environment/SubEnvironment.cpp b/src/storm/environment/SubEnvironment.cpp index 5226fdca4..c42e21c93 100644 --- a/src/storm/environment/SubEnvironment.cpp +++ b/src/storm/environment/SubEnvironment.cpp @@ -8,6 +8,7 @@ #include "storm/environment/solver/EigenSolverEnvironment.h" #include "storm/environment/solver/GmmxxSolverEnvironment.h" #include "storm/environment/solver/NativeSolverEnvironment.h" +#include "storm/environment/solver/LongRunAverageSolverEnvironment.h" #include "storm/environment/solver/MinMaxSolverEnvironment.h" #include "storm/environment/solver/MultiplierEnvironment.h" #include "storm/environment/solver/GameSolverEnvironment.h" @@ -63,6 +64,7 @@ namespace storm { template class SubEnvironment; template class SubEnvironment; template class SubEnvironment; + template class SubEnvironment; template class SubEnvironment; template class SubEnvironment; template class SubEnvironment; diff --git a/src/storm/environment/solver/LongRunAverageSolverEnvironment.cpp b/src/storm/environment/solver/LongRunAverageSolverEnvironment.cpp new file mode 100644 index 000000000..83e9d4f7d --- /dev/null +++ b/src/storm/environment/solver/LongRunAverageSolverEnvironment.cpp @@ -0,0 +1,94 @@ +#include "storm/environment/solver/LongRunAverageSolverEnvironment.h" + +#include "storm/settings/SettingsManager.h" +#include "storm/settings/modules/LongRunAverageSolverSettings.h" +#include "storm/utility/constants.h" +#include "storm/utility/macros.h" + +namespace storm { + + LongRunAverageSolverEnvironment::LongRunAverageSolverEnvironment() { + auto const& lraSettings = storm::settings::getModule(); + detMethod = lraSettings.getDetLraMethod(); + detMethodSetFromDefault = lraSettings.isDetLraMethodSetFromDefaultValue(); + nondetMethod = lraSettings.getNondetLraMethod(); + nondetMethodSetFromDefault = lraSettings.isNondetLraMethodSetFromDefaultValue(); + precision = storm::utility::convertNumber(lraSettings.getPrecision()); + relative = lraSettings.isRelativePrecision(); + if (lraSettings.isMaximalIterationCountSet()) { + maxIters = lraSettings.getMaximalIterationCount(); + } + aperiodicFactor = lraSettings.getAperiodicFactor(); + } + + LongRunAverageSolverEnvironment::~LongRunAverageSolverEnvironment() { + // Intentionally left empty + } + + storm::solver::LraMethod const& LongRunAverageSolverEnvironment::getDetLraMethod() const { + return detMethod; + } + + bool const& LongRunAverageSolverEnvironment::isDetLraMethodSetFromDefault() const { + return detMethodSetFromDefault; + } + + void LongRunAverageSolverEnvironment::setDetLraMethod(storm::solver::LraMethod value, bool isSetFromDefault) { + detMethod = value; + detMethodSetFromDefault = isSetFromDefault; + } + + storm::solver::LraMethod const& LongRunAverageSolverEnvironment::getNondetLraMethod() const { + return nondetMethod; + } + + bool const& LongRunAverageSolverEnvironment::isNondetLraMethodSetFromDefault() const { + return nondetMethodSetFromDefault; + } + + void LongRunAverageSolverEnvironment::setNondetLraMethod(storm::solver::LraMethod value, bool isSetFromDefault) { + nondetMethod = value; + nondetMethodSetFromDefault = isSetFromDefault; + } + + storm::RationalNumber const& LongRunAverageSolverEnvironment::getPrecision() const { + return precision; + } + + void LongRunAverageSolverEnvironment::setPrecision(storm::RationalNumber value) { + precision = value; + } + + bool const& LongRunAverageSolverEnvironment::getRelativeTerminationCriterion() const { + return relative; + } + + void LongRunAverageSolverEnvironment::setRelativeTerminationCriterion(bool value) { + relative = value; + } + + bool LongRunAverageSolverEnvironment::isMaximalIterationCountSet() const { + return maxIters.is_initialized(); + } + + uint64_t LongRunAverageSolverEnvironment::getMaximalIterationCount() const { + return maxIters.get(); + } + + void LongRunAverageSolverEnvironment::setMaximalIterationCount(uint64_t value) { + maxIters = value; + } + + void LongRunAverageSolverEnvironment::unsetMaximalIterationCount() { + maxIters = boost::none; + } + + storm::RationalNumber const& LongRunAverageSolverEnvironment::getAperiodicFactor() const { + return aperiodicFactor; + } + + void LongRunAverageSolverEnvironment::setAperiodicFactor(storm::RationalNumber value) { + aperiodicFactor = value; + } + +} diff --git a/src/storm/environment/solver/LongRunAverageSolverEnvironment.h b/src/storm/environment/solver/LongRunAverageSolverEnvironment.h new file mode 100644 index 000000000..bf21972d4 --- /dev/null +++ b/src/storm/environment/solver/LongRunAverageSolverEnvironment.h @@ -0,0 +1,49 @@ +#pragma once + +#include "storm/environment/solver/SolverEnvironment.h" +#include "storm/solver/SolverSelectionOptions.h" + +namespace storm { + + class LongRunAverageSolverEnvironment { + public: + + LongRunAverageSolverEnvironment(); + ~LongRunAverageSolverEnvironment(); + + storm::solver::LraMethod const& getDetLraMethod() const; + bool const& isDetLraMethodSetFromDefault() const; + void setDetLraMethod(storm::solver::LraMethod value, bool isSetFromDefault = false); + + storm::solver::LraMethod const& getNondetLraMethod() const; + bool const& isNondetLraMethodSetFromDefault() const; + void setNondetLraMethod(storm::solver::LraMethod value, bool isSetFromDefault = false); + + storm::RationalNumber const& getPrecision() const; + void setPrecision(storm::RationalNumber value); + bool const& getRelativeTerminationCriterion() const; + void setRelativeTerminationCriterion(bool value); + + bool isMaximalIterationCountSet() const; + uint64_t getMaximalIterationCount() const; + void setMaximalIterationCount(uint64_t value); + void unsetMaximalIterationCount(); + + storm::RationalNumber const& getAperiodicFactor() const; + void setAperiodicFactor(storm::RationalNumber value); + + private: + storm::solver::LraMethod detMethod; + bool detMethodSetFromDefault; + + storm::solver::LraMethod nondetMethod; + bool nondetMethodSetFromDefault; + + storm::RationalNumber precision; + bool relative; + boost::optional maxIters; + + storm::RationalNumber aperiodicFactor; + }; +} + diff --git a/src/storm/environment/solver/SolverEnvironment.cpp b/src/storm/environment/solver/SolverEnvironment.cpp index 42600adbb..07a02aca1 100644 --- a/src/storm/environment/solver/SolverEnvironment.cpp +++ b/src/storm/environment/solver/SolverEnvironment.cpp @@ -1,5 +1,6 @@ #include "storm/environment/solver/SolverEnvironment.h" +#include "storm/environment/solver/LongRunAverageSolverEnvironment.h" #include "storm/environment/solver/MinMaxSolverEnvironment.h" #include "storm/environment/solver/MultiplierEnvironment.h" #include "storm/environment/solver/EigenSolverEnvironment.h" @@ -29,6 +30,14 @@ namespace storm { // Intentionally left empty } + LongRunAverageSolverEnvironment& SolverEnvironment::lra() { + return longRunAverageSolverEnvironment.get(); + } + + LongRunAverageSolverEnvironment const& SolverEnvironment::lra() const { + return longRunAverageSolverEnvironment.get(); + } + MinMaxSolverEnvironment& SolverEnvironment::minMax() { return minMaxSolverEnvironment.get(); } diff --git a/src/storm/environment/solver/SolverEnvironment.h b/src/storm/environment/solver/SolverEnvironment.h index 5e24652e1..26ac33fea 100644 --- a/src/storm/environment/solver/SolverEnvironment.h +++ b/src/storm/environment/solver/SolverEnvironment.h @@ -14,6 +14,7 @@ namespace storm { class EigenSolverEnvironment; class GmmxxSolverEnvironment; class NativeSolverEnvironment; + class LongRunAverageSolverEnvironment; class MinMaxSolverEnvironment; class MultiplierEnvironment; class GameSolverEnvironment; @@ -31,6 +32,8 @@ namespace storm { GmmxxSolverEnvironment const& gmmxx() const; NativeSolverEnvironment& native(); NativeSolverEnvironment const& native() const; + LongRunAverageSolverEnvironment& lra(); + LongRunAverageSolverEnvironment const& lra() const; MinMaxSolverEnvironment& minMax(); MinMaxSolverEnvironment const& minMax() const; MultiplierEnvironment& multiplier(); @@ -56,6 +59,7 @@ namespace storm { SubEnvironment nativeSolverEnvironment; SubEnvironment gameSolverEnvironment; SubEnvironment topologicalSolverEnvironment; + SubEnvironment longRunAverageSolverEnvironment; SubEnvironment minMaxSolverEnvironment; SubEnvironment multiplierEnvironment; diff --git a/src/storm/solver/SolverSelectionOptions.h b/src/storm/solver/SolverSelectionOptions.h index 448878b72..d21262099 100644 --- a/src/storm/solver/SolverSelectionOptions.h +++ b/src/storm/solver/SolverSelectionOptions.h @@ -9,7 +9,7 @@ namespace storm { ExtendEnumsWithSelectionField(MinMaxMethod, PolicyIteration, ValueIteration, LinearProgramming, Topological, RationalSearch, IntervalIteration, SoundValueIteration, TopologicalCuda, ViToPi) ExtendEnumsWithSelectionField(MultiplierType, Native, Gmmxx) ExtendEnumsWithSelectionField(GameMethod, PolicyIteration, ValueIteration) - ExtendEnumsWithSelectionField(LraMethod, LinearProgramming, ValueIteration) + ExtendEnumsWithSelectionField(LraMethod, LinearProgramming, ValueIteration, GainBiasEquations, LraDistributionEquations) ExtendEnumsWithSelectionField(LpSolverType, Gurobi, Glpk, Z3) ExtendEnumsWithSelectionField(EquationSolverType, Native, Gmmxx, Eigen, Elimination, Topological) From 324eb23cdd2a43ab9bceec6d4e562c2cd6632826 Mon Sep 17 00:00:00 2001 From: Tim Quatmann Date: Tue, 3 Dec 2019 10:56:57 +0100 Subject: [PATCH 04/24] Using new LRA environment --- .../helper/SparseMarkovAutomatonCslHelper.cpp | 50 ++++++++++++++----- .../prctl/helper/SparseMdpPrctlHelper.cpp | 33 ++++++++---- 2 files changed, 62 insertions(+), 21 deletions(-) diff --git a/src/storm/modelchecker/csl/helper/SparseMarkovAutomatonCslHelper.cpp b/src/storm/modelchecker/csl/helper/SparseMarkovAutomatonCslHelper.cpp index 14e133014..2a39daf3c 100644 --- a/src/storm/modelchecker/csl/helper/SparseMarkovAutomatonCslHelper.cpp +++ b/src/storm/modelchecker/csl/helper/SparseMarkovAutomatonCslHelper.cpp @@ -13,6 +13,9 @@ #include "storm/environment/Environment.h" #include "storm/environment/solver/MinMaxSolverEnvironment.h" +#include "storm/environment/solver/TopologicalSolverEnvironment.h" +#include "storm/environment/solver/LongRunAverageSolverEnvironment.h" +#include "storm/environment/solver/EigenSolverEnvironment.h" #include "storm/utility/macros.h" #include "storm/utility/vector.h" @@ -703,7 +706,9 @@ namespace storm { auto underlyingSolverEnvironment = env; if (env.solver().isForceSoundness()) { // For sound computations, the error in the MECS plus the error in the remaining system should be less then the user defined precsion. - underlyingSolverEnvironment.solver().minMax().setPrecision(env.solver().minMax().getPrecision() / storm::utility::convertNumber(2)); + underlyingSolverEnvironment.solver().minMax().setPrecision(env.solver().lra().getPrecision() / storm::utility::convertNumber(2)); + underlyingSolverEnvironment.solver().minMax().setRelativeTerminationCriterion(env.solver().lra().getRelativeTerminationCriterion()); + underlyingSolverEnvironment.solver().lra().setPrecision(env.solver().lra().getPrecision() / storm::utility::convertNumber(2)); } for (uint64_t currentMecIndex = 0; currentMecIndex < mecDecomposition.size(); ++currentMecIndex) { @@ -877,11 +882,11 @@ namespace storm { } // Solve MEC with the method specified in the settings - storm::solver::LraMethod method = env.solver().minMax().getLraMethod(); - if (storm::NumberTraits::IsExact && env.solver().minMax().isLraMethodSetFromDefault() && method != storm::solver::LraMethod::LinearProgramming) { + storm::solver::LraMethod method = env.solver().lra().getNondetLraMethod(); + if (storm::NumberTraits::IsExact && env.solver().lra().isNondetLraMethodSetFromDefault() && method != storm::solver::LraMethod::LinearProgramming) { STORM_LOG_INFO("Selecting 'LP' as the solution technique for long-run properties to guarantee exact results. If you want to override this, please explicitly specify a different LRA method."); method = storm::solver::LraMethod::LinearProgramming; - } else if (env.solver().isForceSoundness() && env.solver().minMax().isLraMethodSetFromDefault() && method != storm::solver::LraMethod::ValueIteration) { + } else if (env.solver().isForceSoundness() && env.solver().lra().isNondetLraMethodSetFromDefault() && method != storm::solver::LraMethod::ValueIteration) { STORM_LOG_INFO("Selecting 'VI' as the solution technique for long-run properties to guarantee sound results. If you want to override this, please explicitly specify a different LRA method."); method = storm::solver::LraMethod::ValueIteration; } @@ -982,7 +987,7 @@ namespace storm { ValueType uniformizationRate = storm::utility::vector::max_if(exitRateVector, markovianMecStates); // To ensure that the model is aperiodic, we need to make sure that every Markovian state gets a self loop. // Hence, we increase the uniformization rate a little. - uniformizationRate += storm::utility::one(); // Todo: try other values such as *=1.01 + uniformizationRate *= (storm::utility::one() + storm::utility::convertNumber(env.solver().lra().getAperiodicFactor())); // Get the transitions of the submodel, that is // * a matrix aMarkovian with all (uniformized) transitions from Markovian mec states to all Markovian mec states. @@ -1042,24 +1047,35 @@ namespace storm { // start the iterations - ValueType precision = storm::utility::convertNumber(env.solver().minMax().getPrecision()) / uniformizationRate; - bool relative = env.solver().minMax().getRelativeTerminationCriterion(); + ValueType precision = storm::utility::convertNumber(env.solver().lra().getPrecision()) / uniformizationRate; + bool relative = env.solver().lra().getRelativeTerminationCriterion(); std::vector v(aMarkovian.getRowCount(), storm::utility::zero()); std::vector w = v; std::vector x, b; + auto solverEnv = env; std::unique_ptr> solver; if (hasProbabilisticStates) { + if (env.solver().isForceSoundness()) { + // To get correct results, the inner equation systems are solved exactly. + // TODO investigate how an error would propagate + solverEnv.solver().minMax().setMethod(storm::solver::MinMaxMethod::Topological); + solverEnv.solver().topological().setUnderlyingMinMaxMethod(storm::solver::MinMaxMethod::PolicyIteration); + solverEnv.solver().setLinearEquationSolverType(storm::solver::EquationSolverType::Topological); + solverEnv.solver().topological().setUnderlyingEquationSolverType(storm::solver::EquationSolverType::Eigen); + solverEnv.solver().eigen().setMethod(storm::solver::EigenLinearEquationSolverMethod::SparseLU); + } + x.resize(aProbabilistic.getRowGroupCount(), storm::utility::zero()); b = probabilisticChoiceRewards; // Check for requirements of the solver. // The solution is unique as we assume non-zeno MAs. storm::solver::GeneralMinMaxLinearEquationSolverFactory minMaxLinearEquationSolverFactory; - storm::solver::MinMaxLinearEquationSolverRequirements requirements = minMaxLinearEquationSolverFactory.getRequirements(env, true, true, dir); + storm::solver::MinMaxLinearEquationSolverRequirements requirements = minMaxLinearEquationSolverFactory.getRequirements(solverEnv, true, true, dir); requirements.clearLowerBounds(); STORM_LOG_THROW(!requirements.hasEnabledCriticalRequirement(), storm::exceptions::UncheckedRequirementException, "Solver requirements " + requirements.getEnabledRequirementsAsString() + " not checked."); - solver = minMaxLinearEquationSolverFactory.create(env, std::move(aProbabilistic)); + solver = minMaxLinearEquationSolverFactory.create(solverEnv, std::move(aProbabilistic)); solver->setLowerBound(storm::utility::zero()); solver->setHasUniqueSolution(true); solver->setHasNoEndComponents(true); @@ -1067,10 +1083,16 @@ namespace storm { solver->setCachingEnabled(true); } - while (true) { + uint64_t iter = 0; + boost::optional maxIter; + if (env.solver().lra().isMaximalIterationCountSet()) { + maxIter = env.solver().lra().getMaximalIterationCount(); + } + while (!maxIter.is_initialized() || iter < maxIter.get()) { + ++iter; // Compute the expected total rewards for the probabilistic states if (hasProbabilisticStates) { - solver->solveEquations(env, dir, x, b); + solver->solveEquations(solverEnv, dir, x, b); } // now compute the values for the markovian states. We also keep track of the maximal and minimal difference between two values (for convergence checking) auto vIt = v.begin(); @@ -1106,8 +1128,12 @@ namespace storm { storm::utility::vector::addVectors(b, probabilisticChoiceRewards, b); } } + if (maxIter.is_initialized() && iter == maxIter.get()) { + STORM_LOG_WARN("LRA computation did not converge within " << iter << " iterations."); + } else { + STORM_LOG_TRACE("LRA computation converged after " << iter << " iterations."); + } return v.front() * uniformizationRate; - } template std::vector SparseMarkovAutomatonCslHelper::computeBoundedUntilProbabilities(Environment const& env, OptimizationDirection dir, storm::storage::SparseMatrix const& transitionMatrix, std::vector const& exitRateVector, storm::storage::BitVector const& markovianStates, storm::storage::BitVector const& psiStates, std::pair const& boundsPair); diff --git a/src/storm/modelchecker/prctl/helper/SparseMdpPrctlHelper.cpp b/src/storm/modelchecker/prctl/helper/SparseMdpPrctlHelper.cpp index f970244d1..0bcc31726 100644 --- a/src/storm/modelchecker/prctl/helper/SparseMdpPrctlHelper.cpp +++ b/src/storm/modelchecker/prctl/helper/SparseMdpPrctlHelper.cpp @@ -38,6 +38,7 @@ #include "storm/transformer/EndComponentEliminator.h" #include "storm/environment/solver/MinMaxSolverEnvironment.h" +#include "storm/environment/solver/LongRunAverageSolverEnvironment.h" #include "storm/exceptions/InvalidStateException.h" #include "storm/exceptions/InvalidPropertyException.h" @@ -1251,8 +1252,10 @@ namespace storm { auto underlyingSolverEnvironment = env; if (env.solver().isForceSoundness()) { - // For sound computations, the error in the MECS plus the error in the remaining system should be less then the user defined precsion. - underlyingSolverEnvironment.solver().minMax().setPrecision(env.solver().minMax().getPrecision() / storm::utility::convertNumber(2)); + // For sound computations, the error in the MECS plus the error in the remaining system should be less than the user defined precsion. + underlyingSolverEnvironment.solver().lra().setPrecision(env.solver().lra().getPrecision() / storm::utility::convertNumber(2)); + underlyingSolverEnvironment.solver().minMax().setPrecision(env.solver().lra().getPrecision() / storm::utility::convertNumber(2)); + underlyingSolverEnvironment.solver().minMax().setRelativeTerminationCriterion(env.solver().lra().getRelativeTerminationCriterion()); } for (uint_fast64_t currentMecIndex = 0; currentMecIndex < mecDecomposition.size(); ++currentMecIndex) { @@ -1473,11 +1476,11 @@ namespace storm { } // Solve MEC with the method specified in the settings - storm::solver::LraMethod method = env.solver().minMax().getLraMethod(); - if (storm::NumberTraits::IsExact && env.solver().minMax().isLraMethodSetFromDefault() && method != storm::solver::LraMethod::LinearProgramming) { + storm::solver::LraMethod method = env.solver().lra().getNondetLraMethod(); + if (storm::NumberTraits::IsExact && env.solver().lra().isNondetLraMethodSetFromDefault() && method != storm::solver::LraMethod::LinearProgramming) { STORM_LOG_INFO("Selecting 'LP' as the solution technique for long-run properties to guarantee exact results. If you want to override this, please explicitly specify a different LRA method."); method = storm::solver::LraMethod::LinearProgramming; - } else if (env.solver().isForceSoundness() && env.solver().minMax().isLraMethodSetFromDefault() && method != storm::solver::LraMethod::ValueIteration) { + } else if (env.solver().isForceSoundness() && env.solver().lra().isNondetLraMethodSetFromDefault() && method != storm::solver::LraMethod::ValueIteration) { STORM_LOG_INFO("Selecting 'VI' as the solution technique for long-run properties to guarantee sound results. If you want to override this, please explicitly specify a different LRA method."); method = storm::solver::LraMethod::ValueIteration; } @@ -1518,7 +1521,7 @@ namespace storm { std::vector choiceRewards; choiceRewards.reserve(mecChoices.getNumberOfSetBits()); uint64_t currRow = 0; - ValueType selfLoopProb = storm::utility::convertNumber(0.1); // todo try other values + ValueType selfLoopProb = storm::utility::convertNumber(env.solver().lra().getAperiodicFactor()); ValueType scalingFactor = storm::utility::one() - selfLoopProb; for (auto const& mecState : mecStates) { mecTransitionBuilder.newRowGroup(currRow); @@ -1553,14 +1556,21 @@ namespace storm { STORM_LOG_ASSERT(mecTransitions.isProbabilistic(), "The MEC-Matrix is not probabilistic."); // start the iterations - ValueType precision = storm::utility::convertNumber(env.solver().minMax().getPrecision()) / scalingFactor; - bool relative = env.solver().minMax().getRelativeTerminationCriterion(); + ValueType precision = storm::utility::convertNumber(env.solver().lra().getPrecision()) / scalingFactor; + bool relative = env.solver().lra().getRelativeTerminationCriterion(); std::vector x(mecTransitions.getRowGroupCount(), storm::utility::zero()); std::vector xPrime = x; auto multiplier = storm::solver::MultiplierFactory().create(env, mecTransitions); ValueType maxDiff, minDiff; - while (true) { + + uint64_t iter = 0; + boost::optional maxIter; + if (env.solver().lra().isMaximalIterationCountSet()) { + maxIter = env.solver().lra().getMaximalIterationCount(); + } + while (!maxIter.is_initialized() || iter < maxIter.get()) { + ++iter; // Compute the obtained rewards for the next step multiplier->multiplyAndReduce(env, dir, x, &choiceRewards, x); @@ -1585,6 +1595,11 @@ namespace storm { break; } } + if (maxIter.is_initialized() && iter == maxIter.get()) { + STORM_LOG_WARN("LRA computation did not converge within " << iter << " iterations."); + } else { + STORM_LOG_TRACE("LRA computation converged after " << iter << " iterations."); + } if (scheduler) { std::vector localMecChoices(mecTransitions.getRowGroupCount(), 0); From 068c1b3ea67b40510386b02f685c1d7ad47fc1cb Mon Sep 17 00:00:00 2001 From: Tim Quatmann Date: Tue, 3 Dec 2019 10:58:25 +0100 Subject: [PATCH 05/24] Removed obsolete settings --- .../solver/MinMaxSolverEnvironment.cpp | 17 ----------------- .../solver/MinMaxSolverEnvironment.h | 5 ----- .../modules/MinMaxEquationSolverSettings.cpp | 19 ------------------- .../modules/MinMaxEquationSolverSettings.h | 13 ------------- .../prctl/mdp/LraMdpPrctlModelCheckerTest.cpp | 12 ++++++------ 5 files changed, 6 insertions(+), 60 deletions(-) diff --git a/src/storm/environment/solver/MinMaxSolverEnvironment.cpp b/src/storm/environment/solver/MinMaxSolverEnvironment.cpp index 63f29573e..63cf5c698 100644 --- a/src/storm/environment/solver/MinMaxSolverEnvironment.cpp +++ b/src/storm/environment/solver/MinMaxSolverEnvironment.cpp @@ -22,8 +22,6 @@ namespace storm { STORM_LOG_ASSERT(considerRelativeTerminationCriterion || minMaxSettings.getConvergenceCriterion() == storm::settings::modules::MinMaxEquationSolverSettings::ConvergenceCriterion::Absolute, "Unknown convergence criterion"); multiplicationStyle = minMaxSettings.getValueIterationMultiplicationStyle(); symmetricUpdates = minMaxSettings.isForceIntervalIterationSymmetricUpdatesSet(); - lraMethod = minMaxSettings.getLraMethod(); - lraMethodSetFromDefault = minMaxSettings.isLraMethodSetFromDefaultValue(); } MinMaxSolverEnvironment::~MinMaxSolverEnvironment() { @@ -83,19 +81,4 @@ namespace storm { symmetricUpdates = value; } - storm::solver::LraMethod const& MinMaxSolverEnvironment::getLraMethod() const { - return lraMethod; - } - - bool const& MinMaxSolverEnvironment::isLraMethodSetFromDefault() const { - return lraMethodSetFromDefault; - } - - void MinMaxSolverEnvironment::setLraMethod(storm::solver::LraMethod value, bool isSetFromDefault) { - lraMethod = value; - lraMethodSetFromDefault = isSetFromDefault; - } - - - } diff --git a/src/storm/environment/solver/MinMaxSolverEnvironment.h b/src/storm/environment/solver/MinMaxSolverEnvironment.h index 104eee739..b0a19d8c4 100644 --- a/src/storm/environment/solver/MinMaxSolverEnvironment.h +++ b/src/storm/environment/solver/MinMaxSolverEnvironment.h @@ -27,9 +27,6 @@ namespace storm { void setMultiplicationStyle(storm::solver::MultiplicationStyle value); bool isSymmetricUpdatesSet() const; void setSymmetricUpdates(bool value); - storm::solver::LraMethod const& getLraMethod() const; - bool const& isLraMethodSetFromDefault() const; - void setLraMethod(storm::solver::LraMethod value, bool isSetFromDefault = false); private: storm::solver::MinMaxMethod minMaxMethod; @@ -39,8 +36,6 @@ namespace storm { bool considerRelativeTerminationCriterion; storm::solver::MultiplicationStyle multiplicationStyle; bool symmetricUpdates; - storm::solver::LraMethod lraMethod; - bool lraMethodSetFromDefault; }; } diff --git a/src/storm/settings/modules/MinMaxEquationSolverSettings.cpp b/src/storm/settings/modules/MinMaxEquationSolverSettings.cpp index eb429d272..c8a4f7141 100644 --- a/src/storm/settings/modules/MinMaxEquationSolverSettings.cpp +++ b/src/storm/settings/modules/MinMaxEquationSolverSettings.cpp @@ -17,7 +17,6 @@ namespace storm { const std::string MinMaxEquationSolverSettings::maximalIterationsOptionShortName = "i"; const std::string MinMaxEquationSolverSettings::precisionOptionName = "precision"; const std::string MinMaxEquationSolverSettings::absoluteOptionName = "absolute"; - const std::string MinMaxEquationSolverSettings::lraMethodOptionName = "lramethod"; const std::string MinMaxEquationSolverSettings::markovAutomatonBoundedReachabilityMethodOptionName = "mamethod"; const std::string MinMaxEquationSolverSettings::valueIterationMultiplicationStyleOptionName = "vimult"; const std::string MinMaxEquationSolverSettings::intervalIterationSymmetricUpdatesOptionName = "symmetricupdates"; @@ -33,10 +32,6 @@ namespace storm { this->addOption(storm::settings::OptionBuilder(moduleName, absoluteOptionName, false, "Sets whether the relative or the absolute error is considered for detecting convergence.").setIsAdvanced().build()); - std::vector lraMethods = {"vi", "value-iteration", "linear-programming", "lp"}; - this->addOption(storm::settings::OptionBuilder(moduleName, lraMethodOptionName, false, "Sets which method is preferred for computing long run averages.").setIsAdvanced() - .addArgument(storm::settings::ArgumentBuilder::createStringArgument("name", "The name of a long run average computation method.").addValidatorString(ArgumentValidatorFactory::createMultipleChoiceValidator(lraMethods)).setDefaultValueString("vi").build()).build()); - std::vector maMethods = {"imca", "unifplus"}; this->addOption(storm::settings::OptionBuilder(moduleName, markovAutomatonBoundedReachabilityMethodOptionName, false, "The method to use to solve bounded reachability queries on MAs.").setIsAdvanced().addArgument(storm::settings::ArgumentBuilder::createStringArgument("name", "The name of the method to use.").addValidatorString(ArgumentValidatorFactory::createMultipleChoiceValidator(maMethods)).setDefaultValueString("unifplus").build()).build()); @@ -103,20 +98,6 @@ namespace storm { return this->getOption(absoluteOptionName).getHasOptionBeenSet() ? MinMaxEquationSolverSettings::ConvergenceCriterion::Absolute : MinMaxEquationSolverSettings::ConvergenceCriterion::Relative; } - storm::solver::LraMethod MinMaxEquationSolverSettings::getLraMethod() const { - std::string lraMethodString = this->getOption(lraMethodOptionName).getArgumentByName("name").getValueAsString(); - if (lraMethodString == "value-iteration" || lraMethodString == "vi") { - return storm::solver::LraMethod::ValueIteration; - } else if (lraMethodString == "linear-programming" || lraMethodString == "lp") { - return storm::solver::LraMethod::LinearProgramming; - } - STORM_LOG_THROW(false, storm::exceptions::IllegalArgumentValueException, "Unknown lra solving technique '" << lraMethodString << "'."); - } - - bool MinMaxEquationSolverSettings::isLraMethodSetFromDefaultValue() const { - return !this->getOption(lraMethodOptionName).getArgumentByName("name").getHasBeenSet() || this->getOption(lraMethodOptionName).getArgumentByName("name").wasSetFromDefaultValue(); - } - MinMaxEquationSolverSettings::MarkovAutomatonBoundedReachabilityMethod MinMaxEquationSolverSettings::getMarkovAutomatonBoundedReachabilityMethod() const { std::string techniqueAsString = this->getOption(markovAutomatonBoundedReachabilityMethodOptionName).getArgumentByName("name").getValueAsString(); if (techniqueAsString == "imca") { diff --git a/src/storm/settings/modules/MinMaxEquationSolverSettings.h b/src/storm/settings/modules/MinMaxEquationSolverSettings.h index 2b958d374..5e7bf9b98 100644 --- a/src/storm/settings/modules/MinMaxEquationSolverSettings.h +++ b/src/storm/settings/modules/MinMaxEquationSolverSettings.h @@ -86,18 +86,6 @@ namespace storm { */ ConvergenceCriterion getConvergenceCriterion() const; - /*! - * Retrieves the selected long run average method. - * - * @return The selected long run average method. - */ - storm::solver::LraMethod getLraMethod() const; - - /*! - * Retrieves whether the LraMethod was set from a default value. - */ - bool isLraMethodSetFromDefaultValue() const; - /*! * Retrieves the method to be used for bounded reachability on MAs. * @@ -126,7 +114,6 @@ namespace storm { static const std::string maximalIterationsOptionShortName; static const std::string precisionOptionName; static const std::string absoluteOptionName; - static const std::string lraMethodOptionName; static const std::string markovAutomatonBoundedReachabilityMethodOptionName; static const std::string valueIterationMultiplicationStyleOptionName; static const std::string intervalIterationSymmetricUpdatesOptionName; diff --git a/src/test/storm/modelchecker/prctl/mdp/LraMdpPrctlModelCheckerTest.cpp b/src/test/storm/modelchecker/prctl/mdp/LraMdpPrctlModelCheckerTest.cpp index 0a6ed97a3..ff5140626 100755 --- a/src/test/storm/modelchecker/prctl/mdp/LraMdpPrctlModelCheckerTest.cpp +++ b/src/test/storm/modelchecker/prctl/mdp/LraMdpPrctlModelCheckerTest.cpp @@ -22,7 +22,7 @@ #include "storm/settings/modules/NativeEquationSolverSettings.h" #include "storm-parsers/parser/AutoParser.h" -#include "storm/environment/solver/MinMaxSolverEnvironment.h" +#include "storm/environment/solver/LongRunAverageSolverEnvironment.h" namespace { @@ -34,8 +34,8 @@ namespace { typedef storm::models::sparse::Mdp ModelType; static storm::Environment createEnvironment() { storm::Environment env; - env.solver().minMax().setLraMethod(storm::solver::LraMethod::ValueIteration); - env.solver().minMax().setPrecision(storm::utility::convertNumber(1e-10)); + env.solver().lra().setNondetLraMethod(storm::solver::LraMethod::ValueIteration); + env.solver().lra().setPrecision(storm::utility::convertNumber(1e-10)); return env; } }; @@ -47,8 +47,8 @@ namespace { typedef storm::models::sparse::Mdp ModelType; static storm::Environment createEnvironment() { storm::Environment env; - env.solver().minMax().setLraMethod(storm::solver::LraMethod::LinearProgramming); - env.solver().minMax().setPrecision(storm::utility::convertNumber(1e-10)); + env.solver().lra().setNondetLraMethod(storm::solver::LraMethod::LinearProgramming); + env.solver().lra().setPrecision(storm::utility::convertNumber(1e-10)); return env; } }; @@ -60,7 +60,7 @@ namespace { typedef storm::models::sparse::Mdp ModelType; static storm::Environment createEnvironment() { storm::Environment env; - env.solver().minMax().setLraMethod(storm::solver::LraMethod::LinearProgramming); + env.solver().lra().setNondetLraMethod(storm::solver::LraMethod::LinearProgramming); return env; } }; From 6110a677f5796f32db182051fb251918de83d38a Mon Sep 17 00:00:00 2001 From: Tim Quatmann Date: Tue, 3 Dec 2019 13:46:24 +0100 Subject: [PATCH 06/24] More environments checked in Lra Dtmc test. --- .../dtmc/LraDtmcPrctlModelCheckerTest.cpp | 91 +++++++++++++++++-- 1 file changed, 81 insertions(+), 10 deletions(-) diff --git a/src/test/storm/modelchecker/prctl/dtmc/LraDtmcPrctlModelCheckerTest.cpp b/src/test/storm/modelchecker/prctl/dtmc/LraDtmcPrctlModelCheckerTest.cpp index 6c63a1d21..3282607d4 100755 --- a/src/test/storm/modelchecker/prctl/dtmc/LraDtmcPrctlModelCheckerTest.cpp +++ b/src/test/storm/modelchecker/prctl/dtmc/LraDtmcPrctlModelCheckerTest.cpp @@ -16,80 +16,147 @@ #include "storm/environment/solver/NativeSolverEnvironment.h" #include "storm/environment/solver/EigenSolverEnvironment.h" #include "storm/environment/solver/GmmxxSolverEnvironment.h" +#include "storm/environment/solver/LongRunAverageSolverEnvironment.h" #include "storm-parsers/parser/AutoParser.h" #include "storm/builder/ExplicitModelBuilder.h" namespace { - class GmmxxDoubleGmresEnvironment { + class GBGmmxxDoubleGmresEnvironment { public: typedef double ValueType; static const bool isExact = false; static storm::Environment createEnvironment() { storm::Environment env; + env.solver().lra().setDetLraMethod(storm::solver::LraMethod::GainBiasEquations); env.solver().setLinearEquationSolverType(storm::solver::EquationSolverType::Gmmxx); env.solver().gmmxx().setMethod(storm::solver::GmmxxLinearEquationSolverMethod::Gmres); env.solver().gmmxx().setPreconditioner(storm::solver::GmmxxLinearEquationSolverPreconditioner::Ilu); env.solver().gmmxx().setPrecision(storm::utility::convertNumber(1e-8)); + env.solver().lra().setPrecision(storm::utility::convertNumber(1e-8)); return env; } }; - class EigenDoubleDGmresEnvironment { + class GBEigenDoubleDGmresEnvironment { public: typedef double ValueType; static const bool isExact = false; static storm::Environment createEnvironment() { storm::Environment env; + env.solver().lra().setDetLraMethod(storm::solver::LraMethod::GainBiasEquations); env.solver().setLinearEquationSolverType(storm::solver::EquationSolverType::Eigen); env.solver().eigen().setMethod(storm::solver::EigenLinearEquationSolverMethod::DGmres); env.solver().eigen().setPreconditioner(storm::solver::EigenLinearEquationSolverPreconditioner::Ilu); + env.solver().lra().setPrecision(storm::utility::convertNumber(1e-8)); env.solver().eigen().setPrecision(storm::utility::convertNumber(1e-8)); return env; } }; - class EigenRationalLUEnvironment { + class GBEigenRationalLUEnvironment { public: typedef storm::RationalNumber ValueType; static const bool isExact = true; static storm::Environment createEnvironment() { storm::Environment env; + env.solver().lra().setDetLraMethod(storm::solver::LraMethod::GainBiasEquations); env.solver().setLinearEquationSolverType(storm::solver::EquationSolverType::Eigen); env.solver().eigen().setMethod(storm::solver::EigenLinearEquationSolverMethod::SparseLU); return env; } }; - class NativeSorEnvironment { + class GBNativeSorEnvironment { public: typedef double ValueType; static const bool isExact = false; static storm::Environment createEnvironment() { storm::Environment env; + env.solver().lra().setDetLraMethod(storm::solver::LraMethod::GainBiasEquations); env.solver().setLinearEquationSolverType(storm::solver::EquationSolverType::Native); env.solver().native().setMethod(storm::solver::NativeLinearEquationSolverMethod::SOR); env.solver().native().setSorOmega(storm::utility::convertNumber(0.8)); // A test fails if this is set to 0.9... env.solver().native().setPrecision(storm::utility::convertNumber(1e-8)); + env.solver().lra().setPrecision(storm::utility::convertNumber(1e-8)); return env; } }; - class NativeWalkerChaeEnvironment { + class GBNativeWalkerChaeEnvironment { public: typedef double ValueType; static const bool isExact = false; static storm::Environment createEnvironment() { storm::Environment env; + env.solver().lra().setDetLraMethod(storm::solver::LraMethod::GainBiasEquations); env.solver().setLinearEquationSolverType(storm::solver::EquationSolverType::Native); env.solver().native().setMethod(storm::solver::NativeLinearEquationSolverMethod::WalkerChae); env.solver().native().setPrecision(storm::utility::convertNumber(1e-8)); + env.solver().lra().setPrecision(storm::utility::convertNumber(1e-8)); env.solver().native().setMaximalNumberOfIterations(50000); return env; } }; + class DistrGmmxxDoubleGmresEnvironment { + public: + typedef double ValueType; + static const bool isExact = false; + static storm::Environment createEnvironment() { + storm::Environment env; + env.solver().lra().setDetLraMethod(storm::solver::LraMethod::LraDistributionEquations); + env.solver().setLinearEquationSolverType(storm::solver::EquationSolverType::Gmmxx); + env.solver().gmmxx().setMethod(storm::solver::GmmxxLinearEquationSolverMethod::Gmres); + env.solver().gmmxx().setPreconditioner(storm::solver::GmmxxLinearEquationSolverPreconditioner::Ilu); + env.solver().gmmxx().setPrecision(storm::utility::convertNumber(1e-8)); + env.solver().lra().setPrecision(storm::utility::convertNumber(1e-8)); + return env; + } + }; + + class DistrEigenRationalLUEnvironment { + public: + typedef storm::RationalNumber ValueType; + static const bool isExact = true; + static storm::Environment createEnvironment() { + storm::Environment env; + env.solver().lra().setDetLraMethod(storm::solver::LraMethod::LraDistributionEquations); + env.solver().setLinearEquationSolverType(storm::solver::EquationSolverType::Eigen); + env.solver().eigen().setMethod(storm::solver::EigenLinearEquationSolverMethod::SparseLU); + return env; + } + }; + + class DistrNativeWalkerChaeEnvironment { + public: + typedef double ValueType; + static const bool isExact = false; + static storm::Environment createEnvironment() { + storm::Environment env; + env.solver().lra().setDetLraMethod(storm::solver::LraMethod::GainBiasEquations); + env.solver().setLinearEquationSolverType(storm::solver::EquationSolverType::Native); + env.solver().native().setMethod(storm::solver::NativeLinearEquationSolverMethod::WalkerChae); + env.solver().lra().setPrecision(storm::utility::convertNumber(1e-8)); + env.solver().native().setPrecision(storm::utility::convertNumber(1e-8)); + env.solver().native().setMaximalNumberOfIterations(50000); + return env; + } + }; + + class ValueIterationEnvironment { + public: + typedef double ValueType; + static const bool isExact = false; + static storm::Environment createEnvironment() { + storm::Environment env; + env.solver().lra().setDetLraMethod(storm::solver::LraMethod::ValueIteration); + env.solver().lra().setPrecision(storm::utility::convertNumber(1e-8)); + return env; + } + }; + template class LraDtmcPrctlModelCheckerTest : public ::testing::Test { public: @@ -104,11 +171,15 @@ namespace { }; typedef ::testing::Types< - GmmxxDoubleGmresEnvironment, - EigenDoubleDGmresEnvironment, - EigenRationalLUEnvironment, - NativeSorEnvironment, - NativeWalkerChaeEnvironment + GBGmmxxDoubleGmresEnvironment, + GBEigenDoubleDGmresEnvironment, + GBEigenRationalLUEnvironment, + GBNativeSorEnvironment, + GBNativeWalkerChaeEnvironment, + DistrGmmxxDoubleGmresEnvironment, + DistrEigenRationalLUEnvironment, + DistrNativeWalkerChaeEnvironment, + ValueIterationEnvironment > TestingTypes; TYPED_TEST_SUITE(LraDtmcPrctlModelCheckerTest, TestingTypes,); From f9f845bb792dafee8f70603306b6e24f77e9f461 Mon Sep 17 00:00:00 2001 From: Tim Quatmann Date: Tue, 3 Dec 2019 13:46:51 +0100 Subject: [PATCH 07/24] Separated LRA tests from CTMC tests and added a testcase for LRA Rewards --- .../examples/testfiles/ctmc/lrarewards.sm | 11 + .../csl/CtmcCslModelCheckerTest.cpp | 14 - .../csl/LraCtmcCslModelCheckerTest.cpp | 400 ++++++++++++++++++ 3 files changed, 411 insertions(+), 14 deletions(-) create mode 100644 resources/examples/testfiles/ctmc/lrarewards.sm create mode 100755 src/test/storm/modelchecker/csl/LraCtmcCslModelCheckerTest.cpp diff --git a/resources/examples/testfiles/ctmc/lrarewards.sm b/resources/examples/testfiles/ctmc/lrarewards.sm new file mode 100644 index 000000000..cdbcd2ea9 --- /dev/null +++ b/resources/examples/testfiles/ctmc/lrarewards.sm @@ -0,0 +1,11 @@ +ctmc + +module main + s : [0..3]; + [] s < 3 -> 3/2: (s'=s+1); + [] s > 0 -> 3: (s'=s-1); +endmodule + +rewards + true: s; +endrewards diff --git a/src/test/storm/modelchecker/csl/CtmcCslModelCheckerTest.cpp b/src/test/storm/modelchecker/csl/CtmcCslModelCheckerTest.cpp index 98a91f945..1adc011b5 100755 --- a/src/test/storm/modelchecker/csl/CtmcCslModelCheckerTest.cpp +++ b/src/test/storm/modelchecker/csl/CtmcCslModelCheckerTest.cpp @@ -302,7 +302,6 @@ namespace { formulasString += "; P=? [ !\"minimum\" U>=1 \"minimum\"]"; formulasString += "; P=? [ \"minimum\" U>=1 !\"minimum\"]"; formulasString += "; R=? [C<=100]"; - formulasString += "; LRA=? [\"minimum\"]"; auto modelFormulas = this->buildModelFormulas(STORM_TEST_RESOURCES_DIR "/ctmc/cluster2.sm", formulasString); auto model = std::move(modelFormulas.first); @@ -333,9 +332,6 @@ namespace { result = checker->check(this->env(), tasks[6]); EXPECT_NEAR(this->parseNumber("0.8602815057967503"), this->getQuantitativeResultAtInitialState(model, result), this->precision()); - - result = checker->check(this->env(), tasks[7]); - EXPECT_NEAR(this->parseNumber("0.99999766034263426"), this->getQuantitativeResultAtInitialState(model, result), this->precision()); } @@ -345,7 +341,6 @@ namespace { formulasString += "; P=? [ !\"down\" U<=10000 \"fail_io\"]"; formulasString += "; P=? [ !\"down\" U<=10000 \"fail_sensors\"]"; formulasString += "; R{\"up\"}=? [C<=10000]"; - formulasString += "; LRA=? [\"fail_sensors\"]"; auto modelFormulas = this->buildModelFormulas(STORM_TEST_RESOURCES_DIR "/ctmc/embedded2.sm", formulasString); auto model = std::move(modelFormulas.first); @@ -371,13 +366,10 @@ namespace { result = checker->check(this->env(), tasks[4]); EXPECT_NEAR(this->parseNumber("2.7745274082080154"), this->getQuantitativeResultAtInitialState(model, result), this->precision()); - result = checker->check(this->env(), tasks[5]); - EXPECT_NEAR(this->parseNumber("0.93458866427696596"), this->getQuantitativeResultAtInitialState(model, result), this->precision()); } TYPED_TEST(CtmcCslModelCheckerTest, Polling) { std::string formulasString = "P=?[ F<=10 \"target\"]"; - formulasString += "; LRA=?[\"target\"]"; auto modelFormulas = this->buildModelFormulas(STORM_TEST_RESOURCES_DIR "/ctmc/polling2.sm", formulasString); auto model = std::move(modelFormulas.first); @@ -391,8 +383,6 @@ namespace { result = checker->check(this->env(), tasks[0]); EXPECT_NEAR(this->parseNumber("1"), this->getQuantitativeResultAtInitialState(model, result), this->precision()); - result = checker->check(this->env(), tasks[1]); - EXPECT_NEAR(this->parseNumber("0.20079750055570736"), this->getQuantitativeResultAtInitialState(model, result), this->precision()); } @@ -403,7 +393,6 @@ namespace { formulasString += "; R=? [I=10]"; formulasString += "; R=? [C<=10]"; formulasString += "; R=? [F \"first_queue_full\"&\"second_queue_full\"]"; - formulasString += "; LRA=? [\"first_queue_full\"]"; auto modelFormulas = this->buildModelFormulas(STORM_TEST_RESOURCES_DIR "/ctmc/tandem5.sm", formulasString); auto model = std::move(modelFormulas.first); @@ -431,9 +420,6 @@ namespace { result = checker->check(this->env(), tasks[5]); EXPECT_NEAR(this->parseNumber("262.85103824276308"), this->getQuantitativeResultAtInitialState(model, result), this->precision()); - - result = checker->check(this->env(), tasks[6]); - EXPECT_NEAR(this->parseNumber("0.9100373532"), this->getQuantitativeResultAtInitialState(model, result), this->precision()); } TYPED_TEST(CtmcCslModelCheckerTest, simple2) { diff --git a/src/test/storm/modelchecker/csl/LraCtmcCslModelCheckerTest.cpp b/src/test/storm/modelchecker/csl/LraCtmcCslModelCheckerTest.cpp new file mode 100755 index 000000000..629684560 --- /dev/null +++ b/src/test/storm/modelchecker/csl/LraCtmcCslModelCheckerTest.cpp @@ -0,0 +1,400 @@ +#include "test/storm_gtest.h" +#include "test/storm_gtest.h" +#include "storm-config.h" + +#include "storm/api/builder.h" +#include "storm-parsers/api/model_descriptions.h" +#include "storm/api/properties.h" +#include "storm-conv/api/storm-conv.h" +#include "storm-parsers/api/properties.h" +#include "storm-parsers/parser/FormulaParser.h" +#include "storm/logic/Formulas.h" +#include "storm/models/sparse/StandardRewardModel.h" +#include "storm/models/symbolic/StandardRewardModel.h" +#include "storm/models/sparse/Ctmc.h" +#include "storm/models/symbolic/Ctmc.h" +#include "storm/modelchecker/csl/SparseCtmcCslModelChecker.h" +#include "storm/modelchecker/csl/HybridCtmcCslModelChecker.h" +#include "storm/modelchecker/results/QuantitativeCheckResult.h" +#include "storm/modelchecker/results/ExplicitQualitativeCheckResult.h" +#include "storm/modelchecker/results/SymbolicQualitativeCheckResult.h" +#include "storm/modelchecker/results/QualitativeCheckResult.h" +#include "storm-parsers/parser/PrismParser.h" +#include "storm/storage/expressions/ExpressionManager.h" +#include "storm/settings/modules/CoreSettings.h" +#include "storm/environment/solver/NativeSolverEnvironment.h" +#include "storm/environment/solver/GmmxxSolverEnvironment.h" +#include "storm/environment/solver/EigenSolverEnvironment.h" +#include "storm/environment/solver/LongRunAverageSolverEnvironment.h" + +namespace { + + enum class CtmcEngine {PrismSparse, JaniSparse, JaniHybrid}; + + class GBSparseGmmxxGmresIluEnvironment { + public: + static const storm::dd::DdType ddType = storm::dd::DdType::Sylvan; // unused for sparse models + static const CtmcEngine engine = CtmcEngine::PrismSparse; + static const bool isExact = false; + typedef double ValueType; + typedef storm::models::sparse::Ctmc ModelType; + static storm::Environment createEnvironment() { + storm::Environment env; + env.solver().setLinearEquationSolverType(storm::solver::EquationSolverType::Gmmxx); + env.solver().gmmxx().setMethod(storm::solver::GmmxxLinearEquationSolverMethod::Gmres); + env.solver().gmmxx().setPreconditioner(storm::solver::GmmxxLinearEquationSolverPreconditioner::Ilu); + env.solver().lra().setDetLraMethod(storm::solver::LraMethod::GainBiasEquations); + return env; + } + }; + + class GBJaniSparseGmmxxGmresIluEnvironment { + public: + static const storm::dd::DdType ddType = storm::dd::DdType::Sylvan; // unused for sparse models + static const CtmcEngine engine = CtmcEngine::JaniSparse; + static const bool isExact = false; + typedef double ValueType; + typedef storm::models::sparse::Ctmc ModelType; + static storm::Environment createEnvironment() { + storm::Environment env; + env.solver().setLinearEquationSolverType(storm::solver::EquationSolverType::Gmmxx); + env.solver().gmmxx().setMethod(storm::solver::GmmxxLinearEquationSolverMethod::Gmres); + env.solver().gmmxx().setPreconditioner(storm::solver::GmmxxLinearEquationSolverPreconditioner::Ilu); + env.solver().lra().setDetLraMethod(storm::solver::LraMethod::GainBiasEquations); + return env; + } + }; + + class GBJaniHybridCuddGmmxxGmresEnvironment { + public: + static const storm::dd::DdType ddType = storm::dd::DdType::CUDD; + static const CtmcEngine engine = CtmcEngine::JaniHybrid; + static const bool isExact = false; + typedef double ValueType; + typedef storm::models::symbolic::Ctmc ModelType; + static storm::Environment createEnvironment() { + storm::Environment env; + env.solver().setLinearEquationSolverType(storm::solver::EquationSolverType::Gmmxx); + env.solver().gmmxx().setMethod(storm::solver::GmmxxLinearEquationSolverMethod::Gmres); + env.solver().lra().setDetLraMethod(storm::solver::LraMethod::GainBiasEquations); + env.solver().lra().setPrecision(storm::utility::convertNumber(1e-8)); // Need to increase precision because eq sys yields incorrect results + return env; + } + }; + + class GBJaniHybridSylvanGmmxxGmresEnvironment { + public: + static const storm::dd::DdType ddType = storm::dd::DdType::Sylvan; + static const CtmcEngine engine = CtmcEngine::JaniHybrid; + static const bool isExact = false; + typedef double ValueType; + typedef storm::models::symbolic::Ctmc ModelType; + static storm::Environment createEnvironment() { + storm::Environment env; + env.solver().setLinearEquationSolverType(storm::solver::EquationSolverType::Gmmxx); + env.solver().gmmxx().setMethod(storm::solver::GmmxxLinearEquationSolverMethod::Gmres); + env.solver().lra().setDetLraMethod(storm::solver::LraMethod::GainBiasEquations); + env.solver().lra().setPrecision(storm::utility::convertNumber(1e-8)); // Need to increase precision because eq sys yields incorrect results + return env; + } + }; + + class GBSparseEigenDGmresEnvironment { + public: + static const storm::dd::DdType ddType = storm::dd::DdType::Sylvan; // unused for sparse models + static const CtmcEngine engine = CtmcEngine::PrismSparse; + static const bool isExact = false; + typedef double ValueType; + typedef storm::models::sparse::Ctmc ModelType; + static storm::Environment createEnvironment() { + storm::Environment env; + env.solver().setLinearEquationSolverType(storm::solver::EquationSolverType::Eigen); + env.solver().eigen().setMethod(storm::solver::EigenLinearEquationSolverMethod::DGmres); + env.solver().eigen().setPreconditioner(storm::solver::EigenLinearEquationSolverPreconditioner::Ilu); + env.solver().lra().setDetLraMethod(storm::solver::LraMethod::GainBiasEquations); + return env; + } + }; + + class GBSparseEigenRationalLUEnvironment { + public: + static const storm::dd::DdType ddType = storm::dd::DdType::Sylvan; // unused for sparse models + static const CtmcEngine engine = CtmcEngine::PrismSparse; + static const bool isExact = false; + typedef double ValueType; + typedef storm::models::sparse::Ctmc ModelType; + static storm::Environment createEnvironment() { + storm::Environment env; + env.solver().setLinearEquationSolverType(storm::solver::EquationSolverType::Eigen); + env.solver().eigen().setMethod(storm::solver::EigenLinearEquationSolverMethod::SparseLU); + env.solver().lra().setDetLraMethod(storm::solver::LraMethod::GainBiasEquations); + return env; + } + }; + + class GBSparseNativeSorEnvironment { + public: + static const storm::dd::DdType ddType = storm::dd::DdType::Sylvan; // unused for sparse models + static const CtmcEngine engine = CtmcEngine::PrismSparse; + static const bool isExact = false; + typedef double ValueType; + typedef storm::models::sparse::Ctmc ModelType; + static storm::Environment createEnvironment() { + storm::Environment env; + env.solver().setLinearEquationSolverType(storm::solver::EquationSolverType::Native); + env.solver().native().setMethod(storm::solver::NativeLinearEquationSolverMethod::SOR); + env.solver().native().setSorOmega(storm::utility::convertNumber(0.7)); // LRA computation fails for 0.9 + env.solver().native().setPrecision(storm::utility::convertNumber(1e-9)); + env.solver().native().setRelativeTerminationCriterion(false); + env.solver().lra().setDetLraMethod(storm::solver::LraMethod::GainBiasEquations); + return env; + } + }; + + class DistrSparseGmmxxGmresIluEnvironment { + public: + static const storm::dd::DdType ddType = storm::dd::DdType::Sylvan; // unused for sparse models + static const CtmcEngine engine = CtmcEngine::PrismSparse; + static const bool isExact = false; + typedef double ValueType; + typedef storm::models::sparse::Ctmc ModelType; + static storm::Environment createEnvironment() { + storm::Environment env; + env.solver().setLinearEquationSolverType(storm::solver::EquationSolverType::Gmmxx); + env.solver().gmmxx().setMethod(storm::solver::GmmxxLinearEquationSolverMethod::Gmres); + env.solver().gmmxx().setPreconditioner(storm::solver::GmmxxLinearEquationSolverPreconditioner::Ilu); + env.solver().lra().setDetLraMethod(storm::solver::LraMethod::LraDistributionEquations); + env.solver().lra().setPrecision(storm::utility::convertNumber(1e-8)); // Need to increase precision because eq sys yields incorrect results + return env; + } + }; + + class DistrSparseEigenDoubleLUEnvironment { + public: + static const storm::dd::DdType ddType = storm::dd::DdType::Sylvan; // unused for sparse models + static const CtmcEngine engine = CtmcEngine::PrismSparse; + static const bool isExact = false; + typedef double ValueType; + typedef storm::models::sparse::Ctmc ModelType; + static storm::Environment createEnvironment() { + storm::Environment env; + env.solver().setLinearEquationSolverType(storm::solver::EquationSolverType::Eigen); + env.solver().eigen().setMethod(storm::solver::EigenLinearEquationSolverMethod::DGmres); + env.solver().eigen().setPreconditioner(storm::solver::EigenLinearEquationSolverPreconditioner::Ilu); + env.solver().eigen().setPrecision(storm::utility::convertNumber(1e-8)); + env.solver().lra().setDetLraMethod(storm::solver::LraMethod::LraDistributionEquations); + return env; + } + }; + + class ValueIterationSparseEnvironment { + public: + static const storm::dd::DdType ddType = storm::dd::DdType::Sylvan; // unused for sparse models + static const CtmcEngine engine = CtmcEngine::PrismSparse; + static const bool isExact = false; + typedef double ValueType; + typedef storm::models::sparse::Ctmc ModelType; + static storm::Environment createEnvironment() { + storm::Environment env; + env.solver().lra().setDetLraMethod(storm::solver::LraMethod::ValueIteration); + return env; + } + }; + + template + class LraCtmcCslModelCheckerTest : public ::testing::Test { + public: + typedef typename TestType::ValueType ValueType; + typedef typename storm::models::sparse::Ctmc SparseModelType; + typedef typename storm::models::symbolic::Ctmc SymbolicModelType; + + LraCtmcCslModelCheckerTest() : _environment(TestType::createEnvironment()) {} + storm::Environment const& env() const { return _environment; } + ValueType parseNumber(std::string const& input) const { return storm::utility::convertNumber(input);} + ValueType precision() const { return TestType::isExact ? parseNumber("0") : parseNumber("1e-6");} + bool isSparseModel() const { return std::is_same::value; } + bool isSymbolicModel() const { return std::is_same::value; } + CtmcEngine getEngine() const { return TestType::engine; } + + template + typename std::enable_if::value, std::pair, std::vector>>>::type + buildModelFormulas(std::string const& pathToPrismFile, std::string const& formulasAsString, std::string const& constantDefinitionString = "") const { + std::pair, std::vector>> result; + storm::prism::Program program = storm::api::parseProgram(pathToPrismFile, true); + program = storm::utility::prism::preprocess(program, constantDefinitionString); + if (TestType::engine == CtmcEngine::PrismSparse) { + result.second = storm::api::extractFormulasFromProperties(storm::api::parsePropertiesForPrismProgram(formulasAsString, program)); + result.first = storm::api::buildSparseModel(program, result.second)->template as(); + } else if (TestType::engine == CtmcEngine::JaniSparse) { + auto janiData = storm::api::convertPrismToJani(program, storm::api::parsePropertiesForPrismProgram(formulasAsString, program)); + janiData.first.substituteFunctions(); + result.second = storm::api::extractFormulasFromProperties(janiData.second); + result.first = storm::api::buildSparseModel(janiData.first, result.second)->template as(); + } + return result; + } + + template + typename std::enable_if::value, std::pair, std::vector>>>::type + buildModelFormulas(std::string const& pathToPrismFile, std::string const& formulasAsString, std::string const& constantDefinitionString = "") const { + std::pair, std::vector>> result; + storm::prism::Program program = storm::api::parseProgram(pathToPrismFile, true); + program = storm::utility::prism::preprocess(program, constantDefinitionString); + if (TestType::engine == CtmcEngine::JaniHybrid) { + auto janiData = storm::api::convertPrismToJani(program, storm::api::parsePropertiesForPrismProgram(formulasAsString, program)); + janiData.first.substituteFunctions(); + result.second = storm::api::extractFormulasFromProperties(janiData.second); + result.first = storm::api::buildSymbolicModel(janiData.first, result.second)->template as(); + } + return result; + } + + std::vector> getTasks(std::vector> const& formulas) const { + std::vector> result; + for (auto const& f : formulas) { + result.emplace_back(*f); + } + return result; + } + + template + typename std::enable_if::value, std::shared_ptr>>::type + createModelChecker(std::shared_ptr const& model) const { + if (TestType::engine == CtmcEngine::PrismSparse || TestType::engine == CtmcEngine::JaniSparse) { + return std::make_shared>(*model); + } + return nullptr; + } + + template + typename std::enable_if::value, std::shared_ptr>>::type + createModelChecker(std::shared_ptr const& model) const { + if (TestType::engine == CtmcEngine::JaniHybrid) { + return std::make_shared>(*model); + } + return nullptr; + } + + bool getQualitativeResultAtInitialState(std::shared_ptr> const& model, std::unique_ptr& result) { + auto filter = getInitialStateFilter(model); + result->filter(*filter); + return result->asQualitativeCheckResult().forallTrue(); + } + + ValueType getQuantitativeResultAtInitialState(std::shared_ptr> const& model, std::unique_ptr& result) { + auto filter = getInitialStateFilter(model); + result->filter(*filter); + return result->asQuantitativeCheckResult().getMin(); + } + + private: + storm::Environment _environment; + + std::unique_ptr getInitialStateFilter(std::shared_ptr> const& model) const { + if (isSparseModel()) { + return std::make_unique(model->template as()->getInitialStates()); + } else { + return std::make_unique>(model->template as()->getReachableStates(), model->template as()->getInitialStates()); + } + } + }; + + typedef ::testing::Types< + GBSparseGmmxxGmresIluEnvironment, + GBJaniSparseGmmxxGmresIluEnvironment, + GBJaniHybridCuddGmmxxGmresEnvironment, + GBJaniHybridSylvanGmmxxGmresEnvironment, + GBSparseEigenDGmresEnvironment, + GBSparseEigenRationalLUEnvironment, + GBSparseNativeSorEnvironment, + DistrSparseGmmxxGmresIluEnvironment, + DistrSparseEigenDoubleLUEnvironment, + ValueIterationSparseEnvironment + > TestingTypes; + + TYPED_TEST_SUITE(LraCtmcCslModelCheckerTest, TestingTypes,); + + TYPED_TEST(LraCtmcCslModelCheckerTest, Cluster) { + std::string formulasString = "LRA=? [\"minimum\"]"; + + auto modelFormulas = this->buildModelFormulas(STORM_TEST_RESOURCES_DIR "/ctmc/cluster2.sm", formulasString); + auto model = std::move(modelFormulas.first); + auto tasks = this->getTasks(modelFormulas.second); + EXPECT_EQ(276ul, model->getNumberOfStates()); + EXPECT_EQ(1120ul, model->getNumberOfTransitions()); + ASSERT_EQ(model->getType(), storm::models::ModelType::Ctmc); + auto checker = this->createModelChecker(model); + std::unique_ptr result; + + result = checker->check(this->env(), tasks[0]); + EXPECT_NEAR(this->parseNumber("0.99999766034263426"), this->getQuantitativeResultAtInitialState(model, result), this->precision()); + + } + + TYPED_TEST(LraCtmcCslModelCheckerTest, Embedded) { + std::string formulasString = "LRA=? [\"fail_sensors\"]"; + + auto modelFormulas = this->buildModelFormulas(STORM_TEST_RESOURCES_DIR "/ctmc/embedded2.sm", formulasString); + auto model = std::move(modelFormulas.first); + auto tasks = this->getTasks(modelFormulas.second); + EXPECT_EQ(3478ul, model->getNumberOfStates()); + EXPECT_EQ(14639ul, model->getNumberOfTransitions()); + ASSERT_EQ(model->getType(), storm::models::ModelType::Ctmc); + auto checker = this->createModelChecker(model); + std::unique_ptr result; + + result = checker->check(this->env(), tasks[0]); + EXPECT_NEAR(this->parseNumber("0.93458866427696596"), this->getQuantitativeResultAtInitialState(model, result), this->precision()); + } + + TYPED_TEST(LraCtmcCslModelCheckerTest, Polling) { + std::string formulasString = "LRA=?[\"target\"]"; + + auto modelFormulas = this->buildModelFormulas(STORM_TEST_RESOURCES_DIR "/ctmc/polling2.sm", formulasString); + auto model = std::move(modelFormulas.first); + auto tasks = this->getTasks(modelFormulas.second); + EXPECT_EQ(12ul, model->getNumberOfStates()); + EXPECT_EQ(22ul, model->getNumberOfTransitions()); + ASSERT_EQ(model->getType(), storm::models::ModelType::Ctmc); + auto checker = this->createModelChecker(model); + std::unique_ptr result; + + result = checker->check(this->env(), tasks[0]); + EXPECT_NEAR(this->parseNumber("0.20079750055570736"), this->getQuantitativeResultAtInitialState(model, result), this->precision()); + + } + + TYPED_TEST(LraCtmcCslModelCheckerTest, Tandem) { + std::string formulasString = "LRA=? [\"first_queue_full\"]"; + + auto modelFormulas = this->buildModelFormulas(STORM_TEST_RESOURCES_DIR "/ctmc/tandem5.sm", formulasString); + auto model = std::move(modelFormulas.first); + auto tasks = this->getTasks(modelFormulas.second); + EXPECT_EQ(66ul, model->getNumberOfStates()); + EXPECT_EQ(189ul, model->getNumberOfTransitions()); + ASSERT_EQ(model->getType(), storm::models::ModelType::Ctmc); + auto checker = this->createModelChecker(model); + std::unique_ptr result; + + result = checker->check(this->env(), tasks[0]); + EXPECT_NEAR(this->parseNumber("0.9100373532"), this->getQuantitativeResultAtInitialState(model, result), this->precision()); + } + + + TYPED_TEST(LraCtmcCslModelCheckerTest, Rewards) { + std::string formulasString = "R=? [ LRA ]"; + + auto modelFormulas = this->buildModelFormulas(STORM_TEST_RESOURCES_DIR "/ctmc/lrarewards.sm", formulasString); + auto model = std::move(modelFormulas.first); + auto tasks = this->getTasks(modelFormulas.second); + EXPECT_EQ(4ul, model->getNumberOfStates()); + EXPECT_EQ(6ul, model->getNumberOfTransitions()); + ASSERT_EQ(model->getType(), storm::models::ModelType::Ctmc); + auto checker = this->createModelChecker(model); + std::unique_ptr result; + + result = checker->check(this->env(), tasks[0]); + EXPECT_NEAR(this->parseNumber("11/15"), this->getQuantitativeResultAtInitialState(model, result), this->precision()); + } + +} From 86506e2b2500c365c6c4fde82ac5d87c621a7be8 Mon Sep 17 00:00:00 2001 From: Tim Quatmann Date: Tue, 3 Dec 2019 13:55:58 +0100 Subject: [PATCH 08/24] Added LRA distribution equation system for computing the LRA of a BSCC. Fixes in the gain/bias characterization. --- .../csl/helper/SparseCtmcCslHelper.cpp | 252 +++++++++++++----- .../csl/helper/SparseCtmcCslHelper.h | 8 +- 2 files changed, 189 insertions(+), 71 deletions(-) diff --git a/src/storm/modelchecker/csl/helper/SparseCtmcCslHelper.cpp b/src/storm/modelchecker/csl/helper/SparseCtmcCslHelper.cpp index 4d381e885..528182ea5 100644 --- a/src/storm/modelchecker/csl/helper/SparseCtmcCslHelper.cpp +++ b/src/storm/modelchecker/csl/helper/SparseCtmcCslHelper.cpp @@ -451,15 +451,23 @@ namespace storm { std::vector bsccLra; bsccLra.reserve(bsccDecomposition.size()); + auto underlyingSolverEnvironment = env; + auto precision = env.solver().lra().getPrecision(); + if (env.solver().isForceSoundness()) { + // For sound computations, the error in the MECS plus the error in the remaining system should be less then the user defined precision. + precision /= storm::utility::convertNumber(2); + underlyingSolverEnvironment.solver().lra().setPrecision(precision); + } + underlyingSolverEnvironment.solver().setLinearEquationSolverPrecision(precision, env.solver().lra().getRelativeTerminationCriterion()); + // Keep track of the maximal and minimal value occuring in one of the BSCCs ValueType maxValue, minValue; storm::storage::BitVector statesInBsccs(numberOfStates); - auto backwardTransitions = rateMatrix.transpose(); for (auto const& bscc : bsccDecomposition) { for (auto const& state : bscc) { statesInBsccs.set(state); } - bsccLra.push_back(computeLongRunAveragesForBscc(env, bscc, rateMatrix, backwardTransitions, valueGetter, exitRateVector)); + bsccLra.push_back(computeLongRunAveragesForBscc(underlyingSolverEnvironment, bscc, rateMatrix, valueGetter, exitRateVector)); if (bsccLra.size() == 1) { maxValue = bsccLra.back(); minValue = bsccLra.back(); @@ -505,7 +513,7 @@ namespace storm { // Compute reachability rewards storm::solver::GeneralLinearEquationSolverFactory linearEquationSolverFactory; - bool isEqSysFormat = linearEquationSolverFactory.getEquationProblemFormat(env) == storm::solver::LinearEquationSolverProblemFormat::EquationSystem; + bool isEqSysFormat = linearEquationSolverFactory.getEquationProblemFormat(underlyingSolverEnvironment) == storm::solver::LinearEquationSolverProblemFormat::EquationSystem; storm::storage::SparseMatrix rewardEquationSystemMatrix = rateMatrix.getSubmatrix(false, statesNotInBsccs, statesNotInBsccs, isEqSysFormat); if (exitRateVector) { uint64_t localRow = 0; @@ -520,14 +528,12 @@ namespace storm { rewardEquationSystemMatrix.convertToEquationSystem(); } rewardSolution = std::vector(rewardEquationSystemMatrix.getColumnCount(), (maxValue + minValue) / storm::utility::convertNumber(2)); - // std::cout << rewardEquationSystemMatrix << std::endl; - // std::cout << storm::utility::vector::toString(rewardRightSide) << std::endl; - std::unique_ptr> solver = linearEquationSolverFactory.create(env, std::move(rewardEquationSystemMatrix)); + std::unique_ptr> solver = linearEquationSolverFactory.create(underlyingSolverEnvironment, std::move(rewardEquationSystemMatrix)); solver->setBounds(minValue, maxValue); // Check solver requirements - auto requirements = solver->getRequirements(env); + auto requirements = solver->getRequirements(underlyingSolverEnvironment); STORM_LOG_THROW(!requirements.hasEnabledCriticalRequirement(), storm::exceptions::UncheckedRequirementException, "Solver requirements " + requirements.getEnabledRequirementsAsString() + " not checked."); - solver->solveEquations(env, rewardSolution, rewardRightSide); + solver->solveEquations(underlyingSolverEnvironment, rewardSolution, rewardRightSide); } // Fill the result vector. @@ -553,13 +559,7 @@ namespace storm { } template - ValueType SparseCtmcCslHelper::computeLongRunAveragesForBscc(Environment const& env, storm::storage::StronglyConnectedComponent const& bscc, storm::storage::SparseMatrix const& rateMatrix, storm::storage::SparseMatrix const& backwardTransitions, std::function const& valueGetter, std::vector const* exitRateVector) { - std::cout << std::endl << "========== BSCC =======" << std::endl; - storm::storage::BitVector bsccStates(rateMatrix.getRowCount()); - for (auto const& s : bscc) { bsccStates.set(s); std::cout << s << ": " << valueGetter(s) << "\t";} - auto subm = rateMatrix.getSubmatrix(false,bsccStates,bsccStates); - std::cout << std::endl << subm << std::endl; - + ValueType SparseCtmcCslHelper::computeLongRunAveragesForBscc(Environment const& env, storm::storage::StronglyConnectedComponent const& bscc, storm::storage::SparseMatrix const& rateMatrix, std::function const& valueGetter, std::vector const* exitRateVector) { // Catch the case where all values are the same (this includes the special case where the bscc is of size 1) auto it = bscc.begin(); @@ -571,28 +571,28 @@ namespace storm { } if (it == bscc.end()) { // All entries have the same LRA - std::cout << std::endl << "==== TRIVIAL: RESULT=" << val << "=======" << std::endl; return val; } - std::cout << std::endl << "========== SOLVING =======" << std::endl; storm::solver::LraMethod method = env.solver().lra().getDetLraMethod(); if (method == storm::solver::LraMethod::ValueIteration) { - return computeLongRunAveragesForBsccVi(env, bscc, rateMatrix, backwardTransitions, valueGetter, exitRateVector); + return computeLongRunAveragesForBsccVi(env, bscc, rateMatrix, valueGetter, exitRateVector); } else if (method == storm::solver::LraMethod::LraDistributionEquations) { - //return computeLongRunAveragesForBsccLraDistrEqSys(env, bscc, rateMatrix, backwardTransitions, valueGetter, exitRateVector); + // We only need the first element of the pair as the lra distribution is not relevant at this point. + return computeLongRunAveragesForBsccLraDistr(env, bscc, rateMatrix, valueGetter, exitRateVector).first; } STORM_LOG_WARN_COND(method == storm::solver::LraMethod::GainBiasEquations, "Unsupported lra method selected. Defaulting to gain-bias-equations."); - return computeLongRunAveragesForBsccEqSys(env, bscc, rateMatrix, backwardTransitions, valueGetter, exitRateVector); + // We don't need the bias values + return computeLongRunAveragesForBsccGainBias(env, bscc, rateMatrix, valueGetter, exitRateVector).first; } template <> - storm::RationalFunction SparseCtmcCslHelper::computeLongRunAveragesForBsccVi(Environment const&, storm::storage::StronglyConnectedComponent const&, storm::storage::SparseMatrix const&, storm::storage::SparseMatrix const&, std::function const&, std::vector const*) { + storm::RationalFunction SparseCtmcCslHelper::computeLongRunAveragesForBsccVi(Environment const&, storm::storage::StronglyConnectedComponent const&, storm::storage::SparseMatrix const&, std::function const&, std::vector const*) { STORM_LOG_THROW(false, storm::exceptions::NotSupportedException, "The requested Method for LRA computation is not supported for parametric models."); } template - ValueType SparseCtmcCslHelper::computeLongRunAveragesForBsccVi(Environment const& env, storm::storage::StronglyConnectedComponent const& bscc, storm::storage::SparseMatrix const& rateMatrix, storm::storage::SparseMatrix const& backwardTransitions, std::function const& valueGetter, std::vector const* exitRateVector) { + ValueType SparseCtmcCslHelper::computeLongRunAveragesForBsccVi(Environment const& env, storm::storage::StronglyConnectedComponent const& bscc, storm::storage::SparseMatrix const& rateMatrix, std::function const& valueGetter, std::vector const* exitRateVector) { // Initialize data about the bscc storm::storage::BitVector bsccStates(rateMatrix.getRowGroupCount(), false); @@ -678,24 +678,24 @@ namespace storm { if ((maxDiff - minDiff) <= (relative ? (precision * minDiff) : precision)) { break; } - + } + if (maxIter.is_initialized() && iter == maxIter.get()) { + STORM_LOG_WARN("LRA computation did not converge within " << iter << " iterations."); + } else { + STORM_LOG_TRACE("LRA computation converged after " << iter << " iterations."); } return (maxDiff + minDiff) * uniformizationRate / (storm::utility::convertNumber(2.0)); } - template <> - storm::RationalFunction SparseCtmcCslHelper::computeLongRunAveragesForBsccEqSys(Environment const&, storm::storage::StronglyConnectedComponent const&, storm::storage::SparseMatrix const&, storm::storage::SparseMatrix const&, std::function const&, std::vector const*) { - STORM_LOG_THROW(false, storm::exceptions::NotSupportedException, "The requested Method for LRA computation is not supported for parametric models."); - // This is not possible due to uniformization - } - template - ValueType SparseCtmcCslHelper::computeLongRunAveragesForBsccEqSys(Environment const& env, storm::storage::StronglyConnectedComponent const& bscc, storm::storage::SparseMatrix const& rateMatrix, storm::storage::SparseMatrix const& backwardTransitions, std::function const& valueGetter, std::vector const* exitRateVector) { + std::pair> SparseCtmcCslHelper::computeLongRunAveragesForBsccGainBias(Environment const& env, storm::storage::StronglyConnectedComponent const& bscc, storm::storage::SparseMatrix const& rateMatrix, std::function const& valueGetter, std::vector const* exitRateVector) { // We build the equation system as in Line 3 of Algorithm 3 from // Kretinsky, Meggendorfer: Efficient Strategy Iteration for Mean Payoff in Markov Decision Processes (ATVA 2017) // The first variable corresponds to the gain of the bscc whereas the subsequent variables yield the bias for each state s_1, s_2, .... // No bias variable for s_0 is needed since it is always set to zero, yielding an nxn equation system matrix - // To properly deal with CMTC, we also need to uniformize the transitions. + // To make this work for CTMC, we could uniformize the model. This preserves LRA and ensures that we can compute the + // LRA as for a DTMC (the soujourn time in each state is the same). If we then multiply the equations with the uniformization rate, + // the uniformization rate cancels out. Hence, we obtain the equation system below. // Get a mapping from global state indices to local ones. std::unordered_map toLocalIndexMap; @@ -705,12 +705,6 @@ namespace storm { ++localIndex; } - // Get the uniformization rate - ValueType uniformizationRate = storm::utility::one(); - if (exitRateVector) { - uniformizationRate = *std::max_element(exitRateVector->begin(), exitRateVector->end()); - } - // Build the equation system matrix and vector. storm::solver::GeneralLinearEquationSolverFactory linearEquationSolverFactory; bool isEquationSystemFormat = linearEquationSolverFactory.getEquationProblemFormat(env) == storm::solver::LinearEquationSolverProblemFormat::EquationSystem; @@ -730,14 +724,10 @@ namespace storm { builder.addNextValue(row, 0, -storm::utility::one()); } // Compute weighted sum over successor state. As this is a BSCC, each successor state will again be in the BSCC. - bool needDiagonalEntry = ((exitRateVector && (*exitRateVector)[globalState] != uniformizationRate) || isEquationSystemFormat) && (row > 0); + bool needDiagonalEntry = isEquationSystemFormat && (row > 0); ValueType diagonalValue; if (needDiagonalEntry) { - if (isEquationSystemFormat) { - diagonalValue = (*exitRateVector)[globalState] / uniformizationRate; // 1 - (1 - r(s)/unifrate) - } else { - diagonalValue = storm::utility::one() - (*exitRateVector)[globalState] / uniformizationRate; - } + diagonalValue = exitRateVector ? (*exitRateVector)[globalState] : storm::utility::one(); } for (auto const& entry : rateMatrix.getRow(globalState)) { uint64_t col = toLocalIndexMap[entry.getColumn()]; @@ -746,9 +736,6 @@ namespace storm { continue; } entryValue = entry.getValue(); - if (exitRateVector) { - // entryValue /= uniformizationRate; - } if (isEquationSystemFormat) { entryValue = -entryValue; } @@ -765,31 +752,160 @@ namespace storm { if (needDiagonalEntry) { builder.addNextValue(row, row, diagonalValue); } - if (exitRateVector) { - eqsysVector.push_back(valueGetter(globalState) / uniformizationRate); - } else { - eqsysVector.push_back(valueGetter(globalState)); - } + eqsysVector.push_back(valueGetter(globalState)); ++row; } // Create a linear equation solver - auto matrix = builder.build(); - std::cout << matrix << std::endl; - std::cout << "RHS: " << storm::utility::vector::toString(eqsysVector) << std::endl; - auto solver = linearEquationSolverFactory.create(env, std::move(matrix)); + auto subEnv = env; + subEnv.solver().setLinearEquationSolverPrecision(env.solver().lra().getPrecision(), env.solver().lra().getRelativeTerminationCriterion()); + auto solver = linearEquationSolverFactory.create(subEnv, builder.build()); // Check solver requirements. - auto requirements = solver->getRequirements(env); + auto requirements = solver->getRequirements(subEnv); STORM_LOG_THROW(!requirements.hasEnabledCriticalRequirement(), storm::exceptions::UncheckedRequirementException, "Solver requirements " + requirements.getEnabledRequirementsAsString() + " not checked."); // Todo: Find bounds on the bias variable. Just inserting the maximal value from the vector probably does not work. std::vector eqSysSol(bscc.size(), storm::utility::zero()); // Take the mean of the rewards as an initial guess for the gain eqSysSol.front() = std::accumulate(eqsysVector.begin(), eqsysVector.end(), storm::utility::zero()) / storm::utility::convertNumber(bscc.size()); - solver->solveEquations(env, eqSysSol, eqsysVector); + solver->solveEquations(subEnv, eqSysSol, eqsysVector); + + ValueType gain = eqSysSol.front(); + // insert bias value for state 0 + eqSysSol.front() = storm::utility::zero(); + // Return the gain and the bias values + return std::pair>(std::move(gain), std::move(eqSysSol)); + } + + template + std::pair> SparseCtmcCslHelper::computeLongRunAveragesForBsccLraDistr(Environment const& env, storm::storage::StronglyConnectedComponent const& bscc, storm::storage::SparseMatrix const& rateMatrix, std::function const& valueGetter, std::vector const* exitRateVector) { + + // Let A be ab auxiliary Matrix with A[s,s] = R(s,s) - r(s) & A[s,s'] = R(s,s') for s,s' in BSCC and s!=s'. + // We build and solve the equation system for + // x*A=0 & x_0+...+x_n=1 <=> A^t*x=0=x-x & x_0+...+x_n=1 <=> (1+A^t)*x = x & 1-x_0-...-x_n-1=x_n + // Then, x[i] will be the fraction of the time we are in state i. + + // This method assumes that this BSCC consist of more than one state + if (bscc.size() == 1) { + return { storm::utility::one(), {storm::utility::one()} }; + } - std::cout << std::endl << "========== Result = " << eqSysSol.front() << "*" << uniformizationRate << " =======" << std::endl; - return eqSysSol.front();// * uniformizationRate; + // Get a mapping from global state indices to local ones as well as a bitvector containing states within the BSCC. + std::unordered_map toLocalIndexMap; + storm::storage::BitVector bsccStates(rateMatrix.getRowCount(), false); + uint64_t localIndex = 0; + for (auto const& globalIndex : bscc) { + bsccStates.set(globalIndex, true); + toLocalIndexMap[globalIndex] = localIndex; + ++localIndex; + } + + // Build the auxiliary Matrix A. + auto auxMatrix = rateMatrix.getSubmatrix(false, bsccStates, bsccStates, true); // add diagonal entries! + uint64_t row = 0; + for (auto const& globalIndex : bscc) { + for (auto& entry : auxMatrix.getRow(row)) { + if (entry.getColumn() == row) { + // This value is non-zero since we have a BSCC with more than one state + if (exitRateVector) { + entry.setValue(entry.getValue() - (*exitRateVector)[globalIndex]); + } else { + entry.setValue(entry.getValue() - storm::utility::one()); + } + } + } + ++row; + } + assert(row == auxMatrix.getRowCount()); + + // We need to consider A^t. This will not delete diagonal entries since they are non-zero. + auxMatrix = auxMatrix.transpose(); + + // Check whether we need the fixpoint characterization + storm::solver::GeneralLinearEquationSolverFactory linearEquationSolverFactory; + bool isFixpointFormat = linearEquationSolverFactory.getEquationProblemFormat(env) == storm::solver::LinearEquationSolverProblemFormat::FixedPointSystem; + + if (isFixpointFormat) { + // Add a 1 on the diagonal + for (row = 0; row < auxMatrix.getRowCount(); ++row) { + for (auto& entry : auxMatrix) { + if (entry.getColumn() == row) { + entry.setValue(storm::utility::one() + entry.getValue()); + } + } + } + } + + // We now build the equation system matrix. + // We can drop the last row of A and add ones in this row instead to assert that the variables sum up to one + // Phase 1: replace the existing entries of the last row with ones + uint64_t col = 0; + uint64_t lastRow = auxMatrix.getRowCount() - 1; + for (auto& entry : auxMatrix.getRow(lastRow)) { + entry.setColumn(col); + if (isFixpointFormat) { + if (col == lastRow) { + entry.setValue(storm::utility::zero()); + } else { + entry.setValue(-storm::utility::one()); + } + } else { + entry.setValue(storm::utility::one()); + } + ++col; + } + storm::storage::SparseMatrixBuilder builder(std::move(auxMatrix)); + for (; col <= lastRow; ++col) { + if (isFixpointFormat) { + if (col != lastRow) { + builder.addNextValue(lastRow, col, -storm::utility::one()); + } + } else { + builder.addNextValue(lastRow, col, storm::utility::one()); + } + } + + std::vector bsccEquationSystemRightSide(bscc.size(), storm::utility::zero()); + bsccEquationSystemRightSide.back() = storm::utility::one(); + + // Create a linear equation solver + auto subEnv = env; + subEnv.solver().setLinearEquationSolverPrecision(env.solver().lra().getPrecision(), env.solver().lra().getRelativeTerminationCriterion()); + auto solver = linearEquationSolverFactory.create(subEnv, builder.build()); + solver->setBounds(storm::utility::zero(), storm::utility::one()); + // Check solver requirements. + auto requirements = solver->getRequirements(subEnv); + STORM_LOG_THROW(!requirements.hasEnabledCriticalRequirement(), storm::exceptions::UncheckedRequirementException, "Solver requirements " + requirements.getEnabledRequirementsAsString() + " not checked."); + + std::vector lraDistr(bscc.size(), storm::utility::one() / storm::utility::convertNumber(bscc.size())); + solver->solveEquations(subEnv, lraDistr, bsccEquationSystemRightSide); + // If exit rates were given, we need to 'fix' the results to also account for the timing behaviour. + if (false && exitRateVector != nullptr) { + ValueType totalValue = storm::utility::zero(); + auto solIt = lraDistr.begin(); + for (auto const& globalState : bscc) { + totalValue += (*solIt) * (storm::utility::one() / (*exitRateVector)[globalState]); + ++solIt; + } + assert(solIt == lraDistr.end()); + solIt = lraDistr.begin(); + for (auto const& globalState : bscc) { + *solIt = ((*solIt) * (storm::utility::one() / (*exitRateVector)[globalState])) / totalValue; + ++solIt; + } + assert(solIt == lraDistr.end()); + } + + // Calculate final LRA Value + ValueType result = storm::utility::zero(); + auto solIt = lraDistr.begin(); + for (auto const& globalState : bscc) { + result += valueGetter(globalState) * (*solIt); + ++solIt; + } + assert(solIt == lraDistr.end()); + + return std::pair>(std::move(result), std::move(lraDistr)); } template ::SupportsExponential, int>::type> @@ -1012,9 +1128,9 @@ namespace storm { template std::vector SparseCtmcCslHelper::computeTotalRewards(Environment const& env, storm::solver::SolveGoal&& goal, storm::storage::SparseMatrix const& rateMatrix, storm::storage::SparseMatrix const& backwardTransitions, std::vector const& exitRateVector, storm::models::sparse::StandardRewardModel const& rewardModel, bool qualitative); - template std::vector SparseCtmcCslHelper::computeLongRunAverageProbabilities(Environment const& env, storm::solver::SolveGoal&& goal, storm::storage::SparseMatrix const& probabilityMatrix, storm::storage::BitVector const& psiStates, std::vector const* exitRateVector); - template std::vector SparseCtmcCslHelper::computeLongRunAverageRewards(Environment const& env, storm::solver::SolveGoal&& goal, storm::storage::SparseMatrix const& probabilityMatrix, storm::models::sparse::StandardRewardModel const& rewardModel, std::vector const* exitRateVector); - template std::vector SparseCtmcCslHelper::computeLongRunAverageRewards(Environment const& env, storm::solver::SolveGoal&& goal, storm::storage::SparseMatrix const& probabilityMatrix, std::vector const& stateRewardVector, std::vector const* exitRateVector); + template std::vector SparseCtmcCslHelper::computeLongRunAverageProbabilities(Environment const& env, storm::solver::SolveGoal&& goal, storm::storage::SparseMatrix const& rateMatrix, storm::storage::BitVector const& psiStates, std::vector const* exitRateVector); + template std::vector SparseCtmcCslHelper::computeLongRunAverageRewards(Environment const& env, storm::solver::SolveGoal&& goal, storm::storage::SparseMatrix const& rateMatrix, storm::models::sparse::StandardRewardModel const& rewardModel, std::vector const* exitRateVector); + template std::vector SparseCtmcCslHelper::computeLongRunAverageRewards(Environment const& env, storm::solver::SolveGoal&& goal, storm::storage::SparseMatrix const& rateMatrix, std::vector const& stateRewardVector, std::vector const* exitRateVector); template std::vector SparseCtmcCslHelper::computeCumulativeRewards(Environment const& env, storm::solver::SolveGoal&& goal, storm::storage::SparseMatrix const& rateMatrix, std::vector const& exitRateVector, storm::models::sparse::StandardRewardModel const& rewardModel, double timeBound); @@ -1049,14 +1165,14 @@ namespace storm { template std::vector SparseCtmcCslHelper::computeTotalRewards(Environment const& env, storm::solver::SolveGoal&& goal, storm::storage::SparseMatrix const& rateMatrix, storm::storage::SparseMatrix const& backwardTransitions, std::vector const& exitRateVector, storm::models::sparse::StandardRewardModel const& rewardModel, bool qualitative); template std::vector SparseCtmcCslHelper::computeTotalRewards(Environment const& env, storm::solver::SolveGoal&& goal, storm::storage::SparseMatrix const& rateMatrix, storm::storage::SparseMatrix const& backwardTransitions, std::vector const& exitRateVector, storm::models::sparse::StandardRewardModel const& rewardModel, bool qualitative); - template std::vector SparseCtmcCslHelper::computeLongRunAverageProbabilities(Environment const& env, storm::solver::SolveGoal&& goal, storm::storage::SparseMatrix const& probabilityMatrix, storm::storage::BitVector const& psiStates, std::vector const* exitRateVector); - template std::vector SparseCtmcCslHelper::computeLongRunAverageProbabilities(Environment const& env, storm::solver::SolveGoal&& goal, storm::storage::SparseMatrix const& probabilityMatrix, storm::storage::BitVector const& psiStates, std::vector const* exitRateVector); + template std::vector SparseCtmcCslHelper::computeLongRunAverageProbabilities(Environment const& env, storm::solver::SolveGoal&& goal, storm::storage::SparseMatrix const& rateMatrix, storm::storage::BitVector const& psiStates, std::vector const* exitRateVector); + template std::vector SparseCtmcCslHelper::computeLongRunAverageProbabilities(Environment const& env, storm::solver::SolveGoal&& goal, storm::storage::SparseMatrix const& rateMatrix, storm::storage::BitVector const& psiStates, std::vector const* exitRateVector); - template std::vector SparseCtmcCslHelper::computeLongRunAverageRewards(Environment const& env, storm::solver::SolveGoal&& goal, storm::storage::SparseMatrix const& probabilityMatrix, storm::models::sparse::StandardRewardModel const& rewardModel, std::vector const* exitRateVector); - template std::vector SparseCtmcCslHelper::computeLongRunAverageRewards(Environment const& env, storm::solver::SolveGoal&& goal, storm::storage::SparseMatrix const& probabilityMatrix, storm::models::sparse::StandardRewardModel const& rewardModel, std::vector const* exitRateVector); + template std::vector SparseCtmcCslHelper::computeLongRunAverageRewards(Environment const& env, storm::solver::SolveGoal&& goal, storm::storage::SparseMatrix const& rateMatrix, storm::models::sparse::StandardRewardModel const& rewardModel, std::vector const* exitRateVector); + template std::vector SparseCtmcCslHelper::computeLongRunAverageRewards(Environment const& env, storm::solver::SolveGoal&& goal, storm::storage::SparseMatrix const& rateMatrix, storm::models::sparse::StandardRewardModel const& rewardModel, std::vector const* exitRateVector); - template std::vector SparseCtmcCslHelper::computeLongRunAverageRewards(Environment const& env, storm::solver::SolveGoal&& goal, storm::storage::SparseMatrix const& probabilityMatrix, std::vector const& stateRewardVector, std::vector const* exitRateVector); - template std::vector SparseCtmcCslHelper::computeLongRunAverageRewards(Environment const& env, storm::solver::SolveGoal&& goal, storm::storage::SparseMatrix const& probabilityMatrix, std::vector const& stateRewardVector, std::vector const* exitRateVector); + template std::vector SparseCtmcCslHelper::computeLongRunAverageRewards(Environment const& env, storm::solver::SolveGoal&& goal, storm::storage::SparseMatrix const& rateMatrix, std::vector const& stateRewardVector, std::vector const* exitRateVector); + template std::vector SparseCtmcCslHelper::computeLongRunAverageRewards(Environment const& env, storm::solver::SolveGoal&& goal, storm::storage::SparseMatrix const& rateMatrix, std::vector const& stateRewardVector, std::vector const* exitRateVector); template std::vector SparseCtmcCslHelper::computeCumulativeRewards(Environment const& env, storm::solver::SolveGoal&& goal, storm::storage::SparseMatrix const& rateMatrix, std::vector const& exitRateVector, storm::models::sparse::StandardRewardModel const& rewardModel, double timeBound); template std::vector SparseCtmcCslHelper::computeCumulativeRewards(Environment const& env, storm::solver::SolveGoal&& goal, storm::storage::SparseMatrix const& rateMatrix, std::vector const& exitRateVector, storm::models::sparse::StandardRewardModel const& rewardModel, double timeBound); diff --git a/src/storm/modelchecker/csl/helper/SparseCtmcCslHelper.h b/src/storm/modelchecker/csl/helper/SparseCtmcCslHelper.h index 8a45ea75f..95840ae19 100644 --- a/src/storm/modelchecker/csl/helper/SparseCtmcCslHelper.h +++ b/src/storm/modelchecker/csl/helper/SparseCtmcCslHelper.h @@ -125,11 +125,13 @@ namespace storm { template static std::vector computeLongRunAverages(Environment const& env, storm::solver::SolveGoal&& goal, storm::storage::SparseMatrix const& rateMatrix, std::function const& valueGetter, std::vector const* exitRateVector); template - static ValueType computeLongRunAveragesForBscc(Environment const& env, storm::storage::StronglyConnectedComponent const& bscc, storm::storage::SparseMatrix const& rateMatrix, storm::storage::SparseMatrix const& backwardTransitions, std::function const& valueGetter, std::vector const* exitRateVector); + static ValueType computeLongRunAveragesForBscc(Environment const& env, storm::storage::StronglyConnectedComponent const& bscc, storm::storage::SparseMatrix const& rateMatrix, std::function const& valueGetter, std::vector const* exitRateVector); template - static ValueType computeLongRunAveragesForBsccVi(Environment const& env, storm::storage::StronglyConnectedComponent const& bscc, storm::storage::SparseMatrix const& rateMatrix, storm::storage::SparseMatrix const& backwardTransitions, std::function const& valueGetter, std::vector const* exitRateVector); + static ValueType computeLongRunAveragesForBsccVi(Environment const& env, storm::storage::StronglyConnectedComponent const& bscc, storm::storage::SparseMatrix const& rateMatrix, std::function const& valueGetter, std::vector const* exitRateVector); template - static ValueType computeLongRunAveragesForBsccEqSys(Environment const& env, storm::storage::StronglyConnectedComponent const& bscc, storm::storage::SparseMatrix const& rateMatrix, storm::storage::SparseMatrix const& backwardTransitions, std::function const& valueGetter, std::vector const* exitRateVector); + static std::pair> computeLongRunAveragesForBsccGainBias(Environment const& env, storm::storage::StronglyConnectedComponent const& bscc, storm::storage::SparseMatrix const& rateMatrix, std::function const& valueGetter, std::vector const* exitRateVector); + template + static std::pair> computeLongRunAveragesForBsccLraDistr(Environment const& env, storm::storage::StronglyConnectedComponent const& bscc, storm::storage::SparseMatrix const& rateMatrix, std::function const& valueGetter, std::vector const* exitRateVector); }; } From ea04f6dcd2a32fbf75164de1c6a999b8d4e0cdff Mon Sep 17 00:00:00 2001 From: Tim Quatmann Date: Tue, 3 Dec 2019 14:29:29 +0100 Subject: [PATCH 09/24] Fixes for LRA computation. --- .../csl/helper/SparseCtmcCslHelper.cpp | 82 +++++++++++-------- .../modules/LongRunAverageSolverSettings.cpp | 2 +- src/storm/solver/SolverSelectionOptions.cpp | 8 +- src/storm/storage/SparseMatrix.cpp | 8 ++ .../csl/LraCtmcCslModelCheckerTest.cpp | 3 +- 5 files changed, 62 insertions(+), 41 deletions(-) diff --git a/src/storm/modelchecker/csl/helper/SparseCtmcCslHelper.cpp b/src/storm/modelchecker/csl/helper/SparseCtmcCslHelper.cpp index 528182ea5..6d82a6d79 100644 --- a/src/storm/modelchecker/csl/helper/SparseCtmcCslHelper.cpp +++ b/src/storm/modelchecker/csl/helper/SparseCtmcCslHelper.cpp @@ -15,6 +15,7 @@ #include "storm/adapters/RationalFunctionAdapter.h" #include "storm/environment/solver/LongRunAverageSolverEnvironment.h" +#include "storm/environment/solver/TopologicalSolverEnvironment.h" #include "storm/utility/macros.h" #include "storm/utility/vector.h" @@ -575,13 +576,21 @@ namespace storm { } storm::solver::LraMethod method = env.solver().lra().getDetLraMethod(); + if (storm::NumberTraits::IsExact && env.solver().lra().isDetLraMethodSetFromDefault() && method == storm::solver::LraMethod::ValueIteration) { + method = storm::solver::LraMethod::GainBiasEquations; + STORM_LOG_INFO("Selecting " << storm::solver::toString(method) << " as the solution technique for long-run properties to guarantee exact results. If you want to override this, please explicitly specify a different LRA method."); + } else if (env.solver().isForceSoundness() && env.solver().lra().isDetLraMethodSetFromDefault() && method != storm::solver::LraMethod::ValueIteration) { + method = storm::solver::LraMethod::ValueIteration; + STORM_LOG_INFO("Selecting " << storm::solver::toString(method) << " as the solution technique for long-run properties to guarantee sound results. If you want to override this, please explicitly specify a different LRA method."); + } + STORM_LOG_TRACE("Computing LRA for BSCC of size " << bscc.size() << " using '" << storm::solver::toString(method) << "'."); if (method == storm::solver::LraMethod::ValueIteration) { return computeLongRunAveragesForBsccVi(env, bscc, rateMatrix, valueGetter, exitRateVector); } else if (method == storm::solver::LraMethod::LraDistributionEquations) { // We only need the first element of the pair as the lra distribution is not relevant at this point. return computeLongRunAveragesForBsccLraDistr(env, bscc, rateMatrix, valueGetter, exitRateVector).first; } - STORM_LOG_WARN_COND(method == storm::solver::LraMethod::GainBiasEquations, "Unsupported lra method selected. Defaulting to gain-bias-equations."); + STORM_LOG_WARN_COND(method == storm::solver::LraMethod::GainBiasEquations, "Unsupported lra method selected. Defaulting to " << storm::solver::toString(storm::solver::LraMethod::GainBiasEquations) << "."); // We don't need the bias values return computeLongRunAveragesForBsccGainBias(env, bscc, rateMatrix, valueGetter, exitRateVector).first; } @@ -705,12 +714,20 @@ namespace storm { ++localIndex; } + // Prepare an environment for the underlying equation solver + auto subEnv = env; + if (subEnv.solver().getLinearEquationSolverType() == storm::solver::EquationSolverType::Topological) { + // Topological solver does not make any sense since the BSCC is connected. + subEnv.solver().setLinearEquationSolverType(subEnv.solver().topological().getUnderlyingEquationSolverType()); + } + subEnv.solver().setLinearEquationSolverPrecision(env.solver().lra().getPrecision(), env.solver().lra().getRelativeTerminationCriterion()); + // Build the equation system matrix and vector. storm::solver::GeneralLinearEquationSolverFactory linearEquationSolverFactory; - bool isEquationSystemFormat = linearEquationSolverFactory.getEquationProblemFormat(env) == storm::solver::LinearEquationSolverProblemFormat::EquationSystem; + bool isEquationSystemFormat = linearEquationSolverFactory.getEquationProblemFormat(subEnv) == storm::solver::LinearEquationSolverProblemFormat::EquationSystem; storm::storage::SparseMatrixBuilder builder(bscc.size(), bscc.size()); - std::vector eqsysVector; - eqsysVector.reserve(bscc.size()); + std::vector eqSysVector; + eqSysVector.reserve(bscc.size()); // The first row asserts that the weighted bias variables and the reward at s_0 sum up to the gain uint64_t row = 0; ValueType entryValue; @@ -724,11 +741,15 @@ namespace storm { builder.addNextValue(row, 0, -storm::utility::one()); } // Compute weighted sum over successor state. As this is a BSCC, each successor state will again be in the BSCC. - bool needDiagonalEntry = isEquationSystemFormat && (row > 0); - ValueType diagonalValue; - if (needDiagonalEntry) { - diagonalValue = exitRateVector ? (*exitRateVector)[globalState] : storm::utility::one(); + auto diagonalValue = storm::utility::zero(); + if (row > 0) { + if (isEquationSystemFormat) { + diagonalValue = exitRateVector ? (*exitRateVector)[globalState] : storm::utility::one(); + } else { + diagonalValue = storm::utility::one() - (exitRateVector ? (*exitRateVector)[globalState] : storm::utility::one()); + } } + bool needDiagonalEntry = !storm::utility::isZero(diagonalValue); for (auto const& entry : rateMatrix.getRow(globalState)) { uint64_t col = toLocalIndexMap[entry.getColumn()]; if (col == 0) { @@ -752,23 +773,21 @@ namespace storm { if (needDiagonalEntry) { builder.addNextValue(row, row, diagonalValue); } - eqsysVector.push_back(valueGetter(globalState)); + eqSysVector.push_back(valueGetter(globalState)); ++row; } // Create a linear equation solver - auto subEnv = env; - subEnv.solver().setLinearEquationSolverPrecision(env.solver().lra().getPrecision(), env.solver().lra().getRelativeTerminationCriterion()); auto solver = linearEquationSolverFactory.create(subEnv, builder.build()); // Check solver requirements. auto requirements = solver->getRequirements(subEnv); STORM_LOG_THROW(!requirements.hasEnabledCriticalRequirement(), storm::exceptions::UncheckedRequirementException, "Solver requirements " + requirements.getEnabledRequirementsAsString() + " not checked."); - // Todo: Find bounds on the bias variable. Just inserting the maximal value from the vector probably does not work. + // Todo: Find bounds on the bias variables. Just inserting the maximal value from the vector probably does not work. std::vector eqSysSol(bscc.size(), storm::utility::zero()); // Take the mean of the rewards as an initial guess for the gain - eqSysSol.front() = std::accumulate(eqsysVector.begin(), eqsysVector.end(), storm::utility::zero()) / storm::utility::convertNumber(bscc.size()); - solver->solveEquations(subEnv, eqSysSol, eqsysVector); + //eqSysSol.front() = std::accumulate(eqSysVector.begin(), eqSysVector.end(), storm::utility::zero()) / storm::utility::convertNumber(bscc.size()); + solver->solveEquations(subEnv, eqSysSol, eqSysVector); ValueType gain = eqSysSol.front(); // insert bias value for state 0 @@ -787,9 +806,17 @@ namespace storm { // This method assumes that this BSCC consist of more than one state if (bscc.size() == 1) { - return { storm::utility::one(), {storm::utility::one()} }; + return { valueGetter(*bscc.begin()), {storm::utility::one()} }; } + // Prepare an environment for the underlying linear equation solver + auto subEnv = env; + if (subEnv.solver().getLinearEquationSolverType() == storm::solver::EquationSolverType::Topological) { + // Topological solver does not make any sense since the BSCC is connected. + subEnv.solver().setLinearEquationSolverType(subEnv.solver().topological().getUnderlyingEquationSolverType()); + } + subEnv.solver().setLinearEquationSolverPrecision(env.solver().lra().getPrecision(), env.solver().lra().getRelativeTerminationCriterion()); + // Get a mapping from global state indices to local ones as well as a bitvector containing states within the BSCC. std::unordered_map toLocalIndexMap; storm::storage::BitVector bsccStates(rateMatrix.getRowCount(), false); @@ -823,12 +850,11 @@ namespace storm { // Check whether we need the fixpoint characterization storm::solver::GeneralLinearEquationSolverFactory linearEquationSolverFactory; - bool isFixpointFormat = linearEquationSolverFactory.getEquationProblemFormat(env) == storm::solver::LinearEquationSolverProblemFormat::FixedPointSystem; - + bool isFixpointFormat = linearEquationSolverFactory.getEquationProblemFormat(subEnv) == storm::solver::LinearEquationSolverProblemFormat::FixedPointSystem; if (isFixpointFormat) { // Add a 1 on the diagonal for (row = 0; row < auxMatrix.getRowCount(); ++row) { - for (auto& entry : auxMatrix) { + for (auto& entry : auxMatrix.getRow(row)) { if (entry.getColumn() == row) { entry.setValue(storm::utility::one() + entry.getValue()); } @@ -869,32 +895,16 @@ namespace storm { bsccEquationSystemRightSide.back() = storm::utility::one(); // Create a linear equation solver - auto subEnv = env; - subEnv.solver().setLinearEquationSolverPrecision(env.solver().lra().getPrecision(), env.solver().lra().getRelativeTerminationCriterion()); auto solver = linearEquationSolverFactory.create(subEnv, builder.build()); solver->setBounds(storm::utility::zero(), storm::utility::one()); // Check solver requirements. auto requirements = solver->getRequirements(subEnv); + requirements.clearLowerBounds(); + requirements.clearUpperBounds(); STORM_LOG_THROW(!requirements.hasEnabledCriticalRequirement(), storm::exceptions::UncheckedRequirementException, "Solver requirements " + requirements.getEnabledRequirementsAsString() + " not checked."); std::vector lraDistr(bscc.size(), storm::utility::one() / storm::utility::convertNumber(bscc.size())); solver->solveEquations(subEnv, lraDistr, bsccEquationSystemRightSide); - // If exit rates were given, we need to 'fix' the results to also account for the timing behaviour. - if (false && exitRateVector != nullptr) { - ValueType totalValue = storm::utility::zero(); - auto solIt = lraDistr.begin(); - for (auto const& globalState : bscc) { - totalValue += (*solIt) * (storm::utility::one() / (*exitRateVector)[globalState]); - ++solIt; - } - assert(solIt == lraDistr.end()); - solIt = lraDistr.begin(); - for (auto const& globalState : bscc) { - *solIt = ((*solIt) * (storm::utility::one() / (*exitRateVector)[globalState])) / totalValue; - ++solIt; - } - assert(solIt == lraDistr.end()); - } // Calculate final LRA Value ValueType result = storm::utility::zero(); diff --git a/src/storm/settings/modules/LongRunAverageSolverSettings.cpp b/src/storm/settings/modules/LongRunAverageSolverSettings.cpp index 9826bdca8..8337bb779 100644 --- a/src/storm/settings/modules/LongRunAverageSolverSettings.cpp +++ b/src/storm/settings/modules/LongRunAverageSolverSettings.cpp @@ -41,7 +41,7 @@ namespace storm { } storm::solver::LraMethod LongRunAverageSolverSettings::getDetLraMethod() const { - std::string lraMethodString = this->getOption(nondetLraMethodOptionName).getArgumentByName("name").getValueAsString(); + std::string lraMethodString = this->getOption(detLraMethodOptionName).getArgumentByName("name").getValueAsString(); if (lraMethodString == "gain-bias-equations" || lraMethodString == "gb") { return storm::solver::LraMethod::GainBiasEquations; } diff --git a/src/storm/solver/SolverSelectionOptions.cpp b/src/storm/solver/SolverSelectionOptions.cpp index 26001226f..f635b5a38 100644 --- a/src/storm/solver/SolverSelectionOptions.cpp +++ b/src/storm/solver/SolverSelectionOptions.cpp @@ -49,9 +49,13 @@ namespace storm { std::string toString(LraMethod m) { switch(m) { case LraMethod::LinearProgramming: - return "linearprogramming"; + return "linear-programming"; case LraMethod::ValueIteration: - return "valueiteration"; + return "value-iteration"; + case LraMethod::LraDistributionEquations: + return "lra-distribution-equations"; + case LraMethod::GainBiasEquations: + return "gain-bias-equations"; } return "invalid"; } diff --git a/src/storm/storage/SparseMatrix.cpp b/src/storm/storage/SparseMatrix.cpp index 857a8d478..d37177b6c 100644 --- a/src/storm/storage/SparseMatrix.cpp +++ b/src/storm/storage/SparseMatrix.cpp @@ -2301,6 +2301,7 @@ namespace storm { template class SparseMatrixBuilder; template class SparseMatrix; template std::ostream& operator<<(std::ostream& out, SparseMatrix const& matrix); + template double SparseMatrix::getPointwiseProductRowSum(storm::storage::SparseMatrix const& otherMatrix, typename SparseMatrix::index_type const& row) const; template std::vector SparseMatrix::getPointwiseProductRowSumVector(storm::storage::SparseMatrix const& otherMatrix) const; template bool SparseMatrix::isSubmatrixOf(SparseMatrix const& matrix) const; @@ -2313,6 +2314,7 @@ namespace storm { template class SparseMatrixBuilder; template class SparseMatrix; template std::ostream& operator<<(std::ostream& out, SparseMatrix const& matrix); + template float SparseMatrix::getPointwiseProductRowSum(storm::storage::SparseMatrix const& otherMatrix, typename SparseMatrix::index_type const& row) const; template std::vector SparseMatrix::getPointwiseProductRowSumVector(storm::storage::SparseMatrix const& otherMatrix) const; template bool SparseMatrix::isSubmatrixOf(SparseMatrix const& matrix) const; @@ -2341,6 +2343,7 @@ namespace storm { template class SparseMatrixBuilder; template class SparseMatrix; template std::ostream& operator<<(std::ostream& out, SparseMatrix const& matrix); + template storm::ClnRationalNumber SparseMatrix::getPointwiseProductRowSum(storm::storage::SparseMatrix const& otherMatrix, typename SparseMatrix::index_type const& row) const; template std::vector SparseMatrix::getPointwiseProductRowSumVector(storm::storage::SparseMatrix const& otherMatrix) const; template bool SparseMatrix::isSubmatrixOf(SparseMatrix const& matrix) const; #endif @@ -2351,6 +2354,7 @@ namespace storm { template class SparseMatrixBuilder; template class SparseMatrix; template std::ostream& operator<<(std::ostream& out, SparseMatrix const& matrix); + template storm::GmpRationalNumber SparseMatrix::getPointwiseProductRowSum(storm::storage::SparseMatrix const& otherMatrix, typename SparseMatrix::index_type const& row) const; template std::vector SparseMatrix::getPointwiseProductRowSumVector(storm::storage::SparseMatrix const& otherMatrix) const; template bool SparseMatrix::isSubmatrixOf(SparseMatrix const& matrix) const; #endif @@ -2361,6 +2365,10 @@ namespace storm { template class SparseMatrixBuilder; template class SparseMatrix; template std::ostream& operator<<(std::ostream& out, SparseMatrix const& matrix); + template storm::RationalFunction SparseMatrix::getPointwiseProductRowSum(storm::storage::SparseMatrix const& otherMatrix, typename SparseMatrix::index_type const& row) const; + template storm::RationalFunction SparseMatrix::getPointwiseProductRowSum(storm::storage::SparseMatrix const& otherMatrix, typename SparseMatrix::index_type const& row) const; + template storm::RationalFunction SparseMatrix::getPointwiseProductRowSum(storm::storage::SparseMatrix const& otherMatrix, typename SparseMatrix::index_type const& row) const; + template storm::RationalFunction SparseMatrix::getPointwiseProductRowSum(storm::storage::SparseMatrix const& otherMatrix, typename SparseMatrix::index_type const& row) const; template std::vector SparseMatrix::getPointwiseProductRowSumVector(storm::storage::SparseMatrix const& otherMatrix) const; template std::vector SparseMatrix::getPointwiseProductRowSumVector(storm::storage::SparseMatrix const& otherMatrix) const; template std::vector SparseMatrix::getPointwiseProductRowSumVector(storm::storage::SparseMatrix const& otherMatrix) const; diff --git a/src/test/storm/modelchecker/csl/LraCtmcCslModelCheckerTest.cpp b/src/test/storm/modelchecker/csl/LraCtmcCslModelCheckerTest.cpp index 629684560..b8d3b1fc7 100755 --- a/src/test/storm/modelchecker/csl/LraCtmcCslModelCheckerTest.cpp +++ b/src/test/storm/modelchecker/csl/LraCtmcCslModelCheckerTest.cpp @@ -144,9 +144,8 @@ namespace { env.solver().setLinearEquationSolverType(storm::solver::EquationSolverType::Native); env.solver().native().setMethod(storm::solver::NativeLinearEquationSolverMethod::SOR); env.solver().native().setSorOmega(storm::utility::convertNumber(0.7)); // LRA computation fails for 0.9 - env.solver().native().setPrecision(storm::utility::convertNumber(1e-9)); - env.solver().native().setRelativeTerminationCriterion(false); env.solver().lra().setDetLraMethod(storm::solver::LraMethod::GainBiasEquations); + env.solver().lra().setPrecision(storm::utility::convertNumber(1e-9)); return env; } }; From 62dc50035c543efd9379ed04fa941ea03eb631ac Mon Sep 17 00:00:00 2001 From: Tim Quatmann Date: Tue, 3 Dec 2019 18:05:20 +0100 Subject: [PATCH 10/24] Removed a test case that is not relevant anymore. --- .../csl/CtmcCslModelCheckerTest.cpp | 18 ------------------ 1 file changed, 18 deletions(-) diff --git a/src/test/storm/modelchecker/csl/CtmcCslModelCheckerTest.cpp b/src/test/storm/modelchecker/csl/CtmcCslModelCheckerTest.cpp index 1adc011b5..83414e72f 100755 --- a/src/test/storm/modelchecker/csl/CtmcCslModelCheckerTest.cpp +++ b/src/test/storm/modelchecker/csl/CtmcCslModelCheckerTest.cpp @@ -368,24 +368,6 @@ namespace { } - TYPED_TEST(CtmcCslModelCheckerTest, Polling) { - std::string formulasString = "P=?[ F<=10 \"target\"]"; - - auto modelFormulas = this->buildModelFormulas(STORM_TEST_RESOURCES_DIR "/ctmc/polling2.sm", formulasString); - auto model = std::move(modelFormulas.first); - auto tasks = this->getTasks(modelFormulas.second); - EXPECT_EQ(12ul, model->getNumberOfStates()); - EXPECT_EQ(22ul, model->getNumberOfTransitions()); - ASSERT_EQ(model->getType(), storm::models::ModelType::Ctmc); - auto checker = this->createModelChecker(model); - std::unique_ptr result; - - result = checker->check(this->env(), tasks[0]); - EXPECT_NEAR(this->parseNumber("1"), this->getQuantitativeResultAtInitialState(model, result), this->precision()); - - - } - TYPED_TEST(CtmcCslModelCheckerTest, Tandem) { std::string formulasString = "P=? [ F<=10 \"network_full\" ]"; formulasString += "; P=? [ F<=10 \"first_queue_full\" ]"; From 3912d59a3b295927f42a93264c5dc6a9f4d6eb6d Mon Sep 17 00:00:00 2001 From: Tim Quatmann Date: Tue, 3 Dec 2019 18:05:40 +0100 Subject: [PATCH 11/24] Added kanban model for LRA test --- .../examples/testfiles/ctmc/kanban.prism | 111 ++++++++++++++++++ .../csl/LraCtmcCslModelCheckerTest.cpp | 17 +++ 2 files changed, 128 insertions(+) create mode 100644 resources/examples/testfiles/ctmc/kanban.prism diff --git a/resources/examples/testfiles/ctmc/kanban.prism b/resources/examples/testfiles/ctmc/kanban.prism new file mode 100644 index 000000000..8f5191efb --- /dev/null +++ b/resources/examples/testfiles/ctmc/kanban.prism @@ -0,0 +1,111 @@ +// kanban manufacturing example [CT96] +// dxp/gxn 3/2/00 + +ctmc + +// number of tokens +const int t = 1; + +// rates +const double in1 = 1.0; +const double out4 = 0.9; +const double synch123 = 0.4; +const double synch234 = 0.5; +const double back = 0.3; +const double redo1 = 0.36; +const double redo2 = 0.42; +const double redo3 = 0.39; +const double redo4 = 0.33; +const double ok1 = 0.84; +const double ok2 = 0.98; +const double ok3 = 0.91; +const double ok4 = 0.77; + +module k1 + + w1 : [0..t]; + x1 : [0..t]; + y1 : [0..t]; + z1 : [0..t]; + + [in] (w1 in1 : (w1'=w1+1) & (x1'=x1+1); + [] (x1>0) & (y1 redo1 : (x1'=x1-1) & (y1'=y1+1); + [] (x1>0) & (z1 ok1 : (x1'=x1-1) & (z1'=z1+1); + [] (y1>0) & (x1 back : (y1'=y1-1) & (x1'=x1+1); + [s1] (z1>0) & (w1>0) -> synch123 : (z1'=z1-1) & (w1'=w1-1); + +endmodule + +module k2 + + w2 : [0..t]; + x2 : [0..t]; + y2 : [0..t]; + z2 : [0..t]; + + [s1] (w2 1 : (w2'=w2+1) & (x2'=x2+1); + [] (x2>0) & (y2 redo2 : (x2'=x2-1) & (y2'=y2+1); + [] (x2>0) & (z2 ok2 : (x2'=x2-1) & (z2'=z2+1); + [] (y2>0) & (x2 back : (y2'=y2-1) & (x2'=x2+1); + [s2] (z2>0) & (w2>0) -> 1 : (z2'=z2-1) & (w2'=w2-1); + +endmodule + +module k3 + + w3 : [0..t]; + x3 : [0..t]; + y3 : [0..t]; + z3 : [0..t]; + + [s1] (w3 1 : (w3'=w3+1) & (x3'=x3+1); + [] (x3>0) & (y3 redo3 : (x3'=x3-1) & (y3'=y3+1); + [] (x3>0) & (z3 ok3 : (x3'=x3-1) & (z3'=z3+1); + [] (y3>0) & (x3 back : (y3'=y3-1) & (x3'=x3+1); + [s2] (z3>0) & (w3>0) -> 1 : (z3'=z3-1) & (w3'=w3-1); + + +endmodule + +module k4 + + w4 : [0..t]; + x4 : [0..t]; + y4 : [0..t]; + z4 : [0..t]; + + [s2] (w4 synch234 : (w4'=w4+1) & (x4'=x4+1); + [] (x4>0) & (y4 redo4 : (x4'=x4-1) & (y4'=y4+1); + [] (x4>0) & (z4 ok4 : (x4'=x4-1) & (z4'=z4+1); + [] (y4>0) & (x4 back : (y4'=y4-1) & (x4'=x4+1); + [] (z4>0) & (w4>0) -> out4 : (z4'=z4-1) & (w4'=w4-1); + +endmodule + +// reward structures + +// tokens in cell1 +rewards "tokens_cell1" + true : x1+y1+z1; +endrewards + +// tokens in cell2 +rewards "tokens_cell2" + true : x2+y2+z2; +endrewards + +// tokens in cell3 +rewards "tokens_cell3" + true : x3+y3+z3; +endrewards + +// tokens in cell4 +rewards "tokens_cell4" + true : x4+y4+z4; +endrewards + +// throughput of the system +rewards "throughput" + [in] true : 1; +endrewards + diff --git a/src/test/storm/modelchecker/csl/LraCtmcCslModelCheckerTest.cpp b/src/test/storm/modelchecker/csl/LraCtmcCslModelCheckerTest.cpp index b8d3b1fc7..f142a764c 100755 --- a/src/test/storm/modelchecker/csl/LraCtmcCslModelCheckerTest.cpp +++ b/src/test/storm/modelchecker/csl/LraCtmcCslModelCheckerTest.cpp @@ -396,4 +396,21 @@ namespace { EXPECT_NEAR(this->parseNumber("11/15"), this->getQuantitativeResultAtInitialState(model, result), this->precision()); } + + TYPED_TEST(LraCtmcCslModelCheckerTest, kanban) { + std::string formulasString = "R{\"throughput\"}=? [ LRA ]"; + + auto modelFormulas = this->buildModelFormulas(STORM_TEST_RESOURCES_DIR "/ctmc/kanban.prism", formulasString); + auto model = std::move(modelFormulas.first); + auto tasks = this->getTasks(modelFormulas.second); + EXPECT_EQ(160ul, model->getNumberOfStates()); + EXPECT_EQ(616ul, model->getNumberOfTransitions()); + ASSERT_EQ(model->getType(), storm::models::ModelType::Ctmc); + auto checker = this->createModelChecker(model); + std::unique_ptr result; + + result = checker->check(this->env(), tasks[0]); + EXPECT_NEAR(this->parseNumber("113237255213395163953677015242972426399989689654967642609491830216061334090202313396984106738516704120069048184391587092670711590526535239899047608853509681914074220789038015289373871985431257486278/1223067474012215838745994023624181143448435715858923491568712969277184634645504346456117443333632535902774869182827010789201972713368729205674432492059242349591780604188152950845769793378621446766887"), this->getQuantitativeResultAtInitialState(model, result), this->precision()); + } + } From 267291ea900b0911d05a31d90190ffa17bbcdc7e Mon Sep 17 00:00:00 2001 From: Tim Quatmann Date: Tue, 3 Dec 2019 18:05:58 +0100 Subject: [PATCH 12/24] Updated Changelog --- CHANGELOG.md | 3 ++- 1 file changed, 2 insertions(+), 1 deletion(-) diff --git a/CHANGELOG.md b/CHANGELOG.md index 150df6c4b..5647623f5 100644 --- a/CHANGELOG.md +++ b/CHANGELOG.md @@ -8,7 +8,8 @@ Version 1.4.x ------------- ### Version 1.4.1 (under development) - +- Implemented Long run average computation for DTMCs/CTMCs via value iteration and via gain/bias equations. +- Added several LRA related settings in a new settings module. Note that `--minmax:lramethod` has been replaced by `--lra:nondetmethod`. ### Version 1.4.0 (2019/11) - Added support for multi-dimensional quantile queries. From 1de76b36c842460b4e2ca75cdf35f3481eab6539 Mon Sep 17 00:00:00 2001 From: Matthias Volk Date: Tue, 3 Dec 2019 17:27:59 +0100 Subject: [PATCH 13/24] Updated steps for new release --- doc/checklist_new_release.md | 16 ++++++++++++---- 1 file changed, 12 insertions(+), 4 deletions(-) diff --git a/doc/checklist_new_release.md b/doc/checklist_new_release.md index e14f02667..553ea68d9 100644 --- a/doc/checklist_new_release.md +++ b/doc/checklist_new_release.md @@ -1,11 +1,12 @@ The following steps should be performed before releasing a new storm version. Note that in most cases a simultaneous release of [carl](https://github.com/smtrat/carl), [storm](https://github.com/moves-rwth/storm), [pycarl](https://github.com/moves-rwth/pycarl/) and [stormpy](https://github.com/moves-rwth/stormpy/) is preferred. -1. Update `CHANGELOG.md`. - To get all the commits from an author since the last tag execute: +1. Update `CHANGELOG.md`: + * To get all the commits from an author since the last tag execute: ```console git log last_tag..HEAD --author "author_name" ``` + * Set release month 2. Update used carl version: * Update `GIT_TAG` in `resources/3rdparty/carl/CMakeLists.txt` @@ -32,6 +33,13 @@ Note that in most cases a simultaneous release of [carl](https://github.com/smtr 6. [Add new release](https://github.com/moves-rwth/storm/releases/new) in GitHub. -7. Update [Homebrew formula](https://github.com/moves-rwth/homebrew-storm). +7. Update `stable` branch: + ```console + git checkout stable + git merge master + git push origin stable + ``` + +8. Update [Homebrew formula](https://github.com/moves-rwth/homebrew-storm). -8. Announce new storm version on [website](http://www.stormchecker.org/news.html). +9. Announce new storm version on [website](http://www.stormchecker.org/news.html). From 6fb9f7e74387c136e4a637b1553fcca4ea51902d Mon Sep 17 00:00:00 2001 From: Matthias Volk Date: Tue, 3 Dec 2019 17:28:33 +0100 Subject: [PATCH 14/24] Warn if property could not be checked on DFT --- .../modelchecker/dft/DFTModelChecker.cpp | 16 ++++++++++------ 1 file changed, 10 insertions(+), 6 deletions(-) diff --git a/src/storm-dft/modelchecker/dft/DFTModelChecker.cpp b/src/storm-dft/modelchecker/dft/DFTModelChecker.cpp index 5bc4b0f4c..83887c503 100644 --- a/src/storm-dft/modelchecker/dft/DFTModelChecker.cpp +++ b/src/storm-dft/modelchecker/dft/DFTModelChecker.cpp @@ -466,13 +466,17 @@ namespace storm { singleModelCheckingTimer.start(); //STORM_PRINT_AND_LOG("Model checking property " << *property << " ..." << std::endl); std::unique_ptr result( - storm::api::verifyWithSparseEngine(model, storm::api::createTask(property, - true))); - STORM_LOG_ASSERT(result, "Result does not exist."); - result->filter(storm::modelchecker::ExplicitQualitativeCheckResult(model->getInitialStates())); - ValueType resultValue = result->asExplicitQuantitativeCheckResult().getValueMap().begin()->second; + storm::api::verifyWithSparseEngine(model, storm::api::createTask(property,true))); + + if (result) { + result->filter(storm::modelchecker::ExplicitQualitativeCheckResult(model->getInitialStates())); + ValueType resultValue = result->asExplicitQuantitativeCheckResult().getValueMap().begin()->second; + results.push_back(resultValue); + } else { + STORM_LOG_WARN("The property '" << *property << "' could not be checked with the current settings."); + results.push_back(-storm::utility::one()); + } //STORM_PRINT_AND_LOG("Result (initial states): " << resultValue << std::endl); - results.push_back(resultValue); singleModelCheckingTimer.stop(); //STORM_PRINT_AND_LOG("Time for model checking: " << singleModelCheckingTimer << "." << std::endl); } From bc85e6742e11ea33dab2ae2dde35a0cb81616049 Mon Sep 17 00:00:00 2001 From: Matthias Volk Date: Tue, 3 Dec 2019 18:25:02 +0100 Subject: [PATCH 15/24] Fixed parsing of RationalFunctions if no parameters are given --- src/storm-parsers/parser/ValueParser.h | 2 ++ 1 file changed, 2 insertions(+) diff --git a/src/storm-parsers/parser/ValueParser.h b/src/storm-parsers/parser/ValueParser.h index 768be44ad..777db0679 100644 --- a/src/storm-parsers/parser/ValueParser.h +++ b/src/storm-parsers/parser/ValueParser.h @@ -19,6 +19,8 @@ namespace storm { * Constructor. */ ValueParser() : manager(new storm::expressions::ExpressionManager()), parser(*manager), evaluator(*manager) { + // Set empty mapping to enable expression creation even without parameters + parser.setIdentifierMapping(identifierMapping); } /*! From 419013025bb44467a741260da7e6726e868ac348 Mon Sep 17 00:00:00 2001 From: TimQu Date: Thu, 5 Dec 2019 12:28:58 +0100 Subject: [PATCH 16/24] Fixed reduction to state-based rewards for CTMCs in Dd engine. When action rewards are reduced to state rewards, they have to be multiplied with the exit rate. --- src/storm/models/symbolic/Ctmc.cpp | 11 +++++++++++ src/storm/models/symbolic/Ctmc.h | 2 ++ src/storm/models/symbolic/DeterministicModel.cpp | 7 ------- src/storm/models/symbolic/DeterministicModel.h | 1 - src/storm/models/symbolic/Dtmc.cpp | 7 +++++++ src/storm/models/symbolic/Dtmc.h | 2 ++ 6 files changed, 22 insertions(+), 8 deletions(-) diff --git a/src/storm/models/symbolic/Ctmc.cpp b/src/storm/models/symbolic/Ctmc.cpp index 4b05372be..d54b3156d 100644 --- a/src/storm/models/symbolic/Ctmc.cpp +++ b/src/storm/models/symbolic/Ctmc.cpp @@ -84,6 +84,17 @@ namespace storm { return exitRates.get(); } + template + void Ctmc::reduceToStateBasedRewards() { + for (auto& rewardModel : this->getRewardModels()) { + if (rewardModel.second.hasStateActionRewards()) { + rewardModel.second.getStateActionRewardVector() *= getExitRateVector(); + } + rewardModel.second.reduceToStateBasedRewards(this->getTransitionMatrix(), this->getRowVariables(), this->getColumnVariables(), true); + } + } + + template template std::shared_ptr> Ctmc::toValueType() const { diff --git a/src/storm/models/symbolic/Ctmc.h b/src/storm/models/symbolic/Ctmc.h index 0035e3ffb..c73340946 100644 --- a/src/storm/models/symbolic/Ctmc.h +++ b/src/storm/models/symbolic/Ctmc.h @@ -141,6 +141,8 @@ namespace storm { */ storm::dd::Add const& getExitRateVector() const; + virtual void reduceToStateBasedRewards() override; + template std::shared_ptr> toValueType() const; diff --git a/src/storm/models/symbolic/DeterministicModel.cpp b/src/storm/models/symbolic/DeterministicModel.cpp index 5eba1bbc9..3555dcaff 100644 --- a/src/storm/models/symbolic/DeterministicModel.cpp +++ b/src/storm/models/symbolic/DeterministicModel.cpp @@ -45,13 +45,6 @@ namespace storm { // Intentionally left empty. } - template - void DeterministicModel::reduceToStateBasedRewards() { - for (auto& rewardModel : this->getRewardModels()) { - rewardModel.second.reduceToStateBasedRewards(this->getTransitionMatrix(), this->getRowVariables(), this->getColumnVariables(), true); - } - } - // Explicitly instantiate the template class. template class DeterministicModel; template class DeterministicModel; diff --git a/src/storm/models/symbolic/DeterministicModel.h b/src/storm/models/symbolic/DeterministicModel.h index a8fde85b4..4612354b7 100644 --- a/src/storm/models/symbolic/DeterministicModel.h +++ b/src/storm/models/symbolic/DeterministicModel.h @@ -81,7 +81,6 @@ namespace storm { std::map> labelToBddMap = std::map>(), std::unordered_map const& rewardModels = std::unordered_map()); - virtual void reduceToStateBasedRewards() override; }; } // namespace symbolic diff --git a/src/storm/models/symbolic/Dtmc.cpp b/src/storm/models/symbolic/Dtmc.cpp index 2b4adb149..a4b422360 100644 --- a/src/storm/models/symbolic/Dtmc.cpp +++ b/src/storm/models/symbolic/Dtmc.cpp @@ -43,6 +43,13 @@ namespace storm { // Intentionally left empty. } + template + void Dtmc::reduceToStateBasedRewards() { + for (auto& rewardModel : this->getRewardModels()) { + rewardModel.second.reduceToStateBasedRewards(this->getTransitionMatrix(), this->getRowVariables(), this->getColumnVariables(), true); + } + } + template template std::shared_ptr> Dtmc::toValueType() const { diff --git a/src/storm/models/symbolic/Dtmc.h b/src/storm/models/symbolic/Dtmc.h index 76a3e81ca..be98cf34d 100644 --- a/src/storm/models/symbolic/Dtmc.h +++ b/src/storm/models/symbolic/Dtmc.h @@ -77,6 +77,8 @@ namespace storm { std::map> labelToBddMap = std::map>(), std::unordered_map const& rewardModels = std::unordered_map()); + virtual void reduceToStateBasedRewards() override; + template std::shared_ptr> toValueType() const; From be3cffe8ba9dd579d26498c0a6356f0443d0a9ff Mon Sep 17 00:00:00 2001 From: Jip Spel Date: Fri, 6 Dec 2019 11:27:28 +0100 Subject: [PATCH 17/24] Write output monotonicity checking to user-specified file --- src/storm-pars-cli/storm-pars.cpp | 9 ++++++++- src/storm-pars/analysis/MonotonicityChecker.cpp | 8 +++----- src/storm-pars/analysis/MonotonicityChecker.h | 8 ++------ src/storm/settings/modules/IOSettings.cpp | 10 ++++++++++ src/storm/settings/modules/IOSettings.h | 11 +++++++++++ .../storm-pars/analysis/MonotonicityCheckerTest.cpp | 4 ++-- 6 files changed, 36 insertions(+), 14 deletions(-) diff --git a/src/storm-pars-cli/storm-pars.cpp b/src/storm-pars-cli/storm-pars.cpp index 75dba4492..de8542e25 100644 --- a/src/storm-pars-cli/storm-pars.cpp +++ b/src/storm-pars-cli/storm-pars.cpp @@ -744,7 +744,14 @@ namespace storm { STORM_LOG_THROW(regions.size() <= 1, storm::exceptions::InvalidArgumentException, "Monotonicity analysis only allowed on single region"); storm::analysis::MonotonicityChecker monotonicityChecker = storm::analysis::MonotonicityChecker(model, formulas, regions, monSettings.isValidateAssumptionsSet(), monSettings.getNumberOfSamples(), monSettings.getMonotonicityAnalysisPrecision()); - monotonicityChecker.checkMonotonicity(); + if (ioSettings.isExportMonotonicitySet()) { + std::ofstream outfile; + utility::openFile(ioSettings.getExportMonotonicityFilename(), outfile); + monotonicityChecker.checkMonotonicity(outfile); + utility::closeFile(outfile); + } else { + monotonicityChecker.checkMonotonicity(std::cout); + } monotonicityWatch.stop(); STORM_PRINT(std::endl << "Total time for monotonicity checking: " << monotonicityWatch << "." << std::endl << std::endl); diff --git a/src/storm-pars/analysis/MonotonicityChecker.cpp b/src/storm-pars/analysis/MonotonicityChecker.cpp index c53609220..f8ba84c7e 100644 --- a/src/storm-pars/analysis/MonotonicityChecker.cpp +++ b/src/storm-pars/analysis/MonotonicityChecker.cpp @@ -67,19 +67,18 @@ namespace storm { } template - std::map::type, std::pair>> MonotonicityChecker::checkMonotonicity() { + std::map::type, std::pair>> MonotonicityChecker::checkMonotonicity(std::ostream& outfile) { auto map = createOrder(); std::shared_ptr> sparseModel = model->as>(); auto matrix = sparseModel->getTransitionMatrix(); - return checkMonotonicity(map, matrix); + return checkMonotonicity(outfile, map, matrix); } template - std::map::type, std::pair>> MonotonicityChecker::checkMonotonicity(std::map>> map, storm::storage::SparseMatrix matrix) { + std::map::type, std::pair>> MonotonicityChecker::checkMonotonicity(std::ostream& outfile, std::map>> map, storm::storage::SparseMatrix matrix) { storm::utility::Stopwatch monotonicityCheckWatch(true); std::map::type, std::pair>> result; - outfile.open(filename, std::ios_base::app); if (map.size() == 0) { // Nothing is known @@ -166,7 +165,6 @@ namespace storm { outfile << ", "; monotonicityCheckWatch.stop(); - outfile.close(); return result; } diff --git a/src/storm-pars/analysis/MonotonicityChecker.h b/src/storm-pars/analysis/MonotonicityChecker.h index 5d32ecf70..31fc0e193 100644 --- a/src/storm-pars/analysis/MonotonicityChecker.h +++ b/src/storm-pars/analysis/MonotonicityChecker.h @@ -41,7 +41,7 @@ namespace storm { /*! * Checks for model and formula as provided in constructor for monotonicity */ - std::map::type, std::pair>> checkMonotonicity(); + std::map::type, std::pair>> checkMonotonicity(std::ostream& outfile); /*! * Checks if monotonicity can be found in this order. Unordered states are not checked @@ -113,7 +113,7 @@ namespace storm { } private: - std::map::type, std::pair>> checkMonotonicity(std::map>> map, storm::storage::SparseMatrix matrix); + std::map::type, std::pair>> checkMonotonicity(std::ostream& outfile, std::map>> map, storm::storage::SparseMatrix matrix); std::map::type, std::pair> analyseMonotonicity(uint_fast64_t i, Order* order, storm::storage::SparseMatrix matrix) ; @@ -141,10 +141,6 @@ namespace storm { OrderExtender *extender; - std::ofstream outfile; - - std::string filename = "monotonicity.txt"; - double precision; storm::storage::ParameterRegion region; diff --git a/src/storm/settings/modules/IOSettings.cpp b/src/storm/settings/modules/IOSettings.cpp index 86866e7e6..7f78e93cb 100644 --- a/src/storm/settings/modules/IOSettings.cpp +++ b/src/storm/settings/modules/IOSettings.cpp @@ -25,6 +25,7 @@ namespace storm { const std::string IOSettings::exportCdfOptionName = "exportcdf"; const std::string IOSettings::exportCdfOptionShortName = "cdf"; const std::string IOSettings::exportSchedulerOptionName = "exportscheduler"; + const std::string IOSettings::exportMonotonicityName = "exportmonotonicity"; const std::string IOSettings::explicitOptionName = "explicit"; const std::string IOSettings::explicitOptionShortName = "exp"; const std::string IOSettings::explicitDrnOptionName = "explicit-drn"; @@ -59,6 +60,7 @@ namespace storm { .addArgument(storm::settings::ArgumentBuilder::createStringArgument("filename", "The name of the file to which the model is to be written.").build()).build()); this->addOption(storm::settings::OptionBuilder(moduleName, exportCdfOptionName, false, "Exports the cumulative density function for reward bounded properties into a .csv file.").setIsAdvanced().setShortName(exportCdfOptionShortName).addArgument(storm::settings::ArgumentBuilder::createStringArgument("directory", "A path to an existing directory where the cdf files will be stored.").build()).build()); this->addOption(storm::settings::OptionBuilder(moduleName, exportSchedulerOptionName, false, "Exports the choices of an optimal scheduler to the given file (if supported by engine).").setIsAdvanced().addArgument(storm::settings::ArgumentBuilder::createStringArgument("filename", "The output file.").build()).build()); + this->addOption(storm::settings::OptionBuilder(moduleName, exportMonotonicityName, false, "Exports the result of monotonicity checking to the given file.").setIsAdvanced().addArgument(storm::settings::ArgumentBuilder::createStringArgument("filename", "The output file.").build()).build()); this->addOption(storm::settings::OptionBuilder(moduleName, exportExplicitOptionName, "", "If given, the loaded model will be written to the specified file in the drn format.") .addArgument(storm::settings::ArgumentBuilder::createStringArgument("filename", "the name of the file to which the model is to be writen.").build()).build()); this->addOption(storm::settings::OptionBuilder(moduleName, exportDdOptionName, "", "If given, the loaded model will be written to the specified file in the drdd format.") @@ -162,6 +164,14 @@ namespace storm { std::string IOSettings::getExportSchedulerFilename() const { return this->getOption(exportSchedulerOptionName).getArgumentByName("filename").getValueAsString(); } + + bool IOSettings::isExportMonotonicitySet() const { + return this->getOption(exportMonotonicityName).getHasOptionBeenSet(); + } + + std::string IOSettings::getExportMonotonicityFilename() const { + return this->getOption(exportMonotonicityName).getArgumentByName("filename").getValueAsString(); + } bool IOSettings::isExplicitSet() const { return this->getOption(explicitOptionName).getHasOptionBeenSet(); diff --git a/src/storm/settings/modules/IOSettings.h b/src/storm/settings/modules/IOSettings.h index 91cc3b765..2cbdf65bf 100644 --- a/src/storm/settings/modules/IOSettings.h +++ b/src/storm/settings/modules/IOSettings.h @@ -105,6 +105,16 @@ namespace storm { * Retrieves a filename to which an optimal scheduler will be exported. */ std::string getExportSchedulerFilename() const; + + /*! + * Retrieves whether an optimal scheduler is to be exported + */ + bool isExportMonotonicitySet() const; + + /*! + * Retrieves a filename to which an optimal scheduler will be exported. + */ + std::string getExportMonotonicityFilename() const; /*! * Retrieves whether the explicit option was set. @@ -339,6 +349,7 @@ namespace storm { static const std::string exportCdfOptionName; static const std::string exportCdfOptionShortName; static const std::string exportSchedulerOptionName; + static const std::string exportMonotonicityName; static const std::string explicitOptionName; static const std::string explicitOptionShortName; static const std::string explicitDrnOptionName; diff --git a/src/test/storm-pars/analysis/MonotonicityCheckerTest.cpp b/src/test/storm-pars/analysis/MonotonicityCheckerTest.cpp index 0eda2d56d..14973cfc2 100644 --- a/src/test/storm-pars/analysis/MonotonicityCheckerTest.cpp +++ b/src/test/storm-pars/analysis/MonotonicityCheckerTest.cpp @@ -149,7 +149,7 @@ std::vector> regions = ASSERT_EQ(dtmc->getNumberOfTransitions(), 195ull); storm::analysis::MonotonicityChecker monotonicityChecker = storm::analysis::MonotonicityChecker(dtmc, formulas, regions, true); - auto result = monotonicityChecker.checkMonotonicity(); + auto result = monotonicityChecker.checkMonotonicity(std::cout); EXPECT_EQ(1ul, result.size()); EXPECT_EQ(2ul, result.begin()->second.size()); auto monotone = result.begin()->second.begin(); @@ -197,7 +197,7 @@ TEST(MonotonicityCheckerTest, Brp_with_bisimulation_samples) { ASSERT_EQ(dtmc->getNumberOfTransitions(), 195ull); auto monotonicityChecker = storm::analysis::MonotonicityChecker(dtmc, formulas, regions, true, 50); - auto result = monotonicityChecker.checkMonotonicity(); + auto result = monotonicityChecker.checkMonotonicity(std::cout); EXPECT_EQ(1ul, result.size()); EXPECT_EQ(2ul, result.begin()->second.size()); auto monotone = result.begin()->second.begin(); From b1bb7872fd08df9f9b5feddfe10686dd05fc5167 Mon Sep 17 00:00:00 2001 From: TimQu Date: Tue, 10 Dec 2019 14:51:34 +0100 Subject: [PATCH 18/24] Parsing integers as long long ints (instead of ints) to fix github issue #60. --- src/storm-parsers/parser/ExpressionParser.cpp | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/storm-parsers/parser/ExpressionParser.cpp b/src/storm-parsers/parser/ExpressionParser.cpp index e88068ea7..9c371e4ac 100644 --- a/src/storm-parsers/parser/ExpressionParser.cpp +++ b/src/storm-parsers/parser/ExpressionParser.cpp @@ -81,7 +81,7 @@ namespace storm { literalExpression = qi::lit("true")[qi::_val = phoenix::bind(&ExpressionCreator::createBooleanLiteralExpression, phoenix::ref(*expressionCreator), true, qi::_pass)] | qi::lit("false")[qi::_val = phoenix::bind(&ExpressionCreator::createBooleanLiteralExpression, phoenix::ref(*expressionCreator), false, qi::_pass)] | rationalLiteral_[qi::_val = phoenix::bind(&ExpressionCreator::createRationalLiteralExpression, phoenix::ref(*expressionCreator), qi::_1, qi::_pass)] - | qi::int_[qi::_val = phoenix::bind(&ExpressionCreator::createIntegerLiteralExpression, phoenix::ref(*expressionCreator), qi::_1, qi::_pass)]; + | qi::long_long[qi::_val = phoenix::bind(&ExpressionCreator::createIntegerLiteralExpression, phoenix::ref(*expressionCreator), qi::_1, qi::_pass)]; literalExpression.name("literal expression"); atomicExpression = floorCeilExpression | roundExpression | prefixPowerModuloExpression | minMaxExpression | (qi::lit("(") >> expression >> qi::lit(")")) | identifierExpression | literalExpression; From cf8337583b6384228286caa5315f9d80a0123355 Mon Sep 17 00:00:00 2001 From: Matthias Volk Date: Fri, 13 Dec 2019 14:29:18 +0100 Subject: [PATCH 19/24] Storm version 1.4.1 --- CHANGELOG.md | 6 ++++-- version.cmake | 2 +- 2 files changed, 5 insertions(+), 3 deletions(-) diff --git a/CHANGELOG.md b/CHANGELOG.md index 5647623f5..53253488f 100644 --- a/CHANGELOG.md +++ b/CHANGELOG.md @@ -7,8 +7,8 @@ The releases of major and minor versions contain an overview of changes since th Version 1.4.x ------------- -### Version 1.4.1 (under development) -- Implemented Long run average computation for DTMCs/CTMCs via value iteration and via gain/bias equations. +### Version 1.4.1 (2019/12) +- Implemented long run average (LRA) computation for DTMCs/CTMCs via value iteration and via gain/bias equations. - Added several LRA related settings in a new settings module. Note that `--minmax:lramethod` has been replaced by `--lra:nondetmethod`. ### Version 1.4.0 (2019/11) @@ -45,6 +45,7 @@ Version 1.4.x - Fixed compilation for macOS Mojave and higher. - Several bug fixes. + Version 1.3.x ------------- @@ -80,6 +81,7 @@ Version 1.3.x - `storm-gspn` improved - Sound value iteration + Version 1.2.x ------------- diff --git a/version.cmake b/version.cmake index 6800eea7d..5f396abc3 100644 --- a/version.cmake +++ b/version.cmake @@ -1,4 +1,4 @@ set(STORM_VERSION_MAJOR 1) set(STORM_VERSION_MINOR 4) -set(STORM_VERSION_PATCH 0) +set(STORM_VERSION_PATCH 1) From 32371c44d4b9953c5e73268850482aa3733f0200 Mon Sep 17 00:00:00 2001 From: Matthias Volk Date: Fri, 13 Dec 2019 15:35:21 +0100 Subject: [PATCH 20/24] Updated documentation for new release --- CHANGELOG.md | 3 +++ doc/checklist_new_release.md | 2 ++ 2 files changed, 5 insertions(+) diff --git a/CHANGELOG.md b/CHANGELOG.md index 53253488f..de69b3b1d 100644 --- a/CHANGELOG.md +++ b/CHANGELOG.md @@ -7,6 +7,9 @@ The releases of major and minor versions contain an overview of changes since th Version 1.4.x ------------- +## Version 1.4.2 (under development) +- n.a. + ### Version 1.4.1 (2019/12) - Implemented long run average (LRA) computation for DTMCs/CTMCs via value iteration and via gain/bias equations. - Added several LRA related settings in a new settings module. Note that `--minmax:lramethod` has been replaced by `--lra:nondetmethod`. diff --git a/doc/checklist_new_release.md b/doc/checklist_new_release.md index 553ea68d9..dcb1351a1 100644 --- a/doc/checklist_new_release.md +++ b/doc/checklist_new_release.md @@ -43,3 +43,5 @@ Note that in most cases a simultaneous release of [carl](https://github.com/smtr 8. Update [Homebrew formula](https://github.com/moves-rwth/homebrew-storm). 9. Announce new storm version on [website](http://www.stormchecker.org/news.html). + +10. Create [Docker containers](https://hub.docker.com/r/movesrwth/storm) for new version. From 3f717202cd9bdf0a4e6c7003d52cd99c1a442637 Mon Sep 17 00:00:00 2001 From: Matthias Volk Date: Fri, 13 Dec 2019 16:01:00 +0100 Subject: [PATCH 21/24] Fixed handling of optional arguments. Missing an optional argument now throws an exception. --- src/storm/settings/Argument.cpp | 16 ++++++++++++---- src/storm/settings/SettingsManager.cpp | 2 +- 2 files changed, 13 insertions(+), 5 deletions(-) diff --git a/src/storm/settings/Argument.cpp b/src/storm/settings/Argument.cpp index 555bcdb9d..7e2d64c0b 100644 --- a/src/storm/settings/Argument.cpp +++ b/src/storm/settings/Argument.cpp @@ -171,20 +171,28 @@ namespace storm { out << std::setw(0) << std::left << "<" << this->getName() << ">"; if (!this->validators.empty() || this->hasDefaultValue) { out << " ("; + bool previousEntry = false; + if (this->getIsOptional()) { + out << "optional"; + previousEntry = true; + } if (!this->validators.empty()) { + if (previousEntry) { + out << "; "; + } for (uint64_t i = 0; i < this->validators.size(); ++i) { out << this->validators[i]->toString(); if (i + 1 < this->validators.size()) { out << ", "; } } - - if (this->hasDefaultValue) { - out << "; "; - } + previousEntry = true; } if (this->hasDefaultValue) { + if (previousEntry) { + out << "; "; + } out << "default: "; printValue(out, defaultValue); } diff --git a/src/storm/settings/SettingsManager.cpp b/src/storm/settings/SettingsManager.cpp index 7088f924e..058d92d99 100644 --- a/src/storm/settings/SettingsManager.cpp +++ b/src/storm/settings/SettingsManager.cpp @@ -491,7 +491,7 @@ namespace storm { // In case there are optional arguments that were not set, we set them to their default value. for (uint_fast64_t i = argumentCache.size(); i < option->getArgumentCount(); ++i) { ArgumentBase& argument = option->getArgument(i); - STORM_LOG_THROW(argument.getHasDefaultValue() || argument.getIsOptional(), storm::exceptions::OptionParserException, "Non-optional argument <" << argument.getName() << "> of option:\n" << *option); + STORM_LOG_THROW(argument.getIsOptional(), storm::exceptions::OptionParserException, "Non-optional argument <" << argument.getName() << "> of option:\n" << *option); argument.setFromDefaultValue(); } From ba6358c3fa26f1cc23a034170a640367b7cd5351 Mon Sep 17 00:00:00 2001 From: Matthias Volk Date: Fri, 13 Dec 2019 16:09:00 +0100 Subject: [PATCH 22/24] Set optional arguments for settings --- .../settings/modules/ConversionGeneralSettings.cpp | 2 +- .../settings/modules/ConversionInputSettings.cpp | 4 ++-- src/storm-dft/settings/modules/DftIOSettings.cpp | 2 +- src/storm-pars/settings/modules/RegionSettings.cpp | 4 ++-- src/storm-pgcl/settings/modules/PGCLSettings.cpp | 2 +- src/storm/settings/ArgumentBuilder.h | 7 +++---- src/storm/settings/modules/EliminationSettings.cpp | 2 +- src/storm/settings/modules/GeneralSettings.cpp | 2 +- src/storm/settings/modules/IOSettings.cpp | 4 ++-- 9 files changed, 14 insertions(+), 15 deletions(-) diff --git a/src/storm-conv/settings/modules/ConversionGeneralSettings.cpp b/src/storm-conv/settings/modules/ConversionGeneralSettings.cpp index 19ab31765..ac872ebc9 100644 --- a/src/storm-conv/settings/modules/ConversionGeneralSettings.cpp +++ b/src/storm-conv/settings/modules/ConversionGeneralSettings.cpp @@ -23,7 +23,7 @@ namespace storm { ConversionGeneralSettings::ConversionGeneralSettings() : ModuleSettings(moduleName) { this->addOption(storm::settings::OptionBuilder(moduleName, helpOptionName, false, "Shows available options, arguments and descriptions.").setShortName(helpOptionShortName) - .addArgument(storm::settings::ArgumentBuilder::createStringArgument("filter", "'frequent' for frequently used options, 'all' for the complete help, or a regular expression to show help for all matching entities.").setDefaultValueString("frequent").build()).build()); + .addArgument(storm::settings::ArgumentBuilder::createStringArgument("filter", "'frequent' for frequently used options, 'all' for the complete help, or a regular expression to show help for all matching entities.").setDefaultValueString("frequent").makeOptional().build()).build()); this->addOption(storm::settings::OptionBuilder(moduleName, versionOptionName, false, "Prints the version information.").build()); this->addOption(storm::settings::OptionBuilder(moduleName, verboseOptionName, false, "Enables more verbose output.").setShortName(verboseOptionShortName).build()); this->addOption(storm::settings::OptionBuilder(moduleName, debugOptionName, false, "Enables verbose and debug output.").build()); diff --git a/src/storm-conv/settings/modules/ConversionInputSettings.cpp b/src/storm-conv/settings/modules/ConversionInputSettings.cpp index 62093a7a7..aa6dc27fa 100644 --- a/src/storm-conv/settings/modules/ConversionInputSettings.cpp +++ b/src/storm-conv/settings/modules/ConversionInputSettings.cpp @@ -30,7 +30,7 @@ namespace storm { // General this->addOption(storm::settings::OptionBuilder(moduleName, propertyOptionName, false, "Specifies the properties to be checked on the model.").setShortName(propertyOptionShortName) .addArgument(storm::settings::ArgumentBuilder::createStringArgument("property or filename", "The formula or the file containing the formulas.").build()) - .addArgument(storm::settings::ArgumentBuilder::createStringArgument("filter", "The names of the properties to check.").setDefaultValueString("all").build()) + .addArgument(storm::settings::ArgumentBuilder::createStringArgument("filter", "The names of the properties to check.").setDefaultValueString("all").makeOptional().build()) .build()); this->addOption(storm::settings::OptionBuilder(moduleName, constantsOptionName, false, "Specifies the constant replacements to use in symbolic models. Note that this requires the model to be given as an symbolic model (e.g., via --" + prismInputOptionName + ").").setShortName(constantsOptionShortName) .addArgument(storm::settings::ArgumentBuilder::createStringArgument("values", "A comma separated list of constants and their value, e.g. a=1,b=2,c=3.").setDefaultValueString("").build()).build()); @@ -44,7 +44,7 @@ namespace storm { this->addOption(storm::settings::OptionBuilder(moduleName, janiInputOptionName, false, "Parses the model given in the JANI format.") .addArgument(storm::settings::ArgumentBuilder::createStringArgument("filename", "The name of the file from which to read the JANI input.").addValidatorString(ArgumentValidatorFactory::createExistingFileValidator()).build()).build()); this->addOption(storm::settings::OptionBuilder(moduleName, janiPropertyOptionName, false, "Specifies the properties from the jani model (given by --" + janiInputOptionName + ") to be checked.").setShortName(janiPropertyOptionShortName) - .addArgument(storm::settings::ArgumentBuilder::createStringArgument("values", "A comma separated list of properties to be checked").setDefaultValueString("").build()).build()); + .addArgument(storm::settings::ArgumentBuilder::createStringArgument("values", "A comma separated list of properties to be checked").setDefaultValueString("").makeOptional().makeOptional().build()).build()); } diff --git a/src/storm-dft/settings/modules/DftIOSettings.cpp b/src/storm-dft/settings/modules/DftIOSettings.cpp index 6a039341a..662fd0f34 100644 --- a/src/storm-dft/settings/modules/DftIOSettings.cpp +++ b/src/storm-dft/settings/modules/DftIOSettings.cpp @@ -50,7 +50,7 @@ namespace storm { .addArgument(storm::settings::ArgumentBuilder::createDoubleArgument("endtime", "The timebound to end with.") .addValidatorDouble(ArgumentValidatorFactory::createDoubleGreaterEqualValidator(0.0)).build()) .addArgument(storm::settings::ArgumentBuilder::createDoubleArgument("inc", "The value to increment with to get the next timepoint.") - .addValidatorDouble(ArgumentValidatorFactory::createDoubleGreaterEqualValidator(0.0)).build()) + .addValidatorDouble(ArgumentValidatorFactory::createDoubleGreaterEqualValidator(0.0)).setDefaultValueDouble(1.0).makeOptional().build()) .build()); this->addOption(storm::settings::OptionBuilder(moduleName, minValueOptionName, false, "Compute minimal value in case of non-determinism.").build()); this->addOption(storm::settings::OptionBuilder(moduleName, maxValueOptionName, false, "Compute maximal value in case of non-determinism.").build()); diff --git a/src/storm-pars/settings/modules/RegionSettings.cpp b/src/storm-pars/settings/modules/RegionSettings.cpp index e572560f7..cad9f4abe 100644 --- a/src/storm-pars/settings/modules/RegionSettings.cpp +++ b/src/storm-pars/settings/modules/RegionSettings.cpp @@ -35,12 +35,12 @@ namespace storm { this->addOption(storm::settings::OptionBuilder(moduleName, refineOptionName, false, "Enables region refinement.") .addArgument(storm::settings::ArgumentBuilder::createDoubleArgument("coverage-threshold", "Refinement converges if the fraction of unknown area falls below this threshold.").setDefaultValueDouble(0.05).addValidatorDouble(storm::settings::ArgumentValidatorFactory::createDoubleRangeValidatorIncluding(0.0,1.0)).build()) - .addArgument(storm::settings::ArgumentBuilder::createIntegerArgument("depth-limit", "If given, limits the number of times a region is refined.").setDefaultValueInteger(-1).build()).build()); + .addArgument(storm::settings::ArgumentBuilder::createIntegerArgument("depth-limit", "If given, limits the number of times a region is refined.").setDefaultValueInteger(-1).makeOptional().build()).build()); std::vector directions = {"min", "max"}; this->addOption(storm::settings::OptionBuilder(moduleName, extremumOptionName, false, "Computes the extremum within the region.") .addArgument(storm::settings::ArgumentBuilder::createStringArgument("direction", "The optimization direction").addValidatorString(storm::settings::ArgumentValidatorFactory::createMultipleChoiceValidator(directions)).build()) - .addArgument(storm::settings::ArgumentBuilder::createDoubleArgument("precision", "The desired precision").setDefaultValueDouble(0.05).addValidatorDouble(storm::settings::ArgumentValidatorFactory::createDoubleRangeValidatorIncluding(0.0,1.0)).build()).build()); + .addArgument(storm::settings::ArgumentBuilder::createDoubleArgument("precision", "The desired precision").setDefaultValueDouble(0.05).makeOptional().addValidatorDouble(storm::settings::ArgumentValidatorFactory::createDoubleRangeValidatorIncluding(0.0,1.0)).build()).build()); std::vector engines = {"pl", "exactpl", "validatingpl"}; this->addOption(storm::settings::OptionBuilder(moduleName, checkEngineOptionName, true, "Sets which engine is used for analyzing regions.") diff --git a/src/storm-pgcl/settings/modules/PGCLSettings.cpp b/src/storm-pgcl/settings/modules/PGCLSettings.cpp index 16b9de190..deedb85e7 100644 --- a/src/storm-pgcl/settings/modules/PGCLSettings.cpp +++ b/src/storm-pgcl/settings/modules/PGCLSettings.cpp @@ -32,7 +32,7 @@ namespace storm { this->addOption(storm::settings::OptionBuilder(moduleName, programVariableRestrictionsOptionName, false, "Restrictions of program variables").setShortName(programVariableRestrictionShortOptionName).addArgument(storm::settings::ArgumentBuilder::createStringArgument("description", "description of the variable restrictions").build()).build()); this->addOption(storm::settings::OptionBuilder(moduleName, propertyOptionName, false, "Specifies the properties to be checked on the model.").setShortName(propertyOptionShortName) .addArgument(storm::settings::ArgumentBuilder::createStringArgument("property or filename", "The formula or the file containing the formulas.").build()) - .addArgument(storm::settings::ArgumentBuilder::createStringArgument("filter", "The names of the properties to check.").setDefaultValueString("all").build()) + .addArgument(storm::settings::ArgumentBuilder::createStringArgument("filter", "The names of the properties to check.").setDefaultValueString("all").makeOptional().build()) .build()); } diff --git a/src/storm/settings/ArgumentBuilder.h b/src/storm/settings/ArgumentBuilder.h index 1c1e0e723..08a9b3219 100644 --- a/src/storm/settings/ArgumentBuilder.h +++ b/src/storm/settings/ArgumentBuilder.h @@ -89,13 +89,12 @@ namespace storm { } /*! - * Sets whether the argument is to be optional. + * Make the argument optional. * - * @param isOptional A flag that indicates whether the argument is to be optional. * @return A reference to the argument builder. */ - ArgumentBuilder& setIsOptional(bool isOptional) { - this->isOptional = isOptional; + ArgumentBuilder& makeOptional() { + this->isOptional = true; STORM_LOG_THROW(this->hasDefaultValue, storm::exceptions::IllegalFunctionCallException, "Unable to make argument '" << this->name << "' optional without default value."); return *this; } diff --git a/src/storm/settings/modules/EliminationSettings.cpp b/src/storm/settings/modules/EliminationSettings.cpp index 38b9faa8c..77c0382b5 100644 --- a/src/storm/settings/modules/EliminationSettings.cpp +++ b/src/storm/settings/modules/EliminationSettings.cpp @@ -28,7 +28,7 @@ namespace storm { this->addOption(storm::settings::OptionBuilder(moduleName, entryStatesLastOptionName, true, "Sets whether the entry states are eliminated last.").setIsAdvanced().build()); this->addOption(storm::settings::OptionBuilder(moduleName, maximalSccSizeOptionName, true, "Sets the maximal size of the SCCs for which state elimination is applied.").setIsAdvanced() - .addArgument(storm::settings::ArgumentBuilder::createUnsignedIntegerArgument("maxsize", "The maximal size of an SCC on which state elimination is applied.").setDefaultValueUnsignedInteger(20).setIsOptional(true).build()).build()); + .addArgument(storm::settings::ArgumentBuilder::createUnsignedIntegerArgument("maxsize", "The maximal size of an SCC on which state elimination is applied.").setDefaultValueUnsignedInteger(20).build()).build()); this->addOption(storm::settings::OptionBuilder(moduleName, useDedicatedModelCheckerOptionName, true, "Sets whether to use the dedicated model elimination checker (only DTMCs).").setIsAdvanced().build()); } diff --git a/src/storm/settings/modules/GeneralSettings.cpp b/src/storm/settings/modules/GeneralSettings.cpp index 05fe6c542..67244b163 100644 --- a/src/storm/settings/modules/GeneralSettings.cpp +++ b/src/storm/settings/modules/GeneralSettings.cpp @@ -35,7 +35,7 @@ namespace storm { GeneralSettings::GeneralSettings() : ModuleSettings(moduleName) { this->addOption(storm::settings::OptionBuilder(moduleName, helpOptionName, false, "Shows available options, arguments and descriptions.").setShortName(helpOptionShortName) - .addArgument(storm::settings::ArgumentBuilder::createStringArgument("filter", "'frequent' for frequently used options, 'all' for the complete help, or a regular expression to show help for all matching entities.").setDefaultValueString("frequent").build()).build()); + .addArgument(storm::settings::ArgumentBuilder::createStringArgument("filter", "'frequent' for frequently used options, 'all' for the complete help, or a regular expression to show help for all matching entities.").setDefaultValueString("frequent").makeOptional().build()).build()); this->addOption(storm::settings::OptionBuilder(moduleName, versionOptionName, false, "Prints the version information.").build()); this->addOption(storm::settings::OptionBuilder(moduleName, verboseOptionName, false, "Enables more verbose output.").setShortName(verboseOptionShortName).build()); this->addOption(storm::settings::OptionBuilder(moduleName, showProgressOptionName, false, "Sets when additional information (if available) about the progress is printed.").addArgument(storm::settings::ArgumentBuilder::createUnsignedIntegerArgument("delay", "The delay to wait (in seconds) between emitting information (0 means never print progress).").setDefaultValueUnsignedInteger(5).build()).build()); diff --git a/src/storm/settings/modules/IOSettings.cpp b/src/storm/settings/modules/IOSettings.cpp index 7f78e93cb..b6753760b 100644 --- a/src/storm/settings/modules/IOSettings.cpp +++ b/src/storm/settings/modules/IOSettings.cpp @@ -81,7 +81,7 @@ namespace storm { this->addOption(storm::settings::OptionBuilder(moduleName, prismToJaniOptionName, false, "If set, the input PRISM model is transformed to JANI.").setIsAdvanced().build()); this->addOption(storm::settings::OptionBuilder(moduleName, propertyOptionName, false, "Specifies the properties to be checked on the model.").setShortName(propertyOptionShortName) .addArgument(storm::settings::ArgumentBuilder::createStringArgument("property or filename", "The formula or the file containing the formulas.").build()) - .addArgument(storm::settings::ArgumentBuilder::createStringArgument("filter", "The names of the properties to check.").setDefaultValueString("all").build()) + .addArgument(storm::settings::ArgumentBuilder::createStringArgument("filter", "The names of the properties to check.").setDefaultValueString("all").makeOptional().build()) .build()); this->addOption(storm::settings::OptionBuilder(moduleName, transitionRewardsOptionName, false, "If given, the transition rewards are read from this file and added to the explicit model. Note that this requires the model to be given as an explicit model (i.e., via --" + explicitOptionName + ").") @@ -93,7 +93,7 @@ namespace storm { this->addOption(storm::settings::OptionBuilder(moduleName, constantsOptionName, false, "Specifies the constant replacements to use in symbolic models. Note that this requires the model to be given as an symbolic model (i.e., via --" + prismInputOptionName + " or --" + janiInputOptionName + ").").setShortName(constantsOptionShortName) .addArgument(storm::settings::ArgumentBuilder::createStringArgument("values", "A comma separated list of constants and their value, e.g. a=1,b=2,c=3.").setDefaultValueString("").build()).build()); this->addOption(storm::settings::OptionBuilder(moduleName, janiPropertyOptionName, false, "Specifies the properties from the jani model (given by --" + janiInputOptionName + ") to be checked.").setShortName(janiPropertyOptionShortName) - .addArgument(storm::settings::ArgumentBuilder::createStringArgument("values", "A comma separated list of properties to be checked").setDefaultValueString("").build()).build()); + .addArgument(storm::settings::ArgumentBuilder::createStringArgument("values", "A comma separated list of properties to be checked").setDefaultValueString("").makeOptional().build()).build()); this->addOption(storm::settings::OptionBuilder(moduleName, qvbsInputOptionName, false, "Selects a model from the Quantitative Verification Benchmark Set.").setShortName(qvbsInputOptionShortName) .addArgument(storm::settings::ArgumentBuilder::createStringArgument("model", "The short model name as in the benchmark set.").build()) From add8193dd44e2a1a3fedd2ef398d02f6f886a490 Mon Sep 17 00:00:00 2001 From: Matthias Volk Date: Fri, 13 Dec 2019 16:16:42 +0100 Subject: [PATCH 23/24] Removed duplicate makeOptional() --- src/storm-conv/settings/modules/ConversionInputSettings.cpp | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/storm-conv/settings/modules/ConversionInputSettings.cpp b/src/storm-conv/settings/modules/ConversionInputSettings.cpp index aa6dc27fa..ec8dcdd9f 100644 --- a/src/storm-conv/settings/modules/ConversionInputSettings.cpp +++ b/src/storm-conv/settings/modules/ConversionInputSettings.cpp @@ -44,7 +44,7 @@ namespace storm { this->addOption(storm::settings::OptionBuilder(moduleName, janiInputOptionName, false, "Parses the model given in the JANI format.") .addArgument(storm::settings::ArgumentBuilder::createStringArgument("filename", "The name of the file from which to read the JANI input.").addValidatorString(ArgumentValidatorFactory::createExistingFileValidator()).build()).build()); this->addOption(storm::settings::OptionBuilder(moduleName, janiPropertyOptionName, false, "Specifies the properties from the jani model (given by --" + janiInputOptionName + ") to be checked.").setShortName(janiPropertyOptionShortName) - .addArgument(storm::settings::ArgumentBuilder::createStringArgument("values", "A comma separated list of properties to be checked").setDefaultValueString("").makeOptional().makeOptional().build()).build()); + .addArgument(storm::settings::ArgumentBuilder::createStringArgument("values", "A comma separated list of properties to be checked").setDefaultValueString("").makeOptional().build()).build()); } From 743b01f00aa70182a36385fa5ba1228b7507e12e Mon Sep 17 00:00:00 2001 From: Sebastian Junges Date: Thu, 19 Dec 2019 11:07:31 +0100 Subject: [PATCH 24/24] add badges for zenodo and version --- README.md | 3 +++ 1 file changed, 3 insertions(+) diff --git a/README.md b/README.md index a8713c02c..ffc16915e 100644 --- a/README.md +++ b/README.md @@ -1,6 +1,9 @@ Storm - A Modern Probabilistic Model Checker ============================================ [![Build Status](https://travis-ci.org/moves-rwth/storm.svg?branch=master)](https://travis-ci.org/moves-rwth/storm) +[![GitHub release](https://img.shields.io/github/release/moves-rwth/storm.svg)](https://github.com/moves-rwth/storm/releases/) +[![DOI](https://zenodo.org/badge/DOI/10.5281/zenodo.3574502.svg)](https://doi.org/10.5281/zenodo.3574502) + For more instructions, check out the documentation found in [Getting Started](http://www.stormchecker.org/getting-started.html).