Added construction of the state remapping for elimination of non-Markovian states in MAs

src/storm/models/sparse/MarkovAutomaton.cpp

@@ -1,3 +1,5 @@
+#include <queue>
+
 #include "storm/models/sparse/MarkovAutomaton.h"
 
 #include "storm/adapters/RationalFunctionAdapter.h"
@@ -16,7 +18,7 @@
 namespace storm {
     namespace models {
         namespace sparse {
-            
+
             template <typename ValueType, typename RewardModelType>
             MarkovAutomaton<ValueType, RewardModelType>::MarkovAutomaton(storm::storage::SparseMatrix<ValueType> const& transitionMatrix,
                                                                          storm::models::sparse::StateLabeling const& stateLabeling,
@@ -25,7 +27,7 @@ namespace storm {
                     : MarkovAutomaton<ValueType, RewardModelType>(storm::storage::sparse::ModelComponents<ValueType, RewardModelType>(transitionMatrix, stateLabeling, rewardModels, true, markovianStates)) {
                 // Intentionally left empty
             }
-            
+
             template <typename ValueType, typename RewardModelType>
             MarkovAutomaton<ValueType, RewardModelType>::MarkovAutomaton(storm::storage::SparseMatrix<ValueType>&& transitionMatrix,
                                                                          storm::models::sparse::StateLabeling&& stateLabeling,
@@ -34,23 +36,23 @@ namespace storm {
                     : MarkovAutomaton<ValueType, RewardModelType>(storm::storage::sparse::ModelComponents<ValueType, RewardModelType>(std::move(transitionMatrix), std::move(stateLabeling), std::move(rewardModels), true, std::move(markovianStates))) {
                 // Intentionally left empty
             }
-            
+
             template <typename ValueType, typename RewardModelType>
             MarkovAutomaton<ValueType, RewardModelType>::MarkovAutomaton(storm::storage::sparse::ModelComponents<ValueType, RewardModelType> const& components) : NondeterministicModel<ValueType, RewardModelType>(ModelType::MarkovAutomaton, components), markovianStates(components.markovianStates.get()) {
-                
+
                 if (components.exitRates) {
                     exitRates = components.exitRates.get();
                 }
-                
+
                 if (components.rateTransitions) {
                     this->turnRatesToProbabilities();
                 }
                 closed = this->checkIsClosed();
             }
-            
+
             template <typename ValueType, typename RewardModelType>
             MarkovAutomaton<ValueType, RewardModelType>::MarkovAutomaton(storm::storage::sparse::ModelComponents<ValueType, RewardModelType>&& components) : NondeterministicModel<ValueType, RewardModelType>(ModelType::MarkovAutomaton, std::move(components)), markovianStates(std::move(components.markovianStates.get())) {
-                
+
                 if (components.exitRates) {
                     exitRates = std::move(components.exitRates.get());
                 }
@@ -60,52 +62,52 @@ namespace storm {
                 }
                 closed = this->checkIsClosed();
             }
-            
+
             template <typename ValueType, typename RewardModelType>
             bool MarkovAutomaton<ValueType, RewardModelType>::isClosed() const {
                 return closed;
             }
-            
+
             template <typename ValueType, typename RewardModelType>
             bool MarkovAutomaton<ValueType, RewardModelType>::isHybridState(storm::storage::sparse::state_type state) const {
                 return isMarkovianState(state) && (this->getTransitionMatrix().getRowGroupSize(state) > 1);
             }
-            
+
             template <typename ValueType, typename RewardModelType>
             bool MarkovAutomaton<ValueType, RewardModelType>::isMarkovianState(storm::storage::sparse::state_type state) const {
                 return this->markovianStates.get(state);
             }
-            
+
             template <typename ValueType, typename RewardModelType>
             bool MarkovAutomaton<ValueType, RewardModelType>::isProbabilisticState(storm::storage::sparse::state_type state) const {
                 return !this->markovianStates.get(state);
             }
-            
+
             template <typename ValueType, typename RewardModelType>
             std::vector<ValueType> const& MarkovAutomaton<ValueType, RewardModelType>::getExitRates() const {
                 return this->exitRates;
             }
-            
+
             template <typename ValueType, typename RewardModelType>
             std::vector<ValueType>& MarkovAutomaton<ValueType, RewardModelType>::getExitRates() {
                 return this->exitRates;
             }
-            
+
             template <typename ValueType, typename RewardModelType>
             ValueType const& MarkovAutomaton<ValueType, RewardModelType>::getExitRate(storm::storage::sparse::state_type state) const {
                 return this->exitRates[state];
             }
-            
+
             template <typename ValueType, typename RewardModelType>
             ValueType MarkovAutomaton<ValueType, RewardModelType>::getMaximalExitRate() const {
                 return storm::utility::vector::max_if(this->exitRates, this->markovianStates);
             }
-            
+
             template <typename ValueType, typename RewardModelType>
             storm::storage::BitVector const& MarkovAutomaton<ValueType, RewardModelType>::getMarkovianStates() const {
                 return this->markovianStates;
             }
-            
+
             template <typename ValueType, typename RewardModelType>
             void MarkovAutomaton<ValueType, RewardModelType>::close() {
                 if (!closed) {
@@ -120,16 +122,16 @@ namespace storm {
                             exitRates[state] = storm::utility::zero<ValueType>();
                         }
                     }
-                    
+
                     if (!keptChoices.full()) {
                         *this = std::move(*storm::transformer::buildSubsystem(*this, storm::storage::BitVector(this->getNumberOfStates(), true), keptChoices, false).model->template as<MarkovAutomaton<ValueType, RewardModelType>>());
                     }
-                    
+
                     // Mark the automaton as closed.
                     closed = true;
                 }
             }
-            
+
             template <typename ValueType, typename RewardModelType>
             void MarkovAutomaton<ValueType, RewardModelType>::turnRatesToProbabilities() {
                 bool assertRates = (this->exitRates.size() == this->getNumberOfStates());
@@ -137,7 +139,7 @@ namespace storm {
                     STORM_LOG_THROW(this->exitRates.empty(), storm::exceptions::InvalidArgumentException, "The specified exit rate vector has an unexpected size.");
                     this->exitRates.reserve(this->getNumberOfStates());
                 }
-                
+
                 storm::utility::ConstantsComparator<ValueType> comparator;
                 for (uint_fast64_t state = 0; state< this->getNumberOfStates(); ++state) {
                     uint_fast64_t row = this->getTransitionMatrix().getRowGroupIndices()[state];
@@ -163,12 +165,12 @@ namespace storm {
                     }
                 }
             }
-            
+
             template <typename ValueType, typename RewardModelType>
             bool MarkovAutomaton<ValueType, RewardModelType>::isConvertibleToCtmc() const {
                 return isClosed() && markovianStates.full();
             }
-            
+
             template <typename ValueType, typename RewardModelType>
             bool MarkovAutomaton<ValueType, RewardModelType>::hasOnlyTrivialNondeterminism() const {
                 // Check every state
@@ -185,7 +187,7 @@ namespace storm {
                 }
                 return true;
             }
-            
+
             template <typename ValueType, typename RewardModelType>
             bool MarkovAutomaton<ValueType, RewardModelType>::checkIsClosed() const {
                 for (auto state : markovianStates) {
@@ -195,7 +197,7 @@ namespace storm {
                 }
                 return true;
             }
-            
+
             template <typename ValueType, typename RewardModelType>
             std::shared_ptr<storm::models::sparse::Ctmc<ValueType, RewardModelType>> MarkovAutomaton<ValueType, RewardModelType>::convertToCtmc() const {
                 if (isClosed() && markovianStates.full()) {
@@ -265,7 +267,109 @@ namespace storm {
                 return std::make_shared<storm::models::sparse::Ctmc<ValueType, RewardModelType>>(std::move(rateMatrix), std::move(stateLabeling));
             }
 
-            
+            template<typename ValueType, typename RewardModelType>
+            std::shared_ptr<MarkovAutomaton<ValueType, RewardModelType>>
+            MarkovAutomaton<ValueType, RewardModelType>::eliminateNonmarkovianStates() const {
+                if (isClosed() && markovianStates.full()) {
+                    storm::storage::sparse::ModelComponents<ValueType, RewardModelType> components(
+                            this->getTransitionMatrix(), this->getStateLabeling(), this->getRewardModels(), false);
+                    components.exitRates = this->getExitRates();
+                    if (this->hasChoiceLabeling()) {
+                        components.choiceLabeling = this->getChoiceLabeling();
+                    }
+                    if (this->hasStateValuations()) {
+                        components.stateValuations = this->getStateValuations();
+                    }
+                    if (this->hasChoiceOrigins()) {
+                        components.choiceOrigins = this->getChoiceOrigins();
+                    }
+                    return std::make_shared<MarkovAutomaton<ValueType, RewardModelType>>(std::move(components));
+                }
+
+                std::map<uint_fast64_t, uint_fast64_t> stateRemapping;
+                std::set<uint_fast64_t> statesToKeep;
+                std::queue<uint_fast64_t> changedStates;
+                std::queue<uint_fast64_t> queue;
+
+                storm::storage::SparseMatrix<ValueType> backwards = this->getBackwardTransitions();
+
+                // Determine the state remapping
+                // TODO Consider state labels
+                for (uint_fast64_t base_state = 0; base_state < this->getNumberOfStates(); ++base_state) {
+                    STORM_LOG_ASSERT(!this->isHybridState(base_state), "Base state is hybrid.");
+                    if (this->isMarkovianState(base_state)) {
+                        queue.push(base_state);
+
+                        while (!queue.empty()) {
+                            auto currState = queue.front();
+                            queue.pop();
+                            // Get predecessors from matrix
+                            typename storm::storage::SparseMatrix<ValueType>::rows entriesInRow = backwards.getRow(
+                                    currState);
+                            for (auto entryIt = entriesInRow.begin(), entryIte = entriesInRow.end();
+                                 entryIt != entryIte; ++entryIt) {
+                                uint_fast64_t predecessor = entryIt->getColumn();
+                                if (!this->isMarkovianState(predecessor) && !statesToKeep.count(predecessor)) {
+                                    if (!stateRemapping.count(predecessor)) {
+                                        stateRemapping[predecessor] = base_state;
+                                        queue.push(predecessor);
+                                    } else if (stateRemapping[predecessor] != base_state) {
+                                        stateRemapping.erase(predecessor);
+                                        statesToKeep.insert(predecessor);
+                                        changedStates.push(predecessor);
+                                    }
+                                }
+                            }
+                        }
+                    }
+                }
+
+                // Correct the mapping with the states which have to be kept
+                while (!changedStates.empty()) {
+                    uint_fast64_t base_state = changedStates.front();
+                    queue.push(base_state);
+
+                    while (!queue.empty()) {
+                        auto currState = queue.front();
+                        queue.pop();
+                        // Get predecessors from matrix
+                        typename storm::storage::SparseMatrix<ValueType>::rows entriesInRow = backwards.getRow(
+                                currState);
+                        for (auto entryIt = entriesInRow.begin(), entryIte = entriesInRow.end();
+                             entryIt != entryIte; ++entryIt) {
+                            uint_fast64_t predecessor = entryIt->getColumn();
+                            if (!this->isMarkovianState(predecessor) && !statesToKeep.count(predecessor)) {
+                                if (!stateRemapping.count(predecessor)) {
+                                    stateRemapping[predecessor] = base_state;
+                                    queue.push(predecessor);
+                                } else if (stateRemapping[predecessor] != base_state) {
+                                    stateRemapping.erase(predecessor);
+                                    statesToKeep.insert(predecessor);
+                                    changedStates.push(predecessor);
+                                }
+                            }
+                        }
+                    }
+
+                    changedStates.pop();
+                }
+
+                // At this point, we hopefully have a valid mapping which eliminates a lot of states
+
+                STORM_PRINT("Remapping \n")
+                for (auto entry : stateRemapping) {
+                    STORM_PRINT(std::to_string(entry.first) << " -> " << std::to_string(entry.second) << "\n")
+                }
+                STORM_PRINT("Remapped States: " << stateRemapping.size() << "\n")
+                // TODO test some examples, especially ones containing non-determinism
+
+                // Build the new matrix
+                // TODO
+
+                return nullptr;
+            }
+
+
             template<typename ValueType, typename RewardModelType>
             void MarkovAutomaton<ValueType, RewardModelType>::printModelInformationToStream(std::ostream& out) const {
                 this->printModelInformationHeaderToStream(out);
@@ -274,13 +378,15 @@ namespace storm {
                 out << "Max. Rate.: \t" << this->getMaximalExitRate() << std::endl;
                 this->printModelInformationFooterToStream(out);
             }
-            
-            
+
+
             template class MarkovAutomaton<double>;
 #ifdef STORM_HAVE_CARL
+
             template class MarkovAutomaton<storm::RationalNumber>;
-            
+
             template class MarkovAutomaton<double, storm::models::sparse::StandardRewardModel<storm::Interval>>;
+
             template class MarkovAutomaton<storm::RationalFunction>;
 #endif
         } // namespace sparse

src/storm/models/sparse/MarkovAutomaton.h

@@ -147,6 +147,8 @@ namespace storm {
                  * @return The resulting CTMC.
                  */
                 std::shared_ptr<storm::models::sparse::Ctmc<ValueType, RewardModelType>> convertToCtmc() const;
+
+                std::shared_ptr<MarkovAutomaton<ValueType, RewardModelType>> eliminateNonmarkovianStates() const;
                 
                 virtual void printModelInformationToStream(std::ostream& out) const override;