From b6be0f335619910f625a9b616da52aa79e100412 Mon Sep 17 00:00:00 2001 From: Sebastian Junges Date: Wed, 16 Sep 2020 08:50:19 -0700 Subject: [PATCH] convex reduction --- .../NondeterministicBeliefTracker.cpp | 63 +++++----- .../generator/NondeterministicBeliefTracker.h | 3 +- .../transformer/ObservationTraceUnfolder.cpp | 13 +- .../storage/geometry/ReduceVertexCloud.cpp | 119 ++++++++++++++++++ .../storage/geometry/ReduceVertexCloud.h | 29 +++++ 5 files changed, 185 insertions(+), 42 deletions(-) create mode 100644 src/storm/storage/geometry/ReduceVertexCloud.cpp create mode 100644 src/storm/storage/geometry/ReduceVertexCloud.h diff --git a/src/storm-pomdp/generator/NondeterministicBeliefTracker.cpp b/src/storm-pomdp/generator/NondeterministicBeliefTracker.cpp index e330b7f90..37147f4ce 100644 --- a/src/storm-pomdp/generator/NondeterministicBeliefTracker.cpp +++ b/src/storm-pomdp/generator/NondeterministicBeliefTracker.cpp @@ -2,6 +2,7 @@ #include "storm-pomdp/generator/NondeterministicBeliefTracker.h" #include "storm/utility/ConstantsComparator.h" #include "storm/storage/geometry/nativepolytopeconversion/QuickHull.h" +#include "storm/storage/geometry/ReduceVertexCloud.h" #include "storm/utility/vector.h" namespace storm { @@ -198,6 +199,12 @@ namespace storm { return manager->getNumberOfStates(); } + template + std::map const& SparseBeliefState::getBeliefMap() const { + return belief; + } + + template void SparseBeliefState::setSupport(storm::storage::BitVector& support) const { assert(false); @@ -486,45 +493,33 @@ namespace storm { } // -// template -// void NondeterministicBeliefTracker::reduce() { -// std::cout << "reduce" << std::endl; -// storm::storage::BitVector support(beliefs.begin()->getSupportSize()); -// std::cout << "supp generated" << std::endl; -// for(auto const& belief : beliefs) { -// belief.setSupport(support); -// } -// std::cout << "Support:" << support << std::endl; -// -// if (beliefs.size() <= support.getNumberOfSetBits()) { -// return; -// } -// -// std::vector::EigenVector> points; -// for(auto const& bel : beliefs) { -// std::cout << bel.toEigenVector(support) << std::endl; -// points.push_back(bel.toEigenVector(support)); -// } -// -// storm::storage::geometry::QuickHull qh; -// std::cout << "Run QH (dim=" << points[0].size() << ")" << std::endl; -// qh.generateHalfspacesFromPoints(points, true); -// std::cout << "End QH" << std::endl; -// auto filteredPoints = qh.getRelevantVertices(); -// //if (filteredPoints.size() < points.size()) { -// // for(auto const& fp : filteredPoints) { -// // beliefs.emplace(manager,filteredPoints); -// // } -// //} -// } + template + uint64_t NondeterministicBeliefTracker::reduce() { + std::shared_ptr solverFactory = std::make_shared(); + storm::storage::geometry::ReduceVertexCloud rvc(solverFactory); + std::vector> points; + std::vector::iterator> iterators; + for (auto it = beliefs.begin(); it != beliefs.end(); ++it) { + // TODO get rid of the getBeliefMap function. + points.push_back(it->getBeliefMap()); + iterators.push_back(it); + } + storm::storage::BitVector eliminate = ~rvc.eliminate(points, pomdp.getNumberOfStates()); + + auto selectedIterators = storm::utility::vector::filterVector(iterators, eliminate); + for (auto iter : selectedIterators) { + beliefs.erase(iter); + } + return eliminate.getNumberOfSetBits(); + } template class SparseBeliefState; template bool operator==(SparseBeliefState const&, SparseBeliefState const&); template class NondeterministicBeliefTracker>; - template class ObservationDenseBeliefState; - template bool operator==(ObservationDenseBeliefState const&, ObservationDenseBeliefState const&); - template class NondeterministicBeliefTracker>; + //template class ObservationDenseBeliefState; + //template bool operator==(ObservationDenseBeliefState const&, ObservationDenseBeliefState const&); + //template class NondeterministicBeliefTracker>; } diff --git a/src/storm-pomdp/generator/NondeterministicBeliefTracker.h b/src/storm-pomdp/generator/NondeterministicBeliefTracker.h index 29f315322..ef3b75a03 100644 --- a/src/storm-pomdp/generator/NondeterministicBeliefTracker.h +++ b/src/storm-pomdp/generator/NondeterministicBeliefTracker.h @@ -46,6 +46,7 @@ namespace storm { Eigen::Matrix toEigenVector(storm::storage::BitVector const& support) const; uint64_t getSupportSize() const; void setSupport(storm::storage::BitVector&) const; + std::map const& getBeliefMap() const; friend bool operator==<>(SparseBeliefState const& lhs, SparseBeliefState const& rhs); private: @@ -102,7 +103,7 @@ namespace storm { uint32_t getCurrentObservation() const; ValueType getCurrentRisk(bool max=true); void setRisk(std::vector const& risk); - //void reduce(); + uint64_t reduce(); private: diff --git a/src/storm-pomdp/transformer/ObservationTraceUnfolder.cpp b/src/storm-pomdp/transformer/ObservationTraceUnfolder.cpp index 1d8ca5399..c8f1348d0 100644 --- a/src/storm-pomdp/transformer/ObservationTraceUnfolder.cpp +++ b/src/storm-pomdp/transformer/ObservationTraceUnfolder.cpp @@ -59,15 +59,15 @@ namespace storm { for (auto const& unfoldedToOldEntry : unfoldedToOld) { transitionMatrixBuilder.newRowGroup(newRowGroupStart); - std::cout << "\tconsider new state " << unfoldedToOldEntry.first << std::endl; + //std::cout << "\tconsider new state " << unfoldedToOldEntry.first << std::endl; assert(step == 0 || newRowCount == transitionMatrixBuilder.getLastRow() + 1); svbuilder.addState(unfoldedToOldEntry.first, {}, {static_cast(unfoldedToOldEntry.second)}); uint64_t oldRowIndexStart = model.getNondeterministicChoiceIndices()[unfoldedToOldEntry.second]; uint64_t oldRowIndexEnd = model.getNondeterministicChoiceIndices()[unfoldedToOldEntry.second+1]; for (uint64_t oldRowIndex = oldRowIndexStart; oldRowIndex != oldRowIndexEnd; oldRowIndex++) { - std::cout << "\t\tconsider old action " << oldRowIndex << std::endl; - std::cout << "\t\tconsider new row nr " << newRowCount << std::endl; + //std::cout << "\t\tconsider old action " << oldRowIndex << std::endl; + //std::cout << "\t\tconsider new row nr " << newRowCount << std::endl; ValueType resetProb = storm::utility::zero(); // We first find the reset probability @@ -76,14 +76,14 @@ namespace storm { resetProb += oldRowEntry.getValue(); } } - std::cout << "\t\t\t add reset" << std::endl; + //std::cout << "\t\t\t add reset" << std::endl; // Add the resets if (resetProb != storm::utility::zero()) { transitionMatrixBuilder.addNextValue(newRowCount, 0, resetProb); } - std::cout << "\t\t\t add other transitions..." << std::endl; + //std::cout << "\t\t\t add other transitions..." << std::endl; // Now, we build the outgoing transitions. for (auto const &oldRowEntry : model.getTransitionMatrix().getRow(oldRowIndex)) { @@ -101,7 +101,7 @@ namespace storm { } else { column = entryIt->second; } - std::cout << "\t\t\t\t transition to " << column << std::endl; + //std::cout << "\t\t\t\t transition to " << column << std::endl; transitionMatrixBuilder.addNextValue(newRowCount, column, oldRowEntry.getValue()); } @@ -145,7 +145,6 @@ namespace storm { storm::storage::sparse::ModelComponents components; components.transitionMatrix = transitionMatrixBuilder.build(); - std::cout << components.transitionMatrix << std::endl; STORM_LOG_ASSERT(components.transitionMatrix.getRowGroupCount() == targetState + 1, "Expect row group count (" << components.transitionMatrix.getRowGroupCount() << ") one more as target state index " << targetState << ")"); diff --git a/src/storm/storage/geometry/ReduceVertexCloud.cpp b/src/storm/storage/geometry/ReduceVertexCloud.cpp new file mode 100644 index 000000000..8c092928d --- /dev/null +++ b/src/storm/storage/geometry/ReduceVertexCloud.cpp @@ -0,0 +1,119 @@ +#include "storm/storage/geometry/ReduceVertexCloud.h" +#include "storm/utility/Stopwatch.h" +#undef _DEBUG_REDUCE_VERTEX_CLOUD + + +namespace storm { + namespace storage { + namespace geometry { + + template + std::string toString(std::map const& point) { + std::stringstream sstr; + bool first = true; + for (auto const& entry : point) { + if (first) { + first = false; + } else { + sstr << ", "; + } + sstr << entry.first << " : " << entry.second; + } + return sstr.str(); + } + + template + storm::storage::BitVector ReduceVertexCloud::eliminate(std::vector> const& input, uint64_t maxdimension) { + std::shared_ptr expressionManager = std::make_shared(); + std::vector supports; + std::vector weightVariables; + std::vector weightVariableExpressions; + + for (uint64_t pointIndex = 0; pointIndex < input.size(); ++pointIndex) { + + // Compute the support vectors to quickly determine which input points could be relevant. + supports.push_back(storm::storage::BitVector(maxdimension)); + for (auto const& entry : input[pointIndex]) { + supports.back().set(entry.first, true); + } + // Add a weight variable for each input point + weightVariables.push_back(expressionManager->declareRationalVariable("w_"+ std::to_string(pointIndex))); + // For convenience and performance, obtain the expression. + weightVariableExpressions.push_back(weightVariables.back().getExpression()); + } + + std::unique_ptr smtSolver = smtSolverFactory->create(*expressionManager); + for (auto const& weightVariableExpr : weightVariableExpressions) { + //smtSolver->add((weightVariableExpr == expressionManager->rational(0.0)) || (weightVariableExpr > expressionManager->rational(0.00001))); + smtSolver->add((weightVariableExpr >= expressionManager->rational(0.0))); + smtSolver->add(weightVariableExpr < expressionManager->rational(1.0)); + } + smtSolver->add(storm::expressions::sum(weightVariableExpressions) <= expressionManager->rational(1.0 + wiggle)); + smtSolver->add(storm::expressions::sum(weightVariableExpressions) >= expressionManager->rational(1 - wiggle)); + + storm::utility::Stopwatch solverTime; + storm::utility::Stopwatch totalTime(true); + storm::storage::BitVector vertices(input.size()); + for (uint64_t pointIndex = 0; pointIndex < input.size(); ++pointIndex) { + std::cout << pointIndex << " out of " << input.size() << std::endl; + smtSolver->push(); + std::map> dimensionTerms; + for (auto const& entry : input[pointIndex]) { + dimensionTerms[entry.first] = {expressionManager->rational(-entry.second)}; + } + for (uint64_t potentialSupport = 0; potentialSupport < input.size(); ++potentialSupport) { + if (pointIndex == potentialSupport) { + smtSolver->add(weightVariableExpressions[potentialSupport] == expressionManager->rational(0.0)); + } else if (potentialSupport < pointIndex && !vertices.get(potentialSupport)) { + smtSolver->add(weightVariableExpressions[potentialSupport] == expressionManager->rational(0.0)); + } else if (supports[potentialSupport].isSubsetOf(supports[pointIndex])) { + for (auto const& entry : input[potentialSupport]) { + dimensionTerms[entry.first].push_back(weightVariableExpressions[potentialSupport] * expressionManager->rational(entry.second)); + } + } else { + smtSolver->add(weightVariableExpressions[potentialSupport] == expressionManager->rational(0.0)); + } + } + for (auto const& entry : dimensionTerms) { + smtSolver->add(storm::expressions::sum(entry.second) == expressionManager->rational(0.0)); + } + + solverTime.start(); + auto result = smtSolver->check(); + solverTime.stop(); + if (result == storm::solver::SmtSolver::CheckResult::Unsat) { +#ifdef _DEBUG_REDUCE_VERTEX_CLOUD + if (input[pointIndex].size() == 2) { + std::cout << "point " << toString(input[pointIndex]) << " is a vertex:"; + std::cout << smtSolver->getSmtLibString() << std::endl; + } +#endif + vertices.set(pointIndex, true); + } +#ifdef _DEBUG_REDUCE_VERTEX_CLOUD + else + { + std::cout << "point " << toString(input[pointIndex]) << " is a convex combination of "; + auto val = smtSolver->getModelAsValuation(); + uint64_t varIndex = 0; + for (auto const& wvar : weightVariables) { + if (!storm::utility::isZero(val.getRationalValue(wvar))) { + std::cout << toString(input[varIndex]) << " (weight: " << val.getRationalValue(wvar) << ")"; + } + varIndex++; + } + std::cout << std::endl; + } +#endif + smtSolver->pop(); + std::cout << "Solver time " << solverTime.getTimeInMilliseconds() << std::endl; + std::cout << "Total time " << totalTime.getTimeInMilliseconds() << std::endl; + } + return vertices; + + } + + template class ReduceVertexCloud; + } + } +} \ No newline at end of file diff --git a/src/storm/storage/geometry/ReduceVertexCloud.h b/src/storm/storage/geometry/ReduceVertexCloud.h new file mode 100644 index 000000000..b9132fa98 --- /dev/null +++ b/src/storm/storage/geometry/ReduceVertexCloud.h @@ -0,0 +1,29 @@ +#pragma once + +#include "storm/storage/BitVector.h" +#include "storm/solver/SmtSolver.h" +#include "storm/utility/solver.h" + +namespace storm { + namespace storage { + namespace geometry { + template + class ReduceVertexCloud { + public: + ReduceVertexCloud(std::shared_ptr& smtSolverFactory, ValueType wiggle = storm::utility::convertNumber(0.001)) + : smtSolverFactory(smtSolverFactory), wiggle(wiggle) + { + + } + + storm::storage::BitVector eliminate(std::vector> const& input, uint64_t maxdimension); + + private: + std::shared_ptr& smtSolverFactory; + ValueType wiggle; + + }; + + } + } +} \ No newline at end of file