MFEM  v4.1.0
Finite element discretization library
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Pages
tmatrix.hpp
Go to the documentation of this file.
1 // Copyright (c) 2010-2020, Lawrence Livermore National Security, LLC. Produced
2 // at the Lawrence Livermore National Laboratory. All Rights reserved. See files
3 // LICENSE and NOTICE for details. LLNL-CODE-806117.
4 //
5 // This file is part of the MFEM library. For more information and source code
6 // availability visit https://mfem.org.
7 //
8 // MFEM is free software; you can redistribute it and/or modify it under the
9 // terms of the BSD-3 license. We welcome feedback and contributions, see file
10 // CONTRIBUTING.md for details.
11 
12 #ifndef MFEM_TEMPLATE_MATRIX
13 #define MFEM_TEMPLATE_MATRIX
14 
15 #include "../config/tconfig.hpp"
16 #include "../general/tassign.hpp"
17 #include "../general/cuda.hpp"
18 #include "../general/hip.hpp"
19 
20 namespace mfem
21 {
22 
23 // Matrix-matrix products
24 
25 // C {=|+=} A.B -- simple version (no blocks)
26 template <bool Add,
27  typename A_layout_t, typename A_data_t,
28  typename B_layout_t, typename B_data_t,
29  typename C_layout_t, typename C_data_t>
30 MFEM_ALWAYS_INLINE inline
31 void sMult_AB(const A_layout_t &A_layout, const A_data_t &A_data,
32  const B_layout_t &B_layout, const B_data_t &B_data,
33  const C_layout_t &C_layout, C_data_t &C_data)
34 {
35  MFEM_STATIC_ASSERT(A_layout_t::rank == 2 && B_layout_t::rank == 2 &&
36  C_layout_t::rank == 2, "invalid ranks");
37  const int A1 = A_layout_t::dim_1;
38  const int A2 = A_layout_t::dim_2;
39  const int B1 = B_layout_t::dim_1;
40  const int B2 = B_layout_t::dim_2;
41  const int C1 = C_layout_t::dim_1;
42  const int C2 = C_layout_t::dim_2;
43  MFEM_STATIC_ASSERT(A2 == B1 && A1 == C1 && B2 == C2,
44  "invalid dimensions");
45 
46  MFEM_FLOPS_ADD(Add ? 2*A1*A2*B2 : 2*A1*A2*B2-A1*B2);
47  for (int b2 = 0; b2 < B2; b2++)
48  {
49  for (int s = 0; s < A2; s++)
50  {
51  for (int a1 = 0; a1 < A1; a1++)
52  {
53  if (!Add && s == 0)
54  {
55  // C(a1,b2) = A(a1,s) * B(s,b2);
56  C_data[C_layout.ind(a1,b2)] =
57  A_data[A_layout.ind(a1,s)] * B_data[B_layout.ind(s,b2)];
58  }
59  else
60  {
61  // C(a1,b2) += A(a1,s) * B(s,b2);
62  C_data[C_layout.ind(a1,b2)] +=
63  A_data[A_layout.ind(a1,s)] * B_data[B_layout.ind(s,b2)];
64  }
65  }
66  }
67  }
68 }
69 
70 // C {=|+=} A.B -- block version
71 template <int bA1, int bA2, int bB2, // block sizes
72  bool Add,
73  typename A_layout_t, typename A_data_t,
74  typename B_layout_t, typename B_data_t,
75  typename C_layout_t, typename C_data_t>
76 MFEM_ALWAYS_INLINE inline
77 void bMult_AB(const A_layout_t &A_layout, const A_data_t &A_data,
78  const B_layout_t &B_layout, const B_data_t &B_data,
79  const C_layout_t &C_layout, C_data_t &C_data)
80 {
81  MFEM_STATIC_ASSERT(A_layout_t::rank == 2 && B_layout_t::rank == 2 &&
82  C_layout_t::rank == 2, "invalid ranks");
83  const int A1 = A_layout_t::dim_1;
84  const int A2 = A_layout_t::dim_2;
85  const int B1 = B_layout_t::dim_1;
86  const int B2 = B_layout_t::dim_2;
87  const int C1 = C_layout_t::dim_1;
88  const int C2 = C_layout_t::dim_2;
89  MFEM_STATIC_ASSERT(A2 == B1 && A1 == C1 && B2 == C2,
90  "invalid dimensions");
91 
92  const int rA1 = A1%bA1;
93  const int rA2 = A2%bA2;
94  const int rB2 = B2%bB2;
95 
96  for (int b2_b = 0; b2_b < B2/bB2; b2_b++)
97  {
98  if (A2/bA2 > 0)
99  {
100  // s_b == 0
101  for (int a1_b = 0; a1_b < A1/bA1; a1_b++)
102  {
103  sMult_AB<Add>(
104  A_layout.template sub<bA1,bA2>(a1_b*bA1,0), A_data,
105  B_layout.template sub<bA2,bB2>(0,b2_b*bB2), B_data,
106  C_layout.template sub<bA1,bB2>(a1_b*bA1,b2_b*bB2), C_data);
107  }
108  if (rA1)
109  {
110  sMult_AB<Add>(
111  A_layout.template sub<rA1,bA2>(A1-rA1,0), A_data,
112  B_layout.template sub<bA2,bB2>(0,b2_b*bB2), B_data,
113  C_layout.template sub<rA1,bB2>(A1-rA1,b2_b*bB2), C_data);
114  }
115  for (int s_b = 1; s_b < A2/bA2; s_b++)
116  {
117  for (int a1_b = 0; a1_b < A1/bA1; a1_b++)
118  {
119  sMult_AB<true>(
120  A_layout.template sub<bA1,bA2>(a1_b*bA1,s_b*bA2), A_data,
121  B_layout.template sub<bA2,bB2>(s_b*bA2,b2_b*bB2), B_data,
122  C_layout.template sub<bA1,bB2>(a1_b*bA1,b2_b*bB2), C_data);
123  }
124  if (rA1)
125  {
126  sMult_AB<true>(
127  A_layout.template sub<rA1,bA2>(A1-rA1,s_b*bA2), A_data,
128  B_layout.template sub<bA2,bB2>(s_b*bA2,b2_b*bB2), B_data,
129  C_layout.template sub<rA1,bB2>(A1-rA1,b2_b*bB2), C_data);
130  }
131  }
132  }
133  if (rA2)
134  {
135  const bool rAdd = Add || (A2/bA2 > 0);
136  for (int a1_b = 0; a1_b < A1/bA1; a1_b++)
137  {
138  sMult_AB<rAdd>(
139  A_layout.template sub<bA1,rA2>(a1_b*bA1,A2-rA2), A_data,
140  B_layout.template sub<rA2,bB2>(A2-rA2,b2_b*bB2), B_data,
141  C_layout.template sub<bA1,bB2>(a1_b*bA1,b2_b*bB2), C_data);
142  }
143  if (rA1)
144  {
145  sMult_AB<rAdd>(
146  A_layout.template sub<rA1,rA2>(A1-rA1,A2-rA2), A_data,
147  B_layout.template sub<rA2,bB2>(A2-rA2,b2_b*bB2), B_data,
148  C_layout.template sub<rA1,bB2>(A1-rA1,b2_b*bB2), C_data);
149  }
150  }
151  }
152  if (rB2)
153  {
154  if (A2/bA2 > 0)
155  {
156  // s_b == 0
157  for (int a1_b = 0; a1_b < A1/bA1; a1_b++)
158  {
159  sMult_AB<Add>(
160  A_layout.template sub<bA1,bA2>(a1_b*bA1,0), A_data,
161  B_layout.template sub<bA2,rB2>(0,B2-rB2), B_data,
162  C_layout.template sub<bA1,rB2>(a1_b*bA1,B2-rB2), C_data);
163  }
164  if (rA1)
165  {
166  sMult_AB<Add>(
167  A_layout.template sub<rA1,bA2>(A1-rA1,0), A_data,
168  B_layout.template sub<bA2,rB2>(0,B2-rB2), B_data,
169  C_layout.template sub<rA1,rB2>(A1-rA1,B2-rB2), C_data);
170  }
171  }
172  if (A2/bA2 > 1)
173  {
174  for (int s_b = 1; s_b < A2/bA2; s_b++)
175  {
176  for (int a1_b = 0; a1_b < A1/bA1; a1_b++)
177  {
178  sMult_AB<true>(
179  A_layout.template sub<bA1,bA2>(a1_b*bA1,s_b*bA2), A_data,
180  B_layout.template sub<bA2,rB2>(s_b*bA2,B2-rB2), B_data,
181  C_layout.template sub<bA1,rB2>(a1_b*bA1,B2-rB2), C_data);
182  }
183  if (rA1)
184  {
185  sMult_AB<true>(
186  A_layout.template sub<rA1,bA2>(A1-rA1,s_b*bA2), A_data,
187  B_layout.template sub<bA2,rB2>(s_b*bA2,B2-rB2), B_data,
188  C_layout.template sub<rA1,rB2>(A1-rA1,B2-rB2), C_data);
189  }
190  }
191  }
192  if (rA2)
193  {
194  const bool rAdd = Add || (A2/bA2 > 0);
195  for (int a1_b = 0; a1_b < A1/bA1; a1_b++)
196  {
197  sMult_AB<rAdd>(
198  A_layout.template sub<bA1,rA2>(a1_b*bA1,A2-rA2), A_data,
199  B_layout.template sub<rA2,rB2>(A2-rA2,B2-rB2), B_data,
200  C_layout.template sub<bA1,rB2>(a1_b*bA1,B2-rB2), C_data);
201  }
202  if (rA1)
203  {
204  sMult_AB<rAdd>(
205  A_layout.template sub<rA1,rA2>(A1-rA1,A2-rA2), A_data,
206  B_layout.template sub<rA2,rB2>(A2-rA2,B2-rB2), B_data,
207  C_layout.template sub<rA1,rB2>(A1-rA1,B2-rB2), C_data);
208  }
209  }
210  }
211 }
212 
213 template <bool Add,
214  typename A_layout_t, typename A_data_t,
215  typename B_layout_t, typename B_data_t,
216  typename C_layout_t, typename C_data_t>
217 MFEM_ALWAYS_INLINE inline
218 void Mult_AB(const A_layout_t &A_layout, const A_data_t &A_data,
219  const B_layout_t &B_layout, const B_data_t &B_data,
220  const C_layout_t &C_layout, C_data_t &C_data)
221 {
222  const int b = MFEM_TEMPLATE_BLOCK_SIZE;
223  bMult_AB<b,b,b,Add>(A_layout, A_data, B_layout, B_data, C_layout, C_data);
224 }
225 
226 
227 // Small matrix operations (determinant, adjugate,...) defined by specialization
228 
229 namespace internal
230 {
231 
232 template <int N1, int N2>
233 struct MatrixOps { };
234 
235 template <>
236 struct MatrixOps<1,1>
237 {
238  // Compute det(A).
239  template <typename scalar_t, typename layout_t, typename data_t>
240  static inline scalar_t Det(const layout_t &a, const data_t &A)
241  {
242  return A[a.ind(0,0)];
243  }
244 
245  // Compute det(A). Batched version: D[i] {=,+=,*=} det(A[i,*,*])
246  template <AssignOp::Type Op, typename A_layout_t, typename A_data_t,
247  typename D_data_t>
248  static inline void Det(const A_layout_t &a, const A_data_t &A, D_data_t &D)
249  {
250  const int M = A_layout_t::dim_1;
251  for (int i = 0; i < M; i++)
252  {
253  Assign<Op>(D[i], A[a.ind(i,0,0)]);
254  }
255  }
256 
257  // Compute B = adj(A).
258  template <typename scalar_t,
259  typename A_layout_t, typename A_data_t,
260  typename B_layout_t, typename B_data_t>
261  static inline void Adjugate(const A_layout_t &a, const A_data_t &A,
262  const B_layout_t &b, B_data_t &B)
263  {
264  B[b.ind(0,0)] = scalar_t(1);
265  }
266 
267  // Compute adj(A) and det(A).
268  template <typename scalar_t,
269  typename A_layout_t, typename A_data_t,
270  typename B_layout_t, typename B_data_t>
271  MFEM_HOST_DEVICE
272  static inline scalar_t AdjDet(const A_layout_t &a, const A_data_t &A,
273  const B_layout_t &b, B_data_t &B)
274  {
275  Adjugate<scalar_t>(a, A, b, B);
276  return Det<scalar_t>(a, A);
277  }
278 };
279 
280 template <>
281 struct MatrixOps<2,2>
282 {
283  // Compute det(A).
284  template <typename scalar_t, typename layout_t, typename data_t>
285  MFEM_HOST_DEVICE
286  static inline scalar_t Det(const layout_t &a, const data_t &A)
287  {
288  MFEM_FLOPS_ADD(3);
289  return (A[a.ind(0,0)]*A[a.ind(1,1)] -
290  A[a.ind(1,0)]*A[a.ind(0,1)]);
291  }
292 
293  // Compute det(A). Batched version: D[i] {=,+=,*=} det(A[i,*,*])
294  template <AssignOp::Type Op, typename A_layout_t, typename A_data_t,
295  typename D_data_t>
296  static inline void Det(const A_layout_t &a, const A_data_t &A, D_data_t &D)
297  {
298  const int M = A_layout_t::dim_1;
299  MFEM_FLOPS_ADD(3*M);
300  for (int i = 0; i < M; i++)
301  {
302  Assign<Op>(D[i], (A[a.ind(i,0,0)]*A[a.ind(i,1,1)] -
303  A[a.ind(i,1,0)]*A[a.ind(i,0,1)]));
304  }
305  }
306 
307  // Compute B = adj(A).
308  template <typename scalar_t,
309  typename A_layout_t, typename A_data_t,
310  typename B_layout_t, typename B_data_t>
311  MFEM_HOST_DEVICE
312  static inline void Adjugate(const A_layout_t &a, const A_data_t &A,
313  const B_layout_t &b, B_data_t &B)
314  {
315  B[b.ind(0,0)] = A[a.ind(1,1)];
316  B[b.ind(0,1)] = -A[a.ind(0,1)];
317  B[b.ind(1,0)] = -A[a.ind(1,0)];
318  B[b.ind(1,1)] = A[a.ind(0,0)];
319  }
320 
321  // Compute adj(A) and det(A).
322  template <typename scalar_t,
323  typename A_layout_t, typename A_data_t,
324  typename B_layout_t, typename B_data_t>
325  MFEM_HOST_DEVICE
326  static inline scalar_t AdjDet(const A_layout_t &a, const A_data_t &A,
327  const B_layout_t &b, B_data_t &B)
328  {
329  Adjugate<scalar_t>(a, A, b, B);
330  return Det<scalar_t>(a, A);
331  }
332 
333  template <bool symm> struct Symm;
334 };
335 
336 template <>
337 struct MatrixOps<2,2>::Symm<true>
338 {
339  template <typename A_layout_t, typename A_data_t, typename scalar_t>
340  static inline MFEM_ALWAYS_INLINE
341  void Set(const A_layout_t &a, A_data_t &A,
342  const scalar_t a11, const scalar_t a21, const scalar_t a22)
343  {
344  A[a.ind(0)] = a11;
345  A[a.ind(1)] = a21;
346  A[a.ind(2)] = a22;
347  }
348 };
349 
350 template <>
351 struct MatrixOps<2,2>::Symm<false>
352 {
353  template <typename A_layout_t, typename A_data_t, typename scalar_t>
354  static inline MFEM_ALWAYS_INLINE
355  void Set(const A_layout_t &a, A_data_t &A,
356  const scalar_t a11, const scalar_t a21, const scalar_t a22)
357  {
358  A[a.ind(0,0)] = a11;
359  A[a.ind(1,0)] = a21;
360  A[a.ind(0,1)] = a21;
361  A[a.ind(1,1)] = a22;
362  }
363 };
364 
365 template <>
366 struct MatrixOps<3,3>
367 {
368  // Compute det(A).
369  template <typename scalar_t, typename layout_t, typename data_t>
370  MFEM_HOST_DEVICE
371  static inline scalar_t Det(const layout_t &a, const data_t &A)
372  {
373  MFEM_FLOPS_ADD(14);
374  return (A[a.ind(0,0)]*(A[a.ind(1,1)]*A[a.ind(2,2)] -
375  A[a.ind(2,1)]*A[a.ind(1,2)]) -
376  A[a.ind(1,0)]*(A[a.ind(0,1)]*A[a.ind(2,2)] -
377  A[a.ind(2,1)]*A[a.ind(0,2)]) +
378  A[a.ind(2,0)]*(A[a.ind(0,1)]*A[a.ind(1,2)] -
379  A[a.ind(1,1)]*A[a.ind(0,2)]));
380  }
381 
382  // Compute det(A). Batched version: D[i] {=,+=,*=} det(A[i,*,*])
383  template <AssignOp::Type Op, typename A_layout_t, typename A_data_t,
384  typename D_data_t>
385  static inline void Det(const A_layout_t &a, const A_data_t &A, D_data_t &D)
386  {
387  const int M = A_layout_t::dim_1;
388  MFEM_FLOPS_ADD(14*M);
389  for (int i = 0; i < M; i++)
390  {
391  Assign<Op>(
392  D[i],
393  A[a.ind(i,0,0)]*(A[a.ind(i,1,1)]*A[a.ind(i,2,2)] -
394  A[a.ind(i,2,1)]*A[a.ind(i,1,2)]) -
395  A[a.ind(i,1,0)]*(A[a.ind(i,0,1)]*A[a.ind(i,2,2)] -
396  A[a.ind(i,2,1)]*A[a.ind(i,0,2)]) +
397  A[a.ind(i,2,0)]*(A[a.ind(i,0,1)]*A[a.ind(i,1,2)] -
398  A[a.ind(i,1,1)]*A[a.ind(i,0,2)]));
399  }
400  }
401 
402  // Compute B = adj(A).
403  template <typename scalar_t,
404  typename A_layout_t, typename A_data_t,
405  typename B_layout_t, typename B_data_t>
406  MFEM_HOST_DEVICE
407  static inline void Adjugate(const A_layout_t &a, const A_data_t &A,
408  const B_layout_t &b, B_data_t &B)
409  {
410  MFEM_FLOPS_ADD(27);
411  B[b.ind(0,0)] = A[a.ind(1,1)]*A[a.ind(2,2)] - A[a.ind(1,2)]*A[a.ind(2,1)];
412  B[b.ind(0,1)] = A[a.ind(0,2)]*A[a.ind(2,1)] - A[a.ind(0,1)]*A[a.ind(2,2)];
413  B[b.ind(0,2)] = A[a.ind(0,1)]*A[a.ind(1,2)] - A[a.ind(0,2)]*A[a.ind(1,1)];
414  B[b.ind(1,0)] = A[a.ind(1,2)]*A[a.ind(2,0)] - A[a.ind(1,0)]*A[a.ind(2,2)];
415  B[b.ind(1,1)] = A[a.ind(0,0)]*A[a.ind(2,2)] - A[a.ind(0,2)]*A[a.ind(2,0)];
416  B[b.ind(1,2)] = A[a.ind(0,2)]*A[a.ind(1,0)] - A[a.ind(0,0)]*A[a.ind(1,2)];
417  B[b.ind(2,0)] = A[a.ind(1,0)]*A[a.ind(2,1)] - A[a.ind(1,1)]*A[a.ind(2,0)];
418  B[b.ind(2,1)] = A[a.ind(0,1)]*A[a.ind(2,0)] - A[a.ind(0,0)]*A[a.ind(2,1)];
419  B[b.ind(2,2)] = A[a.ind(0,0)]*A[a.ind(1,1)] - A[a.ind(0,1)]*A[a.ind(1,0)];
420  }
421 
422  // Compute adj(A) and det(A).
423  template <typename scalar_t,
424  typename A_layout_t, typename A_data_t,
425  typename B_layout_t, typename B_data_t>
426  MFEM_HOST_DEVICE
427  static inline scalar_t AdjDet(const A_layout_t &a, const A_data_t &A,
428  const B_layout_t &b, B_data_t &B)
429  {
430  MFEM_FLOPS_ADD(5);
431  Adjugate<scalar_t>(a, A, b, B);
432  return (A[a.ind(0,0)]*B[b.ind(0,0)] +
433  A[a.ind(1,0)]*B[b.ind(0,1)] +
434  A[a.ind(2,0)]*B[b.ind(0,2)]);
435  }
436 
437  template <bool symm> struct Symm;
438 };
439 
440 template <>
441 struct MatrixOps<3,3>::Symm<true>
442 {
443  template <typename A_layout_t, typename A_data_t, typename scalar_t>
444  static inline MFEM_ALWAYS_INLINE
445  void Set(const A_layout_t &a, A_data_t &A,
446  const scalar_t a11, const scalar_t a21, const scalar_t a31,
447  const scalar_t a22, const scalar_t a32, const scalar_t a33)
448  {
449  A[a.ind(0)] = a11;
450  A[a.ind(1)] = a21;
451  A[a.ind(2)] = a31;
452  A[a.ind(3)] = a22;
453  A[a.ind(4)] = a32;
454  A[a.ind(5)] = a33;
455  }
456 };
457 
458 template <>
459 struct MatrixOps<3,3>::Symm<false>
460 {
461  template <typename A_layout_t, typename A_data_t, typename scalar_t>
462  static inline MFEM_ALWAYS_INLINE
463  void Set(const A_layout_t &a, A_data_t &A,
464  const scalar_t a11, const scalar_t a21, const scalar_t a31,
465  const scalar_t a22, const scalar_t a32, const scalar_t a33)
466  {
467  A[a.ind(0,0)] = a11;
468  A[a.ind(1,0)] = a21;
469  A[a.ind(2,0)] = a31;
470  A[a.ind(0,1)] = a21;
471  A[a.ind(1,1)] = a22;
472  A[a.ind(2,1)] = a32;
473  A[a.ind(0,2)] = a31;
474  A[a.ind(1,2)] = a32;
475  A[a.ind(2,2)] = a33;
476  }
477 };
478 
479 } // namespace mfem::internal
480 
481 // Compute the determinant of a (small) matrix: det(A).
482 template <typename scalar_t, typename layout_t, typename data_t>
483 MFEM_HOST_DEVICE
484 inline scalar_t TDet(const layout_t &a, const data_t &A)
485 {
486  MFEM_STATIC_ASSERT(layout_t::rank == 2, "invalid rank");
487 #if !defined(__xlC__) || (__xlC__ >= 0x0d00)
488  return internal::MatrixOps<layout_t::dim_1,layout_t::dim_2>::
489  template Det<scalar_t>(a, A);
490 #else
491  return internal::MatrixOps<layout_t::dim_1,layout_t::dim_2>::
492  Det<scalar_t>(a, A);
493 #endif
494 }
495 
496 // Compute the determinants of a set of (small) matrices: D[i] = det(A[i,*,*]).
497 // The layout of A is (M x N1 x N2) and the size of D is M.
498 template <AssignOp::Type Op, typename A_layout_t, typename A_data_t,
499  typename D_data_t>
500 MFEM_HOST_DEVICE
501 inline void TDet(const A_layout_t &a, const A_data_t &A, D_data_t &D)
502 {
503  MFEM_STATIC_ASSERT(A_layout_t::rank == 3, "invalid rank");
504 #if !defined(__xlC__) || (__xlC__ >= 0x0d00)
505  internal::MatrixOps<A_layout_t::dim_2,A_layout_t::dim_3>::
506  template Det<Op>(a, A, D);
507 #else
508  internal::MatrixOps<A_layout_t::dim_2,A_layout_t::dim_3>::
509  Det<Op>(a, A, D);
510 #endif
511 }
512 
513 // Compute the adjugate matrix of a (small) matrix: B = adj(A).
514 template <typename scalar_t,
515  typename A_layout_t, typename A_data_t,
516  typename B_layout_t, typename B_data_t>
517 inline void TAdjugate(const A_layout_t &a, const A_data_t &A,
518  const B_layout_t &b, B_data_t &B)
519 {
520  MFEM_STATIC_ASSERT(A_layout_t::rank == 2 && B_layout_t::rank == 2,
521  "invalid ranks");
522  internal::MatrixOps<A_layout_t::dim_1,A_layout_t::dim_2>::
523  template Adjugate<scalar_t>(a, A, b, B);
524 }
525 
526 // Compute the adjugate and the determinant of a (small) matrix: B = adj(A),
527 // return det(A).
528 template <typename scalar_t,
529  typename A_layout_t, typename A_data_t,
530  typename B_layout_t, typename B_data_t>
531 MFEM_HOST_DEVICE
532 inline scalar_t TAdjDet(const A_layout_t &a, const A_data_t &A,
533  const B_layout_t &b, B_data_t &B)
534 {
535  MFEM_STATIC_ASSERT(A_layout_t::rank == 2 && B_layout_t::rank == 2,
536  "invalid ranks");
537  return internal::MatrixOps<A_layout_t::dim_1,A_layout_t::dim_2>::
538  template AdjDet<scalar_t>(a, A, b, B);
539 }
540 
541 } // namespace mfem
542 
543 #endif // MFEM_TEMPLATE_MATRIX
mfem::kernels::Det
MFEM_HOST_DEVICE T Det(const T *data)
Compute the determinant of a square matrix of size dim with given data.
Definition: kernels.hpp:121
mfem::TDet
MFEM_HOST_DEVICE scalar_t TDet(const layout_t &a, const data_t &A)
Definition: tmatrix.hpp:484
mfem::TAdjDet
MFEM_HOST_DEVICE scalar_t TAdjDet(const A_layout_t &a, const A_data_t &A, const B_layout_t &b, B_data_t &B)
Definition: tmatrix.hpp:532
mfem::bMult_AB
MFEM_ALWAYS_INLINE void bMult_AB(const A_layout_t &A_layout, const A_data_t &A_data, const B_layout_t &B_layout, const B_data_t &B_data, const C_layout_t &C_layout, C_data_t &C_data)
Definition: tmatrix.hpp:77
mfem::Add
void Add(const DenseMatrix &A, const DenseMatrix &B, double alpha, DenseMatrix &C)
C = A + alpha*B.
Definition: densemat.cpp:1926
b
double b
Definition: lissajous.cpp:42
mfem::TAdjugate
void TAdjugate(const A_layout_t &a, const A_data_t &A, const B_layout_t &b, B_data_t &B)
Definition: tmatrix.hpp:517
mfem::Mult_AB
MFEM_ALWAYS_INLINE void Mult_AB(const A_layout_t &A_layout, const A_data_t &A_data, const B_layout_t &B_layout, const B_data_t &B_data, const C_layout_t &C_layout, C_data_t &C_data)
Definition: tmatrix.hpp:218
a
double a
Definition: lissajous.cpp:41
mfem::sMult_AB
MFEM_ALWAYS_INLINE void sMult_AB(const A_layout_t &A_layout, const A_data_t &A_data, const B_layout_t &B_layout, const B_data_t &B_data, const C_layout_t &C_layout, C_data_t &C_data)
Definition: tmatrix.hpp:31
Generated on Fri Feb 18 2022 18:54:18 for MFEM by   doxygen 1.8.5