diff --git a/src/storm/modelchecker/prctl/helper/SparseMdpPrctlHelper.cpp b/src/storm/modelchecker/prctl/helper/SparseMdpPrctlHelper.cpp index 17bef0ff8..9193c8f26 100644 --- a/src/storm/modelchecker/prctl/helper/SparseMdpPrctlHelper.cpp +++ b/src/storm/modelchecker/prctl/helper/SparseMdpPrctlHelper.cpp @@ -52,19 +52,19 @@ namespace storm { namespace modelchecker { namespace helper { - + template std::map SparseMdpPrctlHelper::computeRewardBoundedValues(Environment const& env, OptimizationDirection dir, rewardbounded::MultiDimensionalRewardUnfolding& rewardUnfolding, storm::storage::BitVector const& initialStates) { storm::utility::Stopwatch swAll(true), swBuild, swCheck; - + // Get lower and upper bounds for the solution. auto lowerBound = rewardUnfolding.getLowerObjectiveBound(); auto upperBound = rewardUnfolding.getUpperObjectiveBound(); - + // Initialize epoch models auto initEpoch = rewardUnfolding.getStartEpoch(); auto epochOrder = rewardUnfolding.getEpochComputationOrder(initEpoch); - + // initialize data that will be needed for each epoch std::vector x, b; std::unique_ptr> minMaxSolver; @@ -72,7 +72,7 @@ namespace storm { ValueType precision = rewardUnfolding.getRequiredEpochModelPrecision(initEpoch, storm::utility::convertNumber(storm::settings::getModule().getPrecision())); Environment preciseEnv = env; preciseEnv.solver().minMax().setPrecision(storm::utility::convertNumber(precision)); - + // In case of cdf export we store the necessary data. std::vector> cdfData; @@ -101,14 +101,14 @@ namespace storm { break; } } - + std::map result; for (auto initState : initialStates) { result[initState] = rewardUnfolding.getInitialStateResult(initEpoch, initState); } - + swAll.stop(); - + if (storm::settings::getModule().isExportCdfSet()) { std::vector headers; for (uint64_t i = 0; i < rewardUnfolding.getEpochManager().getDimensionCount(); ++i) { @@ -118,7 +118,7 @@ namespace storm { storm::utility::exportDataToCSVFile(storm::settings::getModule().getExportCdfDirectory() + "cdf.csv", cdfData, headers); } - + if (storm::settings::getModule().isShowStatisticsSet()) { STORM_PRINT_AND_LOG("---------------------------------" << std::endl); STORM_PRINT_AND_LOG("Statistics:" << std::endl); @@ -129,7 +129,7 @@ namespace storm { STORM_PRINT_AND_LOG("Epoch Model checking Time: " << swCheck << "." << std::endl); STORM_PRINT_AND_LOG("---------------------------------" << std::endl); } - + return result; } @@ -154,7 +154,7 @@ namespace storm { return MDPSparseModelCheckingHelperReturnType(std::move(result), std::move(allStates), nullptr, std::move(choiceValues)); } - + template std::vector computeValidSchedulerHint(Environment const& env, SolutionType const& type, storm::storage::SparseMatrix const& transitionMatrix, storm::storage::SparseMatrix const& backwardTransitions, storm::storage::BitVector const& maybeStates, storm::storage::BitVector const& filterStates, storm::storage::BitVector const& targetStates) { storm::storage::Scheduler validScheduler(maybeStates.size()); @@ -166,7 +166,7 @@ namespace storm { } else { STORM_LOG_ASSERT(false, "Unexpected equation system type."); } - + // Extract the relevant parts of the scheduler for the solver. std::vector schedulerHint(maybeStates.getNumberOfSetBits()); auto maybeIt = maybeStates.begin(); @@ -176,13 +176,13 @@ namespace storm { } return schedulerHint; } - + template struct SparseMdpHintType { SparseMdpHintType() : eliminateEndComponents(false), computeUpperBounds(false), uniqueSolution(false), noEndComponents(false) { // Intentionally left empty. } - + bool hasSchedulerHint() const { return static_cast(schedulerHint); } @@ -198,11 +198,11 @@ namespace storm { ValueType const& getLowerResultBound() const { return lowerResultBound.get(); } - + bool hasUpperResultBound() const { return static_cast(upperResultBound); } - + bool hasUpperResultBounds() const { return static_cast(upperResultBounds); } @@ -214,19 +214,19 @@ namespace storm { std::vector& getUpperResultBounds() { return upperResultBounds.get(); } - + std::vector const& getUpperResultBounds() const { return upperResultBounds.get(); } - + std::vector& getSchedulerHint() { return schedulerHint.get(); } - + std::vector& getValueHint() { return valueHint.get(); } - + bool getEliminateEndComponents() const { return eliminateEndComponents; } @@ -238,11 +238,11 @@ namespace storm { bool hasUniqueSolution() const { return uniqueSolution; } - + bool hasNoEndComponents() const { return noEndComponents; } - + boost::optional> schedulerHint; boost::optional> valueHint; boost::optional lowerResultBound; @@ -253,10 +253,10 @@ namespace storm { bool uniqueSolution; bool noEndComponents; }; - + template void extractValueAndSchedulerHint(SparseMdpHintType& hintStorage, storm::storage::SparseMatrix const& transitionMatrix, storm::storage::SparseMatrix const& backwardTransitions, storm::storage::BitVector const& maybeStates, boost::optional const& selectedChoices, ModelCheckerHint const& hint, bool skipECWithinMaybeStatesCheck) { - + // Deal with scheduler hint. if (hint.isExplicitModelCheckerHint() && hint.template asExplicitModelCheckerHint().hasSchedulerHint()) { if (hintStorage.hasSchedulerHint()) { @@ -276,7 +276,7 @@ namespace storm { } else { hintApplicable = true; } - + if (hintApplicable) { // Compute the hint w.r.t. the given subsystem. hintChoices.clear(); @@ -297,13 +297,13 @@ namespace storm { } } } - + // Deal with solution value hint. Only applicable if there are no End Components consisting of maybe states. if (hint.isExplicitModelCheckerHint() && hint.template asExplicitModelCheckerHint().hasResultHint() && (skipECWithinMaybeStatesCheck || hintStorage.hasSchedulerHint() || storm::utility::graph::performProb1A(transitionMatrix, transitionMatrix.getRowGroupIndices(), backwardTransitions, maybeStates, ~maybeStates).full())) { hintStorage.valueHint = storm::utility::vector::filterVector(hint.template asExplicitModelCheckerHint().getResultHint(), maybeStates); } } - + template SparseMdpHintType computeHints(Environment const& env, SolutionType const& type, ModelCheckerHint const& hint, storm::OptimizationDirection const& dir, storm::storage::SparseMatrix const& transitionMatrix, storm::storage::SparseMatrix const& backwardTransitions, storm::storage::BitVector const& maybeStates, storm::storage::BitVector const& phiStates, storm::storage::BitVector const& targetStates, bool produceScheduler, boost::optional const& selectedChoices = boost::none) { SparseMdpHintType result; @@ -313,11 +313,11 @@ namespace storm { result.noEndComponents = (dir == storm::solver::OptimizationDirection::Minimize && type == SolutionType::UntilProbabilities) || (dir == storm::solver::OptimizationDirection::Maximize && type == SolutionType::ExpectedRewards) || (hint.isExplicitModelCheckerHint() && hint.asExplicitModelCheckerHint().getNoEndComponentsInMaybeStates()); - + // If there are no end components, the solution is unique. (Note that the other direction does not hold, // e.g., end components in which infinite reward is collected. result.uniqueSolution = result.hasNoEndComponents(); - + // Check for requirements of the solver. bool hasSchedulerHint = hint.isExplicitModelCheckerHint() && hint.template asExplicitModelCheckerHint().hasSchedulerHint(); storm::solver::GeneralMinMaxLinearEquationSolverFactory minMaxLinearEquationSolverFactory; @@ -335,14 +335,14 @@ namespace storm { // Note that in the case of minimizing expected rewards there might still be end components in which reward is collected. result.noEndComponents = (type == SolutionType::UntilProbabilities); } - + // If the solver requires an initial scheduler, compute one now. Note that any scheduler is valid if there are no end components. if (requirements.validInitialScheduler() && !result.noEndComponents) { STORM_LOG_DEBUG("Computing valid scheduler, because the solver requires it."); result.schedulerHint = computeValidSchedulerHint(env, type, transitionMatrix, backwardTransitions, maybeStates, phiStates, targetStates); requirements.clearValidInitialScheduler(); } - + // Finally, we have information on the bounds depending on the problem type. if (type == SolutionType::UntilProbabilities) { requirements.clearBounds(); @@ -373,7 +373,7 @@ namespace storm { if (!result.hasUpperResultBound() && type == SolutionType::UntilProbabilities) { result.upperResultBound = storm::utility::one(); } - + // If we received an upper bound, we can drop the requirement to compute one. if (result.hasUpperResultBound()) { result.computeUpperBounds = false; @@ -381,35 +381,35 @@ namespace storm { return result; } - + template struct MaybeStateResult { MaybeStateResult(std::vector&& values) : values(std::move(values)) { // Intentionally left empty. } - + bool hasScheduler() const { return static_cast(scheduler); } - + std::vector const& getScheduler() const { return scheduler.get(); } - + std::vector const& getValues() const { return values; } - + std::vector values; boost::optional> scheduler; }; - + template MaybeStateResult computeValuesForMaybeStates(Environment const& env, storm::solver::SolveGoal&& goal, storm::storage::SparseMatrix&& submatrix, std::vector const& b, bool produceScheduler, SparseMdpHintType& hint) { - + // Initialize the solution vector. std::vector x = hint.hasValueHint() ? std::move(hint.getValueHint()) : std::vector(submatrix.getRowGroupCount(), hint.hasLowerResultBound() ? hint.getLowerResultBound() : storm::utility::zero()); - + // Set up the solver. storm::solver::GeneralMinMaxLinearEquationSolverFactory minMaxLinearEquationSolverFactory; std::unique_ptr> solver = storm::solver::configureMinMaxLinearEquationSolver(env, std::move(goal), minMaxLinearEquationSolverFactory, std::move(submatrix)); @@ -429,10 +429,10 @@ namespace storm { solver->setInitialScheduler(std::move(hint.getSchedulerHint())); } solver->setTrackScheduler(produceScheduler); - + // Solve the corresponding system of equations. solver->solveEquations(env, x, b); - + #ifndef NDEBUG // As a sanity check, make sure our local upper bounds were in fact correct. if (solver->hasUpperBound(storm::solver::AbstractEquationSolver::BoundType::Local)) { @@ -443,7 +443,7 @@ namespace storm { } } #endif - + // Create result. MaybeStateResult result(std::move(x)); @@ -453,18 +453,18 @@ namespace storm { } return result; } - + struct QualitativeStateSetsUntilProbabilities { storm::storage::BitVector maybeStates; storm::storage::BitVector statesWithProbability0; storm::storage::BitVector statesWithProbability1; }; - + template QualitativeStateSetsUntilProbabilities getQualitativeStateSetsUntilProbabilitiesFromHint(ModelCheckerHint const& hint) { QualitativeStateSetsUntilProbabilities result; result.maybeStates = hint.template asExplicitModelCheckerHint().getMaybeStates(); - + // Treat the states with probability zero/one. std::vector const& resultsForNonMaybeStates = hint.template asExplicitModelCheckerHint().getResultHint(); result.statesWithProbability1 = storm::storage::BitVector(result.maybeStates.size()); @@ -478,10 +478,10 @@ namespace storm { result.statesWithProbability0.set(state, true); } } - + return result; } - + template QualitativeStateSetsUntilProbabilities computeQualitativeStateSetsUntilProbabilities(storm::solver::SolveGoal const& goal, storm::storage::SparseMatrix const& transitionMatrix, storm::storage::SparseMatrix const& backwardTransitions, storm::storage::BitVector const& phiStates, storm::storage::BitVector const& psiStates) { QualitativeStateSetsUntilProbabilities result; @@ -496,10 +496,10 @@ namespace storm { result.statesWithProbability0 = std::move(statesWithProbability01.first); result.statesWithProbability1 = std::move(statesWithProbability01.second); result.maybeStates = ~(result.statesWithProbability0 | result.statesWithProbability1); - + return result; } - + template QualitativeStateSetsUntilProbabilities getQualitativeStateSetsUntilProbabilities(storm::solver::SolveGoal const& goal, storm::storage::SparseMatrix const& transitionMatrix, storm::storage::SparseMatrix const& backwardTransitions, storm::storage::BitVector const& phiStates, storm::storage::BitVector const& psiStates, ModelCheckerHint const& hint) { if (hint.isExplicitModelCheckerHint() && hint.template asExplicitModelCheckerHint().getComputeOnlyMaybeStates()) { @@ -508,7 +508,7 @@ namespace storm { return computeQualitativeStateSetsUntilProbabilities(goal, transitionMatrix, backwardTransitions, phiStates, psiStates); } } - + template void extractSchedulerChoices(storm::storage::Scheduler& scheduler, std::vector const& subChoices, storm::storage::BitVector const& maybeStates) { auto subChoiceIt = subChoices.begin(); @@ -518,10 +518,10 @@ namespace storm { } assert(subChoiceIt == subChoices.end()); } - + template void extendScheduler(storm::storage::Scheduler& scheduler, storm::solver::SolveGoal const& goal, QualitativeStateSetsUntilProbabilities const& qualitativeStateSets, storm::storage::SparseMatrix const& transitionMatrix, storm::storage::SparseMatrix const& backwardTransitions, storm::storage::BitVector const& phiStates, storm::storage::BitVector const& psiStates) { - + // Finally, if we need to produce a scheduler, we also need to figure out the parts of the scheduler for // the states with probability 1 or 0 (depending on whether we maximize or minimize). // We also need to define some arbitrary choice for the remaining states to obtain a fully defined scheduler. @@ -537,40 +537,40 @@ namespace storm { } } } - + template void computeFixedPointSystemUntilProbabilities(storm::solver::SolveGoal& goal, storm::storage::SparseMatrix const& transitionMatrix, QualitativeStateSetsUntilProbabilities const& qualitativeStateSets, storm::storage::SparseMatrix& submatrix, std::vector& b) { // First, we can eliminate the rows and columns from the original transition probability matrix for states // whose probabilities are already known. submatrix = transitionMatrix.getSubmatrix(true, qualitativeStateSets.maybeStates, qualitativeStateSets.maybeStates, false); - + // Prepare the right-hand side of the equation system. For entry i this corresponds to // the accumulated probability of going from state i to some state that has probability 1. b = transitionMatrix.getConstrainedRowGroupSumVector(qualitativeStateSets.maybeStates, qualitativeStateSets.statesWithProbability1); - + // If the solve goal has relevant values, we need to adjust them. goal.restrictRelevantValues(qualitativeStateSets.maybeStates); } - + template boost::optional> computeFixedPointSystemUntilProbabilitiesEliminateEndComponents(storm::solver::SolveGoal& goal, storm::storage::SparseMatrix const& transitionMatrix, storm::storage::SparseMatrix const& backwardTransitions, QualitativeStateSetsUntilProbabilities const& qualitativeStateSets, storm::storage::SparseMatrix& submatrix, std::vector& b, bool produceScheduler) { - + // Get the set of states that (under some scheduler) can stay in the set of maybestates forever storm::storage::BitVector candidateStates = storm::utility::graph::performProb0E(transitionMatrix, transitionMatrix.getRowGroupIndices(), backwardTransitions, qualitativeStateSets.maybeStates, ~qualitativeStateSets.maybeStates); - + bool doDecomposition = !candidateStates.empty(); - + storm::storage::MaximalEndComponentDecomposition endComponentDecomposition; if (doDecomposition) { // Compute the states that are in MECs. endComponentDecomposition = storm::storage::MaximalEndComponentDecomposition(transitionMatrix, backwardTransitions, candidateStates); } - + // Only do more work if there are actually end-components. if (doDecomposition && !endComponentDecomposition.empty()) { STORM_LOG_DEBUG("Eliminating " << endComponentDecomposition.size() << " EC(s)."); SparseMdpEndComponentInformation result = SparseMdpEndComponentInformation::eliminateEndComponents(endComponentDecomposition, transitionMatrix, qualitativeStateSets.maybeStates, &qualitativeStateSets.statesWithProbability1, nullptr, nullptr, submatrix, &b, nullptr, produceScheduler); - + // If the solve goal has relevant values, we need to adjust them. if (goal.hasRelevantValues()) { storm::storage::BitVector newRelevantValues(submatrix.getRowGroupCount()); @@ -583,38 +583,38 @@ namespace storm { goal.setRelevantValues(std::move(newRelevantValues)); } } - + return result; } else { STORM_LOG_DEBUG("Not eliminating ECs as there are none."); computeFixedPointSystemUntilProbabilities(goal, transitionMatrix, qualitativeStateSets, submatrix, b); - + return boost::none; } } - + template MDPSparseModelCheckingHelperReturnType SparseMdpPrctlHelper::computeUntilProbabilities(Environment const& env, storm::solver::SolveGoal&& goal, storm::storage::SparseMatrix const& transitionMatrix, storm::storage::SparseMatrix const& backwardTransitions, storm::storage::BitVector const& phiStates, storm::storage::BitVector const& psiStates, bool qualitative, bool produceScheduler, ModelCheckerHint const& hint) { STORM_LOG_THROW(!qualitative || !produceScheduler, storm::exceptions::InvalidSettingsException, "Cannot produce scheduler when performing qualitative model checking only."); - + // Prepare resulting vector. std::vector result(transitionMatrix.getRowGroupCount(), storm::utility::zero()); - + // We need to identify the maybe states (states which have a probability for satisfying the until formula // that is strictly between 0 and 1) and the states that satisfy the formula with probablity 1 and 0, respectively. QualitativeStateSetsUntilProbabilities qualitativeStateSets = getQualitativeStateSetsUntilProbabilities(goal, transitionMatrix, backwardTransitions, phiStates, psiStates, hint); STORM_LOG_INFO("Preprocessing: " << qualitativeStateSets.statesWithProbability1.getNumberOfSetBits() << " states with probability 1, " << qualitativeStateSets.statesWithProbability0.getNumberOfSetBits() << " with probability 0 (" << qualitativeStateSets.maybeStates.getNumberOfSetBits() << " states remaining)."); - + // Set values of resulting vector that are known exactly. storm::utility::vector::setVectorValues(result, qualitativeStateSets.statesWithProbability1, storm::utility::one()); - + // If requested, we will produce a scheduler. std::unique_ptr> scheduler; if (produceScheduler) { scheduler = std::make_unique>(transitionMatrix.getRowGroupCount()); } - + // Check if the values of the maybe states are relevant for the SolveGoal bool maybeStatesNotRelevant = goal.hasRelevantValues() && goal.relevantValues().isDisjointFrom(qualitativeStateSets.maybeStates); @@ -629,7 +629,7 @@ namespace storm { } std::vector maybeStateChoiceValues = std::vector(sizeMaybeStateChoiceValues, storm::utility::zero()); - + // TODO: if a scheduler is to be produced and maybestatesNotRelevant is true, we have to set the scheduler for maybsetsates as "unreachable" TODO // Check whether we need to compute exact probabilities for some states. if ((qualitative || maybeStatesNotRelevant) && !goal.isShieldingTask()) { // Set the values for all maybe-states to 0.5 to indicate that their probability values are neither 0 nor 1. @@ -640,7 +640,7 @@ namespace storm { // Obtain proper hint information either from the provided hint or from requirements of the solver. SparseMdpHintType hintInformation = computeHints(env, SolutionType::UntilProbabilities, hint, goal.direction(), transitionMatrix, backwardTransitions, qualitativeStateSets.maybeStates, phiStates, qualitativeStateSets.statesWithProbability1, produceScheduler); - + // Declare the components of the equation system we will solve. storm::storage::SparseMatrix submatrix; std::vector b; @@ -710,7 +710,7 @@ namespace storm { if (produceScheduler) { extendScheduler(*scheduler, goal, qualitativeStateSets, transitionMatrix, backwardTransitions, phiStates, psiStates); } - + // Sanity check for created scheduler. STORM_LOG_ASSERT(!produceScheduler || scheduler, "Expected that a scheduler was obtained."); STORM_LOG_ASSERT((!produceScheduler && !scheduler) || !scheduler->isPartialScheduler(), "Expected a fully defined scheduler"); @@ -745,42 +745,42 @@ namespace storm { return result; } } - + template template std::vector SparseMdpPrctlHelper::computeInstantaneousRewards(Environment const& env, storm::solver::SolveGoal&& goal, storm::storage::SparseMatrix const& transitionMatrix, RewardModelType const& rewardModel, uint_fast64_t stepCount) { // Only compute the result if the model has a state-based reward this->getModel(). STORM_LOG_THROW(rewardModel.hasStateRewards(), storm::exceptions::InvalidPropertyException, "Missing reward model for formula. Skipping formula."); - + // Initialize result to state rewards of the this->getModel(). std::vector result(rewardModel.getStateRewardVector()); - + auto multiplier = storm::solver::MultiplierFactory().create(env, transitionMatrix); multiplier->repeatedMultiplyAndReduce(env, goal.direction(), result, nullptr, stepCount); return result; } - + template template std::vector SparseMdpPrctlHelper::computeCumulativeRewards(Environment const& env, storm::solver::SolveGoal&& goal, storm::storage::SparseMatrix const& transitionMatrix, RewardModelType const& rewardModel, uint_fast64_t stepBound) { // Only compute the result if the model has at least one reward this->getModel(). STORM_LOG_THROW(!rewardModel.empty(), storm::exceptions::InvalidPropertyException, "Missing reward model for formula. Skipping formula."); - + // Compute the reward vector to add in each step based on the available reward models. std::vector totalRewardVector = rewardModel.getTotalRewardVector(transitionMatrix); - + // Initialize result to the zero vector. std::vector result(transitionMatrix.getRowGroupCount(), storm::utility::zero()); - + auto multiplier = storm::solver::MultiplierFactory().create(env, transitionMatrix); multiplier->repeatedMultiplyAndReduce(env, goal.direction(), result, &totalRewardVector, stepBound); - + return result; } - + template template MDPSparseModelCheckingHelperReturnType SparseMdpPrctlHelper::computeTotalRewards(Environment const& env, storm::solver::SolveGoal&& goal, storm::storage::SparseMatrix const& transitionMatrix, storm::storage::SparseMatrix const& backwardTransitions, RewardModelType const& rewardModel, bool qualitative, bool produceScheduler, ModelCheckerHint const& hint) { @@ -804,7 +804,7 @@ namespace storm { storm::storage::BitVector statesWithoutReward = rewardModel.getStatesWithZeroReward(transitionMatrix); storm::storage::BitVector rew0AStates = storm::utility::graph::performProbGreater0E(backwardTransitions, statesWithoutReward, ~statesWithoutReward); rew0AStates.complement(); - + // There might be end components that consists only of states/choices with zero rewards. The reachability reward semantics would assign such // end components reward infinity. To avoid this, we potentially need to eliminate such end components storm::storage::BitVector trueStates(transitionMatrix.getRowGroupCount(), true); @@ -819,7 +819,7 @@ namespace storm { for (auto oldRew0AState : rew0AStates) { newRew0AStates.set(ecElimResult.oldToNewStateMapping[oldRew0AState]); } - + MDPSparseModelCheckingHelperReturnType result = computeReachabilityRewardsHelper(env, std::move(goal), ecElimResult.matrix, ecElimResult.matrix.transpose(true), [&] (uint_fast64_t rowCount, storm::storage::SparseMatrix const& transitionMatrix, storm::storage::BitVector const& maybeStates) { std::vector result; @@ -850,7 +850,7 @@ namespace storm { } return newChoicesWithoutReward; }); - + std::vector resultInEcQuotient = std::move(result.values); result.values.resize(ecElimResult.oldToNewStateMapping.size()); storm::utility::vector::selectVectorValues(result.values, ecElimResult.oldToNewStateMapping, resultInEcQuotient); @@ -858,7 +858,7 @@ namespace storm { } } } - + template template MDPSparseModelCheckingHelperReturnType SparseMdpPrctlHelper::computeReachabilityRewards(Environment const& env, storm::solver::SolveGoal&& goal, storm::storage::SparseMatrix const& transitionMatrix, storm::storage::SparseMatrix const& backwardTransitions, RewardModelType const& rewardModel, storm::storage::BitVector const& targetStates, bool qualitative, bool produceScheduler, ModelCheckerHint const& hint) { @@ -877,7 +877,7 @@ namespace storm { }, hint); } - + template MDPSparseModelCheckingHelperReturnType SparseMdpPrctlHelper::computeReachabilityTimes(Environment const& env, storm::solver::SolveGoal&& goal, storm::storage::SparseMatrix const& transitionMatrix, storm::storage::SparseMatrix const& backwardTransitions, storm::storage::BitVector const& targetStates, bool qualitative, bool produceScheduler, ModelCheckerHint const& hint) { return computeReachabilityRewardsHelper(env, std::move(goal), transitionMatrix, backwardTransitions, @@ -893,7 +893,7 @@ namespace storm { }, hint); } - + #ifdef STORM_HAVE_CARL template std::vector SparseMdpPrctlHelper::computeReachabilityRewards(Environment const& env, storm::solver::SolveGoal&& goal, storm::storage::SparseMatrix const& transitionMatrix, storm::storage::SparseMatrix const& backwardTransitions, storm::models::sparse::StandardRewardModel const& intervalRewardModel, bool lowerBoundOfIntervals, storm::storage::BitVector const& targetStates, bool qualitative) { @@ -917,13 +917,13 @@ namespace storm { return intervalRewardModel.getChoicesWithFilter(transitionMatrix, [&](storm::Interval const& i) {return storm::utility::isZero(lowerBoundOfIntervals ? i.lower() : i.upper());}); }).values; } - + template<> std::vector SparseMdpPrctlHelper::computeReachabilityRewards(Environment const& env, storm::solver::SolveGoal&&, storm::storage::SparseMatrix const&, storm::storage::SparseMatrix const&, storm::models::sparse::StandardRewardModel const&, bool, storm::storage::BitVector const&, bool) { STORM_LOG_THROW(false, storm::exceptions::IllegalFunctionCallException, "Computing reachability rewards is unsupported for this data type."); } #endif - + struct QualitativeStateSetsReachabilityRewards { storm::storage::BitVector maybeStates; storm::storage::BitVector infinityStates; @@ -934,7 +934,7 @@ namespace storm { QualitativeStateSetsReachabilityRewards getQualitativeStateSetsReachabilityRewardsFromHint(ModelCheckerHint const& hint, storm::storage::BitVector const& targetStates) { QualitativeStateSetsReachabilityRewards result; result.maybeStates = hint.template asExplicitModelCheckerHint().getMaybeStates(); - + // Treat the states with reward zero/infinity. std::vector const& resultsForNonMaybeStates = hint.template asExplicitModelCheckerHint().getResultHint(); result.infinityStates = storm::storage::BitVector(result.maybeStates.size()); @@ -950,7 +950,7 @@ namespace storm { } return result; } - + template QualitativeStateSetsReachabilityRewards computeQualitativeStateSetsReachabilityRewards(storm::solver::SolveGoal const& goal, storm::storage::SparseMatrix const& transitionMatrix, storm::storage::SparseMatrix const& backwardTransitions, storm::storage::BitVector const& targetStates, std::function const& zeroRewardStatesGetter, std::function const& zeroRewardChoicesGetter) { QualitativeStateSetsReachabilityRewards result; @@ -961,7 +961,7 @@ namespace storm { result.infinityStates = storm::utility::graph::performProb1A(transitionMatrix, transitionMatrix.getRowGroupIndices(), backwardTransitions, trueStates, targetStates); } result.infinityStates.complement(); - + if (storm::settings::getModule().isFilterRewZeroSet()) { if (goal.minimize()) { result.rewardZeroStates = storm::utility::graph::performProb1E(transitionMatrix, transitionMatrix.getRowGroupIndices(), backwardTransitions, trueStates, targetStates, zeroRewardChoicesGetter()); @@ -974,7 +974,7 @@ namespace storm { result.maybeStates = ~(result.rewardZeroStates | result.infinityStates); return result; } - + template QualitativeStateSetsReachabilityRewards getQualitativeStateSetsReachabilityRewards(storm::solver::SolveGoal const& goal, storm::storage::SparseMatrix const& transitionMatrix, storm::storage::SparseMatrix const& backwardTransitions, storm::storage::BitVector const& targetStates, ModelCheckerHint const& hint, std::function const& zeroRewardStatesGetter, std::function const& zeroRewardChoicesGetter) { if (hint.isExplicitModelCheckerHint() && hint.template asExplicitModelCheckerHint().getComputeOnlyMaybeStates()) { @@ -983,7 +983,7 @@ namespace storm { return computeQualitativeStateSetsReachabilityRewards(goal, transitionMatrix, backwardTransitions, targetStates, zeroRewardStatesGetter, zeroRewardChoicesGetter); } } - + template void extendScheduler(storm::storage::Scheduler& scheduler, storm::solver::SolveGoal const& goal, QualitativeStateSetsReachabilityRewards const& qualitativeStateSets, storm::storage::SparseMatrix const& transitionMatrix, storm::storage::SparseMatrix const& backwardTransitions, storm::storage::BitVector const& targetStates, std::function const& zeroRewardChoicesGetter) { // Finally, if we need to produce a scheduler, we also need to figure out the parts of the scheduler for @@ -1000,7 +1000,7 @@ namespace storm { } } } - + template void extractSchedulerChoices(storm::storage::Scheduler& scheduler, storm::storage::SparseMatrix const& transitionMatrix, std::vector const& subChoices, storm::storage::BitVector const& maybeStates, boost::optional const& selectedChoices) { auto subChoiceIt = subChoices.begin(); @@ -1023,7 +1023,7 @@ namespace storm { } assert(subChoiceIt == subChoices.end()); } - + template void computeFixedPointSystemReachabilityRewards(storm::solver::SolveGoal& goal, storm::storage::SparseMatrix const& transitionMatrix, QualitativeStateSetsReachabilityRewards const& qualitativeStateSets, boost::optional const& selectedChoices, std::function(uint_fast64_t, storm::storage::SparseMatrix const&, storm::storage::BitVector const&)> const& totalStateRewardVectorGetter, storm::storage::SparseMatrix& submatrix, std::vector& b, std::vector* oneStepTargetProbabilities = nullptr) { // Remove rows and columns from the original transition probability matrix for states whose reward values are already known. @@ -1042,20 +1042,20 @@ namespace storm { (*oneStepTargetProbabilities) = transitionMatrix.getConstrainedRowSumVector(*selectedChoices, qualitativeStateSets.rewardZeroStates); } } - + // If the solve goal has relevant values, we need to adjust them. goal.restrictRelevantValues(qualitativeStateSets.maybeStates); } - + template boost::optional> computeFixedPointSystemReachabilityRewardsEliminateEndComponents(storm::solver::SolveGoal& goal, storm::storage::SparseMatrix const& transitionMatrix, storm::storage::SparseMatrix const& backwardTransitions, QualitativeStateSetsReachabilityRewards const& qualitativeStateSets, boost::optional const& selectedChoices, std::function(uint_fast64_t, storm::storage::SparseMatrix const&, storm::storage::BitVector const&)> const& totalStateRewardVectorGetter, storm::storage::SparseMatrix& submatrix, std::vector& b, boost::optional>& oneStepTargetProbabilities, bool produceScheduler) { - + // Start by computing the choices with reward 0, as we only want ECs within this fragment. storm::storage::BitVector zeroRewardChoices(transitionMatrix.getRowCount()); // Get the rewards of all choices. std::vector rewardVector = totalStateRewardVectorGetter(transitionMatrix.getRowCount(), transitionMatrix, storm::storage::BitVector(transitionMatrix.getRowGroupCount(), true)); - + uint64_t index = 0; for (auto const& e : rewardVector) { if (storm::utility::isZero(e)) { @@ -1063,40 +1063,40 @@ namespace storm { } ++index; } - + // Compute the states that have some zero reward choice. storm::storage::BitVector candidateStates(qualitativeStateSets.maybeStates); for (auto state : qualitativeStateSets.maybeStates) { bool keepState = false; - + for (auto row = transitionMatrix.getRowGroupIndices()[state], rowEnd = transitionMatrix.getRowGroupIndices()[state + 1]; row < rowEnd; ++row) { if (zeroRewardChoices.get(row)) { keepState = true; break; } } - + if (!keepState) { candidateStates.set(state, false); } } - + // Only keep the candidate states that (under some scheduler) can stay in the set of candidates forever candidateStates = storm::utility::graph::performProb0E(transitionMatrix, transitionMatrix.getRowGroupIndices(), backwardTransitions, candidateStates, ~candidateStates); - + bool doDecomposition = !candidateStates.empty(); - + storm::storage::MaximalEndComponentDecomposition endComponentDecomposition; if (doDecomposition) { // Then compute the states that are in MECs with zero reward. endComponentDecomposition = storm::storage::MaximalEndComponentDecomposition(transitionMatrix, backwardTransitions, candidateStates, zeroRewardChoices); } - + // Only do more work if there are actually end-components. if (doDecomposition && !endComponentDecomposition.empty()) { STORM_LOG_DEBUG("Eliminating " << endComponentDecomposition.size() << " ECs."); SparseMdpEndComponentInformation result = SparseMdpEndComponentInformation::eliminateEndComponents(endComponentDecomposition, transitionMatrix, qualitativeStateSets.maybeStates, oneStepTargetProbabilities ? &qualitativeStateSets.rewardZeroStates : nullptr, selectedChoices ? &selectedChoices.get() : nullptr, &rewardVector, submatrix, oneStepTargetProbabilities ? &oneStepTargetProbabilities.get() : nullptr, &b, produceScheduler); - + // If the solve goal has relevant values, we need to adjust them. if (goal.hasRelevantValues()) { storm::storage::BitVector newRelevantValues(submatrix.getRowGroupCount()); @@ -1109,7 +1109,7 @@ namespace storm { goal.setRelevantValues(std::move(newRelevantValues)); } } - + return result; } else { STORM_LOG_DEBUG("Not eliminating ECs as there are none."); @@ -1117,10 +1117,10 @@ namespace storm { return boost::none; } } - + template void computeUpperRewardBounds(SparseMdpHintType& hintInformation, storm::OptimizationDirection const& direction, storm::storage::SparseMatrix const& submatrix, std::vector const& choiceRewards, std::vector const& oneStepTargetProbabilities) { - + // For the min-case, we use DS-MPI, for the max-case variant 2 of the Baier et al. paper (CAV'17). if (direction == storm::OptimizationDirection::Minimize) { DsMpiMdpUpperRewardBoundsComputer dsmpi(submatrix, choiceRewards, oneStepTargetProbabilities); @@ -1130,29 +1130,29 @@ namespace storm { hintInformation.upperResultBound = baier.computeUpperBound(); } } - + template MDPSparseModelCheckingHelperReturnType SparseMdpPrctlHelper::computeReachabilityRewardsHelper(Environment const& env, storm::solver::SolveGoal&& goal, storm::storage::SparseMatrix const& transitionMatrix, storm::storage::SparseMatrix const& backwardTransitions, std::function(uint_fast64_t, storm::storage::SparseMatrix const&, storm::storage::BitVector const&)> const& totalStateRewardVectorGetter, storm::storage::BitVector const& targetStates, bool qualitative, bool produceScheduler, std::function const& zeroRewardStatesGetter, std::function const& zeroRewardChoicesGetter, ModelCheckerHint const& hint) { - + // Prepare resulting vector. std::vector result(transitionMatrix.getRowGroupCount(), storm::utility::zero()); // Determine which states have a reward that is infinity or less than infinity. QualitativeStateSetsReachabilityRewards qualitativeStateSets = getQualitativeStateSetsReachabilityRewards(goal, transitionMatrix, backwardTransitions, targetStates, hint, zeroRewardStatesGetter, zeroRewardChoicesGetter); - + STORM_LOG_INFO("Preprocessing: " << qualitativeStateSets.infinityStates.getNumberOfSetBits() << " states with reward infinity, " << qualitativeStateSets.rewardZeroStates.getNumberOfSetBits() << " states with reward zero (" << qualitativeStateSets.maybeStates.getNumberOfSetBits() << " states remaining)."); storm::utility::vector::setVectorValues(result, qualitativeStateSets.infinityStates, storm::utility::infinity()); - + // If requested, we will produce a scheduler. std::unique_ptr> scheduler; if (produceScheduler) { scheduler = std::make_unique>(transitionMatrix.getRowGroupCount()); } - + // Check if the values of the maybe states are relevant for the SolveGoal bool maybeStatesNotRelevant = goal.hasRelevantValues() && goal.relevantValues().isDisjointFrom(qualitativeStateSets.maybeStates); - + // Check whether we need to compute exact rewards for some states. if (qualitative || maybeStatesNotRelevant) { STORM_LOG_INFO("The rewards for the initial states were determined in a preprocessing step. No exact rewards were computed."); @@ -1168,10 +1168,10 @@ namespace storm { if (!qualitativeStateSets.infinityStates.empty()) { selectedChoices = transitionMatrix.getRowFilter(qualitativeStateSets.maybeStates, ~qualitativeStateSets.infinityStates); } - + // Obtain proper hint information either from the provided hint or from requirements of the solver. SparseMdpHintType hintInformation = computeHints(env, SolutionType::ExpectedRewards, hint, goal.direction(), transitionMatrix, backwardTransitions, qualitativeStateSets.maybeStates, ~qualitativeStateSets.rewardZeroStates, qualitativeStateSets.rewardZeroStates, produceScheduler, selectedChoices); - + // Declare the components of the equation system we will solve. storm::storage::SparseMatrix submatrix; std::vector b; @@ -1182,7 +1182,7 @@ namespace storm { if (hintInformation.getComputeUpperBounds()) { oneStepTargetProbabilities = std::vector(); } - + // If the hint information tells us that we have to eliminate MECs, we do so now. boost::optional> ecInformation; if (hintInformation.getEliminateEndComponents()) { @@ -1191,13 +1191,13 @@ namespace storm { // Otherwise, we compute the standard equations. computeFixedPointSystemReachabilityRewards(goal, transitionMatrix, qualitativeStateSets, selectedChoices, totalStateRewardVectorGetter, submatrix, b, oneStepTargetProbabilities ? &oneStepTargetProbabilities.get() : nullptr); } - + // If we need to compute upper bounds, do so now. if (hintInformation.getComputeUpperBounds()) { STORM_LOG_ASSERT(oneStepTargetProbabilities, "Expecting one step target probability vector to be available."); computeUpperRewardBounds(hintInformation, goal.direction(), submatrix, b, oneStepTargetProbabilities.get()); } - + // Now compute the results for the maybe states. MaybeStateResult resultForMaybeStates = computeValuesForMaybeStates(env, std::move(goal), std::move(submatrix), b, produceScheduler, hintInformation); @@ -1216,12 +1216,12 @@ namespace storm { } } } - + // Extend scheduler with choices for the states in the qualitative state sets. if (produceScheduler) { extendScheduler(*scheduler, goal, qualitativeStateSets, transitionMatrix, backwardTransitions, targetStates, zeroRewardChoicesGetter); } - + // Sanity check for created scheduler. STORM_LOG_ASSERT(!produceScheduler || scheduler, "Expected that a scheduler was obtained."); STORM_LOG_ASSERT((!produceScheduler && !scheduler) || !scheduler->isPartialScheduler(), "Expected a fully defined scheduler"); @@ -1232,12 +1232,12 @@ namespace storm { return MDPSparseModelCheckingHelperReturnType(std::move(result), std::move(qualitativeStateSets.maybeStates), std::move(scheduler), std::move(choiceValues)); } - + template std::unique_ptr SparseMdpPrctlHelper::computeConditionalProbabilities(Environment const& env, storm::solver::SolveGoal&& goal, storm::storage::SparseMatrix const& transitionMatrix, storm::storage::SparseMatrix const& backwardTransitions, storm::storage::BitVector const& targetStates, storm::storage::BitVector const& conditionStates) { - + std::chrono::high_resolution_clock::time_point start = std::chrono::high_resolution_clock::now(); - + // For the max-case, we can simply take the given target states. For the min-case, however, we need to // find the MECs of non-target states and make them the new target states. storm::storage::BitVector fixedTargetStates; @@ -1252,18 +1252,18 @@ namespace storm { } } } - + storm::storage::BitVector allStates(fixedTargetStates.size(), true); - + // Extend the target states by computing all states that have probability 1 to go to a target state // under *all* schedulers. fixedTargetStates = storm::utility::graph::performProb1A(transitionMatrix, transitionMatrix.getRowGroupIndices(), backwardTransitions, allStates, fixedTargetStates); - + // We solve the max-case and later adjust the result if the optimization direction was to minimize. storm::storage::BitVector initialStatesBitVector = goal.relevantValues(); STORM_LOG_THROW(initialStatesBitVector.getNumberOfSetBits() == 1, storm::exceptions::NotSupportedException, "Computing conditional probabilities in MDPs is only supported for models with exactly one initial state."); storm::storage::sparse::state_type initialState = *initialStatesBitVector.begin(); - + // Extend the condition states by computing all states that have probability 1 to go to a condition state // under *all* schedulers. storm::storage::BitVector extendedConditionStates = storm::utility::graph::performProb1A(transitionMatrix, transitionMatrix.getRowGroupIndices(), backwardTransitions, allStates, conditionStates); @@ -1273,12 +1273,12 @@ namespace storm { std::vector conditionProbabilities = std::move(computeUntilProbabilities(env, OptimizationDirection::Maximize, transitionMatrix, backwardTransitions, allStates, extendedConditionStates, false, false).values); std::chrono::high_resolution_clock::time_point conditionEnd = std::chrono::high_resolution_clock::now(); STORM_LOG_DEBUG("Computed probabilities to satisfy for condition in " << std::chrono::duration_cast(conditionEnd - conditionStart).count() << "ms."); - + // If the conditional probability is undefined for the initial state, we return directly. if (storm::utility::isZero(conditionProbabilities[initialState])) { return std::unique_ptr(new ExplicitQuantitativeCheckResult(initialState, storm::utility::infinity())); } - + STORM_LOG_DEBUG("Computing probabilities to reach target."); std::chrono::high_resolution_clock::time_point targetStart = std::chrono::high_resolution_clock::now(); std::vector targetProbabilities = std::move(computeUntilProbabilities(env, OptimizationDirection::Maximize, transitionMatrix, backwardTransitions, allStates, fixedTargetStates, false, false).values); @@ -1306,7 +1306,7 @@ namespace storm { storm::storage::sparse::state_type newGoalState = relevantStates.getNumberOfSetBits(); storm::storage::sparse::state_type newStopState = newGoalState + 1; storm::storage::sparse::state_type newFailState = newStopState + 1; - + // Build the transitions of the (relevant) states of the original model. storm::storage::SparseMatrixBuilder builder(0, newFailState + 1, 0, true, true); uint_fast64_t currentRow = 0; @@ -1344,7 +1344,7 @@ namespace storm { } } } - + // Now build the transitions of the newly introduced states. builder.newRowGroup(currentRow); builder.addNextValue(currentRow, newGoalState, storm::utility::one()); @@ -1355,10 +1355,10 @@ namespace storm { builder.newRowGroup(currentRow); builder.addNextValue(currentRow, numberOfStatesBeforeRelevantStates[initialState], storm::utility::one()); ++currentRow; - + std::chrono::high_resolution_clock::time_point end = std::chrono::high_resolution_clock::now(); STORM_LOG_DEBUG("Computed transformed model in " << std::chrono::duration_cast(end - start).count() << "ms."); - + // Finally, build the matrix and dispatch the query as a reachability query. STORM_LOG_DEBUG("Computing conditional probabilties."); storm::storage::BitVector newGoalStates(newFailState + 1); @@ -1366,20 +1366,20 @@ namespace storm { storm::storage::SparseMatrix newTransitionMatrix = builder.build(); STORM_LOG_DEBUG("Transformed model has " << newTransitionMatrix.getRowGroupCount() << " states and " << newTransitionMatrix.getNonzeroEntryCount() << " transitions."); storm::storage::SparseMatrix newBackwardTransitions = newTransitionMatrix.transpose(true); - + storm::solver::OptimizationDirection dir = goal.direction(); if (goal.minimize()) { goal.oneMinus(); } - + std::chrono::high_resolution_clock::time_point conditionalStart = std::chrono::high_resolution_clock::now(); std::vector goalProbabilities = std::move(computeUntilProbabilities(env, std::move(goal), newTransitionMatrix, newBackwardTransitions, storm::storage::BitVector(newFailState + 1, true), newGoalStates, false, false).values); std::chrono::high_resolution_clock::time_point conditionalEnd = std::chrono::high_resolution_clock::now(); STORM_LOG_DEBUG("Computed conditional probabilities in transformed model in " << std::chrono::duration_cast(conditionalEnd - conditionalStart).count() << "ms."); - + return std::unique_ptr(new ExplicitQuantitativeCheckResult(initialState, dir == OptimizationDirection::Maximize ? goalProbabilities[numberOfStatesBeforeRelevantStates[initialState]] : storm::utility::one() - goalProbabilities[numberOfStatesBeforeRelevantStates[initialState]])); } - + template class SparseMdpPrctlHelper; template std::vector SparseMdpPrctlHelper::computeInstantaneousRewards(Environment const& env, storm::solver::SolveGoal&& goal, storm::storage::SparseMatrix const& transitionMatrix, storm::models::sparse::StandardRewardModel const& rewardModel, uint_fast64_t stepCount); template std::vector SparseMdpPrctlHelper::computeCumulativeRewards(Environment const& env, storm::solver::SolveGoal&& goal, storm::storage::SparseMatrix const& transitionMatrix, storm::models::sparse::StandardRewardModel const& rewardModel, uint_fast64_t stepBound); diff --git a/src/storm/storage/Scheduler.h b/src/storm/storage/Scheduler.h index 43c7bb1cb..5107bb1cf 100644 --- a/src/storm/storage/Scheduler.h +++ b/src/storm/storage/Scheduler.h @@ -5,6 +5,7 @@ #include "storm/storage/memorystructure/MemoryStructure.h" #include "storm/storage/SchedulerChoice.h" +#include "storm/storage/BitVector.h" namespace storm { namespace storage { @@ -65,7 +66,7 @@ namespace storm { storm::storage::BitVector computeActionSupport(std::vector const& nondeterministicChoiceIndicies) const; /*! - * Retrieves whether there is a pair of model and memory state for which the choice is undefined. + * Retrieves whether there is a *reachable* pair of model and memory state for which the choice is undefined. */ bool isPartialScheduler() const; @@ -127,8 +128,11 @@ namespace storm { std::vector>> schedulerChoices; bool printUndefinedChoices = false; - uint_fast64_t numOfUndefinedChoices; + + std::vector reachableStates; + uint_fast64_t numOfUndefinedChoices; // Only consider reachable ones uint_fast64_t numOfDeterministicChoices; + uint_fast64_t numOfUnreachableStates; }; } } diff --git a/src/storm/storage/memorystructure/MemoryStructureBuilder.h b/src/storm/storage/memorystructure/MemoryStructureBuilder.h index bb310ccf3..0b3c6e9cf 100644 --- a/src/storm/storage/memorystructure/MemoryStructureBuilder.h +++ b/src/storm/storage/memorystructure/MemoryStructureBuilder.h @@ -19,7 +19,7 @@ namespace storm { * @param numberOfMemoryStates The number of states the resulting memory structure should have */ MemoryStructureBuilder(uint_fast64_t numberOfMemoryStates, storm::models::sparse::Model const& model); - + // TODO: Add variant with a flag: Consider non-initial model states /*! * Initializes a new builder with the data from the provided memory structure */