diff --git a/dune/matrix-vector/transformmatrix.hh b/dune/matrix-vector/transformmatrix.hh index ee851397d334ef3cf692e3f21a31175b85490320..eb74e8b064e06ee00c6a8537d716d77fae9d6a11 100644 --- a/dune/matrix-vector/transformmatrix.hh +++ b/dune/matrix-vector/transformmatrix.hh @@ -22,17 +22,13 @@ namespace MatrixVector { const TransformationMatrix1& T1, const MatrixB& B, const TransformationMatrix2& T2) { - for (size_t i = 0; i < A.N(); ++i) { - for (auto jIt = A[i].begin(); jIt != A[i].end(); ++jIt) { - for (size_t k = 0; k < B.N(); ++k) { - if (T1.exists(k, i)) { + for (size_t i = 0; i < A.N(); ++i) + for (auto jIt = A[i].begin(); jIt != A[i].end(); ++jIt) + for (size_t k = 0; k < B.N(); ++k) + if (T1.exists(k, i)) for (auto lIt = B[k].begin(); lIt != B[k].end(); ++lIt) if (T2.exists(lIt.index(), jIt.index())) *jIt += T1[k][i] * (*lIt) * T2[lIt.index()][jIt.index()]; - } - } - } - } } // static void transformMatrix(MatrixA& A, const @@ -85,20 +81,16 @@ namespace MatrixVector { const TransformationMatrix& T2) { Dune::FieldMatrix<K1, m, n> T1transposedB; T1transposedB = 0; - for (size_t i = 0; i < MatrixA::rows; ++i) { - for (size_t k = 0; k < MatrixB::rows; ++k) { + for (size_t i = 0; i < MatrixA::rows; ++i) + for (size_t k = 0; k < MatrixB::rows; ++k) if (T1[k][i] != 0) for (size_t l = 0; l < MatrixB::cols; ++l) T1transposedB[i][l] += T1[k][i] * B[k][l]; - } - } - for (size_t k = 0; k < TransformationMatrix::rows; ++k) { - for (size_t l = 0; l < TransformationMatrix::cols; ++l) { + for (size_t k = 0; k < TransformationMatrix::rows; ++k) + for (size_t l = 0; l < TransformationMatrix::cols; ++l) if (T2[k][l] != 0) for (size_t i = 0; i < MatrixA::rows; ++i) A[i][l] += T1transposedB[i][k] * T2[k][l]; - } - } } // static void transformMatrix(MatrixA& A, const @@ -137,13 +129,11 @@ namespace MatrixVector { const ScalarB& B, const TransformationMatrix2& T2) { assert(TransformationMatrix1::rows == TransformationMatrix2::rows); - for (size_t k = 0; k < TransformationMatrix1::rows; ++k) { - for (auto Skj = T2[k].begin(); Skj != T2[k].end(); ++Skj) { + for (size_t k = 0; k < TransformationMatrix1::rows; ++k) + for (auto Skj = T2[k].begin(); Skj != T2[k].end(); ++Skj) for (auto Tki = T1[k].begin(); Tki != T1[k].end(); Tki++) if (A.exists(Tki.index(), Skj.index())) A[Tki.index()][Skj.index()] += (*Tki) * B * (*Skj); - } - } } }; @@ -185,14 +175,11 @@ namespace MatrixVector { MatrixA& A, const TransformationMatrix1& T1, const Dune::DiagonalMatrix<FieldType, n>& B, const TransformationMatrix2& T2) { - for (size_t k = 0; k < n; ++k) { - for (auto Skj = T2[k].begin(); Skj != T2[k].end(); ++Skj) { + for (size_t k = 0; k < n; ++k) + for (auto Skj = T2[k].begin(); Skj != T2[k].end(); ++Skj) for (auto Tki = T1[k].begin(); Tki != T1[k].end(); Tki++) - if (A.exists(Tki.index(), Skj.index())) { + if (A.exists(Tki.index(), Skj.index())) A[Tki.index()][Skj.index()] += (*Tki) * B.diagonal(k) * (*Skj); - } - } - } } }; @@ -221,11 +208,10 @@ namespace MatrixVector { Dune::DiagonalMatrix<FieldType, n>& A, const TransformationMatrix1& T1, const Dune::DiagonalMatrix<FieldType, n>& B, const TransformationMatrix2& T2) { - for (size_t k = 0; k < n; k++) { + for (size_t k = 0; k < n; k++) for (auto Tki = T1[k].begin(); Tki != T1[k].end(); ++Tki) A.diagonal(Tki.index()) += (*Tki) * B.diagonal(k) * T2[k][Tki.index()]; - } } };