fixes github issue #24 related to MEC decomposition

resources/examples/testfiles/mdp/prism-mec-example1.nm

@@ -0,0 +1,12 @@
+mdp
+
+module test
+
+	x : [0..2];
+	
+	[] x=0 -> true;
+	[] x=0 -> 0.5 : (x'=1) + 0.5: (x'=2);
+	[] x=1 -> (x'=0);
+	[] x=2 -> true;
+	
+endmodule 

resources/examples/testfiles/mdp/prism-mec-example2.nm

@@ -0,0 +1,13 @@
+mdp
+
+module test
+
+	x : [0..2];
+	
+	[] x=0 -> true;
+	[] x=0 -> 0.5 : (x'=1) + 0.5: (x'=1);
+	[] x=0 -> (x'=2);
+	[] x=1 -> (x'=0);
+	[] x=2 -> true;
+	
+endmodule 

src/storm/storage/MaximalEndComponentDecomposition.cpp

@@ -80,6 +80,19 @@ namespace storm {
                 endComponentStateSets.emplace_back(states.begin(), states.end(), true);
             }
             storm::storage::BitVector statesToCheck(numberOfStates);
+            storm::storage::BitVector includedChoices;
+            if (choices) {
+                includedChoices = *choices;
+            } else if (states) {
+                includedChoices = storm::storage::BitVector(transitionMatrix.getRowCount());
+                for (auto state : *states) {
+                    for (uint_fast64_t choice = nondeterministicChoiceIndices[state]; choice < nondeterministicChoiceIndices[state + 1]; ++choice) {
+                        includedChoices.set(choice, true);
+                    }
+                }
+            } else {
+                includedChoices = storm::storage::BitVector(transitionMatrix.getRowCount(), true);
+            }
                         
             for (std::list<StateBlock>::const_iterator mecIterator = endComponentStateSets.begin(); mecIterator != endComponentStateSets.end();) {
                 StateBlock const& mec = *mecIterator;
@@ -88,7 +101,7 @@ namespace storm {
                 bool mecChanged = false;
                 
                 // Get an SCC decomposition of the current MEC candidate.
-                StronglyConnectedComponentDecomposition<ValueType> sccs(transitionMatrix, mec, true);
+                StronglyConnectedComponentDecomposition<ValueType> sccs(transitionMatrix, mec, includedChoices, true);
                 
                 // We need to do another iteration in case we have either more than once SCC or the SCC is smaller than
                 // the MEC canditate itself.
@@ -105,11 +118,17 @@ namespace storm {
                             bool keepStateInMEC = false;
                             
                             for (uint_fast64_t choice = nondeterministicChoiceIndices[state]; choice < nondeterministicChoiceIndices[state + 1]; ++choice) {
+                                
                                 // If the choice is not part of our subsystem, skip it.
                                 if (choices && !choices->get(choice)) {
                                     continue;
                                 }
 
+                                // If the choice is not included any more, skip it.
+                                if (!includedChoices.get(choice)) {
+                                    continue;
+                                }
+                                
                                 bool choiceContainedInMEC = true;
                                 for (auto const& entry : transitionMatrix.getRow(choice)) {
                                     if (storm::utility::isZero(entry.getValue())) {
@@ -117,6 +136,7 @@ namespace storm {
                                     }
                                         
                                     if (!scc.containsState(entry.getColumn())) {
+                                        includedChoices.set(choice, false);
                                         choiceContainedInMEC = false;
                                         break;
                                     }
@@ -125,7 +145,6 @@ namespace storm {
                                 // If there is at least one choice whose successor states are fully contained in the MEC, we can leave the state in the MEC.
                                 if (choiceContainedInMEC) {
                                     keepStateInMEC = true;
-                                    break;
                                 }
                             }
                             
@@ -185,15 +204,7 @@ namespace storm {
                             continue;
                         }
                         
-                        bool choiceContained = true;
-                        for (auto const& entry : transitionMatrix.getRow(choice)) {
-                            if (!mecStateSet.containsState(entry.getColumn())) {
-                                choiceContained = false;
-                                break;
-                            }
-                        }
-                        
-                        if (choiceContained) {
+                        if (includedChoices.get(choice)) {
                             containedChoices.insert(choice);
                         }
                     }

src/storm/storage/StronglyConnectedComponentDecomposition.cpp

@@ -38,6 +38,11 @@ namespace storm {
             performSccDecomposition(transitionMatrix, &subsystem, nullptr, dropNaiveSccs, onlyBottomSccs);
         }
         
+        template <typename ValueType>
+        StronglyConnectedComponentDecomposition<ValueType>::StronglyConnectedComponentDecomposition(storm::storage::SparseMatrix<ValueType> const& transitionMatrix, StateBlock const& block, storm::storage::BitVector const& choices, bool dropNaiveSccs, bool onlyBottomSccs) {
+            storm::storage::BitVector subsystem(transitionMatrix.getRowGroupCount(), block.begin(), block.end());
+            performSccDecomposition(transitionMatrix, &subsystem, &choices, dropNaiveSccs, onlyBottomSccs);
+        }
         
         template <typename ValueType>
         StronglyConnectedComponentDecomposition<ValueType>::StronglyConnectedComponentDecomposition(storm::storage::SparseMatrix<ValueType> const& transitionMatrix, bool dropNaiveSccs, bool onlyBottomSccs) {

src/storm/storage/StronglyConnectedComponentDecomposition.h

@@ -79,6 +79,20 @@ namespace storm {
              */
             StronglyConnectedComponentDecomposition(storm::storage::SparseMatrix<ValueType> const& transitionMatrix, StateBlock const& block, 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 block The block to decompose 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, StateBlock const& block, storm::storage::BitVector const& choices, bool dropNaiveSccs = false, bool onlyBottomSccs = false);
+
             /*
              * Creates an SCC decomposition of the given system (whose transition relation is given by a sparse matrix).
              *

src/test/storm/storage/MaximalEndComponentDecompositionTest.cpp

@@ -3,7 +3,11 @@
 #include "storm/parser/AutoParser.h"
 #include "storm/storage/MaximalEndComponentDecomposition.h"
 #include "storm/models/sparse/MarkovAutomaton.h"
+#include "storm/models/sparse/Mdp.h"
 #include "storm/models/sparse/StandardRewardModel.h"
+#include "storm/builder/ExplicitModelBuilder.h"
+#include "storm/storage/SymbolicModelDescription.h"
+#include "storm/parser/PrismParser.h"
 
 TEST(MaximalEndComponentDecomposition, FullSystem1) {
     std::shared_ptr<storm::models::sparse::Model<double>> abstractModel = storm::parser::AutoParser<>::parseModel(STORM_TEST_RESOURCES_DIR "/tra/tiny1.tra", STORM_TEST_RESOURCES_DIR "/lab/tiny1.lab", "", "");
@@ -133,3 +137,42 @@ TEST(MaximalEndComponentDecomposition, Subsystem) {
         ASSERT_TRUE(false);
     }
 }
+
+TEST(MaximalEndComponentDecomposition, Example1) {
+    std::string prismModelPath = STORM_TEST_RESOURCES_DIR "/mdp/prism-mec-example1.nm";
+    storm::storage::SymbolicModelDescription modelDescription = storm::parser::PrismParser::parse(prismModelPath);
+    storm::prism::Program program = modelDescription.preprocess().asPrismProgram();
+    
+    std::shared_ptr<storm::models::sparse::Model<double>> model = storm::builder::ExplicitModelBuilder<double>(program).build();
+    std::shared_ptr<storm::models::sparse::Mdp<double>> mdp = model->as<storm::models::sparse::Mdp<double>>();
+    
+    storm::storage::MaximalEndComponentDecomposition<double> mecDecomposition(*mdp);
+    
+    EXPECT_EQ(mecDecomposition.size(), 2);
+    
+    ASSERT_TRUE(mecDecomposition[0].getStateSet() == storm::storage::MaximalEndComponent::set_type{2});
+    EXPECT_TRUE(mecDecomposition[0].getChoicesForState(2) == storm::storage::MaximalEndComponent::set_type{3});
+    
+    ASSERT_TRUE(mecDecomposition[1].getStateSet() == storm::storage::MaximalEndComponent::set_type{0});
+    EXPECT_TRUE(mecDecomposition[1].getChoicesForState(0) == storm::storage::MaximalEndComponent::set_type{0});
+}
+
+TEST(MaximalEndComponentDecomposition, Example2) {
+    std::string prismModelPath = STORM_TEST_RESOURCES_DIR "/mdp/prism-mec-example2.nm";
+    storm::storage::SymbolicModelDescription modelDescription = storm::parser::PrismParser::parse(prismModelPath);
+    storm::prism::Program program = modelDescription.preprocess().asPrismProgram();
+    
+    std::shared_ptr<storm::models::sparse::Model<double>> model = storm::builder::ExplicitModelBuilder<double>(program).build();
+    std::shared_ptr<storm::models::sparse::Mdp<double>> mdp = model->as<storm::models::sparse::Mdp<double>>();
+    
+    storm::storage::MaximalEndComponentDecomposition<double> mecDecomposition(*mdp);
+    
+    EXPECT_EQ(mecDecomposition.size(), 2);
+    
+    ASSERT_TRUE(mecDecomposition[0].getStateSet() == storm::storage::MaximalEndComponent::set_type{2});
+    EXPECT_TRUE(mecDecomposition[0].getChoicesForState(2) == storm::storage::MaximalEndComponent::set_type{4});
+    
+    ASSERT_TRUE((mecDecomposition[1].getStateSet() == storm::storage::MaximalEndComponent::set_type{0, 1}));
+    EXPECT_TRUE((mecDecomposition[1].getChoicesForState(0) == storm::storage::MaximalEndComponent::set_type{0, 1}));
+    EXPECT_TRUE((mecDecomposition[1].getChoicesForState(1) == storm::storage::MaximalEndComponent::set_type{3}));
+}