From 530c0fe57f2d1f24d5b6b3a6d30fb3a88ec900d3 Mon Sep 17 00:00:00 2001
From: Robert Kloefkorn <robert.kloefkorn@iris.no>
Date: Fri, 12 Feb 2016 12:05:08 +0100
Subject: [PATCH] MatrixBlock: cleanup invert implementation to be less
 confusing.

---
 .../NewtonIterationBlackoilInterleaved.cpp    | 45 ++++++++++++++-----
 1 file changed, 34 insertions(+), 11 deletions(-)

diff --git a/opm/autodiff/NewtonIterationBlackoilInterleaved.cpp b/opm/autodiff/NewtonIterationBlackoilInterleaved.cpp
index e9135e5ed..88e16d224 100644
--- a/opm/autodiff/NewtonIterationBlackoilInterleaved.cpp
+++ b/opm/autodiff/NewtonIterationBlackoilInterleaved.cpp
@@ -52,18 +52,42 @@
 
 namespace Dune
 {
-    namespace FMatrixHelp {
 
-    //! invert Matrix without changing the original matrix
-    //! return inverse matrix (specialization for n <= 3 in dune/common/fmatrix.hh)
-    template <typename K, int n>
-    static inline K invertMatrix (const FieldMatrix<K,n,n> &matrix, FieldMatrix<K,n,n> &inverse)
-    {
-        inverse.invert();
-        return K(0);
-    }
+namespace ISTLUtility {
+
+//! invert matrix by calling FMatrixHelp::invert
+template <typename K>
+static inline void invertMatrix (FieldMatrix<K,1,1> &matrix)
+{
+    FieldMatrix<K,1,1> A ( matrix );
+    FMatrixHelp::invertMatrix(A, matrix );
 }
 
+//! invert matrix by calling FMatrixHelp::invert
+template <typename K>
+static inline void invertMatrix (FieldMatrix<K,2,2> &matrix)
+{
+    FieldMatrix<K,2,2> A ( matrix );
+    FMatrixHelp::invertMatrix(A, matrix );
+}
+
+//! invert matrix by calling FMatrixHelp::invert
+template <typename K>
+static inline void invertMatrix (FieldMatrix<K,3,3> &matrix)
+{
+    FieldMatrix<K,3,3> A ( matrix );
+    FMatrixHelp::invertMatrix(A, matrix );
+}
+
+//! invert matrix by calling matrix.invert
+template <typename K, int n>
+static inline void invertMatrix (FieldMatrix<K,n,n> &matrix)
+{
+    matrix.invert();
+}
+
+} // end ISTLUtility
+
 template <class Scalar, int n, int m>
 class MatrixBlock : public Dune::FieldMatrix<Scalar, n, m>
 {
@@ -76,8 +100,7 @@ public:
     explicit MatrixBlock( const Scalar scalar = 0 ) : BaseType( scalar ) {}
     void invert()
     {
-        MatrixBlock matrix( *this );
-        Dune::FMatrixHelp::invertMatrix(matrix, *this);
+        ISTLUtility::invertMatrix( *this );
     }
     const BaseType& asBase() const { return static_cast< const BaseType& > (*this); }
     BaseType& asBase() { return static_cast< BaseType& > (*this); }