Storm 1.10.0.1
A Modern Probabilistic Model Checker
Loading...
Searching...
No Matches
ConditionalHelper.cpp
Go to the documentation of this file.
1#include <algorithm>
2#include <iterator>
3
5
19#include "storm/utility/graph.h"
21
24
25namespace storm::modelchecker {
26
27namespace internal {
28
29template<typename ValueType>
30void eliminateEndComponents(storm::storage::BitVector possibleEcStates, bool addRowAtRepresentativeState, std::optional<uint64_t> representativeRowEntry,
31 storm::storage::SparseMatrix<ValueType>& matrix, uint64_t& initialState, storm::storage::BitVector& rowsWithSum1,
32 std::vector<ValueType>& rowValues1, storm::OptionalRef<std::vector<ValueType>> rowValues2 = {}) {
33 storm::storage::MaximalEndComponentDecomposition<ValueType> ecs(matrix, matrix.transpose(true), possibleEcStates, rowsWithSum1);
34 if (ecs.empty()) {
35 return; // nothing to do
36 }
37
38 storm::storage::BitVector allRowGroups(matrix.getRowGroupCount(), true);
40 matrix, ecs, allRowGroups, addRowAtRepresentativeState ? allRowGroups : ~allRowGroups, representativeRowEntry.has_value());
41
42 // Update matrix
43 matrix = std::move(ecElimResult.matrix);
44 if (addRowAtRepresentativeState && representativeRowEntry) {
45 auto const columnIndex = ecElimResult.oldToNewStateMapping[*representativeRowEntry];
46 for (auto representativeRowIndex : ecElimResult.sinkRows) {
47 auto row = matrix.getRow(representativeRowIndex);
48 STORM_LOG_ASSERT(row.getNumberOfEntries() == 1, "unexpected number of entries in representative row.");
49 auto& entry = *row.begin();
50 entry.setColumn(columnIndex);
51 }
52 }
53
54 // update vectors
55 auto updateRowValue = [&ecElimResult](std::vector<ValueType>& rowValues) {
56 std::vector<ValueType> newRowValues;
57 newRowValues.reserve(ecElimResult.newToOldRowMapping.size());
58 for (auto oldRowIndex : ecElimResult.newToOldRowMapping) {
59 newRowValues.push_back(rowValues[oldRowIndex]);
60 }
61 rowValues = std::move(newRowValues);
63 std::all_of(ecElimResult.sinkRows.begin(), ecElimResult.sinkRows.end(), [&rowValues](auto i) { return storm::utility::isZero(rowValues[i]); }),
64 "Sink rows are expected to have zero value");
65 };
66 updateRowValue(rowValues1);
67 if (rowValues2) {
68 updateRowValue(*rowValues2);
69 }
70
71 // update initial state
72 initialState = ecElimResult.oldToNewStateMapping[initialState];
73
74 // update bitvector
75 storm::storage::BitVector newRowsWithSum1(ecElimResult.newToOldRowMapping.size(), true);
76 uint64_t newRowIndex = 0;
77 for (auto oldRowIndex : ecElimResult.newToOldRowMapping) {
78 if ((addRowAtRepresentativeState && !representativeRowEntry.has_value() && ecElimResult.sinkRows.get(newRowIndex)) || !rowsWithSum1.get(oldRowIndex)) {
79 newRowsWithSum1.set(newRowIndex, false);
80 }
81 ++newRowIndex;
82 }
83 rowsWithSum1 = std::move(newRowsWithSum1);
84}
85
86template<typename ValueType, typename SolutionType = ValueType>
88 std::vector<ValueType> const& rowValues, storm::storage::BitVector const& rowsWithSum1,
89 storm::solver::OptimizationDirection const dir, uint64_t const initialState) {
90 // Initialize the solution vector.
91 std::vector<SolutionType> x(matrix.getRowGroupCount(), storm::utility::zero<ValueType>());
92
93 // Set up the solver.
95 solver->setOptimizationDirection(dir);
96 solver->setRequirementsChecked();
97 solver->setHasUniqueSolution(true);
98 solver->setHasNoEndComponents(true);
99 solver->setLowerBound(storm::utility::zero<ValueType>());
100 solver->setUpperBound(storm::utility::one<ValueType>());
101
102 // Solve the corresponding system of equations.
103 solver->solveEquations(env, x, rowValues);
104 return x[initialState];
105}
106
111template<typename ValueType>
112void computeReachabilityProbabilities(Environment const& env, std::map<uint64_t, ValueType>& nonZeroResults, storm::solver::OptimizationDirection const dir,
113 storm::storage::SparseMatrix<ValueType> const& transitionMatrix, storm::storage::BitVector const& initialStates,
114 storm::storage::BitVector const& allowedStates, storm::storage::BitVector const& targetStates) {
115 if (initialStates.empty()) { // nothing to do
116 return;
117 }
118 auto const reachableStates = storm::utility::graph::getReachableStates(transitionMatrix, initialStates, allowedStates, targetStates);
119 auto const subTargets = targetStates % reachableStates;
120 // Catch the case where no target is reachable from an initial state. In this case, there is nothing to do since all probabilities are zero.
121 if (subTargets.empty()) {
122 return;
123 }
124 auto const subInits = initialStates % reachableStates;
125 auto const submatrix = transitionMatrix.getSubmatrix(true, reachableStates, reachableStates);
127 env, storm::solver::SolveGoal<ValueType>(dir, subInits), submatrix, submatrix.transpose(true), storm::storage::BitVector(subTargets.size(), true),
128 subTargets, false, false);
129 auto origInitIt = initialStates.begin();
130 for (auto subInit : subInits) {
131 auto const& val = subResult.values[subInit];
132 if (!storm::utility::isZero(val)) {
133 nonZeroResults.emplace(*origInitIt, val);
134 }
135 ++origInitIt;
136 }
137}
138
139template<typename ValueType>
141 storm::storage::BitVector const maybeStates; // Those states that can be reached from initial without reaching a terminal state
142 storm::storage::BitVector const terminalStates; // Those states where we already know the probability to reach the condition and the target value
143 storm::storage::BitVector const conditionStates; // Those states where the condition holds almost surely (under all schedulers)
144 storm::storage::BitVector const universalObservationFailureStates; // Those states where the condition is not reachable (under all schedulers)
145 storm::storage::BitVector const existentialObservationFailureStates; // Those states s where a scheduler exists that (i) does not reach the condition from
146 // s and (ii) acts optimal in all terminal states
147 std::map<uint64_t, ValueType> const nonZeroTargetStateValues; // The known non-zero target values. (default is zero)
148 // There are three cases of terminal states:
149 // 1. conditionStates: The condition holds, so the target value is the optimal probability to reach target from there
150 // 2. targetStates: The target is reached, so the target value is the optimal probability to reach a condition from there.
151 // The remaining probability mass is the probability of an observation failure
152 // 3. states that can not reach the condition under any scheduler. The target value is zero.
153
154 // TerminalStates is a superset of conditionStates and dom(nonZeroTargetStateValues).
155 // For a terminalState that is not a conditionState, it is impossible to (reach the condition and not reach the target).
156
157 ValueType getTargetValue(uint64_t state) const {
158 STORM_LOG_ASSERT(terminalStates.get(state), "Tried to get target value for non-terminal state");
159 auto const it = nonZeroTargetStateValues.find(state);
160 return it == nonZeroTargetStateValues.end() ? storm::utility::zero<ValueType>() : it->second;
161 }
162
163 ValueType failProbability(uint64_t state) const {
164 STORM_LOG_ASSERT(terminalStates.get(state), "Tried to get fail probability for non-terminal state");
165 STORM_LOG_ASSERT(!conditionStates.get(state), "Tried to get fail probability for a condition state");
166 // condition states have fail probability zero
167 return storm::utility::one<ValueType>() - getTargetValue(state);
168 }
169};
170
171template<typename ValueType>
173 storm::storage::SparseMatrix<ValueType> const& transitionMatrix,
174 storm::storage::SparseMatrix<ValueType> const& backwardTransitions, storm::storage::BitVector const& relevantStates,
175 storm::storage::BitVector const& targetStates, storm::storage::BitVector const& conditionStates) {
176 storm::storage::BitVector const allStates(transitionMatrix.getRowGroupCount(), true);
177 auto extendedConditionStates =
178 storm::utility::graph::performProb1A(transitionMatrix, transitionMatrix.getRowGroupIndices(), backwardTransitions, allStates, conditionStates);
179 auto universalObservationFailureStates = storm::utility::graph::performProb0A(backwardTransitions, allStates, extendedConditionStates);
180 std::map<uint64_t, ValueType> nonZeroTargetStateValues;
181 auto const extendedTargetStates =
182 storm::utility::graph::performProb1A(transitionMatrix, transitionMatrix.getRowGroupIndices(), backwardTransitions, allStates, targetStates);
183 computeReachabilityProbabilities(env, nonZeroTargetStateValues, dir, transitionMatrix, extendedConditionStates, allStates, extendedTargetStates);
184 auto const targetAndNotCondFailStates = extendedTargetStates & ~(extendedConditionStates | universalObservationFailureStates);
185 computeReachabilityProbabilities(env, nonZeroTargetStateValues, dir, transitionMatrix, targetAndNotCondFailStates, allStates, extendedConditionStates);
186
187 // get states where the optimal policy reaches the condition with positive probability
188 auto terminalStatesThatReachCondition = extendedConditionStates;
189 for (auto state : targetAndNotCondFailStates) {
190 if (nonZeroTargetStateValues.contains(state)) {
191 terminalStatesThatReachCondition.set(state, true);
192 }
193 }
194
195 // get the terminal states following the three cases described above
196 auto terminalStates = extendedConditionStates | extendedTargetStates | universalObservationFailureStates;
197 if (storm::solver::minimize(dir)) {
198 // There can be target states from which (only) the *minimal* probability to reach a condition is zero.
199 // For those states, the optimal policy is to enforce observation failure.
200 // States that can only reach (target states with almost sure observation failure) or observation failure will be treated as terminal states with
201 // targetValue zero and failProbability one.
202 terminalStates |= storm::utility::graph::performProb0A(backwardTransitions, ~terminalStates, terminalStatesThatReachCondition);
203 }
204
205 auto nonTerminalStates = ~terminalStates;
206
207 auto existentialObservationFailureStates = storm::utility::graph::performProb0E(transitionMatrix, transitionMatrix.getRowGroupIndices(),
208 backwardTransitions, nonTerminalStates, terminalStatesThatReachCondition);
209
210 // Restrict non-terminal states to those that are still relevant
211 nonTerminalStates &= storm::utility::graph::getReachableStates(transitionMatrix, relevantStates, nonTerminalStates, terminalStates);
212
213 return NormalFormData<ValueType>{.maybeStates = std::move(nonTerminalStates),
214 .terminalStates = std::move(terminalStates),
215 .conditionStates = std::move(extendedConditionStates),
216 .universalObservationFailureStates = std::move(universalObservationFailureStates),
217 .existentialObservationFailureStates = std::move(existentialObservationFailureStates),
218 .nonZeroTargetStateValues = std::move(nonZeroTargetStateValues)};
219}
220
225template<typename ValueType, typename SolutionType = ValueType>
226SolutionType computeViaRestartMethod(Environment const& env, uint64_t const initialState, storm::solver::OptimizationDirection const dir,
227 storm::storage::SparseMatrix<ValueType> const& transitionMatrix, NormalFormData<ValueType> const& normalForm) {
228 auto const& maybeStates = normalForm.maybeStates;
229 auto const stateToMatrixIndexMap = maybeStates.getNumberOfSetBitsBeforeIndices();
230 auto const numMaybeStates = maybeStates.getNumberOfSetBits();
231 auto const numMaybeChoices = transitionMatrix.getNumRowsInRowGroups(maybeStates);
232
233 // Build the transitions that include a backwards loop to the initial state
234 storm::storage::SparseMatrixBuilder<ValueType> matrixBuilder(numMaybeChoices, numMaybeStates, 0, true, true, numMaybeStates);
235 std::vector<ValueType> rowValues;
236 storm::storage::BitVector rowsWithSum1(numMaybeChoices, true);
237 rowValues.reserve(numMaybeChoices);
238 uint64_t currentRow = 0;
239 for (auto state : maybeStates) {
240 matrixBuilder.newRowGroup(currentRow);
241 for (auto origRowIndex : transitionMatrix.getRowGroupIndices(state)) {
242 // We make two passes over the successors. First, we find out the reset probabilities and target probabilities
243 // Then, we insert the matrix entries in the correct order
244 // This two-phase approach is to avoid a costly out-of-order insertion into the matrix
245 ValueType targetProbability = storm::utility::zero<ValueType>();
246 ValueType restartProbability = storm::utility::zero<ValueType>();
247 bool rowSumIsLess1 = false;
248 for (auto const& entry : transitionMatrix.getRow(origRowIndex)) {
249 if (normalForm.terminalStates.get(entry.getColumn())) {
250 ValueType const targetValue = normalForm.getTargetValue(entry.getColumn());
251 targetProbability += targetValue * entry.getValue();
252 if (normalForm.conditionStates.get(entry.getColumn())) {
253 rowSumIsLess1 = true;
254 } else {
255 if (!storm::utility::isZero(targetValue)) {
256 rowSumIsLess1 = true;
257 }
258 restartProbability += entry.getValue() * normalForm.failProbability(entry.getColumn());
259 }
260 }
261 }
262 if (rowSumIsLess1) {
263 rowsWithSum1.set(currentRow, false);
264 }
265 rowValues.push_back(targetProbability);
266 bool addRestartTransition = !storm::utility::isZero(restartProbability);
267 for (auto const& entry : transitionMatrix.getRow(origRowIndex)) {
268 // Insert backloop probability if we haven't done so yet and are past the initial state index
269 // This is to avoid a costly out-of-order insertion into the matrix
270 if (addRestartTransition && entry.getColumn() > initialState) {
271 matrixBuilder.addNextValue(currentRow, stateToMatrixIndexMap[initialState], restartProbability);
272 addRestartTransition = false;
273 }
274 if (maybeStates.get(entry.getColumn())) {
275 matrixBuilder.addNextValue(currentRow, stateToMatrixIndexMap[entry.getColumn()], entry.getValue());
276 }
277 }
278 // Add the backloop if we haven't done this already
279 if (addRestartTransition) {
280 matrixBuilder.addNextValue(currentRow, stateToMatrixIndexMap[initialState], restartProbability);
281 }
282 ++currentRow;
283 }
284 }
285
286 auto matrix = matrixBuilder.build();
287 auto initStateInMatrix = stateToMatrixIndexMap[initialState];
288
289 // Eliminate end components in two phases
290 // First, we catch all end components that do not contain the initial state. It is possible to stay in those ECs forever
291 // without reaching the condition. This is reflected by a backloop to the initial state.
292 storm::storage::BitVector selectedStatesInMatrix(numMaybeStates, true);
293 selectedStatesInMatrix.set(initStateInMatrix, false);
294 eliminateEndComponents(selectedStatesInMatrix, true, initStateInMatrix, matrix, initStateInMatrix, rowsWithSum1, rowValues);
295 // Second, eliminate the remaining ECs. These must involve the initial state and might have been introduced in the previous step.
296 // A policy selecting such an EC must reach the condition with probability zero and is thus invalid.
297 selectedStatesInMatrix.set(initStateInMatrix, true);
298 eliminateEndComponents(selectedStatesInMatrix, false, std::nullopt, matrix, initStateInMatrix, rowsWithSum1, rowValues);
299
300 STORM_LOG_INFO("Processed model has " << matrix.getRowGroupCount() << " states and " << matrix.getRowGroupCount() << " choices and "
301 << matrix.getEntryCount() << " transitions.");
302 // Finally, solve the equation system
303 return solveMinMaxEquationSystem(env, matrix, rowValues, rowsWithSum1, dir, initStateInMatrix);
304}
305
311template<typename ValueType, typename SolutionType = ValueType>
313 public:
314 WeightedReachabilityHelper(uint64_t const initialState, storm::storage::SparseMatrix<ValueType> const& transitionMatrix,
315 NormalFormData<ValueType> const& normalForm) {
316 // Determine rowgroups (states) and rows (choices) of the submatrix
317 auto subMatrixRowGroups = normalForm.maybeStates;
318 // Identify and eliminate the initial component to enforce that it is eventually exited
319 // The initial component is the largest subset of maybestates C such that
320 // (i) the initial state is contained in C
321 // (ii) each state in C can be reached from the initial state while only playing actions that stay inside C or observation failure and
322 // (iii) for each state in C except the initial state there is a policy that almost surely reaches an observation failure
323 // An optimal scheduler can intuitively pick the best exiting action of C and enforce that all paths that satisfy the condition exit C through that
324 // action. By eliminating the initial component, we ensure that only policies that actually exit C are considered. The remaining policies have
325 // probability zero of satisfying the condition.
326 storm::storage::BitVector initialComponentExitRows(transitionMatrix.getRowCount(), false);
327 subMatrixRowGroups.set(initialState, false); // temporarily unset initial state
328 std::vector<uint64_t> dfsStack = {initialState};
329 while (!dfsStack.empty()) {
330 auto const state = dfsStack.back();
331 dfsStack.pop_back();
332 for (auto rowIndex : transitionMatrix.getRowGroupIndices(state)) {
333 auto const row = transitionMatrix.getRow(rowIndex);
334 if (std::all_of(row.begin(), row.end(),
335 [&normalForm](auto const& entry) { return normalForm.existentialObservationFailureStates.get(entry.getColumn()); })) {
336 for (auto const& entry : row) {
337 auto const& successor = entry.getColumn();
338 if (subMatrixRowGroups.get(successor)) {
339 subMatrixRowGroups.set(successor, false);
340 dfsStack.push_back(successor);
341 }
342 }
343 } else {
344 initialComponentExitRows.set(rowIndex, true);
345 }
346 }
347 }
348 auto const numSubmatrixRows = transitionMatrix.getNumRowsInRowGroups(subMatrixRowGroups) + initialComponentExitRows.getNumberOfSetBits();
349 subMatrixRowGroups.set(initialState, true); // set initial state again, as single representative state for the initial component
350 auto const numSubmatrixRowGroups = subMatrixRowGroups.getNumberOfSetBits();
351
352 // state index mapping and initial state
353 auto stateToMatrixIndexMap = subMatrixRowGroups.getNumberOfSetBitsBeforeIndices();
354 initialStateInSubmatrix = stateToMatrixIndexMap[initialState];
355 auto const eliminatedInitialComponentStates = normalForm.maybeStates & ~subMatrixRowGroups;
356 for (auto state : eliminatedInitialComponentStates) {
357 stateToMatrixIndexMap[state] = initialStateInSubmatrix; // map all eliminated states to the initial state
358 }
359
360 // build matrix, rows that sum up to 1, target values, condition values
361 storm::storage::SparseMatrixBuilder<ValueType> matrixBuilder(numSubmatrixRows, numSubmatrixRowGroups, 0, true, true, numSubmatrixRowGroups);
362 rowsWithSum1 = storm::storage::BitVector(numSubmatrixRows, true);
363 targetRowValues.reserve(numSubmatrixRows);
364 conditionRowValues.reserve(numSubmatrixRows);
365 uint64_t currentRow = 0;
366 for (auto state : subMatrixRowGroups) {
367 matrixBuilder.newRowGroup(currentRow);
368
369 // Put the row processing into a lambda for avoiding code duplications
370 auto processRow = [&](uint64_t origRowIndex) {
371 // We make two passes. First, we find out the probability to reach an eliminated initial component state
372 ValueType const eliminatedInitialComponentProbability = transitionMatrix.getConstrainedRowSum(origRowIndex, eliminatedInitialComponentStates);
373 // Second, we insert the submatrix entries and find out the target and condition probabilities for this row
374 ValueType targetProbability = storm::utility::zero<ValueType>();
375 ValueType conditionProbability = storm::utility::zero<ValueType>();
376 bool rowSumIsLess1 = false;
377 bool initialStateEntryInserted = false;
378 for (auto const& entry : transitionMatrix.getRow(origRowIndex)) {
379 if (normalForm.terminalStates.get(entry.getColumn())) {
380 STORM_LOG_ASSERT(!storm::utility::isZero(entry.getValue()), "Transition probability must be non-zero");
381 rowSumIsLess1 = true;
382 ValueType const scaledTargetValue = normalForm.getTargetValue(entry.getColumn()) * entry.getValue();
383 targetProbability += scaledTargetValue;
384 if (normalForm.conditionStates.get(entry.getColumn())) {
385 conditionProbability += entry.getValue(); // conditionValue of successor is 1
386 } else {
387 conditionProbability += scaledTargetValue; // for terminal, non-condition states, the condition value equals the target value
388 }
389 } else if (!eliminatedInitialComponentStates.get(entry.getColumn())) {
390 auto const columnIndex = stateToMatrixIndexMap[entry.getColumn()];
391 if (!initialStateEntryInserted && columnIndex >= initialStateInSubmatrix) {
392 if (columnIndex == initialStateInSubmatrix) {
393 matrixBuilder.addNextValue(currentRow, initialStateInSubmatrix, eliminatedInitialComponentProbability + entry.getValue());
394 } else {
395 matrixBuilder.addNextValue(currentRow, initialStateInSubmatrix, eliminatedInitialComponentProbability);
396 matrixBuilder.addNextValue(currentRow, columnIndex, entry.getValue());
397 }
398 initialStateEntryInserted = true;
399 } else {
400 matrixBuilder.addNextValue(currentRow, columnIndex, entry.getValue());
401 }
402 }
403 }
404 if (rowSumIsLess1) {
405 rowsWithSum1.set(currentRow, false);
406 }
407 targetRowValues.push_back(targetProbability);
408 conditionRowValues.push_back(conditionProbability);
409 ++currentRow;
410 };
411 // invoke the lambda
412 if (state == initialState) {
413 for (auto origRowIndex : initialComponentExitRows) {
414 processRow(origRowIndex);
415 }
416 } else {
417 for (auto origRowIndex : transitionMatrix.getRowGroupIndices(state)) {
418 processRow(origRowIndex);
419 }
420 }
421 }
422 submatrix = matrixBuilder.build();
423
424 // eliminate ECs if present. We already checked that the initial state can not yield observation failure, so it cannot be part of an EC.
425 // For all remaining ECs, staying in an EC forever is reflected by collecting a value of zero for both, target and condition
426 storm::storage::BitVector allExceptInit(numSubmatrixRowGroups, true);
427 allExceptInit.set(initialStateInSubmatrix, false);
428 eliminateEndComponents<ValueType>(allExceptInit, true, std::nullopt, submatrix, initialStateInSubmatrix, rowsWithSum1, targetRowValues,
429 conditionRowValues);
430 STORM_LOG_INFO("Processed model has " << submatrix.getRowGroupCount() << " states and " << submatrix.getRowGroupCount() << " choices and "
431 << submatrix.getEntryCount() << " transitions.");
432 }
433
434 SolutionType computeWeightedDiff(storm::Environment const& env, storm::OptimizationDirection const dir, ValueType const& targetWeight,
435 ValueType const& conditionWeight) const {
436 auto rowValues = createScaledVector(targetWeight, targetRowValues, conditionWeight, conditionRowValues);
437
438 // Initialize the solution vector.
439 std::vector<SolutionType> x(submatrix.getRowGroupCount(), storm::utility::zero<ValueType>());
440
441 // Set up the solver.
443 solver->setOptimizationDirection(dir);
444 solver->setRequirementsChecked();
445 solver->setHasUniqueSolution(true);
446 solver->setHasNoEndComponents(true);
447 solver->setLowerBound(-storm::utility::one<ValueType>());
448 solver->setUpperBound(storm::utility::one<ValueType>());
449
450 solver->solveEquations(env, x, rowValues);
451 return x[initialStateInSubmatrix];
452 }
453
455 return initialStateInSubmatrix;
456 }
457
458 void evaluateScheduler(storm::Environment const& env, std::vector<uint64_t>& scheduler, std::vector<SolutionType>& targetResults,
459 std::vector<SolutionType>& conditionResults) const {
460 if (scheduler.empty()) {
461 scheduler.resize(submatrix.getRowGroupCount(), 0);
462 }
463 if (targetResults.empty()) {
464 targetResults.resize(submatrix.getRowGroupCount(), storm::utility::zero<ValueType>());
465 }
466 if (conditionResults.empty()) {
467 conditionResults.resize(submatrix.getRowGroupCount(), storm::utility::zero<ValueType>());
468 }
469 // apply the scheduler
471 bool const convertToEquationSystem = factory.getEquationProblemFormat(env) == storm::solver::LinearEquationSolverProblemFormat::EquationSystem;
472 auto scheduledMatrix = submatrix.selectRowsFromRowGroups(scheduler, convertToEquationSystem);
473 if (convertToEquationSystem) {
474 scheduledMatrix.convertToEquationSystem();
475 }
476 auto solver = factory.create(env, std::move(scheduledMatrix));
477 solver->setBounds(storm::utility::zero<ValueType>(), storm::utility::one<ValueType>());
478 solver->setCachingEnabled(true);
479
480 std::vector<ValueType> subB(submatrix.getRowGroupCount());
481 storm::utility::vector::selectVectorValues<ValueType>(subB, scheduler, submatrix.getRowGroupIndices(), targetRowValues);
482 solver->solveEquations(env, targetResults, subB);
483
484 storm::utility::vector::selectVectorValues<ValueType>(subB, scheduler, submatrix.getRowGroupIndices(), conditionRowValues);
485 solver->solveEquations(env, conditionResults, subB);
486 }
487
488 template<OptimizationDirection Dir>
489 bool improveScheduler(std::vector<uint64_t>& scheduler, ValueType const& lambda, std::vector<SolutionType> const& targetResults,
490 std::vector<SolutionType> const& conditionResults) {
491 bool improved = false;
492 for (uint64_t rowGroupIndex = 0; rowGroupIndex < scheduler.size(); ++rowGroupIndex) {
494 uint64_t optimalRowIndex{0};
495 ValueType scheduledValue;
496 for (auto rowIndex : submatrix.getRowGroupIndices(rowGroupIndex)) {
497 ValueType rowValue = targetRowValues[rowIndex] - lambda * conditionRowValues[rowIndex];
498 for (auto const& entry : submatrix.getRow(rowIndex)) {
499 rowValue += entry.getValue() * (targetResults[entry.getColumn()] - lambda * conditionResults[entry.getColumn()]);
500 }
501 if (rowIndex == scheduler[rowGroupIndex] + submatrix.getRowGroupIndices()[rowGroupIndex]) {
502 scheduledValue = rowValue;
503 }
504 if (groupValue &= rowValue) {
505 optimalRowIndex = rowIndex;
506 }
507 }
508 if (scheduledValue != *groupValue) {
509 scheduler[rowGroupIndex] = optimalRowIndex - submatrix.getRowGroupIndices()[rowGroupIndex];
510 improved = true;
511 }
512 }
513 return improved;
514 }
515
516 private:
517 std::vector<ValueType> createScaledVector(ValueType const& w1, std::vector<ValueType> const& v1, ValueType const& w2,
518 std::vector<ValueType> const& v2) const {
519 STORM_LOG_ASSERT(v1.size() == v2.size(), "Vector sizes must match");
520 std::vector<ValueType> result;
521 result.reserve(v1.size());
522 for (size_t i = 0; i < v1.size(); ++i) {
523 result.push_back(w1 * v1[i] + w2 * v2[i]);
524 }
525 return result;
526 }
527
529 storm::storage::BitVector rowsWithSum1;
530 std::vector<ValueType> targetRowValues;
531 std::vector<ValueType> conditionRowValues;
532 uint64_t initialStateInSubmatrix;
533};
534
536template<typename ValueType, typename SolutionType = ValueType>
537SolutionType computeViaBisection(Environment const& env, BisectionMethodBounds boundOption, uint64_t const initialState,
539 NormalFormData<ValueType> const& normalForm) {
540 // We currently handle sound model checking incorrectly: we would need the actual lower/upper bounds of the weightedReachabilityHelper
542 "Bisection method does not adequately handle propagation of errors. Result is not necessarily sound.");
543 SolutionType const precision = [&env, boundOption]() {
546 "Selected bisection method with exact precision in a setting that might not terminate.");
547 return storm::utility::zero<SolutionType>();
548 } else {
549 return storm::utility::convertNumber<SolutionType>(env.solver().minMax().getPrecision());
550 }
551 }();
552 bool const relative = env.solver().minMax().getRelativeTerminationCriterion();
553
554 WeightedReachabilityHelper wrh(initialState, transitionMatrix, normalForm);
555 SolutionType pMin{storm::utility::zero<SolutionType>()};
556 SolutionType pMax{storm::utility::one<SolutionType>()};
557
558 if (boundOption == BisectionMethodBounds::Advanced) {
559 pMin = wrh.computeWeightedDiff(env, storm::OptimizationDirection::Minimize, storm::utility::zero<ValueType>(), storm::utility::one<ValueType>());
560 pMax = wrh.computeWeightedDiff(env, storm::OptimizationDirection::Maximize, storm::utility::zero<ValueType>(), storm::utility::one<ValueType>());
561 STORM_LOG_TRACE("Conditioning event bounds:\n\t Lower bound: " << storm::utility::convertNumber<double>(pMin)
562 << ",\n\t Upper bound: " << storm::utility::convertNumber<double>(pMax));
563 }
564 storm::utility::Extremum<storm::OptimizationDirection::Maximize, SolutionType> lowerBound = storm::utility::zero<ValueType>();
565 storm::utility::Extremum<storm::OptimizationDirection::Minimize, SolutionType> upperBound = storm::utility::one<ValueType>();
566 SolutionType middle = (*lowerBound + *upperBound) / 2;
567 for (uint64_t iterationCount = 1; true; ++iterationCount) {
568 // evaluate the current middle
569 SolutionType const middleValue = wrh.computeWeightedDiff(env, dir, storm::utility::one<ValueType>(), -middle);
570 // update the bounds and new middle value according to the bisection method
571 if (boundOption == BisectionMethodBounds::Simple) {
572 if (middleValue >= storm::utility::zero<ValueType>()) {
573 lowerBound &= middle;
574 }
575 if (middleValue <= storm::utility::zero<ValueType>()) {
576 upperBound &= middle;
577 }
578 middle = (*lowerBound + *upperBound) / 2; // update middle to the average of the bounds
579 } else {
580 STORM_LOG_ASSERT(boundOption == BisectionMethodBounds::Advanced, "Unknown bisection method bounds");
581 if (middleValue >= storm::utility::zero<ValueType>()) {
582 lowerBound &= middle + (middleValue / pMax);
583 upperBound &= middle + (middleValue / pMin);
584 }
585 if (middleValue <= storm::utility::zero<ValueType>()) {
586 lowerBound &= middle + (middleValue / pMin);
587 upperBound &= middle + (middleValue / pMax);
588 }
589 // update middle to the average of the bounds, but scale it according to the middle value (which is in [-1,1])
590 middle = *lowerBound + (storm::utility::one<SolutionType>() + middleValue) * (*upperBound - *lowerBound) / 2;
591 STORM_LOG_ASSERT(middle >= *lowerBound && middle <= *upperBound, "Bisection method bounds are inconsistent.");
592 }
593 // check for convergence
594 SolutionType const boundDiff = *upperBound - *lowerBound;
595 STORM_LOG_TRACE("Iteration #" << iterationCount << ":\n\t Lower bound: " << storm::utility::convertNumber<double>(*lowerBound)
596 << ",\n\t Upper bound: " << storm::utility::convertNumber<double>(*upperBound)
597 << ",\n\t Difference: " << storm::utility::convertNumber<double>(boundDiff)
598 << ",\n\t Middle val: " << storm::utility::convertNumber<double>(middleValue) << ".");
599 if (boundDiff <= (relative ? (precision * *lowerBound) : precision)) {
600 STORM_LOG_INFO("Bisection method converged after " << iterationCount << " iterations. Difference is "
601 << std::setprecision(std::numeric_limits<double>::digits10)
602 << storm::utility::convertNumber<double>(boundDiff) << ".");
603 break;
604 }
605 // check for early termination
607 STORM_LOG_WARN("Bisection solver aborted after " << iterationCount << "iterations. Bound difference is "
608 << storm::utility::convertNumber<double>(boundDiff) << ".");
609 break;
610 }
611 // process the middle value for the next iteration
613 // find a rational number with a concise representation close to middle and within the bounds
614 auto const exactMiddle = middle;
615 auto numDigits = storm::utility::convertNumber<uint64_t>(
616 storm::utility::floor(storm::utility::log10<SolutionType>(storm::utility::one<SolutionType>() / (*upperBound - *lowerBound))));
617 do {
618 ++numDigits;
619 middle = storm::utility::kwek_mehlhorn::sharpen<SolutionType, SolutionType>(numDigits, exactMiddle);
620 } while (middle <= *lowerBound || middle >= *upperBound);
621 }
622 // Since above code never sets 'middle' to exactly zero or one, we check if that could be necessary after a couple of iterations
623 if (iterationCount == 8) { // 8 is just a heuristic value, it could be any number
624 if (storm::utility::isZero(*lowerBound)) {
625 middle = storm::utility::zero<SolutionType>();
626 } else if (storm::utility::isOne(*upperBound)) {
627 middle = storm::utility::one<SolutionType>();
628 }
629 }
630 }
631 return (*lowerBound + *upperBound) / 2;
632}
633
634template<typename ValueType, typename SolutionType = ValueType>
635SolutionType computeViaPolicyIteration(Environment const& env, uint64_t const initialState, storm::solver::OptimizationDirection const dir,
636 storm::storage::SparseMatrix<ValueType> const& transitionMatrix, NormalFormData<ValueType> const& normalForm) {
637 WeightedReachabilityHelper wrh(initialState, transitionMatrix, normalForm);
638
639 std::vector<uint64_t> scheduler;
640 std::vector<SolutionType> targetResults, conditionResults;
641 for (uint64_t iterationCount = 1; true; ++iterationCount) {
642 wrh.evaluateScheduler(env, scheduler, targetResults, conditionResults);
644 targetResults[wrh.getInternalInitialState()] <= conditionResults[wrh.getInternalInitialState()],
645 "Potential numerical issues: the probability to reach the target is greater than the probability to reach the condition. Difference is "
646 << (storm::utility::convertNumber<double, ValueType>(targetResults[wrh.getInternalInitialState()] -
647 conditionResults[wrh.getInternalInitialState()]))
648 << ".");
649 ValueType const lambda = storm::utility::isZero(conditionResults[wrh.getInternalInitialState()])
650 ? storm::utility::zero<ValueType>()
651 : ValueType(targetResults[wrh.getInternalInitialState()] / conditionResults[wrh.getInternalInitialState()]);
652 bool schedulerChanged{false};
653 if (storm::solver::minimize(dir)) {
654 schedulerChanged = wrh.template improveScheduler<storm::OptimizationDirection::Minimize>(scheduler, lambda, targetResults, conditionResults);
655 } else {
656 schedulerChanged = wrh.template improveScheduler<storm::OptimizationDirection::Maximize>(scheduler, lambda, targetResults, conditionResults);
657 }
658 if (!schedulerChanged) {
659 STORM_LOG_INFO("Policy iteration for conditional probabilities converged after " << iterationCount << " iterations.");
660 return lambda;
661 }
663 STORM_LOG_WARN("Policy iteration for conditional probabilities converged aborted after " << iterationCount << "iterations.");
664 return lambda;
665 }
666 }
667}
668
669template<typename ValueType, typename SolutionType = ValueType>
670std::optional<SolutionType> handleTrivialCases(uint64_t const initialState, NormalFormData<ValueType> const& normalForm) {
671 if (normalForm.terminalStates.get(initialState)) {
672 STORM_LOG_DEBUG("Initial state is terminal.");
673 if (normalForm.conditionStates.get(initialState)) {
674 return normalForm.getTargetValue(initialState); // The value is already known, nothing to do.
675 } else {
676 STORM_LOG_THROW(!normalForm.universalObservationFailureStates.get(initialState), storm::exceptions::NotSupportedException,
677 "Trying to compute undefined conditional probability: the condition has probability 0 under all policies.");
678 // The last case for a terminal initial state is that it is already target and the condition is reachable with non-zero probability.
679 // In this case, all schedulers induce a conditional probability of 1 (or do not reach the condition, i.e., have undefined value)
680 return storm::utility::one<SolutionType>();
681 }
682 } else {
683 // Catch the case where all terminal states have value zero
684 if (normalForm.nonZeroTargetStateValues.empty()) {
685 return storm::utility::zero<SolutionType>();
686 };
687 }
688 return std::nullopt; // No trivial case applies, we need to compute the value.
689}
690
691} // namespace internal
692
693template<typename ValueType, typename SolutionType>
695 storm::storage::SparseMatrix<ValueType> const& transitionMatrix,
696 storm::storage::SparseMatrix<ValueType> const& backwardTransitions,
697 storm::storage::BitVector const& targetStates, storm::storage::BitVector const& conditionStates) {
698 // We might require adapting the precision of the solver to counter error propagation (e.g. when computing the normal form).
699 auto normalFormConstructionEnv = env;
700 auto analysisEnv = env;
701 if (env.solver().isForceSoundness()) {
702 // We intuitively have to divide the precision into two parts, one for computations when constructing the normal form and one for the actual analysis.
703 // As the former is usually less numerically challenging, we use a factor of 1/10 for the normal form construction and 9/10 for the analysis.
704 auto const normalFormPrecisionFactor = storm::utility::convertNumber<storm::RationalNumber, std::string>("1/10");
705 normalFormConstructionEnv.solver().minMax().setPrecision(env.solver().minMax().getPrecision() * normalFormPrecisionFactor);
706 analysisEnv.solver().minMax().setPrecision(env.solver().minMax().getPrecision() *
707 (storm::utility::one<storm::RationalNumber>() - normalFormPrecisionFactor));
708 }
709
710 // We first translate the problem into a normal form.
711 // @see doi.org/10.1007/978-3-642-54862-8_43
712 STORM_LOG_THROW(goal.hasRelevantValues(), storm::exceptions::NotSupportedException,
713 "No initial state given. Conditional probabilities can only be computed for models with a single initial state.");
714 STORM_LOG_THROW(goal.relevantValues().hasUniqueSetBit(), storm::exceptions::NotSupportedException,
715 "Only one initial state is supported for conditional probabilities");
716 STORM_LOG_TRACE("Computing conditional probabilities for a model with " << transitionMatrix.getRowGroupCount() << " states and "
717 << transitionMatrix.getEntryCount() << " transitions.");
718 // storm::utility::Stopwatch sw(true);
719 auto normalFormData = internal::obtainNormalForm(normalFormConstructionEnv, goal.direction(), transitionMatrix, backwardTransitions, goal.relevantValues(),
720 targetStates, conditionStates);
721 // sw.stop();
722 // STORM_PRINT_AND_LOG("Time for obtaining the normal form: " << sw << ".\n");
723 // Then, we solve the induced problem using the selected algorithm
724 auto const initialState = *goal.relevantValues().begin();
725 ValueType initialStateValue = -storm::utility::one<ValueType>();
726 if (auto trivialValue = internal::handleTrivialCases<ValueType, SolutionType>(initialState, normalFormData); trivialValue.has_value()) {
727 initialStateValue = *trivialValue;
728 STORM_LOG_DEBUG("Initial state has trivial value " << initialStateValue);
729 } else {
730 STORM_LOG_ASSERT(normalFormData.maybeStates.get(initialState), "Initial state must be a maybe state if it is not a terminal state");
731 auto alg = analysisEnv.modelchecker().getConditionalAlgorithmSetting();
734 }
735 STORM_LOG_INFO("Analyzing normal form with " << normalFormData.maybeStates.getNumberOfSetBits() << " maybe states using algorithm '" << alg << ".");
736 // sw.restart();
737 switch (alg) {
739 initialStateValue = internal::computeViaRestartMethod(analysisEnv, initialState, goal.direction(), transitionMatrix, normalFormData);
740 break;
742 initialStateValue = internal::computeViaBisection(analysisEnv, internal::BisectionMethodBounds::Simple, initialState, goal.direction(),
743 transitionMatrix, normalFormData);
744 break;
746 initialStateValue = internal::computeViaBisection(analysisEnv, internal::BisectionMethodBounds::Advanced, initialState, goal.direction(),
747 transitionMatrix, normalFormData);
748 break;
750 initialStateValue = internal::computeViaPolicyIteration(analysisEnv, initialState, goal.direction(), transitionMatrix, normalFormData);
751 break;
752 default:
753 STORM_LOG_THROW(false, storm::exceptions::NotImplementedException, "Unknown conditional probability algorithm: " << alg);
754 }
755 // sw.stop();
756 // STORM_PRINT_AND_LOG("Time for analyzing the normal form: " << sw << ".\n");
757 }
758 return std::unique_ptr<CheckResult>(new ExplicitQuantitativeCheckResult<ValueType>(initialState, initialStateValue));
759}
760
761template std::unique_ptr<CheckResult> computeConditionalProbabilities(Environment const& env, storm::solver::SolveGoal<double>&& goal,
762 storm::storage::SparseMatrix<double> const& transitionMatrix,
763 storm::storage::SparseMatrix<double> const& backwardTransitions,
764 storm::storage::BitVector const& targetStates,
765 storm::storage::BitVector const& conditionStates);
766
770 storm::storage::BitVector const& targetStates,
771 storm::storage::BitVector const& conditionStates);
772
773} // namespace storm::modelchecker
SolverEnvironment & solver()
storm::RationalNumber const & getPrecision() const
bool const & getRelativeTerminationCriterion() const
Helper class that optionally holds a reference to an object of type T.
Definition OptionalRef.h:48
MinMaxSolverEnvironment & minMax()
static MDPSparseModelCheckingHelperReturnType< SolutionType > computeUntilProbabilities(Environment const &env, storm::solver::SolveGoal< ValueType, SolutionType > &&goal, storm::storage::SparseMatrix< ValueType > const &transitionMatrix, storm::storage::SparseMatrix< ValueType > const &backwardTransitions, storm::storage::BitVector const &phiStates, storm::storage::BitVector const &psiStates, bool qualitative, bool produceScheduler, ModelCheckerHint const &hint=ModelCheckerHint())
A helper class that computes (weighted) reachability probabilities for a given MDP in normal form.
bool improveScheduler(std::vector< uint64_t > &scheduler, ValueType const &lambda, std::vector< SolutionType > const &targetResults, std::vector< SolutionType > const &conditionResults)
SolutionType computeWeightedDiff(storm::Environment const &env, storm::OptimizationDirection const dir, ValueType const &targetWeight, ValueType const &conditionWeight) const
WeightedReachabilityHelper(uint64_t const initialState, storm::storage::SparseMatrix< ValueType > const &transitionMatrix, NormalFormData< ValueType > const &normalForm)
void evaluateScheduler(storm::Environment const &env, std::vector< uint64_t > &scheduler, std::vector< SolutionType > &targetResults, std::vector< SolutionType > &conditionResults) const
virtual std::unique_ptr< LinearEquationSolver< ValueType > > create(Environment const &env) const override
Creates an equation solver with the current settings, but without a matrix.
virtual std::unique_ptr< MinMaxLinearEquationSolver< ValueType, SolutionType > > create(Environment const &env) const override
virtual LinearEquationSolverProblemFormat getEquationProblemFormat(Environment const &env) const
Retrieves the problem format that the solver expects if it was created with the current settings.
A bit vector that is internally represented as a vector of 64-bit values.
Definition BitVector.h:18
bool empty() const
Retrieves whether no bits are set to true in this bit vector.
const_iterator begin() const
Returns an iterator to the indices of the set bits in the bit vector.
void set(uint_fast64_t index, bool value=true)
Sets the given truth value at the given index.
uint_fast64_t getNumberOfSetBits() const
Returns the number of bits that are set to true in this bit vector.
std::vector< uint_fast64_t > getNumberOfSetBitsBeforeIndices() const
Retrieves a vector that holds at position i the number of bits set before index i.
bool get(uint_fast64_t index) const
Retrieves the truth value of the bit at the given index and performs a bound check.
This class represents the decomposition of a nondeterministic model into its maximal end components.
A class that can be used to build a sparse matrix by adding value by value.
void addNextValue(index_type row, index_type column, value_type const &value)
Sets the matrix entry at the given row and column to the given value.
void newRowGroup(index_type startingRow)
Starts a new row group in the matrix.
SparseMatrix< value_type > build(index_type overriddenRowCount=0, index_type overriddenColumnCount=0, index_type overriddenRowGroupCount=0)
A class that holds a possibly non-square matrix in the compressed row storage format.
const_rows getRow(index_type row) const
Returns an object representing the given row.
index_type getEntryCount() const
Returns the number of entries in the matrix.
index_type getNumRowsInRowGroups(storm::storage::BitVector const &groupConstraint) const
Returns the total number of rows that are in one of the specified row groups.
SparseMatrix getSubmatrix(bool useGroups, storm::storage::BitVector const &rowConstraint, storm::storage::BitVector const &columnConstraint, bool insertDiagonalEntries=false, storm::storage::BitVector const &makeZeroColumns=storm::storage::BitVector()) const
Creates a submatrix of the current matrix by dropping all rows and columns whose bits are not set to ...
index_type getRowGroupCount() const
Returns the number of row groups in the matrix.
value_type getConstrainedRowSum(index_type row, storm::storage::BitVector const &columns) const
Sums the entries in the given row and columns.
std::vector< index_type > const & getRowGroupIndices() const
Returns the grouping of rows of this matrix.
storm::storage::SparseMatrix< value_type > transpose(bool joinGroups=false, bool keepZeros=false) const
Transposes the matrix.
index_type getRowCount() const
Returns the number of rows of the matrix.
static EndComponentEliminatorReturnType transform(storm::storage::SparseMatrix< ValueType > const &originalMatrix, storm::storage::MaximalEndComponentDecomposition< ValueType > ecs, storm::storage::BitVector const &subsystemStates, storm::storage::BitVector const &addSinkRowStates, bool addSelfLoopAtSinkStates=false)
Stores and manages an extremal (maximal or minimal) value.
Definition Extremum.h:15
#define STORM_LOG_INFO(message)
Definition logging.h:29
#define STORM_LOG_WARN(message)
Definition logging.h:30
#define STORM_LOG_DEBUG(message)
Definition logging.h:23
#define STORM_LOG_TRACE(message)
Definition logging.h:17
#define STORM_LOG_ASSERT(cond, message)
Definition macros.h:11
#define STORM_LOG_WARN_COND(cond, message)
Definition macros.h:38
#define STORM_LOG_THROW(cond, exception, message)
Definition macros.h:30
void computeReachabilityProbabilities(Environment const &env, std::map< uint64_t, ValueType > &nonZeroResults, storm::solver::OptimizationDirection const dir, storm::storage::SparseMatrix< ValueType > const &transitionMatrix, storm::storage::BitVector const &initialStates, storm::storage::BitVector const &allowedStates, storm::storage::BitVector const &targetStates)
Computes the reachability probabilities for the given target states and inserts all non-zero values i...
void eliminateEndComponents(storm::storage::BitVector possibleEcStates, bool addRowAtRepresentativeState, std::optional< uint64_t > representativeRowEntry, storm::storage::SparseMatrix< ValueType > &matrix, uint64_t &initialState, storm::storage::BitVector &rowsWithSum1, std::vector< ValueType > &rowValues1, storm::OptionalRef< std::vector< ValueType > > rowValues2={})
SolutionType computeViaPolicyIteration(Environment const &env, uint64_t const initialState, storm::solver::OptimizationDirection const dir, storm::storage::SparseMatrix< ValueType > const &transitionMatrix, NormalFormData< ValueType > const &normalForm)
NormalFormData< ValueType > obtainNormalForm(Environment const &env, storm::solver::OptimizationDirection const dir, storm::storage::SparseMatrix< ValueType > const &transitionMatrix, storm::storage::SparseMatrix< ValueType > const &backwardTransitions, storm::storage::BitVector const &relevantStates, storm::storage::BitVector const &targetStates, storm::storage::BitVector const &conditionStates)
std::optional< SolutionType > handleTrivialCases(uint64_t const initialState, NormalFormData< ValueType > const &normalForm)
SolutionType solveMinMaxEquationSystem(storm::Environment const &env, storm::storage::SparseMatrix< ValueType > const &matrix, std::vector< ValueType > const &rowValues, storm::storage::BitVector const &rowsWithSum1, storm::solver::OptimizationDirection const dir, uint64_t const initialState)
SolutionType computeViaBisection(Environment const &env, BisectionMethodBounds boundOption, uint64_t const initialState, storm::solver::OptimizationDirection const dir, storm::storage::SparseMatrix< ValueType > const &transitionMatrix, NormalFormData< ValueType > const &normalForm)
SolutionType computeViaRestartMethod(Environment const &env, uint64_t const initialState, storm::solver::OptimizationDirection const dir, storm::storage::SparseMatrix< ValueType > const &transitionMatrix, NormalFormData< ValueType > const &normalForm)
Uses the restart method by Baier et al.
std::unique_ptr< CheckResult > computeConditionalProbabilities(Environment const &env, storm::solver::SolveGoal< ValueType, SolutionType > &&goal, storm::storage::SparseMatrix< ValueType > const &transitionMatrix, storm::storage::SparseMatrix< ValueType > const &backwardTransitions, storm::storage::BitVector const &targetStates, storm::storage::BitVector const &conditionStates)
bool constexpr minimize(OptimizationDirection d)
storm::storage::BitVector getReachableStates(storm::storage::SparseMatrix< T > const &transitionMatrix, storm::storage::BitVector const &initialStates, storm::storage::BitVector const &constraintStates, storm::storage::BitVector const &targetStates, bool useStepBound, uint_fast64_t maximalSteps, boost::optional< storm::storage::BitVector > const &choiceFilter)
Performs a forward depth-first search through the underlying graph structure to identify the states t...
Definition graph.cpp:48
storm::storage::BitVector performProb1A(storm::models::sparse::NondeterministicModel< T, RM > const &model, storm::storage::SparseMatrix< T > const &backwardTransitions, storm::storage::BitVector const &phiStates, storm::storage::BitVector const &psiStates)
Computes the sets of states that have probability 1 of satisfying phi until psi under all possible re...
Definition graph.cpp:988
storm::storage::BitVector performProb0A(storm::storage::SparseMatrix< T > const &backwardTransitions, storm::storage::BitVector const &phiStates, storm::storage::BitVector const &psiStates)
Definition graph.cpp:740
storm::storage::BitVector performProb0E(storm::models::sparse::NondeterministicModel< T, RM > const &model, storm::storage::SparseMatrix< T > const &backwardTransitions, storm::storage::BitVector const &phiStates, storm::storage::BitVector const &psiStates)
Computes the sets of states that have probability 0 of satisfying phi until psi under at least one po...
Definition graph.cpp:967
bool isTerminate()
Check whether the program should terminate (due to some abort signal).
bool isOne(ValueType const &a)
Definition constants.cpp:36
ValueType floor(ValueType const &number)
bool isZero(ValueType const &a)
Definition constants.cpp:41
storm::storage::BitVector const existentialObservationFailureStates
ValueType failProbability(uint64_t state) const
std::map< uint64_t, ValueType > const nonZeroTargetStateValues
storm::storage::BitVector const conditionStates
storm::storage::BitVector const universalObservationFailureStates