21#ifndef OPM_SMALL_DENSE_MATRIX_UTILS_HEADER_INCLUDED
22#define OPM_SMALL_DENSE_MATRIX_UTILS_HEADER_INCLUDED
24#include <dune/common/dynmatrix.hh>
31 template<
class TA,
class TB,
class TC,
class PositiveSign >
37 using size_type =
typename TA :: size_type;
38 using K =
typename TA :: field_type;
39 assert( A.N() == B.N() );
40 assert( A.M() == ret.N() );
41 assert( B.M() == ret.M() );
43 const size_type n = A.N();
44 const size_type m = ret.N();
45 const size_type p = B.M();
46 for( size_type i = 0; i < m; ++i )
48 for( size_type j = 0; j < p; ++j )
51 for( size_type k = 0; k < n; ++k )
53 sum += A[ i ][ k ] * B[ k ][ j ];
56 ret[ i ][ j ] = PositiveSign::value ? sum : -sum;
64 template<
class TA,
class TB,
class TC,
class PositiveSign >
70 using size_type =
typename TA :: size_type;
71 using K =
typename TA :: field_type;
72 assert( A.N() == B.N() );
73 assert( A.M() == ret.N() );
74 assert( B.M() == ret.M() );
76 const size_type n = A.N();
77 const size_type m = ret.N();
78 const size_type p = B.M();
79 for( size_type i = 0; i < m; ++i )
81 for( size_type j = 0; j < p; ++j )
84 for( size_type k = 0; k < n; ++k )
86 sum += A[ k ][ i ] * B[ k ][ j ];
89 ret[ i ][ j ] = PositiveSign::value ? sum : -sum;
95 template <
class DenseMatrixA,
class DenseMatrixB,
class DenseMatrixC>
97 const DenseMatrixB& B,
104 template <
class DenseMatrixA,
class DenseMatrixB,
class DenseMatrixC>
106 const DenseMatrixB& B,
114 static inline void multMatrix(
const Dune::DynamicMatrix<K>& A,
115 const Dune::DynamicMatrix<K>& B,
116 Dune::DynamicMatrix<K>& ret )
118 using size_type =
typename Dune::DynamicMatrix<K> :: size_type;
120 const size_type m = A.rows();
121 const size_type n = A.cols();
123 assert(n == B.rows() );
125 const size_type p = B.cols();
129 for( size_type i = 0; i < m; ++i )
131 for( size_type j = 0; j < p; ++j )
133 ret[ i ][ j ] = K( 0 );
134 for( size_type k = 0; k < n; ++k )
135 ret[ i ][ j ] += A[ i ][ k ] * B[ k ][ j ];
static void multMatrix(const Dune::DynamicMatrix< K > &A, const Dune::DynamicMatrix< K > &B, Dune::DynamicMatrix< K > &ret)
calculates ret = A * B
Definition: SmallDenseMatrixUtils.hpp:114
static void multMatrixTransposed(const DenseMatrixA &A, const DenseMatrixB &B, DenseMatrixC &ret)
calculates ret = A^T * B
Definition: SmallDenseMatrixUtils.hpp:96
static void multMatrixTransposedImpl(const TA &A, const TB &B, TC &ret, const PositiveSign)
Definition: SmallDenseMatrixUtils.hpp:65
static void negativeMultMatrixTransposed(const DenseMatrixA &A, const DenseMatrixB &B, DenseMatrixC &ret)
calculates ret = -A^T * B
Definition: SmallDenseMatrixUtils.hpp:105
static void multMatrixImpl(const TA &A, const TB &B, TC &ret, const PositiveSign)
calculates ret = A * B
Definition: SmallDenseMatrixUtils.hpp:32
Definition: blackoilboundaryratevector.hh:37