AMDiS  2.10
The Adaptive Multi-Dimensional Simulation Toolbox
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 
12 namespace 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 
27 namespace 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 
209  // some arithmetic operations with FieldVector and FieldMatrix
210 
211  template <class A,
212  std::enable_if_t<MatVec::IsMatVec<A>::value, int> = 0>
213  auto operator-(A const& a)
214  {
215  return MatVec::negate(MatVec::as_scalar(a));
216  }
217 
218  template <class A, class B,
219  std::enable_if_t<MatVec::IsMatVec<A>::value || MatVec::IsMatVec<B>::value, int> = 0>
220  auto operator+(A const& a, B const& b)
221  {
222  return MatVec::plus(MatVec::as_scalar(a), MatVec::as_scalar(b));
223  }
224 
225  template <class A, class B,
226  std::enable_if_t<MatVec::IsMatVec<A>::value || MatVec::IsMatVec<B>::value, int> = 0>
227  auto operator-(A const& a, B const& b)
228  {
229  return MatVec::minus(MatVec::as_scalar(a), MatVec::as_scalar(b));
230  }
231 
232  template <class A, class B,
233  std::enable_if_t<MatVec::IsMatVec<A>::value || MatVec::IsMatVec<B>::value, int> = 0>
234  auto operator*(A const& a, B const& b)
235  {
236  return MatVec::multiplies(MatVec::as_scalar(a), MatVec::as_scalar(b));
237  }
238 
239  template <class A, class B,
240  std::enable_if_t<MatVec::IsMatVec<A>::value || MatVec::IsMatVec<B>::value, int> = 0>
241  auto operator/(A const& a, B const& b)
242  {
243  return MatVec::divides(MatVec::as_scalar(a), MatVec::as_scalar(b));
244  }
245 
246  // ----------------------------------------------------------------------------
247 
249  template <class T>
250  FieldVector<T, 2> cross(FieldVector<T, 2> const& a);
251 
253  template <class T>
254  FieldVector<T, 3> cross(FieldVector<T, 3> const& a, FieldVector<T, 3> const& b);
255 
257  template <class T, class S, int N>
258  auto dot(FieldVector<T,N> const& vec1, FieldVector<S,N> const& vec2);
259 
260  template <class T, class S, int N>
261  auto dot(FieldMatrix<T,1,N> const& vec1, FieldMatrix<S,1,N> const& vec2);
262 
263  // ----------------------------------------------------------------------------
264 
266  template <class T, int N>
267  T sum(FieldVector<T, N> const& x);
268 
269  template <class T, int N>
270  T sum(FieldMatrix<T, 1, N> const& x);
271 
272 
274  template <class T,
275  std::enable_if_t<Dune::IsNumber<T>::value, int> = 0>
276  auto unary_dot(T const& x);
277 
278  template <class T, int N>
279  auto unary_dot(FieldVector<T, N> const& x);
280 
281  template <class T, int N>
282  auto unary_dot(FieldMatrix<T, 1, N> const& x);
283 
285  template <class T, int N>
286  auto max(FieldVector<T, N> const& x);
287 
288  template <class T, int N>
289  auto max(FieldMatrix<T, 1, N> const& x);
290 
292  template <class T, int N>
293  auto min(FieldVector<T, N> const& x);
294 
295  template <class T, int N>
296  auto min(FieldMatrix<T, 1, N> const& x);
297 
299  template <class T, int N>
300  auto abs_max(FieldVector<T, N> const& x);
301 
302  template <class T, int N>
303  auto abs_max(FieldMatrix<T, 1, N> const& x);
304 
306  template <class T, int N>
307  auto abs_min(FieldVector<T, N> const& x);
308 
309  template <class T, int N>
310  auto abs_min(FieldMatrix<T, 1, N> const& x);
311 
312  // ----------------------------------------------------------------------------
313 
317  template <class T, int N>
318  auto one_norm(FieldVector<T, N> const& x);
319 
320  template <class T, int N>
321  auto one_norm(FieldMatrix<T, 1, N> const& x);
322 
326  template <class T,
327  std::enable_if_t<Dune::IsNumber<T>::value, int> = 0>
328  auto two_norm(T const& x);
329 
330  template <class T, int N>
331  auto two_norm(FieldVector<T, N> const& x);
332 
333  template <class T, int N>
334  auto two_norm(FieldMatrix<T, 1, N> const& x);
335 
339  template <int p, class T, int N>
340  auto p_norm(FieldVector<T, N> const& x);
341 
342  template <int p, class T, int N>
343  auto p_norm(FieldMatrix<T, 1, N> const& x);
344 
348  template <class T, int N>
349  auto infty_norm(FieldVector<T, N> const& x);
350 
351  template <class T, int N>
352  auto infty_norm(FieldMatrix<T, 1, N> const& x);
353 
354  // ----------------------------------------------------------------------------
355 
357  template <class T,
358  std::enable_if_t<Dune::IsNumber<T>::value, int> = 0>
359  T distance(T const& lhs, T const& rhs);
360 
361  template <class T, int N>
362  T distance(FieldVector<T, N> const& lhs, FieldVector<T, N> const& rhs);
363 
364  // ----------------------------------------------------------------------------
365 
367  template <class T, class S, int N, int M, int K>
368  auto outer(FieldMatrix<T,N,K> const& vec1, FieldMatrix<S,M,K> const& vec2);
369 
370  template <class T, class S, int N, int M>
371  auto outer(FieldVector<T,N> const& vec1, FieldVector<S,M> const& vec2);
372 
373  // ----------------------------------------------------------------------------
374 
375  template <class T>
376  T det(FieldMatrix<T, 0, 0> const& /*mat*/);
377 
379  template <class T>
380  T det(FieldMatrix<T, 1, 1> const& mat);
381 
383  template <class T>
384  T det(FieldMatrix<T, 2, 2> const& mat);
385 
387  template <class T>
388  T det(FieldMatrix<T, 3, 3> const& mat);
389 
391  template <class T, int N>
392  T det(FieldMatrix<T, N, N> const& mat);
393 
394 
396  template <class T, int N>
397  auto inv(FieldMatrix<T, N, N> mat);
398 
400  template <class T, int N>
401  void solve(FieldMatrix<T, N, N> const& A, FieldVector<T, N>& x, FieldVector<T, N> const& b);
402 
403 
405  template <class T, int N, int M>
406  T gramian(FieldMatrix<T,N,M> const& DF);
407 
409  template <class T, int M>
410  T gramian(FieldMatrix<T, 1, M> const& DF);
411 
412  // ----------------------------------------------------------------------------
413  // some arithmetic operations with FieldMatrix
414 
415  template <class T, int M, int N>
416  FieldMatrix<T,N,M> trans(FieldMatrix<T, M, N> const& A);
417 
418  template <class T, int N>
419  DiagonalMatrix<T,N> const& trans(DiagonalMatrix<T,N> const& A)
420  {
421  return A;
422  }
423 
424  // -----------------------------------------------------------------------------
425 
426  template <class T1, class T2, int M, int N, int L>
427  FieldMatrix<std::common_type_t<T1,T2>,M,N> multiplies_AtB(FieldMatrix<T1, L, M> const& A, FieldMatrix<T2, N, L> const& B);
428 
429  template <class T1, class T2, int M, int N, int L>
430  FieldMatrix<std::common_type_t<T1,T2>,M,N> multiplies_ABt(FieldMatrix<T1, M, L> const& A, FieldMatrix<T2, N, L> const& B);
431 
432 
433  template <class T1, class T2, class T3, int M, int N, int L>
434  FieldMatrix<T3,M,N>& multiplies_ABt(FieldMatrix<T1, M, L> const& A, FieldMatrix<T2, N, L> const& B, FieldMatrix<T3,M,N>& C);
435 
436  template <class T1, class T2, class T3, int N, int L>
437  FieldVector<T3,N>& multiplies_ABt(FieldMatrix<T1, 1, L> const& A, FieldMatrix<T2, N, L> const& B, FieldVector<T3,N>& C);
438 
439  template <class T1, class T2, class T3, int N, int L>
440  FieldVector<T3,N>& multiplies_ABt(FieldVector<T1, L> const& A, FieldMatrix<T2, N, L> const& B, FieldVector<T3,N>& C);
441 
442  template <class T1, class T2, class T3, int N, int L>
443  FieldMatrix<T3,1,N>& multiplies_ABt(FieldVector<T1, L> const& A, FieldMatrix<T2, N, L> const& B, FieldMatrix<T3,1,N>& C);
444 
445 
446  template <class T1, class T2, class T3, int M, int N>
447  FieldMatrix<T3,M,N>& multiplies_ABt(FieldMatrix<T1, M, N> const& A, DiagonalMatrix<T2, N> const& B, FieldMatrix<T3,M,N>& C);
448 
449  template <class T1, class T2, class T3, int N>
450  FieldVector<T3,N>& multiplies_ABt(FieldMatrix<T1, 1, N> const& A, DiagonalMatrix<T2, N> const& B, FieldVector<T3,N>& C);
451 
452  template <class T1, class T2, class T3, int N>
453  FieldVector<T3,N>& multiplies_ABt(FieldVector<T1, N> const& A, DiagonalMatrix<T2, N> const& B, FieldVector<T3,N>& C);
454 
455  template <class T1, class T2, class T3, int N>
456  FieldMatrix<T3,1,N>& multiplies_ABt(FieldVector<T1, N> const& A, DiagonalMatrix<T2, N> const& B, FieldMatrix<T3,1,N>& C);
457 
458  // -----------------------------------------------------------------------------
459 
460  template <class T, int N>
461  T const& at(FieldMatrix<T,N,1> const& vec, std::size_t i);
462 
463  template <class T, int M>
464  T const& at(FieldMatrix<T,1,M> const& vec, std::size_t i);
465 
466  // necessary specialization to resolve ambiguous calls
467  template <class T>
468  T const& at(FieldMatrix<T,1,1> const& vec, std::size_t i);
469 
470  template <class T, int N>
471  T const& at(FieldVector<T,N> const& vec, std::size_t i);
472 
473 } // end namespace Dune
474 
475 namespace AMDiS
476 {
477  using Dune::FieldMatrix;
478  using Dune::FieldVector;
479 }
480 
481 #include "FieldMatVec.inc.hpp"
Definition: AdaptiveGrid.hpp:393
Definition: FieldMatVec.hpp:12
Definition: AdaptBase.hpp:6