diff --git a/dune/solvers/common/staticmatrixtools.hh b/dune/solvers/common/staticmatrixtools.hh index 522af0eb76eb9c568ba82c7668165d9f3dfbd2a5..733407ce05722122cfa3cd8a334873c41dc189d1 100644 --- a/dune/solvers/common/staticmatrixtools.hh +++ b/dune/solvers/common/staticmatrixtools.hh @@ -9,6 +9,7 @@ #include "dune/istl/bcrsmatrix.hh" #include "dune/istl/matrixindexset.hh" +#include <dune/matrix-vector/addtodiagonal.hh> #include <dune/matrix-vector/axpy.hh> #include <dune/matrix-vector/scalartraits.hh> @@ -58,21 +59,7 @@ namespace StaticMatrix template<class Matrix> static void addToDiagonal(Matrix& x, const typename Matrix::field_type a) { - for(int i=0; i<Matrix::rows; ++i) - x[i][i] += a; - } - - template <typename FieldType, int n> -// static void addToDiagonal(Dune::ScaledIdentityMatrix<FieldType,n>& x, const FieldType a) - // the commented line should NOT be used as it may lead to ambiguities in which case the general method above will be used. - // an example is the line taken from massassembler.hh (line 83 at the time): - // StaticMatrix::addToDiagonal(localMatrix[i][i], values[i] * zi); - // where values[i] is of type FieldVector<FieldType,1> and zi is a double. This call then does not exactly fit this specialization (without the implicit cast of FV<FT,1>) - // and hence some wild template voodoo decides which of the templates is to be taken - in this case with a gcc 4.4.5 that was the one above leading to a wrong factor n - // in the diagonal value - static void addToDiagonal(Dune::ScaledIdentityMatrix<FieldType,n>& x, const typename Dune::ScaledIdentityMatrix<FieldType,n>::field_type a) - { - x.scalar() += a; + Dune::MatrixVector::addToDiagonal(x, a); } // add transformed matrix A += T1^t*B*T2 ******************************************************