diff --git a/src/storage/BisimulationDecomposition2.cpp b/src/storage/BisimulationDecomposition2.cpp index 35bd1c703..abda69f66 100644 --- a/src/storage/BisimulationDecomposition2.cpp +++ b/src/storage/BisimulationDecomposition2.cpp @@ -8,23 +8,30 @@ namespace storm { namespace storage { - static int globalId = 0; + static std::size_t globalId = 0; template - BisimulationDecomposition2::Block::Block(storm::storage::sparse::state_type begin, storm::storage::sparse::state_type end, Block* prev, Block* next) : originalBegin(begin), begin(begin), end(end), prev(prev), next(next), numberOfStates(end - begin), isMarked(false), myId(globalId++) { - // Intentionally left empty. - std::cout << "created new block from " << begin << " to " << end << std::endl; + BisimulationDecomposition2::Block::Block(storm::storage::sparse::state_type begin, storm::storage::sparse::state_type end, Block* prev, Block* next) : next(next), prev(prev), begin(begin), end(end), markedAsSplitter(false), markedPosition(begin), id(globalId++) { + if (next != nullptr) { + next->prev = this; + } + if (prev != nullptr) { + prev->next = this; + } } template void BisimulationDecomposition2::Block::print(Partition const& partition) const { - std::cout << "block " << this->myId << " and ptr " << this << std::endl; - std::cout << "begin: " << this->begin << " and end: " << this->end << " (number of states: " << this->numberOfStates << ")" << std::endl; - std::cout << "next: " << this->next << " and prev " << this->prev << std::endl; + std::cout << "block " << this->getId() << " with marked position " << this->getMarkedPosition() << std::endl; + std::cout << "begin: " << this->begin << " and end: " << this->end << " (number of states: " << this->getNumberOfStates() << ")" << std::endl; std::cout << "states:" << std::endl; for (storm::storage::sparse::state_type index = this->begin; index < this->end; ++index) { std::cout << partition.states[index] << " "; } + std::cout << std::endl << "original: " << std::endl; + for (storm::storage::sparse::state_type index = this->getOriginalBegin(); index < this->end; ++index) { + std::cout << partition.states[index] << " "; + } std::cout << std::endl << "values:" << std::endl; for (storm::storage::sparse::state_type index = this->begin; index < this->end; ++index) { std::cout << std::setprecision(3) << partition.values[index] << " "; @@ -35,13 +42,187 @@ namespace storm { template void BisimulationDecomposition2::Block::setBegin(storm::storage::sparse::state_type begin) { this->begin = begin; - this->originalBegin = begin; + this->markedPosition = begin; + } + + template + void BisimulationDecomposition2::Block::setEnd(storm::storage::sparse::state_type end) { + this->end = end; + } + + template + void BisimulationDecomposition2::Block::incrementBegin() { + ++this->begin; + std::cout << "incremented begin to " << this->begin << std::endl; + } + + template + void BisimulationDecomposition2::Block::decrementEnd() { + ++this->begin; + } + + template + storm::storage::sparse::state_type BisimulationDecomposition2::Block::getBegin() const { + return this->begin; + } + + template + storm::storage::sparse::state_type BisimulationDecomposition2::Block::getOriginalBegin() const { + if (this->hasPreviousBlock()) { + return this->getPreviousBlock().getEnd(); + } else { + return 0; + } + } + + template + storm::storage::sparse::state_type BisimulationDecomposition2::Block::getEnd() const { + return this->end; + } + + template + void BisimulationDecomposition2::Block::setIterator(const_iterator it) { + this->selfIt = it; + } + + template + typename BisimulationDecomposition2::Block::const_iterator BisimulationDecomposition2::Block::getIterator() const { + return this->selfIt; + } + + template + typename BisimulationDecomposition2::Block::const_iterator BisimulationDecomposition2::Block::getNextIterator() const { + return this->getNextBlock().getIterator(); + } + + template + typename BisimulationDecomposition2::Block::const_iterator BisimulationDecomposition2::Block::getPreviousIterator() const { + return this->getPreviousBlock().getIterator(); + } + + template + typename BisimulationDecomposition2::Block& BisimulationDecomposition2::Block::getNextBlock() { + return *this->next; + } + + template + typename BisimulationDecomposition2::Block const& BisimulationDecomposition2::Block::getNextBlock() const { + return *this->next; + } + + template + bool BisimulationDecomposition2::Block::hasNextBlock() const { + return this->next != nullptr; + } + + template + typename BisimulationDecomposition2::Block& BisimulationDecomposition2::Block::getPreviousBlock() { + return *this->prev; + } + + template + typename BisimulationDecomposition2::Block* BisimulationDecomposition2::Block::getPreviousBlockPointer() { + return this->prev; + } + + template + typename BisimulationDecomposition2::Block* BisimulationDecomposition2::Block::getNextBlockPointer() { + return this->next; + } + + template + typename BisimulationDecomposition2::Block const& BisimulationDecomposition2::Block::getPreviousBlock() const { + return *this->prev; + } + + template + bool BisimulationDecomposition2::Block::hasPreviousBlock() const { + return this->prev != nullptr; + } + + template + bool BisimulationDecomposition2::Block::check() const { + if (this->begin >= this->end) { + std::cout << "beg: " << this->begin << " and end " << this->end << std::endl; + assert(false); + } + if (this->prev != nullptr) { + assert (this->prev->next == this); + } + if (this->next != nullptr) { + assert (this->next->prev == this); + } + return true; + } + + template + std::size_t BisimulationDecomposition2::Block::getNumberOfStates() const { + // We need to take the original begin here, because the begin is temporarily moved. + return (this->end - this->getOriginalBegin()); + } + + template + bool BisimulationDecomposition2::Block::isMarkedAsSplitter() const { + return this->markedAsSplitter; + } + + template + void BisimulationDecomposition2::Block::markAsSplitter() { + this->markedAsSplitter = true; + } + + template + void BisimulationDecomposition2::Block::unmarkAsSplitter() { + this->markedAsSplitter = false; + } + + template + std::size_t BisimulationDecomposition2::Block::getId() const { + return this->id; + } + + template + storm::storage::sparse::state_type BisimulationDecomposition2::Block::getMarkedPosition() const { + return this->markedPosition; + } + + template + void BisimulationDecomposition2::Block::setMarkedPosition(storm::storage::sparse::state_type position) { + this->markedPosition = position; + } + + template + void BisimulationDecomposition2::Block::resetMarkedPosition() { + this->markedPosition = this->begin; + } + + template + void BisimulationDecomposition2::Block::incrementMarkedPosition() { + ++this->markedPosition; + } + + template + void BisimulationDecomposition2::Block::markAsPredecessorBlock() { + this->markedAsPredecessorBlock = true; + } + + template + void BisimulationDecomposition2::Block::unmarkAsPredecessorBlock() { + this->markedAsPredecessorBlock = false; + } + + template + bool BisimulationDecomposition2::Block::isMarkedAsPredecessor() const { + return markedAsPredecessorBlock; } template BisimulationDecomposition2::Partition::Partition(std::size_t numberOfStates) : stateToBlockMapping(numberOfStates), states(numberOfStates), positions(numberOfStates), values(numberOfStates) { + // Create the block and give it an iterator to itself. typename std::list::iterator it = blocks.emplace(this->blocks.end(), 0, numberOfStates, nullptr, nullptr); - it->itToSelf = it; + it->setIterator(it); + + // Set up the different parts of the internal structure. for (storm::storage::sparse::state_type state = 0; state < numberOfStates; ++state) { states[state] = state; positions[state] = state; @@ -49,40 +230,192 @@ namespace storm { } } + template + void BisimulationDecomposition2::Partition::swapStates(storm::storage::sparse::state_type state1, storm::storage::sparse::state_type state2) { + std::cout << "swapping states " << state1 << " and " << state2 << std::endl; + std::swap(this->states[this->positions[state1]], this->states[this->positions[state2]]); + std::swap(this->values[this->positions[state1]], this->values[this->positions[state2]]); + std::swap(this->positions[state1], this->positions[state2]); + assert(this->check()); + } + + template + void BisimulationDecomposition2::Partition::swapStatesAtPositions(storm::storage::sparse::state_type position1, storm::storage::sparse::state_type position2) { + storm::storage::sparse::state_type state1 = this->states[position1]; + storm::storage::sparse::state_type state2 = this->states[position2]; + std::cout << "swapping states " << state1 << " and " << state2 << " at positions " << position1 << " and " << position2 << std::endl; + + std::swap(this->states[position1], this->states[position2]); + std::swap(this->values[position1], this->values[position2]); + + this->positions[state1] = position2; + this->positions[state2] = position1; + std::cout << "pos of " << state1 << " is now " << position2 << " and pos of " << state2 << " is now " << position1 << std::endl; + std::cout << this->states[position1] << " vs " << state2 << " and " << this->states[position2] << " vs " << state1 << std::endl; + assert(this->check()); + } + + template + storm::storage::sparse::state_type const& BisimulationDecomposition2::Partition::getPosition(storm::storage::sparse::state_type state) const { + return this->positions[state]; + } + + template + void BisimulationDecomposition2::Partition::setPosition(storm::storage::sparse::state_type state, storm::storage::sparse::state_type position) { + this->positions[state] = position; + } + + template + storm::storage::sparse::state_type const& BisimulationDecomposition2::Partition::getState(storm::storage::sparse::state_type position) const { + return this->states[position]; + } + + template + ValueType const& BisimulationDecomposition2::Partition::getValue(storm::storage::sparse::state_type state) const { + return this->values[this->positions[state]]; + } + + template + ValueType const& BisimulationDecomposition2::Partition::getValueAtPosition(storm::storage::sparse::state_type position) const { + return this->values[position]; + } + + template + void BisimulationDecomposition2::Partition::setValue(storm::storage::sparse::state_type state, ValueType value) { + this->values[this->positions[state]] = value; + } + + template + std::vector& BisimulationDecomposition2::Partition::getValues() { + return this->values; + } + + template + void BisimulationDecomposition2::Partition::increaseValue(storm::storage::sparse::state_type state, ValueType value) { + this->values[this->positions[state]] += value; + } + + template + void BisimulationDecomposition2::Partition::updateBlockMapping(Block& block, std::vector::iterator first, std::vector::iterator last) { + for (; first != last; ++first) { + this->stateToBlockMapping[*first] = █ + } + } + + template + typename BisimulationDecomposition2::Block& BisimulationDecomposition2::Partition::getBlock(storm::storage::sparse::state_type state) { + return *this->stateToBlockMapping[state]; + } + + template + std::vector::iterator BisimulationDecomposition2::Partition::getBeginOfStates(Block const& block) { + return this->states.begin() + block.getBegin(); + } + + template + std::vector::iterator BisimulationDecomposition2::Partition::getEndOfStates(Block const& block) { + return this->states.begin() + block.getEnd(); + } + + template + std::vector::const_iterator BisimulationDecomposition2::Partition::getBeginOfStates(Block const& block) const { + return this->states.begin() + block.getBegin(); + } + + template + std::vector::const_iterator BisimulationDecomposition2::Partition::getEndOfStates(Block const& block) const { + return this->states.begin() + block.getEnd(); + } + + template + typename std::vector::iterator BisimulationDecomposition2::Partition::getBeginOfValues(Block const& block) { + return this->values.begin() + block.getBegin(); + } + + template + typename std::vector::iterator BisimulationDecomposition2::Partition::getEndOfValues(Block const& block) { + return this->values.begin() + block.getEnd(); + } + + template + typename BisimulationDecomposition2::Block& BisimulationDecomposition2::Partition::splitBlock(Block& block, storm::storage::sparse::state_type position) { + // FIXME: this could be improved by splitting off the smaller of the two resulting blocks. + + std::cout << "splitting block (" << block.getBegin() << "," << block.getEnd() << ") at position " << position << std::endl; + block.print(*this); + + // In case one of the resulting blocks would be empty, we simply return the current block and do not create + // a new one. + if (position == block.getBegin() || position == block.getEnd()) { + return block; + } + + // Actually create the new block and insert it at the correct position. + typename std::list::iterator selfIt = this->blocks.emplace(block.getIterator(), block.getBegin(), position, block.getPreviousBlockPointer(), &block); + std::cout << "created new block from " << block.getBegin() << " to " << position << std::endl; + selfIt->setIterator(selfIt); + Block& newBlock = *selfIt; + + // Make the current block end at the given position. + block.setBegin(position); + + std::cout << "old block: " << std::endl; + block.print(*this); + + // Update the mapping of the states in the newly created block. + for (auto it = this->getBeginOfStates(newBlock), ite = this->getEndOfStates(newBlock); it != ite; ++it) { + stateToBlockMapping[*it] = &newBlock; + } + + return newBlock; + } + + template + typename BisimulationDecomposition2::Block& BisimulationDecomposition2::Partition::insertBlock(Block& block) { + // Find the beginning of the new block. + storm::storage::sparse::state_type begin; + if (block.hasPreviousBlock()) { + begin = block.getPreviousBlock().getEnd(); + std::cout << "previous block ended at " << begin << std::endl; + block.getPreviousBlock().print(*this); + } else { + begin = 0; + } + + // Actually insert the new block. + typename std::list::iterator it = this->blocks.emplace(block.getIterator(), begin, block.getBegin(), block.getPreviousBlockPointer(), &block); + Block& newBlock = *it; + newBlock.setIterator(it); + + // Update the mapping of the states in the newly created block. + for (auto it = this->getBeginOfStates(newBlock), ite = this->getEndOfStates(newBlock); it != ite; ++it) { + stateToBlockMapping[*it] = &newBlock; + } + + return *it; + } + template void BisimulationDecomposition2::Partition::splitLabel(storm::storage::BitVector const& statesWithLabel) { for (auto blockIterator = this->blocks.begin(), ite = this->blocks.end(); blockIterator != ite; ) { // The update of the loop was intentionally moved to the bottom of the loop. Block& block = *blockIterator; // Sort the range of the block such that all states that have the label are moved to the front. - std::sort(this->states.begin() + block.begin, this->states.begin() + block.end, [&statesWithLabel] (storm::storage::sparse::state_type const& a, storm::storage::sparse::state_type const& b) { return statesWithLabel.get(a) && !statesWithLabel.get(b); } ); + std::sort(this->getBeginOfStates(block), this->getEndOfStates(block), [&statesWithLabel] (storm::storage::sparse::state_type const& a, storm::storage::sparse::state_type const& b) { return statesWithLabel.get(a) && !statesWithLabel.get(b); } ); // Update the positions vector. - storm::storage::sparse::state_type position = block.begin; - for (auto stateIt = this->states.begin() + block.begin, stateIte = this->states.begin() + block.end; stateIt != stateIte; ++stateIt, ++position) { + storm::storage::sparse::state_type position = block.getBegin(); + for (auto stateIt = this->getBeginOfStates(block), stateIte = this->getEndOfStates(block); stateIt != stateIte; ++stateIt, ++position) { this->positions[*stateIt] = position; } // Now we can find the first position in the block that does not have the label and create new blocks. - std::vector::iterator it = std::find_if(this->states.begin() + block.begin, this->states.begin() + block.end, [&] (storm::storage::sparse::state_type const& a) { return !statesWithLabel.get(a); }); + std::vector::iterator it = std::find_if(this->getBeginOfStates(block), this->getEndOfStates(block), [&] (storm::storage::sparse::state_type const& a) { return !statesWithLabel.get(a); }); // If not all the states agreed on the validity of the label, we need to split the block. - if (it != this->states.begin() + block.begin && it != this->states.begin() + block.end) { + if (it != this->getBeginOfStates(block) && it != this->getEndOfStates(block)) { auto cutPoint = std::distance(this->states.begin(), it); - - ++blockIterator; - auto newBlockIterator = this->blocks.emplace(blockIterator, cutPoint, block.end, &block, block.next); - newBlockIterator->itToSelf = newBlockIterator; - - // Make the old block end at the cut position and insert a new block after it. - block.end = cutPoint; - block.next = &(*newBlockIterator); - block.numberOfStates = block.end - block.begin; - - // Update the block mapping for all states that we just removed from the block. - for (auto it = this->states.begin() + newBlockIterator->begin, ite = this->states.begin() + newBlockIterator->end; it != ite; ++it) { - stateToBlockMapping[*it] = &(*newBlockIterator); - } + this->splitBlock(block, cutPoint); } else { // Otherwise, we simply proceed to the next block. ++blockIterator; @@ -90,6 +423,40 @@ namespace storm { } } + template + std::list::Block> const& BisimulationDecomposition2::Partition::getBlocks() const { + return this->blocks; + } + + template + std::vector& BisimulationDecomposition2::Partition::getStates() { + return this->states; + } + + template + std::list::Block>& BisimulationDecomposition2::Partition::getBlocks() { + return this->blocks; + } + + template + bool BisimulationDecomposition2::Partition::check() const { + for (uint_fast64_t state = 0; state < this->positions.size(); ++state) { + if (this->states[this->positions[state]] != state) { + assert(false); + } + } + for (auto const& block : this->blocks) { + assert(block.check()); + for (auto stateIt = this->getBeginOfStates(block), stateIte = this->getEndOfStates(block); stateIt != stateIte; ++stateIt) { + if (this->stateToBlockMapping[*stateIt] != &block) { + std::cout << "state " << *stateIt << " has wrong block mapping " << this->stateToBlockMapping[*stateIt] << " (should have " << &block << ")" << std::endl; + assert(false); + } + } + } + return true; + } + template void BisimulationDecomposition2::Partition::print() const { for (auto const& block : this->blocks) { @@ -109,6 +476,7 @@ namespace storm { } std::cout << std::endl; std::cout << "size: " << this->blocks.size() << std::endl; + assert(this->check()); } template @@ -129,7 +497,7 @@ namespace storm { template void BisimulationDecomposition2::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. Partition partition(dtmc.getNumberOfStates()); partition.print(); @@ -142,10 +510,11 @@ namespace storm { std::cout << "initial partition:" << std::endl; partition.print(); + assert(partition.check()); // Initially, all blocks are potential splitter, so we insert them in the splitterQueue. std::deque splitterQueue; - std::for_each(partition.blocks.begin(), partition.blocks.end(), [&] (Block& a) { splitterQueue.push_back(&a); }); + std::for_each(partition.getBlocks().begin(), partition.getBlocks().end(), [&] (Block& a) { splitterQueue.push_back(&a); }); storm::storage::SparseMatrix backwardTransitions = dtmc.getBackwardTransitions(); @@ -166,318 +535,213 @@ namespace storm { } template - std::size_t BisimulationDecomposition2::splitBlockProbabilities(Block* block, Partition& partition, std::deque& splitterQueue) { + std::size_t BisimulationDecomposition2::splitBlockProbabilities(Block& block, Partition& partition, std::deque& splitterQueue) { std::cout << "part before split prob" << std::endl; partition.print(); - Block& currentBlock = *block; - // Sort the states in the block based on their probabilities. - std::sort(partition.states.begin() + currentBlock.begin, partition.states.begin() + currentBlock.end, [&partition] (storm::storage::sparse::state_type const& a, storm::storage::sparse::state_type const& b) { return partition.values[a] < partition.values[b]; } ); + std::sort(partition.getBeginOfStates(block), partition.getEndOfStates(block), [&partition] (storm::storage::sparse::state_type const& a, storm::storage::sparse::state_type const& b) { return partition.getValue(a) < partition.getValue(b); } ); // FIXME: This can probably be done more efficiently. - std::sort(partition.values.begin() + currentBlock.begin, partition.values.begin() + currentBlock.end); + std::sort(partition.getBeginOfValues(block), partition.getEndOfValues(block)); // Update the positions vector. - storm::storage::sparse::state_type position = currentBlock.begin; - for (auto stateIt = partition.states.begin() + currentBlock.begin, stateIte = partition.states.begin() + currentBlock.end; stateIt != stateIte; ++stateIt, ++position) { - partition.positions[*stateIt] = position; + storm::storage::sparse::state_type position = block.getBegin(); + for (auto stateIt = partition.getBeginOfStates(block), stateIte = partition.getEndOfStates(block); stateIt != stateIte; ++stateIt, ++position) { + partition.setPosition(*stateIt, position); } - + // Finally, we need to scan the ranges of states that agree on the probability. - storm::storage::sparse::state_type beginIndex = currentBlock.begin; + storm::storage::sparse::state_type beginIndex = block.getBegin(); storm::storage::sparse::state_type currentIndex = beginIndex; - storm::storage::sparse::state_type endIndex = currentBlock.end; + storm::storage::sparse::state_type endIndex = block.getEnd(); // 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. - // Note that this check requires the block to be non-empty. - if (std::abs(partition.values[beginIndex] - partition.values[endIndex - 1]) < 1e-6) { - return 1; - } - - // 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 currentValue = partition.values[beginIndex]; - ValueType* valueIterator = &(partition.values[beginIndex]) + 1; - ++currentIndex; - while (std::abs(*valueIterator - currentValue) < 1e-6) { - std::cout << "prob: " << *valueIterator << std::endl; - ++currentIndex; - ++valueIterator; - } - - // Now resize the block. - std::cout << "resizing block from " << std::endl; - currentBlock.print(partition); - currentBlock.end = currentIndex; - currentBlock.numberOfStates = currentBlock.end - currentBlock.begin; - std::cout << " to " << std::endl; - currentBlock.print(partition); - beginIndex = currentIndex; - - // If it is not already a splitter, we mark it as such. - if (!currentBlock.isMarked) { - currentBlock.isMarked = true; - } - - Block* prevBlock = ¤tBlock; - - // Now scan for new blocks. - std::list createdBlocks; - std::cout << currentIndex << " < " << endIndex << std::endl; - typename std::list::const_iterator insertPosition = currentBlock.itToSelf; - ++insertPosition; - while (currentIndex < endIndex) { - ValueType currentValue = *(partition.values.begin() + currentIndex); + std::size_t createdBlocks = 0; + while (std::abs(partition.getValueAtPosition(beginIndex) - partition.getValueAtPosition(endIndex)) >= 1e-6) { + // 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 = partition.getValueAtPosition(beginIndex); + typename std::vector::const_iterator valueIterator = partition.getValues().begin() + beginIndex + 1; ++currentIndex; - ValueType* nextValuePtr = ¤tValue; - while (currentIndex < endIndex && std::abs(currentValue - *nextValuePtr) < 1e-6) { + while (currentIndex < endIndex && std::abs(*valueIterator - currentValue) <= 1e-6) { + ++valueIterator; ++currentIndex; - ++nextValuePtr; } - // Create a new block from the states that agree on the values. - ++insertPosition; - typename std::list::iterator newBlockIterator = partition.blocks.emplace(insertPosition, beginIndex, currentIndex, prevBlock, prevBlock->next); - newBlockIterator->itToSelf = newBlockIterator; - Block* newBlock = &(*newBlockIterator); - prevBlock->next = newBlock; - prevBlock = newBlock; - if (prevBlock->numberOfStates > 1) { - createdBlocks.emplace_back(prevBlock); + // Now we split the block and mark it as a potential splitter. + ++createdBlocks; + Block& newBlock = partition.splitBlock(block, currentIndex); + if (!newBlock.isMarkedAsPredecessor()) { + newBlock.markAsSplitter(); + splitterQueue.push_back(&newBlock); } beginIndex = currentIndex; - - // Update the block mapping for the moved states. - for (auto it = partition.states.begin() + newBlock->begin, ite = partition.states.begin() + newBlock->end; it != ite; ++it) { - partition.stateToBlockMapping[*it] = newBlock; - } - } - - for (auto block : createdBlocks) { - splitterQueue.push_back(block); - block->isMarked = true; } - std::cout << "created " << createdBlocks.size() << " blocks" << std::endl; + assert(partition.check()); - return createdBlocks.size(); + return createdBlocks; } template - std::size_t BisimulationDecomposition2::splitPartition(storm::storage::SparseMatrix const& backwardTransitions, Block const& splitter, Partition& partition, std::deque& splitterQueue) { - std::cout << "getting block " << splitter.myId << " as splitter" << std::endl; + std::size_t BisimulationDecomposition2::splitPartition(storm::storage::SparseMatrix const& backwardTransitions, Block& splitter, Partition& partition, std::deque& splitterQueue) { + std::cout << "treating block " << splitter.getId() << " as splitter" << std::endl; + splitter.print(partition); + std::list predecessorBlocks; // Iterate over all states of the splitter and check its predecessors. - for (auto stateIterator = partition.states.begin() + splitter.begin, stateIte = partition.states.begin() + splitter.end; stateIterator != stateIte; ++stateIterator) { - storm::storage::sparse::state_type& state = *stateIterator; + storm::storage::sparse::state_type currentPosition = splitter.getBegin(); + for (auto stateIterator = partition.getBeginOfStates(splitter), stateIte = partition.getEndOfStates(splitter); stateIterator != stateIte; ++stateIterator, ++currentPosition) { + std::cout << "states -----" << std::endl; + for (auto stateIterator = partition.getStates().begin() + splitter.getOriginalBegin(), stateIte = partition.getStates().begin() + splitter.getEnd(); stateIterator != stateIte; ++stateIterator) { + std::cout << *stateIterator << " "; + } + std::cout << std::endl; + + storm::storage::sparse::state_type currentState = *stateIterator; + std::cout << "current state " << currentState << " at pos " << currentPosition << std::endl; - for (auto const& predecessorEntry : backwardTransitions.getRow(state)) { + uint_fast64_t elementsToSkip = 0; + for (auto const& predecessorEntry : backwardTransitions.getRow(currentState)) { storm::storage::sparse::state_type predecessor = predecessorEntry.getColumn(); - std::cout << "found predecessor " << predecessor << " of state " << *stateIterator << std::endl; - Block* predecessorBlock = partition.stateToBlockMapping[predecessor]; - std::cout << "predecessor block " << predecessorBlock->myId << std::endl; + std::cout << "found predecessor " << predecessor << " of state " << currentState << std::endl; + + Block& predecessorBlock = partition.getBlock(predecessor); + std::cout << "predecessor block " << predecessorBlock.getId() << std::endl; // If the predecessor block has just one state, there is no point in splitting it. - if (predecessorBlock->numberOfStates <= 1) { + if (predecessorBlock.getNumberOfStates() <= 1) { continue; } - storm::storage::sparse::state_type predecessorPosition = partition.positions[predecessor]; + storm::storage::sparse::state_type predecessorPosition = partition.getPosition(predecessor); // If we have not seen this predecessor before, we move it to a part near the beginning of the block. - std::cout << "predecessor position: " << predecessorPosition << " and begin " << predecessorBlock->begin << std::endl; - if (predecessorPosition >= predecessorBlock->begin) { - -// We cannot directly move the states at this point, otherwise we would consider states multiple -// times while others are skipped. -// std::swap(partition.states[predecessorPosition], partition.states[predecessorBlock->begin]); -// std::cout << "swapping positions of " << predecessor << " and " << partition.states[predecessorPosition] << std::endl; - - // Instead, we only virtually move the states by setting their positions and move them in place - // later. - storm::storage::sparse::state_type tmp = partition.positions[predecessor]; - partition.positions[predecessor] = predecessorBlock->begin; - std::cout << "moved " << predecessor << " to pos " << predecessorBlock->begin << std::endl; - partition.positions[partition.states[predecessorBlock->begin]] = tmp; - -// storm::storage::sparse::state_type tmp = partition.positions[partition.states[predecessorPosition]]; -// partition.positions[partition.states[predecessorPosition]] = partition.positions[predecessor]; -// partition.positions[predecessor] = tmp; - -// std::swap(partition.positions[predecessor], partition.positions[predecessorBlock->begin]); - - ++predecessorBlock->begin; - std::cout << "incrementing begin, setting probability at " << partition.positions[predecessor] << " to " << predecessorEntry.getValue() << " ... " << std::endl; - partition.values[partition.positions[predecessor]] = predecessorEntry.getValue(); - assert(partition.values[partition.positions[predecessor]] <= 1); + if (predecessorPosition >= predecessorBlock.getBegin()) { + if (&predecessorBlock == &splitter) { + // If the predecessor we just found was already processed (in terms of visiting its predecessors), + // we swap it with the state that is currently at the beginning of the block and move the start + // of the block one step further. + if (predecessorPosition <= currentPosition) { + partition.swapStates(predecessor, partition.getState(predecessorBlock.getBegin())); + predecessorBlock.incrementBegin(); + } else { + std::cout << "current position is " << currentPosition << std::endl; + // Otherwise, we need to move the predecessor, but we need to make sure that we explore its + // predecessors later. + if (predecessorBlock.getMarkedPosition() == predecessorBlock.getBegin()) { + partition.swapStatesAtPositions(predecessorBlock.getMarkedPosition(), predecessorPosition); + partition.swapStatesAtPositions(predecessorPosition, currentPosition + elementsToSkip + 1); + } else { + partition.swapStatesAtPositions(predecessorBlock.getMarkedPosition(), predecessorPosition); + partition.swapStatesAtPositions(predecessorPosition, predecessorBlock.getBegin()); + partition.swapStatesAtPositions(predecessorPosition, currentPosition + elementsToSkip + 1); + } + + ++elementsToSkip; + predecessorBlock.incrementMarkedPosition(); + predecessorBlock.incrementBegin(); + } + } else { + partition.swapStates(predecessor, partition.getState(predecessorBlock.getBegin())); + predecessorBlock.incrementBegin(); + } + partition.setValue(predecessor, predecessorEntry.getValue()); } else { // Otherwise, we just need to update the probability for this predecessor. - partition.values[predecessorPosition] += predecessorEntry.getValue(); - std::cout << "updating probability " << predecessorPosition << " to " << std::setprecision(10) << partition.values[predecessorPosition] << std::endl; - assert(partition.values[predecessorPosition] <= 1 + 1e-6); + partition.increaseValue(predecessor, predecessorEntry.getValue()); } - if (!predecessorBlock->isMarked) { - predecessorBlocks.emplace_back(predecessorBlock); - predecessorBlock->isMarked = true; + if (!predecessorBlock.isMarkedAsPredecessor()) { + predecessorBlocks.emplace_back(&predecessorBlock); + predecessorBlock.markAsPredecessorBlock(); } } + + // If we had to move some elements beyond the current element, we may have to skip them. + if (elementsToSkip > 0) { + std::cout << "skipping " << elementsToSkip << " elements" << std::endl; + stateIterator += elementsToSkip; + currentPosition += elementsToSkip; + } } - - // Now that we have computed the new positions of the states we want to move, we need to actually move them. - for (auto block : predecessorBlocks) { - std::cout << "moving block " << block->myId << std::endl; - for (auto stateIterator = partition.states.begin() + block->originalBegin, stateIte = partition.states.begin() + block->end; stateIterator != stateIte; ++stateIterator) { - std::cout << "swapping " << *stateIterator << " to " << partition.positions[*stateIterator] << std::endl; - std::swap(partition.states[partition.positions[*stateIterator]], *stateIterator); + + // Now we can traverse the list of states of the splitter whose predecessors we have not yet explored. + for (auto stateIterator = partition.getStates().begin() + splitter.getOriginalBegin(), stateIte = partition.getStates().begin() + splitter.getMarkedPosition(); stateIterator != stateIte; ++stateIterator) { + storm::storage::sparse::state_type currentState = *stateIterator; + + for (auto const& predecessorEntry : backwardTransitions.getRow(currentState)) { + storm::storage::sparse::state_type predecessor = predecessorEntry.getColumn(); + Block& predecessorBlock = partition.getBlock(predecessor); + storm::storage::sparse::state_type predecessorPosition = partition.getPosition(predecessor); + + if (predecessorPosition >= predecessorBlock.getBegin()) { + partition.swapStatesAtPositions(predecessorPosition, predecessorBlock.getBegin()); + predecessorBlock.incrementBegin(); + partition.setValue(predecessor, predecessorEntry.getValue()); + } else { + partition.increaseValue(predecessor, predecessorEntry.getValue()); + } + + if (!predecessorBlock.isMarkedAsPredecessor()) { + predecessorBlocks.emplace_back(&predecessorBlock); + predecessorBlock.markAsPredecessorBlock(); + } } - // Set the original begin and begin to the same value. - block->setBegin(block->begin); } + splitter.setMarkedPosition(splitter.getBegin()); + std::list blocksToSplit; + std::cout << "having " << predecessorBlocks.size() << " pred blocks " << std::endl; + // 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 block : predecessorBlocks) { - block->isMarked = false; + 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->begin != block->end) { - std::cout << "moved begin to " << block->begin << " and end to " << block->end << std::endl; - storm::storage::sparse::state_type tmpBegin = block->begin; - storm::storage::sparse::state_type tmpEnd = block->end; - - block->setBegin(block->prev != nullptr ? block->prev->end : 0); - block->end = tmpBegin; - block->numberOfStates = block->end - block->begin; - - // Create a new block that holds all states that do not have a successor in the current splitter. - typename std::list::iterator it = partition.blocks.emplace(block->next != nullptr ? block->next->itToSelf : partition.blocks.end(), tmpBegin, tmpEnd, block, block->next); - Block* newBlock = &(*it); - newBlock->itToSelf = it; - if (block->next != nullptr) { - block->next->prev = newBlock; - } - block->next = newBlock; + if (block.getBegin() != block.getEnd()) { + std::cout << "moved begin of block " << block.getId() << " to " << block.getBegin() << " and end to " << block.getEnd() << std::endl; + Block& newBlock = partition.insertBlock(block); std::cout << "created new block " << std::endl; - newBlock->print(partition); - - // Update the block mapping in the partition. - for (auto it = partition.states.begin() + newBlock->begin, ite = partition.states.begin() + newBlock->end; it != ite; ++it) { - partition.stateToBlockMapping[*it] = newBlock; - } - - // Mark the half of the block that can be further refined using the probability information. - blocksToSplit.emplace_back(block); - block->print(partition); - - splitterQueue.push_back(newBlock); + newBlock.print(partition); + splitterQueue.push_back(&newBlock); } else { - std::cout << "found block to split" << std::endl; - + std::cout << "all states are predecessors" << std::endl; + // In this case, we can keep the block by setting its begin to the old value. - block->setBegin((block->prev != nullptr) ? block->prev->end : 0); - blocksToSplit.emplace_back(block); + block.setBegin(block.getOriginalBegin()); + blocksToSplit.emplace_back(&block); } - - block->setBegin(block->begin); } + assert(partition.check()); + // Finally, we walk through the blocks that have a transition to the splitter and split them using // probabilistic information. - for (auto block : blocksToSplit) { - if (block->numberOfStates <= 1) { + for (auto blockPtr : blocksToSplit) { + std::cout << "block to split: " << blockPtr->getId() << std::endl; + if (blockPtr->getNumberOfStates() <= 1) { continue; } - splitBlockProbabilities(block, partition, splitterQueue); + splitBlockProbabilities(*blockPtr, partition, splitterQueue); } return 0; } -// template -// std::size_t BisimulationDecomposition2::splitPartition(storm::models::Dtmc const& dtmc, storm::storage::SparseMatrix const& backwardTransitions, std::size_t const& blockId, std::vector& stateToBlockMapping, storm::storage::BitVector& blocksInSplitterQueue, std::deque& splitterQueue, bool weakBisimulation) { -// std::chrono::high_resolution_clock::time_point totalStart = std::chrono::high_resolution_clock::now(); -// std::unordered_map, typename BisimulationDecomposition2::block_type> distributionToNewBlocks; -// -// // Traverse all states of the block and check whether they have different distributions. -// std::chrono::high_resolution_clock::time_point gatherStart = std::chrono::high_resolution_clock::now(); -// 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 we are requested to compute a weak bisimulation, we need to scale the distribution with the -// // self-loop probability. -// if (weakBisimulation) { -// distribution.scale(blockId); -// } -// -// // 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); -// } -// } -// -// std::chrono::high_resolution_clock::duration gatherTime = std::chrono::high_resolution_clock::now() - gatherStart; -// std::cout << "time to iterate over all states was " << std::chrono::duration_cast(gatherTime).count() << "ms." << std::endl; -// -// // 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 BisimulationDecomposition2::block_type tmpBlock; -// -// auto distributionIterator = distributionToNewBlocks.begin(); -// 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)); -// -// // 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; -// } -// } -// -// // 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 BisimulationDecomposition2; } } \ No newline at end of file diff --git a/src/storage/BisimulationDecomposition2.h b/src/storage/BisimulationDecomposition2.h index 2cc4a98d8..50153d4b5 100644 --- a/src/storage/BisimulationDecomposition2.h +++ b/src/storage/BisimulationDecomposition2.h @@ -30,6 +30,8 @@ namespace storm { class Block { public: + typedef typename std::list::const_iterator const_iterator; + Block(storm::storage::sparse::state_type begin, storm::storage::sparse::state_type end, Block* prev, Block* next); // Prints the block. @@ -37,31 +39,130 @@ namespace storm { // Sets the beginning index of the block. void setBegin(storm::storage::sparse::state_type begin); + + // Moves the beginning index of the block one step further. + void incrementBegin(); + + // Sets the end index of the block. + void setEnd(storm::storage::sparse::state_type end); + + // Moves the end index of the block one step to the front. + void decrementEnd(); + + // Returns the beginning index of the block. + storm::storage::sparse::state_type getBegin() const; + + // Returns the beginning index of the block. + storm::storage::sparse::state_type getEnd() const; + + // Retrieves the original beginning index of the block in case the begin index has been moved. + storm::storage::sparse::state_type getOriginalBegin() const; + + // Returns the iterator the block in the list of overall blocks. + const_iterator getIterator() const; + + // Returns the iterator the block in the list of overall blocks. + void setIterator(const_iterator it); + + // Returns the iterator the next block in the list of overall blocks if it exists. + const_iterator getNextIterator() const; + + // Returns the iterator the next block in the list of overall blocks if it exists. + const_iterator getPreviousIterator() const; + + // Gets the next block (if there is one). + Block& getNextBlock(); + + // Gets the next block (if there is one). + Block const& getNextBlock() const; + + // Gets a pointer to the next block (if there is one). + Block* getNextBlockPointer(); + + // Retrieves whether the block as a successor block. + bool hasNextBlock() const; + + // Gets the previous block (if there is one). + Block& getPreviousBlock(); + + // Gets a pointer to the previous block (if there is one). + Block* getPreviousBlockPointer(); + + // Gets the next block (if there is one). + Block const& getPreviousBlock() const; + + // Retrieves whether the block as a successor block. + bool hasPreviousBlock() const; + + // Checks consistency of the information in the block. + bool check() const; + + // Retrieves the number of states in this block. + std::size_t getNumberOfStates() const; + + // Checks whether the block is marked as a splitter. + bool isMarkedAsSplitter() const; + + // Marks the block as being a splitter. + void markAsSplitter(); + + // Removes the mark. + void unmarkAsSplitter(); + + // Retrieves the ID of the block. + std::size_t getId() const; + + // Retrieves the marked position in the block. + storm::storage::sparse::state_type getMarkedPosition() const; + + // Sets the marked position to the given value.. + void setMarkedPosition(storm::storage::sparse::state_type position); + + // Increases the marked position by one. + void incrementMarkedPosition(); + + // Resets the marked position to the begin of the block. + void resetMarkedPosition(); + + // Retrieves whether the block is marked as a predecessor. + bool isMarkedAsPredecessor() const; + + // Marks the block as being a predecessor block. + void markAsPredecessorBlock(); + + // Removes the marking. + void unmarkAsPredecessorBlock(); + private: // An iterator to itself. This is needed to conveniently insert elements in the overall list of blocks // kept by the partition. - typename std::list::const_iterator itToSelf; + const_iterator selfIt; + + // Pointers to the next and previous block. + Block* next; + Block* prev; // The begin and end indices of the block in terms of the state vector of the partition. - storm::storage::sparse::state_type originalBegin; storm::storage::sparse::state_type begin; storm::storage::sparse::state_type end; - // The block before and after the current one. - Block* prev; - Block* next; + // A field that can be used for marking the block. + bool markedAsSplitter; - // The number of states in the block. - std::size_t numberOfStates; + // A field that can be used for marking the block as a predecessor block. + bool markedAsPredecessorBlock; - // A field that can be used for marking the block. - bool isMarked; + // A position that can be used to store a certain position within the block. + storm::storage::sparse::state_type markedPosition; - int myId; + // The ID of the block. This is only used for debugging purposes. + std::size_t id; }; class Partition { public: + friend class Block; + /*! * Creates a partition with one block consisting of all the states. */ @@ -73,6 +174,87 @@ namespace storm { */ void splitLabel(storm::storage::BitVector const& statesWithLabel); + // Retrieves the size of the partition, i.e. the number of blocks. + std::size_t size() const; + + // Prints the partition to the standard output. + void print() const; + + // Splits the block at the given position and inserts a new block after the current one holding the rest + // of the states. + Block& splitBlock(Block& block, storm::storage::sparse::state_type position); + + // Inserts a block before the given block. The new block will cover all states between the beginning + // of the given block and the end of the previous block. + Block& insertBlock(Block& block); + + // Retrieves the blocks of the partition. + std::list const& getBlocks() const; + + // Retrieves the blocks of the partition. + std::list& getBlocks(); + + // Retrieves the vector of all the states. + std::vector& getStates(); + + // Checks the partition for internal consistency. + bool check() const; + + // Returns an iterator to the beginning of the states of the given block. + std::vector::iterator getBeginOfStates(Block const& block); + + // Returns an iterator to the beginning of the states of the given block. + std::vector::iterator getEndOfStates(Block const& block); + + // Returns an iterator to the beginning of the states of the given block. + std::vector::const_iterator getBeginOfStates(Block const& block) const; + + // Returns an iterator to the beginning of the states of the given block. + std::vector::const_iterator getEndOfStates(Block const& block) const; + + // Returns an iterator to the beginning of the states of the given block. + typename std::vector::iterator getBeginOfValues(Block const& block); + + // Returns an iterator to the beginning of the states of the given block. + typename std::vector::iterator getEndOfValues(Block const& block); + + // Swaps the positions of the two given states. + void swapStates(storm::storage::sparse::state_type state1, storm::storage::sparse::state_type state2); + + // Swaps the positions of the two states given by their positions. + void swapStatesAtPositions(storm::storage::sparse::state_type position1, storm::storage::sparse::state_type position2); + + // Retrieves the block of the given state. + Block& getBlock(storm::storage::sparse::state_type state); + + // Retrieves the position of the given state. + storm::storage::sparse::state_type const& getPosition(storm::storage::sparse::state_type state) const; + + // Retrieves the position of the given state. + void setPosition(storm::storage::sparse::state_type state, storm::storage::sparse::state_type position); + + // Sets the position of the state to the given position. + storm::storage::sparse::state_type const& getState(storm::storage::sparse::state_type position) const; + + // Retrieves the value for the given state. + ValueType const& getValue(storm::storage::sparse::state_type state) const; + + // Retrieves the value at the given position. + ValueType const& getValueAtPosition(storm::storage::sparse::state_type position) const; + + // Sets the given value for the given state. + void setValue(storm::storage::sparse::state_type state, ValueType value); + + // Retrieves the vector with the probabilities going into the current splitter. + std::vector& getValues(); + + // Increases the value for the given state by the specified amount. + 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, std::vector::iterator first, std::vector::iterator end); + + private: // The list of blocks in the partition. std::list blocks; @@ -88,17 +270,13 @@ namespace storm { // This vector stores the probabilities of going to the current splitter. std::vector values; - - std::size_t size() const; - - void print() const; }; void computeBisimulationEquivalenceClasses(storm::models::Dtmc const& model, bool weak); - std::size_t splitPartition(storm::storage::SparseMatrix const& backwardTransitions, Block const& splitter, Partition& partition, std::deque& splitterQueue); + std::size_t splitPartition(storm::storage::SparseMatrix const& backwardTransitions, Block& splitter, Partition& partition, std::deque& splitterQueue); - std::size_t splitBlockProbabilities(Block* block, Partition& partition, std::deque& splitterQueue); + std::size_t splitBlockProbabilities(Block& block, Partition& partition, std::deque& splitterQueue); }; } }