From 999fd0752cae53775fd04fa2ff52bee1f1657dd8 Mon Sep 17 00:00:00 2001
From: TimQu <tim.quatmann@cs.rwth-aachen.de>
Date: Tue, 5 Sep 2017 13:49:12 +0200
Subject: [PATCH] The model memory product can now retrieve the reachable
 states without actually building the product

---
 .../MultiDimensionalRewardUnfolding.cpp       |   4 +-
 .../models/sparse/NondeterministicModel.cpp   |  19 +-
 .../SparseModelMemoryProduct.cpp              | 166 ++++++++++--------
 .../SparseModelMemoryProduct.h                |  45 +++--
 4 files changed, 132 insertions(+), 102 deletions(-)

diff --git a/src/storm/modelchecker/multiobjective/rewardbounded/MultiDimensionalRewardUnfolding.cpp b/src/storm/modelchecker/multiobjective/rewardbounded/MultiDimensionalRewardUnfolding.cpp
index 479286ad2..04c099e6e 100644
--- a/src/storm/modelchecker/multiobjective/rewardbounded/MultiDimensionalRewardUnfolding.cpp
+++ b/src/storm/modelchecker/multiobjective/rewardbounded/MultiDimensionalRewardUnfolding.cpp
@@ -484,8 +484,8 @@ namespace storm {
                 productToMemoryStateMap.resize(numProductStates, std::numeric_limits<uint64_t>::max());
                 for (uint64_t modelState = 0; modelState < numModelStates; ++modelState) {
                     for (uint64_t memoryState = 0; memoryState < numMemoryStates; ++memoryState) {
-                        uint64_t productState = productBuilder.getResultState(modelState, memoryState);
-                        if (productState < numProductStates) {
+                        if (productBuilder.isStateReachable(modelState, memoryState)) {
+                            uint64_t productState = productBuilder.getResultState(modelState, memoryState);
                             modelMemoryToProductStateMap[modelState * numMemoryStates + memoryState] = productState;
                             productToModelStateMap[productState] = modelState;
                             productToMemoryStateMap[productState] = memoryState;
diff --git a/src/storm/models/sparse/NondeterministicModel.cpp b/src/storm/models/sparse/NondeterministicModel.cpp
index acdbd3885..cd67e384f 100644
--- a/src/storm/models/sparse/NondeterministicModel.cpp
+++ b/src/storm/models/sparse/NondeterministicModel.cpp
@@ -51,22 +51,11 @@ namespace storm {
             
             template<typename ValueType, typename RewardModelType>
             std::shared_ptr<storm::models::sparse::Model<ValueType, RewardModelType>> NondeterministicModel<ValueType, RewardModelType>::applyScheduler(storm::storage::Scheduler<ValueType> const& scheduler, bool dropUnreachableStates) {
-                if (scheduler.isMemorylessScheduler()) {
-                    auto memStruct = storm::storage::MemoryStructureBuilder<ValueType, RewardModelType>::buildTrivialMemoryStructure(*this);
-                    auto memoryProduct = memStruct.product(*this);
-                    if (!dropUnreachableStates) {
-                        memoryProduct.setBuildFullProduct();
-                    }
-                    return memoryProduct.build(scheduler);
-                } else {
-                    boost::optional<storm::storage::MemoryStructure> const& memStruct = scheduler.getMemoryStructure();
-                    STORM_LOG_ASSERT(memStruct, "Memoryless scheduler without memory structure.");
-                    auto memoryProduct = memStruct->product(*this);
-                    if (!dropUnreachableStates) {
-                        memoryProduct.setBuildFullProduct();
-                    }
-                    return memoryProduct.build(scheduler);
+                storm::storage::SparseModelMemoryProduct<ValueType, RewardModelType> memoryProduct(*this, scheduler);
+                if (!dropUnreachableStates) {
+                    memoryProduct.setBuildFullProduct();
                 }
+                return memoryProduct.build();
             }
             
             template<typename ValueType, typename RewardModelType>
diff --git a/src/storm/storage/memorystructure/SparseModelMemoryProduct.cpp b/src/storm/storage/memorystructure/SparseModelMemoryProduct.cpp
index 49d69e240..3907f0625 100644
--- a/src/storm/storage/memorystructure/SparseModelMemoryProduct.cpp
+++ b/src/storm/storage/memorystructure/SparseModelMemoryProduct.cpp
@@ -9,6 +9,7 @@
 #include "storm/modelchecker/propositional/SparsePropositionalModelChecker.h"
 #include "storm/modelchecker/results/ExplicitQualitativeCheckResult.h"
 #include "storm/utility/constants.h"
+#include "storm/storage/memorystructure/MemoryStructureBuilder.h"
 #include "storm/utility/macros.h"
 #include "storm/utility/builder.h"
 
@@ -18,91 +19,113 @@ namespace storm {
     namespace storage {
 
         template <typename ValueType, typename RewardModelType>
-        SparseModelMemoryProduct<ValueType, RewardModelType>::SparseModelMemoryProduct(storm::models::sparse::Model<ValueType, RewardModelType> const& sparseModel, storm::storage::MemoryStructure const& memoryStructure) : memoryStateCount(memoryStructure.getNumberOfStates()), model(sparseModel), memory(memoryStructure) {
+        SparseModelMemoryProduct<ValueType, RewardModelType>::SparseModelMemoryProduct(storm::models::sparse::Model<ValueType, RewardModelType> const& sparseModel, storm::storage::MemoryStructure const& memoryStructure) : isInitialized(false), memoryStateCount(memoryStructure.getNumberOfStates()), model(sparseModel), memory(memoryStructure), scheduler(boost::none) {
             reachableStates = storm::storage::BitVector(model.getNumberOfStates() * memoryStateCount, false);
         }
         
+        template <typename ValueType, typename RewardModelType>
+        SparseModelMemoryProduct<ValueType, RewardModelType>::SparseModelMemoryProduct(storm::models::sparse::Model<ValueType, RewardModelType> const& sparseModel, storm::storage::Scheduler<ValueType> const& scheduler) : isInitialized(false), memoryStateCount(scheduler.getNumberOfMemoryStates()), model(sparseModel), localMemory(scheduler.getMemoryStructure() ? boost::optional<MemoryStructure>() : boost::optional<MemoryStructure>(storm::storage::MemoryStructureBuilder<ValueType, RewardModelType>::buildTrivialMemoryStructure(model))), memory(scheduler.getMemoryStructure() ? scheduler.getMemoryStructure().get() : localMemory.get()), scheduler(scheduler) {
+            reachableStates = storm::storage::BitVector(model.getNumberOfStates() * memoryStateCount, false);
+        }
+        
+        template <typename ValueType, typename RewardModelType>
+        void SparseModelMemoryProduct<ValueType, RewardModelType>::initialize() {
+            if (!isInitialized) {
+                uint64_t modelStateCount = model.getNumberOfStates();
+            
+                computeMemorySuccessors();
+                
+                // Get the initial states and reachable states. A stateIndex s corresponds to the model state (s / memoryStateCount) and memory state (s % memoryStateCount)
+                storm::storage::BitVector initialStates(modelStateCount * memoryStateCount, false);
+                auto memoryInitIt = memory.getInitialMemoryStates().begin();
+                for (auto const& modelInit : model.getInitialStates()) {
+                    initialStates.set(modelInit * memoryStateCount + *memoryInitIt, true);
+                    ++memoryInitIt;
+                }
+                STORM_LOG_ASSERT(memoryInitIt == memory.getInitialMemoryStates().end(), "Unexpected number of initial states.");
+                
+                computeReachableStates(initialStates);
+                
+                // Compute the mapping to the states of the result
+                uint64_t reachableStateCount = 0;
+                toResultStateMapping = std::vector<uint64_t> (model.getNumberOfStates() * memoryStateCount, std::numeric_limits<uint64_t>::max());
+                for (auto const& reachableState : reachableStates) {
+                    toResultStateMapping[reachableState] = reachableStateCount;
+                    ++reachableStateCount;
+                }
+                
+                isInitialized = true;
+            }
+        }
+        
         template <typename ValueType, typename RewardModelType>
         void SparseModelMemoryProduct<ValueType, RewardModelType>::addReachableState(uint64_t const& modelState, uint64_t const& memoryState) {
+            isInitialized = false;
             reachableStates.set(modelState * memoryStateCount + memoryState, true);
         }
 
        template <typename ValueType, typename RewardModelType>
        void SparseModelMemoryProduct<ValueType, RewardModelType>::setBuildFullProduct() {
+            isInitialized = false;
             reachableStates.fill();
        }
         
         template <typename ValueType, typename RewardModelType>
-        std::shared_ptr<storm::models::sparse::Model<ValueType, RewardModelType>> SparseModelMemoryProduct<ValueType, RewardModelType>::build(boost::optional<storm::storage::Scheduler<ValueType>> const& scheduler) {
-            
-            uint64_t modelStateCount = model.getNumberOfStates();
+        std::shared_ptr<storm::models::sparse::Model<ValueType, RewardModelType>> SparseModelMemoryProduct<ValueType, RewardModelType>::build() {
             
-            std::vector<uint64_t> memorySuccessors = computeMemorySuccessors();
+            initialize();
             
-            // Get the initial states and reachable states. A stateIndex s corresponds to the model state (s / memoryStateCount) and memory state (s % memoryStateCount)
-            storm::storage::BitVector initialStates(modelStateCount * memoryStateCount, false);
-            auto memoryInitIt = memory.getInitialMemoryStates().begin();
-            for (auto const& modelInit : model.getInitialStates()) {
-                initialStates.set(modelInit * memoryStateCount + *memoryInitIt, true);
-                ++memoryInitIt;
-            }
-            STORM_LOG_ASSERT(memoryInitIt == memory.getInitialMemoryStates().end(), "Unexpected number of initial states.");
-            
-            computeReachableStates(memorySuccessors, initialStates, scheduler);
-            
-            // Compute the mapping to the states of the result
-            uint64_t reachableStateCount = 0;
-            toResultStateMapping = std::vector<uint64_t> (model.getNumberOfStates() * memoryStateCount, std::numeric_limits<uint64_t>::max());
-            for (auto const& reachableState : reachableStates) {
-                toResultStateMapping[reachableState] = reachableStateCount;
-                ++reachableStateCount;
-            }
-                
             // Build the model components
             storm::storage::SparseMatrix<ValueType> transitionMatrix;
             if (scheduler) {
-                transitionMatrix = buildTransitionMatrixForScheduler(memorySuccessors, scheduler.get());
+                transitionMatrix = buildTransitionMatrixForScheduler();
             } else if (model.getTransitionMatrix().hasTrivialRowGrouping()) {
-                transitionMatrix = buildDeterministicTransitionMatrix(memorySuccessors);
+                transitionMatrix = buildDeterministicTransitionMatrix();
             } else {
-                transitionMatrix = buildNondeterministicTransitionMatrix(memorySuccessors);
+                transitionMatrix = buildNondeterministicTransitionMatrix();
             }
             storm::models::sparse::StateLabeling labeling = buildStateLabeling(transitionMatrix);
-            std::unordered_map<std::string, RewardModelType> rewardModels = buildRewardModels(transitionMatrix, memorySuccessors, scheduler);
+            std::unordered_map<std::string, RewardModelType> rewardModels = buildRewardModels(transitionMatrix);
 
-            // Add the label for the initial states. We need to translate the state indices w.r.t. the set of reachable states.
-            labeling.addLabel("init", initialStates % reachableStates);
-            
             
-            return buildResult(std::move(transitionMatrix), std::move(labeling), std::move(rewardModels), scheduler);
+            return buildResult(std::move(transitionMatrix), std::move(labeling), std::move(rewardModels));
 
         }
             
         template <typename ValueType, typename RewardModelType>
-        uint64_t const& SparseModelMemoryProduct<ValueType, RewardModelType>::getResultState(uint64_t const& modelState, uint64_t const& memoryState) const {
-                return toResultStateMapping[modelState * memoryStateCount + memoryState];
+        bool SparseModelMemoryProduct<ValueType, RewardModelType>::isStateReachable(uint64_t const& modelState, uint64_t const& memoryState) {
+            STORM_LOG_ASSERT(modelState < getOriginalModel().getNumberOfStates(), "Invalid model state: " << modelState << ".");
+            STORM_LOG_ASSERT(memoryState < memoryStateCount, "Invalid memory state: " << memoryState << ".");
+            initialize();
+            return reachableStates.get(modelState * memoryStateCount + memoryState);
+        }
+        
+        template <typename ValueType, typename RewardModelType>
+        uint64_t const& SparseModelMemoryProduct<ValueType, RewardModelType>::getResultState(uint64_t const& modelState, uint64_t const& memoryState) {
+            initialize();
+            STORM_LOG_ASSERT(isStateReachable(modelState, memoryState), "Tried to get unreachable product state (" << modelState << "," << memoryState << ")");
+            return toResultStateMapping[modelState * memoryStateCount + memoryState];
         }
             
         template <typename ValueType, typename RewardModelType>
-        std::vector<uint64_t> SparseModelMemoryProduct<ValueType, RewardModelType>::computeMemorySuccessors() const {
+        void SparseModelMemoryProduct<ValueType, RewardModelType>::computeMemorySuccessors() {
             uint64_t modelTransitionCount = model.getTransitionMatrix().getEntryCount();
-            std::vector<uint64_t> result(modelTransitionCount * memoryStateCount, std::numeric_limits<uint64_t>::max());
+            memorySuccessors = std::vector<uint64_t>(modelTransitionCount * memoryStateCount, std::numeric_limits<uint64_t>::max());
             
             for (uint64_t memoryState = 0; memoryState < memoryStateCount; ++memoryState) {
                 for (uint64_t transitionGoal = 0; transitionGoal < memoryStateCount; ++transitionGoal) {
                     auto const& memoryTransition = memory.getTransitionMatrix()[memoryState][transitionGoal];
                     if (memoryTransition) {
                         for (auto const& modelTransitionIndex : memoryTransition.get()) {
-                            result[modelTransitionIndex * memoryStateCount + memoryState] = transitionGoal;
+                            memorySuccessors[modelTransitionIndex * memoryStateCount + memoryState] = transitionGoal;
                         }
                     }
                 }
             }
-            return result;
         }
             
         template <typename ValueType, typename RewardModelType>
-        void SparseModelMemoryProduct<ValueType, RewardModelType>::computeReachableStates(std::vector<uint64_t> const& memorySuccessors, storm::storage::BitVector const& initialStates,                     boost::optional<storm::storage::Scheduler<ValueType>> const& scheduler) {
+        void SparseModelMemoryProduct<ValueType, RewardModelType>::computeReachableStates(storm::storage::BitVector const& initialStates) {
             // Explore the reachable states via DFS.
             // A state s on the stack corresponds to the model state (s / memoryStateCount) and memory state (s % memoryStateCount)
             reachableStates |= initialStates;
@@ -153,7 +176,7 @@ namespace storm {
         }
         
         template <typename ValueType, typename RewardModelType>
-        storm::storage::SparseMatrix<ValueType> SparseModelMemoryProduct<ValueType, RewardModelType>::buildDeterministicTransitionMatrix(std::vector<uint64_t> const& memorySuccessors) const {
+        storm::storage::SparseMatrix<ValueType> SparseModelMemoryProduct<ValueType, RewardModelType>::buildDeterministicTransitionMatrix() {
             uint64_t numResStates = reachableStates.getNumberOfSetBits();
             uint64_t numResTransitions = 0;
             for (auto const& stateIndex : reachableStates) {
@@ -168,7 +191,7 @@ namespace storm {
                 auto const& modelRow = model.getTransitionMatrix().getRow(modelState);
                 for (auto entryIt = modelRow.begin(); entryIt != modelRow.end(); ++entryIt) {
                     uint64_t transitionId = entryIt - model.getTransitionMatrix().begin();
-                    uint64_t const& successorMemoryState = memorySuccessors[transitionId * memoryStateCount + memoryState];
+                    uint64_t successorMemoryState = memorySuccessors[transitionId * memoryStateCount + memoryState];
                     builder.addNextValue(currentRow, getResultState(entryIt->getColumn(), successorMemoryState), entryIt->getValue());
                 }
                 ++currentRow;
@@ -178,7 +201,7 @@ namespace storm {
         }
         
         template <typename ValueType, typename RewardModelType>
-        storm::storage::SparseMatrix<ValueType> SparseModelMemoryProduct<ValueType, RewardModelType>::buildNondeterministicTransitionMatrix(std::vector<uint64_t> const& memorySuccessors) const {
+        storm::storage::SparseMatrix<ValueType> SparseModelMemoryProduct<ValueType, RewardModelType>::buildNondeterministicTransitionMatrix() {
             uint64_t numResStates = reachableStates.getNumberOfSetBits();
             uint64_t numResChoices = 0;
             uint64_t numResTransitions = 0;
@@ -200,7 +223,7 @@ namespace storm {
                     auto const& modelRow = model.getTransitionMatrix().getRow(modelRowIndex);
                     for (auto entryIt = modelRow.begin(); entryIt != modelRow.end(); ++entryIt) {
                         uint64_t transitionId = entryIt - model.getTransitionMatrix().begin();
-                        uint64_t const& successorMemoryState = memorySuccessors[transitionId * memoryStateCount + memoryState];
+                        uint64_t successorMemoryState = memorySuccessors[transitionId * memoryStateCount + memoryState];
                         builder.addNextValue(currentRow, getResultState(entryIt->getColumn(), successorMemoryState), entryIt->getValue());
                     }
                     ++currentRow;
@@ -211,7 +234,7 @@ namespace storm {
         }
         
         template <typename ValueType, typename RewardModelType>
-        storm::storage::SparseMatrix<ValueType> SparseModelMemoryProduct<ValueType, RewardModelType>::buildTransitionMatrixForScheduler(std::vector<uint64_t> const& memorySuccessors, storm::storage::Scheduler<ValueType> const& scheduler) const {
+        storm::storage::SparseMatrix<ValueType> SparseModelMemoryProduct<ValueType, RewardModelType>::buildTransitionMatrixForScheduler() {
             uint64_t numResStates = reachableStates.getNumberOfSetBits();
             uint64_t numResChoices = 0;
             uint64_t numResTransitions = 0;
@@ -219,7 +242,7 @@ namespace storm {
             for (auto const& stateIndex : reachableStates) {
                 uint64_t modelState = stateIndex / memoryStateCount;
                 uint64_t memoryState = stateIndex % memoryStateCount;
-                storm::storage::SchedulerChoice<ValueType> choice = scheduler.getChoice(modelState, memoryState);
+                storm::storage::SchedulerChoice<ValueType> choice = scheduler->getChoice(modelState, memoryState);
                 if (choice.isDefined()) {
                     ++numResChoices;
                     if (choice.isDeterministic()) {
@@ -256,14 +279,14 @@ namespace storm {
                 if (!hasTrivialNondeterminism) {
                     builder.newRowGroup(currentRow);
                 }
-                storm::storage::SchedulerChoice<ValueType> choice = scheduler.getChoice(modelState, memoryState);
+                storm::storage::SchedulerChoice<ValueType> choice = scheduler->getChoice(modelState, memoryState);
                 if (choice.isDefined()) {
                     if (choice.isDeterministic()) {
                         uint64_t modelRowIndex = model.getTransitionMatrix().getRowGroupIndices()[modelState] + choice.getDeterministicChoice();
                         auto const& modelRow = model.getTransitionMatrix().getRow(modelRowIndex);
                         for (auto entryIt = modelRow.begin(); entryIt != modelRow.end(); ++entryIt) {
                             uint64_t transitionId = entryIt - model.getTransitionMatrix().begin();
-                            uint64_t const& successorMemoryState = memorySuccessors[transitionId * memoryStateCount + memoryState];
+                            uint64_t successorMemoryState = memorySuccessors[transitionId * memoryStateCount + memoryState];
                             builder.addNextValue(currentRow, getResultState(entryIt->getColumn(), successorMemoryState), entryIt->getValue());
                         }
                     } else {
@@ -274,7 +297,7 @@ namespace storm {
                                 auto const& modelRow = model.getTransitionMatrix().getRow(modelRowIndex);
                                 for (auto entryIt = modelRow.begin(); entryIt != modelRow.end(); ++entryIt) {
                                     uint64_t transitionId = entryIt - model.getTransitionMatrix().begin();
-                                    uint64_t const& successorMemoryState = memorySuccessors[transitionId * memoryStateCount + memoryState];
+                                    uint64_t successorMemoryState = memorySuccessors[transitionId * memoryStateCount + memoryState];
                                     ValueType transitionValue = choiceIndex.second * entryIt->getValue();
                                     auto insertionRes = transitions.insert(std::make_pair(getResultState(entryIt->getColumn(), successorMemoryState), transitionValue));
                                     if (!insertionRes.second) {
@@ -293,7 +316,7 @@ namespace storm {
                         auto const& modelRow = model.getTransitionMatrix().getRow(modelRowIndex);
                         for (auto entryIt = modelRow.begin(); entryIt != modelRow.end(); ++entryIt) {
                             uint64_t transitionId = entryIt - model.getTransitionMatrix().begin();
-                            uint64_t const& successorMemoryState = memorySuccessors[transitionId * memoryStateCount + memoryState];
+                            uint64_t successorMemoryState = memorySuccessors[transitionId * memoryStateCount + memoryState];
                             builder.addNextValue(currentRow, getResultState(entryIt->getColumn(), successorMemoryState), entryIt->getValue());
                         }
                         ++currentRow;
@@ -305,7 +328,7 @@ namespace storm {
         }
         
         template <typename ValueType, typename RewardModelType>
-        storm::models::sparse::StateLabeling SparseModelMemoryProduct<ValueType, RewardModelType>::buildStateLabeling(storm::storage::SparseMatrix<ValueType> const& resultTransitionMatrix) const {
+        storm::models::sparse::StateLabeling SparseModelMemoryProduct<ValueType, RewardModelType>::buildStateLabeling(storm::storage::SparseMatrix<ValueType> const& resultTransitionMatrix) {
             uint64_t modelStateCount = model.getNumberOfStates();
             
             uint64_t numResStates = resultTransitionMatrix.getRowGroupCount();
@@ -316,10 +339,8 @@ namespace storm {
                     storm::storage::BitVector resLabeledStates(numResStates, false);
                     for (auto const& modelState : model.getStateLabeling().getStates(modelLabel)) {
                         for (uint64_t memoryState = 0; memoryState < memoryStateCount; ++memoryState) {
-                            uint64_t const& resState = getResultState(modelState, memoryState);
-                            // Check if the state exists in the result (i.e. if it is reachable)
-                            if (resState < numResStates) {
-                                resLabeledStates.set(resState, true);
+                            if (isStateReachable(modelState, memoryState)) {
+                                resLabeledStates.set(getResultState(modelState, memoryState), true);
                             }
                         }
                     }
@@ -331,20 +352,27 @@ namespace storm {
                 storm::storage::BitVector resLabeledStates(numResStates, false);
                 for (auto const& memoryState : memory.getStateLabeling().getStates(memoryLabel)) {
                     for (uint64_t modelState = 0; modelState < modelStateCount; ++modelState) {
-                        uint64_t const& resState = getResultState(modelState, memoryState);
-                        // Check if the state exists in the result (i.e. if it is reachable)
-                        if (resState < numResStates) {
-                            resLabeledStates.set(resState, true);
+                        if (isStateReachable(modelState, memoryState)) {
+                            resLabeledStates.set(getResultState(modelState, memoryState), true);
                         }
                     }
                 }
                 resultLabeling.addLabel(memoryLabel, std::move(resLabeledStates));
             }
+            
+            storm::storage::BitVector initialStates(numResStates, false);
+            auto memoryInitIt = memory.getInitialMemoryStates().begin();
+            for (auto const& modelInit : model.getInitialStates()) {
+                    initialStates.set(getResultState(modelInit, *memoryInitIt), true);
+                    ++memoryInitIt;
+                }
+            resultLabeling.addLabel("init", std::move(initialStates));
+            
             return resultLabeling;
         }
         
         template <typename ValueType, typename RewardModelType>
-        std::unordered_map<std::string, RewardModelType> SparseModelMemoryProduct<ValueType, RewardModelType>::buildRewardModels(storm::storage::SparseMatrix<ValueType> const& resultTransitionMatrix, std::vector<uint64_t> const& memorySuccessors, boost::optional<storm::storage::Scheduler<ValueType>> const& scheduler) const {
+        std::unordered_map<std::string, RewardModelType> SparseModelMemoryProduct<ValueType, RewardModelType>::buildRewardModels(storm::storage::SparseMatrix<ValueType> const& resultTransitionMatrix) {
             
             typedef typename RewardModelType::ValueType RewardValueType;
             
@@ -359,10 +387,8 @@ namespace storm {
                     for (auto const& modelStateReward : rewardModel.second.getStateRewardVector()) {
                         if (!storm::utility::isZero(modelStateReward)) {
                             for (uint64_t memoryState = 0; memoryState < memoryStateCount; ++memoryState) {
-                                uint64_t const& resState = getResultState(modelState, memoryState);
-                                // Check if the state exists in the result (i.e. if it is reachable)
-                                if (resState < numResStates) {
-                                    stateRewards.get()[resState] = modelStateReward;
+                                if (isStateReachable(modelState, memoryState)) {
+                                    stateRewards.get()[getResultState(modelState, memoryState)] = modelStateReward;
                                 }
                             }
                         }
@@ -381,14 +407,12 @@ namespace storm {
                             }
                             uint64_t rowOffset = modelRow - model.getTransitionMatrix().getRowGroupIndices()[modelState];
                             for (uint64_t memoryState = 0; memoryState < memoryStateCount; ++memoryState) {
-                                uint64_t const& resState = getResultState(modelState, memoryState);
-                                // Check if the state exists in the result (i.e. if it is reachable)
-                                if (resState < numResStates) {
+                                if (isStateReachable(modelState, memoryState)) {
                                     if (scheduler && scheduler->getChoice(modelState, memoryState).isDefined()) {
                                         ValueType factor = scheduler->getChoice(modelState, memoryState).getChoiceAsDistribution().getProbability(rowOffset);
-                                        stateActionRewards.get()[resultTransitionMatrix.getRowGroupIndices()[resState]] = factor * modelStateActionReward;
+                                        stateActionRewards.get()[resultTransitionMatrix.getRowGroupIndices()[getResultState(modelState, memoryState)]] = factor * modelStateActionReward;
                                     } else {
-                                        stateActionRewards.get()[resultTransitionMatrix.getRowGroupIndices()[resState] + rowOffset] = modelStateActionReward;
+                                        stateActionRewards.get()[resultTransitionMatrix.getRowGroupIndices()[getResultState(modelState, memoryState)] + rowOffset] = modelStateActionReward;
                                     }
                                 }
                             }
@@ -416,7 +440,7 @@ namespace storm {
                                             ++transitionEntryIt;
                                         }
                                         uint64_t transitionId = transitionEntryIt - model.getTransitionMatrix().begin();
-                                        uint64_t const& successorMemoryState = memorySuccessors[transitionId * memoryStateCount + memoryState];
+                                        uint64_t successorMemoryState = memorySuccessors[transitionId * memoryStateCount + memoryState];
                                         auto insertionRes = rewards.insert(std::make_pair(getResultState(rewardEntry.getColumn(), successorMemoryState), rewardEntry.getValue()));
                                         if (!insertionRes.second) {
                                             insertionRes.first->second += rewardEntry.getValue();
@@ -438,7 +462,7 @@ namespace storm {
                                             ++transitionEntryIt;
                                         }
                                         uint64_t transitionId = transitionEntryIt - model.getTransitionMatrix().begin();
-                                        uint64_t const& successorMemoryState = memorySuccessors[transitionId * memoryStateCount + memoryState];
+                                        uint64_t successorMemoryState = memorySuccessors[transitionId * memoryStateCount + memoryState];
                                         builder.addNextValue(resRowIndex, getResultState(rewardEntry.getColumn(), successorMemoryState), rewardEntry.getValue());
                                     }
                                 }
@@ -454,7 +478,7 @@ namespace storm {
         }
             
         template <typename ValueType, typename RewardModelType>
-        std::shared_ptr<storm::models::sparse::Model<ValueType, RewardModelType>> SparseModelMemoryProduct<ValueType, RewardModelType>::buildResult(storm::storage::SparseMatrix<ValueType>&& matrix, storm::models::sparse::StateLabeling&& labeling, std::unordered_map<std::string, RewardModelType>&& rewardModels, boost::optional<storm::storage::Scheduler<ValueType>> const& scheduler) const {
+        std::shared_ptr<storm::models::sparse::Model<ValueType, RewardModelType>> SparseModelMemoryProduct<ValueType, RewardModelType>::buildResult(storm::storage::SparseMatrix<ValueType>&& matrix, storm::models::sparse::StateLabeling&& labeling, std::unordered_map<std::string, RewardModelType>&& rewardModels) {
             storm::storage::sparse::ModelComponents<ValueType, RewardModelType> components (std::move(matrix), std::move(labeling), std::move(rewardModels));
             
             if (model.isOfType(storm::models::ModelType::Ctmc)) {
diff --git a/src/storm/storage/memorystructure/SparseModelMemoryProduct.h b/src/storm/storage/memorystructure/SparseModelMemoryProduct.h
index f32f092d1..45f4d2525 100644
--- a/src/storm/storage/memorystructure/SparseModelMemoryProduct.h
+++ b/src/storm/storage/memorystructure/SparseModelMemoryProduct.h
@@ -25,48 +25,63 @@ namespace storm {
         class SparseModelMemoryProduct {
         public:
             
+            // Constructs the product w.r.t. the given model and the given memory structure
             SparseModelMemoryProduct(storm::models::sparse::Model<ValueType, RewardModelType> const& sparseModel, storm::storage::MemoryStructure const& memoryStructure);
             
+            // Constructs the product w.r.t. the given model and the given (finite memory) scheduler.
+            SparseModelMemoryProduct(storm::models::sparse::Model<ValueType, RewardModelType> const& sparseModel, storm::storage::Scheduler<ValueType> const& scheduler);
+            
             // Enforces that the given model and memory state as well as the successor(s) are considered reachable -- even if they are not reachable from an initial state.
             void addReachableState(uint64_t const& modelState, uint64_t const& memoryState);
     
             // Enforces that every state is considered reachable. If this is set, the result has size #modelStates * #memoryStates
             void setBuildFullProduct();
             
-            // Invokes the building of the product under the specified scheduler (if given).
-            std::shared_ptr<storm::models::sparse::Model<ValueType, RewardModelType>> build(boost::optional<storm::storage::Scheduler<ValueType>> const& scheduler = boost::none);
+            // Returns true iff the given model and memory state is reachable in the product
+            bool isStateReachable(uint64_t const& modelState, uint64_t const& memoryState);
             
             // Retrieves the state of the resulting model that represents the given memory and model state.
-            // Should only be called AFTER calling build();
-            // An invalid index is returned, if the specifyied state does not exist (i.e., if it is not part of product).
-            uint64_t const& getResultState(uint64_t const& modelState, uint64_t const& memoryState) const;
+            // This method should only be called if the given state is reachable.
+            uint64_t const& getResultState(uint64_t const& modelState, uint64_t const& memoryState);
+            
+            // Invokes the building of the product under the specified scheduler (if given).
+            std::shared_ptr<storm::models::sparse::Model<ValueType, RewardModelType>> build();
             
             storm::models::sparse::Model<ValueType, RewardModelType> const& getOriginalModel() const;
             storm::storage::MemoryStructure const& getMemory() const;
+            
         private:
             
+            // Initializes auxiliary data for building the product
+            void initialize();
+            
             // Computes for each pair of memory state and model transition index the successor memory state
             // The resulting vector maps (modelTransition * memoryStateCount) + memoryState to the corresponding successor state of the memory structure
-            std::vector<uint64_t> computeMemorySuccessors() const;
+            void computeMemorySuccessors();
             
             // Computes the reachable states of the resulting model
-            void computeReachableStates(std::vector<uint64_t> const& memorySuccessors, storm::storage::BitVector const& initialStates, boost::optional<storm::storage::Scheduler<ValueType>> const& scheduler);
+            void computeReachableStates(storm::storage::BitVector const& initialStates);
             
             // Methods that build the model components
             // Matrix for deterministic models
-            storm::storage::SparseMatrix<ValueType> buildDeterministicTransitionMatrix(std::vector<uint64_t> const& memorySuccessors) const;
+            storm::storage::SparseMatrix<ValueType> buildDeterministicTransitionMatrix();
             // Matrix for nondeterministic models
-            storm::storage::SparseMatrix<ValueType> buildNondeterministicTransitionMatrix(std::vector<uint64_t> const& memorySuccessors) const;
+            storm::storage::SparseMatrix<ValueType> buildNondeterministicTransitionMatrix();
             // Matrix for models that consider a scheduler
-            storm::storage::SparseMatrix<ValueType> buildTransitionMatrixForScheduler(std::vector<uint64_t> const& memorySuccessors, storm::storage::Scheduler<ValueType> const& scheduler) const;
-            // State labeling. Note: DOES NOT ADD A LABEL FOR THE INITIAL STATES
-            storm::models::sparse::StateLabeling buildStateLabeling(storm::storage::SparseMatrix<ValueType> const& resultTransitionMatrix) const;
+            storm::storage::SparseMatrix<ValueType> buildTransitionMatrixForScheduler();
+            // State labeling.
+            storm::models::sparse::StateLabeling buildStateLabeling(storm::storage::SparseMatrix<ValueType> const& resultTransitionMatrix);
             // Reward models
-            std::unordered_map<std::string, RewardModelType> buildRewardModels(storm::storage::SparseMatrix<ValueType> const& resultTransitionMatrix, std::vector<uint64_t> const& memorySuccessors, boost::optional<storm::storage::Scheduler<ValueType>> const& scheduler) const;
+            std::unordered_map<std::string, RewardModelType> buildRewardModels(storm::storage::SparseMatrix<ValueType> const& resultTransitionMatrix);
             
             // Builds the resulting model
-            std::shared_ptr<storm::models::sparse::Model<ValueType, RewardModelType>> buildResult(storm::storage::SparseMatrix<ValueType>&& matrix, storm::models::sparse::StateLabeling&& labeling, std::unordered_map<std::string, RewardModelType>&& rewardModels, boost::optional<storm::storage::Scheduler<ValueType>> const& scheduler) const;
+            std::shared_ptr<storm::models::sparse::Model<ValueType, RewardModelType>> buildResult(storm::storage::SparseMatrix<ValueType>&& matrix, storm::models::sparse::StateLabeling&& labeling, std::unordered_map<std::string, RewardModelType>&& rewardModels);
+            
+            // Stores whether this builder has already been initialized.
+            bool isInitialized;
             
+            // stores the successor memory states for each transition in the product
+            std::vector<uint64_t> memorySuccessors;
             
             // Maps (modelState * memoryStateCount) + memoryState to the state in the result that represents (memoryState,modelState)
             std::vector<uint64_t> toResultStateMapping;
@@ -77,7 +92,9 @@ namespace storm {
             uint64_t const memoryStateCount;
             
             storm::models::sparse::Model<ValueType, RewardModelType> const& model;
+            boost::optional<storm::storage::MemoryStructure> localMemory;
             storm::storage::MemoryStructure const& memory;
+            boost::optional<storm::storage::Scheduler<ValueType> const&>  scheduler;
             
         };
     }