diff --git a/CHANGELOG.md b/CHANGELOG.md index d9b1274da..43b7c547a 100644 --- a/CHANGELOG.md +++ b/CHANGELOG.md @@ -8,7 +8,7 @@ The releases of major and minor versions contain an overview of changes since th Version 1.6.x ------------- -## Version 1.6.1 (??) +## Version 1.6.2 (2020/09) - Prism program simplification improved. - Revamped implementation of long-run-average algorithms, including scheduler export for LRA properties on Markov automata. - Support for step-bounded properties of the form ... [F[x,y] ... ] for DTMCs and MDPs (sparse engine). diff --git a/CMakeLists.txt b/CMakeLists.txt index 3d882f18c..edf9fa3e3 100644 --- a/CMakeLists.txt +++ b/CMakeLists.txt @@ -61,6 +61,9 @@ export_option(STORM_USE_CLN_RF) option(BUILD_SHARED_LIBS "Build the Storm library dynamically" OFF) option(STORM_DEBUG_CUDD "Build CUDD in debug mode." OFF) MARK_AS_ADVANCED(STORM_DEBUG_CUDD) +option(STORM_EXCLUDE_TESTS_FROM_ALL "If set, tests will not be compiled by default" OFF ) +export_option(STORM_EXCLUDE_TESTS_FROM_ALL) +MARK_AS_ADVANCED(STORM_EXCLUDE_TESTS_FROM_ALL) set(BOOST_ROOT "" CACHE STRING "A hint to the root directory of Boost (optional).") set(GUROBI_ROOT "" CACHE STRING "A hint to the root directory of Gurobi (optional).") set(Z3_ROOT "" CACHE STRING "A hint to the root directory of Z3 (optional).") diff --git a/doc/checklist_new_release.md b/doc/checklist_new_release.md index dcb1351a1..2735ca46d 100644 --- a/doc/checklist_new_release.md +++ b/doc/checklist_new_release.md @@ -34,11 +34,18 @@ 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 `stable` branch: + ```console git checkout stable - git merge master + git rebase master git push origin stable ``` + Note: Rebasing might fail if `stable` is ahead of `master` (e.g. because of merge commits). In this case we can do: + ```console + git checkout stable + git reset --hard master + git push --force origin stable + ``` 8. Update [Homebrew formula](https://github.com/moves-rwth/homebrew-storm). diff --git a/resources/3rdparty/CMakeLists.txt b/resources/3rdparty/CMakeLists.txt index d939ab369..3379689bf 100644 --- a/resources/3rdparty/CMakeLists.txt +++ b/resources/3rdparty/CMakeLists.txt @@ -260,14 +260,26 @@ if (NOT STORM_FORCE_SHIPPED_CARL) find_package(carl QUIET) endif() endif() + +set(STORM_SHIPPED_CARL OFF) if(carl_FOUND AND NOT STORM_FORCE_SHIPPED_CARL) get_target_property(carlLOCATION lib_carl LOCATION) if("${carlLOCATION}" STREQUAL "carlLOCATION-NOTFOUND") - message(FATAL_ERROR "Library location for carl is not found, did you build carl?") + if(EXISTS ${STORM_3RDPARTY_BINARY_DIR}/carl) + message(WARNING "Storm - Library for carl location is not found but the directory ${STORM_3RDPARTY_BINARY_DIR}/carl exists. Will (re-)try to build a shipped version of carl.") + set(STORM_SHIPPED_CARL ON) + else() + message(FATAL_ERROR "Library location for carl is not found, did you build carl?") + endif() elseif(EXISTS ${carlLOCATION}) #empty on purpose else() - message(FATAL_ERROR "File ${carlLOCATION} does not exist, did you build carl?") + if(EXISTS ${STORM_3RDPARTY_BINARY_DIR}/carl) + message(WARNING "Storm - File ${carlLOCATION} does not exist but the directory ${STORM_3RDPARTY_BINARY_DIR}/carl exists. Will (re-)try to build a shipped version of carl.") + set(STORM_SHIPPED_CARL ON) + else() + message(FATAL_ERROR "File ${carlLOCATION} does not exist, did you build carl?") + endif() endif() if("${carl_VERSION_MAJOR}" STREQUAL "${CARL_C14VERSION}") message(STATUS "Storm - Found carl using master14 branch.") @@ -276,14 +288,16 @@ if(carl_FOUND AND NOT STORM_FORCE_SHIPPED_CARL) message(FATAL_ERROR "Carl outdated, require ${CARL_MINVERSION}, have ${carl_VERSION}") endif() - set(STORM_SHIPPED_CARL OFF) set(STORM_HAVE_CARL ON) message(STATUS "Storm - Use system version of carl.") message(STATUS "Storm - Linking with preinstalled carl ${carl_VERSION} (include: ${carl_INCLUDE_DIR}, library ${carl_LIBRARIES}, CARL_USE_CLN_NUMBERS: ${CARL_USE_CLN_NUMBERS}, CARL_USE_GINAC: ${CARL_USE_GINAC}).") set(STORM_HAVE_CLN ${CARL_USE_CLN_NUMBERS}) set(STORM_HAVE_GINAC ${CARL_USE_GINAC}) else() - set(STORM_SHIPPED_CARL ON) + set(STORM_SHIPPED_CARL ON) +endif() + +if (STORM_SHIPPED_CARL) # The first external project will be built at *configure stage* message(STATUS "Carl - Start of config process") file(MAKE_DIRECTORY ${STORM_3RDPARTY_BINARY_DIR}/carl_download) @@ -333,7 +347,11 @@ else() message(STATUS "Storm - Linking with shipped carl ${carl_VERSION} (include: ${carl_INCLUDE_DIR}, library ${carl_LIBRARIES}, CARL_USE_CLN_NUMBERS: ${CARL_USE_CLN_NUMBERS}, CARL_USE_GINAC: ${CARL_USE_GINAC}).") # install the carl dynamic library if we built it - install(FILES ${STORM_3RDPARTY_BINARY_DIR}/carl/lib/libcarl.${carl_VERSION}${DYNAMIC_EXT} DESTINATION lib) + if(MACOSX) + install(FILES ${STORM_3RDPARTY_BINARY_DIR}/carl/lib/libcarl.${carl_VERSION}${DYNAMIC_EXT} DESTINATION lib) + else() + install(FILES ${STORM_3RDPARTY_BINARY_DIR}/carl/lib/libcarl${DYNAMIC_EXT}.${carl_VERSION} DESTINATION lib) + endif() endif() if(STORM_USE_CLN_RF AND NOT STORM_HAVE_CLN) diff --git a/src/CMakeLists.txt b/src/CMakeLists.txt index 24a213e76..2f0e8ac34 100644 --- a/src/CMakeLists.txt +++ b/src/CMakeLists.txt @@ -23,7 +23,10 @@ add_subdirectory(storm-pomdp-cli) add_subdirectory(storm-conv) add_subdirectory(storm-conv-cli) - -add_subdirectory(test) +if (STORM_EXCLUDE_TESTS_FROM_ALL) + add_subdirectory(test EXCLUDE_FROM_ALL) +else() + add_subdirectory(test) +endif() set(STORM_TARGETS ${STORM_TARGETS} PARENT_SCOPE) diff --git a/src/storm-pars/modelchecker/instantiation/SparseDtmcInstantiationModelChecker.cpp b/src/storm-pars/modelchecker/instantiation/SparseDtmcInstantiationModelChecker.cpp index 5de3d8fd3..6c743422a 100644 --- a/src/storm-pars/modelchecker/instantiation/SparseDtmcInstantiationModelChecker.cpp +++ b/src/storm-pars/modelchecker/instantiation/SparseDtmcInstantiationModelChecker.cpp @@ -42,8 +42,29 @@ namespace storm { this->currentCheckTask->setHint(std::make_shared>()); } ExplicitModelCheckerHint& hint = this->currentCheckTask->getHint().template asExplicitModelCheckerHint(); - std::unique_ptr result; + if (this->getInstantiationsAreGraphPreserving() && !hint.hasMaybeStates()) { + // Perform purely qualitative analysis once + std::vector qualitativeResult; + if (this->currentCheckTask->getFormula().asOperatorFormula().hasQuantitativeResult()) { + auto newCheckTask = *this->currentCheckTask; + newCheckTask.setQualitative(true); + newCheckTask.setOnlyInitialStatesRelevant(false); + qualitativeResult = modelChecker.check(env, newCheckTask)->template asExplicitQuantitativeCheckResult().getValueVector(); + } else { + auto newCheckTask = this->currentCheckTask->substituteFormula(this->currentCheckTask->getFormula().asOperatorFormula().getSubformula()); + newCheckTask.setQualitative(true); + newCheckTask.setOnlyInitialStatesRelevant(false); + qualitativeResult = modelChecker.computeProbabilities(env, newCheckTask)->template asExplicitQuantitativeCheckResult().getValueVector(); + } + storm::storage::BitVector maybeStates = storm::utility::vector::filter(qualitativeResult, + [] (ConstantType const& value) -> bool { return !(storm::utility::isZero(value) || storm::utility::isOne(value)); }); + hint.setMaybeStates(std::move(maybeStates)); + hint.setResultHint(std::move(qualitativeResult)); + hint.setComputeOnlyMaybeStates(true); + } + + std::unique_ptr result; // Check the formula and store the result as a hint for the next call. // For qualitative properties, we still want a quantitative result hint. Hence we perform the check on the subformula if (this->currentCheckTask->getFormula().asOperatorFormula().hasQuantitativeResult()) { @@ -56,15 +77,6 @@ namespace storm { hint.setResultHint(std::move(quantitativeResult->template asExplicitQuantitativeCheckResult().getValueVector())); } - if (this->getInstantiationsAreGraphPreserving() && !hint.hasMaybeStates()) { - // Extract the maybe states from the current result. - assert(hint.hasResultHint()); - storm::storage::BitVector maybeStates = ~storm::utility::vector::filter(hint.getResultHint(), - [] (ConstantType const& value) -> bool { return storm::utility::isZero(value) || storm::utility::isOne(value); }); - hint.setMaybeStates(std::move(maybeStates)); - hint.setComputeOnlyMaybeStates(true); - } - return result; } @@ -75,6 +87,28 @@ namespace storm { this->currentCheckTask->setHint(std::make_shared>()); } ExplicitModelCheckerHint& hint = this->currentCheckTask->getHint().template asExplicitModelCheckerHint(); + + if (this->getInstantiationsAreGraphPreserving() && !hint.hasMaybeStates()) { + // Perform purely qualitative analysis once + std::vector qualitativeResult; + if (this->currentCheckTask->getFormula().asOperatorFormula().hasQuantitativeResult()) { + auto newCheckTask = *this->currentCheckTask; + newCheckTask.setQualitative(true); + newCheckTask.setOnlyInitialStatesRelevant(false); + qualitativeResult = modelChecker.check(env, newCheckTask)->template asExplicitQuantitativeCheckResult().getValueVector(); + } else { + auto newCheckTask = this->currentCheckTask->substituteFormula(this->currentCheckTask->getFormula().asOperatorFormula().getSubformula()); + newCheckTask.setQualitative(true); + newCheckTask.setOnlyInitialStatesRelevant(false); + qualitativeResult = modelChecker.computeRewards(env, this->currentCheckTask->getFormula().asRewardOperatorFormula().getMeasureType(), newCheckTask)->template asExplicitQuantitativeCheckResult().getValueVector(); + } + storm::storage::BitVector maybeStates = storm::utility::vector::filter(qualitativeResult, + [] (ConstantType const& value) -> bool { return !(storm::utility::isZero(value) || storm::utility::isInfinity(value)); }); + hint.setMaybeStates(std::move(maybeStates)); + hint.setResultHint(std::move(qualitativeResult)); + hint.setComputeOnlyMaybeStates(true); + } + std::unique_ptr result; // Check the formula and store the result as a hint for the next call. @@ -89,18 +123,6 @@ namespace storm { this->currentCheckTask->getHint().template asExplicitModelCheckerHint().setResultHint(std::move(quantitativeResult->template asExplicitQuantitativeCheckResult().getValueVector())); } - if (this->getInstantiationsAreGraphPreserving() && !hint.hasMaybeStates()) { - // Extract the maybe states from the current result. - assert(hint.hasResultHint()); - storm::storage::BitVector maybeStates = ~storm::utility::vector::filterInfinity(hint.getResultHint()); - // We need to exclude the target states from the maybe states. - // Note that we can not consider the states with reward zero since a valuation might set a reward to zero - std::unique_ptr subFormulaResult = modelChecker.check(env, this->currentCheckTask->getFormula().asOperatorFormula().getSubformula().asEventuallyFormula().getSubformula()); - maybeStates = maybeStates & ~(subFormulaResult->asExplicitQualitativeCheckResult().getTruthValuesVector()); - hint.setMaybeStates(std::move(maybeStates)); - hint.setComputeOnlyMaybeStates(true); - } - return result; } diff --git a/src/storm-pars/modelchecker/instantiation/SparseMdpInstantiationModelChecker.cpp b/src/storm-pars/modelchecker/instantiation/SparseMdpInstantiationModelChecker.cpp index 3dc4866d3..6bbbaa158 100644 --- a/src/storm-pars/modelchecker/instantiation/SparseMdpInstantiationModelChecker.cpp +++ b/src/storm-pars/modelchecker/instantiation/SparseMdpInstantiationModelChecker.cpp @@ -47,8 +47,31 @@ namespace storm { this->currentCheckTask->setHint(std::make_shared>()); } ExplicitModelCheckerHint& hint = this->currentCheckTask->getHint().template asExplicitModelCheckerHint(); - std::unique_ptr result; + if (this->getInstantiationsAreGraphPreserving() && !hint.hasMaybeStates()) { + // Perform purely qualitative analysis once + std::vector qualitativeResult; + if (this->currentCheckTask->getFormula().asOperatorFormula().hasQuantitativeResult()) { + auto newCheckTask = *this->currentCheckTask; + newCheckTask.setQualitative(true); + newCheckTask.setOnlyInitialStatesRelevant(false); + newCheckTask.setProduceSchedulers(false); + qualitativeResult = modelChecker.check(env, newCheckTask)->template asExplicitQuantitativeCheckResult().getValueVector(); + } else { + auto newCheckTask = this->currentCheckTask->substituteFormula(this->currentCheckTask->getFormula().asOperatorFormula().getSubformula()); + newCheckTask.setQualitative(true); + newCheckTask.setOnlyInitialStatesRelevant(false); + newCheckTask.setProduceSchedulers(false); + qualitativeResult = modelChecker.computeProbabilities(env, newCheckTask)->template asExplicitQuantitativeCheckResult().getValueVector(); + } + storm::storage::BitVector maybeStates = storm::utility::vector::filter(qualitativeResult, + [] (ConstantType const& value) -> bool { return !(storm::utility::isZero(value) || storm::utility::isOne(value)); }); + hint.setMaybeStates(std::move(maybeStates)); + hint.setResultHint(std::move(qualitativeResult)); + hint.setComputeOnlyMaybeStates(true); + } + + std::unique_ptr result; // Check the formula and store the result as a hint for the next call. // For qualitative properties, we still want a quantitative result hint. Hence we perform the check on the subformula if (this->currentCheckTask->getFormula().asOperatorFormula().hasQuantitativeResult()) { @@ -65,21 +88,6 @@ namespace storm { hint.setSchedulerHint(std::move(dynamic_cast&>(scheduler))); } - if (this->getInstantiationsAreGraphPreserving() && !hint.hasMaybeStates()) { - // Extract the maybe states from the current result. - assert(hint.hasResultHint()); - storm::storage::BitVector maybeStates = ~storm::utility::vector::filter(hint.getResultHint(), - [] (ConstantType const& value) -> bool { return storm::utility::isZero(value) || storm::utility::isOne(value); }); - hint.setMaybeStates(std::move(maybeStates)); - hint.setComputeOnlyMaybeStates(true); - - // Check if there can be end components within the maybestates - if (storm::solver::minimize(this->currentCheckTask->getOptimizationDirection()) || - storm::utility::graph::performProb1A(instantiatedModel.getTransitionMatrix(), instantiatedModel.getTransitionMatrix().getRowGroupIndices(), instantiatedModel.getBackwardTransitions(), hint.getMaybeStates(), ~hint.getMaybeStates()).full()) { - hint.setNoEndComponentsInMaybeStates(true); - } - } - return result; } @@ -92,8 +100,31 @@ namespace storm { this->currentCheckTask->setHint(std::make_shared>()); } ExplicitModelCheckerHint& hint = this->currentCheckTask->getHint().template asExplicitModelCheckerHint(); - std::unique_ptr result; + if (this->getInstantiationsAreGraphPreserving() && !hint.hasMaybeStates()) { + // Perform purely qualitative analysis once + std::vector qualitativeResult; + if (this->currentCheckTask->getFormula().asOperatorFormula().hasQuantitativeResult()) { + auto newCheckTask = *this->currentCheckTask; + newCheckTask.setQualitative(true); + newCheckTask.setOnlyInitialStatesRelevant(false); + newCheckTask.setProduceSchedulers(false); + qualitativeResult = modelChecker.check(env, newCheckTask)->template asExplicitQuantitativeCheckResult().getValueVector(); + } else { + auto newCheckTask = this->currentCheckTask->substituteFormula(this->currentCheckTask->getFormula().asOperatorFormula().getSubformula()); + newCheckTask.setQualitative(true); + newCheckTask.setOnlyInitialStatesRelevant(false); + newCheckTask.setProduceSchedulers(false); + qualitativeResult = modelChecker.computeRewards(env, this->currentCheckTask->getFormula().asRewardOperatorFormula().getMeasureType(), newCheckTask)->template asExplicitQuantitativeCheckResult().getValueVector(); + } + storm::storage::BitVector maybeStates = storm::utility::vector::filter(qualitativeResult, + [] (ConstantType const& value) -> bool { return !(storm::utility::isZero(value) || storm::utility::isInfinity(value)); }); + hint.setMaybeStates(std::move(maybeStates)); + hint.setResultHint(std::move(qualitativeResult)); + hint.setComputeOnlyMaybeStates(true); + } + + std::unique_ptr result; // Check the formula and store the result as a hint for the next call. // For qualitative properties, we still want a quantitative result hint. Hence we perform the check on the subformula if(this->currentCheckTask->getFormula().asOperatorFormula().hasQuantitativeResult()) { @@ -110,23 +141,6 @@ namespace storm { hint.setSchedulerHint(std::move(dynamic_cast&>(scheduler))); } - if (this->getInstantiationsAreGraphPreserving() && !hint.hasMaybeStates()) { - // Extract the maybe states from the current result. - assert(hint.hasResultHint()); - storm::storage::BitVector maybeStates = ~storm::utility::vector::filterInfinity(hint.getResultHint()); - // We need to exclude the target states from the maybe states. - // Note that we can not consider the states with reward zero since a valuation might set a reward to zero - std::unique_ptr subFormulaResult = modelChecker.check(env, this->currentCheckTask->getFormula().asOperatorFormula().getSubformula().asEventuallyFormula().getSubformula()); - maybeStates = maybeStates & ~(subFormulaResult->asExplicitQualitativeCheckResult().getTruthValuesVector()); - hint.setMaybeStates(std::move(maybeStates)); - hint.setComputeOnlyMaybeStates(true); - - // Check if there can be end components within the maybestates - if (storm::solver::maximize(this->currentCheckTask->getOptimizationDirection()) || - storm::utility::graph::performProb1A(instantiatedModel.getTransitionMatrix(), instantiatedModel.getTransitionMatrix().getRowGroupIndices(), instantiatedModel.getBackwardTransitions(), hint.getMaybeStates(), ~hint.getMaybeStates()).full()) { - hint.setNoEndComponentsInMaybeStates(true); - } - } return result; } diff --git a/src/storm-pars/storage/ParameterRegion.cpp b/src/storm-pars/storage/ParameterRegion.cpp index 68d6e94dc..71030ff4d 100644 --- a/src/storm-pars/storage/ParameterRegion.cpp +++ b/src/storm-pars/storage/ParameterRegion.cpp @@ -1,7 +1,10 @@ #include "storm-pars/storage/ParameterRegion.h" +#include + #include "storm/utility/macros.h" #include "storm/exceptions/InvalidArgumentException.h" +#include "storm/exceptions/OutOfRangeException.h" #include "storm/utility/constants.h" namespace storm { @@ -88,6 +91,7 @@ namespace storm { template std::vector::Valuation> ParameterRegion::getVerticesOfRegion(std::set const& consideredVariables) const { std::size_t const numOfVariables = consideredVariables.size(); + STORM_LOG_THROW(numOfVariables <= std::numeric_limits::digits, storm::exceptions::OutOfRangeException, "Number of variables " << numOfVariables << " is too high."); std::size_t const numOfVertices = std::pow(2, numOfVariables); std::vector resultingVector(numOfVertices); diff --git a/src/storm/environment/solver/OviSolverEnvironment.cpp b/src/storm/environment/solver/OviSolverEnvironment.cpp index fe6a21a12..382b27a54 100644 --- a/src/storm/environment/solver/OviSolverEnvironment.cpp +++ b/src/storm/environment/solver/OviSolverEnvironment.cpp @@ -11,7 +11,6 @@ namespace storm { auto const& oviSettings = storm::settings::getModule(); precisionUpdateFactor = storm::utility::convertNumber(oviSettings.getPrecisionUpdateFactor()); maxVerificationIterationFactor = storm::utility::convertNumber(oviSettings.getMaxVerificationIterationFactor()); - relevantValuesForPrecisionUpdate = oviSettings.useRelevantValuesForPrecisionUpdate(); upperBoundGuessingFactor = storm::utility::convertNumber(oviSettings.getUpperBoundGuessingFactor()); upperBoundOnlyIterations = oviSettings.getUpperBoundOnlyIterations(); noTerminationGuaranteeMinimumMethod = oviSettings.useNoTerminationGuaranteeMinimumMethod(); @@ -29,10 +28,6 @@ namespace storm { return maxVerificationIterationFactor; } - bool OviSolverEnvironment::useRelevantValuesForPrecisionUpdate() const { - return relevantValuesForPrecisionUpdate; - } - storm::RationalNumber OviSolverEnvironment::getUpperBoundGuessingFactor() const { return upperBoundGuessingFactor; } diff --git a/src/storm/environment/solver/OviSolverEnvironment.h b/src/storm/environment/solver/OviSolverEnvironment.h index 0875bad9a..207717982 100644 --- a/src/storm/environment/solver/OviSolverEnvironment.h +++ b/src/storm/environment/solver/OviSolverEnvironment.h @@ -14,7 +14,6 @@ namespace storm { storm::RationalNumber getPrecisionUpdateFactor() const; storm::RationalNumber getMaxVerificationIterationFactor() const; - bool useRelevantValuesForPrecisionUpdate() const; storm::RationalNumber getUpperBoundGuessingFactor() const; uint64_t getUpperBoundOnlyIterations() const; bool useNoTerminationGuaranteeMinimumMethod() const; @@ -22,7 +21,6 @@ namespace storm { private: storm::RationalNumber precisionUpdateFactor; storm::RationalNumber maxVerificationIterationFactor; - bool relevantValuesForPrecisionUpdate; storm::RationalNumber upperBoundGuessingFactor; uint64_t upperBoundOnlyIterations; bool noTerminationGuaranteeMinimumMethod; diff --git a/src/storm/modelchecker/CheckTask.h b/src/storm/modelchecker/CheckTask.h index 45f2a9bc3..898795fda 100644 --- a/src/storm/modelchecker/CheckTask.h +++ b/src/storm/modelchecker/CheckTask.h @@ -207,6 +207,14 @@ namespace storm { bool isQualitativeSet() const { return qualitative; } + + /*! + * sets whether the computation only needs to be performed qualitatively, because the values will only + * be compared to 0/1. + */ + void setQualitative(bool value) { + qualitative = value; + } /*! * Sets whether to produce schedulers (if supported). diff --git a/src/storm/modelchecker/multiobjective/preprocessing/SparseMultiObjectivePreprocessor.cpp b/src/storm/modelchecker/multiobjective/preprocessing/SparseMultiObjectivePreprocessor.cpp index 294a598ae..9f253bc38 100644 --- a/src/storm/modelchecker/multiobjective/preprocessing/SparseMultiObjectivePreprocessor.cpp +++ b/src/storm/modelchecker/multiobjective/preprocessing/SparseMultiObjectivePreprocessor.cpp @@ -343,8 +343,6 @@ namespace storm { template void SparseMultiObjectivePreprocessor::preprocessTimeOperatorFormula(storm::logic::TimeOperatorFormula const& formula, storm::logic::OperatorInformation const& opInfo, PreprocessorData& data) { - // Time formulas are only supported for Markov automata - STORM_LOG_THROW(data.model->isOfType(storm::models::ModelType::MarkovAutomaton), storm::exceptions::InvalidPropertyException, "Time operator formulas are only supported for Markov automata."); data.objectives.back()->lowerResultBound = storm::utility::zero(); @@ -474,11 +472,15 @@ namespace storm { } } data.model->addRewardModel(rewardModelName, std::move(objectiveRewards)); - } else if (formula.isReachabilityTimeFormula() && data.model->isOfType(storm::models::ModelType::MarkovAutomaton)) { + } else if (formula.isReachabilityTimeFormula()) { - // build stateAction reward vector that only gives reward for relevant states + // build state reward vector that only gives reward for relevant states std::vector timeRewards(data.model->getNumberOfStates(), storm::utility::zero()); - storm::utility::vector::setVectorValues(timeRewards, dynamic_cast const&>(*data.model).getMarkovianStates() & reachableFromInit, storm::utility::one()); + if (data.model->isOfType(storm::models::ModelType::MarkovAutomaton)) { + storm::utility::vector::setVectorValues(timeRewards, dynamic_cast const&>(*data.model).getMarkovianStates() & reachableFromInit, storm::utility::one()); + } else { + storm::utility::vector::setVectorValues(timeRewards, reachableFromInit, storm::utility::one()); + } data.model->addRewardModel(rewardModelName, typename SparseModelType::RewardModelType(std::move(timeRewards))); } else { STORM_LOG_THROW(false, storm::exceptions::InvalidPropertyException, "The formula " << formula << " neither considers reachability probabilities nor reachability rewards " << (data.model->isOfType(storm::models::ModelType::MarkovAutomaton) ? "nor reachability time" : "") << ". This is not supported."); diff --git a/src/storm/settings/modules/OviSolverSettings.cpp b/src/storm/settings/modules/OviSolverSettings.cpp index 931aa1af8..4d8311737 100644 --- a/src/storm/settings/modules/OviSolverSettings.cpp +++ b/src/storm/settings/modules/OviSolverSettings.cpp @@ -14,7 +14,6 @@ namespace storm { const std::string OviSolverSettings::moduleName = "ovi"; const std::string OviSolverSettings::precisionUpdateFactorOptionName = "precision-update-factor"; const std::string OviSolverSettings::maxVerificationIterationFactorOptionName = "max-verification-iter-factor"; - const std::string OviSolverSettings::useRelevantValuesForPrecisionUpdateOptionName = "use-relevant-values"; const std::string OviSolverSettings::upperBoundGuessingFactorOptionName = "upper-bound-factor"; const std::string OviSolverSettings::upperBoundOnlyIterationsOptionName = "check-upper-only-iter"; const std::string OviSolverSettings::useNoTerminationGuaranteeMinimumMethodOptionName = "no-termination-guarantee"; @@ -23,9 +22,7 @@ namespace storm { this->addOption(storm::settings::OptionBuilder(moduleName, precisionUpdateFactorOptionName, false, "Sets with which factor the precision of the inner value iteration is updated.").setIsAdvanced().addArgument(storm::settings::ArgumentBuilder::createDoubleArgument("factor", "The factor.").setDefaultValueDouble(0.4).addValidatorDouble(ArgumentValidatorFactory::createDoubleRangeValidatorExcluding(0.0, 1.0)).build()).build()); - this->addOption(storm::settings::OptionBuilder(moduleName, useRelevantValuesForPrecisionUpdateOptionName, false, "Sets whether the precision of the inner value iteration is only based on the relevant values (i.e. initial states).").setIsAdvanced().build()); - - this->addOption(storm::settings::OptionBuilder(moduleName, maxVerificationIterationFactorOptionName, false, "Controls how many verification iterations are performed before guessing a new upper bound.").setIsAdvanced().addArgument(storm::settings::ArgumentBuilder::createDoubleArgument("factor", "The factor.").setDefaultValueDouble(0.1).addValidatorDouble(ArgumentValidatorFactory::createDoubleGreaterValidator(0.0)).build()).build()); + this->addOption(storm::settings::OptionBuilder(moduleName, maxVerificationIterationFactorOptionName, false, "Controls how many verification iterations are performed before guessing a new upper bound.").setIsAdvanced().addArgument(storm::settings::ArgumentBuilder::createDoubleArgument("factor", "The factor.").setDefaultValueDouble(1.0).addValidatorDouble(ArgumentValidatorFactory::createDoubleGreaterValidator(0.0)).build()).build()); this->addOption(storm::settings::OptionBuilder(moduleName, upperBoundGuessingFactorOptionName, false, "Sets with which factor the precision is multiplied to guess the upper bound.").setIsAdvanced().addArgument(storm::settings::ArgumentBuilder::createDoubleArgument("factor", "The factor.").setDefaultValueDouble(1.0).addValidatorDouble(ArgumentValidatorFactory::createDoubleGreaterValidator(0.0)).build()).build()); @@ -42,10 +39,6 @@ namespace storm { return this->getOption(maxVerificationIterationFactorOptionName).getArgumentByName("factor").getValueAsDouble(); } - bool OviSolverSettings::useRelevantValuesForPrecisionUpdate() const { - return this->getOption(useRelevantValuesForPrecisionUpdateOptionName).getHasOptionBeenSet(); - } - double OviSolverSettings::getUpperBoundGuessingFactor() const { return this->getOption(upperBoundGuessingFactorOptionName).getArgumentByName("factor").getValueAsDouble(); } diff --git a/src/storm/settings/modules/OviSolverSettings.h b/src/storm/settings/modules/OviSolverSettings.h index 70c1c4f0f..1b033121a 100644 --- a/src/storm/settings/modules/OviSolverSettings.h +++ b/src/storm/settings/modules/OviSolverSettings.h @@ -19,8 +19,6 @@ namespace storm { double getMaxVerificationIterationFactor() const; - bool useRelevantValuesForPrecisionUpdate() const; - double getUpperBoundGuessingFactor() const; uint64_t getUpperBoundOnlyIterations() const; @@ -33,7 +31,6 @@ namespace storm { private: static const std::string precisionUpdateFactorOptionName; static const std::string maxVerificationIterationFactorOptionName; - static const std::string useRelevantValuesForPrecisionUpdateOptionName; static const std::string upperBoundGuessingFactorOptionName; static const std::string upperBoundOnlyIterationsOptionName; static const std::string useNoTerminationGuaranteeMinimumMethodOptionName; diff --git a/src/storm/solver/AbstractEquationSolver.cpp b/src/storm/solver/AbstractEquationSolver.cpp index d8e3e7b91..c068e5d31 100644 --- a/src/storm/solver/AbstractEquationSolver.cpp +++ b/src/storm/solver/AbstractEquationSolver.cpp @@ -66,6 +66,11 @@ namespace storm { return relevantValues.get(); } + template + boost::optional const& AbstractEquationSolver::getOptionalRelevantValues() const { + return relevantValues; + } + template void AbstractEquationSolver::setRelevantValues(storm::storage::BitVector&& relevantValues) { this->relevantValues = std::move(relevantValues); diff --git a/src/storm/solver/AbstractEquationSolver.h b/src/storm/solver/AbstractEquationSolver.h index 8e303f042..0301ad42e 100644 --- a/src/storm/solver/AbstractEquationSolver.h +++ b/src/storm/solver/AbstractEquationSolver.h @@ -52,6 +52,7 @@ namespace storm { * Retrieves the relevant values (if there are any). */ storm::storage::BitVector const& getRelevantValues() const; + boost::optional const& getOptionalRelevantValues() const; /*! * Sets the relevant values. diff --git a/src/storm/solver/IterativeMinMaxLinearEquationSolver.cpp b/src/storm/solver/IterativeMinMaxLinearEquationSolver.cpp index 089ba266e..5e71c0503 100644 --- a/src/storm/solver/IterativeMinMaxLinearEquationSolver.cpp +++ b/src/storm/solver/IterativeMinMaxLinearEquationSolver.cpp @@ -2,7 +2,6 @@ #include #include "storm/solver/IterativeMinMaxLinearEquationSolver.h" -#include "storm/solver/helper/OptimisticValueIterationHelper.h" #include "storm/environment/solver/MinMaxSolverEnvironment.h" #include "storm/environment/solver/OviSolverEnvironment.h" @@ -369,63 +368,28 @@ namespace storm { return true; } - if (!this->multiplierA) { - this->multiplierA = storm::solver::MultiplierFactory().create(env, *this->A); - } - if (!auxiliaryRowGroupVector) { auxiliaryRowGroupVector = std::make_unique>(this->A->getRowGroupCount()); } - if (!auxiliaryRowGroupVector2) { - auxiliaryRowGroupVector2 = std::make_unique>(this->A->getRowGroupCount()); + if (!optimisticValueIterationHelper) { + optimisticValueIterationHelper = std::make_unique>(*this->A); } - // By default, we can not provide any guarantee - SolverGuarantee guarantee = SolverGuarantee::None; - // Get handle to multiplier. - storm::solver::Multiplier const &multiplier = *this->multiplierA; - // Allow aliased multiplications. - storm::solver::MultiplicationStyle multiplicationStyle = env.solver().minMax().getMultiplicationStyle(); - bool useGaussSeidelMultiplication = multiplicationStyle == storm::solver::MultiplicationStyle::GaussSeidel; - - boost::optional relevantValues; - if (this->hasRelevantValues()) { - relevantValues = this->getRelevantValues(); - } + storm::solver::helper::OptimisticValueIterationHelper helper(*this->A); // x has to start with a lower bound. this->createLowerBoundsVector(x); std::vector* lowerX = &x; std::vector* upperX = auxiliaryRowGroupVector.get(); - std::vector* auxVector = auxiliaryRowGroupVector2.get(); - this->startMeasureProgress(); - storm::solver::helper::OptimisticValueIterationHelper helper(*this->A); - auto statusIters = helper.solveEquationsOptimisticValueIteration(env, lowerX, upperX, auxVector, b, - [&] (std::vector*& y, std::vector*& yPrime, ValueType const& precision, bool const& relative, uint64_t const& i, uint64_t const& maxI) { - this->showProgressIterative(i); - auto result = performValueIteration(env, dir, y, yPrime, b, precision, relative, guarantee, i, maxI, multiplicationStyle); - return std::make_pair(result.iterations, result.status); - }, - [&] (std::vector* y, std::vector* yPrime, uint64_t const& i) { - this->showProgressIterative(i); - if (useGaussSeidelMultiplication) { - // Copy over the current vectors so we can modify them in-place. - // This is necessary as we want to compare the new values with the current ones. - *yPrime = *y; - multiplier.multiplyAndReduceGaussSeidel(env, dir, *y, &b); - } else { - multiplier.multiplyAndReduce(env, dir, *y, &b, *yPrime); - std::swap(y, yPrime); - } - }, - env.solver().minMax().getRelativeTerminationCriterion(), - storm::utility::convertNumber(env.solver().minMax().getPrecision()), - env.solver().minMax().getMaximalNumberOfIterations(), - dir, - relevantValues); + auto statusIters = helper.solveEquations(env, lowerX, upperX, b, + env.solver().minMax().getRelativeTerminationCriterion(), + storm::utility::convertNumber(env.solver().minMax().getPrecision()), + env.solver().minMax().getMaximalNumberOfIterations(), + dir, + this->getOptionalRelevantValues()); auto two = storm::utility::convertNumber(2.0); storm::utility::vector::applyPointwise(*lowerX, *upperX, x, [&two] (ValueType const& a, ValueType const& b) -> ValueType { return (a + b) / two; }); @@ -434,8 +398,7 @@ namespace storm { // If requested, we store the scheduler for retrieval. if (this->isTrackSchedulerSet()) { this->schedulerChoices = std::vector(this->A->getRowGroupCount()); - this->multiplierA->multiplyAndReduce(env, dir, x, &b, *auxiliaryRowGroupVector.get(), &this->schedulerChoices.get()); - this->multiplierA->multiplyAndReduce(env, dir, x, &b, *auxiliaryRowGroupVector.get(), &this->schedulerChoices.get()); + this->A->multiplyAndReduce(dir, this->A->getRowGroupIndices(), x, &b, *auxiliaryRowGroupVector.get(), &this->schedulerChoices.get()); } if (!this->isCachingEnabled()) { @@ -1109,6 +1072,7 @@ namespace storm { auxiliaryRowGroupVector.reset(); auxiliaryRowGroupVector2.reset(); soundValueIterationHelper.reset(); + optimisticValueIterationHelper.reset(); StandardMinMaxLinearEquationSolver::clearCache(); } diff --git a/src/storm/solver/IterativeMinMaxLinearEquationSolver.h b/src/storm/solver/IterativeMinMaxLinearEquationSolver.h index 6582caefa..be268f35f 100644 --- a/src/storm/solver/IterativeMinMaxLinearEquationSolver.h +++ b/src/storm/solver/IterativeMinMaxLinearEquationSolver.h @@ -8,6 +8,7 @@ #include "storm/solver/Multiplier.h" #include "storm/solver/StandardMinMaxLinearEquationSolver.h" #include "storm/solver/helper/SoundValueIterationHelper.h" +#include "storm/solver/helper/OptimisticValueIterationHelper.h" #include "storm/solver/SolverStatus.h" @@ -85,6 +86,7 @@ namespace storm { mutable std::unique_ptr> auxiliaryRowGroupVector; // A.rowGroupCount() entries mutable std::unique_ptr> auxiliaryRowGroupVector2; // A.rowGroupCount() entries mutable std::unique_ptr> soundValueIterationHelper; + mutable std::unique_ptr> optimisticValueIterationHelper; }; diff --git a/src/storm/solver/NativeLinearEquationSolver.cpp b/src/storm/solver/NativeLinearEquationSolver.cpp index 9ca9683fc..82bf7d5ec 100644 --- a/src/storm/solver/NativeLinearEquationSolver.cpp +++ b/src/storm/solver/NativeLinearEquationSolver.cpp @@ -643,62 +643,26 @@ namespace storm { return true; } - if (!this->multiplier) { - this->multiplier = storm::solver::MultiplierFactory().create(env, *this->A); - } - if (!this->cachedRowVector) { this->cachedRowVector = std::make_unique>(this->A->getRowCount()); } - if (!this->cachedRowVector2) { - this->cachedRowVector2 = std::make_unique>(this->A->getRowCount()); + if (!optimisticValueIterationHelper) { + optimisticValueIterationHelper = std::make_unique>(*this->A); } - // By default, we can not provide any guarantee - SolverGuarantee guarantee = SolverGuarantee::None; - // Get handle to multiplier. - storm::solver::Multiplier const &multiplier = *this->multiplier; - // Allow aliased multiplications. - storm::solver::MultiplicationStyle multiplicationStyle = env.solver().native().getPowerMethodMultiplicationStyle(); - bool useGaussSeidelMultiplication = multiplicationStyle == storm::solver::MultiplicationStyle::GaussSeidel; - - boost::optional relevantValues; - if (this->hasRelevantValues()) { - relevantValues = this->getRelevantValues(); - } - // x has to start with a lower bound. this->createLowerBoundsVector(x); std::vector* lowerX = &x; std::vector* upperX = this->cachedRowVector.get(); - std::vector* auxVector = this->cachedRowVector2.get(); - this->startMeasureProgress(); storm::solver::helper::OptimisticValueIterationHelper helper(*this->A); - auto statusIters = helper.solveEquationsOptimisticValueIteration(env, lowerX, upperX, auxVector, b, - [&] (std::vector*& y, std::vector*& yPrime, ValueType const& precision, bool const& relative, uint64_t const& i, uint64_t const& maxI) { - this->showProgressIterative(i); - auto result = performPowerIteration(env, y, yPrime, b, precision, relative, guarantee, i, maxI, multiplicationStyle); - return std::make_pair(result.iterations, result.status); - }, - [&] (std::vector* y, std::vector* yPrime, uint64_t const& i) { - this->showProgressIterative(i); - if (useGaussSeidelMultiplication) { - // Copy over the current vectors so we can modify them in-place. - // This is necessary as we want to compare the new values with the current ones. - *yPrime = *y; - multiplier.multiplyGaussSeidel(env, *y, &b); - } else { - multiplier.multiply(env, *y, &b, *yPrime); - std::swap(y, yPrime); - } - }, - env.solver().native().getRelativeTerminationCriterion(), - storm::utility::convertNumber(env.solver().native().getPrecision()), - env.solver().native().getMaximalNumberOfIterations(), - boost::none, // No optimization dir - relevantValues); + auto statusIters = helper.solveEquations(env, lowerX, upperX, b, + env.solver().native().getRelativeTerminationCriterion(), + storm::utility::convertNumber(env.solver().native().getPrecision()), + env.solver().native().getMaximalNumberOfIterations(), + boost::none, // No optimization dir + this->getOptionalRelevantValues()); auto two = storm::utility::convertNumber(2.0); storm::utility::vector::applyPointwise(*lowerX, *upperX, x, [&two] (ValueType const& a, ValueType const& b) -> ValueType { return (a + b) / two; }); this->logIterations(statusIters.first == SolverStatus::Converged, statusIters.first == SolverStatus::TerminatedEarly, statusIters.second); @@ -1066,6 +1030,7 @@ namespace storm { walkerChaeData.reset(); multiplier.reset(); soundValueIterationHelper.reset(); + optimisticValueIterationHelper.reset(); LinearEquationSolver::clearCache(); } diff --git a/src/storm/solver/NativeLinearEquationSolver.h b/src/storm/solver/NativeLinearEquationSolver.h index d551cd931..93c627a86 100644 --- a/src/storm/solver/NativeLinearEquationSolver.h +++ b/src/storm/solver/NativeLinearEquationSolver.h @@ -9,6 +9,7 @@ #include "storm/solver/NativeMultiplier.h" #include "storm/solver/SolverStatus.h" #include "storm/solver/helper/SoundValueIterationHelper.h" +#include "storm/solver/helper/OptimisticValueIterationHelper.h" #include "storm/utility/NumberTraits.h" @@ -96,6 +97,7 @@ namespace storm { // cached auxiliary data mutable std::unique_ptr> cachedRowVector2; // A.getRowCount() rows mutable std::unique_ptr> soundValueIterationHelper; + mutable std::unique_ptr> optimisticValueIterationHelper; struct JacobiDecomposition { JacobiDecomposition(Environment const& env, storm::storage::SparseMatrix const& A); diff --git a/src/storm/solver/TopologicalLinearEquationSolver.cpp b/src/storm/solver/TopologicalLinearEquationSolver.cpp index 76e4ce3ea..5ae744203 100644 --- a/src/storm/solver/TopologicalLinearEquationSolver.cpp +++ b/src/storm/solver/TopologicalLinearEquationSolver.cpp @@ -144,11 +144,12 @@ namespace storm { } if (hasDiagonalEntry) { - if (storm::utility::isZero(denominator)) { - STORM_LOG_THROW(storm::utility::isZero(xi), storm::exceptions::InvalidOperationException, "The equation system has no solution."); - } else { - xi /= denominator; - } + STORM_LOG_WARN_COND_DEBUG(storm::NumberTraits::IsExact || !storm::utility::isAlmostZero(denominator) || storm::utility::isZero(denominator), "State " << sccState << " has a selfloop with probability '1-(" << denominator << ")'. This could be an indication for numerical issues."); + if (storm::utility::isZero(denominator)) { + STORM_LOG_THROW(storm::utility::isZero(xi), storm::exceptions::InvalidOperationException, "The equation system has no solution."); + } else { + xi /= denominator; + } } return true; } diff --git a/src/storm/solver/TopologicalMinMaxLinearEquationSolver.cpp b/src/storm/solver/TopologicalMinMaxLinearEquationSolver.cpp index 7c366e436..6cf09bd45 100644 --- a/src/storm/solver/TopologicalMinMaxLinearEquationSolver.cpp +++ b/src/storm/solver/TopologicalMinMaxLinearEquationSolver.cpp @@ -155,6 +155,7 @@ namespace storm { } } if (hasDiagonalEntry) { + STORM_LOG_WARN_COND_DEBUG(storm::NumberTraits::IsExact || !storm::utility::isAlmostZero(denominator) || storm::utility::isZero(denominator), "State " << sccState << " has a selfloop with probability '1-(" << denominator << ")'. This could be an indication for numerical issues."); if (storm::utility::isZero(denominator)) { // In this case we have a selfloop on this state. This can never an optimal choice: // When minimizing, we are looking for the largest fixpoint (which will never be attained by this action) diff --git a/src/storm/solver/helper/OptimisticValueIterationHelper.cpp b/src/storm/solver/helper/OptimisticValueIterationHelper.cpp index d52d41d37..296dd2ceb 100644 --- a/src/storm/solver/helper/OptimisticValueIterationHelper.cpp +++ b/src/storm/solver/helper/OptimisticValueIterationHelper.cpp @@ -15,56 +15,9 @@ namespace storm { namespace oviinternal { template - ValueType computeMaxAbsDiff(std::vector const& allOldValues, std::vector const& allNewValues, storm::storage::BitVector const& relevantValues) { - ValueType result = storm::utility::zero(); - for (auto value : relevantValues) { - result = storm::utility::max(result, storm::utility::abs(allNewValues[value] - allOldValues[value])); - } - return result; - } - - template - ValueType computeMaxAbsDiff(std::vector const& allOldValues, std::vector const& allNewValues) { - ValueType result = storm::utility::zero(); - for (uint64_t i = 0; i < allOldValues.size(); ++i) { - result = storm::utility::max(result, storm::utility::abs(allNewValues[i] - allOldValues[i])); - } - return result; - } - - template - ValueType computeMaxRelDiff(std::vector const& allOldValues, std::vector const& allNewValues, storm::storage::BitVector const& relevantValues) { - ValueType result = storm::utility::zero(); - for (auto const& i : relevantValues) { - STORM_LOG_ASSERT(!storm::utility::isZero(allNewValues[i]) || storm::utility::isZero(allOldValues[i]), "Unexpected entry in iteration vector."); - if (!storm::utility::isZero(allNewValues[i])) { - result = storm::utility::max(result, storm::utility::abs(allNewValues[i] - allOldValues[i]) / allNewValues[i]); - } - } - return result; - } - - template - ValueType computeMaxRelDiff(std::vector const& allOldValues, std::vector const& allNewValues) { - ValueType result = storm::utility::zero(); - for (uint64_t i = 0; i < allOldValues.size(); ++i) { - STORM_LOG_ASSERT(!storm::utility::isZero(allNewValues[i]) || storm::utility::isZero(allOldValues[i]), "Unexpected entry in iteration vector."); - if (!storm::utility::isZero(allNewValues[i])) { - result = storm::utility::max(result, storm::utility::abs(allNewValues[i] - allOldValues[i]) / allNewValues[i]); - } - } - return result; - } - - template - ValueType updateIterationPrecision(storm::Environment const& env, std::vector const& currentX, std::vector const& newX, bool const& relative, boost::optional const& relevantValues) { + ValueType updateIterationPrecision(storm::Environment const& env, ValueType const& diff) { auto factor = storm::utility::convertNumber(env.solver().ovi().getPrecisionUpdateFactor()); - bool useRelevant = relevantValues.is_initialized() && env.solver().ovi().useRelevantValuesForPrecisionUpdate(); - if (relative) { - return (useRelevant ? computeMaxRelDiff(newX, currentX, relevantValues.get()) : computeMaxRelDiff(newX, currentX)) * factor; - } else { - return (useRelevant ? computeMaxAbsDiff(newX, currentX, relevantValues.get()) : computeMaxAbsDiff(newX, currentX)) * factor; - } + return factor * diff; } template @@ -78,7 +31,7 @@ namespace storm { } template - UpperBoundIterator::UpperBoundIterator(storm::storage::SparseMatrix const& matrix) { + IterationHelper::IterationHelper(storm::storage::SparseMatrix const& matrix) { STORM_LOG_THROW(static_cast(std::numeric_limits::max()) > matrix.getRowCount() + 1, storm::exceptions::NotSupportedException, "Matrix dimensions too large."); STORM_LOG_THROW(static_cast(std::numeric_limits::max()) > matrix.getEntryCount(), storm::exceptions::NotSupportedException, "Matrix dimensions too large."); matrixValues.reserve(matrix.getNonzeroEntryCount()); @@ -98,22 +51,118 @@ namespace storm { } template - typename UpperBoundIterator::IterateResult UpperBoundIterator::iterate(std::vector& x, std::vector const& b, bool takeMinOfOldAndNew) { - return iterateInternal(x, b, takeMinOfOldAndNew); + ValueType IterationHelper::singleIterationWithDiff(std::vector& x, std::vector const& b, bool computeRelativeDiff) { + return singleIterationWithDiffInternal(x, b, computeRelativeDiff); + } + + template + ValueType IterationHelper::singleIterationWithDiff(storm::solver::OptimizationDirection const& dir, std::vector& x, std::vector const& b, bool computeRelativeDiff) { + if (minimize(dir)) { + return singleIterationWithDiffInternal(x, b, computeRelativeDiff); + } else { + return singleIterationWithDiffInternal(x, b, computeRelativeDiff); + } + } + + template + template + ValueType IterationHelper::singleIterationWithDiffInternal(std::vector& x, std::vector const& b, bool computeRelativeDiff) { + STORM_LOG_ASSERT(x.size() > 0, "Empty equation system not expected."); + ValueType diff = storm::utility::zero(); + + IndexType i = x.size(); + while (i > 0) { + --i; + ValueType newXi = HasRowGroups ? multiplyRowGroup(i, b, x) : multiplyRow(i, b[i], x); + ValueType& oldXi = x[i]; + if (computeRelativeDiff) { + if (storm::utility::isZero(newXi)) { + if (storm::utility::isZero(oldXi)) { + // this operation has no effect: + // diff = std::max(diff, storm::utility::zero()); + } else { + diff = std::max(diff, storm::utility::one()); + } + } else { + diff = std::max(diff, storm::utility::abs((newXi - oldXi) / newXi)); + } + } else { + diff = std::max(diff, storm::utility::abs(newXi - oldXi)); + } + oldXi = std::move(newXi); + } + return diff; } template - typename UpperBoundIterator::IterateResult UpperBoundIterator::iterate(storm::solver::OptimizationDirection const& dir, std::vector& x, std::vector const& b, bool takeMinOfOldAndNew) { + uint64_t IterationHelper::repeatedIterate(storm::solver::OptimizationDirection const& dir, std::vector& x, std::vector const& b, ValueType precision, bool relative) { if (minimize(dir)) { - return iterateInternal(x, b, takeMinOfOldAndNew); + return repeatedIterateInternal(x, b, precision, relative); } else { - return iterateInternal(x, b, takeMinOfOldAndNew); + return repeatedIterateInternal(x, b, precision, relative); } } + template + uint64_t IterationHelper::repeatedIterate(std::vector& x, const std::vector& b, ValueType precision, bool relative) { + return repeatedIterateInternal(x, b, precision, relative); + } + template template - typename UpperBoundIterator::IterateResult UpperBoundIterator::iterateInternal(std::vector& x, std::vector const& b, bool takeMinOfOldAndNew) { + uint64_t IterationHelper::repeatedIterateInternal(std::vector& x, std::vector const& b, ValueType precision, bool relative) { + // Do a backwards gauss-seidel style iteration + bool convergence = true; + IndexType i = x.size(); + while (i > 0) { + --i; + ValueType newXi = HasRowGroups ? multiplyRowGroup(i, b, x) : multiplyRow(i, b[i], x); + ValueType& oldXi = x[i]; + // Check if we converged + if (relative) { + if (storm::utility::isZero(oldXi)) { + if (!storm::utility::isZero(newXi)) { + convergence = false; + break; + } + } else if (storm::utility::abs((newXi - oldXi) / oldXi) > precision) { + convergence = false; + break; + } + } else { + if (storm::utility::abs((newXi - oldXi)) > precision) { + convergence = false; + break; + } + } + } + if (!convergence) { + // we now know that we did not converge. We still need to set the remaining values + while (i > 0) { + --i; + x[i] = HasRowGroups ? multiplyRowGroup(i, b, x) : multiplyRow(i, b[i], x); + } + } + return convergence; + } + + template + typename IterationHelper::IterateResult IterationHelper::iterateUpper(storm::solver::OptimizationDirection const& dir, std::vector& x, std::vector const& b, bool takeMinOfOldAndNew) { + if (minimize(dir)) { + return iterateUpperInternal(x, b, takeMinOfOldAndNew); + } else { + return iterateUpperInternal(x, b, takeMinOfOldAndNew); + } + } + + template + typename IterationHelper::IterateResult IterationHelper::iterateUpper(std::vector& x, std::vector const& b, bool takeMinOfOldAndNew) { + return iterateUpperInternal(x, b, takeMinOfOldAndNew); + } + + template + template + typename IterationHelper::IterateResult IterationHelper::iterateUpperInternal(std::vector& x, std::vector const& b, bool takeMinOfOldAndNew) { // For each row compare the new upper bound candidate with the old one bool newUpperBoundAlwaysHigherEqual = true; bool newUpperBoundAlwaysLowerEqual = true; @@ -150,7 +199,7 @@ namespace storm { } template - ValueType UpperBoundIterator::multiplyRow(IndexType const& rowIndex, ValueType const& bi, std::vector const& x) { + ValueType IterationHelper::multiplyRow(IndexType const& rowIndex, ValueType const& bi, std::vector const& x) { assert(rowIndex < rowIndications.size()); ValueType xRes = bi; @@ -165,7 +214,7 @@ namespace storm { template template - ValueType UpperBoundIterator::multiplyRowGroup(IndexType const& rowGroupIndex, std::vector const& b, std::vector const& x) { + ValueType IterationHelper::multiplyRowGroup(IndexType const& rowGroupIndex, std::vector const& b, std::vector const& x) { STORM_LOG_ASSERT(rowGroupIndices != nullptr, "No row group indices available."); auto row = (*rowGroupIndices)[rowGroupIndex]; auto const& groupEnd = (*rowGroupIndices)[rowGroupIndex + 1]; @@ -180,24 +229,21 @@ namespace storm { } template - OptimisticValueIterationHelper::OptimisticValueIterationHelper(storm::storage::SparseMatrix const& matrix) : upperBoundIterator(matrix) { + OptimisticValueIterationHelper::OptimisticValueIterationHelper(storm::storage::SparseMatrix const& matrix) : iterationHelper(matrix) { // Intentionally left empty. - } + } template - std::pair OptimisticValueIterationHelper::solveEquationsOptimisticValueIteration(Environment const& env, std::vector* lowerX, std::vector* upperX, std::vector* auxVector, std::vector const& b, ValueIterationCallBackType const& valueIterationCallback, SingleIterationCallBackType const& singleIterationCallback, bool relative, ValueType precision, uint64_t maxOverallIterations, boost::optional dir, boost::optional relevantValues) { + std::pair OptimisticValueIterationHelper::solveEquations(Environment const& env, std::vector* lowerX, std::vector* upperX, std::vector const& b, bool relative, ValueType precision, uint64_t maxOverallIterations, boost::optional dir, boost::optional const& relevantValues) { STORM_LOG_ASSERT(lowerX->size() == upperX->size(), "Dimension missmatch."); - STORM_LOG_ASSERT(lowerX->size() == auxVector->size(), "Dimension missmatch."); // As we will shuffle pointers around, let's store the original positions here. std::vector* initLowerX = lowerX; std::vector* initUpperX = upperX; - std::vector* initAux = auxVector; uint64_t overallIterations = 0; uint64_t lastValueIterationIterations = 0; uint64_t currentVerificationIterations = 0; - uint64_t valueIterationInvocations = 0; // Get some parameters for the algorithm // 2 @@ -215,135 +261,101 @@ namespace storm { SolverStatus status = SolverStatus::InProgress; while (status == SolverStatus::InProgress && overallIterations < maxOverallIterations) { - // Perform value iteration until convergence - ++valueIterationInvocations; - auto result = valueIterationCallback(lowerX, auxVector, iterationPrecision, relative, overallIterations, maxOverallIterations); - lastValueIterationIterations = result.first; - overallIterations += result.first; + lastValueIterationIterations = dir ? iterationHelper.repeatedIterate(dir.get(), *lowerX, b, iterationPrecision, relative) : iterationHelper.repeatedIterate(*lowerX, b, iterationPrecision, relative); + overallIterations += lastValueIterationIterations; - if (result.second != SolverStatus::Converged) { - status = result.second; - } else { - bool intervalIterationNeeded = false; - currentVerificationIterations = 0; + bool intervalIterationNeeded = false; + currentVerificationIterations = 0; - if (relative) { - oviinternal::guessUpperBoundRelative(*lowerX, *upperX, relativeBoundGuessingScaler); - } else { - oviinternal::guessUpperBoundAbsolute(*lowerX, *upperX, precision); - } + if (relative) { + oviinternal::guessUpperBoundRelative(*lowerX, *upperX, relativeBoundGuessingScaler); + } else { + oviinternal::guessUpperBoundAbsolute(*lowerX, *upperX, precision); + } - bool cancelGuess = false; - while (status == SolverStatus::InProgress && overallIterations < maxOverallIterations) { - ++overallIterations; - ++currentVerificationIterations; - // Perform value iteration stepwise for lower bound and guessed upper bound + bool cancelGuess = false; + while (status == SolverStatus::InProgress && overallIterations < maxOverallIterations) { + ++overallIterations; + ++currentVerificationIterations; + // Perform value iteration stepwise for lower bound and guessed upper bound - // Upper bound iteration - auto upperBoundIterResult = dir ? upperBoundIterator.iterate(dir.get(), *upperX, b, !noTerminationGuarantee) : upperBoundIterator.iterate(*upperX, b, !noTerminationGuarantee); + // Upper bound iteration + auto upperBoundIterResult = dir ? iterationHelper.iterateUpper(dir.get(), *upperX, b, !noTerminationGuarantee) : iterationHelper.iterateUpper(*upperX, b, !noTerminationGuarantee); - if (upperBoundIterResult == oviinternal::UpperBoundIterator::IterateResult::AlwaysHigherOrEqual) { - // All values moved up (and did not stay the same) - // That means the guess for an upper bound is actually a lower bound - iterationPrecision = oviinternal::updateIterationPrecision(env, *auxVector, *upperX, relative, relevantValues); - // We assume to have a single fixed point. We can thus safely set the new lower bound, to the wrongly guessed upper bound - // Set lowerX to the upper bound candidate - std::swap(lowerX, upperX); - break; - } else if (upperBoundIterResult == oviinternal::UpperBoundIterator::IterateResult::AlwaysLowerOrEqual) { - // All values moved down (and stayed not the same) - // This is a valid upper bound. We still need to check the precision. - // We can safely use twice the requested precision, as we calculate the center of both vectors - bool reachedPrecision; - if (relevantValues) { - reachedPrecision = storm::utility::vector::equalModuloPrecision(*lowerX, *upperX, relevantValues.get(), doublePrecision, relative); - } else { - reachedPrecision = storm::utility::vector::equalModuloPrecision(*lowerX, *upperX, doublePrecision, relative); - } - if (reachedPrecision) { - status = SolverStatus::Converged; - break; - } else { - // From now on, we keep updating both bounds - intervalIterationNeeded = true; - } - // The following case below covers that both vectors (old and new) are equal. - // Theoretically, this means that the precise fixpoint has been reached. However, numerical instabilities can be tricky and this detection might be incorrect (see the haddad-monmege model). - // We therefore disable it. It is very unlikely that we guessed the right fixpoint anyway. - //} else if (upperBoundIterResult == oviinternal::UpperBoundIterator::IterateResult::Equal) { - // In this case, the guessed upper bound is the precise fixpoint - // status = SolverStatus::Converged; - // std::swap(lowerX, auxVector); - // break; + if (upperBoundIterResult == oviinternal::IterationHelper::IterateResult::AlwaysHigherOrEqual) { + // All values moved up (and did not stay the same) + // That means the guess for an upper bound is actually a lower bound + auto diff = dir ? iterationHelper.singleIterationWithDiff(dir.get(), *upperX, b, relative) : iterationHelper.singleIterationWithDiff(*upperX, b, relative); + iterationPrecision = oviinternal::updateIterationPrecision(env, diff); + // We assume to have a single fixed point. We can thus safely set the new lower bound, to the wrongly guessed upper bound + // Set lowerX to the upper bound candidate + std::swap(lowerX, upperX); + break; + } else if (upperBoundIterResult == oviinternal::IterationHelper::IterateResult::AlwaysLowerOrEqual) { + // All values moved down (and stayed not the same) + // This is a valid upper bound. We still need to check the precision. + // We can safely use twice the requested precision, as we calculate the center of both vectors + bool reachedPrecision; + if (relevantValues) { + reachedPrecision = storm::utility::vector::equalModuloPrecision(*lowerX, *upperX, relevantValues.get(), doublePrecision, relative); + } else { + reachedPrecision = storm::utility::vector::equalModuloPrecision(*lowerX, *upperX, doublePrecision, relative); } - - // Check whether we tried this guess for too long - ValueType scaledIterationCount = storm::utility::convertNumber(currentVerificationIterations) * storm::utility::convertNumber(env.solver().ovi().getMaxVerificationIterationFactor()); - if (!intervalIterationNeeded && scaledIterationCount >= storm::utility::convertNumber(lastValueIterationIterations)) { - cancelGuess = true; - // In this case we will make one more iteration on the lower bound (mainly to obtain a new iterationPrecision) + if (reachedPrecision) { + status = SolverStatus::Converged; + break; + } else { + // From now on, we keep updating both bounds + intervalIterationNeeded = true; } + // The following case below covers that both vectors (old and new) are equal. + // Theoretically, this means that the precise fixpoint has been reached. However, numerical instabilities can be tricky and this detection might be incorrect (see the haddad-monmege model). + // We therefore disable it. It is very unlikely that we guessed the right fixpoint anyway. + //} else if (upperBoundIterResult == oviinternal::IterationHelper::IterateResult::Equal) { + // In this case, the guessed upper bound is the precise fixpoint + // status = SolverStatus::Converged; + // break; + } - // Lower bound iteration (only if needed) - if (cancelGuess || intervalIterationNeeded || currentVerificationIterations > upperBoundOnlyIterations) { - singleIterationCallback(lowerX, auxVector, overallIterations); - // At this point, auxVector contains the old values for the lower bound whereas lowerX contains the new ones. + // Check whether we tried this guess for too long + ValueType scaledIterationCount = storm::utility::convertNumber(currentVerificationIterations) * storm::utility::convertNumber(env.solver().ovi().getMaxVerificationIterationFactor()); + if (!intervalIterationNeeded && scaledIterationCount * iterationPrecision >= storm::utility::one()) { + cancelGuess = true; + // In this case we will make one more iteration on the lower bound (mainly to obtain a new iterationPrecision) + } - // Check whether the upper and lower bounds have crossed, i.e., the upper bound is smaller than the lower bound. - bool valuesCrossed = false; - for (uint64_t i = 0; i < lowerX->size(); ++i) { - if ((*upperX)[i] < (*lowerX)[i]) { - valuesCrossed = true; - break; - } - } + // Lower bound iteration (only if needed) + if (cancelGuess || intervalIterationNeeded || currentVerificationIterations > upperBoundOnlyIterations) { + auto diff = dir ? iterationHelper.singleIterationWithDiff(dir.get(), *lowerX, b, relative) : iterationHelper.singleIterationWithDiff(*lowerX, b, relative); - if (cancelGuess || valuesCrossed) { - // A new guess is needed. - iterationPrecision = oviinternal::updateIterationPrecision(env, *auxVector, *lowerX, relative, relevantValues); + // Check whether the upper and lower bounds have crossed, i.e., the upper bound is smaller than the lower bound. + bool valuesCrossed = false; + for (uint64_t i = 0; i < lowerX->size(); ++i) { + if ((*upperX)[i] < (*lowerX)[i]) { + valuesCrossed = true; break; } } - } - if (storm::utility::resources::isTerminate()) { - status = SolverStatus::Aborted; + + if (cancelGuess || valuesCrossed) { + // A new guess is needed. + iterationPrecision = oviinternal::updateIterationPrecision(env, diff); + break; + } } } + if (storm::utility::resources::isTerminate()) { + status = SolverStatus::Aborted; + } } // end while - // Swap the results into the output vectors. - if (initLowerX == lowerX) { - // lowerX is already at the correct position. We still have to care for upperX - if (initUpperX != upperX) { - // UpperX is not at the correct position. It has to be at the auxVector - assert(initAux == upperX); - std::swap(*initUpperX, *initAux); - } - } else if (initUpperX == upperX) { - // UpperX is already at the correct position. - // We already know that lowerX is at the wrong position. It has to be at the auxVector - assert(initAux == lowerX); - std::swap(*initLowerX, *initAux); - } else if (initAux == auxVector) { - // We know that upperX and lowerX are swapped. - assert(initLowerX == upperX); + // Swap the results into the output vectors (if necessary). + assert(initLowerX != lowerX || (initLowerX == lowerX && initUpperX == upperX)); + if (initLowerX != lowerX) { assert(initUpperX == lowerX); - std::swap(*initUpperX, *initLowerX); - } else { - // Now we know that all vectors are at the wrong position. There are only two possibilities left - if (initLowerX == upperX) { - assert(initUpperX == auxVector); - assert(initAux == lowerX); - std::swap(*initLowerX, *initAux); - std::swap(*initUpperX, *initAux); - } else { - assert(initLowerX == auxVector); - assert(initUpperX == lowerX); - assert (initAux == upperX); - std::swap(*initUpperX, *initAux); - std::swap(*initLowerX, *initAux); - } + assert(initLowerX == upperX); + lowerX->swap(*upperX); } if (overallIterations > maxOverallIterations) { diff --git a/src/storm/solver/helper/OptimisticValueIterationHelper.h b/src/storm/solver/helper/OptimisticValueIterationHelper.h index abf67dd3c..d559d508d 100644 --- a/src/storm/solver/helper/OptimisticValueIterationHelper.h +++ b/src/storm/solver/helper/OptimisticValueIterationHelper.h @@ -18,10 +18,10 @@ namespace storm { namespace oviinternal { template - class UpperBoundIterator { + class IterationHelper { public: typedef uint32_t IndexType; - UpperBoundIterator(storm::storage::SparseMatrix const& matrix); + IterationHelper(storm::storage::SparseMatrix const& matrix); enum class IterateResult { AlwaysHigherOrEqual, AlwaysLowerOrEqual, @@ -29,13 +29,26 @@ namespace storm { Incomparable }; - IterateResult iterate(std::vector& x, std::vector const& b, bool takeMinOfOldAndNew); - IterateResult iterate(storm::solver::OptimizationDirection const& dir, std::vector& x, std::vector const& b, bool takeMinOfOldAndNew); + /// performs a single iteration and returns the maximal difference between the iterations as well as the index where this difference happened + ValueType singleIterationWithDiff(std::vector& x, std::vector const& b, bool computeRelativeDiff); + ValueType singleIterationWithDiff(storm::solver::OptimizationDirection const& dir, std::vector& x, std::vector const& b, bool computeRelativeDiff); + + /// returns the number of performed iterations + uint64_t repeatedIterate(std::vector& x, std::vector const& b, ValueType precision, bool relative); + uint64_t repeatedIterate(storm::solver::OptimizationDirection const& dir, std::vector& x, std::vector const& b, ValueType precision, bool relative); + + /// Performs a single iteration for the upper bound + IterateResult iterateUpper(std::vector& x, std::vector const& b, bool takeMinOfOldAndNew); + IterateResult iterateUpper(storm::solver::OptimizationDirection const& dir, std::vector& x, std::vector const& b, bool takeMinOfOldAndNew); private: template - IterateResult iterateInternal(std::vector& x, std::vector const& b, bool takeMinOfOldAndNew); + ValueType singleIterationWithDiffInternal(std::vector& x, std::vector const& b, bool computeRelativeDiff); + template + uint64_t repeatedIterateInternal(std::vector& x, std::vector const& b, ValueType precision, bool relative); + template + IterateResult iterateUpperInternal(std::vector& x, std::vector const& b, bool takeMinOfOldAndNew); ValueType multiplyRow(IndexType const& rowIndex, ValueType const& bi, std::vector const& x); template ValueType multiplyRowGroup(IndexType const& rowGroupIndex, std::vector const& b, std::vector const& x); @@ -56,49 +69,20 @@ namespace storm { public: OptimisticValueIterationHelper(storm::storage::SparseMatrix const& matrix); - /*! - * Return type of value iteration callback. - * The first component shall be the number of performed iterations, the second component is the final convergence status. - */ - typedef std::pair ValueIterationReturnType; - - /*! - * Value iteration callback. Performs conventional value iteration for the given input. - * @pre y points to a vector that contains starting values - * @post y points to a vector that contains resulting values - * @param yPrime points to auxiliary storage - * @param precision is the target precision - * @param relative sets whether relative precision is considered - * @param maxI the maximal number of iterations to perform - */ - typedef std::function*& y, std::vector*& yPrime, ValueType const& precision, bool const& relative, uint64_t const& i, uint64_t const& maxI)> ValueIterationCallBackType; - - /*! - * Should perform a single value iteration step (using conventional value iteration). - * @pre y points to a vector that contains starting values - * @post y points to a vector that contains resulting values - * @param yPrime points to auxiliary storage - * @param currI the current iteration (can be used to display progress within the callback) - */ - typedef std::function* y, std::vector* yPrime, uint64_t const& currI)> SingleIterationCallBackType; - /*! * @param env * @param lowerX Needs to be some arbitrary lower bound on the actual values initially * @param upperX Does not need to be an upper bound initially - * @param auxVector auxiliary storage (same size as lowerX and upperX) * @param b the values added to each matrix row (the b in A*x+b) - * @param valueIterationCallback Function that should perform standard value iteration on the input vector - * @param singleIterationCallback Function that should perform a single value iteration step on the input vector e.g. ( x' = min/max(A*x + b)) * @param dir The optimization direction * @param relevantValues If given, we only check the precision at the states with the given indices. * @return The status upon termination as well as the number of iterations Also, the maximum (relative/absolute) difference between lowerX and upperX will be 2*epsilon * with the provided precision parameters. */ - std::pair solveEquationsOptimisticValueIteration(Environment const& env, std::vector* lowerX, std::vector* upperX, std::vector* auxVector, std::vector const& b, ValueIterationCallBackType const& valueIterationCallback, SingleIterationCallBackType const& singleIterationCallback, bool relative, ValueType precision, uint64_t maxOverallIterations, boost::optional dir, boost::optional relevantValues = boost::none); + std::pair solveEquations(Environment const& env, std::vector* lowerX, std::vector* upperX, std::vector const& b, bool relative, ValueType precision, uint64_t maxOverallIterations, boost::optional dir, boost::optional const& relevantValues); private: - oviinternal::UpperBoundIterator upperBoundIterator; + oviinternal::IterationHelper iterationHelper; }; } diff --git a/src/storm/storage/jani/EdgeContainer.cpp b/src/storm/storage/jani/EdgeContainer.cpp index f2e3b9855..65230dcbb 100644 --- a/src/storm/storage/jani/EdgeContainer.cpp +++ b/src/storm/storage/jani/EdgeContainer.cpp @@ -67,9 +67,13 @@ namespace storm { e.setTemplateEdge(map[e.getTemplateEdge()]); } } - - - + } + + EdgeContainer& EdgeContainer::operator=(EdgeContainer const& other) { + EdgeContainer otherCpy(other); + this->templates = std::move(otherCpy.templates); + this->edges = std::move(otherCpy.edges); + return *this; } void EdgeContainer::finalize(Model const& containingModel) { diff --git a/src/storm/storage/jani/EdgeContainer.h b/src/storm/storage/jani/EdgeContainer.h index bf1e0cbe9..dc2cb8b69 100644 --- a/src/storm/storage/jani/EdgeContainer.h +++ b/src/storm/storage/jani/EdgeContainer.h @@ -91,6 +91,7 @@ namespace storm { EdgeContainer() = default; EdgeContainer(EdgeContainer const& other); + EdgeContainer& operator=(EdgeContainer const& other); void clearConcreteEdges(); std::vector const& getConcreteEdges() const; diff --git a/src/storm/storage/jani/TemplateEdgeContainer.cpp b/src/storm/storage/jani/TemplateEdgeContainer.cpp index 81ce8f874..66d93a4fb 100644 --- a/src/storm/storage/jani/TemplateEdgeContainer.cpp +++ b/src/storm/storage/jani/TemplateEdgeContainer.cpp @@ -3,10 +3,18 @@ namespace storm { namespace jani { - TemplateEdgeContainer::TemplateEdgeContainer(TemplateEdgeContainer const &other) { + TemplateEdgeContainer::TemplateEdgeContainer(TemplateEdgeContainer const &other) : std::unordered_set>() { for (auto const& te : other) { this->insert(std::make_shared(*te)); } } + + TemplateEdgeContainer& TemplateEdgeContainer::operator=(const TemplateEdgeContainer& other) { + this->clear(); + for (auto const& te : other) { + this->insert(std::make_shared(*te)); + } + return *this; + } } -} \ No newline at end of file +} diff --git a/src/storm/storage/jani/TemplateEdgeContainer.h b/src/storm/storage/jani/TemplateEdgeContainer.h index 293302594..12dd6050f 100644 --- a/src/storm/storage/jani/TemplateEdgeContainer.h +++ b/src/storm/storage/jani/TemplateEdgeContainer.h @@ -11,6 +11,7 @@ namespace storm { struct TemplateEdgeContainer : public std::unordered_set> { TemplateEdgeContainer() = default; TemplateEdgeContainer(TemplateEdgeContainer const& other); + TemplateEdgeContainer& operator=(TemplateEdgeContainer const& other); }; } } \ No newline at end of file diff --git a/src/storm/utility/constants.cpp b/src/storm/utility/constants.cpp index 1b2bdff97..9172807ad 100644 --- a/src/storm/utility/constants.cpp +++ b/src/storm/utility/constants.cpp @@ -55,12 +55,14 @@ namespace storm { return std::isnan(value); } - bool isAlmostZero(double const& a) { - return a < 1e-12 && a > -1e-12; + template + bool isAlmostZero(ValueType const& a) { + return a < convertNumber(1e-12) && a > -convertNumber(1e-12); } - bool isAlmostOne(double const& a) { - return a < (1.0 + 1e-12) && a > (1.0 - 1e-12); + template + bool isAlmostOne(ValueType const& a) { + return a < convertNumber(1.0 + 1e-12) && a > convertNumber(1.0 - 1e-12); } template @@ -818,6 +820,16 @@ namespace storm { return std::move(value); } + template<> + bool isAlmostZero(storm::RationalFunction const& a) { + return a.isConstant() && isAlmostZero(convertNumber(a)); + } + + template<> + bool isAlmostOne(storm::RationalFunction const& a) { + return a.isConstant() && isAlmostOne(convertNumber(a)); + } + template<> std::pair minmax(std::vector const&) { STORM_LOG_THROW(false, storm::exceptions::InvalidArgumentException, "Minimum/maximum for rational functions is not defined."); @@ -885,6 +897,8 @@ namespace storm { template double infinity(); template bool isOne(double const& value); template bool isZero(double const& value); + template bool isAlmostZero(double const& value); + template bool isAlmostOne(double const& value); template bool isConstant(double const& value); template bool isInfinity(double const& value); template bool isInteger(double const& number); @@ -982,6 +996,8 @@ namespace storm { template bool isConstant(storm::ClnRationalNumber const& value); template bool isInfinity(storm::ClnRationalNumber const& value); template bool isNan(storm::ClnRationalNumber const& value); + template bool isAlmostZero(storm::ClnRationalNumber const& value); + template bool isAlmostOne(storm::ClnRationalNumber const& value); template storm::NumberTraits::IntegerType convertNumber(storm::NumberTraits::IntegerType const& number); template storm::ClnRationalNumber convertNumber(storm::ClnRationalNumber const& number); template storm::ClnRationalNumber simplify(storm::ClnRationalNumber value); @@ -1008,6 +1024,8 @@ namespace storm { template bool isConstant(storm::GmpRationalNumber const& value); template bool isInfinity(storm::GmpRationalNumber const& value); template bool isNan(storm::GmpRationalNumber const& value); + template bool isAlmostZero(storm::GmpRationalNumber const& value); + template bool isAlmostOne(storm::GmpRationalNumber const& value); template storm::NumberTraits::IntegerType convertNumber(storm::NumberTraits::IntegerType const& number); template storm::GmpRationalNumber convertNumber(storm::GmpRationalNumber const& number); template storm::GmpRationalNumber simplify(storm::GmpRationalNumber value); diff --git a/src/storm/utility/constants.h b/src/storm/utility/constants.h index ddd125edb..787ead089 100644 --- a/src/storm/utility/constants.h +++ b/src/storm/utility/constants.h @@ -88,9 +88,11 @@ namespace storm { template bool isNan(ValueType const& a); - bool isAlmostZero(double const& a); + template + bool isAlmostZero(ValueType const& a); - bool isAlmostOne(double const& a); + template + bool isAlmostOne(ValueType const& a); template bool isConstant(ValueType const& a); diff --git a/src/storm/utility/shortestPaths.h b/src/storm/utility/shortestPaths.h index 60b206bca..ab3396a3f 100644 --- a/src/storm/utility/shortestPaths.h +++ b/src/storm/utility/shortestPaths.h @@ -47,7 +47,7 @@ namespace storm { if (predecessorNode != rhs.predecessorNode) { return predecessorNode < rhs.predecessorNode; } - return predecessorK < predecessorK; + return predecessorK < rhs.predecessorK; } bool operator==(const Path& rhs) const { diff --git a/version.cmake b/version.cmake index 801d8403c..4fcf2a030 100644 --- a/version.cmake +++ b/version.cmake @@ -1,4 +1,4 @@ set(STORM_VERSION_MAJOR 1) set(STORM_VERSION_MINOR 6) -set(STORM_VERSION_PATCH 0) +set(STORM_VERSION_PATCH 2)