Added bisimulation settings. Further work on weak bisimulation.

Former-commit-id: c04759575a

src/settings/SettingsManager.cpp

@@ -24,6 +24,7 @@ namespace storm {
             this->addModule(std::unique_ptr<modules::ModuleSettings>(new modules::CuddSettings(*this)));
             this->addModule(std::unique_ptr<modules::ModuleSettings>(new modules::GmmxxEquationSolverSettings(*this)));
             this->addModule(std::unique_ptr<modules::ModuleSettings>(new modules::NativeEquationSolverSettings(*this)));
+            this->addModule(std::unique_ptr<modules::ModuleSettings>(new modules::BisimulationSettings(*this)));
             this->addModule(std::unique_ptr<modules::ModuleSettings>(new modules::GlpkSettings(*this)));
             this->addModule(std::unique_ptr<modules::ModuleSettings>(new modules::GurobiSettings(*this)));
         }
@@ -495,6 +496,10 @@ namespace storm {
             return dynamic_cast<storm::settings::modules::NativeEquationSolverSettings const&>(manager().getModule(storm::settings::modules::NativeEquationSolverSettings::moduleName));
         }
         
+        storm::settings::modules::BisimulationSettings const& bisimulationSettings() {
+            return dynamic_cast<storm::settings::modules::BisimulationSettings const&>(manager().getModule(storm::settings::modules::BisimulationSettings::moduleName));
+        }
+        
         storm::settings::modules::GlpkSettings const& glpkSettings() {
             return dynamic_cast<storm::settings::modules::GlpkSettings const&>(manager().getModule(storm::settings::modules::GlpkSettings::moduleName));
         }

src/settings/SettingsManager.h

@@ -25,6 +25,7 @@
 #include "src/settings/modules/CuddSettings.h"
 #include "src/settings/modules/GmmxxEquationSolverSettings.h"
 #include "src/settings/modules/NativeEquationSolverSettings.h"
+#include "src/settings/modules/BisimulationSettings.h"
 #include "src/settings/modules/GlpkSettings.h"
 #include "src/settings/modules/GurobiSettings.h"
 
@@ -286,6 +287,13 @@ namespace storm {
          * @return An object that allows accessing the settings of the native equation solver.
          */
         storm::settings::modules::NativeEquationSolverSettings const& nativeEquationSolverSettings();
+
+        /*!
+         * Retrieves the settings of the native equation solver.
+         *
+         * @return An object that allows accessing the settings of the native equation solver.
+         */
+        storm::settings::modules::BisimulationSettings const& bisimulationSettings();
         
         /*!
          * Retrieves the settings of glpk.

src/settings/modules/BisimulationSettings.cpp

@@ -0,0 +1,33 @@
+#include "src/settings/modules/BisimulationSettings.h"
+
+#include "src/settings/SettingsManager.h"
+
+namespace storm {
+    namespace settings {
+        namespace modules {
+            
+            const std::string BisimulationSettings::moduleName = "bisimulation";
+            const std::string BisimulationSettings::typeOptionName = "type";
+            
+            BisimulationSettings::BisimulationSettings(storm::settings::SettingsManager& settingsManager) : ModuleSettings(settingsManager, moduleName) {
+                std::vector<std::string> types = { "strong", "weak" };
+                this->addOption(storm::settings::OptionBuilder(moduleName, typeOptionName, true, "Sets the kind of bisimulation quotienting used. Available are: { strong, weak }.").addArgument(storm::settings::ArgumentBuilder::createStringArgument("name", "The name of the type to use.").addValidationFunctionString(storm::settings::ArgumentValidators::stringInListValidator(types)).setDefaultValueString("strong").build()).build());
+            }
+            
+            bool BisimulationSettings::isStrongBisimulationSet() const {
+                if (this->getOption(typeOptionName).getArgumentByName("name").getValueAsString() == "strong") {
+                    return true;
+                }
+                return false;
+            }
+            
+            bool BisimulationSettings::isWeakBisimulationSet() const {
+                if (this->getOption(typeOptionName).getArgumentByName("name").getValueAsString() == "weak") {
+                    return true;
+                }
+                return false;
+            }
+            
+        } // namespace modules
+    } // namespace settings
+} // namespace storm

src/settings/modules/BisimulationSettings.h

@@ -0,0 +1,50 @@
+#ifndef STORM_SETTINGS_MODULES_BISIMULATIONSETTINGS_H_
+#define STORM_SETTINGS_MODULES_BISIMULATIONSETTINGS_H_
+
+#include "src/settings/modules/ModuleSettings.h"
+
+namespace storm {
+    namespace settings {
+        namespace modules {
+            
+            /*!
+             * This class represents the bisimulation settings.
+             */
+            class BisimulationSettings : public ModuleSettings {
+            public:
+                // An enumeration of all available bisimulation types.
+                enum class BisimulationType { Strong, Weak };
+                
+                /*!
+                 * Creates a new set of bisimulation settings that is managed by the given manager.
+                 *
+                 * @param settingsManager The responsible manager.
+                 */
+                BisimulationSettings(storm::settings::SettingsManager& settingsManager);
+                
+                /*!
+                 * Retrieves whether strong bisimulation is to be used.
+                 *
+                 * @return True iff strong bisimulation is to be used.
+                 */
+                bool isStrongBisimulationSet() const;
+
+                /*!
+                 * Retrieves whether weak bisimulation is to be used.
+                 *
+                 * @return True iff weak bisimulation is to be used.
+                 */
+                bool isWeakBisimulationSet() const;
+
+                // The name of the module.
+                static const std::string moduleName;
+                
+            private:
+                // Define the string names of the options as constants.
+                static const std::string typeOptionName;
+            };
+        } // namespace modules
+    } // namespace settings
+} // namespace storm
+
+#endif /* STORM_SETTINGS_MODULES_BISIMULATIONSETTINGS_H_ */

src/storage/BitVector.cpp

@@ -1,4 +1,5 @@
 #include <boost/container/flat_set.hpp>
+#include <iostream>
 
 #include "src/storage/BitVector.h"
 #include "src/exceptions/InvalidArgumentException.h"
@@ -96,6 +97,28 @@ namespace storm {
             return *this;
         }
         
+        bool BitVector::operator<(BitVector const& other) const {
+            if (this->size() < other.size()) {
+                return true;
+            } else if (this->size() > other.size()) {
+                return false;
+            }
+            
+            std::vector<uint64_t>::const_iterator first1 = this->bucketVector.begin();
+            std::vector<uint64_t>::const_iterator last1 = this->bucketVector.end();
+            std::vector<uint64_t>::const_iterator first2 = other.bucketVector.begin();
+            std::vector<uint64_t>::const_iterator last2 = other.bucketVector.end();
+            
+            for (; first1 != last1; ++first1, ++first2) {
+                if (*first1 < *first2) {
+                    return true;
+                } else if (*first1 > *first2) {
+                    return false;
+                }
+            }
+            return false;
+        }
+        
         BitVector& BitVector::operator=(BitVector&& other) {
             // Only perform the assignment if the source and target are not identical.
             if (this != &other) {
@@ -122,8 +145,8 @@ namespace storm {
         }
         
         void BitVector::set(uint_fast64_t index, bool value) {
-            uint_fast64_t bucket = index >> 6;
             if (index >= bitCount) throw storm::exceptions::OutOfRangeException() << "Invalid call to BitVector::set: written index " << index << " out of bounds.";
+            uint_fast64_t bucket = index >> 6;
             
             uint_fast64_t mask = static_cast<uint_fast64_t>(1) << (index & mod64mask);
             if (value) {

src/storage/BitVector.h

@@ -155,7 +155,15 @@ namespace storm {
              * into it.
              */
             BitVector& operator=(BitVector&& other);
-                        
+            
+            /*!
+             * Retrieves whether the current bit vector is (in some order) smaller than the given one.
+             *
+             * @param other The other bit vector.
+             * @return True iff the current bit vector is smaller than the given one.
+             */
+            bool operator<(BitVector const& other) const;
+            
             /*!
              * Sets the given truth value at the given index.
              *

src/storage/DeterministicModelStrongBisimulationDecomposition.cpp

@@ -149,15 +149,9 @@ namespace storm {
         
         template<typename ValueType>
         bool DeterministicModelStrongBisimulationDecomposition<ValueType>::Block::check() const {
-            if (this->begin >= this->end) {
-                assert(false);
-            }
-            if (this->prev != nullptr) {
-                assert (this->prev->next == this);
-            }
-            if (this->next != nullptr) {
-                assert (this->next->prev == this);
-            }
+            assert(this->begin < this->end);
+            assert(this->prev == nullptr || this->prev->next == this);
+            assert(this->next == nullptr || this->next->prev == this);
             return true;
         }
         
@@ -260,6 +254,11 @@ namespace storm {
                 positions[state] = state;
                 stateToBlockMapping[state] = &blocks.back();
             }
+            
+            // If we are requested to store silent probabilities, we need to prepare the vector.
+            if (this->keepSilentProbabilities) {
+                silentProbabilities = std::vector<ValueType>(numberOfStates);
+            }
         }
         
         template<typename ValueType>
@@ -300,6 +299,11 @@ namespace storm {
                 stateToBlockMapping[state] = &thirdBlock;
                 ++position;
             }
+            
+            // If we are requested to store silent probabilities, we need to prepare the vector.
+            if (this->keepSilentProbabilities) {
+                silentProbabilities = std::vector<ValueType>(numberOfStates);
+            }
         }
         
         template<typename ValueType>
@@ -469,6 +473,46 @@ namespace storm {
             }
         }
         
+        template<typename ValueType>
+        bool DeterministicModelStrongBisimulationDecomposition<ValueType>::Partition::isSilent(storm::storage::sparse::state_type state, storm::utility::ConstantsComparator<ValueType> const& comparator) const {
+            STORM_LOG_ASSERT(this->keepSilentProbabilities, "Unable to retrieve silentness of state, because silent probabilities are not tracked.");
+            return comparator.isOne(this->silentProbabilities[state]);
+        }
+        
+        template<typename ValueType>
+        bool DeterministicModelStrongBisimulationDecomposition<ValueType>::Partition::hasSilentProbability(storm::storage::sparse::state_type state, storm::utility::ConstantsComparator<ValueType> const& comparator) const {
+            STORM_LOG_ASSERT(this->keepSilentProbabilities, "Unable to retrieve silentness of state, because silent probabilities are not tracked.");
+            return !comparator.isZero(this->silentProbabilities[state]);
+        }
+        
+        template<typename ValueType>
+        ValueType const& DeterministicModelStrongBisimulationDecomposition<ValueType>::Partition::getSilentProbability(storm::storage::sparse::state_type state) const {
+            STORM_LOG_ASSERT(this->keepSilentProbabilities, "Unable to retrieve silent probability of state, because silent probabilities are not tracked.");
+            return this->silentProbabilities[state];
+        }
+        
+        template<typename ValueType>
+        void DeterministicModelStrongBisimulationDecomposition<ValueType>::Partition::setSilentProbabilities(typename std::vector<std::pair<storm::storage::sparse::state_type, ValueType>>::iterator first, typename std::vector<std::pair<storm::storage::sparse::state_type, ValueType>>::iterator last) {
+            STORM_LOG_ASSERT(this->keepSilentProbabilities, "Unable to set silent probability of state, because silent probabilities are not tracked.");
+            for (; first != last; ++first) {
+                this->silentProbabilities[first->first] = first->second;
+            }
+        }
+        
+        template<typename ValueType>
+        void DeterministicModelStrongBisimulationDecomposition<ValueType>::Partition::setSilentProbabilitiesToZero(typename std::vector<std::pair<storm::storage::sparse::state_type, ValueType>>::iterator first, typename std::vector<std::pair<storm::storage::sparse::state_type, ValueType>>::iterator last) {
+            STORM_LOG_ASSERT(this->keepSilentProbabilities, "Unable to set silent probability of state, because silent probabilities are not tracked.");
+            for (; first != last; ++first) {
+                this->silentProbabilities[first->first] = storm::utility::zero<ValueType>();
+            }
+        }
+        
+        template<typename ValueType>
+        void DeterministicModelStrongBisimulationDecomposition<ValueType>::Partition::setSilentProbability(storm::storage::sparse::state_type state, ValueType const& value) {
+            STORM_LOG_ASSERT(this->keepSilentProbabilities, "Unable to set silent probability of state, because silent probabilities are not tracked.");
+            this->silentProbabilities[state] = value;
+        }
+        
         template<typename ValueType>
         std::list<typename DeterministicModelStrongBisimulationDecomposition<ValueType>::Block> const& DeterministicModelStrongBisimulationDecomposition<ValueType>::Partition::getBlocks() const {
             return this->blocks;
@@ -530,17 +574,23 @@ namespace storm {
         }
         
         template<typename ValueType>
-        DeterministicModelStrongBisimulationDecomposition<ValueType>::DeterministicModelStrongBisimulationDecomposition(storm::models::Dtmc<ValueType> const& model, bool buildQuotient) {
+        DeterministicModelStrongBisimulationDecomposition<ValueType>::DeterministicModelStrongBisimulationDecomposition(storm::models::Dtmc<ValueType> const& model, bool weak, bool buildQuotient) {
             STORM_LOG_THROW(!model.hasStateRewards() && !model.hasTransitionRewards(), storm::exceptions::IllegalFunctionCallException, "Bisimulation is currently only supported for models without reward structures.");
-            Partition initialPartition = getLabelBasedInitialPartition(model);
-            partitionRefinement(model, model.getBackwardTransitions(), initialPartition, buildQuotient);
+            Partition initialPartition = getLabelBasedInitialPartition(model, weak);
+            if (weak) {
+                this->initializeSilentProbabilities(model, initialPartition);
+            }
+            partitionRefinement(model, model.getBackwardTransitions(), initialPartition, weak, buildQuotient);
         }
 
         template<typename ValueType>
-        DeterministicModelStrongBisimulationDecomposition<ValueType>::DeterministicModelStrongBisimulationDecomposition(storm::models::Ctmc<ValueType> const& model, bool buildQuotient) {
+        DeterministicModelStrongBisimulationDecomposition<ValueType>::DeterministicModelStrongBisimulationDecomposition(storm::models::Ctmc<ValueType> const& model, bool weak, bool buildQuotient) {
             STORM_LOG_THROW(!model.hasStateRewards() && !model.hasTransitionRewards(), storm::exceptions::IllegalFunctionCallException, "Bisimulation is currently only supported for models without reward structures.");
-            Partition initialPartition = getLabelBasedInitialPartition(model);
-            partitionRefinement(model, model.getBackwardTransitions(), initialPartition, buildQuotient);
+            Partition initialPartition = getLabelBasedInitialPartition(model, weak);
+            if (weak) {
+                this->initializeSilentProbabilities(model, initialPartition);
+            }
+            partitionRefinement(model, model.getBackwardTransitions(), initialPartition, weak, buildQuotient);
         }
         
         template<typename ValueType>
@@ -548,7 +598,7 @@ namespace storm {
             STORM_LOG_THROW(!model.hasStateRewards() && !model.hasTransitionRewards(), storm::exceptions::IllegalFunctionCallException, "Bisimulation is currently only supported for models without reward structures.");
             storm::storage::SparseMatrix<ValueType> backwardTransitions = model.getBackwardTransitions();
             Partition initialPartition = getMeasureDrivenInitialPartition(model, backwardTransitions, phiLabel, psiLabel, bounded);
-            partitionRefinement(model, model.getBackwardTransitions(), initialPartition, buildQuotient);
+            partitionRefinement(model, model.getBackwardTransitions(), initialPartition, false, buildQuotient);
         }
         
         template<typename ValueType>
@@ -556,7 +606,7 @@ namespace storm {
             STORM_LOG_THROW(!model.hasStateRewards() && !model.hasTransitionRewards(), storm::exceptions::IllegalFunctionCallException, "Bisimulation is currently only supported for models without reward structures.");
             storm::storage::SparseMatrix<ValueType> backwardTransitions = model.getBackwardTransitions();
             Partition initialPartition = getMeasureDrivenInitialPartition(model, backwardTransitions, phiLabel, psiLabel, bounded);
-            partitionRefinement(model, model.getBackwardTransitions(), initialPartition, buildQuotient);
+            partitionRefinement(model, model.getBackwardTransitions(), initialPartition, false, buildQuotient);
         }
         
         template<typename ValueType>
@@ -648,7 +698,7 @@ namespace storm {
         
         template<typename ValueType>
         template<typename ModelType>
-        void DeterministicModelStrongBisimulationDecomposition<ValueType>::partitionRefinement(ModelType const& model, storm::storage::SparseMatrix<ValueType> const& backwardTransitions, Partition& partition, bool buildQuotient) {
+        void DeterministicModelStrongBisimulationDecomposition<ValueType>::partitionRefinement(ModelType const& model, storm::storage::SparseMatrix<ValueType> const& backwardTransitions, Partition& partition, bool weak, bool buildQuotient) {
             std::chrono::high_resolution_clock::time_point totalStart = std::chrono::high_resolution_clock::now();
 
             // Initially, all blocks are potential splitter, so we insert them in the splitterQueue.
@@ -659,7 +709,7 @@ namespace storm {
             // Then perform the actual splitting until there are no more splitters.
             while (!splitterQueue.empty()) {
                 std::sort(splitterQueue.begin(), splitterQueue.end(), [] (Block const* b1, Block const* b2) { return b1->getNumberOfStates() < b2->getNumberOfStates(); } );
-                refinePartition(backwardTransitions, *splitterQueue.front(), partition, splitterQueue);
+                refinePartition(backwardTransitions, *splitterQueue.front(), partition, weak, splitterQueue);
                 splitterQueue.pop_front();
             }
             std::chrono::high_resolution_clock::duration refinementTime = std::chrono::high_resolution_clock::now() - refinementStart;
@@ -679,6 +729,7 @@ namespace storm {
 
             // If we are required to build the quotient model, do so now.
             if (buildQuotient) {
+                // FIXME: this needs to do a bit more work for weak bisimulation.
                 this->buildQuotient(model, partition);
             }
 
@@ -693,6 +744,7 @@ namespace storm {
                 std::cout << "------------------------------------------" << std::endl;
                 std::cout << "    * total time: " << std::chrono::duration_cast<std::chrono::milliseconds>(totalTime).count() << "ms" << std::endl;
                 std::cout << std::endl;
+                std::cout << "Number of equivalence classes: " << this->size() << std::endl;
             }
         }
         
@@ -713,9 +765,9 @@ namespace storm {
             typename std::vector<std::pair<storm::storage::sparse::state_type, ValueType>>::const_iterator end = partition.getEnd(block) - 1;
             storm::storage::sparse::state_type currentIndex = block.getBegin();
             
-            
             // Now we can check whether the block needs to be split, which is the case iff the probabilities for the
             // first and the last state are different.
+            bool blockSplit = !comparator.isEqual(begin->second, end->second);
             while (!comparator.isEqual(begin->second, end->second)) {
                 // Now we scan for the first state in the block that disagrees on the probability value.
                 // Note that we do not have to check currentIndex for staying within bounds, because we know the matching
@@ -736,13 +788,168 @@ namespace storm {
                     newBlock.markAsSplitter();
                 }
             }
+            
+            // If the block was split, we also need to insert itself into the splitter queue.
+            if (blockSplit) {
+                if (!block.isMarkedAsSplitter()) {
+                    splitterQueue.push_back(&block);
+                    block.markAsSplitter();
+                }
+            }
+        }
+        
+        template<typename ValueType>
+        void DeterministicModelStrongBisimulationDecomposition<ValueType>::refineBlockWeak(Block& block, Partition& partition, storm::storage::SparseMatrix<ValueType> const& backwardTransitions, std::deque<Block*>& splitterQueue) {
+            std::vector<uint_fast64_t> splitPoints = getSplitPointsWeak(block, partition);
+            
+            // Restore the original begin of the block.
+            block.setBegin(block.getOriginalBegin());
+            
+            // Now that we have the split points of the non-silent states, we perform a backward search from
+            // each non-silent state and label the predecessors with the class of the non-silent state.
+            std::cout << "creating labels " << block.getEnd() << " // " << block.getBegin() << std::endl;
+            std::cout << "split points: " << std::endl;
+            for (auto const& point : splitPoints) {
+                std::cout << point << std::endl;
+            }
+            block.print(partition);
+            std::vector<storm::storage::BitVector> stateLabels(block.getEnd() - block.getBegin(), storm::storage::BitVector(splitPoints.size() - 1));
+            
+            std::vector<storm::storage::sparse::state_type> stateStack;
+            stateStack.reserve(block.getEnd() - block.getBegin());
+            for (uint_fast64_t stateClassIndex = 0; stateClassIndex < splitPoints.size() - 1; ++stateClassIndex) {
+                for (auto stateIt = partition.getStatesAndValues().begin() + splitPoints[stateClassIndex], stateIte = partition.getStatesAndValues().begin() + splitPoints[stateClassIndex + 1]; stateIt != stateIte; ++stateIt) {
+                    
+                    stateStack.push_back(stateIt->first);
+                    stateLabels[partition.getPosition(stateIt->first) - block.getBegin()].set(stateClassIndex);
+                    while (!stateStack.empty()) {
+                        storm::storage::sparse::state_type currentState = stateStack.back();
+                        stateStack.pop_back();
+                        
+                        for (auto const& predecessorEntry : backwardTransitions.getRow(currentState)) {
+                            if (comparator.isZero(predecessorEntry.getValue())) {
+                                continue;
+                            }
+                            
+                            storm::storage::sparse::state_type predecessor = predecessorEntry.getColumn();
+                            
+                            // Only if the state is in the same block, is a silent state and it has not yet been
+                            // labeled with the current label.
+                            if (&partition.getBlock(predecessor) == &block && partition.isSilent(predecessor, comparator) && !stateLabels[partition.getPosition(predecessor) - block.getBegin()].get(stateClassIndex)) {
+                                stateStack.push_back(predecessor);
+                                stateLabels[partition.getPosition(predecessor) - block.getBegin()].set(stateClassIndex);
+                            }
+                        }
+                    }
+                }
+            }
+            
+            // Now that all silent states were appropriately labeled, we can sort the states according to their
+            // labels and then scan for ranges that agree on the label.
+            std::sort(partition.getBegin(block), partition.getEnd(block), [&] (std::pair<storm::storage::sparse::state_type, ValueType> const& a, std::pair<storm::storage::sparse::state_type, ValueType> const& b) { return stateLabels[partition.getPosition(a.first) - block.getBegin()] < stateLabels[partition.getPosition(b.first) - block.getBegin()]; });
+            
+            // Update the positions vector.
+            storm::storage::sparse::state_type position = block.getBegin();
+            for (auto stateIt = partition.getBegin(block), stateIte = partition.getEnd(block); stateIt != stateIte; ++stateIt, ++position) {
+                partition.setPosition(stateIt->first, position);
+            }
+            
+            // Now we have everything in place to actually split the block by just scanning for ranges of equal label.
+            typename std::vector<std::pair<storm::storage::sparse::state_type, ValueType>>::const_iterator begin = partition.getBegin(block);
+            typename std::vector<std::pair<storm::storage::sparse::state_type, ValueType>>::const_iterator current = begin;
+            typename std::vector<std::pair<storm::storage::sparse::state_type, ValueType>>::const_iterator end = partition.getEnd(block) - 1;
+            storm::storage::sparse::state_type currentIndex = block.getBegin();
+            
+            // Now we can check whether the block needs to be split, which is the case iff the probabilities for the
+            // first and the last state are different.
+            bool blockSplit = !comparator.isEqual(begin->second, end->second);
+            while (!comparator.isEqual(begin->second, end->second)) {
+                // Now we scan for the first state in the block that disagrees on the labeling value.
+                // Note that we do not have to check currentIndex for staying within bounds, because we know the matching
+                // state is within bounds.
+                ValueType const& currentValue = begin->second;
+                
+                ++begin;
+                ++currentIndex;
+                while (begin != end && comparator.isEqual(begin->second, currentValue)) {
+                    ++begin;
+                    ++currentIndex;
+                }
+                
+                // Now we split the block and mark it as a potential splitter.
+                Block& newBlock = partition.splitBlock(block, currentIndex);
+                if (!newBlock.isMarkedAsSplitter()) {
+                    splitterQueue.push_back(&newBlock);
+                    newBlock.markAsSplitter();
+                }
+            }
+            
+            // If the block was split, we also need to insert itself into the splitter queue.
+            if (blockSplit) {
+                if (!block.isMarkedAsSplitter()) {
+                    splitterQueue.push_back(&block);
+                    block.markAsSplitter();
+                }
+            }
         }
         
         template<typename ValueType>
-        void DeterministicModelStrongBisimulationDecomposition<ValueType>::refinePartition(storm::storage::SparseMatrix<ValueType> const& backwardTransitions, Block& splitter, Partition& partition, std::deque<Block*>& splitterQueue) {
+        std::vector<uint_fast64_t> DeterministicModelStrongBisimulationDecomposition<ValueType>::getSplitPointsWeak(Block& block, Partition& partition) {
+            std::vector<uint_fast64_t> result;
+            // We first scale all probabilities with (1-p[s]) where p[s] is the silent probability of state s.
+            std::for_each(partition.getStatesAndValues().begin() + block.getOriginalBegin(), partition.getStatesAndValues().begin() + block.getBegin(), [&] (std::pair<storm::storage::sparse::state_type, ValueType>& stateValuePair) {
+                ValueType const& silentProbability = partition.getSilentProbability(stateValuePair.first);
+                if (!comparator.isOne(silentProbability)) {
+                    stateValuePair.second /= storm::utility::one<ValueType>() - silentProbability;
+                }
+            });
+            
+            // Now sort the states based on their probabilities.
+            std::sort(partition.getStatesAndValues().begin() + block.getOriginalBegin(), partition.getStatesAndValues().begin() + block.getBegin(), [&partition] (std::pair<storm::storage::sparse::state_type, ValueType> const& a, std::pair<storm::storage::sparse::state_type, ValueType> const& b) { return a.second < b.second; } );
+            
+            // Update the positions vector.
+            storm::storage::sparse::state_type position = block.getOriginalBegin();
+            for (auto stateIt = partition.getStatesAndValues().begin() + block.getOriginalBegin(), stateIte = partition.getStatesAndValues().begin() + block.getBegin(); stateIt != stateIte; ++stateIt, ++position) {
+                partition.setPosition(stateIt->first, position);
+            }
+            
+            // Then, we scan for the ranges of states that agree on the probability.
+            typename std::vector<std::pair<storm::storage::sparse::state_type, ValueType>>::const_iterator begin = partition.getStatesAndValues().begin() + block.getOriginalBegin();
+            typename std::vector<std::pair<storm::storage::sparse::state_type, ValueType>>::const_iterator current = begin;
+            typename std::vector<std::pair<storm::storage::sparse::state_type, ValueType>>::const_iterator end = partition.getStatesAndValues().begin() + block.getBegin() - 1;
+            storm::storage::sparse::state_type currentIndex = block.getOriginalBegin();
+            result.push_back(currentIndex);
+
+            // Now we can check whether the block needs to be split, which is the case iff the probabilities for the
+            // first and the last state are different.
+            while (!comparator.isEqual(begin->second, end->second)) {
+                // Now we scan for the first state in the block that disagrees on the probability value.
+                // Note that we do not have to check currentIndex for staying within bounds, because we know the matching
+                // state is within bounds.
+                ValueType const& currentValue = begin->second;
+                
+                ++begin;
+                ++currentIndex;
+                while (begin != end && comparator.isEqual(begin->second, currentValue)) {
+                    ++begin;
+                    ++currentIndex;
+                }
+
+                // Remember the index at which the probabilities were different.
+                result.push_back(currentIndex);
+            }
+            
+            // Push a sentinel element and return result.
+            result.push_back(block.getEnd());
+            return result;
+        }
+        
+        template<typename ValueType>
+        void DeterministicModelStrongBisimulationDecomposition<ValueType>::refinePartition(storm::storage::SparseMatrix<ValueType> const& backwardTransitions, Block& splitter, Partition& partition, bool weak, std::deque<Block*>& splitterQueue) {
             std::list<Block*> predecessorBlocks;
             
             // Iterate over all states of the splitter and check its predecessors.
+            bool splitterIsPredecessor = false;
             storm::storage::sparse::state_type currentPosition = splitter.getBegin();
             for (auto stateIterator = partition.getBegin(splitter), stateIte = partition.getEnd(splitter); stateIterator != stateIte; ++stateIterator, ++currentPosition) {
                 storm::storage::sparse::state_type currentState = stateIterator->first;
@@ -751,7 +958,11 @@ namespace storm {
                 for (auto const& predecessorEntry : backwardTransitions.getRow(currentState)) {
                     storm::storage::sparse::state_type predecessor = predecessorEntry.getColumn();
                     
+                    // Get predecessor block and remember if the splitter was a predecessor of itself.
                     Block& predecessorBlock = partition.getBlock(predecessor);
+                    if (&predecessorBlock == &splitter) {
+                        splitterIsPredecessor = true;
+                    }
                     
                     // If the predecessor block has just one state, there is no point in splitting it.
                     if (predecessorBlock.getNumberOfStates() <= 1 || predecessorBlock.isAbsorbing()) {
@@ -831,45 +1042,68 @@ namespace storm {
                     }
                 }
             }
-            
-            // Reset the marked position of the splitter to the begin.
-            splitter.setMarkedPosition(splitter.getBegin());
-            
-            std::list<Block*> blocksToSplit;
-            
-            // Now, we can iterate over the predecessor blocks and see whether we have to create a new block for
-            // predecessors of the splitter.
-            for (auto blockPtr : predecessorBlocks) {
-                Block& block = *blockPtr;
+
+            if (!weak) {
+                std::list<Block*> blocksToSplit;
                 
-                block.unmarkAsPredecessorBlock();
-                block.resetMarkedPosition();
+                // Now, we can iterate over the predecessor blocks and see whether we have to create a new block for
+                // predecessors of the splitter.
+                for (auto blockPtr : predecessorBlocks) {
+                    Block& block = *blockPtr;
+                    
+                    block.unmarkAsPredecessorBlock();
+                    block.resetMarkedPosition();
+                    
+                    // If we have moved the begin of the block to somewhere in the middle of the block, we need to split it.
+                    if (block.getBegin() != block.getEnd()) {
+                        Block& newBlock = partition.insertBlock(block);
+                        if (!newBlock.isMarkedAsSplitter()) {
+                            splitterQueue.push_back(&newBlock);
+                            newBlock.markAsSplitter();
+                        }
+                        
+                        // Schedule the block of predecessors for refinement based on probabilities.
+                        blocksToSplit.emplace_back(&newBlock);
+                    } else {
+                        // In this case, we can keep the block by setting its begin to the old value.
+                        block.setBegin(block.getOriginalBegin());
+                        blocksToSplit.emplace_back(&block);
+                    }
+                }
                 
-                // If we have moved the begin of the block to somewhere in the middle of the block, we need to split it.
-                if (block.getBegin() != block.getEnd()) {
-                    Block& newBlock = partition.insertBlock(block);
-                    if (!newBlock.isMarkedAsSplitter()) {
-                        splitterQueue.push_back(&newBlock);
-                        newBlock.markAsSplitter();
+                // Finally, we walk through the blocks that have a transition to the splitter and split them using
+                // probabilistic information.
+                for (auto blockPtr : blocksToSplit) {
+                    if (blockPtr->getNumberOfStates() <= 1) {
+                        continue;
                     }
                     
-                    // Schedule the block of predecessors for refinement based on probabilities.
-                    blocksToSplit.emplace_back(&newBlock);
-                } else {
-                    // In this case, we can keep the block by setting its begin to the old value.
-                    block.setBegin(block.getOriginalBegin());
-                    blocksToSplit.emplace_back(&block);
+                    refineBlockProbabilities(*blockPtr, partition, splitterQueue);
                 }
-            }
-            
-            // Finally, we walk through the blocks that have a transition to the splitter and split them using
-            // probabilistic information.
-            for (auto blockPtr : blocksToSplit) {
-                if (blockPtr->getNumberOfStates() <= 1) {
-                    continue;
+            } else {
+                // If the splitter was a predecessor of itself and we are computing a weak bisimulation, we need to update
+                // the silent probabilities.
+                if (splitterIsPredecessor) {
+                    partition.setSilentProbabilities(partition.getStatesAndValues().begin() + splitter.getOriginalBegin(), partition.getStatesAndValues().begin() + splitter.getBegin());
+                    partition.setSilentProbabilitiesToZero(partition.getStatesAndValues().begin() + splitter.getBegin(), partition.getStatesAndValues().begin() + splitter.getEnd());
+                }
+
+                // Now refine all predecessor blocks in a weak manner. That is, we split according to the criterion of
+                // weak bisimulation.
+                for (auto blockPtr : predecessorBlocks) {
+                    Block& block = *blockPtr;
+                    
+                    // If the splitter is also the predecessor block, we must not refine it at this point.
+                    if (&block != &splitter) {
+                        refineBlockWeak(block, partition, backwardTransitions, splitterQueue);
+                    } else {
+                        // Restore the begin of the block.
+                        block.setBegin(block.getOriginalBegin());
+                    }
+
+                    block.unmarkAsPredecessorBlock();
+                    block.resetMarkedPosition();
                 }
-                
-                refineBlockProbabilities(*blockPtr, partition, splitterQueue);
             }
         }
         
@@ -883,8 +1117,8 @@ namespace storm {
         
         template<typename ValueType>
         template<typename ModelType>
-        typename DeterministicModelStrongBisimulationDecomposition<ValueType>::Partition DeterministicModelStrongBisimulationDecomposition<ValueType>::getLabelBasedInitialPartition(ModelType const& model) {
-            Partition partition(model.getNumberOfStates());
+        typename DeterministicModelStrongBisimulationDecomposition<ValueType>::Partition DeterministicModelStrongBisimulationDecomposition<ValueType>::getLabelBasedInitialPartition(ModelType const& model, bool weak) {
+            Partition partition(model.getNumberOfStates(), weak);
             for (auto const& label : model.getStateLabeling().getAtomicPropositions()) {
                 if (label == "init") {
                     continue;
@@ -894,6 +1128,24 @@ namespace storm {
             return partition;
         }
         
+        template<typename ValueType>
+        template<typename ModelType>
+        void DeterministicModelStrongBisimulationDecomposition<ValueType>::initializeSilentProbabilities(ModelType const& model, Partition& partition) {
+            for (auto const& block : partition.getBlocks()) {
+                for (auto stateIt = partition.getBegin(block), stateIte = partition.getEnd(block); stateIt != stateIte; ++stateIt) {
+                    ValueType silentProbability = storm::utility::zero<ValueType>();
+                    
+                    for (auto const& successorEntry : model.getRows(stateIt->first)) {
+                        if (&partition.getBlock(successorEntry.getColumn()) == &block) {
+                            silentProbability += successorEntry.getValue();
+                        }
+                    }
+                    
+                    partition.setSilentProbability(stateIt->first, silentProbability);
+                }
+            }
+        }
+        
         template class DeterministicModelStrongBisimulationDecomposition<double>;
     }
 }

src/storage/DeterministicModelStrongBisimulationDecomposition.h

@@ -22,20 +22,20 @@ namespace storm {
         class DeterministicModelStrongBisimulationDecomposition : public Decomposition<StateBlock> {
         public:
             /*!
-             * Decomposes the given DTMC into equivalence classes under strong bisimulation.
+             * Decomposes the given DTMC into equivalence classes under weak or strong bisimulation.
              *
              * @param model The model to decompose.
              * @param buildQuotient Sets whether or not the quotient model is to be built.
              */
-            DeterministicModelStrongBisimulationDecomposition(storm::models::Dtmc<ValueType> const& model, bool buildQuotient = false);
+            DeterministicModelStrongBisimulationDecomposition(storm::models::Dtmc<ValueType> const& model, bool weak = false, bool buildQuotient = false);
 
             /*!
-             * Decomposes the given CTMC into equivalence classes under strong bisimulation.
+             * Decomposes the given CTMC into equivalence classes under weak or strong bisimulation.
              *
              * @param model The model to decompose.
              * @param buildQuotient Sets whether or not the quotient model is to be built.
              */
-            DeterministicModelStrongBisimulationDecomposition(storm::models::Ctmc<ValueType> const& model, bool buildQuotient = false);
+            DeterministicModelStrongBisimulationDecomposition(storm::models::Ctmc<ValueType> const& model, bool weak = false, bool buildQuotient = false);
             
             /*!
              * Decomposes the given DTMC into equivalence classes under strong bisimulation in a way that onle safely
@@ -340,11 +340,29 @@ namespace storm {
                 void increaseValue(storm::storage::sparse::state_type state, ValueType value);
                 
                 // Updates the block mapping for the given range of states to the specified block.
-                void updateBlockMapping(Block& block, typename std::vector<std::pair<storm::storage::sparse::state_type, ValueType>>::iterator first, typename std::vector<std::pair<storm::storage::sparse::state_type, ValueType>>::iterator end);
+                void updateBlockMapping(Block& block, typename std::vector<std::pair<storm::storage::sparse::state_type, ValueType>>::iterator first, typename std::vector<std::pair<storm::storage::sparse::state_type, ValueType>>::iterator last);
                 
                 // Retrieves the first block of the partition.
                 Block& getFirstBlock();
                 
+                // Retrieves whether the given state is fully silent (only in case the silent probabilities are tracked).
+                bool isSilent(storm::storage::sparse::state_type state, storm::utility::ConstantsComparator<ValueType> const& comparator) const;
+                
+                // Retrieves whether the given state has a non-zero silent probability.
+                bool hasSilentProbability(storm::storage::sparse::state_type state, storm::utility::ConstantsComparator<ValueType> const& comparator) const;
+                
+                // Retrieves the silent probability (i.e. the probability to stay within the own equivalence class).
+                ValueType const& getSilentProbability(storm::storage::sparse::state_type state) const;
+                
+                // Sets the silent probabilities for all the states in the range to their values in the range.
+                void setSilentProbabilities(typename std::vector<std::pair<storm::storage::sparse::state_type, ValueType>>::iterator first, typename std::vector<std::pair<storm::storage::sparse::state_type, ValueType>>::iterator last);
+                
+                // Sets the silent probabilities for all states in the range to zero.
+                void setSilentProbabilitiesToZero(typename std::vector<std::pair<storm::storage::sparse::state_type, ValueType>>::iterator first, typename std::vector<std::pair<storm::storage::sparse::state_type, ValueType>>::iterator last);
+
+                // Sets the silent probability for the given state to the given value.
+                void setSilentProbability(storm::storage::sparse::state_type state, ValueType const& value);
+
             private:
                 // The list of blocks in the partition.
                 std::list<Block> blocks;
@@ -375,11 +393,12 @@ namespace storm {
              * @param model The model on whose state space to compute the coarses strong bisimulation relation.
              * @param backwardTransitions The backward transitions of the model.
              * @param The initial partition.
+             * @param weak A flag indicating whether a weak or a strong bisimulation is to be computed.
              * @param buildQuotient If set, the quotient model is built and may be retrieved using the getQuotient()
              * method.
              */
             template<typename ModelType>
-            void partitionRefinement(ModelType const& model, storm::storage::SparseMatrix<ValueType> const& backwardTransitions, Partition& partition, bool buildQuotient);
+            void partitionRefinement(ModelType const& model, storm::storage::SparseMatrix<ValueType> const& backwardTransitions, Partition& partition, bool weak, bool buildQuotient);
             
             /*!
              * Refines the partition based on the provided splitter. After calling this method all blocks are stable
@@ -389,10 +408,11 @@ namespace storm {
              * probabilities).
              * @param splitter The splitter to use.
              * @param partition The partition to split.
+             * @param weak A flag indicating whether a weak or a strong bisimulation is to be computed.
              * @param splitterQueue A queue into which all blocks that were split are inserted so they can be treated
              * as splitters in the future.
              */
-            void refinePartition(storm::storage::SparseMatrix<ValueType> const& backwardTransitions, Block& splitter, Partition& partition, std::deque<Block*>& splitterQueue);
+            void refinePartition(storm::storage::SparseMatrix<ValueType> const& backwardTransitions, Block& splitter, Partition& partition, bool weak, std::deque<Block*>& splitterQueue);
             
             /*!
              * Refines the block based on their probability values (leading into the splitter).
@@ -404,6 +424,16 @@ namespace storm {
              */
             void refineBlockProbabilities(Block& block, Partition& partition, std::deque<Block*>& splitterQueue);
             
+            void refineBlockWeak(Block& block, Partition& partition, storm::storage::SparseMatrix<ValueType> const& backwardTransitions, std::deque<Block*>& splitterQueue);
+            
+            /*!
+             * Determines the split offsets in the given block.
+             *
+             * @param block The block that is to be analyzed for splits.
+             * @param partition The partition that contains the block.
+             */
+            std::vector<uint_fast64_t> getSplitPointsWeak(Block& block, Partition& partition);
+            
             /*!
              * Builds the quotient model based on the previously computed equivalence classes (stored in the blocks
              * of the decomposition.
@@ -435,10 +465,21 @@ namespace storm {
              * Creates the initial partition based on all the labels in the given model.
              *
              * @param model The model whose state space is partitioned based on its labels.
+             * @param weak A flag indicating whether a weak bisimulation is to be computed.
              * @return The resulting partition.
              */
             template<typename ModelType>
-            Partition getLabelBasedInitialPartition(ModelType const& model);
+            Partition getLabelBasedInitialPartition(ModelType const& model, bool weak);
+            
+            /*!
+             * Initializes the silent probabilities by traversing all blocks and adding the probability of going to
+             * the very own equivalence class for each state.
+             *
+             * @param model The model from which to look-up the probabilities.
+             * @param partition The partition that holds the silent probabilities.
+             */
+            template<typename ModelType>
+            void initializeSilentProbabilities(ModelType const& model, Partition& partition);
             
             // If required, a quotient model is built and stored in this member.
             std::shared_ptr<storm::models::AbstractDeterministicModel<ValueType>> quotient;

src/utility/cli.h

@@ -267,7 +267,7 @@ namespace storm {
                     std::shared_ptr<storm::models::Dtmc<double>> dtmc = result->template as<storm::models::Dtmc<double>>();
                     
                     STORM_PRINT(std::endl << "Performing bisimulation minimization..." << std::endl);
-                    storm::storage::DeterministicModelStrongBisimulationDecomposition<double> bisimulationDecomposition(*dtmc, true);
+                    storm::storage::DeterministicModelStrongBisimulationDecomposition<double> bisimulationDecomposition(*dtmc, storm::settings::bisimulationSettings().isWeakBisimulationSet(), true);
                     
                     result = bisimulationDecomposition.getQuotient();