diff --git a/dune/matrix-vector/transformmatrix.hh b/dune/matrix-vector/transformmatrix.hh index 35df78280e82f98eec1b277dfeb97f54f19e9ab9..216b4b1fc03dae293d4d4dbec2e69005dbba2bb7 100644 --- a/dune/matrix-vector/transformmatrix.hh +++ b/dune/matrix-vector/transformmatrix.hh @@ -271,6 +271,23 @@ namespace MatrixVector { isScalar<TransformationMatrix2>()>::addTransformedMatrix(A, T1, B, T2); } + template <class MatrixBlockA, class TransformationMatrix1, class MatrixB, + class TransformationMatrix2> + static void transformMatrixPattern(Dune::BCRSMatrix<MatrixBlockA>& A, + const TransformationMatrix1& T1, + const MatrixB& B, + const TransformationMatrix2& T2) { + Dune::MatrixIndexSet indices(T1.M(), T2.M()); + for (size_t k = 0; k < B.N(); ++k) + for (auto BIt = B[k].begin(); BIt != B[k].end(); ++BIt) { + size_t l = BIt.index(); + for (auto T1kIt = T1[k].begin(); T1kIt != T1[k].end(); ++T1kIt) + for (auto T2lIt = T2[l].begin(); T2lIt != T2[l].end(); ++T2lIt) + indices.add(T1kIt.index(), T2lIt.index()); + } + indices.exportIdx(A); + } + template <class MatrixBlockA, class TransformationMatrix1, class MatrixB, class TransformationMatrix2> static void transformMatrix(Dune::BCRSMatrix<MatrixBlockA>& A, @@ -289,23 +306,6 @@ namespace MatrixVector { } } } - - template <class MatrixBlockA, class TransformationMatrix1, class MatrixB, - class TransformationMatrix2> - static void transformMatrixPattern(Dune::BCRSMatrix<MatrixBlockA>& A, - const TransformationMatrix1& T1, - const MatrixB& B, - const TransformationMatrix2& T2) { - Dune::MatrixIndexSet indices(T1.M(), T2.M()); - for (size_t k = 0; k < B.N(); ++k) - for (auto BIt = B[k].begin(); BIt != B[k].end(); ++BIt) { - size_t l = BIt.index(); - for (auto T1kIt = T1[k].begin(); T1kIt != T1[k].end(); ++T1kIt) - for (auto T2lIt = T2[l].begin(); T2lIt != T2[l].end(); ++T2lIt) - indices.add(T1kIt.index(), T2lIt.index()); - } - indices.exportIdx(A); - } } }