https://mooseframework.inl.gov
KokkosDatum.h
Go to the documentation of this file.
1 //* This file is part of the MOOSE framework
2 //* https://www.mooseframework.org
3 //*
4 //* All rights reserved, see COPYRIGHT for full restrictions
5 //* https://github.com/idaholab/moose/blob/master/COPYRIGHT
6 //*
7 //* Licensed under LGPL 2.1, please see LICENSE for details
8 //* https://www.gnu.org/licenses/lgpl-2.1.html
9 
10 #pragma once
11 
12 #include "KokkosTypes.h"
13 #include "KokkosAssembly.h"
14 #include "KokkosSystem.h"
15 #include "KokkosVariable.h"
16 
17 namespace Moose
18 {
19 namespace Kokkos
20 {
21 
25 class Datum
26 {
27 public:
35  KOKKOS_FUNCTION
37  const unsigned int side,
38  const Assembly & assembly,
39  const Array<System> & systems)
41  _systems(systems),
42  _elem(assembly.kokkosMesh().getElementInfo(elem)),
43  _side(side),
44  _neighbor(_side == libMesh::invalid_uint ? libMesh::DofObject::invalid_id
45  : assembly.kokkosMesh().getNeighbor(_elem.id, side)),
46  _n_qps(side == libMesh::invalid_uint ? assembly.getNumQps(_elem)
47  : assembly.getNumFaceQps(_elem, side)),
49  : assembly.getQpFaceOffset(_elem, side))
50  {
51  }
58  KOKKOS_FUNCTION
59  Datum(const ContiguousNodeID node, const Assembly & assembly, const Array<System> & systems)
60  : _assembly(assembly), _systems(systems), _node(node)
61  {
62  }
63 
68  KOKKOS_FUNCTION const Assembly & assembly() const { return _assembly; }
74  KOKKOS_FUNCTION const System & system(unsigned int sys) const { return _systems[sys]; }
75 
80  KOKKOS_FUNCTION const ElementInfo & elem() const { return _elem; }
85  KOKKOS_FUNCTION ContiguousSubdomainID subdomain() const { return _elem.subdomain; }
90  KOKKOS_FUNCTION unsigned int side() const { return _side; }
95  KOKKOS_FUNCTION ContiguousNodeID node() const { return _node; }
100  KOKKOS_FUNCTION unsigned int n_qps() const { return _n_qps; }
105  KOKKOS_FUNCTION dof_id_type qpOffset() const { return _qp_offset; }
110  KOKKOS_FUNCTION bool hasNeighbor() const { return _neighbor != libMesh::DofObject::invalid_id; }
115  KOKKOS_FUNCTION bool isNodal() const { return _node != libMesh::DofObject::invalid_id; }
116 
125  KOKKOS_FUNCTION const Real33 & J(const unsigned int qp)
126  {
127  if (!isNodal())
128  reinitTransform(qp);
129  else
131 
132  return _J;
133  }
139  KOKKOS_FUNCTION Real JxW(const unsigned int qp)
140  {
141  if (!isNodal())
142  reinitTransform(qp);
143  else
144  _JxW = 1;
145 
146  return _JxW;
147  }
153  KOKKOS_FUNCTION Real3 q_point(const unsigned int qp)
154  {
155  if (!isNodal())
156  reinitTransform(qp);
157  else
159 
160  return _xyz;
161  }
162 
166  KOKKOS_FUNCTION void reinit() { _transform_reinit = false; }
167 
168 protected:
184  const unsigned int _side = libMesh::invalid_uint;
196  const unsigned int _n_qps = 1;
201 
202 private:
207  KOKKOS_FUNCTION void reinitTransform(const unsigned int qp);
208 
212  bool _transform_reinit = false;
216  Real33 _J;
221 };
222 
223 KOKKOS_FUNCTION inline void
224 Datum::reinitTransform(const unsigned int qp)
225 {
226  if (_transform_reinit)
227  return;
228 
230  {
231  _J = _assembly.getJacobian(_elem, qp);
232  _JxW = _assembly.getJxW(_elem, qp);
233  _xyz = _assembly.getQPoint(_elem, qp);
234  }
235  else
237 
238  _transform_reinit = true;
239 }
240 
244 class AssemblyDatum : public Datum
245 {
246 public:
257  KOKKOS_FUNCTION
259  const unsigned int side,
260  const Assembly & assembly,
261  const Array<System> & systems,
262  const Variable & ivar,
263  const unsigned int jvar,
264  const unsigned int comp = 0)
265  : Datum(elem, side, assembly, systems),
266  _tag(ivar.tag()),
267  _ivar(ivar.var(comp)),
268  _jvar(jvar),
269  _ife(systems[ivar.sys(comp)].getFETypeID(_ivar)),
270  _jfe(systems[ivar.sys(comp)].getFETypeID(_jvar)),
271  _n_idofs(assembly.getNumDofs(_elem.type, _ife)),
272  _n_jdofs(assembly.getNumDofs(_elem.type, _jfe))
273  {
274  }
284  KOKKOS_FUNCTION
286  const Assembly & assembly,
287  const Array<System> & systems,
288  const Variable & ivar,
289  const unsigned int jvar,
290  const unsigned int comp = 0)
291  : Datum(node, assembly, systems),
292  _tag(ivar.tag()),
293  _ivar(ivar.var(comp)),
294  _jvar(jvar),
295  _ife(systems[ivar.sys(comp)].getFETypeID(_ivar)),
296  _jfe(systems[ivar.sys(comp)].getFETypeID(_jvar))
297  {
298  }
299 
304  KOKKOS_FUNCTION unsigned int n_dofs() const { return _n_idofs; }
309  KOKKOS_FUNCTION unsigned int n_idofs() const { return _n_idofs; }
314  KOKKOS_FUNCTION unsigned int n_jdofs() const { return _n_jdofs; }
319  KOKKOS_FUNCTION unsigned int var() const { return _ivar; }
324  KOKKOS_FUNCTION unsigned int ivar() const { return _ivar; }
329  KOKKOS_FUNCTION unsigned int jvar() const { return _jvar; }
334  KOKKOS_FUNCTION unsigned int fe() const { return _ife; }
339  KOKKOS_FUNCTION unsigned int ife() const { return _ife; }
344  KOKKOS_FUNCTION unsigned int jfe() const { return _jfe; }
345 
346 protected:
350  const TagID _tag;
354  const unsigned int _ivar, _jvar;
358  const unsigned int _ife, _jfe;
362  const unsigned int _n_idofs = 1, _n_jdofs = 1;
363 };
364 
365 } // namespace Kokkos
366 } // namespace Moose
367 
The Kokkos array class.
Definition: KokkosArray.h:56
The Kokkos assembly class.
The Kokkos object that contains the information of an element The IDs used in Kokkos are different fr...
Definition: KokkosMesh.h:35
The Kokkos object that holds thread-private data in the parallel operations of any Kokkos object...
Definition: KokkosDatum.h:25
Moose::Kokkos::AssemblyDatum AssemblyDatum
Definition: KokkosDatum.h:369
const unsigned int invalid_uint
const unsigned int _n_jdofs
Definition: KokkosDatum.h:362
const Assembly & _assembly
Reference of the Kokkos assembly.
Definition: KokkosDatum.h:172
unsigned int TagID
Definition: MooseTypes.h:210
KOKKOS_INLINE_FUNCTION void identity(const unsigned int dim=3)
Definition: KokkosTypes.h:58
dof_id_type ContiguousElementID
Definition: KokkosMesh.h:20
const unsigned int _jvar
Definition: KokkosDatum.h:354
Moose::Kokkos::Datum Datum
Definition: KokkosDatum.h:368
KOKKOS_FUNCTION unsigned int getDimension() const
Get the mesh dimension.
const unsigned int _side
Current side index.
Definition: KokkosDatum.h:184
KOKKOS_FUNCTION unsigned int jvar() const
Get the coupled variable number.
Definition: KokkosDatum.h:329
KOKKOS_FUNCTION Real3 q_point(const unsigned int qp)
Get the physical quadrature point coordinate.
Definition: KokkosDatum.h:153
const unsigned int _n_qps
Number of local quadrature points.
Definition: KokkosDatum.h:196
The following methods are specializations for using the libMesh::Parallel::packed_range_* routines fo...
KOKKOS_FUNCTION unsigned int n_jdofs() const
Get the number of local DOFs for the coupled variable.
Definition: KokkosDatum.h:314
The Kokkos system class.
Definition: KokkosSystem.h:35
KOKKOS_FUNCTION unsigned int ivar() const
Get the variable number.
Definition: KokkosDatum.h:324
const ElementInfo _elem
Current element information object.
Definition: KokkosDatum.h:180
KOKKOS_FUNCTION const Real33 & J(const unsigned int qp)
Get the inverse of Jacobian matrix | dxi/dx deta/dx dzeta/dx | | dxi/dy deta/dy dzeta/dy | | dxi/dz d...
Definition: KokkosDatum.h:125
KOKKOS_FUNCTION unsigned int fe() const
Get the variable FE type ID.
Definition: KokkosDatum.h:334
KOKKOS_FUNCTION unsigned int side() const
Get the side index.
Definition: KokkosDatum.h:90
const unsigned int _jfe
Definition: KokkosDatum.h:358
KOKKOS_FUNCTION unsigned int jfe() const
Get the coupled variable FE type ID.
Definition: KokkosDatum.h:344
static constexpr dof_id_type invalid_id
dof_id_type ContiguousNodeID
Definition: KokkosMesh.h:21
KOKKOS_FUNCTION unsigned int n_idofs() const
Get the number of local DOFs.
Definition: KokkosDatum.h:309
const unsigned int _n_idofs
Number of local DOFs.
Definition: KokkosDatum.h:362
KOKKOS_FUNCTION unsigned int n_qps() const
Get the number of local quadrature points.
Definition: KokkosDatum.h:100
KOKKOS_FUNCTION AssemblyDatum(const ContiguousElementID elem, const unsigned int side, const Assembly &assembly, const Array< System > &systems, const Variable &ivar, const unsigned int jvar, const unsigned int comp=0)
Constructor for element and side data.
Definition: KokkosDatum.h:258
KOKKOS_FUNCTION const Assembly & assembly() const
Get the Kokkos assembly.
Definition: KokkosDatum.h:68
const TagID _tag
Solution tag ID.
Definition: KokkosDatum.h:350
KOKKOS_FUNCTION ContiguousSubdomainID subdomain() const
Get the contiguous subdomain ID.
Definition: KokkosDatum.h:85
KOKKOS_FUNCTION AssemblyDatum(const ContiguousNodeID node, const Assembly &assembly, const Array< System > &systems, const Variable &ivar, const unsigned int jvar, const unsigned int comp=0)
Constructor for node data.
Definition: KokkosDatum.h:285
KOKKOS_FUNCTION void reinitTransform(const unsigned int qp)
Compute and cache the physical transformation data.
Definition: KokkosDatum.h:224
KOKKOS_FUNCTION unsigned int var() const
Get the variable number.
Definition: KokkosDatum.h:319
KOKKOS_FUNCTION Real3 getNodePoint(ContiguousNodeID node) const
Get the coordinate of a node.
Definition: KokkosMesh.h:215
KOKKOS_FUNCTION Real getJxW(ElementInfo info, unsigned int qp) const
Get the transformed Jacobian weight of an element quadrature point.
KOKKOS_FUNCTION const Mesh & kokkosMesh() const
Get the const reference of the Kokkos mesh.
Definition: KokkosMesh.h:360
ContiguousSubdomainID subdomain
Contiguous subdomain ID.
Definition: KokkosMesh.h:48
KOKKOS_FUNCTION ContiguousNodeID node() const
Get the contiguous node ID.
Definition: KokkosDatum.h:95
const dof_id_type _qp_offset
Starting offset into the global quadrature point index.
Definition: KokkosDatum.h:200
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
bool _transform_reinit
Flag whether the physical transformation data was cached.
Definition: KokkosDatum.h:212
KOKKOS_FUNCTION bool isNodal() const
Get whether the current datum is on a node.
Definition: KokkosDatum.h:115
The Kokkos variable object that carries the coupled variable and tag information. ...
KOKKOS_FUNCTION const ElementInfo & elem() const
Get the element information object.
Definition: KokkosDatum.h:80
KOKKOS_FUNCTION unsigned int ife() const
Get the variable FE type ID.
Definition: KokkosDatum.h:339
The Kokkos object that holds thread-private data in the parallel operations of Kokkos kernels...
Definition: KokkosDatum.h:244
KOKKOS_FUNCTION Datum(const ContiguousNodeID node, const Assembly &assembly, const Array< System > &systems)
Constructor for node data.
Definition: KokkosDatum.h:59
const ContiguousElementID _neighbor
Current contiguous element ID of neighbor.
Definition: KokkosDatum.h:192
KOKKOS_FUNCTION Real JxW(const unsigned int qp)
Get the transformed Jacobian weight.
Definition: KokkosDatum.h:139
KOKKOS_FUNCTION Datum(const ContiguousElementID elem, const unsigned int side, const Assembly &assembly, const Array< System > &systems)
Constructor for element and side data.
Definition: KokkosDatum.h:36
MOOSE now contains C++17 code, so give a reasonable error message stating what the user can do to add...
KOKKOS_FUNCTION unsigned int n_dofs() const
Get the number of local DOFs.
Definition: KokkosDatum.h:304
KOKKOS_FUNCTION const System & system(unsigned int sys) const
Get the Kokkos system.
Definition: KokkosDatum.h:74
KOKKOS_FUNCTION dof_id_type qpOffset() const
Get the starting offset into the global quadrature point index.
Definition: KokkosDatum.h:105
const Array< System > & _systems
Reference of the Kokkos systems.
Definition: KokkosDatum.h:176
KOKKOS_FUNCTION void reinit()
Reset the reinit flag.
Definition: KokkosDatum.h:166
const ContiguousNodeID _node
Current contiguous node ID.
Definition: KokkosDatum.h:188
const unsigned int _ivar
Variable numbers.
Definition: KokkosDatum.h:354
KOKKOS_FUNCTION Real33 getJacobian(ElementInfo info, unsigned int qp) const
Get the inverse of Jacobian matrix of an element quadrature point.
const unsigned int _ife
FE type IDs of variables.
Definition: KokkosDatum.h:358
KOKKOS_FUNCTION bool hasNeighbor() const
Get whether the current side has a neighbor.
Definition: KokkosDatum.h:110
uint8_t dof_id_type
Real33 _J
Cached physical transformation data.
Definition: KokkosDatum.h:217
KOKKOS_FUNCTION void computePhysicalMap(const ElementInfo info, const unsigned int qp, Real33 *const jacobian, Real *const JxW, Real3 *const q_points) const
Compute physical transformation data for an element.
KOKKOS_FUNCTION Real3 getQPoint(ElementInfo info, unsigned int qp) const
Get the coordinate of an element quadrature point.