diff --git a/src/settings/modules/GeneralSettings.cpp b/src/settings/modules/GeneralSettings.cpp index b10204e29..f251aad38 100644 --- a/src/settings/modules/GeneralSettings.cpp +++ b/src/settings/modules/GeneralSettings.cpp @@ -74,7 +74,7 @@ namespace storm { .addArgument(storm::settings::ArgumentBuilder::createStringArgument("filename", "The file from which to read the LTL formulas.").addValidationFunctionString(storm::settings::ArgumentValidators::existingReadableFileValidator()).build()).build()); this->addOption(storm::settings::OptionBuilder(moduleName, counterexampleOptionName, false, "Generates a counterexample for the given PRCTL formulas if not satisfied by the model") .addArgument(storm::settings::ArgumentBuilder::createStringArgument("filename", "The name of the file to which the counterexample is to be written.").setDefaultValueString("-").setIsOptional(true).build()).setShortName(counterexampleOptionShortName).build()); - this->addOption(storm::settings::OptionBuilder(moduleName, bisimulationOptionName, true, "Sets whether to perform bisimulation minimization.").setShortName(bisimulationOptionShortName).build()); + this->addOption(storm::settings::OptionBuilder(moduleName, bisimulationOptionName, false, "Sets whether to perform bisimulation minimization.").setShortName(bisimulationOptionShortName).build()); this->addOption(storm::settings::OptionBuilder(moduleName, transitionRewardsOptionName, "", "If given, the transition rewards are read from this file and added to the explicit model. Note that this requires the model to be given as an explicit model (i.e., via --" + explicitOptionName + ").") .addArgument(storm::settings::ArgumentBuilder::createStringArgument("filename", "The file from which to read the transition rewards.").addValidationFunctionString(storm::settings::ArgumentValidators::existingReadableFileValidator()).build()).build()); this->addOption(storm::settings::OptionBuilder(moduleName, stateRewardsOptionName, false, "If given, the state rewards are read from this file and added to the explicit model. Note that this requires the model to be given as an explicit model (i.e., via --" + explicitOptionName + ").") @@ -92,7 +92,7 @@ namespace storm { .addArgument(storm::settings::ArgumentBuilder::createStringArgument("name", "The name of an LP solver. Available are: gurobi and glpk.").addValidationFunctionString(storm::settings::ArgumentValidators::stringInListValidator(lpSolvers)).setDefaultValueString("glpk").build()).build()); this->addOption(storm::settings::OptionBuilder(moduleName, constantsOptionName, false, "Specifies the constant replacements to use in symbolic models. Note that Note that this requires the model to be given as an symbolic model (i.e., via --" + symbolicOptionName + ").").setShortName(constantsOptionShortName) .addArgument(storm::settings::ArgumentBuilder::createStringArgument("values", "A comma separated list of constants and their value, e.g. a=1,b=2,c=3.").setDefaultValueString("").build()).build()); - this->addOption(storm::settings::OptionBuilder(moduleName, statisticsOptionName, true, "Sets whether to display statistics if available.").setShortName(statisticsOptionShortName).build()); + this->addOption(storm::settings::OptionBuilder(moduleName, statisticsOptionName, false, "Sets whether to display statistics if available.").setShortName(statisticsOptionShortName).build()); } bool GeneralSettings::isHelpSet() const { diff --git a/src/storage/BisimulationDecomposition.cpp b/src/storage/BisimulationDecomposition.cpp index ba2e4111b..53459fb38 100644 --- a/src/storage/BisimulationDecomposition.cpp +++ b/src/storage/BisimulationDecomposition.cpp @@ -1,48 +1,189 @@ #include "src/storage/BisimulationDecomposition.h" -#include +#include +#include namespace storm { namespace storage { template - BisimulationDecomposition::BisimulationDecomposition(storm::models::Dtmc const& model, bool weak) { - computeBisimulationEquivalenceClasses(model, weak); + BisimulationDecomposition::BisimulationDecomposition(storm::models::Dtmc const& dtmc, bool weak) { + computeBisimulationEquivalenceClasses(dtmc, weak); } template - void BisimulationDecomposition::computeBisimulationEquivalenceClasses(storm::models::Dtmc const& model, bool weak) { + void BisimulationDecomposition::computeBisimulationEquivalenceClasses(storm::models::Dtmc const& dtmc, bool weak) { + std::chrono::high_resolution_clock::time_point totalStart = std::chrono::high_resolution_clock::now(); // We start by computing the initial partition. In particular, we also keep a mapping of states to their blocks. - std::vector stateToBlockMapping(model.getNumberOfStates()); - storm::storage::BitVector labeledStates = model.getLabeledStates("one"); + std::vector stateToBlockMapping(dtmc.getNumberOfStates()); + storm::storage::BitVector labeledStates = dtmc.getLabeledStates("observe0Greater1"); this->blocks.emplace_back(labeledStates.begin(), labeledStates.end()); std::for_each(labeledStates.begin(), labeledStates.end(), [&] (storm::storage::sparse::state_type const& state) { stateToBlockMapping[state] = 0; } ); labeledStates.complement(); this->blocks.emplace_back(labeledStates.begin(), labeledStates.end()); std::for_each(labeledStates.begin(), labeledStates.end(), [&] (storm::storage::sparse::state_type const& state) { stateToBlockMapping[state] = 1; } ); + // Create empty distributions for the two initial blocks. + std::vector> distributions(2); + // Retrieve the backward transitions to allow for better checking of states that need to be re-examined. - storm::storage::SparseMatrix const& backwardTransitions = model.getBackwardTransitions(); + storm::storage::SparseMatrix const& backwardTransitions = dtmc.getBackwardTransitions(); - // Initially, both blocks are potential splitters. - std::queue splitterQueue; - splitterQueue.push(0); - splitterQueue.push(1); - - // As long as there is a splitter, we keep refining the current partition. - while (!splitterQueue.empty()) { + // Initially, both blocks are potential splitters. A splitter is marked as a pair in which the first entry + // is the ID of the parent block of the splitter and the second entry is the block ID of the splitter itself. + std::deque refinementQueue; + storm::storage::BitVector blocksInRefinementQueue(this->size()); + refinementQueue.push_back(0); + refinementQueue.push_back(1); + + // As long as there are blocks to refine, well, refine them. + uint_fast64_t iteration = 0; + while (!refinementQueue.empty()) { + ++iteration; + + // Optionally sort the queue to potentially speed up the convergence. + // std::sort(refinementQueue.begin(), refinementQueue.end(), [=] (std::size_t const& a, std::size_t const& b) { return this->blocks[a].size() > this->blocks[b].size(); }); + std::size_t currentBlock = refinementQueue.front(); + refinementQueue.pop_front(); + blocksInRefinementQueue.set(currentBlock, false); + + splitBlock(dtmc, backwardTransitions, currentBlock, stateToBlockMapping, distributions, blocksInRefinementQueue, refinementQueue); } - // While there is a splitter... - // - // check the predecessors of the splitter: - // * if they still have the same signature as before, then they remain unsplit - // * otherwise, split off the states that now behave differently - // and mark the smaller block as a splitter + std::chrono::high_resolution_clock::duration totalTime = std::chrono::high_resolution_clock::now() - totalStart; - // + std::cout << "Bisimulation quotient has " << this->blocks.size() << " blocks and took " << iteration << " iterations and " << std::chrono::duration_cast(totalTime).count() << "ms." << std::endl; + } + + template + std::size_t BisimulationDecomposition::splitPredecessorsGraphBased(storm::models::Dtmc const& dtmc, storm::storage::SparseMatrix const& backwardTransitions, std::size_t const& blockId, std::vector& stateToBlockMapping, std::vector>& distributions, storm::storage::BitVector& blocksInRefinementQueue, std::deque& graphRefinementQueue) { + + // This method tries to split the blocks of predecessors of the provided block by checking whether there is + // a transition into the current block or not. + std::unordered_map::block_type> predecessorBlockToNewBlock; + + // Now for each predecessor block which state could actually reach the current block. + for (auto const& state : this->blocks[blockId]) { + for (auto const& predecessorEntry : backwardTransitions.getRow(state)) { + storm::storage::sparse::state_type predecessor = predecessorEntry.getColumn(); + predecessorBlockToNewBlock[stateToBlockMapping[predecessor]].insert(predecessor); + } + } + + // Now, we can check for each predecessor block whether it needs to be split. + for (auto const& blockNewBlockPair : predecessorBlockToNewBlock) { + if (this->blocks[blockNewBlockPair.first].size() > blockNewBlockPair.second.size()) { + // Add the states which have a successor in the current block to a totally new block. + this->blocks.emplace_back(std::move(blockNewBlockPair.second)); + + // Compute the set of states that remains in the old block; + typename BisimulationDecomposition::block_type newBlock; + std::set_difference(this->blocks[blockId].begin(), this->blocks[blockId].end(), this->blocks.back().begin(), this->blocks.back().end(), std::inserter(newBlock, newBlock.begin())); + this->blocks[blockNewBlockPair.first] = std::move(newBlock); + + blocksInRefinementQueue.resize(this->blocks.size()); + + // Add the smaller part of the old block to the queue. + std::size_t blockToAddToQueue = this->blocks.back().size() < this->blocks[blockNewBlockPair.first].size() ? this->blocks.size() - 1 : blockNewBlockPair.first; + if (!blocksInRefinementQueue.get(blockToAddToQueue)) { + graphRefinementQueue.push_back(blockToAddToQueue); + blocksInRefinementQueue.set(blockToAddToQueue, true); + } + } + } + + // FIXME + return 0; + } + + template + std::size_t BisimulationDecomposition::splitBlock(storm::models::Dtmc const& dtmc, storm::storage::SparseMatrix const& backwardTransitions, std::size_t const& blockId, std::vector& stateToBlockMapping, std::vector>& distributions, storm::storage::BitVector& blocksInRefinementQueue, std::deque& refinementQueue) { + std::chrono::high_resolution_clock::time_point totalStart = std::chrono::high_resolution_clock::now(); + std::unordered_map, typename BisimulationDecomposition::block_type> distributionToNewBlocks; + + // Traverse all states of the block and check whether they have different distributions. + for (auto const& state : this->blocks[blockId]) { + // Now construct the distribution of this state wrt. to the current partitioning. + storm::storage::Distribution distribution; + for (auto const& successorEntry : dtmc.getTransitionMatrix().getRow(state)) { + distribution.addProbability(static_cast(stateToBlockMapping[successorEntry.getColumn()]), successorEntry.getValue()); + } + + // If the distribution already exists, we simply add the state. Otherwise, we open a new block. + auto distributionIterator = distributionToNewBlocks.find(distribution); + if (distributionIterator != distributionToNewBlocks.end()) { + distributionIterator->second.insert(state); + } else { + distributionToNewBlocks[distribution].insert(state); + } + } + + // Now we are ready to split the block. + if (distributionToNewBlocks.size() == 1) { + // If there is just one behavior, we just set the distribution as the new one for this block. + distributions[blockId] = std::move(distributionToNewBlocks.begin()->first); + } else { + // In this case, we need to split the block. + typename BisimulationDecomposition::block_type tmpBlock; + + auto distributionIterator = distributionToNewBlocks.begin(); + distributions[blockId] = std::move(distributionIterator->first); + tmpBlock = std::move(distributionIterator->second); + std::swap(this->blocks[blockId], tmpBlock); + ++distributionIterator; + + // Remember the number of blocks prior to splitting for later use. + std::size_t beforeNumberOfBlocks = this->blocks.size(); + + for (; distributionIterator != distributionToNewBlocks.end(); ++distributionIterator) { + // In this case, we need to move the newly created block to the end of the list of actual blocks. + this->blocks.emplace_back(std::move(distributionIterator->second)); + distributions.emplace_back(std::move(distributionIterator->first)); + + // Update the mapping of states to their blocks. + std::size_t newBlockId = this->blocks.size() - 1; + for (auto const& state : this->blocks.back()) { + stateToBlockMapping[state] = newBlockId; + } + } + + // Now that we have split the block, we try to trigger a chain of graph-based splittings. That is, we + // try to split the predecessors of the current block by checking whether they have some transition + // to one given sub-block of the current block. + std::deque localRefinementQueue; + storm::storage::BitVector blocksInLocalRefinementQueue(this->size()); + localRefinementQueue.push_back(blockId); + for (std::size_t i = beforeNumberOfBlocks; i < this->blocks.size(); ++i) { + localRefinementQueue.push_back(i); + } + + while (!localRefinementQueue.empty()) { + std::size_t currentBlock = localRefinementQueue.front(); + localRefinementQueue.pop_front(); + blocksInLocalRefinementQueue.set(currentBlock, false); + + splitPredecessorsGraphBased(dtmc, backwardTransitions, blockId, stateToBlockMapping, distributions, blocksInLocalRefinementQueue, localRefinementQueue); + } + + // Since we created some new blocks, we need to extend the bit vector storing the blocks in the + // refinement queue. + blocksInRefinementQueue.resize(blocksInRefinementQueue.size() + (distributionToNewBlocks.size() - 1)); + + // Insert blocks that possibly need a refinement into the queue. + for (auto const& state : tmpBlock) { + for (auto const& predecessor : backwardTransitions.getRow(state)) { + if (!blocksInRefinementQueue.get(stateToBlockMapping[predecessor.getColumn()])) { + blocksInRefinementQueue.set(stateToBlockMapping[predecessor.getColumn()]); + refinementQueue.push_back(stateToBlockMapping[predecessor.getColumn()]); + } + } + } + } + std::chrono::high_resolution_clock::duration totalTime = std::chrono::high_resolution_clock::now() - totalStart; + std::cout << "refinement of block " << blockId << " took " << std::chrono::duration_cast(totalTime).count() << "ms." << std::endl; + return distributionToNewBlocks.size(); } template class BisimulationDecomposition; diff --git a/src/storage/BisimulationDecomposition.h b/src/storage/BisimulationDecomposition.h index 1d2535b2e..d02d02bc1 100644 --- a/src/storage/BisimulationDecomposition.h +++ b/src/storage/BisimulationDecomposition.h @@ -1,8 +1,12 @@ #ifndef STORM_STORAGE_BISIMULATIONDECOMPOSITION_H_ #define STORM_STORAGE_BISIMULATIONDECOMPOSITION_H_ +#include +#include + #include "src/storage/Decomposition.h" #include "src/models/Dtmc.h" +#include "src/storage/Distribution.h" namespace storm { namespace storage { @@ -13,15 +17,18 @@ namespace storm { template class BisimulationDecomposition : public Decomposition { public: + BisimulationDecomposition() = default; + /*! * Decomposes the given DTMC into equivalence classes under weak or strong bisimulation. */ - BisimulationDecomposition(storm::models::Dtmc const& model, bool weak); + BisimulationDecomposition(storm::models::Dtmc const& model, bool weak = false); private: void computeBisimulationEquivalenceClasses(storm::models::Dtmc const& model, bool weak); + std::size_t splitPredecessorsGraphBased(storm::models::Dtmc const& dtmc, storm::storage::SparseMatrix const& backwardTransitions, std::size_t const& block, std::vector& stateToBlockMapping, std::vector>& distributions, storm::storage::BitVector& blocksInRefinementQueue, std::deque& refinementQueue); + std::size_t splitBlock(storm::models::Dtmc const& dtmc, storm::storage::SparseMatrix const& backwardTransitions, std::size_t const& block, std::vector& stateToBlockMapping, std::vector>& distributions, storm::storage::BitVector& blocksInRefinementQueue, std::deque& refinementQueue); }; - } } diff --git a/src/storage/Decomposition.cpp b/src/storage/Decomposition.cpp index 06656b483..39c49bf4f 100644 --- a/src/storage/Decomposition.cpp +++ b/src/storage/Decomposition.cpp @@ -89,7 +89,10 @@ namespace storm { out << "]"; return out; } - + + template class Decomposition; + template std::ostream& operator<<(std::ostream& out, Decomposition const& decomposition); + template class Decomposition; template std::ostream& operator<<(std::ostream& out, Decomposition const& decomposition); diff --git a/src/storage/Distribution.cpp b/src/storage/Distribution.cpp index 84bda0abf..55de31d51 100644 --- a/src/storage/Distribution.cpp +++ b/src/storage/Distribution.cpp @@ -1,13 +1,48 @@ #include "src/storage/Distribution.h" #include +#include + +#include "src/settings/SettingsManager.h" namespace storm { namespace storage { + + template + Distribution::Distribution() : hash(0) { + // Intentionally left empty. + } + + template + bool Distribution::operator==(Distribution const& other) const { + // We need to check equality by ourselves, because we need to account for epsilon differences. + if (this->distribution.size() != other.distribution.size() || this->getHash() != other.getHash()) { + return false; + } + + auto first1 = this->distribution.begin(); + auto last1 = this->distribution.end(); + auto first2 = other.distribution.begin(); + auto last2 = other.distribution.end(); + + for (; first1 != last1; ++first1, ++first2) { + if (first1->first != first2->first) { + return false; + } + if (std::abs(first1->second - first2->second) > 1e-6) { + return false; + } + } + + return true; + } + template void Distribution::addProbability(storm::storage::sparse::state_type const& state, ValueType const& probability) { + if (this->distribution.find(state) == this->distribution.end()) { + this->hash += static_cast(state); + } this->distribution[state] += probability; - this->hash += static_cast(state); } template @@ -30,6 +65,23 @@ namespace storm { return this->distribution.end(); } + template + std::size_t Distribution::getHash() const { + return this->hash ^ (this->distribution.size() << 8); + } + + template + std::ostream& operator<<(std::ostream& out, Distribution const& distribution) { + out << "{"; + for (auto const& entry : distribution) { + out << "[" << entry.second << ": " << entry.first << "], "; + } + out << "}"; + + return out; + } + template class Distribution; + template std::ostream& operator<<(std::ostream& out, Distribution const& distribution); } -} \ No newline at end of file +} diff --git a/src/storage/Distribution.h b/src/storage/Distribution.h index 97639501c..a888a88ed 100644 --- a/src/storage/Distribution.h +++ b/src/storage/Distribution.h @@ -2,6 +2,7 @@ #define STORM_STORAGE_DISTRIBUTION_H_ #include +#include #include #include "src/storage/sparse/StateType.h" @@ -19,7 +20,15 @@ namespace storm { /*! * Creates an empty distribution. */ - Distribution() = default; + Distribution(); + + /*! + * Checks whether the two distributions specify the same probabilities to go to the same states. + * + * @param other The distribution with which the current distribution is to be compared. + * @return True iff the two distributions are equal. + */ + bool operator==(Distribution const& other) const; /*! * Assigns the given state the given probability under this distribution. @@ -57,6 +66,13 @@ namespace storm { */ const_iterator end() const; + /*! + * Retrieves the hash value of the distribution. + * + * @return The hash value of the distribution. + */ + std::size_t getHash() const; + private: // A list of states and the probabilities that are assigned to them. container_type distribution; @@ -65,7 +81,21 @@ namespace storm { std::size_t hash; }; + template + std::ostream& operator<<(std::ostream& out, Distribution const& distribution); } } +namespace std { + + template + struct hash> { + std::size_t operator()(storm::storage::Distribution const& distribution) const { + return (distribution.getHash()); + } + }; + +} + + #endif /* STORM_STORAGE_DISTRIBUTION_H_ */ \ No newline at end of file diff --git a/src/storage/StateBlock.cpp b/src/storage/StateBlock.cpp index fb4f21133..bc270872d 100644 --- a/src/storage/StateBlock.cpp +++ b/src/storage/StateBlock.cpp @@ -30,6 +30,10 @@ namespace storm { states.insert(state); } + StateBlock::iterator StateBlock::insert(container_type::const_iterator iterator, value_type const& state) { + return states.insert(iterator, state); + } + void StateBlock::erase(value_type const& state) { states.erase(state); } @@ -55,5 +59,6 @@ namespace storm { out << block.getStates(); return out; } + } } \ No newline at end of file diff --git a/src/storage/StateBlock.h b/src/storage/StateBlock.h index dcd516277..8497f3a2f 100644 --- a/src/storage/StateBlock.h +++ b/src/storage/StateBlock.h @@ -101,6 +101,13 @@ namespace storm { * @param state The state to add to this SCC. */ void insert(value_type const& state); + + /*! + * Inserts the given element into this SCC. + * + * @param state The state to add to this SCC. + */ + iterator insert(container_type::const_iterator iterator, value_type const& state); /*! * Removes the given element from this SCC. diff --git a/src/utility/cli.h b/src/utility/cli.h index 4ceb64a95..0a416e5c3 100644 --- a/src/utility/cli.h +++ b/src/utility/cli.h @@ -43,6 +43,9 @@ log4cplus::Logger logger; // Headers of adapters. #include "src/adapters/ExplicitModelAdapter.h" +// Headers for model processing. +#include "src/storage/BisimulationDecomposition.h" + // Headers for counterexample generation. #include "src/counterexamples/MILPMinimalLabelSetGenerator.h" #include "src/counterexamples/SMTMinimalCommandSetGenerator.h" @@ -254,6 +257,12 @@ namespace storm { STORM_LOG_THROW(false, storm::exceptions::InvalidSettingsException, "No input model."); } + if (settings.isBisimulationSet()) { + STORM_LOG_THROW(result->getType() == storm::models::DTMC, storm::exceptions::InvalidSettingsException, "Bisimulation minimization is currently only compatible with DTMCs."); + std::shared_ptr> dtmc = result->template as>(); + storm::storage::BisimulationDecomposition bisimulationDecomposition(*dtmc); + } + return result; } @@ -298,6 +307,7 @@ namespace storm { std::shared_ptr> filterFormula = storm::parser::PrctlParser::parsePrctlFormula(settings.getPctlProperty()); generateCounterexample(model, filterFormula->getChild()); } + } }