From ac1d74bb47cc20aaea096de204757a09ad7aca10 Mon Sep 17 00:00:00 2001
From: Oliver Sander <sander@igpm.rwth-aachen.de>
Date: Fri, 13 Aug 2010 10:33:48 +0000
Subject: [PATCH] Introduce a global operator* for matrix-matrix
multiplications and use it.
That makes the code more readable, which is more important to me than
speed. The stuff is difficult enough as it is.
[[Imported from SVN: r6214]]
---
dune/gfe/localgeodesicfefunction.hh | 28 +++++++++++++++++++++-------
1 file changed, 21 insertions(+), 7 deletions(-)
diff --git a/dune/gfe/localgeodesicfefunction.hh b/dune/gfe/localgeodesicfefunction.hh
index c9a57ceb..f46dde15 100644
--- a/dune/gfe/localgeodesicfefunction.hh
+++ b/dune/gfe/localgeodesicfefunction.hh
@@ -11,6 +11,24 @@
#include <dune/gfe/svd.hh>
+//! calculates ret = A * B
+template< class K, int m, int n, int p >
+Dune::FieldMatrix< K, m, p > operator* ( const Dune::FieldMatrix< K, m, n > &A, const Dune::FieldMatrix< K, n, p > &B)
+{
+ typedef typename Dune::FieldMatrix< K, m, p > :: size_type size_type;
+ Dune::FieldMatrix< K, m, p > ret;
+
+ for( size_type i = 0; i < m; ++i ) {
+
+ for( size_type j = 0; j < p; ++j ) {
+ ret[ i ][ j ] = K( 0 );
+ for( size_type k = 0; k < n; ++k )
+ ret[ i ][ j ] += A[ i ][ k ] * B[ k ][ j ];
+ }
+ }
+ return ret;
+}
+
/** \brief A function defined by simplicial geodesic interpolation
from the reference element to a Riemannian manifold.
@@ -87,12 +105,9 @@ private:
UT[i] = 0;
}
- Dune::FieldMatrix<double,targetDim,targetDim> pseudoInv;
- Dune::FMatrixHelp::multMatrix(V,UT,pseudoInv);
-
- return pseudoInv;
+ return V*UT;
}
-
+
/** \brief The coefficient vector */
std::vector<TargetSpace> coefficients_;
@@ -164,8 +179,7 @@ evaluateDerivative(const Dune::FieldVector<ctype, dim>& local)
dFdw *= -1;
// multiply the two previous matrices: the result is the right hand side
- Dune::FieldMatrix<ctype,targetDim,dim> RHS;
- Dune::FMatrixHelp::multMatrix(dFdw,B, RHS);
+ Dune::FieldMatrix<ctype,targetDim,dim> RHS = dFdw * B;
// the actual system matrix
std::vector<ctype> w = barycentricCoordinates(local);
--
GitLab