From f5fb67d10e3f0553ed3d990d3941a843e8381e06 Mon Sep 17 00:00:00 2001 From: Max Kahnt <max.kahnt@fu-berlin.de> Date: Tue, 26 Jul 2016 13:44:47 +0200 Subject: [PATCH] Treat ignore nodes correctly (esp. in transpose case). --- dune/matrix-vector/triangularsolve.hh | 16 ++++++++-------- 1 file changed, 8 insertions(+), 8 deletions(-) diff --git a/dune/matrix-vector/triangularsolve.hh b/dune/matrix-vector/triangularsolve.hh index a9ab4c8..9289992 100644 --- a/dune/matrix-vector/triangularsolve.hh +++ b/dune/matrix-vector/triangularsolve.hh @@ -15,16 +15,16 @@ namespace MatrixVector { if (transpose) { for (auto it = L.begin(); it != L.end(); ++it) { const size_t i = it.index(); - if (ignore != nullptr and (*ignore)[i].all()) - continue; auto cIt = it->begin(); assert(cIt.index() == i); - x[i] = b[i] / *cIt; + if (ignore == nullptr or (*ignore)[cIt.index()].none()) + x[i] = b[i] / *cIt; cIt++; for (; cIt != it->end(); ++cIt) { const size_t j = cIt.index(); assert(j > i); - b[j] -= x[i] * *cIt; + if (ignore == nullptr or (*ignore)[j].none()) + b[j] -= x[i] * *cIt; } } } else { @@ -57,16 +57,16 @@ namespace MatrixVector { if (transpose) { for (auto it = U.beforeEnd(); it != U.beforeBegin(); --it) { const size_t i = it.index(); - if (ignore != nullptr and (*ignore)[i].all()) - continue; auto cIt = it->beforeEnd(); assert(cIt.index() == i); - x[i] = b[i] / *cIt; + if (ignore == nullptr or (*ignore)[cIt.index()].none()) + x[i] = b[i] / *cIt; cIt--; for (; cIt != it->beforeBegin(); --cIt) { const size_t j = cIt.index(); assert(j < i); - b[j] -= *cIt * x[i]; + if (ignore == nullptr or (*ignore)[j].none()) + b[j] -= *cIt * x[i]; } } } else { -- GitLab