diff --git a/src/storm/modelchecker/multiobjective/deterministicScheds/DeterministicSchedsLpChecker.cpp b/src/storm/modelchecker/multiobjective/deterministicScheds/DeterministicSchedsLpChecker.cpp index 1a354c0c2..0f91bfc7e 100644 --- a/src/storm/modelchecker/multiobjective/deterministicScheds/DeterministicSchedsLpChecker.cpp +++ b/src/storm/modelchecker/multiobjective/deterministicScheds/DeterministicSchedsLpChecker.cpp @@ -45,7 +45,9 @@ namespace storm { } bool inOutEncoding() { // + 8 - return encodingSettings().get(60); + bool result = encodingSettings().get(60); + STORM_LOG_ERROR_COND(!result || !isMinNegativeEncoding(), "inout-encoding only works without minnegative encoding."); + return result; } bool assertBottomStateSum() { // + 16 @@ -60,27 +62,26 @@ namespace storm { } template - DeterministicSchedsLpChecker::DeterministicSchedsLpChecker(Environment const& env, ModelType const& model, std::vector> const& objectiveHelper) : model(model) , objectiveHelper(objectiveHelper) { + DeterministicSchedsLpChecker::DeterministicSchedsLpChecker(Environment const& env, ModelType const& model, std::vector> const& objectiveHelper) : model(model) , objectiveHelper(objectiveHelper), numLpQueries(0) { swInit.start(); initializeLpModel(env); swInit.stop(); } template - DeterministicSchedsLpChecker::~DeterministicSchedsLpChecker() { - std::cout << "Deterministic Scheds LP CHECKER STATISTICS: " << std::endl; - std::cout << "\t" << swInit << " seconds for initialization" << std::endl; - std::cout << "\t" << swCheck << " seconds for checking, including" << std::endl; - std::cout << "\t\t" << swLpBuild << " seconds for LP building" << std::endl; - std::cout << "\t\t" << swLpSolve << " seconds for LP solving, including" << std::endl; - std::cout << "\t\t\t" << swCheckVertices << " seconds for finding the vertices of the convex hull." << std::endl; - std::cout << "\t" << swAux << " seconds for aux stuff" << std::endl; + std::string DeterministicSchedsLpChecker::getStatistics(std::string const& prefix) const { + std::stringstream out; + out << prefix << swInit << " seconds for LP initialization" << std::endl; + out << prefix << swCheckWeightVectors << " seconds for checking weight vectors" << std::endl; + out << prefix << swCheckAreas << " seconds for checking areas" << std::endl; + out << prefix << swValidate << " seconds for validating LP solutions" << std::endl; + out << prefix << numLpQueries << " calls to LP optimization" << std::endl; + return out.str(); } template void DeterministicSchedsLpChecker::setCurrentWeightVector(std::vector const& weightVector) { STORM_LOG_ASSERT(weightVector.size() == objectiveHelper.size(), "Setting a weight vector with invalid number of entries."); - swLpBuild.start(); if (!currentWeightVector.empty()) { // Pop information of the current weight vector. lpModel->pop(); @@ -101,25 +102,22 @@ namespace storm { } } lpModel->update(); - swLpBuild.stop(); } template - boost::optional> DeterministicSchedsLpChecker::check(storm::Environment const& env, Polytope area) { + boost::optional> DeterministicSchedsLpChecker::check(storm::Environment const& env, Polytope overapproximation) { STORM_LOG_ASSERT(!currentWeightVector.empty(), "Checking invoked before specifying a weight vector."); STORM_LOG_TRACE("Checking a vertex..."); - swCheck.start(); - swLpBuild.start(); lpModel->push(); - auto areaConstraints = area->getConstraints(lpModel->getManager(), currentObjectiveVariables); + auto areaConstraints = overapproximation->getConstraints(lpModel->getManager(), currentObjectiveVariables); for (auto const& c : areaConstraints) { lpModel->addConstraint("", c); } lpModel->update(); - swLpBuild.stop(); - swLpSolve.start(); swCheckVertices.start(); + swCheckWeightVectors.start(); lpModel->optimize(); - swCheckVertices.stop(); swLpSolve.stop(); + swCheckWeightVectors.stop(); + ++numLpQueries; STORM_LOG_TRACE("\t Done checking a vertex..."); boost::optional result; if (!lpModel->isInfeasible()) { @@ -130,14 +128,12 @@ namespace storm { } } lpModel->pop(); - swCheck.stop(); return result; } template std::pair>, std::vector>>> DeterministicSchedsLpChecker::check(storm::Environment const& env, storm::storage::geometry::PolytopeTree& polytopeTree, Point const& eps) { STORM_LOG_INFO("Checking " << polytopeTree.toString()); - swCheck.start(); STORM_LOG_ASSERT(!currentWeightVector.empty(), "Checking invoked before specifying a weight vector."); if (polytopeTree.isEmpty()) { return {{}, {}}; @@ -156,7 +152,6 @@ namespace storm { std::vector infeasableAreas; checkRecursive(env, polytopeTree, eps, foundPoints, infeasableAreas, 0); - swCheck.stop(); return {foundPoints, infeasableAreas}; } @@ -295,7 +290,6 @@ namespace storm { template bool DeterministicSchedsLpChecker::processEndComponents(std::vector>& ecVars) { - bool hasEndComponents = false; uint64_t ecCounter = 0; auto backwardTransitions = model.getBackwardTransitions(); @@ -362,7 +356,7 @@ namespace storm { } } } - hasEndComponents = ecCounter > 0 || storm::utility::graph::checkIfECWithChoiceExists(model.getTransitionMatrix(), backwardTransitions, storm::storage::BitVector(model.getNumberOfStates(), true), ~choicesWithValueZero); + bool hasEndComponents = ecCounter > 0 || storm::utility::graph::checkIfECWithChoiceExists(model.getTransitionMatrix(), backwardTransitions, storm::storage::BitVector(model.getNumberOfStates(), true), ~choicesWithValueZero); STORM_LOG_WARN_COND(!hasEndComponents, "Processed " << ecCounter << " End components."); return hasEndComponents; } @@ -428,7 +422,8 @@ namespace storm { // Compute upper bounds for each state std::vector visitingTimesUpperBounds = DeterministicSchedsObjectiveHelper::computeUpperBoundOnExpectedVisitingTimes(model.getTransitionMatrix(), bottomStates, nonBottomStates, hasEndComponents); - std::cout << "maximal visiting times upper bound is " << *std::max_element(visitingTimesUpperBounds.begin(), visitingTimesUpperBounds.end()) << std::endl; + ValueType largestUpperBound = *std::max_element(visitingTimesUpperBounds.begin(), visitingTimesUpperBounds.end()); + STORM_LOG_WARN_COND(largestUpperBound < storm::utility::convertNumber(1e5), "Found a large upper bound '" << storm::utility::convertNumber(largestUpperBound) << "' in the LP encoding. This might trigger numerical instabilities."); // create choiceValue variables and assert deterministic ones. std::vector choiceValVars(model.getNumberOfChoices()); for (auto const& state : nonBottomStates) { @@ -507,6 +502,7 @@ namespace storm { } else { // 'classic' backward encoding. for (uint64_t objIndex = 0; objIndex < objectiveHelper.size(); ++objIndex) { + STORM_LOG_WARN_COND(objectiveHelper[objIndex].getLargestUpperBound(env) < storm::utility::convertNumber(1e5), "Found a large upper value bound '" << storm::utility::convertNumber(objectiveHelper[objIndex].getLargestUpperBound(env)) << "' in the LP encoding. This might trigger numerical instabilities."); auto const& schedulerIndependentStates = objectiveHelper[objIndex].getSchedulerIndependentStateValues(); // Create state variables and store variables of ecs which contain a state with a scheduler independent value std::vector stateVars; @@ -633,13 +629,10 @@ namespace storm { template void DeterministicSchedsLpChecker::checkRecursive(Environment const& env, storm::storage::geometry::PolytopeTree& polytopeTree, Point const& eps, std::vector& foundPoints, std::vector& infeasableAreas, uint64_t const& depth) { - std::cout << "."; - std::cout.flush(); STORM_LOG_ASSERT(!polytopeTree.isEmpty(), "Tree node is empty"); STORM_LOG_ASSERT(!polytopeTree.getPolytope()->isEmpty(), "Tree node is empty."); STORM_LOG_TRACE("Checking at depth " << depth << ": " << polytopeTree.toString()); - swLpBuild.start(); lpModel->push(); // Assert the constraints of the current polytope auto nodeConstraints = polytopeTree.getPolytope()->getConstraints(lpModel->getManager(), currentObjectiveVariables); @@ -647,20 +640,19 @@ namespace storm { lpModel->addConstraint("", constr); } lpModel->update(); - swLpBuild.stop(); if (polytopeTree.getChildren().empty()) { // At leaf nodes we need to perform the actual check. - swLpSolve.start(); STORM_LOG_TRACE("\tSolving MILP..."); + swCheckAreas.start(); lpModel->optimize(); + swCheckAreas.stop(); + ++numLpQueries; STORM_LOG_TRACE("\tDone solving MILP..."); - swLpSolve.stop(); -#ifndef NDEBUG - STORM_PRINT_AND_LOG("Writing model to file '" << polytopeTree.toId() << ".lp'" << std::endl;); - lpModel->writeModelToFile(polytopeTree.toId() + ".lp"); -#endif + // STORM_PRINT_AND_LOG("Writing model to file '" << polytopeTree.toId() << ".lp'" << std::endl;); + // lpModel->writeModelToFile(polytopeTree.toId() + ".lp"); + if (lpModel->isInfeasible()) { infeasableAreas.push_back(polytopeTree.getPolytope()); polytopeTree.clear(); @@ -731,10 +723,8 @@ namespace storm { } STORM_LOG_TRACE("Checking DONE at depth " << depth << " with node " << polytopeTree.toString()); - swLpBuild.start(); lpModel->pop(); lpModel->update(); - swLpBuild.stop(); } template diff --git a/src/storm/modelchecker/multiobjective/deterministicScheds/DeterministicSchedsObjectiveHelper.cpp b/src/storm/modelchecker/multiobjective/deterministicScheds/DeterministicSchedsObjectiveHelper.cpp index 698ef3c9e..c51930a89 100644 --- a/src/storm/modelchecker/multiobjective/deterministicScheds/DeterministicSchedsObjectiveHelper.cpp +++ b/src/storm/modelchecker/multiobjective/deterministicScheds/DeterministicSchedsObjectiveHelper.cpp @@ -221,38 +221,14 @@ namespace storm { } template - typename ModelType::ValueType const& DeterministicSchedsObjectiveHelper::getUpperValueBoundAtState(Environment const& env, uint64_t state) const{ - if (!upperResultBounds) { - upperResultBounds = computeValueBounds(env, false, model, *objective.formula); - auto upperResultBound = objective.upperResultBound; - if (!upperResultBound.is_initialized() && storm::utility::vector::hasInfinityEntry(upperResultBounds.get())) { - STORM_LOG_THROW(objective.formula->isRewardOperatorFormula(), storm::exceptions::NotSupportedException, "The upper bound for objective " << *objective.originalFormula << " is infinity at some state. This is only supported for reachability rewards / total rewards."); - STORM_LOG_THROW(objective.formula->getSubformula().isTotalRewardFormula() || objective.formula->getSubformula().isEventuallyFormula(), storm::exceptions::NotSupportedException, "The upper bound for objective " << *objective.originalFormula << " is infinity at some state. This is only supported for reachability rewards / total rewards."); - auto rewards = getTotalRewardVector(model, *objective.formula); - auto zeroValuedStates = storm::utility::vector::filterZero(upperResultBounds.get()); - auto expVisits = computeUpperBoundOnExpectedVisitingTimes(model.getTransitionMatrix(), zeroValuedStates, ~zeroValuedStates, true); - ValueType upperBound = storm::utility::zero(); - for (uint64_t state = 0; state < expVisits.size(); ++state) { - ValueType maxReward = storm::utility::zero(); - for (auto row = model.getTransitionMatrix().getRowGroupIndices()[state], endRow = model.getTransitionMatrix().getRowGroupIndices()[state + 1]; row < endRow; ++row) { - maxReward = std::max(maxReward, rewards[row]); - } - upperBound += expVisits[state] * maxReward; - } - upperResultBound = upperBound; - } - storm::utility::vector::clip(upperResultBounds.get(), objective.lowerResultBound, upperResultBound); - } + typename ModelType::ValueType const& DeterministicSchedsObjectiveHelper::getUpperValueBoundAtState(Environment const& env, uint64_t state) const { + computeUpperBounds(env); return upperResultBounds.get()[state]; } template - typename ModelType::ValueType const& DeterministicSchedsObjectiveHelper::getLowerValueBoundAtState(Environment const& env, uint64_t state) const{ - if (!lowerResultBounds) { - lowerResultBounds = computeValueBounds(env, true, model, *objective.formula); - storm::utility::vector::clip(lowerResultBounds.get(), objective.lowerResultBound, objective.upperResultBound); - STORM_LOG_THROW(!storm::utility::vector::hasInfinityEntry(lowerResultBounds.get()), storm::exceptions::NotSupportedException, "The lower bound for objective " << *objective.originalFormula << " is infinity at some state. This is not supported."); - } + typename ModelType::ValueType const& DeterministicSchedsObjectiveHelper::getLowerValueBoundAtState(Environment const& env, uint64_t state) const { + computeLowerBounds(env); return lowerResultBounds.get()[state]; } @@ -511,10 +487,53 @@ namespace storm { return visitingTimesUpperBounds; } + template + typename ModelType::ValueType const& DeterministicSchedsObjectiveHelper::getLargestUpperBound(Environment const& env) const { + computeUpperBounds(env); + return *std::max_element(upperResultBounds->begin(), upperResultBounds->end()); + } + + + template + void DeterministicSchedsObjectiveHelper::computeUpperBounds(Environment const& env) const { + if (!upperResultBounds) { + upperResultBounds = computeValueBounds(env, false, model, *objective.formula); + auto upperResultBound = objective.upperResultBound; + if (!upperResultBound.is_initialized() && storm::utility::vector::hasInfinityEntry(upperResultBounds.get())) { + STORM_LOG_THROW(objective.formula->isRewardOperatorFormula(), storm::exceptions::NotSupportedException, "The upper bound for objective " << *objective.originalFormula << " is infinity at some state. This is only supported for reachability rewards / total rewards."); + STORM_LOG_THROW(objective.formula->getSubformula().isTotalRewardFormula() || objective.formula->getSubformula().isEventuallyFormula(), storm::exceptions::NotSupportedException, "The upper bound for objective " << *objective.originalFormula << " is infinity at some state. This is only supported for reachability rewards / total rewards."); + auto rewards = getTotalRewardVector(model, *objective.formula); + auto zeroValuedStates = storm::utility::vector::filterZero(upperResultBounds.get()); + auto expVisits = computeUpperBoundOnExpectedVisitingTimes(model.getTransitionMatrix(), zeroValuedStates, ~zeroValuedStates, true); + ValueType upperBound = storm::utility::zero(); + for (uint64_t state = 0; state < expVisits.size(); ++state) { + ValueType maxReward = storm::utility::zero(); + for (auto row = model.getTransitionMatrix().getRowGroupIndices()[state], endRow = model.getTransitionMatrix().getRowGroupIndices()[state + 1]; row < endRow; ++row) { + maxReward = std::max(maxReward, rewards[row]); + } + upperBound += expVisits[state] * maxReward; + } + upperResultBound = upperBound; + } + storm::utility::vector::clip(upperResultBounds.get(), objective.lowerResultBound, upperResultBound); + } + } + + template + void DeterministicSchedsObjectiveHelper::computeLowerBounds(Environment const& env) const { + if (!lowerResultBounds) { + lowerResultBounds = computeValueBounds(env, true, model, *objective.formula); + storm::utility::vector::clip(lowerResultBounds.get(), objective.lowerResultBound, objective.upperResultBound); + STORM_LOG_THROW(!storm::utility::vector::hasInfinityEntry(lowerResultBounds.get()), storm::exceptions::NotSupportedException, "The lower bound for objective " << *objective.originalFormula << " is infinity at some state. This is not supported."); + } + } + template typename ModelType::ValueType DeterministicSchedsObjectiveHelper::evaluateOnModel(Environment const& env, ModelType const& evaluatedModel) const { return evaluateOperatorFormula(env, evaluatedModel, *objective.formula)[*evaluatedModel.getInitialStates().begin()]; } + + template class DeterministicSchedsObjectiveHelper>; template class DeterministicSchedsObjectiveHelper>; diff --git a/src/storm/modelchecker/multiobjective/deterministicScheds/DeterministicSchedsObjectiveHelper.h b/src/storm/modelchecker/multiobjective/deterministicScheds/DeterministicSchedsObjectiveHelper.h index e50664e40..207122999 100644 --- a/src/storm/modelchecker/multiobjective/deterministicScheds/DeterministicSchedsObjectiveHelper.h +++ b/src/storm/modelchecker/multiobjective/deterministicScheds/DeterministicSchedsObjectiveHelper.h @@ -37,6 +37,8 @@ namespace storm { ValueType const& getUpperValueBoundAtState(Environment const& env, uint64_t state) const; ValueType const& getLowerValueBoundAtState(Environment const& env, uint64_t state) const; + ValueType const& getLargestUpperBound(Environment const& env) const; + bool minimizing() const; /*! @@ -50,6 +52,9 @@ namespace storm { private: + void computeUpperBounds(Environment const& env) const; + void computeLowerBounds(Environment const& env) const; + mutable boost::optional> schedulerIndependentStateValues; mutable boost::optional> choiceValueOffsets; mutable boost::optional> lowerResultBounds;