5#ifndef DUNE_ALBERTA_ALGEBRA_HH
6#define DUNE_ALBERTA_ALGEBRA_HH
18 inline static FieldVector< K, 3 >
19 vectorProduct (
const FieldVector< K, 3 > &u,
const FieldVector< K, 3 > &v )
21 FieldVector< K, 3 > w;
22 w[ 0 ] = u[ 1 ] * v[ 2 ] - u[ 2 ] * v[ 1 ];
23 w[ 1 ] = u[ 2 ] * v[ 0 ] - u[ 0 ] * v[ 2 ];
24 w[ 2 ] = u[ 0 ] * v[ 1 ] - u[ 1 ] * v[ 0 ];
29 template<
class K,
int m >
30 inline static K determinant ( [[maybe_unused]]
const FieldMatrix< K, 0, m > &matrix )
36 inline static K determinant (
const FieldMatrix< K, 1, 1 > &matrix )
38 return matrix[ 0 ][ 0 ];
41 template<
class K,
int m >
42 inline static K determinant (
const FieldMatrix< K, 1, m > &matrix )
45 K sum = matrix[ 0 ][ 0 ] * matrix[ 0 ][ 0 ];
46 for(
int i = 1; i < m; ++i )
47 sum += matrix[ 0 ][ i ] * matrix[ 0 ][ i ];
52 inline static K determinant (
const FieldMatrix< K, 2, 2 > &matrix )
54 return matrix[ 0 ][ 0 ] * matrix[ 1 ][ 1 ] - matrix[ 0 ][ 1 ] * matrix[ 1 ][ 0 ];
58 inline static K determinant (
const FieldMatrix< K, 2, 3 > &matrix )
60 return vectorProduct( matrix[ 0 ], matrix[ 1 ] ).two_norm();
63 template<
class K,
int m >
64 inline static K determinant (
const FieldMatrix< K, 2, m > &matrix )
67 const K tmpA = matrix[ 0 ].two_norm2();
68 const K tmpB = matrix[ 1 ].two_norm2();
69 const K tmpC = matrix[ 0 ] * matrix[ 1 ];
70 return sqrt( tmpA * tmpB - tmpC * tmpC );
74 inline static K determinant (
const FieldMatrix< K, 3, 3 > &matrix )
76 return matrix[ 0 ] * vectorProduct( matrix[ 1 ], matrix[ 2 ] );
80 template<
class K,
int m >
81 inline static K invert ( [[maybe_unused]]
const FieldMatrix< K, 0, m > &matrix,
82 [[maybe_unused]] FieldMatrix< K, m, 0 > &inverse )
88 inline static K invert (
const FieldMatrix< K, 1, 1 > &matrix,
89 FieldMatrix< K, 1, 1 > &inverse )
91 inverse[ 0 ][ 0 ] = K( 1 ) / matrix[ 0 ][ 0 ];
92 return matrix[ 0 ][ 0 ];
95 template<
class K,
int m >
96 inline static K invert (
const FieldMatrix< K, 1, m > &matrix,
97 FieldMatrix< K, m, 1 > &inverse )
100 K detSqr = matrix[ 0 ].two_norm2();
101 K invDetSqr = K( 1 ) / detSqr;
102 for(
int i = 0; i < m; ++i )
103 inverse[ i ][ 0 ] = invDetSqr * matrix[ 0 ][ i ];
104 return sqrt( detSqr );
108 inline static K invert (
const FieldMatrix< K, 2, 2 > &matrix,
109 FieldMatrix< K, 2, 2 > &inverse )
111 K det = determinant( matrix );
112 K invDet = K( 1 ) / det;
113 inverse[ 0 ][ 0 ] = invDet * matrix[ 1 ][ 1 ];
114 inverse[ 0 ][ 1 ] = - invDet * matrix[ 0 ][ 1 ];
115 inverse[ 1 ][ 0 ] = - invDet * matrix[ 1 ][ 0 ];
116 inverse[ 1 ][ 1 ] = invDet * matrix[ 0 ][ 0 ];
120 template<
class K,
int m >
121 inline static K invert (
const FieldMatrix< K, 2, m > &matrix,
122 FieldMatrix< K, m, 2 > &inverse )
125 const K tmpA = matrix[ 0 ].two_norm2();
126 const K tmpB = matrix[ 1 ].two_norm2();
127 const K tmpC = matrix[ 0 ] * matrix[ 1 ];
128 const K detSqr = tmpA * tmpB - tmpC * tmpC;
129 const K invDetSqr = K( 1 ) / detSqr;
130 for(
int i = 0; i < m; ++i )
132 inverse[ i ][ 0 ] = invDetSqr * (tmpB * matrix[ 0 ][ i ] - tmpC * matrix[ 1 ][ i ]);
133 inverse[ i ][ 1 ] = invDetSqr * (tmpA * matrix[ 1 ][ i ] - tmpC * matrix[ 0 ][ i ]);
135 return sqrt( detSqr );
139 inline static K invert (
const FieldMatrix< K, 3, 3 > &matrix,
140 FieldMatrix< K, 3, 3 > &inverse )
142 return FMatrixHelp::invertMatrix( matrix, inverse );
Implements a matrix constructed from a given type representing a field and compile-time given number ...
Implements a vector constructed from a given type representing a field and a compile-time given size.
Dune namespace.
Definition: alignedallocator.hh:13