AMDiS 2.11-git
The Adaptive Multi-Dimensional Simulation Toolbox
 
Loading...
Searching...
No Matches
FieldMatVec.hpp
1#pragma once
2
3#include <type_traits>
4
5#include <dune/common/diagonalmatrix.hh>
6#include <dune/common/fmatrix.hh>
7#include <dune/common/fvector.hh>
8#include <dune/common/typetraits.hh>
9
10#include <amdis/common/TypeTraits.hpp>
11
12namespace std
13{
14 template <class T, int N>
15 struct common_type<Dune::FieldVector<T,N>, T>
16 {
17 using type = T;
18 };
19
20 template <class T, int N, int M>
21 struct common_type<Dune::FieldMatrix<T,N,M>, T>
22 {
23 using type = T;
24 };
25}
26
27namespace Dune
28{
29 namespace MatVec
30 {
33 template <class T>
34 struct IsMatrix : std::false_type {};
35
36 template <class T, int M, int N>
37 struct IsMatrix<FieldMatrix<T,M,N>> : std::true_type {};
38
39 template <class T, int N>
40 struct IsMatrix<DiagonalMatrix<T,N>> : std::true_type {};
41
42 template <class T>
43 struct IsScalarMatrix : std::false_type {};
44
45 template <class T>
46 struct IsScalarMatrix<FieldMatrix<T,1,1>> : std::true_type {};
47
48 template <class T>
49 struct IsScalarMatrix<DiagonalMatrix<T,1>> : std::true_type {};
50
51
52 template <class T>
53 struct IsVector : std::false_type {};
54
55 template <class T, int N>
56 struct IsVector<FieldVector<T,N>> : std::true_type {};
57
58 template <class T>
59 struct IsScalarVector : std::false_type {};
60
61 template <class T>
62 struct IsScalarVector<FieldVector<T,1>> : std::true_type {};
63
64
65
66 template <class T>
67 struct IsMatVec
68 : std::integral_constant<bool, IsMatrix<T>::value || IsVector<T>::value> {};
70
73 template <class A, class S>
74 struct MakeMatVec
75 {
76 using type = S;
77 };
78
79 template <class T, int M, int N, class S>
80 struct MakeMatVec<FieldMatrix<T,M,N>,S>
81 {
82 using type = FieldMatrix<S,M,N>;
83 };
84
85 template <class T, int N, class S>
86 struct MakeMatVec<DiagonalMatrix<T,N>,S>
87 {
88 using type = DiagonalMatrix<S,N>;
89 };
90
91 template <class T, int N, class S>
92 struct MakeMatVec<FieldVector<T,N>,S>
93 {
94 using type = FieldVector<S,N>;
95 };
97
100 template <class T>
101 decltype(auto) as_scalar(T&& t)
102 {
103 return FWD(t);
104 }
105
106 template <class T>
107 T as_scalar(FieldVector<T,1> const& t)
108 {
109 return t[0];
110 }
111
112 template <class T>
113 T as_scalar(FieldMatrix<T,1,1> const& t)
114 {
115 return t[0][0];
116 }
117
118 template <class T>
119 T as_scalar(DiagonalMatrix<T,1> const& t)
120 {
121 return t.diagonal(0);
122 }
124
127 template <class T>
128 decltype(auto) as_vector(T&& t)
129 {
130 return FWD(t);
131 }
132
133 template <class T, int N>
134 FieldVector<T,N> const& as_vector(FieldMatrix<T,1,N> const& t)
135 {
136 return t[0];
137 }
138
139 template <class T, int N>
140 FieldVector<T,N>& as_vector(FieldMatrix<T,1,N>& t)
141 {
142 return t[0];
143 }
145
146
149 template <class T>
150 decltype(auto) as_matrix(T&& t)
151 {
152 return FWD(t);
153 }
154
155 template <class Mat>
156 class MatrixView;
157
158 template <class T, int N>
159 MatrixView<DiagonalMatrix<T,N>> as_matrix(DiagonalMatrix<T,N> const& mat)
160 {
161 return {mat};
162 }
164
165 // returns -a
166 template <class A>
167 auto negate(A const& a);
168
169 // returns a+b
170 template <class A, class B>
171 auto plus(A const& a, B const& b);
172
173 // returns a-b
174 template <class A, class B>
175 auto minus(A const& a, B const& b);
176
177 // returns a*b
178 template <class A, class B,
179 std::enable_if_t<IsNumber<A>::value && IsNumber<B>::value, int> = 0>
180 auto multiplies(A const& a, B const& b);
181
182 template <class A, class B,
183 std::enable_if_t<IsNumber<A>::value != IsNumber<B>::value, int> = 0>
184 auto multiplies(A const& a, B const& b);
185
186 template <class T, int N, class S>
187 auto multiplies(FieldVector<T,N> const& a, FieldVector<S,N> const& b);
188
189 template <class Mat, class Vec,
190 std::enable_if_t<IsMatrix<Mat>::value && IsVector<Vec>::value, int> = 0>
191 auto multiplies(Mat const& mat, Vec const& vec);
192
193 template <class Vec, class Mat,
194 std::enable_if_t<IsVector<Vec>::value && IsMatrix<Mat>::value, int> = 0>
195 auto multiplies(Vec const& vec, Mat const& mat);
196
197 template <class T, int L, int M, int N, class S>
198 auto multiplies(FieldMatrix<T,L,M> const& a, FieldMatrix<S,M,N> const& b);
199
200 // return a/b
201 template <class A, class B>
202 auto divides(A a, B const& b)
203 {
204 return a /= b;
205 }
206
207 } // end namespace MatVec
208
210 template <class K1, int N, int M, class K2>
211 auto operator*(Dune::FieldMatrix<K1, N, M> const& a, Dune::FieldVector<K2,M> const& b)
212 {
213 return MatVec::multiplies(a, b);
214 }
215
217 template <class K1, int N, int M, class K2>
218 auto operator*(Dune::FieldVector<K1, N> const& a, Dune::FieldMatrix<K2,N,M> const& b)
219 {
220 return MatVec::multiplies(a, b);
221 }
222
224 template <class T>
225 FieldVector<T, 2> cross(FieldVector<T, 2> const& a);
226
228 template <class T>
229 FieldVector<T, 3> cross(FieldVector<T, 3> const& a, FieldVector<T, 3> const& b);
230
232 template <class T, class S, int N>
233 auto dot(FieldVector<T,N> const& vec1, FieldVector<S,N> const& vec2);
234
235 template <class T, class S, int N>
236 auto dot(FieldMatrix<T,1,N> const& vec1, FieldMatrix<S,1,N> const& vec2);
237
238 // ----------------------------------------------------------------------------
239
241 template <class T, int N>
242 T sum(FieldVector<T, N> const& x);
243
244 template <class T, int N>
245 T sum(FieldMatrix<T, 1, N> const& x);
246
247
249 template <class T,
250 std::enable_if_t<Dune::IsNumber<T>::value, int> = 0>
251 auto unary_dot(T const& x);
252
253 template <class T, int N>
254 auto unary_dot(FieldVector<T, N> const& x);
255
256 template <class T, int N>
257 auto unary_dot(FieldMatrix<T, 1, N> const& x);
258
260 template <class T, int N>
261 auto max(FieldVector<T, N> const& x);
262
263 template <class T, int N>
264 auto max(FieldMatrix<T, 1, N> const& x);
265
267 template <class T, int N>
268 auto min(FieldVector<T, N> const& x);
269
270 template <class T, int N>
271 auto min(FieldMatrix<T, 1, N> const& x);
272
274 template <class T, int N>
275 auto abs_max(FieldVector<T, N> const& x);
276
277 template <class T, int N>
278 auto abs_max(FieldMatrix<T, 1, N> const& x);
279
281 template <class T, int N>
282 auto abs_min(FieldVector<T, N> const& x);
283
284 template <class T, int N>
285 auto abs_min(FieldMatrix<T, 1, N> const& x);
286
287 // ----------------------------------------------------------------------------
288
292 template <class T, int N>
293 auto one_norm(FieldVector<T, N> const& x);
294
295 template <class T, int N>
296 auto one_norm(FieldMatrix<T, 1, N> const& x);
297
301 template <class T,
302 std::enable_if_t<Dune::IsNumber<T>::value, int> = 0>
303 auto two_norm(T const& x);
304
305 template <class T, int N>
306 auto two_norm(FieldVector<T, N> const& x);
307
308 template <class T, int N>
309 auto two_norm(FieldMatrix<T, 1, N> const& x);
310
314 template <int p, class T, int N>
315 auto p_norm(FieldVector<T, N> const& x);
316
317 template <int p, class T, int N>
318 auto p_norm(FieldMatrix<T, 1, N> const& x);
319
323 template <class T, int N>
324 auto infty_norm(FieldVector<T, N> const& x);
325
326 template <class T, int N>
327 auto infty_norm(FieldMatrix<T, 1, N> const& x);
328
329 // ----------------------------------------------------------------------------
330
332 template <class T,
333 std::enable_if_t<Dune::IsNumber<T>::value, int> = 0>
334 T distance(T const& lhs, T const& rhs);
335
336 template <class T, int N>
337 T distance(FieldVector<T, N> const& lhs, FieldVector<T, N> const& rhs);
338
339 // ----------------------------------------------------------------------------
340
342 template <class T, class S, int N, int M, int K>
343 auto outer(FieldMatrix<T,N,K> const& vec1, FieldMatrix<S,M,K> const& vec2);
344
345 template <class T, class S, int N, int M>
346 auto outer(FieldVector<T,N> const& vec1, FieldVector<S,M> const& vec2);
347
348 // ----------------------------------------------------------------------------
349
350 template <class T>
351 T det(FieldMatrix<T, 0, 0> const& /*mat*/);
352
354 template <class T>
355 T det(FieldMatrix<T, 1, 1> const& mat);
356
358 template <class T>
359 T det(FieldMatrix<T, 2, 2> const& mat);
360
362 template <class T>
363 T det(FieldMatrix<T, 3, 3> const& mat);
364
366 template <class T, int N>
367 T det(FieldMatrix<T, N, N> const& mat);
368
369
371 template <class T, int N>
372 auto inv(FieldMatrix<T, N, N> mat);
373
375 template <class T, int N>
376 void solve(FieldMatrix<T, N, N> const& A, FieldVector<T, N>& x, FieldVector<T, N> const& b);
377
378
380 template <class T, int N, int M>
381 T gramian(FieldMatrix<T,N,M> const& DF);
382
384 template <class T, int M>
385 T gramian(FieldMatrix<T, 1, M> const& DF);
386
387 // ----------------------------------------------------------------------------
388 // some arithmetic operations with FieldMatrix
389
390 template <class T, int M, int N>
391 FieldMatrix<T,N,M> trans(FieldMatrix<T, M, N> const& A);
392
393 template <class T, int N>
394 DiagonalMatrix<T,N> const& trans(DiagonalMatrix<T,N> const& A)
395 {
396 return A;
397 }
398
399 // -----------------------------------------------------------------------------
400
401 template <class T1, class T2, int M, int N, int L>
402 FieldMatrix<std::common_type_t<T1,T2>,M,N> multiplies_AtB(FieldMatrix<T1, L, M> const& A, FieldMatrix<T2, N, L> const& B);
403
404 template <class T1, class T2, int M, int N, int L>
405 FieldMatrix<std::common_type_t<T1,T2>,M,N> multiplies_ABt(FieldMatrix<T1, M, L> const& A, FieldMatrix<T2, N, L> const& B);
406
407
408 template <class T1, class T2, class T3, int M, int N, int L>
409 FieldMatrix<T3,M,N>& multiplies_ABt(FieldMatrix<T1, M, L> const& A, FieldMatrix<T2, N, L> const& B, FieldMatrix<T3,M,N>& C);
410
411 template <class T1, class T2, class T3, int N, int L>
412 FieldVector<T3,N>& multiplies_ABt(FieldMatrix<T1, 1, L> const& A, FieldMatrix<T2, N, L> const& B, FieldVector<T3,N>& C);
413
414 template <class T1, class T2, class T3, int N, int L>
415 FieldVector<T3,N>& multiplies_ABt(FieldVector<T1, L> const& A, FieldMatrix<T2, N, L> const& B, FieldVector<T3,N>& C);
416
417 template <class T1, class T2, class T3, int N, int L>
418 FieldMatrix<T3,1,N>& multiplies_ABt(FieldVector<T1, L> const& A, FieldMatrix<T2, N, L> const& B, FieldMatrix<T3,1,N>& C);
419
420
421 template <class T1, class T2, class T3, int M, int N>
422 FieldMatrix<T3,M,N>& multiplies_ABt(FieldMatrix<T1, M, N> const& A, DiagonalMatrix<T2, N> const& B, FieldMatrix<T3,M,N>& C);
423
424 template <class T1, class T2, class T3, int N>
425 FieldVector<T3,N>& multiplies_ABt(FieldMatrix<T1, 1, N> const& A, DiagonalMatrix<T2, N> const& B, FieldVector<T3,N>& C);
426
427 template <class T1, class T2, class T3, int N>
428 FieldVector<T3,N>& multiplies_ABt(FieldVector<T1, N> const& A, DiagonalMatrix<T2, N> const& B, FieldVector<T3,N>& C);
429
430 template <class T1, class T2, class T3, int N>
431 FieldMatrix<T3,1,N>& multiplies_ABt(FieldVector<T1, N> const& A, DiagonalMatrix<T2, N> const& B, FieldMatrix<T3,1,N>& C);
432
433 // -----------------------------------------------------------------------------
434
435 template <class T, int N>
436 T const& at(FieldMatrix<T,N,1> const& vec, std::size_t i);
437
438 template <class T, int M>
439 T const& at(FieldMatrix<T,1,M> const& vec, std::size_t i);
440
441 // necessary specialization to resolve ambiguous calls
442 template <class T>
443 T const& at(FieldMatrix<T,1,1> const& vec, std::size_t i);
444
445 template <class T, int N>
446 T const& at(FieldVector<T,N> const& vec, std::size_t i);
447
448} // end namespace Dune
449
450namespace AMDiS
451{
452 using Dune::FieldMatrix;
453 using Dune::FieldVector;
454}
455
456#include "FieldMatVec.inc.hpp"