diff --git a/src/storm/modelchecker/multiobjective/pcaa/SparseMdpRewardBoundedPcaaWeightVectorChecker.cpp b/src/storm/modelchecker/multiobjective/pcaa/SparseMdpRewardBoundedPcaaWeightVectorChecker.cpp index b898262b6..dc782524f 100644 --- a/src/storm/modelchecker/multiobjective/pcaa/SparseMdpRewardBoundedPcaaWeightVectorChecker.cpp +++ b/src/storm/modelchecker/multiobjective/pcaa/SparseMdpRewardBoundedPcaaWeightVectorChecker.cpp @@ -25,14 +25,19 @@ namespace storm { STORM_LOG_THROW(preprocessorResult.rewardFinitenessType == SparseMultiObjectivePreprocessorResult::RewardFinitenessType::AllFinite, storm::exceptions::NotSupportedException, "There is a scheduler that yields infinite reward for one objective. This is not supported."); STORM_LOG_THROW(preprocessorResult.preprocessedModel->getInitialStates().getNumberOfSetBits() == 1, storm::exceptions::NotSupportedException, "The model has multiple initial states."); + + numSchedChanges = 0; + numCheckedEpochs = 0; } template void SparseMdpRewardBoundedPcaaWeightVectorChecker::check(std::vector const& weightVector) { auto initEpoch = rewardUnfolding.getStartEpoch(); auto epochOrder = rewardUnfolding.getEpochComputationOrder(initEpoch); + + EpochCheckingData cachedData; for (auto const& epoch : epochOrder) { - computeEpochSolution(epoch, weightVector); + computeEpochSolution(epoch, weightVector, cachedData); } auto solution = rewardUnfolding.getInitialStateResult(initEpoch); @@ -45,73 +50,96 @@ namespace storm { } template - void SparseMdpRewardBoundedPcaaWeightVectorChecker::computeEpochSolution(typename MultiDimensionalRewardUnfolding::Epoch const& epoch, std::vector const& weightVector) { + void SparseMdpRewardBoundedPcaaWeightVectorChecker::computeEpochSolution(typename MultiDimensionalRewardUnfolding::Epoch const& epoch, std::vector const& weightVector, EpochCheckingData& cachedData) { auto& epochModel = rewardUnfolding.setCurrentEpoch(epoch); + updateCachedData(epochModel, cachedData); + ++numCheckedEpochs; swEqBuilding.start(); - uint64_t stateSolutionSize = this->objectives.size() + 1; + + //TODO result can now be set from the cacheData auto& result = epochModel.inStateSolutions; - result.resize(epochModel.epochInStates.getNumberOfSetBits(), typename MultiDimensionalRewardUnfolding::SolutionType(stateSolutionSize)); + result.resize(epochModel.epochInStates.getNumberOfSetBits(), typename MultiDimensionalRewardUnfolding::SolutionType(this->objectives.size() + 1)); // 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()); + swAux1.start(); + assert(cachedData.bMinMax.capacity() >= epochModel.epochMatrix.getRowCount()); + assert(cachedData.xMinMax.size() == epochModel.epochMatrix.getRowGroupCount()); + cachedData.bMinMax.assign(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]; + cachedData.bMinMax[choice] += weight * objectiveReward[choice]; } } } auto stepSolutionIt = epochModel.stepSolutions.begin(); for (auto const& choice : epochModel.stepChoices) { - b[choice] += stepSolutionIt->front(); + cachedData.bMinMax[choice] += stepSolutionIt->front(); ++stepSolutionIt; } + swAux1.stop(); // 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(epochModel.epochMatrix.getRowGroupCount(), storm::utility::zero()); swEqBuilding.stop(); swMinMaxSolving.start(); - minMaxSolver->solveEquations(x, b); + cachedData.minMaxSolver->solveEquations(cachedData.xMinMax, cachedData.bMinMax); swMinMaxSolving.stop(); swEqBuilding.start(); auto resultIt = result.begin(); for (auto const& state : epochModel.epochInStates) { - resultIt->front() = x[state]; + resultIt->front() = cachedData.xMinMax[state]; ++resultIt; } + // Check whether the linear equation solver needs to be updated + auto const& choices = cachedData.minMaxSolver->getSchedulerChoices(); + if (cachedData.schedulerChoices != choices) { + swAux2.start(); + ++numSchedChanges; + cachedData.schedulerChoices = choices; + storm::storage::SparseMatrix subMatrix = epochModel.epochMatrix.selectRowsFromRowGroups(choices, true); + subMatrix.convertToEquationSystem(); + storm::solver::GeneralLinearEquationSolverFactory linEqSolverFactory; + cachedData.linEqSolver = linEqSolverFactory.create(std::move(subMatrix)); + cachedData.linEqSolver->setCachingEnabled(true); + swAux2.stop(); + } + // 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()); + swAux3.start(); + assert(cachedData.bLinEq.size() == 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]; + auto rowGroupIndexIt = epochModel.epochMatrix.getRowGroupIndices().begin(); + auto choiceIt = choices.begin(); + auto stepChoiceIt = epochModel.stepChoices.begin(); + auto stepSolutionIt = epochModel.stepSolutions.begin(); + for (auto& b_i : cachedData.bLinEq) { + uint64_t i = *rowGroupIndexIt + *choiceIt; + if (epochModel.objectiveRewardFilter[objIndex].get(i)) { + b_i = objectiveReward[i]; } else { - b[state] = storm::utility::zero(); + b_i = storm::utility::zero(); + } + while (*stepChoiceIt < i) { + ++stepChoiceIt; + ++stepSolutionIt; } - if (epochModel.stepChoices.get(choice)) { - b[state] += epochModel.stepSolutions[epochModel.stepChoices.getNumberOfSetBitsBeforeIndex(choice)][objIndex + 1]; + if (i == *stepChoiceIt) { + b_i += (*stepSolutionIt)[objIndex + 1]; + ++stepChoiceIt; + ++stepSolutionIt; } + ++rowGroupIndexIt; + ++choiceIt; } + std::vector& x = cachedData.xLinEq[objIndex]; + assert(x.size() == choices.size()); swEqBuilding.stop(); swLinEqSolving.start(); - linEqSolver->solveEquations(x, b); + cachedData.linEqSolver->solveEquations(x, cachedData.bLinEq); swLinEqSolving.stop(); swEqBuilding.start(); resultIt = result.begin(); @@ -121,10 +149,37 @@ namespace storm { } } swEqBuilding.stop(); - + swAux3.stop(); rewardUnfolding.setSolutionForCurrentEpoch(); } + template + void SparseMdpRewardBoundedPcaaWeightVectorChecker::updateCachedData(typename MultiDimensionalRewardUnfolding::EpochModel const& epochModel, EpochCheckingData& cachedData) { + if (epochModel.epochMatrixChanged) { + swDataUpdate.start(); + + // Update the cached MinMaxSolver data + cachedData.bMinMax.resize(epochModel.epochMatrix.getRowCount()); + cachedData.xMinMax.assign(epochModel.epochMatrix.getRowGroupCount(), storm::utility::zero()); + storm::solver::GeneralMinMaxLinearEquationSolverFactory minMaxSolverFactory; + cachedData.minMaxSolver = minMaxSolverFactory.create(epochModel.epochMatrix); + cachedData.minMaxSolver->setOptimizationDirection(storm::solver::OptimizationDirection::Maximize); + cachedData.minMaxSolver->setTrackScheduler(true); + cachedData.minMaxSolver->setCachingEnabled(true); + + // Set the scheduler choices to invalid choice indices so that an update of the linEqSolver is enforced + cachedData.schedulerChoices.assign(epochModel.epochMatrix.getRowGroupCount(), std::numeric_limits::max()); + + // Update data for linear equation solving + cachedData.bLinEq.resize(epochModel.epochMatrix.getRowGroupCount()); + cachedData.xLinEq.resize(this->objectives.size()); + for (auto& x_o : cachedData.xLinEq) { + x_o.assign(epochModel.epochMatrix.getRowGroupCount(), storm::utility::zero()); + } + swDataUpdate.stop(); + } + } + 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."); diff --git a/src/storm/modelchecker/multiobjective/pcaa/SparseMdpRewardBoundedPcaaWeightVectorChecker.h b/src/storm/modelchecker/multiobjective/pcaa/SparseMdpRewardBoundedPcaaWeightVectorChecker.h index 5f05fd242..55bc30da6 100644 --- a/src/storm/modelchecker/multiobjective/pcaa/SparseMdpRewardBoundedPcaaWeightVectorChecker.h +++ b/src/storm/modelchecker/multiobjective/pcaa/SparseMdpRewardBoundedPcaaWeightVectorChecker.h @@ -4,6 +4,8 @@ #include "storm/modelchecker/multiobjective/pcaa/PcaaWeightVectorChecker.h" #include "storm/modelchecker/multiobjective/rewardbounded/MultiDimensionalRewardUnfolding.h" +#include "storm/solver/MinMaxLinearEquationSolver.h" +#include "storm/solver/LinearEquationSolver.h" #include "storm/utility/Stopwatch.h" @@ -30,9 +32,14 @@ namespace storm { std::cout << "WVC statistics: " << std::endl; std::cout << " overall: " << swAll << " seconds." << std::endl; std::cout << "---------------------------------------------" << std::endl; + std::cout << " Sched changed " << numSchedChanges << "/" << numCheckedEpochs << " times." << std::endl; + std::cout << " dataUpdate: " << swDataUpdate << " seconds." << std::endl; std::cout << " eqSysBuilding: " << swEqBuilding << " seconds." << std::endl; std::cout << " MinMaxSolving: " << swMinMaxSolving << " seconds." << std::endl; std::cout << " LinEqSolving: " << swLinEqSolving << " seconds." << std::endl; + std::cout << " Aux1StopWatch: " << swAux1 << " seconds." << std::endl; + std::cout << " Aux2StopWatch: " << swAux2 << " seconds." << std::endl; + std::cout << " Aux3StopWatch: " << swAux3 << " seconds." << std::endl; } @@ -55,9 +62,25 @@ namespace storm { private: - void computeEpochSolution(typename MultiDimensionalRewardUnfolding::Epoch const& epoch, std::vector const& weightVector); + struct EpochCheckingData { - storm::utility::Stopwatch swAll, swEqBuilding, swLinEqSolving, swMinMaxSolving; + std::vector bMinMax; + std::vector xMinMax; + std::unique_ptr> minMaxSolver; + + std::vector schedulerChoices; + + std::vector bLinEq; + std::vector> xLinEq; + std::unique_ptr> linEqSolver; + }; + + void computeEpochSolution(typename MultiDimensionalRewardUnfolding::Epoch const& epoch, std::vector const& weightVector, EpochCheckingData& cachedData); + + void updateCachedData(typename MultiDimensionalRewardUnfolding::EpochModel const& epochModel, EpochCheckingData& cachedData); + + storm::utility::Stopwatch swAll, swDataUpdate, swEqBuilding, swLinEqSolving, swMinMaxSolving, swAux1, swAux2, swAux3; + uint64_t numSchedChanges, numCheckedEpochs; MultiDimensionalRewardUnfolding rewardUnfolding; diff --git a/src/storm/modelchecker/multiobjective/rewardbounded/MultiDimensionalRewardUnfolding.cpp b/src/storm/modelchecker/multiobjective/rewardbounded/MultiDimensionalRewardUnfolding.cpp index e08866a3f..d6fc54c38 100644 --- a/src/storm/modelchecker/multiobjective/rewardbounded/MultiDimensionalRewardUnfolding.cpp +++ b/src/storm/modelchecker/multiobjective/rewardbounded/MultiDimensionalRewardUnfolding.cpp @@ -210,6 +210,9 @@ namespace storm { // Check if we need to update the current epoch class if (!currentEpoch || !sameEpochClass(epoch, currentEpoch.get())) { setCurrentEpochClass(epoch); + epochModel.epochMatrixChanged = true; + } else { + epochModel.epochMatrixChanged = false; } swSetEpoch.start(); diff --git a/src/storm/modelchecker/multiobjective/rewardbounded/MultiDimensionalRewardUnfolding.h b/src/storm/modelchecker/multiobjective/rewardbounded/MultiDimensionalRewardUnfolding.h index 31e519608..00627e42e 100644 --- a/src/storm/modelchecker/multiobjective/rewardbounded/MultiDimensionalRewardUnfolding.h +++ b/src/storm/modelchecker/multiobjective/rewardbounded/MultiDimensionalRewardUnfolding.h @@ -27,6 +27,7 @@ namespace storm { typedef typename std::conditional>::type SolutionType; struct EpochModel { + bool epochMatrixChanged; storm::storage::SparseMatrix epochMatrix; storm::storage::BitVector stepChoices; std::vector stepSolutions;