From fc97c1fc9dc53ac58a93413948a736c38427d142 Mon Sep 17 00:00:00 2001
From: TimQu <tim.quatmann@cs.rwth-aachen.de>
Date: Tue, 25 Apr 2017 17:35:57 +0200
Subject: [PATCH] introduced memory structure

---
 src/storm/storage/SparseMatrix.cpp            |   6 +-
 .../memorystructure/MemoryStructure.cpp       |  96 +++++
 .../storage/memorystructure/MemoryStructure.h |  65 ++++
 .../MemoryStructureBuilder.cpp                |  36 ++
 .../memorystructure/MemoryStructureBuilder.h  |  49 +++
 .../SparseModelMemoryProduct.cpp              | 327 ++++++++++++++++++
 .../SparseModelMemoryProduct.h                |  67 ++++
 7 files changed, 642 insertions(+), 4 deletions(-)
 create mode 100644 src/storm/storage/memorystructure/MemoryStructure.cpp
 create mode 100644 src/storm/storage/memorystructure/MemoryStructure.h
 create mode 100644 src/storm/storage/memorystructure/MemoryStructureBuilder.cpp
 create mode 100644 src/storm/storage/memorystructure/MemoryStructureBuilder.h
 create mode 100644 src/storm/storage/memorystructure/SparseModelMemoryProduct.cpp
 create mode 100644 src/storm/storage/memorystructure/SparseModelMemoryProduct.h

diff --git a/src/storm/storage/SparseMatrix.cpp b/src/storm/storage/SparseMatrix.cpp
index d71fa7ba4..da12ce7f8 100644
--- a/src/storm/storage/SparseMatrix.cpp
+++ b/src/storm/storage/SparseMatrix.cpp
@@ -14,6 +14,7 @@
 #include "storm/storage/BitVector.h"
 #include "storm/utility/constants.h"
 #include "storm/utility/ConstantsComparator.h"
+#include "storm/utility/vector.h"
 
 #include "storm/exceptions/InvalidStateException.h"
 #include "storm/exceptions/NotImplementedException.h"
@@ -546,10 +547,7 @@ namespace storm {
             // If there is no current row grouping, we need to create it.
             if (!this->rowGroupIndices) {
                 STORM_LOG_ASSERT(trivialRowGrouping, "Only trivial row-groupings can be constructed on-the-fly.");
-                this->rowGroupIndices = std::vector<index_type>(this->getRowCount() + 1);
-                for (uint_fast64_t group = 0; group <= this->getRowCount(); ++group) {
-                    this->rowGroupIndices.get()[group] = group;
-                }
+                this->rowGroupIndices = storm::utility::vector::buildVectorForRange(0, this->getRowGroupCount() + 1);
             }
             return rowGroupIndices.get();
         }
diff --git a/src/storm/storage/memorystructure/MemoryStructure.cpp b/src/storm/storage/memorystructure/MemoryStructure.cpp
new file mode 100644
index 000000000..bbf3137f5
--- /dev/null
+++ b/src/storm/storage/memorystructure/MemoryStructure.cpp
@@ -0,0 +1,96 @@
+#include "storm/storage/memorystructure/MemoryStructure.h"
+
+#include "storm/logic/Formulas.h"
+#include "storm/utility/macros.h"
+#include "storm/storage/memorystructure/SparseModelMemoryProduct.h"
+
+#include "storm/exceptions/InvalidOperationException.h"
+
+namespace storm {
+    namespace storage {
+        
+        MemoryStructure::MemoryStructure(TransitionMatrix const& transitionMatrix, storm::models::sparse::StateLabeling const& memoryStateLabeling) : transitions(transitionMatrix), stateLabeling(memoryStateLabeling) {
+            // intentionally left empty
+        }
+        
+        MemoryStructure::MemoryStructure(TransitionMatrix&& transitionMatrix, storm::models::sparse::StateLabeling&& memoryStateLabeling) : transitions(std::move(transitionMatrix)), stateLabeling(std::move(memoryStateLabeling)) {
+            // intentionally left empty
+        }
+            
+        MemoryStructure::TransitionMatrix const& MemoryStructure::getTransitionMatrix() const {
+            return this->transitions;
+        }
+        
+        storm::models::sparse::StateLabeling const& MemoryStructure::getStateLabeling() const {
+            return stateLabeling;
+        }
+            
+        MemoryStructure MemoryStructure::product(MemoryStructure const& rhs) const {
+            uint_fast64_t lhsNumStates = this->getTransitionMatrix().size();
+            uint_fast64_t rhsNumStates = rhs.getTransitionMatrix().size();
+            uint_fast64_t resNumStates = lhsNumStates * rhsNumStates;
+                
+            // Transition matrix
+            TransitionMatrix resultTransitions(resNumStates, std::vector<std::shared_ptr<storm::logic::Formula const>>(resNumStates));
+            uint_fast64_t resState = 0;
+            for (uint_fast64_t lhsState = 0; lhsState < lhsNumStates; ++lhsState) {
+                for (uint_fast64_t rhsState = 0; rhsState < rhsNumStates; ++rhsState) {
+                    assert (resState == (lhsState * rhsNumStates) + rhsState);
+                    auto resStateTransitions = resultTransitions[resState];
+                    for (uint_fast64_t lhsTransitionTarget = 0; lhsTransitionTarget < lhsNumStates; ++lhsTransitionTarget) {
+                        auto& lhsTransition = this->getTransitionMatrix()[lhsState][lhsTransitionTarget];
+                        if (lhsTransition) {
+                            for (uint_fast64_t rhsTransitionTarget = 0; rhsTransitionTarget < rhsNumStates; ++rhsTransitionTarget) {
+                                auto& rhsTransition = this->getTransitionMatrix()[rhsState][rhsTransitionTarget];
+                                if (rhsTransition) {
+                                    uint_fast64_t resTransitionTarget = (lhsTransitionTarget * rhsNumStates) + rhsTransitionTarget;
+                                    resStateTransitions[resTransitionTarget] = std::make_shared<storm::logic::BinaryBooleanStateFormula const>(storm::logic::BinaryBooleanStateFormula::OperatorType::And, lhsTransition,rhsTransition);
+                                }
+                            }
+                        }
+                    }
+                    ++resState;
+                }
+            }
+                
+            // State Labels
+            storm::models::sparse::StateLabeling resultLabeling(resNumStates);
+            for (std::string lhsLabel : this->getStateLabeling().getLabels()) {
+                storm::storage::BitVector const& lhsLabeledStates = this->getStateLabeling().getStates(lhsLabel);
+                storm::storage::BitVector resLabeledStates(resNumStates, false);
+                for (auto const& lhsState : lhsLabeledStates) {
+                    for (uint_fast64_t rhsState = 0; rhsState < rhsNumStates; ++rhsState) {
+                        resState = (lhsState * rhsNumStates) + rhsState;
+                        resLabeledStates.set(resState, true);
+                    }
+                }
+                resultLabeling.addLabel(lhsLabel, std::move(resLabeledStates));
+            }
+            for (std::string rhsLabel : rhs.getStateLabeling().getLabels()) {
+                STORM_LOG_THROW(!resultLabeling.containsLabel(rhsLabel), storm::exceptions::InvalidOperationException, "Failed to build the product of two memory structures: State labelings are not disjoint as both structures contain the label " << rhsLabel << ".");
+                storm::storage::BitVector const& rhsLabeledStates = rhs.getStateLabeling().getStates(rhsLabel);
+                storm::storage::BitVector resLabeledStates(resNumStates, false);
+                for (auto const& rhsState : rhsLabeledStates) {
+                    for (uint_fast64_t lhsState = 0; lhsState < lhsNumStates; ++lhsState) {
+                        resState = (lhsState * rhsNumStates) + rhsState;
+                        resLabeledStates.set(resState, true);
+                    }
+                }
+                resultLabeling.addLabel(rhsLabel, std::move(resLabeledStates));
+            }
+            return MemoryStructure(std::move(resultTransitions), std::move(resultLabeling));
+        }
+            
+        template <typename ValueType>
+        SparseModelMemoryProduct<ValueType> MemoryStructure::product(storm::models::sparse::Model<ValueType> const& sparseModel) const {
+            return SparseModelMemoryProduct<ValueType>(sparseModel, *this);
+        }
+        
+        template SparseModelMemoryProduct<double> MemoryStructure::product(storm::models::sparse::Model<double> const& sparseModel) const;
+        template SparseModelMemoryProduct<storm::RationalNumber> MemoryStructure::product(storm::models::sparse::Model<storm::RationalNumber> const& sparseModel) const;
+
+
+    }
+}
+
+
diff --git a/src/storm/storage/memorystructure/MemoryStructure.h b/src/storm/storage/memorystructure/MemoryStructure.h
new file mode 100644
index 000000000..651a54de7
--- /dev/null
+++ b/src/storm/storage/memorystructure/MemoryStructure.h
@@ -0,0 +1,65 @@
+#pragma once
+
+#include <vector>
+#include <memory>
+
+#include "storm/logic/Formula.h"
+#include "storm/models/sparse/StateLabeling.h"
+#include "storm/models/sparse/Model.h"
+
+namespace storm {
+    namespace storage {
+        
+        template <typename ValueType>
+        class SparseModelMemoryProduct;
+        
+        /*!
+         * This class represents a (deterministic) memory structure that can be used to encode certain events
+         * (such as reaching a set of goal states) into the state space of the model.
+         */
+        class MemoryStructure {
+        public:
+            
+            typedef std::vector<std::vector<std::shared_ptr<storm::logic::Formula const>>> TransitionMatrix;
+            
+            /*!
+             * Creates a memory structure with the given transition matrix and the given memory state labeling.
+             * The initial state is always the state with index 0.
+             * The transition matrix is assumed to contain propositional state formulas. The entry
+             * transitionMatrix[m][n] specifies the set of model states which trigger a transition from memory
+             * state m to memory state n.
+             * Transitions are assumed to be deterministic and complete, i.e., the formulas in
+             * transitionMatrix[m] form a partition of the state space of the considered model.
+             *
+             * @param transitionMatrix The transition matrix
+             * @param memoryStateLabeling A labeling of the memory states to specify, e.g., accepting states
+             */
+            MemoryStructure(TransitionMatrix const& transitionMatrix, storm::models::sparse::StateLabeling const& memoryStateLabeling);
+            MemoryStructure(TransitionMatrix&& transitionMatrix, storm::models::sparse::StateLabeling&& memoryStateLabeling);
+            
+            TransitionMatrix const& getTransitionMatrix() const;
+            storm::models::sparse::StateLabeling const& getStateLabeling() const;
+            
+            /*!
+             * Builds the product of this memory structure and the given memory structure.
+             * The resulting memory structure will have the state labels of both given structures.
+             * Throws an exception if the state labelings are not disjoint.
+             */
+            MemoryStructure product(MemoryStructure const& rhs) const;
+            
+            /*!
+             * Builds the product of this memory structure and the given sparse model.
+             * An exception is thrown if the state labelings of this memory structure and the given model are not disjoint.
+             */
+            template <typename ValueType>
+            SparseModelMemoryProduct<ValueType> product(storm::models::sparse::Model<ValueType> const& sparseModel) const;
+
+        private:
+            TransitionMatrix transitions;
+            storm::models::sparse::StateLabeling stateLabeling;
+        };
+        
+    }
+}
+
+
diff --git a/src/storm/storage/memorystructure/MemoryStructureBuilder.cpp b/src/storm/storage/memorystructure/MemoryStructureBuilder.cpp
new file mode 100644
index 000000000..b2b654a93
--- /dev/null
+++ b/src/storm/storage/memorystructure/MemoryStructureBuilder.cpp
@@ -0,0 +1,36 @@
+#include "storm/storage/memorystructure/MemoryStructureBuilder.h"
+
+#include "storm/logic/FragmentSpecification.h"
+#include "storm/utility/macros.h"
+
+#include "storm/exceptions/InvalidOperationException.h"
+
+namespace storm {
+    namespace storage {
+        
+        MemoryStructureBuilder::MemoryStructureBuilder(uint_fast64_t const& numberOfStates) : transitions(numberOfStates, std::vector<std::shared_ptr<storm::logic::Formula const>>(numberOfStates)), stateLabeling(numberOfStates) {
+            // Intentionally left empty
+        }
+            
+        void MemoryStructureBuilder::setTransition(uint_fast64_t const& startState, uint_fast64_t const& goalState, std::shared_ptr<storm::logic::Formula const> formula) {
+            STORM_LOG_THROW(startState < transitions.size(), storm::exceptions::InvalidOperationException, "Invalid index of start state: " << startState << ". There are only " << transitions.size() << " states in this memory structure.");
+            STORM_LOG_THROW(goalState < transitions.size(), storm::exceptions::InvalidOperationException, "Invalid index of goal state: " << startState << ". There are only " << transitions.size() << " states in this memory structure.");
+            STORM_LOG_THROW(formula->isInFragment(storm::logic::propositional()), storm::exceptions::InvalidOperationException, "The formula '" << *formula << "' is not propositional");
+            
+            transitions[startState][goalState] = formula;
+        }
+                        
+        void MemoryStructureBuilder::setLabel(uint_fast64_t const& state, std::string const& label) {
+            STORM_LOG_THROW(state < transitions.size(), storm::exceptions::InvalidOperationException, "Can not add label to state with index " << state << ". There are only " << transitions.size() << " states in this memory structure.");
+            if (!stateLabeling.containsLabel(label)) {
+                stateLabeling.addLabel(label);
+            }
+            stateLabeling.addLabelToState(label, state);
+        }
+            
+        MemoryStructure MemoryStructureBuilder::build() {
+            return MemoryStructure(std::move(transitions), std::move(stateLabeling));
+        }
+        
+    }
+}
diff --git a/src/storm/storage/memorystructure/MemoryStructureBuilder.h b/src/storm/storage/memorystructure/MemoryStructureBuilder.h
new file mode 100644
index 000000000..e87992d56
--- /dev/null
+++ b/src/storm/storage/memorystructure/MemoryStructureBuilder.h
@@ -0,0 +1,49 @@
+#pragma once
+
+#include <vector>
+#include <memory>
+
+#include "storm/logic/Formula.h"
+#include "storm/models/sparse/StateLabeling.h"
+#include "storm/storage/memorystructure/MemoryStructure.h"
+
+namespace storm {
+    namespace storage {
+        
+        
+        class MemoryStructureBuilder {
+        public:
+            
+            /*!
+             * Initializes a new builder for a memory structure
+             * @param numberOfStates The number of states the resulting memory structure should have
+             */
+            MemoryStructureBuilder(uint_fast64_t const& numberOfStates);
+            
+            /*!
+             * Specifies a transition. The formula should be a propositional formula
+             */
+            void setTransition(uint_fast64_t const& startState, uint_fast64_t const& goalState, std::shared_ptr<storm::logic::Formula const> formula);
+                        
+            /*!
+             * Sets a label to the given state.
+             */
+            void setLabel(uint_fast64_t const& state, std::string const& label);
+            
+            /*!
+             * Builds the memory structure.
+             * @note Calling this invalidates this builder.
+             * @note When calling this method, the specified transitions should be deterministic and complete. This is not checked.
+             */
+            MemoryStructure build();
+            
+
+        private:
+            MemoryStructure::TransitionMatrix transitions;
+            storm::models::sparse::StateLabeling stateLabeling;
+        };
+        
+    }
+}
+
+
diff --git a/src/storm/storage/memorystructure/SparseModelMemoryProduct.cpp b/src/storm/storage/memorystructure/SparseModelMemoryProduct.cpp
new file mode 100644
index 000000000..949781c92
--- /dev/null
+++ b/src/storm/storage/memorystructure/SparseModelMemoryProduct.cpp
@@ -0,0 +1,327 @@
+#include "storm/storage/memorystructure/SparseModelMemoryProduct.h"
+
+#include <boost/optional.hpp>
+
+#include "storm/models/sparse/Dtmc.h"
+#include "storm/models/sparse/Mdp.h"
+#include "storm/models/sparse/Ctmc.h"
+#include "storm/models/sparse/MarkovAutomaton.h"
+#include "storm/modelchecker/propositional/SparsePropositionalModelChecker.h"
+#include "storm/modelchecker/results/ExplicitQualitativeCheckResult.h"
+#include "storm/utility/constants.h"
+#include "storm/utility/macros.h"
+
+#include "storm/exceptions/InvalidOperationException.h"
+
+namespace storm {
+    namespace storage {
+
+        template <typename ValueType>
+        SparseModelMemoryProduct<ValueType>::SparseModelMemoryProduct(storm::models::sparse::Model<ValueType> const& sparseModel, storm::storage::MemoryStructure const& memoryStructure) : model(sparseModel), memory(memoryStructure) {
+            // intentionally left empty
+        }
+            
+        template <typename ValueType>
+        std::shared_ptr<storm::models::sparse::Model<ValueType>> SparseModelMemoryProduct<ValueType>::build() {
+            
+            std::vector<uint_fast64_t> memorySuccessors = computeMemorySuccessors();
+            storm::storage::BitVector reachableStates = computeReachableStates(memorySuccessors);
+            
+            // Compute the mapping to the states of the result
+            uint_fast64_t reachableStateCount = 0;
+            toResultStateMapping = std::vector<uint_fast64_t> (model.getNumberOfStates() * memory.getTransitionMatrix().size(), std::numeric_limits<uint_fast64_t>::max());
+            for (auto const& reachableState : reachableStates) {
+                toResultStateMapping[reachableState] = reachableStateCount;
+                ++reachableStateCount;
+            }
+                
+            // Build the model components
+            storm::storage::SparseMatrix<ValueType> transitionMatrix = model.getTransitionMatrix().hasTrivialRowGrouping() ?
+                                                                       buildDeterministicTransitionMatrix(reachableStates, memorySuccessors) :
+                                                                       buildNondeterministicTransitionMatrix(reachableStates, memorySuccessors);
+            storm::models::sparse::StateLabeling labeling = buildStateLabeling(transitionMatrix);
+            std::unordered_map<std::string, storm::models::sparse::StandardRewardModel<ValueType>> rewardModels = buildRewardModels(transitionMatrix, memorySuccessors);
+
+            return buildResult(std::move(transitionMatrix), std::move(labeling), std::move(rewardModels));
+
+        }
+            
+        template <typename ValueType>
+        uint_fast64_t const& SparseModelMemoryProduct<ValueType>::getResultState(uint_fast64_t const& modelState, uint_fast64_t const& memoryState) {
+                return toResultStateMapping[modelState * memory.getTransitionMatrix().size() + memoryState];
+        }
+            
+            
+        template <typename ValueType>
+        std::vector<uint_fast64_t> SparseModelMemoryProduct<ValueType>::computeMemorySuccessors() const {
+            uint_fast64_t modelStateCount = model.getNumberOfStates();
+            uint_fast64_t memoryStateCount = memory.getTransitionMatrix().size();
+            std::vector<uint_fast64_t> result(modelStateCount * memoryStateCount, std::numeric_limits<uint_fast64_t>::max());
+            
+            storm::modelchecker::SparsePropositionalModelChecker<storm::models::sparse::Model<ValueType>> mc(model);
+            for (uint_fast64_t memoryState = 0; memoryState < memoryStateCount; ++memoryState) {
+                for (uint_fast64_t transitionGoal = 0; transitionGoal < memoryStateCount; ++transitionGoal) {
+                    auto const& transition = memory.getTransitionMatrix()[memoryState][transitionGoal];
+                    if (transition) {
+                        storm::storage::BitVector const& modelStates = mc.check(*transition)->asExplicitQualitativeCheckResult().getTruthValuesVector();
+                        for (auto const& modelState : modelStates) {
+                            result[modelState * memoryStateCount + memoryState] = transitionGoal;
+                        }
+                    }
+                }
+            }
+            return result;
+        }
+            
+        template <typename ValueType>
+        storm::storage::BitVector SparseModelMemoryProduct<ValueType>::computeReachableStates(std::vector<uint_fast64_t> const& memorySuccessors) const {
+            uint_fast64_t modelStateCount = model.getNumberOfStates();
+            uint_fast64_t memoryStateCount = memory.getTransitionMatrix().size();
+            // Explore the reachable states via with DFS.
+            // A state s on the stack corresponds to the model state (s / memoryStateCount) and memory state (s % memoryStateCount)
+            storm::storage::BitVector reachableStates(modelStateCount * memoryStateCount, false);
+            std::vector<uint_fast64_t> stack;
+            for (auto const& modelInit : model.getInitialStates()) {
+                // Note: 0 is always the initial state of a memory structure.
+                uint_fast64_t stateIndex = modelInit * memoryStateCount + memorySuccessors[modelInit * memoryStateCount];
+                reachableStates.set(stateIndex, true);
+                stack.push_back(stateIndex);
+            }
+            while (!stack.empty()) {
+                uint_fast64_t stateIndex = stack.back();
+                stack.pop_back();
+                uint_fast64_t modelState = stateIndex / memoryStateCount;
+                uint_fast64_t memoryState = stateIndex % memoryStateCount;
+                
+                for (auto const& modelTransition : model.getTransitionMatrix().getRowGroup(modelState)) {
+                    if (!storm::utility::isZero(modelTransition.getValue())) {
+                        uint_fast64_t successorModelState = modelTransition.getColumn();
+                        uint_fast64_t successorMemoryState = memorySuccessors[successorModelState * memoryStateCount + memoryState];
+                        uint_fast64_t successorStateIndex = successorModelState * memoryStateCount + successorMemoryState;
+                        if (!reachableStates.get(successorStateIndex)) {
+                            reachableStates.set(successorStateIndex, true);
+                            stack.push_back(successorStateIndex);
+                        }
+                    }
+                }
+            }
+            return reachableStates;
+        }
+            
+        template <typename ValueType>
+        storm::storage::SparseMatrix<ValueType> SparseModelMemoryProduct<ValueType>::buildDeterministicTransitionMatrix(storm::storage::BitVector const& reachableStates, std::vector<uint_fast64_t> const& memorySuccessors) const {
+            uint_fast64_t memoryStateCount = memory.getTransitionMatrix().size();
+            uint_fast64_t numResStates = reachableStates.getNumberOfSetBits();
+            uint_fast64_t numResTransitions = 0;
+            for (auto const& stateIndex : reachableStates) {
+                numResTransitions += model.getTransitionMatrix().getRow(stateIndex / memoryStateCount).getNumberOfEntries();
+            }
+            
+            storm::storage::SparseMatrixBuilder<ValueType> builder(numResStates, numResStates, numResTransitions, true);
+            uint_fast64_t currentRow = 0;
+            for (auto const& stateIndex : reachableStates) {
+                uint_fast64_t modelState = stateIndex / memoryStateCount;
+                uint_fast64_t memoryState = stateIndex % memoryStateCount;
+                for (auto const& entry : model.getTransitionMatrix().getRow(modelState)) {
+                    uint_fast64_t const& successorMemoryState = memorySuccessors[entry.getColumn() * memoryStateCount + memoryState];
+                    uint_fast64_t transitionTargetStateIndex = entry.getColumn() * memoryStateCount + successorMemoryState;
+                    builder.addNextValue(currentRow, toResultStateMapping[transitionTargetStateIndex], entry.getValue());
+                }
+                ++currentRow;
+            }
+            
+            return builder.build();
+        }
+        
+        template <typename ValueType>
+        storm::storage::SparseMatrix<ValueType> SparseModelMemoryProduct<ValueType>::buildNondeterministicTransitionMatrix(storm::storage::BitVector const& reachableStates, std::vector<uint_fast64_t> const& memorySuccessors) const {
+            uint_fast64_t memoryStateCount = memory.getTransitionMatrix().size();
+            uint_fast64_t numResStates = reachableStates.getNumberOfSetBits();
+            uint_fast64_t numResChoices = 0;
+            uint_fast64_t numResTransitions = 0;
+            for (auto const& stateIndex : reachableStates) {
+                uint_fast64_t modelState = stateIndex / memoryStateCount;
+                for (uint_fast64_t modelRow = model.getTransitionMatrix().getRowGroupIndices()[modelState]; modelRow < model.getTransitionMatrix().getRowGroupIndices()[modelState + 1]; ++modelRow) {
+                    ++numResChoices;
+                    numResTransitions += model.getTransitionMatrix().getRow(modelRow).getNumberOfEntries();
+                }
+            }
+            
+            storm::storage::SparseMatrixBuilder<ValueType> builder(numResChoices, numResStates, numResTransitions, true, true, numResStates);
+            uint_fast64_t currentRow = 0;
+            for (auto const& stateIndex : reachableStates) {
+                uint_fast64_t modelState = stateIndex / memoryStateCount;
+                uint_fast64_t memoryState = stateIndex % memoryStateCount;
+                builder.newRowGroup(currentRow);
+                for (uint_fast64_t modelRow = model.getTransitionMatrix().getRowGroupIndices()[modelState]; modelRow < model.getTransitionMatrix().getRowGroupIndices()[modelState + 1]; ++modelRow) {
+                    for (auto const& entry : model.getTransitionMatrix().getRow(modelRow)) {
+                        uint_fast64_t const& successorMemoryState = memorySuccessors[entry.getColumn() * memoryStateCount + memoryState];
+                        uint_fast64_t transitionTargetStateIndex = entry.getColumn() * memoryStateCount + successorMemoryState;
+                        builder.addNextValue(currentRow, toResultStateMapping[transitionTargetStateIndex], entry.getValue());
+                    }
+                    ++currentRow;
+                }
+            }
+            
+            return builder.build();
+        }
+        
+        template <typename ValueType>
+        storm::models::sparse::StateLabeling SparseModelMemoryProduct<ValueType>::buildStateLabeling(storm::storage::SparseMatrix<ValueType> const& resultTransitionMatrix) const {
+            uint_fast64_t modelStateCount = model.getNumberOfStates();
+            uint_fast64_t memoryStateCount = memory.getTransitionMatrix().size();
+            
+            uint_fast64_t numResStates = resultTransitionMatrix.getRowGroupCount();
+            storm::models::sparse::StateLabeling resultLabeling(numResStates);
+            
+            for (std::string modelLabel : model.getStateLabeling().getLabels()) {
+                storm::storage::BitVector const& modelLabeledStates = model.getStateLabeling().getStates(modelLabel);
+                storm::storage::BitVector resLabeledStates(numResStates, false);
+                for (auto const& modelState : modelLabeledStates) {
+                    for (uint_fast64_t memoryState = 0; memoryState < memoryStateCount; ++memoryState) {
+                        uint_fast64_t const& resState = toResultStateMapping[modelState * memoryStateCount + memoryState];
+                        // Check if the state exists in the result (i.e. if it is reachable)
+                        if (resState < numResStates) {
+                            resLabeledStates.set(resState, true);
+                        }
+                    }
+                }
+                resultLabeling.addLabel(modelLabel, std::move(resLabeledStates));
+            }
+            for (std::string memoryLabel : memory.getStateLabeling().getLabels()) {
+                STORM_LOG_THROW(!resultLabeling.containsLabel(memoryLabel), storm::exceptions::InvalidOperationException, "Failed to build the product of model and memory structure: State labelings are not disjoint as both structures contain the label " << memoryLabel << ".");
+                storm::storage::BitVector const& memoryLabeledStates = memory.getStateLabeling().getStates(memoryLabel);
+                storm::storage::BitVector resLabeledStates(numResStates, false);
+                for (auto const& memoryState : memoryLabeledStates) {
+                    for (uint_fast64_t modelState = 0; modelState < modelStateCount; ++modelState) {
+                        uint_fast64_t const& resState = toResultStateMapping[modelState * memoryStateCount + memoryState];
+                        // Check if the state exists in the result (i.e. if it is reachable)
+                        if (resState < numResStates) {
+                            resLabeledStates.set(resState, true);
+                        }
+                    }
+                }
+                resultLabeling.addLabel(memoryLabel, std::move(resLabeledStates));
+            }
+            return resultLabeling;
+        }
+        
+        template <typename ValueType>
+        std::unordered_map<std::string, storm::models::sparse::StandardRewardModel<ValueType>> SparseModelMemoryProduct<ValueType>::buildRewardModels(storm::storage::SparseMatrix<ValueType> const& resultTransitionMatrix, std::vector<uint_fast64_t> const& memorySuccessors) const {
+            std::unordered_map<std::string, storm::models::sparse::StandardRewardModel<ValueType>> result;
+            uint_fast64_t memoryStateCount = memory.getTransitionMatrix().size();
+            uint_fast64_t numResStates = resultTransitionMatrix.getRowGroupCount();
+
+            for (auto const& rewardModel : model.getRewardModels()) {
+                boost::optional<std::vector<ValueType>> stateRewards;
+                if (rewardModel.second.hasStateRewards()) {
+                    stateRewards = std::vector<ValueType>(numResStates, storm::utility::zero<ValueType>());
+                    uint_fast64_t modelState = 0;
+                    for (auto const& modelStateReward : rewardModel.second.getStateRewardVector()) {
+                        if (!storm::utility::isZero(modelStateReward)) {
+                            for (uint_fast64_t memoryState = 0; memoryState < memoryStateCount; ++memoryState) {
+                                uint_fast64_t const& resState = toResultStateMapping[modelState * memoryStateCount + memoryState];
+                                // Check if the state exists in the result (i.e. if it is reachable)
+                                if (resState < numResStates) {
+                                    stateRewards.get()[resState] = modelStateReward;
+                                }
+                            }
+                        }
+                        ++modelState;
+                    }
+                }
+                boost::optional<std::vector<ValueType>> stateActionRewards;
+                if (rewardModel.second.hasStateActionRewards()) {
+                    stateActionRewards = std::vector<ValueType>(resultTransitionMatrix.getRowCount(), storm::utility::zero<ValueType>());
+                    uint_fast64_t modelState = 0;
+                    uint_fast64_t modelRow = 0;
+                    for (auto const& modelStateActionReward : rewardModel.second.getStateActionRewardVector()) {
+                        if (!storm::utility::isZero(modelStateActionReward)) {
+                            while (modelRow >= model.getTransitionMatrix().getRowGroupIndices()[modelState + 1]) {
+                                ++modelState;
+                            }
+                            uint_fast64_t rowOffset = modelRow - model.getTransitionMatrix().getRowGroupIndices()[modelState];
+                            for (uint_fast64_t memoryState = 0; memoryState < memoryStateCount; ++memoryState) {
+                                uint_fast64_t const& resState = toResultStateMapping[modelState * memoryStateCount + memoryState];
+                                // Check if the state exists in the result (i.e. if it is reachable)
+                                if (resState < numResStates) {
+                                    stateActionRewards.get()[resultTransitionMatrix.getRowGroupIndices()[resState] + rowOffset] = modelStateActionReward;
+                                }
+                            }
+                        }
+                        ++modelRow;
+                    }
+                }
+                boost::optional<storm::storage::SparseMatrix<ValueType>> transitionRewards;
+                if (rewardModel.second.hasTransitionRewards()) {
+                    storm::storage::SparseMatrixBuilder<ValueType> builder(resultTransitionMatrix.getRowCount(), resultTransitionMatrix.getColumnCount());
+                    uint_fast64_t stateIndex = 0;
+                    for (auto const& resState : toResultStateMapping) {
+                        if (resState < numResStates) {
+                            uint_fast64_t modelState = stateIndex / memoryStateCount;
+                            uint_fast64_t memoryState = stateIndex % memoryStateCount;
+                            uint_fast64_t rowGroupSize = resultTransitionMatrix.getRowGroupSize(resState);
+                            for (uint_fast64_t rowOffset = 0; rowOffset < rowGroupSize; ++rowOffset) {
+                                uint_fast64_t resRowIndex = resultTransitionMatrix.getRowGroupIndices()[resState] + rowOffset;
+                                uint_fast64_t modelRowIndex = model.getTransitionMatrix().getRowGroupIndices()[modelState] + rowOffset;
+                                for (auto const& entry : rewardModel.second.getTransitionRewardMatrix().getRow(modelRowIndex)) {
+                                    uint_fast64_t const& successorMemoryState = memorySuccessors[entry.getColumn() * memoryStateCount + memoryState];
+                                    uint_fast64_t transitionTargetStateIndex = entry.getColumn() * memoryStateCount + successorMemoryState;
+                                    builder.addNextValue(resRowIndex, toResultStateMapping[transitionTargetStateIndex], entry.getValue());
+                                }
+                            }
+                        }
+                        ++stateIndex;
+                    }
+                    transitionRewards = builder.build();
+                }
+                result.insert(std::make_pair(rewardModel.first, storm::models::sparse::StandardRewardModel<ValueType>(std::move(stateRewards), std::move(stateActionRewards), std::move(transitionRewards))));
+            }
+            return result;
+        }
+            
+        template <typename ValueType>
+        std::shared_ptr<storm::models::sparse::Model<ValueType>> SparseModelMemoryProduct<ValueType>::buildResult(storm::storage::SparseMatrix<ValueType>&& matrix, storm::models::sparse::StateLabeling&& labeling, std::unordered_map<std::string, storm::models::sparse::StandardRewardModel<ValueType>>&& rewardModels) const {
+            switch (model.getType()) {
+                case storm::models::ModelType::Dtmc:
+                    return std::make_shared<storm::models::sparse::Dtmc<ValueType>> (std::move(matrix), std::move(labeling), std::move(rewardModels));
+                case storm::models::ModelType::Mdp:
+                    return std::make_shared<storm::models::sparse::Mdp<ValueType>> (std::move(matrix), std::move(labeling), std::move(rewardModels));
+                case storm::models::ModelType::Ctmc:
+                    return std::make_shared<storm::models::sparse::Ctmc<ValueType>> (std::move(matrix), std::move(labeling), std::move(rewardModels));
+                case storm::models::ModelType::MarkovAutomaton:
+                {
+                    // We also need to translate the exit rates and the Markovian states
+                    uint_fast64_t numResStates = matrix.getRowGroupCount();
+                    uint_fast64_t memoryStateCount = memory.getTransitionMatrix().size();
+                    std::vector<ValueType> resultExitRates;
+                    resultExitRates.reserve(matrix.getRowGroupCount());
+                    storm::storage::BitVector resultMarkovianStates(numResStates, false);
+                    auto const& modelExitRates = dynamic_cast<storm::models::sparse::MarkovAutomaton<ValueType> const&>(model).getExitRates();
+                    auto const& modelMarkovianStates = dynamic_cast<storm::models::sparse::MarkovAutomaton<ValueType> const&>(model).getMarkovianStates();
+                    
+                    uint_fast64_t stateIndex = 0;
+                    for (auto const& resState : toResultStateMapping) {
+                        if (resState < numResStates) {
+                            assert(resState == resultExitRates.size());
+                            uint_fast64_t modelState = stateIndex / memoryStateCount;
+                            resultExitRates.push_back(modelExitRates[modelState]);
+                            if (modelMarkovianStates.get(modelState)) {
+                                resultMarkovianStates.set(resState, true);
+                            }
+                        }
+                    }
+                    ++stateIndex;
+                    return std::make_shared<storm::models::sparse::MarkovAutomaton<ValueType>> (std::move(matrix), std::move(labeling), std::move(resultMarkovianStates), std::move(resultExitRates), true, std::move(rewardModels));
+                }
+            }
+        }
+        
+        template  class SparseModelMemoryProduct<double>;
+        template  class SparseModelMemoryProduct<storm::RationalNumber>;
+        template  class SparseModelMemoryProduct<storm::RationalFunction>;
+            
+    }
+}
+
+
diff --git a/src/storm/storage/memorystructure/SparseModelMemoryProduct.h b/src/storm/storage/memorystructure/SparseModelMemoryProduct.h
new file mode 100644
index 000000000..92baea200
--- /dev/null
+++ b/src/storm/storage/memorystructure/SparseModelMemoryProduct.h
@@ -0,0 +1,67 @@
+#pragma once
+
+#include <memory>
+#include <string>
+#include <unordered_map>
+
+#include "storm/storage/BitVector.h"
+#include "storm/storage/memorystructure/MemoryStructure.h"
+#include "storm/models/sparse/Model.h"
+#include "storm/models/sparse/StandardRewardModel.h"
+
+namespace storm {
+    namespace storage {
+        /*!
+         * This class builds the product of the given sparse model and the given memory structure.
+         * This is similar to the well-known product of a model with a deterministic rabin automaton.
+         * Note: we already do one memory-structure-step for the initial state, i.e., if s is an initial state of
+         * the given model and s satisfies memoryStructure.getTransitionMatrix[0][n], then (s,n) will be the corresponding
+         * initial state of the resulting model.
+         *
+         * The states of the resulting sparse model will have the original state labels plus the labels of this
+         * memory structure.
+         * An exception is thrown if the state labelings are not disjoint.
+         */
+        template <typename ValueType>
+        class SparseModelMemoryProduct {
+        public:
+            
+            SparseModelMemoryProduct(storm::models::sparse::Model<ValueType> const& sparseModel, storm::storage::MemoryStructure const& memoryStructure);
+            
+            // Invokes the building of the product
+            std::shared_ptr<storm::models::sparse::Model<ValueType>> build();
+            
+            // Retrieves the state of the resulting model that represents the given memory and model state
+            uint_fast64_t const& getResultState(uint_fast64_t const& modelState, uint_fast64_t const& memoryState);
+            
+        private:
+            
+            // Computes for each pair of memory and model state the successor memory state
+            // The resulting vector maps (modelState * memoryStateCount) + memoryState to the corresponding successor state of the memory structure
+            std::vector<uint_fast64_t> computeMemorySuccessors() const;
+            
+            // Computes the reachable states of the resulting model
+            storm::storage::BitVector computeReachableStates(std::vector<uint_fast64_t> const& memorySuccessors) const;
+            
+            // Methods that build the model components
+            storm::storage::SparseMatrix<ValueType> buildDeterministicTransitionMatrix(storm::storage::BitVector const& reachableStates, std::vector<uint_fast64_t> const& memorySuccessors) const;
+            storm::storage::SparseMatrix<ValueType> buildNondeterministicTransitionMatrix(storm::storage::BitVector const& reachableStates, std::vector<uint_fast64_t> const& memorySuccessors) const;
+            storm::models::sparse::StateLabeling buildStateLabeling(storm::storage::SparseMatrix<ValueType> const& resultTransitionMatrix) const;
+            std::unordered_map<std::string, storm::models::sparse::StandardRewardModel<ValueType>> buildRewardModels(storm::storage::SparseMatrix<ValueType> const& resultTransitionMatrix, std::vector<uint_fast64_t> const& memorySuccessors) const;
+            
+            // Builds the resulting model
+            std::shared_ptr<storm::models::sparse::Model<ValueType>> buildResult(storm::storage::SparseMatrix<ValueType>&& matrix, storm::models::sparse::StateLabeling&& labeling, std::unordered_map<std::string, storm::models::sparse::StandardRewardModel<ValueType>>&& rewardModels) const;
+            
+            
+            // Maps (modelState * memoryStateCount) + memoryState to the state in the result that represents (memoryState,modelState)
+            std::vector<uint_fast64_t> toResultStateMapping;
+            
+
+            storm::models::sparse::Model<ValueType> const& model;
+            storm::storage::MemoryStructure const& memory;
+            
+        };
+    }
+}
+
+