MFEM  v4.6.0
Finite element discretization library
teltrans.hpp
Go to the documentation of this file.
1 // Copyright (c) 2010-2023, 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_ELEMENT_TRANSFORMATION
13 #define MFEM_TEMPLATE_ELEMENT_TRANSFORMATION
14 
15 #include "../config/tconfig.hpp"
16 #include "tevaluator.hpp"
17 #include "../mesh/element.hpp"
18 
19 namespace mfem
20 {
21 
22 // Templated element transformation classes, cf. eltrans.?pp
23 
24 /** @brief Element transformation class, templated on a mesh type and an
25  integration rule.
26  It is constructed from a mesh (e.g. class TMesh) and shape evaluator
27  (e.g. class ShapeEvaluator) objects. Allows computation of physical
28  coordinates and Jacobian matrices corresponding to the reference integration
29  points. The desired result is specified through the template subclass Result
30  and stored in an object of the same type.
31 */
32 template <typename Mesh_t, typename IR, typename real_t = double>
34 {
35 public:
36  typedef real_t real_type;
37  typedef typename Mesh_t::FE_type FE_type;
38  typedef typename Mesh_t::FESpace_type FESpace_type;
39  typedef typename Mesh_t::nodeLayout_type nodeLayout_type;
41 
43 
44  /// Enumeration for the result type of the TElementTransformation::Eval()
45  /// method. The types can obtained by summing constants from this enumeration
46  /// and used as a template parameter in struct Result.
48  {
49  EvalNone = 0,
54  };
55 
56  /// Determines at compile-time the operations needed for given coefficient
57  /// and kernel
58  template <typename coeff_t, typename kernel_t> struct Get
59  {
60  static const int EvalOps =
61  (EvalCoordinates * coeff_t::uses_coordinates +
62  EvalJacobians * coeff_t::uses_Jacobians +
63  LoadAttributes * coeff_t::uses_attributes +
64  LoadElementIdxs * coeff_t::uses_element_idxs) |
65  (EvalJacobians * kernel_t::uses_Jacobians);
66  };
67 
68  /** @brief Templated struct Result, used to specify the type result that is
69  computed by the TElementTransformation::Eval() method and stored in this
70  structure.
71  @tparam EvalOps is a sum (bitwise or) of constants from the enum EvalOperations
72  @tparam NE is the number of elements to be processed in the Eval() method.
73  @tparam impl_traits_t specifies additional parameters and types to be used by the Eval() method
74  */
75  template<int EvalOps, typename impl_traits_t> struct Result;
76 
77  static const int dim = Mesh_t::dim;
78  static const int sdim = Mesh_t::space_dim;
79  static const int dofs = FE_type::dofs;
80  static const int qpts = IR::qpts;
81 
82 protected:
83 #ifdef MFEM_TEMPLATE_ELTRANS_HAS_NODE_DOFS
85 #endif
86 
90  const real_t *nodes;
91 
92  const Element* const *elements;
93 
94  template <typename vint_t, int NE>
95  inline MFEM_ALWAYS_INLINE
96  void SetAttributes(int el, vint_t (&attrib)[NE]) const
97  {
98  const int vsize = sizeof(vint_t)/sizeof(attrib[0][0]);
99  for (int i = 0; i < NE; i++)
100  {
101  for (int j = 0; j < vsize; i++)
102  {
103  attrib[i][j] = elements[el+j+i*vsize]->GetAttribute();
104  }
105  }
106  }
107 
108 public:
109  // Constructor.
110  TElementTransformation(const Mesh_t &mesh, const ShapeEval &eval)
111  : evaluator(eval),
112  fes(mesh.t_fes),
113  node_layout(mesh.node_layout),
114  nodes(mesh.Nodes.GetData()), // real_t = double
115  elements(mesh.m_mesh.GetElementsArray())
116  { }
117 
118  /// Evaluate coordinates and/or Jacobian matrices at quadrature points.
119  template<int EvalOps, typename impl_traits_t>
120  inline MFEM_ALWAYS_INLINE
122  {
123  F.Eval(el, *this);
124  }
125 
126 #ifdef MFEM_TEMPLATE_ENABLE_SERIALIZE
127  template<int EvalOps, typename impl_traits_t>
128  inline MFEM_ALWAYS_INLINE
129  void EvalSerialized(int el, const typename impl_traits_t::vreal_t *nodeData,
131  {
132  F.EvalSerialized(el, *this, nodeData);
133  }
134 #endif
135 
136  // Specialization of the Result<> class
137 
138  // Case EvalOps = 0 = EvalNone
139  template <typename it_t> struct Result<0,it_t>
140  {
141  static const int ne = it_t::batch_size;
142  typedef typename it_t::vreal_t vreal_t;
143  // x_type x;
144  // Jt_type Jt;
145  // int attrib[NE];
146  // int first_elem_idx;
147  inline MFEM_ALWAYS_INLINE
148  void Eval(int el, T_type &T)
149  {
150  // T.SetAttributes(el, attrib);
151  // first_elem_idx = el;
152  }
153 #ifdef MFEM_TEMPLATE_ENABLE_SERIALIZE
154  inline MFEM_ALWAYS_INLINE
155  void EvalSerialized(int el, T_type &T, const vreal_t *nodeData) { }
156 #endif
157  };
158 
159  // Case EvalOps = 1 = EvalCoordinates
160  template <typename it_t> struct Result<1,it_t>
161  {
162  static const int ne = it_t::batch_size;
163  typedef typename it_t::vreal_t vreal_t;
164 #ifdef MFEM_TEMPLATE_ELTRANS_RESULT_HAS_NODES
166 #else
167  typedef TTensor3<qpts,sdim,ne,vreal_t/*,true*/> x_type;
168 #endif
170 
172 #ifdef MFEM_TEMPLATE_ELTRANS_RESULT_HAS_NODES
174 #endif
175 
176  inline MFEM_ALWAYS_INLINE
177  void Eval(int el, T_type &T)
178  {
179 #ifdef MFEM_TEMPLATE_ELTRANS_HAS_NODE_DOFS
180  MFEM_STATIC_ASSERT(ne == 1, "only ne == 1 is supported");
182 #elif !defined(MFEM_TEMPLATE_ELTRANS_RESULT_HAS_NODES)
184 #endif
185  T.fes.SetElement(el);
186  T.fes.VectorExtract(T.node_layout, T.nodes,
188  T.evaluator.Calc(nodes_dof.layout.merge_23(), nodes_dof,
189  x.layout.merge_23(), x);
190  }
191 
192 #ifdef MFEM_TEMPLATE_ENABLE_SERIALIZE
193  inline MFEM_ALWAYS_INLINE
194  void EvalSerialized(int el, T_type &T, const vreal_t *nodeData)
195  {
196  const int SS = sizeof(nodeData[0])/sizeof(nodeData[0][0]);
197  MFEM_ASSERT(el % (SS*ne) == 0, "invalid element index: " << el);
198  T.evaluator.Calc(nodes_dof_t::layout.merge_23(),
199  &nodeData[el/SS*nodes_dof_t::size],
200  x.layout.merge_23(), x);
201  }
202 #endif
203  };
204 
205  // Case EvalOps = 2 = EvalJacobians
206  template <typename it_t> struct Result<2,it_t>
207  {
208  static const int ne = it_t::batch_size;
209  typedef typename it_t::vreal_t vreal_t;
210 #ifdef MFEM_TEMPLATE_ELTRANS_RESULT_HAS_NODES
212 #else
213  typedef TTensor4<qpts,dim,sdim,ne,vreal_t/*,true*/> Jt_type;
214 #endif
216 
218 #ifdef MFEM_TEMPLATE_ELTRANS_RESULT_HAS_NODES
220 #endif
221 
222  inline MFEM_ALWAYS_INLINE
223  void Eval(int el, T_type &T)
224  {
225 #ifdef MFEM_TEMPLATE_ELTRANS_HAS_NODE_DOFS
226  MFEM_STATIC_ASSERT(ne == 1, "only ne == 1 is supported");
228 #elif !defined(MFEM_TEMPLATE_ELTRANS_RESULT_HAS_NODES)
230 #endif
231  T.fes.SetElement(el);
232  T.fes.VectorExtract(T.node_layout, T.nodes,
234  T.evaluator.CalcGrad(nodes_dof.layout.merge_23(), nodes_dof,
235  Jt.layout.merge_34(), Jt);
236  }
237 
238 #ifdef MFEM_TEMPLATE_ENABLE_SERIALIZE
239  inline MFEM_ALWAYS_INLINE
240  void EvalSerialized(int el, T_type &T, const vreal_t *nodeData)
241  {
242  const int SS = sizeof(nodeData[0])/sizeof(nodeData[0][0]);
243  MFEM_ASSERT(el % (SS*ne) == 0, "invalid element index: " << el);
244  T.evaluator.CalcGrad(nodes_dof_t::layout.merge_23(),
245  &nodeData[el/SS*nodes_dof_t::size],
246  Jt.layout.merge_34(), Jt);
247  }
248 #endif
249  };
250 
251  // Case EvalOps = 3 = EvalCoordinates|EvalJacobians
252  template <typename it_t> struct Result<3,it_t>
253  {
254  static const int ne = it_t::batch_size;
255  typedef typename it_t::vreal_t vreal_t;
258 #ifdef MFEM_TEMPLATE_ELTRANS_RESULT_HAS_NODES
260 #else
261  typedef TTensor4<qpts,dim,sdim,ne,vreal_t/*,true*/> Jt_type;
262 #endif
264 
266 #ifdef MFEM_TEMPLATE_ELTRANS_RESULT_HAS_NODES
268 #endif
269 
270  inline MFEM_ALWAYS_INLINE
271  void Eval(int el, T_type &T)
272  {
273 #ifdef MFEM_TEMPLATE_ELTRANS_HAS_NODE_DOFS
274  MFEM_STATIC_ASSERT(ne == 1, "only ne == 1 is supported");
276 #elif !defined(MFEM_TEMPLATE_ELTRANS_RESULT_HAS_NODES)
278 #endif
279  T.fes.SetElement(el);
280  T.fes.VectorExtract(T.node_layout, T.nodes,
282  T.evaluator.Calc(nodes_dof.layout.merge_23(), nodes_dof,
283  x.layout.merge_23(), x);
284  T.evaluator.CalcGrad(nodes_dof.layout.merge_23(), nodes_dof,
285  Jt.layout.merge_34(), Jt);
286  }
287 
288 #ifdef MFEM_TEMPLATE_ENABLE_SERIALIZE
289  inline MFEM_ALWAYS_INLINE
290  void EvalSerialized(int el, T_type &T, const vreal_t *nodeData)
291  {
292  const int SS = sizeof(nodeData[0])/sizeof(nodeData[0][0]);
293  MFEM_ASSERT(el % (SS*ne) == 0, "invalid element index: " << el);
294  T.evaluator.Calc(nodes_dof_t::layout.merge_23(),
295  &nodeData[el/SS*nodes_dof_t::size],
296  x.layout.merge_23(), x);
297  T.evaluator.CalcGrad(nodes_dof_t::layout.merge_23(),
298  &nodeData[el/SS*nodes_dof_t::size],
299  Jt.layout.merge_34(), Jt);
300  }
301 #endif
302  };
303 
304  // Case EvalOps = 6 = EvalJacobians|LoadAttributes
305  template <typename it_t> struct Result<6,it_t>
306  {
307  static const int ne = it_t::batch_size;
308  typedef typename it_t::vreal_t vreal_t;
309  typedef typename it_t::vint_t vint_t;
310 #ifdef MFEM_TEMPLATE_ELTRANS_RESULT_HAS_NODES
312 #else
313  typedef TTensor4<qpts,dim,sdim,ne,vreal_t/*,true*/> Jt_type;
314 #endif
316 
318 #ifdef MFEM_TEMPLATE_ELTRANS_RESULT_HAS_NODES
320 #endif
321  vint_t attrib[ne];
322 
323  inline MFEM_ALWAYS_INLINE
324  void Eval(int el, T_type &T)
325  {
326 #ifdef MFEM_TEMPLATE_ELTRANS_HAS_NODE_DOFS
327  MFEM_STATIC_ASSERT(ne == 1, "only ne == 1 is supported");
329 #elif !defined(MFEM_TEMPLATE_ELTRANS_RESULT_HAS_NODES)
331 #endif
332  T.fes.SetElement(el);
333  T.fes.VectorExtract(T.node_layout, T.nodes,
335  T.evaluator.CalcGrad(nodes_dof.layout.merge_23(), nodes_dof,
336  Jt.layout.merge_34(), Jt);
337  T.SetAttributes(el, attrib);
338  }
339 
340 #ifdef MFEM_TEMPLATE_ENABLE_SERIALIZE
341  inline MFEM_ALWAYS_INLINE
342  void EvalSerialized(int el, T_type &T, const vreal_t *nodeData)
343  {
344  const int SS = sizeof(nodeData[0])/sizeof(nodeData[0][0]);
345  MFEM_ASSERT(el % (SS*ne) == 0, "invalid element index: " << el);
346  T.evaluator.CalcGrad(nodes_dof_t::layout.merge_23(),
347  &nodeData[el/SS*nodes_dof_t::size],
348  Jt.layout.merge_34(), Jt);
349  T.SetAttributes(el, attrib);
350  }
351 #endif
352  };
353 
354  // Case EvalOps = 10 = EvalJacobians|LoadElementIdxs
355  template <typename it_t> struct Result<10,it_t>
356  {
357  static const int ne = it_t::batch_size;
358  typedef typename it_t::vreal_t vreal_t;
359 #ifdef MFEM_TEMPLATE_ELTRANS_RESULT_HAS_NODES
361 #else
362  typedef TTensor4<qpts,dim,sdim,ne,vreal_t/*,true*/> Jt_type;
363 #endif
365 
367 #ifdef MFEM_TEMPLATE_ELTRANS_RESULT_HAS_NODES
369 #endif
371 
372  inline MFEM_ALWAYS_INLINE
373  void Eval(int el, T_type &T)
374  {
375 #ifdef MFEM_TEMPLATE_ELTRANS_HAS_NODE_DOFS
376  MFEM_STATIC_ASSERT(ne == 1, "only ne == 1 is supported");
378 #elif !defined(MFEM_TEMPLATE_ELTRANS_RESULT_HAS_NODES)
380 #endif
381  T.fes.SetElement(el);
382  T.fes.VectorExtract(T.node_layout, T.nodes,
384  T.evaluator.CalcGrad(nodes_dof.layout.merge_23(), nodes_dof,
385  Jt.layout.merge_34(), Jt);
386  first_elem_idx = el;
387  }
388 
389 #ifdef MFEM_TEMPLATE_ENABLE_SERIALIZE
390  inline MFEM_ALWAYS_INLINE
391  void EvalSerialized(int el, T_type &T, const vreal_t *nodeData)
392  {
393  const int SS = sizeof(nodeData[0])/sizeof(nodeData[0][0]);
394  MFEM_ASSERT(el % (SS*ne) == 0, "invalid element index: " << el);
395  T.evaluator.CalcGrad(nodes_dof_t::layout.merge_23(),
396  &nodeData[el/SS*nodes_dof_t::size],
397  Jt.layout.merge_34(), Jt);
398  first_elem_idx = el;
399  }
400 #endif
401  };
402 };
403 
404 } // namespace mfem
405 
406 #endif // MFEM_TEMPLATE_ELEMENT_TRANSFORMATION
TTensor4< qpts, dim, sdim, ne, vreal_t > Jt_type
Definition: teltrans.hpp:213
TTensor3< qpts, sdim, ne, vreal_t > x_type
Definition: teltrans.hpp:167
TTensor4< qpts, dim, sdim, ne, vreal_t > Jt_type
Definition: teltrans.hpp:313
Mesh_t::nodeLayout_type nodeLayout_type
Definition: teltrans.hpp:39
TTensor3< dofs, sdim, ne, vreal_t > nodes_dof_t
Definition: teltrans.hpp:217
MFEM_ALWAYS_INLINE void Eval(int el, T_type &T)
Definition: teltrans.hpp:271
int GetAttribute() const
Return element&#39;s attribute.
Definition: element.hpp:55
TTensor3< dofs, sdim, ne, vreal_t > nodes_dof_t
Definition: teltrans.hpp:265
MFEM_ALWAYS_INLINE void Eval(int el, T_type &T)
Definition: teltrans.hpp:373
TTensor3< qpts, sdim, NE, vreal_t, true > x_type
Definition: teltrans.hpp:165
Templated struct Result, used to specify the type result that is computed by the TElementTransformati...
Definition: teltrans.hpp:75
static const int dofs
Definition: teltrans.hpp:79
MFEM_ALWAYS_INLINE void SetAttributes(int el, vint_t(&attrib)[NE]) const
Definition: teltrans.hpp:96
MFEM_ALWAYS_INLINE void Eval(int el, T_type &T)
Definition: teltrans.hpp:177
static const layout_type layout
Definition: ttensor.hpp:392
Element transformation class, templated on a mesh type and an integration rule. It is constructed fro...
Definition: teltrans.hpp:33
TElementTransformation(const Mesh_t &mesh, const ShapeEval &eval)
Definition: teltrans.hpp:110
MFEM_ALWAYS_INLINE void EvalSerialized(int el, T_type &T, const vreal_t *nodeData)
Definition: teltrans.hpp:342
static const layout_type layout
Definition: ttensor.hpp:414
static const int qpts
Definition: teltrans.hpp:80
TTensor3< dofs, sdim, ne, vreal_t > nodes_dof_t
Definition: teltrans.hpp:366
TTensor4< qpts, dim, sdim, ne, vreal_t, true > Jt_type
Definition: teltrans.hpp:311
MFEM_ALWAYS_INLINE void EvalSerialized(int el, T_type &T, const vreal_t *nodeData)
Definition: teltrans.hpp:391
TTensor3< dofs, sdim, 1, real_t > nodes_dof
Definition: teltrans.hpp:84
TTensor4< qpts, dim, sdim, ne, vreal_t, true > Jt_type
Definition: teltrans.hpp:211
TTensor4< qpts, dim, sdim, ne, vreal_t > Jt_type
Definition: teltrans.hpp:261
ShapeEvaluator< FE_type, IR, real_t > ShapeEval
Definition: teltrans.hpp:40
TTensor3< dofs, sdim, ne, vreal_t > nodes_dof_t
Definition: teltrans.hpp:317
MFEM_ALWAYS_INLINE void EvalSerialized(int el, T_type &T, const vreal_t *nodeData)
Definition: teltrans.hpp:155
MFEM_ALWAYS_INLINE void Eval(int el, T_type &T)
Definition: teltrans.hpp:324
const Element *const * elements
Definition: teltrans.hpp:92
static const int dim
Definition: teltrans.hpp:77
TTensor3< dofs, sdim, ne, vreal_t > nodes_dof_t
Definition: teltrans.hpp:171
MFEM_ALWAYS_INLINE void Eval(int el, Result< EvalOps, impl_traits_t > &F)
Evaluate coordinates and/or Jacobian matrices at quadrature points.
Definition: teltrans.hpp:121
int dim
Definition: ex24.cpp:53
TTensor3< qpts, sdim, ne, vreal_t, true > x_type
Definition: teltrans.hpp:256
nodeLayout_type node_layout
Definition: teltrans.hpp:89
MFEM_ALWAYS_INLINE void EvalSerialized(int el, T_type &T, const vreal_t *nodeData)
Definition: teltrans.hpp:290
TTensor4< qpts, dim, sdim, ne, vreal_t > Jt_type
Definition: teltrans.hpp:362
MFEM_ALWAYS_INLINE void Eval(int el, T_type &T)
Definition: teltrans.hpp:223
TTensor4< qpts, dim, sdim, ne, vreal_t, true > Jt_type
Definition: teltrans.hpp:360
MFEM_ALWAYS_INLINE void EvalSerialized(int el, const typename impl_traits_t::vreal_t *nodeData, Result< EvalOps, impl_traits_t > &F)
Definition: teltrans.hpp:129
TElementTransformation< Mesh_t, IR, real_t > T_type
Definition: teltrans.hpp:42
TTensor4< qpts, dim, sdim, ne, vreal_t, true > Jt_type
Definition: teltrans.hpp:259
Abstract data type element.
Definition: element.hpp:28
Mesh_t::FE_type FE_type
Definition: teltrans.hpp:37
MFEM_ALWAYS_INLINE void EvalSerialized(int el, T_type &T, const vreal_t *nodeData)
Definition: teltrans.hpp:240
MFEM_ALWAYS_INLINE void Eval(int el, T_type &T)
Definition: teltrans.hpp:148
static const int sdim
Definition: teltrans.hpp:78
Mesh_t::FESpace_type FESpace_type
Definition: teltrans.hpp:38
MFEM_ALWAYS_INLINE void EvalSerialized(int el, T_type &T, const vreal_t *nodeData)
Definition: teltrans.hpp:194