From 52e9a7e1ed88bcb9d3403da9a0a7e8085eb43c05 Mon Sep 17 00:00:00 2001 From: Elias Pipping <elias.pipping@fu-berlin.de> Date: Mon, 25 Jul 2016 12:57:18 +0200 Subject: [PATCH] Revert "Add (lower|upper)TriangularSolve" This reverts commit 0e1afee43e74c505335e068373f41dae922c40f7. Adding these functions to dune-matrix-vector directly instead --- dune/solvers/common/staticmatrixtools.hh | 84 +----------------------- 1 file changed, 1 insertion(+), 83 deletions(-) diff --git a/dune/solvers/common/staticmatrixtools.hh b/dune/solvers/common/staticmatrixtools.hh index f95c5791..70bd60ca 100644 --- a/dune/solvers/common/staticmatrixtools.hh +++ b/dune/solvers/common/staticmatrixtools.hh @@ -3,8 +3,6 @@ #ifndef STATIC_MATRIX_TOOL_HH #define STATIC_MATRIX_TOOL_HH -#include <cassert> - #include "dune/common/diagonalmatrix.hh" #include "dune/common/fmatrix.hh" #include "dune/istl/scaledidmatrix.hh" @@ -460,87 +458,7 @@ class StaticMatrix } } - template <class Matrix, class Vector, class BitVector> - static void lowerTriangularSolve(Matrix const& A, Vector const& b, - Vector& x, BitVector const* ignore, - bool transpose = false) { - static_assert( - Matrix::block_type::rows == 1 and Matrix::block_type::cols == 1, - "Only implemented for scalar problems"); - x = 0; - Vector r = b; - if (transpose) { - for (auto it = A.begin(); it != A.end(); ++it) { - const size_t i = it.index(); - if (ignore != nullptr and (*ignore)[i].all()) - continue; - auto cIt = it->begin(); - assert(cIt.index() == it.index()); - x[i] = r[i] / *cIt; - for (; cIt != it->end(); ++cIt) { - const size_t j = cIt.index(); - r[j] -= x[i] * *cIt; - } - } - } else { - for (auto it = A.begin(); it != A.end(); ++it) { - const size_t i = it.index(); - if (ignore != nullptr and (*ignore)[i].all()) - continue; - for (auto cIt = it->begin(); cIt != it->end(); ++cIt) { - const size_t j = cIt.index(); - if (i == j) { - x[i] = r[i] / *cIt; - break; - } - assert(j < i); - r[i] -= *cIt * x[j]; - } - } - } - } - - template <class Matrix, class Vector, class BitVector> - static void upperTriangularSolve(Matrix const& U, Vector const& b, - Vector& x, BitVector const* ignore, - bool transpose = false) { - static_assert( - Matrix::block_type::rows == 1 and Matrix::block_type::cols == 1, - "Only implemented for scalar problems"); - x = 0; - Vector r = b; - 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] = r[i] / *cIt; - cIt--; - for (; cIt != it->beforeBegin(); --cIt) { - const size_t j = cIt.index(); - assert(j < i); - r[j] -= *cIt * x[i]; - } - } - } else { - 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->begin(); - assert(cIt.index() == i); - auto diagonal = *cIt; - cIt++; - for (; cIt != it->end(); ++cIt) { - const size_t j = cIt.index(); - r[i] -= *cIt * x[j]; - } - x[i] = r[i] / diagonal; - } - } - } }; #endif + -- GitLab