Merge pull request 'Shields for MDPs' (#24) from mdp_shields into main

Reviewed-on: https://git.pranger.xyz/TEMPEST/tempest-devel/pulls/24

src/storm/modelchecker/helper/finitehorizon/SparseNondeterministicStepBoundedHorizonHelper.cpp

@@ -17,6 +17,30 @@ namespace storm {
     namespace modelchecker {
         namespace helper {
 
+            template<typename ValueType>
+            void SparseNondeterministicStepBoundedHorizonHelper<ValueType>::getMaybeStatesRowGroupSizes(storm::storage::SparseMatrix<ValueType> const& transitionMatrix, storm::storage::BitVector const& maybeStates, std::vector<uint64_t>& maybeStatesRowGroupSizes, uint& choiceValuesCounter) {
+                std::vector<uint64_t> rowGroupIndices = transitionMatrix.getRowGroupIndices();
+                for(uint counter = 0; counter < maybeStates.size(); counter++) {
+                    if(maybeStates.get(counter)) {
+                        maybeStatesRowGroupSizes.push_back(rowGroupIndices.at(counter));
+                        choiceValuesCounter += transitionMatrix.getRowGroupSize(counter);
+                    }
+                }
+            }
+
+            template<typename ValueType>
+            void SparseNondeterministicStepBoundedHorizonHelper<ValueType>::getChoiceValues(storm::storage::SparseMatrix<ValueType> const& transitionMatrix, storm::storage::BitVector const& maybeStates, std::vector<ValueType> const& choiceValues, std::vector<ValueType>& allChoices) {
+                auto choice_it = choiceValues.begin();
+                for(uint state = 0; state < transitionMatrix.getRowGroupIndices().size() - 1; state++) {
+                    uint rowGroupSize = transitionMatrix.getRowGroupSize(state);
+                    if (maybeStates.get(state)) {
+                        for(uint choice = 0; choice < rowGroupSize; choice++, choice_it++) {
+                            allChoices.at(transitionMatrix.getRowGroupIndices().at(state) + choice) = *choice_it;
+                        }
+                    }
+                }
+            }
+
             template<typename ValueType>
             SparseNondeterministicStepBoundedHorizonHelper<ValueType>::SparseNondeterministicStepBoundedHorizonHelper(/*storm::storage::SparseMatrix<ValueType> const& transitionMatrix, storm::storage::SparseMatrix<ValueType> const& backwardTransitions*/)
             //transitionMatrix(transitionMatrix), backwardTransitions(backwardTransitions)
@@ -25,7 +49,7 @@ namespace storm {
             }
 
             template<typename ValueType>
-            std::vector<ValueType> SparseNondeterministicStepBoundedHorizonHelper<ValueType>::compute(Environment const& env, storm::solver::SolveGoal<ValueType>&& goal, storm::storage::SparseMatrix<ValueType> const& transitionMatrix, storm::storage::SparseMatrix<ValueType> const& backwardTransitions, storm::storage::BitVector const& phiStates, storm::storage::BitVector const& psiStates,   uint64_t lowerBound, uint64_t upperBound,  ModelCheckerHint const& hint)
+            std::vector<ValueType> SparseNondeterministicStepBoundedHorizonHelper<ValueType>::compute(Environment const& env, storm::solver::SolveGoal<ValueType>&& goal, storm::storage::SparseMatrix<ValueType> const& transitionMatrix, storm::storage::SparseMatrix<ValueType> const& backwardTransitions, storm::storage::BitVector const& phiStates, storm::storage::BitVector const& psiStates, uint64_t lowerBound, uint64_t upperBound, ModelCheckerHint const& hint, storm::storage::BitVector& resultMaybeStates, std::vector<ValueType>& choiceValues)
             {
                 std::vector<ValueType> result(transitionMatrix.getRowGroupCount(), storm::utility::zero<ValueType>());
                 storm::storage::BitVector makeZeroColumns;
@@ -60,8 +84,38 @@ namespace storm {
                     std::vector<ValueType> subresult(maybeStates.getNumberOfSetBits());
 
                     auto multiplier = storm::solver::MultiplierFactory<ValueType>().create(env, submatrix);
+
+                    std::vector<uint64_t> rowGroupIndices = transitionMatrix.getRowGroupIndices();
+                    std::vector<uint64_t> maybeStatesRowGroupSizes;
+                    uint choiceValuesCounter;
+                    getMaybeStatesRowGroupSizes(transitionMatrix, maybeStates, maybeStatesRowGroupSizes, choiceValuesCounter);
+                    choiceValues = std::vector<ValueType>(choiceValuesCounter, storm::utility::zero<ValueType>());
+
                     if (lowerBound == 0) {
+                        if(goal.isShieldingTask())
+                        {
+                            multiplier->repeatedMultiplyAndReduceWithChoices(env, goal.direction(), subresult, &b, upperBound, nullptr, choiceValues, maybeStatesRowGroupSizes);
+
+                            // fill up choicesValues for shields
+                            std::vector<ValueType> allChoices = std::vector<ValueType>(transitionMatrix.getRowCount(), storm::utility::zero<ValueType>());
+                            getChoiceValues(transitionMatrix, maybeStates, choiceValues, allChoices);
+                            choiceValues = allChoices;
+                        } else {
                             multiplier->repeatedMultiplyAndReduce(env, goal.direction(), subresult, &b, upperBound);
+                        }
+                    } else {
+                        if(goal.isShieldingTask())
+                        {
+                            multiplier->repeatedMultiplyAndReduceWithChoices(env, goal.direction(), subresult, &b, upperBound - lowerBound + 1, nullptr, choiceValues, maybeStatesRowGroupSizes);
+                            storm::storage::SparseMatrix<ValueType> submatrix = transitionMatrix.getSubmatrix(true, maybeStates, maybeStates, false);
+                            auto multiplier = storm::solver::MultiplierFactory<ValueType>().create(env, submatrix);
+                            b = std::vector<ValueType>(b.size(), storm::utility::zero<ValueType>());
+                            multiplier->repeatedMultiplyAndReduceWithChoices(env, goal.direction(), subresult, &b, lowerBound - 1, nullptr, choiceValues, maybeStatesRowGroupSizes);
+
+                            // fill up choicesValues for shields
+                            std::vector<ValueType> allChoices = std::vector<ValueType>(transitionMatrix.getRowCount(), storm::utility::zero<ValueType>());
+                            getChoiceValues(transitionMatrix, maybeStates, choiceValues, allChoices);
+                            choiceValues = allChoices;
                         } else {
                             multiplier->repeatedMultiplyAndReduce(env, goal.direction(), subresult, &b, upperBound - lowerBound + 1);
                             storm::storage::SparseMatrix<ValueType> submatrix = transitionMatrix.getSubmatrix(true, maybeStates, maybeStates, false);
@@ -69,16 +123,18 @@ namespace storm {
                             b = std::vector<ValueType>(b.size(), storm::utility::zero<ValueType>());
                             multiplier->repeatedMultiplyAndReduce(env, goal.direction(), subresult, &b, lowerBound - 1);
                         }
+                    }
                     // Set the values of the resulting vector accordingly.
                     storm::utility::vector::setVectorValues(result, maybeStates, subresult);
                 }
                 if (lowerBound == 0) {
                     storm::utility::vector::setVectorValues(result, psiStates, storm::utility::one<ValueType>());
                 }
+
+                resultMaybeStates = maybeStates;
                 return result;
             }
 
-
             template class SparseNondeterministicStepBoundedHorizonHelper<double>;
             template class SparseNondeterministicStepBoundedHorizonHelper<storm::RationalNumber>;
         }

src/storm/modelchecker/helper/finitehorizon/SparseNondeterministicStepBoundedHorizonHelper.h

@@ -4,6 +4,7 @@
 #include "storm/modelchecker/hints/ModelCheckerHint.h"
 #include "storm/modelchecker/prctl/helper/SolutionType.h"
 #include "storm/storage/SparseMatrix.h"
+#include "storm/storage/BitVector.h"
 #include "storm/utility/solver.h"
 #include "storm/solver/SolveGoal.h"
 
@@ -15,7 +16,10 @@ namespace storm {
             class SparseNondeterministicStepBoundedHorizonHelper {
             public:
                 SparseNondeterministicStepBoundedHorizonHelper();
-                std::vector<ValueType> compute(Environment const& env, storm::solver::SolveGoal<ValueType>&& goal, storm::storage::SparseMatrix<ValueType> const& transitionMatrix, storm::storage::SparseMatrix<ValueType> const& backwardTransitions,  storm::storage::BitVector const& phiStates, storm::storage::BitVector const& psiStates, uint64_t lowerBound, uint64_t upperBound, ModelCheckerHint const& hint = ModelCheckerHint());
+
+                void getMaybeStatesRowGroupSizes(storm::storage::SparseMatrix<ValueType> const& transitionMatrix, storm::storage::BitVector const& maybeStates, std::vector<uint64_t>& maybeStatesRowGroupSizes, uint& choiceValuesCounter);
+                void getChoiceValues(storm::storage::SparseMatrix<ValueType> const& transitionMatrix, storm::storage::BitVector const& maybeStates, std::vector<ValueType> const& choiceValues, std::vector<ValueType>& allChoices);
+                std::vector<ValueType> compute(Environment const& env, storm::solver::SolveGoal<ValueType>&& goal, storm::storage::SparseMatrix<ValueType> const& transitionMatrix, storm::storage::SparseMatrix<ValueType> const& backwardTransitions, storm::storage::BitVector const& phiStates, storm::storage::BitVector const& psiStates, uint64_t lowerBound, uint64_t upperBound, ModelCheckerHint const& hint = ModelCheckerHint(), storm::storage::BitVector& resultMaybeStates = nullptr, std::vector<ValueType>& choiceValues = nullptr);
             private:
             };
 

src/storm/modelchecker/prctl/SparseMdpPrctlModelChecker.cpp

@@ -24,6 +24,10 @@
 
 #include "storm/solver/SolveGoal.h"
 
+#include "storm/storage/BitVector.h"
+
+#include "storm/shields/ShieldHandling.h"
+
 #include "storm/settings/modules/GeneralSettings.h"
 
 #include "storm/exceptions/InvalidStateException.h"
@@ -92,7 +96,16 @@ namespace storm {
                 ExplicitQualitativeCheckResult const& leftResult = leftResultPointer->asExplicitQualitativeCheckResult();
                 ExplicitQualitativeCheckResult const& rightResult = rightResultPointer->asExplicitQualitativeCheckResult();
                 storm::modelchecker::helper::SparseNondeterministicStepBoundedHorizonHelper<ValueType> helper;
-                std::vector<ValueType> numericResult = helper.compute(env, storm::solver::SolveGoal<ValueType>(this->getModel(), checkTask), this->getModel().getTransitionMatrix(), this->getModel().getBackwardTransitions(), leftResult.getTruthValuesVector(), rightResult.getTruthValuesVector(), pathFormula.getNonStrictLowerBound<uint64_t>(), pathFormula.getNonStrictUpperBound<uint64_t>(), checkTask.getHint());
+                std::vector<ValueType> numericResult;
+
+                //This works only with empty vectors, no nullptr
+                storm::storage::BitVector resultMaybeStates;
+                std::vector<ValueType> choiceValues;
+
+                numericResult = helper.compute(env, storm::solver::SolveGoal<ValueType>(this->getModel(), checkTask), this->getModel().getTransitionMatrix(), this->getModel().getBackwardTransitions(), leftResult.getTruthValuesVector(), rightResult.getTruthValuesVector(), pathFormula.getNonStrictLowerBound<uint64_t>(), pathFormula.getNonStrictUpperBound<uint64_t>(), checkTask.getHint(), resultMaybeStates, choiceValues);
+                if(checkTask.isShieldingTask()) {
+                    tempest::shields::createShield<ValueType>(std::make_shared<storm::models::sparse::Mdp<ValueType>>(this->getModel()), std::move(choiceValues), checkTask.getShieldingExpression(), checkTask.getOptimizationDirection(), std::move(resultMaybeStates), storm::storage::BitVector(resultMaybeStates.size(), true));
+                }
                 return std::unique_ptr<CheckResult>(new ExplicitQuantitativeCheckResult<ValueType>(std::move(numericResult)));
             }
         }
@@ -103,8 +116,14 @@ namespace storm {
             STORM_LOG_THROW(checkTask.isOptimizationDirectionSet(), storm::exceptions::InvalidPropertyException, "Formula needs to specify whether minimal or maximal values are to be computed on nondeterministic model.");
             std::unique_ptr<CheckResult> subResultPointer = this->check(env, pathFormula.getSubformula());
             ExplicitQualitativeCheckResult const& subResult = subResultPointer->asExplicitQualitativeCheckResult();
-            std::vector<ValueType> numericResult = storm::modelchecker::helper::SparseMdpPrctlHelper<ValueType>::computeNextProbabilities(env, checkTask.getOptimizationDirection(), this->getModel().getTransitionMatrix(), subResult.getTruthValuesVector());
-            return std::unique_ptr<CheckResult>(new ExplicitQuantitativeCheckResult<ValueType>(std::move(numericResult)));
+            auto ret = storm::modelchecker::helper::SparseMdpPrctlHelper<ValueType>::computeNextProbabilities(env, storm::solver::SolveGoal<ValueType>(this->getModel(), checkTask), checkTask.getOptimizationDirection(), this->getModel().getTransitionMatrix(), subResult.getTruthValuesVector());
+            std::unique_ptr<CheckResult> result(new ExplicitQuantitativeCheckResult<ValueType>(std::move(ret.values)));
+            if(checkTask.isShieldingTask()) {
+                tempest::shields::createShield<ValueType>(std::make_shared<storm::models::sparse::Mdp<ValueType>>(this->getModel()), std::move(ret.choiceValues), checkTask.getShieldingExpression(), checkTask.getOptimizationDirection(), std::move(ret.maybeStates), storm::storage::BitVector(ret.maybeStates.size(), true));
+            } else if (checkTask.isProduceSchedulersSet() && ret.scheduler) {
+                result->asExplicitQuantitativeCheckResult<ValueType>().setScheduler(std::move(ret.scheduler));
+            }
+            return result;
         }
 
         template<typename SparseMdpModelType>
@@ -117,7 +136,9 @@ namespace storm {
             ExplicitQualitativeCheckResult const& rightResult = rightResultPointer->asExplicitQualitativeCheckResult();
             auto ret = storm::modelchecker::helper::SparseMdpPrctlHelper<ValueType>::computeUntilProbabilities(env, storm::solver::SolveGoal<ValueType>(this->getModel(), checkTask), this->getModel().getTransitionMatrix(), this->getModel().getBackwardTransitions(), leftResult.getTruthValuesVector(), rightResult.getTruthValuesVector(), checkTask.isQualitativeSet(), checkTask.isProduceSchedulersSet(), checkTask.getHint());
             std::unique_ptr<CheckResult> result(new ExplicitQuantitativeCheckResult<ValueType>(std::move(ret.values)));
-            if (checkTask.isProduceSchedulersSet() && ret.scheduler) {
+            if(checkTask.isShieldingTask()) {
+                tempest::shields::createShield<ValueType>(std::make_shared<storm::models::sparse::Mdp<ValueType>>(this->getModel()), std::move(ret.choiceValues), checkTask.getShieldingExpression(), checkTask.getOptimizationDirection(), std::move(ret.maybeStates), storm::storage::BitVector(ret.maybeStates.size(), true));
+            } else if (checkTask.isProduceSchedulersSet() && ret.scheduler) {
                 result->asExplicitQuantitativeCheckResult<ValueType>().setScheduler(std::move(ret.scheduler));
             }
             return result;
@@ -131,7 +152,9 @@ namespace storm {
             ExplicitQualitativeCheckResult const& subResult = subResultPointer->asExplicitQualitativeCheckResult();
             auto ret = storm::modelchecker::helper::SparseMdpPrctlHelper<ValueType>::computeGloballyProbabilities(env, storm::solver::SolveGoal<ValueType>(this->getModel(), checkTask), this->getModel().getTransitionMatrix(), this->getModel().getBackwardTransitions(), subResult.getTruthValuesVector(), checkTask.isQualitativeSet(), checkTask.isProduceSchedulersSet());
             std::unique_ptr<CheckResult> result(new ExplicitQuantitativeCheckResult<ValueType>(std::move(ret.values)));
-            if (checkTask.isProduceSchedulersSet() && ret.scheduler) {
+            if(checkTask.isShieldingTask()) {
+                tempest::shields::createShield<ValueType>(std::make_shared<storm::models::sparse::Mdp<ValueType>>(this->getModel()), std::move(ret.choiceValues), checkTask.getShieldingExpression(), checkTask.getOptimizationDirection(),subResult.getTruthValuesVector(), storm::storage::BitVector(ret.maybeStates.size(), true));
+            } else if (checkTask.isProduceSchedulersSet() && ret.scheduler) {
                 result->asExplicitQuantitativeCheckResult<ValueType>().setScheduler(std::move(ret.scheduler));
             }
             return result;

src/storm/modelchecker/prctl/helper/MDPModelCheckingHelperReturnType.h

@@ -18,7 +18,7 @@ namespace storm {
                 MDPSparseModelCheckingHelperReturnType(MDPSparseModelCheckingHelperReturnType const&) = delete;
                 MDPSparseModelCheckingHelperReturnType(MDPSparseModelCheckingHelperReturnType&&) = default;
 
-                MDPSparseModelCheckingHelperReturnType(std::vector<ValueType>&& values, std::unique_ptr<storm::storage::Scheduler<ValueType>>&& scheduler = nullptr) : values(std::move(values)), scheduler(std::move(scheduler)) {
+                MDPSparseModelCheckingHelperReturnType(std::vector<ValueType>&& values, storm::storage::BitVector&& maybeStates = nullptr, std::unique_ptr<storm::storage::Scheduler<ValueType>>&& scheduler = nullptr, std::vector<ValueType>&& choiceValues = nullptr) : values(std::move(values)), maybeStates(maybeStates), scheduler(std::move(scheduler)), choiceValues(std::move(choiceValues)) {
                     // Intentionally left empty.
                 }
 
@@ -29,8 +29,14 @@ namespace storm {
                 // The values computed for the states.
                 std::vector<ValueType> values;
 
+                // The maybe states of the model
+                storm::storage::BitVector maybeStates;
+
                 // A scheduler, if it was computed.
                 std::unique_ptr<storm::storage::Scheduler<ValueType>> scheduler;
+
+                // The values computed for the available choices.
+                std::vector<ValueType> choiceValues;
             };
         }
 

src/storm/modelchecker/prctl/helper/SparseMdpPrctlHelper.cpp

@@ -134,16 +134,25 @@ namespace storm {
             }
 
             template<typename ValueType>
-            std::vector<ValueType> SparseMdpPrctlHelper<ValueType>::computeNextProbabilities(Environment const& env, OptimizationDirection dir, storm::storage::SparseMatrix<ValueType> const& transitionMatrix, storm::storage::BitVector const& nextStates) {
+            MDPSparseModelCheckingHelperReturnType<ValueType> SparseMdpPrctlHelper<ValueType>::computeNextProbabilities(Environment const& env, storm::solver::SolveGoal<ValueType>&& goal, OptimizationDirection dir, storm::storage::SparseMatrix<ValueType> const& transitionMatrix, storm::storage::BitVector const& nextStates) {
 
                 // Create the vector with which to multiply and initialize it correctly.
                 std::vector<ValueType> result(transitionMatrix.getRowGroupCount());
                 storm::utility::vector::setVectorValues(result, nextStates, storm::utility::one<ValueType>());
+                std::vector<ValueType> choiceValues = std::vector<ValueType>(transitionMatrix.getRowCount(), storm::utility::zero<ValueType>());
+                storm::storage::BitVector allStates = storm::storage::BitVector(transitionMatrix.getRowGroupCount(), true);
 
                 auto multiplier = storm::solver::MultiplierFactory<ValueType>().create(env, transitionMatrix);
+
+                if(goal.isShieldingTask()) {
+                    multiplier->multiply(env, result, nullptr, choiceValues);
+                    multiplier->reduce(env, goal.direction(), choiceValues, transitionMatrix.getRowGroupIndices(), result, nullptr);
+                }
+                else {
                     multiplier->multiplyAndReduce(env, dir, result, nullptr, result);
+                }
 
-                return result;
+                return MDPSparseModelCheckingHelperReturnType<ValueType>(std::move(result), std::move(allStates), nullptr, std::move(choiceValues));
             }
             
             template<typename ValueType>
@@ -609,8 +618,20 @@ namespace storm {
                 // Check if the values of the maybe states are relevant for the SolveGoal
                 bool maybeStatesNotRelevant = goal.hasRelevantValues() && goal.relevantValues().isDisjointFrom(qualitativeStateSets.maybeStates);
 
+                // create multiplier and execute the calculation for 1 additional step
+                auto multiplier = storm::solver::MultiplierFactory<ValueType>().create(env, transitionMatrix);
+
+                uint sizeMaybeStateChoiceValues = 0;
+                for(uint counter = 0; counter < qualitativeStateSets.maybeStates.size(); counter++) {
+                    if(qualitativeStateSets.maybeStates.get(counter)) {
+                        sizeMaybeStateChoiceValues += transitionMatrix.getRowGroupSize(counter);
+                    }
+                }
+
+                std::vector<ValueType> maybeStateChoiceValues = std::vector<ValueType>(sizeMaybeStateChoiceValues, storm::utility::zero<ValueType>());
+                
                 // Check whether we need to compute exact probabilities for some states.
-                if (qualitative || maybeStatesNotRelevant) {
+                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.
                     storm::utility::vector::setVectorValues<ValueType>(result, qualitativeStateSets.maybeStates, storm::utility::convertNumber<ValueType>(0.5));
                 } else {
@@ -649,6 +670,39 @@ namespace storm {
                                 extractSchedulerChoices(*scheduler, resultForMaybeStates.getScheduler(), qualitativeStateSets.maybeStates);
                             }
                         }
+                        if (goal.isShieldingTask()) {
+                            std::vector<ValueType> subResult;
+                            uint sizeChoiceValues = 0;
+                            for(uint counter = 0; counter < qualitativeStateSets.maybeStates.size(); counter++) {
+                                if(qualitativeStateSets.maybeStates.get(counter)) {
+                                    subResult.push_back(result.at(counter));
+                                }
+                            }
+
+                            submatrix = transitionMatrix.getSubmatrix(true, qualitativeStateSets.maybeStates, qualitativeStateSets.maybeStates, false);
+                            auto sub_multiplier = storm::solver::MultiplierFactory<ValueType>().create(env, submatrix);
+                            sub_multiplier->multiply(env, subResult, &b, maybeStateChoiceValues);
+
+                        }
+                    }
+                }
+
+                std::vector<ValueType> choiceValues = std::vector<ValueType>(transitionMatrix.getRowGroupIndices().at(transitionMatrix.getRowGroupIndices().size() - 1), storm::utility::zero<ValueType>());
+                auto choice_it = maybeStateChoiceValues.begin();
+                for(uint state = 0; state < transitionMatrix.getRowGroupIndices().size() - 1; state++) {
+                    uint rowGroupSize = transitionMatrix.getRowGroupIndices().at(state + 1) - transitionMatrix.getRowGroupIndices().at(state);
+                    if (qualitativeStateSets.maybeStates.get(state)) {
+                        for(uint choice = 0; choice < rowGroupSize; choice++, choice_it++) {
+                            choiceValues.at(transitionMatrix.getRowGroupIndices().at(state) + choice) = *choice_it;
+                        }
+                    } else if (qualitativeStateSets.statesWithProbability0.get(state)) {
+                        for(uint choice = 0; choice < rowGroupSize; choice++) {
+                            choiceValues.at(transitionMatrix.getRowGroupIndices().at(state) + choice) = 0;
+                        }
+                    } else if (qualitativeStateSets.statesWithProbability1.get(state)) {
+                        for(uint choice = 0; choice < rowGroupSize; choice++) {
+                            choiceValues.at(transitionMatrix.getRowGroupIndices().at(state) + choice) = 1;
+                        }
                     }
                 }
 
@@ -664,7 +718,7 @@ namespace storm {
                 STORM_LOG_ASSERT((!produceScheduler && !scheduler) || scheduler->isMemorylessScheduler(), "Expected a memoryless scheduler");
 
                 // Return result.
-                return MDPSparseModelCheckingHelperReturnType<ValueType>(std::move(result), std::move(scheduler));
+                return MDPSparseModelCheckingHelperReturnType<ValueType>(std::move(result), std::move(qualitativeStateSets.maybeStates), std::move(scheduler), std::move(choiceValues));
             }
 
             template<typename ValueType>
@@ -678,7 +732,6 @@ namespace storm {
                             statesInPsiMecs.set(stateActionsPair.first, true);
                         }
                     }
-                    
                     return computeUntilProbabilities(env, std::move(goal), transitionMatrix, backwardTransitions, psiStates, statesInPsiMecs, qualitative, produceScheduler);
                 } else {
                     goal.oneMinus();
@@ -686,6 +739,9 @@ namespace storm {
                     for (auto& element : result.values) {
                         element = storm::utility::one<ValueType>() - element;
                     }
+                    for (auto& choice : result.choiceValues) {
+                        choice = storm::utility::one<ValueType>() - choice;
+                    }
                     return result;
                 }
             }
@@ -1172,8 +1228,9 @@ namespace storm {
                 STORM_LOG_ASSERT((!produceScheduler && !scheduler) || scheduler->isDeterministicScheduler(), "Expected a deterministic scheduler");
                 STORM_LOG_ASSERT((!produceScheduler && !scheduler) || scheduler->isMemorylessScheduler(), "Expected a memoryless scheduler");
 
+                std::vector<ValueType> choiceValues;
 
-                return MDPSparseModelCheckingHelperReturnType<ValueType>(std::move(result), std::move(scheduler));
+                return MDPSparseModelCheckingHelperReturnType<ValueType>(std::move(result), std::move(qualitativeStateSets.maybeStates), std::move(scheduler), std::move(choiceValues));
             }
             
             template<typename ValueType>

src/storm/modelchecker/prctl/helper/SparseMdpPrctlHelper.h

@@ -41,7 +41,7 @@ namespace storm {
 
                 static std::map<storm::storage::sparse::state_type, ValueType> computeRewardBoundedValues(Environment const& env, OptimizationDirection dir, rewardbounded::MultiDimensionalRewardUnfolding<ValueType, true>& rewardUnfolding, storm::storage::BitVector const& initialStates);
                 
-                static std::vector<ValueType> computeNextProbabilities(Environment const& env, OptimizationDirection dir, storm::storage::SparseMatrix<ValueType> const& transitionMatrix, storm::storage::BitVector const& nextStates);
+                static MDPSparseModelCheckingHelperReturnType<ValueType> computeNextProbabilities(Environment const& env, storm::solver::SolveGoal<ValueType>&& goal, OptimizationDirection dir, storm::storage::SparseMatrix<ValueType> const& transitionMatrix, storm::storage::BitVector const& nextStates);
 
                 static MDPSparseModelCheckingHelperReturnType<ValueType> computeUntilProbabilities(Environment const& env, storm::solver::SolveGoal<ValueType>&& 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());
                 

src/storm/modelchecker/rpatl/SparseSmgRpatlModelChecker.cpp

@@ -24,7 +24,7 @@
 
 #include "storm/models/sparse/StandardRewardModel.h"
 
-#include "storm/shields/shield-handling.h"
+#include "storm/shields/ShieldHandling.h"
 
 #include "storm/settings/modules/GeneralSettings.h"
 
@@ -212,7 +212,7 @@ namespace storm {
 
             std::unique_ptr<CheckResult> result(new ExplicitQuantitativeCheckResult<ValueType>(std::move(values)));
             if(checkTask.isShieldingTask()) {
-                tempest::shields::createOptimalShield<ValueType>(std::make_shared<storm::models::sparse::Smg<ValueType>>(this->getModel()), helper.getProducedOptimalChoices(), checkTask.getShieldingExpression(), checkTask.getOptimizationDirection(), statesOfCoalition, statesOfCoalition);
+                tempest::shields::createQuantitativeShield<ValueType>(std::make_shared<storm::models::sparse::Smg<ValueType>>(this->getModel()), helper.getProducedOptimalChoices(), checkTask.getShieldingExpression(), checkTask.getOptimizationDirection(), statesOfCoalition, statesOfCoalition);
             } else if (checkTask.isProduceSchedulersSet()) {
                 result->asExplicitQuantitativeCheckResult<ValueType>().setScheduler(std::make_unique<storm::storage::Scheduler<ValueType>>(helper.extractScheduler()));
             }

src/storm/modelchecker/rpatl/helper/SparseSmgRpatlHelper.cpp

@@ -41,8 +41,7 @@ namespace storm {
                 }
 
                 viHelper.performValueIteration(env, x, b, goal.direction());
-                //if(goal.isShieldingTask()) {
-                if(true) {
+                if(goal.isShieldingTask()) {
                     viHelper.getChoiceValues(env, x, constrainedChoiceValues);
                 }
                 viHelper.fillResultVector(x, relevantStates, psiStates);
@@ -107,13 +106,12 @@ namespace storm {
                 // create multiplier and execute the calculation for 1 step
                 auto multiplier = storm::solver::MultiplierFactory<ValueType>().create(env, transitionMatrix);
                 std::vector<ValueType> choiceValues = std::vector<ValueType>(transitionMatrix.getRowCount(), storm::utility::zero<ValueType>());
-
-                //if(goal.isShieldingTask()) {
-                if (true) {
+                if (goal.isShieldingTask()) {
                     multiplier->multiply(env, x, &b, choiceValues);
-                }
+                    multiplier->reduce(env, goal.direction(), choiceValues, transitionMatrix.getRowGroupIndices(), x, &statesOfCoalition);
+                } else {
                     multiplier->multiplyAndReduce(env, goal.direction(), x, &b, x, nullptr, &statesOfCoalition);
-
+                }
                 return SMGSparseModelCheckingHelperReturnType<ValueType>(std::move(x), std::move(allStates), nullptr, std::move(choiceValues));
             }
 

src/storm/shields/ShieldHandling.cpp

@@ -0,0 +1,47 @@
+#include "ShieldHandling.h"
+
+namespace tempest {
+    namespace shields {
+        std::string shieldFilename(std::shared_ptr<storm::logic::ShieldExpression const> const& shieldingExpression) {
+            return shieldingExpression->getFilename() + ".shield";
+        }
+
+        template<typename ValueType, typename IndexType>
+        void createShield(std::shared_ptr<storm::models::sparse::Model<ValueType>> model, std::vector<ValueType> const& choiceValues, std::shared_ptr<storm::logic::ShieldExpression const> const& shieldingExpression, storm::OptimizationDirection optimizationDirection, storm::storage::BitVector relevantStates, boost::optional<storm::storage::BitVector> coalitionStates) {
+            std::ofstream stream;
+            storm::utility::openFile(shieldFilename(shieldingExpression), stream);
+            if(shieldingExpression->isPreSafetyShield()) {
+                PreSafetyShield<ValueType, IndexType> shield(model->getTransitionMatrix().getRowGroupIndices(), choiceValues, shieldingExpression, optimizationDirection, relevantStates, coalitionStates);
+                shield.construct().printToStream(stream, shieldingExpression, model);
+            } else if(shieldingExpression->isPostSafetyShield()) {
+                PostSafetyShield<ValueType, IndexType> shield(model->getTransitionMatrix().getRowGroupIndices(), choiceValues, shieldingExpression, optimizationDirection, relevantStates, coalitionStates);
+                shield.construct().printToStream(stream, shieldingExpression, model);
+            } else {
+                STORM_LOG_THROW(false, storm::exceptions::InvalidArgumentException, "Unknown Shielding Type: " + shieldingExpression->typeToString());
+                storm::utility::closeFile(stream);
+            }
+            storm::utility::closeFile(stream);
+        }
+
+        template<typename ValueType, typename IndexType>
+        void createQuantitativeShield(std::shared_ptr<storm::models::sparse::Model<ValueType>> model, std::vector<uint64_t> const& precomputedChoices, std::shared_ptr<storm::logic::ShieldExpression const> const& shieldingExpression, storm::OptimizationDirection optimizationDirection, storm::storage::BitVector relevantStates, boost::optional<storm::storage::BitVector> coalitionStates) {
+            std::ofstream stream;
+            storm::utility::openFile(shieldFilename(shieldingExpression), stream);
+            if(shieldingExpression->isOptimalShield()) {
+                OptimalShield<ValueType, IndexType> shield(model->getTransitionMatrix().getRowGroupIndices(), precomputedChoices, shieldingExpression, optimizationDirection, relevantStates, coalitionStates);
+                shield.construct().printToStream(stream, shieldingExpression, model);
+            } else {
+                STORM_LOG_THROW(false, storm::exceptions::InvalidArgumentException, "Unknown Shielding Type: " + shieldingExpression->typeToString());
+                storm::utility::closeFile(stream);
+            }
+            storm::utility::closeFile(stream);
+        }
+        // Explicitly instantiate appropriate
+        template void createShield<double, typename storm::storage::SparseMatrix<double>::index_type>(std::shared_ptr<storm::models::sparse::Model<double>> model, std::vector<double> const& choiceValues, std::shared_ptr<storm::logic::ShieldExpression const> const& shieldingExpression, storm::OptimizationDirection optimizationDirection, storm::storage::BitVector relevantStates, boost::optional<storm::storage::BitVector> coalitionStates);
+        template void createQuantitativeShield<double, typename storm::storage::SparseMatrix<double>::index_type>(std::shared_ptr<storm::models::sparse::Model<double>> model, std::vector<uint64_t> const& precomputedChoices, std::shared_ptr<storm::logic::ShieldExpression const> const& shieldingExpression, storm::OptimizationDirection optimizationDirection, storm::storage::BitVector relevantStates, boost::optional<storm::storage::BitVector> coalitionStates);
+#ifdef STORM_HAVE_CARL
+        template void createShield<storm::RationalNumber, typename storm::storage::SparseMatrix<storm::RationalNumber>::index_type>(std::shared_ptr<storm::models::sparse::Model<storm::RationalNumber>> model, std::vector<storm::RationalNumber> const& choiceValues, std::shared_ptr<storm::logic::ShieldExpression const> const& shieldingExpression, storm::OptimizationDirection optimizationDirection, storm::storage::BitVector relevantStates, boost::optional<storm::storage::BitVector> coalitionStates);
+        template void createQuantitativeShield<storm::RationalNumber, typename storm::storage::SparseMatrix<storm::RationalNumber>::index_type>(std::shared_ptr<storm::models::sparse::Model<storm::RationalNumber>> model, std::vector<uint64_t> const& precomputedChoices, std::shared_ptr<storm::logic::ShieldExpression const> const& shieldingExpression, storm::OptimizationDirection optimizationDirection, storm::storage::BitVector relevantStates, boost::optional<storm::storage::BitVector> coalitionStates);
+#endif
+    }
+}

src/storm/shields/ShieldHandling.h

@@ -0,0 +1,32 @@
+#pragma once
+
+#include <iostream>
+#include <boost/optional.hpp>
+#include <memory>
+
+#include "storm/storage/Scheduler.h"
+#include "storm/storage/BitVector.h"
+
+#include "storm/logic/ShieldExpression.h"
+
+#include "storm/shields/AbstractShield.h"
+#include "storm/shields/PreSafetyShield.h"
+#include "storm/shields/PostSafetyShield.h"
+#include "storm/shields/OptimalShield.h"
+
+#include "storm/io/file.h"
+#include "storm/utility/macros.h"
+
+#include "storm/exceptions/InvalidArgumentException.h"
+
+namespace tempest {
+    namespace shields {
+        std::string shieldFilename(std::shared_ptr<storm::logic::ShieldExpression const> const& shieldingExpression);
+
+        template<typename ValueType, typename IndexType = storm::storage::sparse::state_type>
+        void createShield(std::shared_ptr<storm::models::sparse::Model<ValueType>> model, std::vector<ValueType> const& choiceValues, std::shared_ptr<storm::logic::ShieldExpression const> const& shieldingExpression, storm::OptimizationDirection optimizationDirection, storm::storage::BitVector relevantStates, boost::optional<storm::storage::BitVector> coalitionStates);
+
+        template<typename ValueType, typename IndexType = storm::storage::sparse::state_type>
+        void createQuantitativeShield(std::shared_ptr<storm::models::sparse::Model<ValueType>> model, std::vector<uint64_t> const& precomputedChoices, std::shared_ptr<storm::logic::ShieldExpression const> const& shieldingExpression, storm::OptimizationDirection optimizationDirection, storm::storage::BitVector relevantStates, boost::optional<storm::storage::BitVector> coalitionStates);
+    }
+}

src/storm/shields/shield-handling.h

@@ -1,60 +0,0 @@
-#pragma once
-
-#include <iostream>
-#include <boost/optional.hpp>
-#include <memory>
-
-#include "storm/storage/Scheduler.h"
-#include "storm/storage/BitVector.h"
-
-#include "storm/logic/ShieldExpression.h"
-
-#include "storm/shields/AbstractShield.h"
-#include "storm/shields/PreSafetyShield.h"
-#include "storm/shields/PostSafetyShield.h"
-#include "storm/shields/OptimalShield.h"
-
-#include "storm/io/file.h"
-#include "storm/utility/macros.h"
-
-#include "storm/exceptions/InvalidArgumentException.h"
-
-
-namespace tempest {
-    namespace shields {
-        std::string shieldFilename(std::shared_ptr<storm::logic::ShieldExpression const> const& shieldingExpression) {
-            return shieldingExpression->getFilename() + ".shield";
-        }
-
-        template<typename ValueType, typename IndexType = storm::storage::sparse::state_type>
-        void createShield(std::shared_ptr<storm::models::sparse::Model<ValueType>> model, std::vector<ValueType> const& choiceValues, std::shared_ptr<storm::logic::ShieldExpression const> const& shieldingExpression, storm::OptimizationDirection optimizationDirection, storm::storage::BitVector relevantStates, boost::optional<storm::storage::BitVector> coalitionStates) {
-            std::ofstream stream;
-            storm::utility::openFile(shieldFilename(shieldingExpression), stream);
-            if(shieldingExpression->isPreSafetyShield()) {
-                PreSafetyShield<ValueType, IndexType> shield(model->getTransitionMatrix().getRowGroupIndices(), choiceValues, shieldingExpression, optimizationDirection, relevantStates, coalitionStates);
-                shield.construct().printToStream(stream, shieldingExpression, model);
-            } else if(shieldingExpression->isPostSafetyShield()) {
-                PostSafetyShield<ValueType, IndexType> shield(model->getTransitionMatrix().getRowGroupIndices(), choiceValues, shieldingExpression, optimizationDirection, relevantStates, coalitionStates);
-                shield.construct().printToStream(stream, shieldingExpression, model);
-            } else {
-                STORM_LOG_THROW(false, storm::exceptions::InvalidArgumentException, "Unknown Shielding Type: " + shieldingExpression->typeToString());
-                storm::utility::closeFile(stream);
-            }
-            storm::utility::closeFile(stream);
-        }
-
-        template<typename ValueType, typename IndexType = storm::storage::sparse::state_type>
-        void createOptimalShield(std::shared_ptr<storm::models::sparse::Model<ValueType>> model, std::vector<uint64_t> const& precomputedChoices, std::shared_ptr<storm::logic::ShieldExpression const> const& shieldingExpression, storm::OptimizationDirection optimizationDirection, storm::storage::BitVector relevantStates, boost::optional<storm::storage::BitVector> coalitionStates) {
-            std::ofstream stream;
-            storm::utility::openFile(shieldFilename(shieldingExpression), stream);
-            if(shieldingExpression->isOptimalShield()) {
-                OptimalShield<ValueType, IndexType> shield(model->getTransitionMatrix().getRowGroupIndices(), precomputedChoices, shieldingExpression, optimizationDirection, relevantStates, coalitionStates);
-                shield.construct().printToStream(stream, shieldingExpression, model);
-            } else {
-                STORM_LOG_THROW(false, storm::exceptions::InvalidArgumentException, "Unknown Shielding Type: " + shieldingExpression->typeToString());
-                storm::utility::closeFile(stream);
-            }
-            storm::utility::closeFile(stream);
-        }
-    }
-}

src/storm/solver/Multiplier.cpp

@@ -68,12 +68,53 @@ namespace storm {
             }
         }
 
+        template<typename ValueType>
+        void Multiplier<ValueType>::repeatedMultiplyAndReduceWithChoices(Environment const& env, OptimizationDirection const& dir, std::vector<ValueType>& x, std::vector<ValueType> const* b, uint64_t n, storm::storage::BitVector const* dirOverride, std::vector<ValueType>& choiceValues, std::vector<storm::storage::SparseMatrix<double>::index_type> rowGroupIndices) const {
+            storm::utility::ProgressMeasurement progress("multiplications");
+            progress.setMaxCount(n);
+            progress.startNewMeasurement(0);
+            for (uint64_t i = 0; i < n; ++i) {
+                multiply(env, x, b, choiceValues);
+                reduce(env, dir, choiceValues, rowGroupIndices, x);
+                if (storm::utility::resources::isTerminate()) {
+                    STORM_LOG_WARN("Aborting after " << i << " of " << n << " multiplications");
+                    break;
+                }
+            }
+        }
+
         template<typename ValueType>
         void Multiplier<ValueType>::multiplyRow2(uint64_t const& rowIndex, std::vector<ValueType> const& x1, ValueType& val1, std::vector<ValueType> const& x2, ValueType& val2) const {
             multiplyRow(rowIndex, x1, val1);
             multiplyRow(rowIndex, x2, val2);
         }
 
+        template<typename ValueType>
+        void Multiplier<ValueType>::reduce(Environment const& env, OptimizationDirection const& dir, std::vector<ValueType> const& choiceValues, std::vector<storm::storage::SparseMatrix<double>::index_type> rowGroupIndices, std::vector<ValueType>& result, storm::storage::BitVector const* dirOverride) const {
+            auto choice_it = choiceValues.begin();
+            for(uint state = 0; state < rowGroupIndices.size() - 1; state++) {
+                uint rowGroupSize = rowGroupIndices[state + 1] - rowGroupIndices[state];
+                if(dirOverride != nullptr) {
+                    if((dir == storm::OptimizationDirection::Minimize && !dirOverride->get(state)) || (dir == storm::OptimizationDirection::Maximize && dirOverride->get(state))) {
+                        result.at(state) = *std::min_element(choice_it, choice_it + rowGroupSize);
+                        choice_it += rowGroupSize;
+                    }
+                    else {
+                        result.at(state) = *std::max_element(choice_it, choice_it + rowGroupSize);
+                        choice_it += rowGroupSize;
+                    }
+                } else {
+                    if(dir == storm::OptimizationDirection::Minimize) {
+                        result.at(state) = *std::min_element(choice_it, choice_it + rowGroupSize);
+                        choice_it += rowGroupSize;
+                    } else {
+                        result.at(state) = *std::max_element(choice_it, choice_it + rowGroupSize);
+                        choice_it += rowGroupSize;
+                    }
+                }
+            }
+        }
+
         template<typename ValueType>
         std::unique_ptr<Multiplier<ValueType>> MultiplierFactory<ValueType>::create(Environment const& env, storm::storage::SparseMatrix<ValueType> const& matrix) {
             auto type = env.solver().multiplier().getType();

src/storm/solver/Multiplier.h

@@ -8,6 +8,8 @@
 
 #include "storm/solver/OptimizationDirection.h"
 #include "storm/solver/MultiplicationStyle.h"
+#include "storm/storage/SparseMatrix.h"
+
 
 namespace storm {
 
@@ -119,6 +121,8 @@ namespace storm {
              */
             void repeatedMultiplyAndReduce(Environment const& env, OptimizationDirection const& dir, std::vector<ValueType>& x, std::vector<ValueType> const* b, uint64_t n, storm::storage::BitVector const* dirOverride = nullptr) const;
 
+            void repeatedMultiplyAndReduceWithChoices(const Environment &env, const OptimizationDirection &dir, std::vector<ValueType> &x, const std::vector<ValueType> *b, uint64_t n, const storage::BitVector *dirOverride, std::vector<ValueType> &choiceValues, std::vector<unsigned long> rowGroupIndices) const;
+
             /*!
              * Multiplies the row with the given index with x and adds the result to the provided value
              * @param rowIndex The index of the considered row
@@ -137,9 +141,12 @@ namespace storm {
              */
             virtual void multiplyRow2(uint64_t const& rowIndex, std::vector<ValueType> const& x1, ValueType& val1, std::vector<ValueType> const& x2, ValueType& val2) const;
 
+            void reduce(Environment const& env, OptimizationDirection const& dir, std::vector<ValueType> const& choiceValues, std::vector<storm::storage::SparseMatrix<double>::index_type> rowGroupIndices, std::vector<ValueType>& result, storm::storage::BitVector const* dirOverride = nullptr) const;
+
         protected:
             mutable std::unique_ptr<std::vector<ValueType>> cachedVector;
             storm::storage::SparseMatrix<ValueType> const& matrix;
+
         };
 
         template<typename ValueType>

src/storm/solver/SolveGoal.cpp

@@ -122,6 +122,11 @@ namespace storm {
             relevantValueVector = std::move(values);
         }
 
+        template<typename ValueType>
+        bool SolveGoal<ValueType>::isShieldingTask() const {
+            return shieldingTask;
+        }
+
         template class SolveGoal<double>;
         
 #ifdef STORM_HAVE_CARL

src/storm/solver/SolveGoal.h

@@ -51,6 +51,7 @@ namespace storm {
                     comparisonType = checkTask.getBoundComparisonType();
                     threshold = checkTask.getBoundThreshold();
                 }
+                shieldingTask = checkTask.isShieldingTask();
             }
             
             SolveGoal(bool minimize);
@@ -84,12 +85,16 @@ namespace storm {
             void restrictRelevantValues(storm::storage::BitVector const& filter);
             void setRelevantValues(storm::storage::BitVector&& values);
 
+            bool isShieldingTask() const;
+
         private:
             boost::optional<OptimizationDirection> optimizationDirection;
             
             boost::optional<storm::logic::ComparisonType> comparisonType;
             boost::optional<ValueType> threshold;
             boost::optional<storm::storage::BitVector> relevantValueVector;
+            // We only want to know if it **is** a shielding task
+            bool shieldingTask;
         };
         
         template<typename ValueType, typename MatrixType>