From 03707f0234f1d6b3f77ca080fbe0149afecccf83 Mon Sep 17 00:00:00 2001
From: dehnert <dehnert@cs.rwth-aachen.de>
Date: Tue, 22 May 2018 21:49:39 +0200
Subject: [PATCH] first step for fixing MEC decomposition: making SCC
 decomposition accept a bit vector of subsystem choices

---
 ...tronglyConnectedComponentDecomposition.cpp | 121 ++++++++++++------
 .../StronglyConnectedComponentDecomposition.h |  24 +++-
 2 files changed, 103 insertions(+), 42 deletions(-)

diff --git a/src/storm/storage/StronglyConnectedComponentDecomposition.cpp b/src/storm/storage/StronglyConnectedComponentDecomposition.cpp
index f6567985c..9eae9d05a 100644
--- a/src/storm/storage/StronglyConnectedComponentDecomposition.cpp
+++ b/src/storm/storage/StronglyConnectedComponentDecomposition.cpp
@@ -23,30 +23,35 @@ namespace storm {
         template <typename RewardModelType>
         StronglyConnectedComponentDecomposition<ValueType>::StronglyConnectedComponentDecomposition(storm::models::sparse::Model<ValueType, RewardModelType> const& model, StateBlock const& block, bool dropNaiveSccs, bool onlyBottomSccs) {
             storm::storage::BitVector subsystem(model.getNumberOfStates(), block.begin(), block.end());
-            performSccDecomposition(model.getTransitionMatrix(), subsystem, dropNaiveSccs, onlyBottomSccs);
+            performSccDecomposition(model.getTransitionMatrix(), &subsystem, nullptr, dropNaiveSccs, onlyBottomSccs);
         }
         
         template <typename ValueType>
         template <typename RewardModelType>
         StronglyConnectedComponentDecomposition<ValueType>::StronglyConnectedComponentDecomposition(storm::models::sparse::Model<ValueType, RewardModelType> const& model, storm::storage::BitVector const& subsystem, bool dropNaiveSccs, bool onlyBottomSccs) {
-            performSccDecomposition(model.getTransitionMatrix(), subsystem, dropNaiveSccs, onlyBottomSccs);
+            performSccDecomposition(model.getTransitionMatrix(), &subsystem, nullptr, dropNaiveSccs, onlyBottomSccs);
         }
         
         template <typename ValueType>
         StronglyConnectedComponentDecomposition<ValueType>::StronglyConnectedComponentDecomposition(storm::storage::SparseMatrix<ValueType> const& transitionMatrix, StateBlock const& block, bool dropNaiveSccs, bool onlyBottomSccs) {
             storm::storage::BitVector subsystem(transitionMatrix.getRowGroupCount(), block.begin(), block.end());
-            performSccDecomposition(transitionMatrix, subsystem, dropNaiveSccs, onlyBottomSccs);
+            performSccDecomposition(transitionMatrix, &subsystem, nullptr, dropNaiveSccs, onlyBottomSccs);
         }
         
         
         template <typename ValueType>
         StronglyConnectedComponentDecomposition<ValueType>::StronglyConnectedComponentDecomposition(storm::storage::SparseMatrix<ValueType> const& transitionMatrix, bool dropNaiveSccs, bool onlyBottomSccs) {
-            performSccDecomposition(transitionMatrix, storm::storage::BitVector(transitionMatrix.getRowGroupCount(), true), dropNaiveSccs, onlyBottomSccs);
+            performSccDecomposition(transitionMatrix, nullptr, nullptr, dropNaiveSccs, onlyBottomSccs);
         }
 
         template <typename ValueType>
         StronglyConnectedComponentDecomposition<ValueType>::StronglyConnectedComponentDecomposition(storm::storage::SparseMatrix<ValueType> const& transitionMatrix, storm::storage::BitVector const& subsystem, bool dropNaiveSccs, bool onlyBottomSccs) {
-            performSccDecomposition(transitionMatrix, subsystem, dropNaiveSccs, onlyBottomSccs);
+            performSccDecomposition(transitionMatrix, &subsystem, nullptr, dropNaiveSccs, onlyBottomSccs);
+        }
+        
+        template <typename ValueType>
+        StronglyConnectedComponentDecomposition<ValueType>::StronglyConnectedComponentDecomposition(storm::storage::SparseMatrix<ValueType> const& transitionMatrix, storm::storage::BitVector const& subsystem, storm::storage::BitVector const& choices, bool dropNaiveSccs, bool onlyBottomSccs) {
+            performSccDecomposition(transitionMatrix, &subsystem, &choices, dropNaiveSccs, onlyBottomSccs);
         }
         
         template <typename ValueType>
@@ -72,7 +77,10 @@ namespace storm {
         }
 
         template <typename ValueType>
-        void StronglyConnectedComponentDecomposition<ValueType>::performSccDecomposition(storm::storage::SparseMatrix<ValueType> const& transitionMatrix, storm::storage::BitVector const& subsystem, bool dropNaiveSccs, bool onlyBottomSccs) {
+        void StronglyConnectedComponentDecomposition<ValueType>::performSccDecomposition(storm::storage::SparseMatrix<ValueType> const& transitionMatrix, storm::storage::BitVector const* subsystem, storm::storage::BitVector const* choices, bool dropNaiveSccs, bool onlyBottomSccs) {
+            
+            STORM_LOG_ASSERT(!choices || subsystem, "Expecting subsystem if choices are given.");
+            
             uint_fast64_t numberOfStates = transitionMatrix.getRowGroupCount();
 
             // Set up the environment of the algorithm.
@@ -94,16 +102,30 @@ namespace storm {
             
             // Start the search for SCCs from every state in the block.
             uint_fast64_t currentIndex = 0;
-            for (auto state : subsystem) {
-                if (!hasPreorderNumber.get(state)) {
-                    performSccDecompositionGCM(transitionMatrix, state, statesWithSelfLoop, subsystem, currentIndex, hasPreorderNumber, preorderNumbers, s, p, stateHasScc, stateToSccMapping, sccCount);
+            if (subsystem) {
+                for (auto state : *subsystem) {
+                    if (!hasPreorderNumber.get(state)) {
+                        performSccDecompositionGCM(transitionMatrix, state, statesWithSelfLoop, subsystem, choices, currentIndex, hasPreorderNumber, preorderNumbers, s, p, stateHasScc, stateToSccMapping, sccCount);
+                    }
+                }
+            } else {
+                for (uint64_t state = 0; state < transitionMatrix.getRowGroupCount(); ++state) {
+                    if (!hasPreorderNumber.get(state)) {
+                        performSccDecompositionGCM(transitionMatrix, state, statesWithSelfLoop, subsystem, choices, currentIndex, hasPreorderNumber, preorderNumbers, s, p, stateHasScc, stateToSccMapping, sccCount);
+                    }
                 }
             }
 
             // After we obtained the state-to-SCC mapping, we build the actual blocks.
             this->blocks.resize(sccCount);
-            for (auto state : subsystem) {
-                this->blocks[stateToSccMapping[state]].insert(state);
+            if (subsystem) {
+                for (auto state : *subsystem) {
+                    this->blocks[stateToSccMapping[state]].insert(state);
+                }
+            } else {
+                for (uint64_t state = 0; state < transitionMatrix.getRowGroupCount(); ++state) {
+                    this->blocks[stateToSccMapping[state]].insert(state);
+                }
             }
             
             // Now flag all trivial SCCs as such.
@@ -132,13 +154,34 @@ namespace storm {
                 
                 // If requested, we need to drop all non-bottom SCCs.
                 if (onlyBottomSccs) {
-                    for (uint_fast64_t state = 0; state < numberOfStates; ++state) {
-                        // If the block of the state is already known to be dropped, we don't need to check the transitions.
-                        if (!blocksToDrop.get(stateToSccMapping[state])) {
-                            for (typename storm::storage::SparseMatrix<ValueType>::const_iterator successorIt = transitionMatrix.getRowGroup(state).begin(), successorIte = transitionMatrix.getRowGroup(state).end(); successorIt != successorIte; ++successorIt) {
-                                if (subsystem.get(successorIt->getColumn()) && stateToSccMapping[state] != stateToSccMapping[successorIt->getColumn()]) {
-                                    blocksToDrop.set(stateToSccMapping[state]);
-                                    break;
+                    if (subsystem) {
+                        for (uint64_t state : *subsystem) {
+                            // If the block of the state is already known to be dropped, we don't need to check the transitions.
+                            if (!blocksToDrop.get(stateToSccMapping[state])) {
+                                for (uint64_t row = transitionMatrix.getRowGroupIndices()[state], endRow = transitionMatrix.getRowGroupIndices()[state + 1]; row != endRow; ++row) {
+                                    if (choices && !choices->get(row)) {
+                                        continue;
+                                    }
+                                    for (auto const& entry : transitionMatrix.getRow(row)) {
+                                        if (subsystem->get(entry.getColumn()) && stateToSccMapping[state] != stateToSccMapping[entry.getColumn()]) {
+                                            blocksToDrop.set(stateToSccMapping[state]);
+                                            break;
+                                        }
+                                    }
+                                }
+                            }
+                        }
+                    } else {
+                        for (uint64_t state = 0; state < transitionMatrix.getRowGroupCount(); ++state) {
+                            // If the block of the state is already known to be dropped, we don't need to check the transitions.
+                            if (!blocksToDrop.get(stateToSccMapping[state])) {
+                                for (uint64_t row = transitionMatrix.getRowGroupIndices()[state], endRow = transitionMatrix.getRowGroupIndices()[state + 1]; row != endRow; ++row) {
+                                    for (auto const& entry : transitionMatrix.getRow(row)) {
+                                        if (stateToSccMapping[state] != stateToSccMapping[entry.getColumn()]) {
+                                            blocksToDrop.set(stateToSccMapping[state]);
+                                            break;
+                                        }
+                                    }
                                 }
                             }
                         }
@@ -163,15 +206,11 @@ namespace storm {
         template <typename ValueType>
         template <typename RewardModelType>
         void StronglyConnectedComponentDecomposition<ValueType>::performSccDecomposition(storm::models::sparse::Model<ValueType, RewardModelType> const& model, bool dropNaiveSccs, bool onlyBottomSccs) {
-            // Prepare a block that contains all states for a call to the other overload of this function.
-            storm::storage::BitVector fullSystem(model.getNumberOfStates(), true);
-            
-            // Call the overloaded function.
-            performSccDecomposition(model.getTransitionMatrix(), fullSystem, dropNaiveSccs, onlyBottomSccs);
+            performSccDecomposition(model.getTransitionMatrix(), nullptr, nullptr, dropNaiveSccs, onlyBottomSccs);
         }
         
         template <typename ValueType>
-        void StronglyConnectedComponentDecomposition<ValueType>::performSccDecompositionGCM(storm::storage::SparseMatrix<ValueType> const& transitionMatrix, uint_fast64_t startState, storm::storage::BitVector& statesWithSelfLoop, storm::storage::BitVector const& subsystem, uint_fast64_t& currentIndex, storm::storage::BitVector& hasPreorderNumber, std::vector<uint_fast64_t>& preorderNumbers, std::vector<uint_fast64_t>& s, std::vector<uint_fast64_t>& p, storm::storage::BitVector& stateHasScc, std::vector<uint_fast64_t>& stateToSccMapping, uint_fast64_t& sccCount) {
+        void StronglyConnectedComponentDecomposition<ValueType>::performSccDecompositionGCM(storm::storage::SparseMatrix<ValueType> const& transitionMatrix, uint_fast64_t startState, storm::storage::BitVector& statesWithSelfLoop, storm::storage::BitVector const* subsystem, storm::storage::BitVector const* choices, uint_fast64_t& currentIndex, storm::storage::BitVector& hasPreorderNumber, std::vector<uint_fast64_t>& preorderNumbers, std::vector<uint_fast64_t>& s, std::vector<uint_fast64_t>& p, storm::storage::BitVector& stateHasScc, std::vector<uint_fast64_t>& stateToSccMapping, uint_fast64_t& sccCount) {
             
             // Prepare the stack used for turning the recursive procedure into an iterative one.
             std::vector<uint_fast64_t> recursionStateStack;
@@ -190,20 +229,26 @@ namespace storm {
                     s.push_back(currentState);
                     p.push_back(currentState);
                 
-                    for (auto const& successor : transitionMatrix.getRowGroup(currentState)) {
-                        if (subsystem.get(successor.getColumn()) && successor.getValue() != storm::utility::zero<ValueType>()) {
-                            if (currentState == successor.getColumn()) {
-                                statesWithSelfLoop.set(currentState);
-                            }
-                            
-                            if (!hasPreorderNumber.get(successor.getColumn())) {
-                                // In this case, we must recursively visit the successor. We therefore push the state
-                                // onto the recursion stack.
-                                recursionStateStack.push_back(successor.getColumn());
-                            } else {
-                                if (!stateHasScc.get(successor.getColumn())) {
-                                    while (preorderNumbers[p.back()] > preorderNumbers[successor.getColumn()]) {
-                                        p.pop_back();
+                    for (uint64_t row = transitionMatrix.getRowGroupIndices()[currentState], rowEnd = transitionMatrix.getRowGroupIndices()[currentState + 1]; row != rowEnd; ++row) {
+                        if (choices && !choices->get(row)) {
+                            continue;
+                        }
+                        
+                        for (auto const& successor : transitionMatrix.getRow(row)) {
+                            if ((!subsystem || subsystem->get(successor.getColumn())) && successor.getValue() != storm::utility::zero<ValueType>()) {
+                                if (currentState == successor.getColumn()) {
+                                    statesWithSelfLoop.set(currentState);
+                                }
+                                
+                                if (!hasPreorderNumber.get(successor.getColumn())) {
+                                    // In this case, we must recursively visit the successor. We therefore push the state
+                                    // onto the recursion stack.
+                                    recursionStateStack.push_back(successor.getColumn());
+                                } else {
+                                    if (!stateHasScc.get(successor.getColumn())) {
+                                        while (preorderNumbers[p.back()] > preorderNumbers[successor.getColumn()]) {
+                                            p.pop_back();
+                                        }
                                     }
                                 }
                             }
diff --git a/src/storm/storage/StronglyConnectedComponentDecomposition.h b/src/storm/storage/StronglyConnectedComponentDecomposition.h
index 04cebff92..5a9fa2cf4 100644
--- a/src/storm/storage/StronglyConnectedComponentDecomposition.h
+++ b/src/storm/storage/StronglyConnectedComponentDecomposition.h
@@ -103,6 +103,20 @@ namespace storm {
              */
             StronglyConnectedComponentDecomposition(storm::storage::SparseMatrix<ValueType> const& transitionMatrix, storm::storage::BitVector const& subsystem, bool dropNaiveSccs = false, bool onlyBottomSccs = false);
             
+            /*
+             * Creates an SCC decomposition of the given subsystem in the given system (whose transition relation is
+             * given by a sparse matrix).
+             *
+             * @param transitionMatrix The transition matrix of the system to decompose.
+             * @param subsystem A bit vector indicating which subsystem to consider for the decomposition into SCCs.
+             * @param choices A bit vector indicating which choices of the states are contained in the subsystem.
+             * @param dropNaiveSccs A flag that indicates whether trivial SCCs (i.e. SCCs consisting of just one state
+             * without a self-loop) are to be kept in the decomposition.
+             * @param onlyBottomSccs If set to true, only bottom SCCs, i.e. SCCs in which all states have no way of
+             * leaving the SCC), are kept.
+             */
+            StronglyConnectedComponentDecomposition(storm::storage::SparseMatrix<ValueType> const& transitionMatrix, storm::storage::BitVector const& subsystem, storm::storage::BitVector const& choices, bool dropNaiveSccs = false, bool onlyBottomSccs = false);
+            
             /*!
              * Creates an SCC decomposition by copying the given SCC decomposition.
              *
@@ -158,13 +172,14 @@ namespace storm {
              * the vector of blocks of the decomposition.
              *
              * @param transitionMatrix The transition matrix of the system to decompose.
-             * @param subsystem A bit vector indicating which subsystem to consider for the decomposition into SCCs.
+             * @param subsystem An optional bit vector indicating which subsystem to consider.
+             * @param choices An optional bit vector indicating which choices belong to the subsystem.
              * @param dropNaiveSccs A flag that indicates whether trivial SCCs (i.e. SCCs consisting of just one state
              * without a self-loop) are to be kept in the decomposition.
              * @param onlyBottomSccs If set to true, only bottom SCCs, i.e. SCCs in which all states have no way of
              * leaving the SCC), are kept.
              */
-            void performSccDecomposition(storm::storage::SparseMatrix<ValueType> const& transitionMatrix, storm::storage::BitVector const& subsystem, bool dropNaiveSccs, bool onlyBottomSccs);
+            void performSccDecomposition(storm::storage::SparseMatrix<ValueType> const& transitionMatrix, storm::storage::BitVector const* subsystem, storm::storage::BitVector const* choices, bool dropNaiveSccs, bool onlyBottomSccs);
             
             /*!
              * Uses the algorithm by Gabow/Cheriyan/Mehlhorn ("Path-based strongly connected component algorithm") to
@@ -175,7 +190,8 @@ namespace storm {
              * @param startState The starting state for the search of Tarjan's algorithm.
              * @param statesWithSelfLoop A bit vector that is to be filled with all states that have a self-loop. This
              * is later needed for identification of the naive SCCs.
-             * @param subsystem The subsystem to search.
+             * @param subsystem An optional bit vector indicating which subsystem to consider.
+             * @param choices An optional bit vector indicating which choices belong to the subsystem.
              * @param currentIndex The next free index that can be assigned to states.
              * @param hasPreorderNumber A bit that is used to keep track of the states that already have a preorder number.
              * @param preorderNumbers A vector storing the preorder number for each state.
@@ -187,7 +203,7 @@ namespace storm {
              * @param sccCount The number of SCCs that have been computed. As a side effect of this function, this count
              * is increased.
              */
-            void performSccDecompositionGCM(storm::storage::SparseMatrix<ValueType> const& transitionMatrix, uint_fast64_t startState, storm::storage::BitVector& statesWithSelfLoop, storm::storage::BitVector const& subsystem, uint_fast64_t& currentIndex, storm::storage::BitVector& hasPreorderNumber, std::vector<uint_fast64_t>& preorderNumbers, std::vector<uint_fast64_t>& s, std::vector<uint_fast64_t>& p, storm::storage::BitVector& stateHasScc, std::vector<uint_fast64_t>& stateToSccMapping, uint_fast64_t& sccCount);
+            void performSccDecompositionGCM(storm::storage::SparseMatrix<ValueType> const& transitionMatrix, uint_fast64_t startState, storm::storage::BitVector& statesWithSelfLoop, storm::storage::BitVector const* subsystem, storm::storage::BitVector const* choices, uint_fast64_t& currentIndex, storm::storage::BitVector& hasPreorderNumber, std::vector<uint_fast64_t>& preorderNumbers, std::vector<uint_fast64_t>& s, std::vector<uint_fast64_t>& p, storm::storage::BitVector& stateHasScc, std::vector<uint_fast64_t>& stateToSccMapping, uint_fast64_t& sccCount);
         };
     }
 }