diff --git a/src/storm/modelchecker/multiobjective/pcaa/PcaaWeightVectorChecker.cpp b/src/storm/modelchecker/multiobjective/pcaa/PcaaWeightVectorChecker.cpp index e201b595c..d134d8ff0 100644 --- a/src/storm/modelchecker/multiobjective/pcaa/PcaaWeightVectorChecker.cpp +++ b/src/storm/modelchecker/multiobjective/pcaa/PcaaWeightVectorChecker.cpp @@ -2,6 +2,7 @@ #include "storm/modelchecker/multiobjective/pcaa/SparseMaPcaaWeightVectorChecker.h" #include "storm/modelchecker/multiobjective/pcaa/SparseMdpPcaaWeightVectorChecker.h" +#include "storm/modelchecker/multiobjective/pcaa/SparseMdpRewardBoundedPcaaWeightVectorChecker.h" #include "storm/utility/macros.h" #include "storm/exceptions/NotSupportedException.h" @@ -33,7 +34,11 @@ namespace storm { template template>::value, int>::type> std::unique_ptr> WeightVectorCheckerFactory::create(SparseMultiObjectivePreprocessorResult const& preprocessorResult) { - return std::make_unique>(preprocessorResult); + if (preprocessorResult.containsOnlyRewardObjectives()) { + return std::make_unique>(preprocessorResult); + } else { + return std::make_unique>(preprocessorResult); + } } template diff --git a/src/storm/modelchecker/multiobjective/pcaa/SparseMdpRewardBoundedPcaaWeightVectorChecker.cpp b/src/storm/modelchecker/multiobjective/pcaa/SparseMdpRewardBoundedPcaaWeightVectorChecker.cpp new file mode 100644 index 000000000..16262a795 --- /dev/null +++ b/src/storm/modelchecker/multiobjective/pcaa/SparseMdpRewardBoundedPcaaWeightVectorChecker.cpp @@ -0,0 +1,129 @@ +#include "storm/modelchecker/multiobjective/pcaa/SparseMdpRewardBoundedPcaaWeightVectorChecker.h" + +#include "storm/adapters/RationalFunctionAdapter.h" +#include "storm/models/sparse/Mdp.h" +#include "storm/models/sparse/StandardRewardModel.h" +#include "storm/utility/macros.h" +#include "storm/utility/vector.h" +#include "storm/logic/Formulas.h" +#include "storm/solver/MinMaxLinearEquationSolver.h" +#include "storm/solver/LinearEquationSolver.h" + + +#include "storm/exceptions/InvalidPropertyException.h" +#include "storm/exceptions/InvalidOperationException.h" +#include "storm/exceptions/IllegalArgumentException.h" +#include "storm/exceptions/NotSupportedException.h" +#include "storm/exceptions/UnexpectedException.h" + + +namespace storm { + namespace modelchecker { + namespace multiobjective { + + template + SparseMdpRewardBoundedPcaaWeightVectorChecker::SparseMdpRewardBoundedPcaaWeightVectorChecker(SparseMultiObjectivePreprocessorResult const& preprocessorResult) : PcaaWeightVectorChecker(preprocessorResult.objectives), rewardUnfolding(*preprocessorResult.preprocessedModel, this->objectives, storm::storage::BitVector(preprocessorResult.preprocessedModel->getNumberOfChoices(), true), preprocessorResult.reward0EStates) { + + } + + template + void SparseMdpRewardBoundedPcaaWeightVectorChecker::check(std::vector const& weightVector) { + auto initEpoch = rewardUnfolding.getStartEpoch(); + auto epochOrder = rewardUnfolding.getEpochComputationOrder(initEpoch); + for (auto const& epoch : epochOrder) { + computeEpochSolution(epoch, weightVector); + } + + auto solution = rewardUnfolding.getInitialStateResult(initEpoch); + // Todo: we currently assume precise results... + underApproxResult = solution.objectiveValues; + overApproxResult = solution.objectiveValues; + + } + + template + void SparseMdpRewardBoundedPcaaWeightVectorChecker::computeEpochSolution(typename MultiDimensionalRewardUnfolding::Epoch const& epoch, std::vector const& weightVector) { + auto const& epochModel = rewardUnfolding.setCurrentEpoch(epoch); + std::vector::SolutionType> result(epochModel.epochMatrix.getRowGroupCount()); + + + // Formulate a min-max equation system max(A*x+b)=x for the weighted sum of the objectives + std::vector b(epochModel.epochMatrix.getRowCount(), storm::utility::zero()); + for (uint64_t objIndex = 0; objIndex < this->objectives.size(); ++objIndex) { + ValueType weight = storm::solver::minimize(this->objectives[objIndex].formula->getOptimalityType()) ? -weightVector[objIndex] : weightVector[objIndex]; + if (!storm::utility::isZero(weight)) { + std::vector const& objectiveReward = epochModel.objectiveRewards[objIndex]; + for (auto const& choice : epochModel.objectiveRewardFilter[objIndex]) { + b[choice] += weight * objectiveReward[choice]; + } + } + } + auto stepSolutionIt = epochModel.stepSolutions.begin(); + for (auto const& choice : epochModel.stepChoices) { + b[choice] += stepSolutionIt->weightedValue; + ++stepSolutionIt; + } + + // Invoke the min max solver + storm::solver::GeneralMinMaxLinearEquationSolverFactory minMaxSolverFactory; + auto minMaxSolver = minMaxSolverFactory.create(epochModel.epochMatrix); + minMaxSolver->setOptimizationDirection(storm::solver::OptimizationDirection::Maximize); + minMaxSolver->setTrackScheduler(true); + //minMaxSolver->setCachingEnabled(true); + std::vector x(result.size(), storm::utility::zero()); + minMaxSolver->solveEquations(x, b); + for (uint64_t state = 0; state < x.size(); ++state) { + result[state].weightedValue = x[state]; + } + + // Formulate for each objective the linear equation system induced by the performed choices + auto const& choices = minMaxSolver->getSchedulerChoices(); + storm::storage::SparseMatrix subMatrix = epochModel.epochMatrix.selectRowsFromRowGroups(choices, true); + subMatrix.convertToEquationSystem(); + storm::solver::GeneralLinearEquationSolverFactory linEqSolverFactory; + auto linEqSolver = linEqSolverFactory.create(std::move(subMatrix)); + b.resize(choices.size()); + // TODO: start with a better initial guess + x.resize(choices.size()); + for (uint64_t objIndex = 0; objIndex < this->objectives.size(); ++objIndex) { + std::vector const& objectiveReward = epochModel.objectiveRewards[objIndex]; + for (uint64_t state = 0; state < choices.size(); ++state) { + uint64_t choice = epochModel.epochMatrix.getRowGroupIndices()[state] + choices[state]; + if (epochModel.objectiveRewardFilter[objIndex].get(choice)) { + b[state] = objectiveReward[choice]; + } else { + b[state] = storm::utility::zero(); + } + if (epochModel.stepChoices.get(choice)) { + b[state] += epochModel.stepSolutions[epochModel.stepChoices.getNumberOfSetBitsBeforeIndex(choice)].objectiveValues[objIndex]; + } + } + linEqSolver->solveEquations(x, b); + for (uint64_t state = 0; state < choices.size(); ++state) { + result[state].objectiveValues.push_back(x[state]); + } + } + + rewardUnfolding.setSolutionForCurrentEpoch(result); + } + + template + std::vector::ValueType> SparseMdpRewardBoundedPcaaWeightVectorChecker::getUnderApproximationOfInitialStateResults() const { + STORM_LOG_THROW(underApproxResult, storm::exceptions::InvalidOperationException, "Tried to retrieve results but check(..) has not been called before."); + return underApproxResult.get(); + } + + template + std::vector::ValueType> SparseMdpRewardBoundedPcaaWeightVectorChecker::getOverApproximationOfInitialStateResults() const { + STORM_LOG_THROW(overApproxResult, storm::exceptions::InvalidOperationException, "Tried to retrieve results but check(..) has not been called before."); + return overApproxResult.get(); + } + + template class SparseMdpRewardBoundedPcaaWeightVectorChecker>; +#ifdef STORM_HAVE_CARL + template class SparseMdpRewardBoundedPcaaWeightVectorChecker>; +#endif + + } + } +} diff --git a/src/storm/modelchecker/multiobjective/pcaa/SparseMdpRewardBoundedPcaaWeightVectorChecker.h b/src/storm/modelchecker/multiobjective/pcaa/SparseMdpRewardBoundedPcaaWeightVectorChecker.h new file mode 100644 index 000000000..55eaa0067 --- /dev/null +++ b/src/storm/modelchecker/multiobjective/pcaa/SparseMdpRewardBoundedPcaaWeightVectorChecker.h @@ -0,0 +1,59 @@ +#pragma once + +#include + +#include "storm/modelchecker/multiobjective/pcaa/PcaaWeightVectorChecker.h" +#include "storm/modelchecker/multiobjective/rewardbounded/MultiDimensionalRewardUnfolding.h" + +namespace storm { + namespace modelchecker { + namespace multiobjective { + + /*! + * Helper Class that takes preprocessed Pcaa data and a weight vector and ... + * - computes the maximal expected reward w.r.t. the weighted sum of the rewards of the individual objectives + * - extracts the scheduler that induces this maximum + * - computes for each objective the value induced by this scheduler + */ + template + class SparseMdpRewardBoundedPcaaWeightVectorChecker : public PcaaWeightVectorChecker { + public: + typedef typename SparseMdpModelType::ValueType ValueType; + typedef typename SparseMdpModelType::RewardModelType RewardModelType; + + SparseMdpRewardBoundedPcaaWeightVectorChecker(SparseMultiObjectivePreprocessorResult const& preprocessorResult); + + virtual ~SparseMdpRewardBoundedPcaaWeightVectorChecker() = default; + + /*! + * - computes the optimal expected reward w.r.t. the weighted sum of the rewards of the individual objectives + * - extracts the scheduler that induces this optimum + * - computes for each objective the value induced by this scheduler + */ + virtual void check(std::vector const& weightVector) override; + + /*! + * Retrieves the results of the individual objectives at the initial state of the given model. + * Note that check(..) has to be called before retrieving results. Otherwise, an exception is thrown. + * Also note that there is no guarantee that the under/over approximation is in fact correct + * as long as the underlying solution methods are unsound (e.g., standard value iteration). + */ + virtual std::vector getUnderApproximationOfInitialStateResults() const override; + virtual std::vector getOverApproximationOfInitialStateResults() const override; + + private: + + void computeEpochSolution(typename MultiDimensionalRewardUnfolding::Epoch const& epoch, std::vector const& weightVector); + + + MultiDimensionalRewardUnfolding rewardUnfolding; + + boost::optional> underApproxResult; + boost::optional> overApproxResult; + + + }; + + } + } +}