From 2693915e66fcb513f2bd9cd6df8dfa0db48fefbe Mon Sep 17 00:00:00 2001
From: Elias Pipping <elias.pipping@fu-berlin.de>
Date: Mon, 25 Jul 2016 12:53:45 +0200
Subject: [PATCH] Incorporate addToDiagonal from StaticMatrix

---
 dune/matrix-vector/CMakeLists.txt   |  1 +
 dune/matrix-vector/addtodiagonal.hh | 44 +++++++++++++++++++++++++++++
 2 files changed, 45 insertions(+)
 create mode 100644 dune/matrix-vector/addtodiagonal.hh

diff --git a/dune/matrix-vector/CMakeLists.txt b/dune/matrix-vector/CMakeLists.txt
index 8d64790..fa409c1 100644
--- a/dune/matrix-vector/CMakeLists.txt
+++ b/dune/matrix-vector/CMakeLists.txt
@@ -1,5 +1,6 @@
 #install headers
 install(FILES
+  addtodiagonal.hh
   axpy.hh
   componentwisematrixmap.hh
   genericvectortools.hh
diff --git a/dune/matrix-vector/addtodiagonal.hh b/dune/matrix-vector/addtodiagonal.hh
new file mode 100644
index 0000000..90ad152
--- /dev/null
+++ b/dune/matrix-vector/addtodiagonal.hh
@@ -0,0 +1,44 @@
+#ifndef DUNE_MATRIX_VECTOR_ADDTODIAGONAL_HH
+#define DUNE_MATRIX_VECTOR_ADDTODIAGONAL_HH
+
+#include <dune/istl/scaledidmatrix.hh>
+
+namespace Dune {
+namespace MatrixVector {
+  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;
+  }
+
+  /*
+    the line
+
+    template <typename FieldType, int n> static void
+    addToDiagonal(Dune::ScaledIdentityMatrix<FieldType,n>& x, const
+    FieldType a)
+
+    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
+  */
+
+  template <typename FieldType, int n>
+  static void addToDiagonal(
+      Dune::ScaledIdentityMatrix<FieldType, n>& x,
+      const typename Dune::ScaledIdentityMatrix<FieldType, n>::field_type a) {
+    x.scalar() += a;
+  }
+}
+}
+
+#endif
-- 
GitLab