From db7cd1053f0388a844715d730c0901a49cda8ec3 Mon Sep 17 00:00:00 2001
From: =?UTF-8?q?Atgeirr=20Fl=C3=B8=20Rasmussen?= <atgeirr@sintef.no>
Date: Thu, 29 Oct 2020 12:56:32 +0100
Subject: [PATCH] Add the multMatrixImpl() utility.

---
 opm/simulators/linalg/MatrixBlock.hpp | 31 +++++++++++++++++++++++++++
 1 file changed, 31 insertions(+)

diff --git a/opm/simulators/linalg/MatrixBlock.hpp b/opm/simulators/linalg/MatrixBlock.hpp
index 1e20a5b96..c2e35bda8 100644
--- a/opm/simulators/linalg/MatrixBlock.hpp
+++ b/opm/simulators/linalg/MatrixBlock.hpp
@@ -458,6 +458,37 @@ namespace Opm
 {
 namespace Detail
 {
+    //! calculates ret = A * B
+    template< class TA, class TB, class TC, class PositiveSign >
+    static inline void multMatrixImpl( const TA &A, // n x m
+                                       const TB &B, // n x p
+                                       TC &ret,     // m x p
+                                       const PositiveSign )
+    {
+        typedef typename TA :: size_type size_type;
+        typedef typename TA :: field_type K;
+        assert( A.N() == B.N() );
+        assert( A.M() == ret.N() );
+        assert( B.M() == ret.M() );
+
+        const size_type n = A.N();
+        const size_type m = ret.N();
+        const size_type p = B.M();
+        for( size_type i = 0; i < m; ++i )
+        {
+            for( size_type j = 0; j < p; ++j )
+            {
+                K sum = 0;
+                for( size_type k = 0; k < n; ++k )
+                {
+                    sum += A[ i ][ k ] * B[ k ][ j ];
+                }
+                // set value depending on given sign
+                ret[ i ][ j ] = PositiveSign::value ? sum : -sum;
+            }
+        }
+    }
+
     //! calculates ret = sign * (A^T * B)
     //! TA, TB, and TC are not necessarily FieldMatrix, but those should
     //! follow the Dune::DenseMatrix interface.