From 26362ed36f02b5222875b25c7dfbdf933fbe4172 Mon Sep 17 00:00:00 2001
From: TimQu <tim.quatmann@cs.rwth-aachen.de>
Date: Fri, 19 Jan 2018 15:58:15 +0100
Subject: [PATCH] trying an alternative implementation of qvi
---
.../IterativeMinMaxLinearEquationSolver.cpp | 82 ++++++++++++++++++-
1 file changed, 81 insertions(+), 1 deletion(-)
diff --git a/src/storm/solver/IterativeMinMaxLinearEquationSolver.cpp b/src/storm/solver/IterativeMinMaxLinearEquationSolver.cpp
index fcb715819..38313a57f 100644
--- a/src/storm/solver/IterativeMinMaxLinearEquationSolver.cpp
+++ b/src/storm/solver/IterativeMinMaxLinearEquationSolver.cpp
@@ -679,7 +679,7 @@ namespace storm {
template<OptimizationDirection dir>
void performIterationStep(storm::storage::SparseMatrix<ValueType> const& A, std::vector<ValueType> const& b) {
if (!decisionValueBlocks) {
- performIterationStepUpdateDecisionValue<dir>(A, b);
+ performIterationStepUpdateDecisionValue2<dir>(A, b);
} else {
assert(decisionValue == getPrimaryBound<dir>());
auto xIt = x.rbegin();
@@ -800,6 +800,86 @@ namespace storm {
}
}
+ template<OptimizationDirection dir>
+ void performIterationStepUpdateDecisionValue2(storm::storage::SparseMatrix<ValueType> const& A, std::vector<ValueType> const& b) {
+ auto const& groupIndices = A.getRowGroupIndices();
+ uint64_t group = groupIndices.size();
+ while (group != 0) {
+ uint64_t const& groupEnd = groupIndices[group];
+ --group;
+ uint64_t row = groupIndices[group];
+ ValueType xBest, yBest;
+ multiplyRow(row, A, b[row], xBest, yBest);
+ ++row;
+
+ // Only do more work if there are still rows in this row group
+ if (row != groupEnd) {
+ ValueType xi, yi;
+ uint64_t xyTmpIndex = 0;
+ if (hasPrimaryBound<dir>()) {
+ ValueType bestValue = xBest + yBest * getPrimaryBound<dir>();
+ for (;row < groupEnd; ++row) {
+ // Get the multiplication results
+ multiplyRow(row, A, b[row], xi, yi);
+ ValueType currentValue = xi + yi * getPrimaryBound<dir>();
+ // Check if the current row is better then the previously found one
+ if (better<dir>(currentValue, bestValue)) {
+ if (yBest < yi) {
+ // We need to store the 'old' best value as it might be relevant for the decision value
+ xTmp[xyTmpIndex] = std::move(xBest);
+ yTmp[xyTmpIndex] = std::move(yBest);
+ ++xyTmpIndex;
+ }
+ xBest = std::move(xi);
+ yBest = std::move(yi);
+ bestValue = std::move(currentValue);
+ } else if (yBest > yi) {
+ // If the value for this row is not strictly better, it might still be equal and have a better y value
+ if (currentValue == bestValue) {
+ xBest = std::move(xi);
+ yBest = std::move(yi);
+ } else {
+ xTmp[xyTmpIndex] = std::move(xi);
+ yTmp[xyTmpIndex] = std::move(yi);
+ ++xyTmpIndex;
+ }
+ }
+ }
+ } else {
+ for (;row < groupEnd; ++row) {
+ multiplyRow(row, A, b[row], xi, yi);
+ // Update the best choice
+ if (yi > yBest || (yi == yBest && better<dir>(xi, xBest))) {
+ xTmp[xyTmpIndex] = std::move(xBest);
+ yTmp[xyTmpIndex] = std::move(yBest);
+ ++xyTmpIndex;
+ xBest = std::move(xi);
+ yBest = std::move(yi);
+ } else {
+ xTmp[xyTmpIndex] = std::move(xi);
+ yTmp[xyTmpIndex] = std::move(yi);
+ ++xyTmpIndex;
+ }
+ }
+ }
+
+ // Update the decision value
+ for (uint64_t i = 0; i < xyTmpIndex; ++i) {
+ ValueType deltaY = yBest - yTmp[i];
+ if (deltaY > storm::utility::zero<ValueType>()) {
+ ValueType newDecisionValue = (xTmp[i] - xBest) / deltaY;
+ if (!hasDecisionValue || better<dir>(newDecisionValue, decisionValue)) {
+ decisionValue = std::move(newDecisionValue);
+ hasDecisionValue = true;
+ }
+ }
+ }
+ }
+ x[group] = std::move(xBest);
+ y[group] = std::move(yBest);
+ }
+ }
+
bool checkRestartCriterion() {
return false;
// iterations <= restartMaxIterations && (minimize(dir) ? restartThreshold * improvedPrimaryBound > primaryBound : restartThreshold * primaryBound > improvedPrimaryBound