https://mooseframework.inl.gov
Public Member Functions | Private Member Functions | Private Attributes | List of all members
Moose::Kokkos::Assembly Class Reference

The Kokkos assembly class. More...

#include <KokkosAssembly.h>

Inheritance diagram for Moose::Kokkos::Assembly:
[legend]

Public Member Functions

 Assembly (FEProblemBase &problem)
 Constructor. More...
 
void init ()
 Initialize assembly. More...
 
unsigned int getFETypeID (FEType type) const
 Get the FE type ID. More...
 
KOKKOS_FUNCTION unsigned int getMaxQpsPerElem () const
 Get the maximum number of quadrature points per element in the current partition. More...
 
KOKKOS_FUNCTION dof_id_type getNumQps (ContiguousSubdomainID subdomain) const
 Get the total number of elemental quadrature points in a subdomain. More...
 
KOKKOS_FUNCTION unsigned int getNumQps (ElementInfo info) const
 Get the number of quadrature points of an element. More...
 
KOKKOS_FUNCTION dof_id_type getNumFaceQps (ContiguousSubdomainID subdomain) const
 Get the total number of facial quadrature points in a subdomain Note: this number does not represent the real number of facial quadrature points but only counts the facial quadrature points that need global caching, such as face material properties. More...
 
KOKKOS_FUNCTION unsigned int getNumFaceQps (ElementInfo info, unsigned int side) const
 Get the number of quadrature points of a side of an element. More...
 
KOKKOS_FUNCTION dof_id_type getQpOffset (ElementInfo info) const
 Get the starting offset of quadrature points of an element into the global quadrature point index. More...
 
KOKKOS_FUNCTION dof_id_type getQpFaceOffset (ElementInfo info, unsigned int side) const
 Get the starting offset of quadrature points of a side of an element into the global quadrature point index. More...
 
KOKKOS_FUNCTION unsigned int getNumDofs (unsigned int elem_type, unsigned int fe_type) const
 Get the number of DOFs of a FE type for an element type. More...
 
KOKKOS_FUNCTION const auto & getPhi (ContiguousSubdomainID subdomain, unsigned int elem_type, unsigned int fe_type) const
 Get the shape functions of a FE type for an element type and subdomain. More...
 
KOKKOS_FUNCTION const auto & getPhiFace (ContiguousSubdomainID subdomain, unsigned int elem_type, unsigned int fe_type) const
 Get the face shape functions of a FE type for an element type and subdomain. More...
 
KOKKOS_FUNCTION const auto & getGradPhi (ContiguousSubdomainID subdomain, unsigned int elem_type, unsigned int fe_type) const
 Get the gradient of shape functions of a FE type for an element type and subdomain. More...
 
KOKKOS_FUNCTION const auto & getGradPhiFace (ContiguousSubdomainID subdomain, unsigned int elem_type, unsigned int fe_type) const
 Get the gradient of face shape functions of a FE type for an element type and subdomain. More...
 
KOKKOS_FUNCTION Real33 getJacobian (ElementInfo info, unsigned int qp) const
 Get the inverse of Jacobian matrix of an element quadrature point. More...
 
KOKKOS_FUNCTION Real getJxW (ElementInfo info, unsigned int qp) const
 Get the transformed Jacobian weight of an element quadrature point. More...
 
KOKKOS_FUNCTION Real3 getQPoint (ElementInfo info, unsigned int qp) const
 Get the coordinate of an element quadrature point. More...
 
KOKKOS_FUNCTION Real coordTransformFactor (const ContiguousSubdomainID subdomain, const Real3 point) const
 Get the coordinate transform factor for a point in a subdomain. More...
 
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. More...
 
KOKKOS_FUNCTION void computePhysicalMap (const ElementInfo info, const unsigned int side, const unsigned int qp, Real33 *const jacobian, Real *const JxW, Real3 *const q_points) const
 Compute physical transformation data for a side. More...
 
KOKKOS_FUNCTION void operator() (const ThreadID tid) const
 Kokkos function for caching physical maps on element quadrature points. More...
 
const auto & getMaterialBoundaries () const
 Get the list of boundaries to cache face material properties. More...
 
KOKKOS_FUNCTION const MeshkokkosMesh () const
 Get the const reference of the Kokkos mesh. More...
 

Private Member Functions

void initQuadrature ()
 Initialize quadrature data. More...
 
void initShape ()
 Initialize shape data. More...
 
void cachePhysicalMap ()
 Cache physical maps on element quadrature points. More...
 

Private Attributes

FEProblemBase_problem
 Reference of the MOOSE problem. More...
 
MooseMesh_mesh
 Reference of the MOOSE mesh. More...
 
std::map< FEType, unsigned int_fe_type_map
 FE type ID map. More...
 
const unsigned int _dimension
 Mesh dimension. More...
 
Array< Moose::CoordinateSystemType_coord_type
 Coordinate system type of each subdomain. More...
 
unsigned int _rz_radial_coord = libMesh::invalid_uint
 Radial coordinate index in cylindrical coordinate system. More...
 
Array< Pair< Real3, Real3 > > _rz_axis
 General axisymmetric axis of each subdomain in cylindrical coordinate system. More...
 
std::set< BoundaryID_material_boundaries
 Boundaries to cache face material properties. More...
 
Array< dof_id_type_qp_offset
 Starting offset into the global quadrature point index. More...
 
Array2D< dof_id_type_qp_offset_face
 
Array< unsigned int_n_qps
 Number of quadrature points. More...
 
Array2D< unsigned int_n_qps_face
 
unsigned int _max_qps_per_elem = 0
 
Array< dof_id_type_n_subdomain_qps
 
Array< dof_id_type_n_subdomain_qps_face
 
Array2D< Array< Real3 > > _q_points
 Quadrature points and weights for reference elements. More...
 
Array2D< Array< Array< Real3 > > > _q_points_face
 
Array2D< Array< Real > > _weights
 
Array2D< Array< Array< Real > > > _weights_face
 
Array3D< Array2D< Real > > _phi
 Shape functions for reference elements. More...
 
Array3D< Array< Array2D< Real > > > _phi_face
 
Array3D< Array2D< Real3 > > _grad_phi
 
Array3D< Array< Array2D< Real3 > > > _grad_phi_face
 
Array2D< unsigned int_n_dofs
 
Array2D< Array2D< Real > > _map_phi
 Shape functions for computing reference-to-physical maps. More...
 
Array2D< Array< Array2D< Real > > > _map_phi_face
 
Array2D< Array< Array2D< Real > > > _map_psi_face
 
Array2D< Array2D< Real3 > > _map_grad_phi
 
Array2D< Array< Array2D< Real3 > > > _map_grad_phi_face
 
Array2D< Array< Array2D< Real3 > > > _map_grad_psi_face
 
Array< Array< Real33 > > _jacobian
 Cached physical maps on element quadrature points. More...
 
Array< Array< Real > > _jxw
 
Array< Array< Real3 > > _xyz
 

Detailed Description

The Kokkos assembly class.

Definition at line 30 of file KokkosAssembly.h.

Constructor & Destructor Documentation

◆ Assembly()

Moose::Kokkos::Assembly::Assembly ( FEProblemBase problem)

Constructor.

Parameters
problemThe MOOSE problem

Member Function Documentation

◆ cachePhysicalMap()

void Moose::Kokkos::Assembly::cachePhysicalMap ( )
private

Cache physical maps on element quadrature points.

◆ computePhysicalMap() [1/2]

KOKKOS_FUNCTION void Assembly::computePhysicalMap ( const ElementInfo  info,
const unsigned int  qp,
Real33 *const  jacobian,
Real *const  JxW,
Real3 *const  q_points 
) const
inline

Compute physical transformation data for an element.

Parameters
infoThe element information object
qpThe local quadrature point index
jacobianThe pointer to store the inverse of Jacobian matrix
JxWThe pointer to store transformed Jacobian weight
q_pointsThe pointer to store physical quadrature point coordinate

Definition at line 378 of file KokkosAssembly.h.

Referenced by Moose::Kokkos::Datum::reinitTransform().

383 {
384  auto sid = info.subdomain;
385  auto eid = info.id;
386  auto elem_type = info.type;
387  auto num_nodes = kokkosMesh().getNumNodes(elem_type);
388 
389  auto & phi = _map_phi(sid, elem_type);
390  auto & grad_phi = _map_grad_phi(sid, elem_type);
391 
392  Real33 J;
393  Real3 xyz;
394 
395  for (unsigned int node = 0; node < num_nodes; ++node)
396  {
397  auto points = kokkosMesh().getNodePoint(kokkosMesh().getContiguousNodeID(eid, node));
398 
399  if (jacobian || JxW)
400  J += grad_phi(node, qp).cartesian_product(points);
401 
402  xyz += phi(node, qp) * points;
403  }
404 
405  if (jacobian)
406  *jacobian = J.inverse(_dimension);
407  if (JxW)
408  *JxW =
409  J.determinant(_dimension) * _weights(sid, elem_type)[qp] * coordTransformFactor(sid, xyz);
410  if (q_points)
411  *q_points = xyz;
412 }
Array2D< Array< Real > > _weights
MPI_Info info
const unsigned int _dimension
Mesh dimension.
KOKKOS_FUNCTION Real coordTransformFactor(const ContiguousSubdomainID subdomain, const Real3 point) const
Get the coordinate transform factor for a point in a subdomain.
Array2D< Array2D< Real > > _map_phi
Shape functions for computing reference-to-physical maps.
KOKKOS_FUNCTION unsigned int getNumNodes(unsigned int elem_type) const
Get the number of nodes of an element type.
Definition: KokkosMesh.h:172
KOKKOS_INLINE_FUNCTION Real33 cartesian_product(const Real3 vector)
Definition: KokkosTypes.h:227
KOKKOS_FUNCTION Real3 getNodePoint(ContiguousNodeID node) const
Get the coordinate of a node.
Definition: KokkosMesh.h:215
KOKKOS_FUNCTION const Mesh & kokkosMesh() const
Get the const reference of the Kokkos mesh.
Definition: KokkosMesh.h:360
KOKKOS_INLINE_FUNCTION Real determinant(const unsigned int dim=3)
Definition: KokkosTypes.h:58
Array2D< Array2D< Real3 > > _map_grad_phi
KOKKOS_INLINE_FUNCTION Real33 inverse(const unsigned int dim=3)
Definition: KokkosTypes.h:75

◆ computePhysicalMap() [2/2]

KOKKOS_FUNCTION void Assembly::computePhysicalMap ( const ElementInfo  info,
const unsigned int  side,
const unsigned int  qp,
Real33 *const  jacobian,
Real *const  JxW,
Real3 *const  q_points 
) const
inline

Compute physical transformation data for a side.

Parameters
infoThe element information object
sideThe side index
qpThe local quadrature point index
jacobianThe pointer to store the inverse of Jacobian matrix
JxWThe pointer to store transformed Jacobian weight
q_pointsThe pointer to store physical quadrature point coordinate

Definition at line 415 of file KokkosAssembly.h.

421 {
422  auto sid = info.subdomain;
423  auto eid = info.id;
424  auto elem_type = info.type;
425  auto num_nodes = kokkosMesh().getNumNodes(elem_type);
426  auto num_side_nodes = kokkosMesh().getNumNodes(elem_type, side);
427 
428  auto & phi = _map_phi_face(sid, elem_type)(side);
429  auto & grad_phi = _map_grad_phi_face(sid, elem_type)(side);
430 
431  Real33 J;
432  Real3 xyz;
433 
434  for (unsigned int node = 0; node < num_nodes; ++node)
435  {
436  auto points = kokkosMesh().getNodePoint(kokkosMesh().getContiguousNodeID(eid, node));
437 
438  if (jacobian)
439  J += grad_phi(node, qp).cartesian_product(points);
440 
441  if (JxW || q_points)
442  xyz += phi(node, qp) * points;
443  }
444 
445  if (jacobian)
446  *jacobian = J.inverse(_dimension);
447  if (q_points)
448  *q_points = xyz;
449 
450  if (JxW)
451  {
452  Real33 J;
453 
454  auto & grad_psi = _map_grad_psi_face(sid, elem_type)(side);
455 
456  for (unsigned int node = 0; node < num_side_nodes; ++node)
457  {
458  auto points = kokkosMesh().getNodePoint(kokkosMesh().getContiguousNodeID(eid, side, node));
459 
460  J += grad_psi(node, qp).cartesian_product(points);
461  }
462 
463  *JxW = ::Kokkos::sqrt((J * J.transpose()).determinant(_dimension - 1)) *
464  _weights_face(sid, elem_type)[side][qp] * coordTransformFactor(sid, xyz);
465  }
466 }
Array2D< Array< Array2D< Real3 > > > _map_grad_psi_face
Array2D< Array< Array2D< Real > > > _map_phi_face
MPI_Info info
const unsigned int _dimension
Mesh dimension.
Array2D< Array< Array< Real > > > _weights_face
KOKKOS_FUNCTION Real coordTransformFactor(const ContiguousSubdomainID subdomain, const Real3 point) const
Get the coordinate transform factor for a point in a subdomain.
Array2D< Array< Array2D< Real3 > > > _map_grad_phi_face
KOKKOS_FUNCTION unsigned int getNumNodes(unsigned int elem_type) const
Get the number of nodes of an element type.
Definition: KokkosMesh.h:172
KOKKOS_INLINE_FUNCTION Real33 cartesian_product(const Real3 vector)
Definition: KokkosTypes.h:227
KOKKOS_FUNCTION Real3 getNodePoint(ContiguousNodeID node) const
Get the coordinate of a node.
Definition: KokkosMesh.h:215
KOKKOS_FUNCTION const Mesh & kokkosMesh() const
Get the const reference of the Kokkos mesh.
Definition: KokkosMesh.h:360
CTSub CT_OPERATOR_BINARY CTMul CTCompareLess CTCompareGreater CTCompareEqual _arg template * sqrt(_arg)) *_arg.template D< dtag >()) CT_SIMPLE_UNARY_FUNCTION(tanh
KOKKOS_INLINE_FUNCTION Real33 inverse(const unsigned int dim=3)
Definition: KokkosTypes.h:75

◆ coordTransformFactor()

KOKKOS_FUNCTION Real Assembly::coordTransformFactor ( const ContiguousSubdomainID  subdomain,
const Real3  point 
) const
inline

Get the coordinate transform factor for a point in a subdomain.

Parameters
subdomainThe contiguous subdomain ID
pointThe point coordinate
Returns
The coordinate transform factor

Definition at line 358 of file KokkosAssembly.h.

Referenced by computePhysicalMap().

359 {
360  switch (_coord_type[subdomain])
361  {
362  case Moose::COORD_XYZ:
363  return 1;
364  case Moose::COORD_RZ:
366  return 2 * M_PI *
367  (point - _rz_axis[subdomain].first).cross_product(_rz_axis[subdomain].second).norm();
368  else
369  return 2 * M_PI * point(_rz_radial_coord);
371  return 4 * M_PI * point(0) * point(0);
372  default:
373  return 0;
374  }
375 }
unsigned int _rz_radial_coord
Radial coordinate index in cylindrical coordinate system.
const unsigned int invalid_uint
Array< Pair< Real3, Real3 > > _rz_axis
General axisymmetric axis of each subdomain in cylindrical coordinate system.
Array< Moose::CoordinateSystemType > _coord_type
Coordinate system type of each subdomain.

◆ getFETypeID()

unsigned int Moose::Kokkos::Assembly::getFETypeID ( FEType  type) const
inline

Get the FE type ID.

Parameters
typeThe libMesh FEType object
Returns
The FE type ID

Definition at line 49 of file KokkosAssembly.h.

49 { return libmesh_map_find(_fe_type_map, type); }
std::map< FEType, unsigned int > _fe_type_map
FE type ID map.

◆ getGradPhi()

KOKKOS_FUNCTION const auto& Moose::Kokkos::Assembly::getGradPhi ( ContiguousSubdomainID  subdomain,
unsigned int  elem_type,
unsigned int  fe_type 
) const
inline

Get the gradient of shape functions of a FE type for an element type and subdomain.

Parameters
subdomainThe contiguous subdomain ID
elem_typeThe element type ID
fe_typeThe FE type ID
Returns
The gradient of shape functions at quadrature points

Definition at line 151 of file KokkosAssembly.h.

Referenced by Moose::Kokkos::VariablePhiGradient::operator()(), and Moose::Kokkos::VariableTestGradient::operator()().

152  {
153  return _grad_phi(subdomain, elem_type, fe_type);
154  }
Array3D< Array2D< Real3 > > _grad_phi

◆ getGradPhiFace()

KOKKOS_FUNCTION const auto& Moose::Kokkos::Assembly::getGradPhiFace ( ContiguousSubdomainID  subdomain,
unsigned int  elem_type,
unsigned int  fe_type 
) const
inline

Get the gradient of face shape functions of a FE type for an element type and subdomain.

Parameters
subdomainThe contiguous subdomain ID
elem_typeThe element type ID
fe_typeThe FE type ID
Returns
The gradient of shape functions of all sides at quadrature points

Definition at line 162 of file KokkosAssembly.h.

Referenced by Moose::Kokkos::System::getVectorQpGradFace(), Moose::Kokkos::VariablePhiGradient::operator()(), and Moose::Kokkos::VariableTestGradient::operator()().

165  {
166  return _grad_phi_face(subdomain, elem_type, fe_type);
167  }
Array3D< Array< Array2D< Real3 > > > _grad_phi_face

◆ getJacobian()

KOKKOS_FUNCTION Real33 Moose::Kokkos::Assembly::getJacobian ( ElementInfo  info,
unsigned int  qp 
) const
inline

Get the inverse of Jacobian matrix of an element quadrature point.

Parameters
infoThe element information object
qpThe local quadrature point index
Returns
The inverse of Jacobian matrix

Definition at line 174 of file KokkosAssembly.h.

Referenced by Moose::Kokkos::Datum::reinitTransform().

175  {
176  return _jacobian[info.subdomain][getQpOffset(info) + qp];
177  }
KOKKOS_FUNCTION dof_id_type getQpOffset(ElementInfo info) const
Get the starting offset of quadrature points of an element into the global quadrature point index...
MPI_Info info
Array< Array< Real33 > > _jacobian
Cached physical maps on element quadrature points.

◆ getJxW()

KOKKOS_FUNCTION Real Moose::Kokkos::Assembly::getJxW ( ElementInfo  info,
unsigned int  qp 
) const
inline

Get the transformed Jacobian weight of an element quadrature point.

Parameters
infoThe element information object
qpThe local quadrature point index
Returns
The inverse of Jacobian matrix

Definition at line 184 of file KokkosAssembly.h.

Referenced by Moose::Kokkos::Datum::reinitTransform().

185  {
186  return _jxw[info.subdomain][getQpOffset(info) + qp];
187  }
KOKKOS_FUNCTION dof_id_type getQpOffset(ElementInfo info) const
Get the starting offset of quadrature points of an element into the global quadrature point index...
MPI_Info info
Array< Array< Real > > _jxw

◆ getMaterialBoundaries()

const auto& Moose::Kokkos::Assembly::getMaterialBoundaries ( ) const
inline

Get the list of boundaries to cache face material properties.

Returns
The list of boundaries

Definition at line 245 of file KokkosAssembly.h.

245 { return _material_boundaries; }
std::set< BoundaryID > _material_boundaries
Boundaries to cache face material properties.

◆ getMaxQpsPerElem()

KOKKOS_FUNCTION unsigned int Moose::Kokkos::Assembly::getMaxQpsPerElem ( ) const
inline

Get the maximum number of quadrature points per element in the current partition.

Returns
The maximum number of quadrature points per element

Definition at line 54 of file KokkosAssembly.h.

54 { return _max_qps_per_elem; }
unsigned int _max_qps_per_elem

◆ getNumDofs()

KOKKOS_FUNCTION unsigned int Moose::Kokkos::Assembly::getNumDofs ( unsigned int  elem_type,
unsigned int  fe_type 
) const
inline

Get the number of DOFs of a FE type for an element type.

Parameters
elem_typeThe element type ID
fe_typeThe FE type ID
Returns
The number of DOFs

Definition at line 115 of file KokkosAssembly.h.

Referenced by Moose::Kokkos::System::getVectorQpGradFace(), and Moose::Kokkos::System::getVectorQpValueFace().

116  {
117  return _n_dofs(elem_type, fe_type);
118  }
Array2D< unsigned int > _n_dofs

◆ getNumFaceQps() [1/2]

KOKKOS_FUNCTION dof_id_type Moose::Kokkos::Assembly::getNumFaceQps ( ContiguousSubdomainID  subdomain) const
inline

Get the total number of facial quadrature points in a subdomain Note: this number does not represent the real number of facial quadrature points but only counts the facial quadrature points that need global caching, such as face material properties.

Parameters
subdomainThe contiguous subdomain ID
Returns
The number of quadrature points

Definition at line 77 of file KokkosAssembly.h.

Referenced by Moose::Kokkos::MaterialProperty< T, dimension >::allocate().

78  {
79  return _n_subdomain_qps_face[subdomain];
80  }
Array< dof_id_type > _n_subdomain_qps_face

◆ getNumFaceQps() [2/2]

KOKKOS_FUNCTION unsigned int Moose::Kokkos::Assembly::getNumFaceQps ( ElementInfo  info,
unsigned int  side 
) const
inline

Get the number of quadrature points of a side of an element.

Parameters
infoThe element information object
sideThe side index
Returns
The number of quadrature points

Definition at line 87 of file KokkosAssembly.h.

88  {
89  return _n_qps_face(side, info.id);
90  }
Array2D< unsigned int > _n_qps_face
MPI_Info info

◆ getNumQps() [1/2]

KOKKOS_FUNCTION dof_id_type Moose::Kokkos::Assembly::getNumQps ( ContiguousSubdomainID  subdomain) const
inline

Get the total number of elemental quadrature points in a subdomain.

Parameters
subdomainThe contiguous subdomain ID
Returns
The number of quadrature points

Definition at line 60 of file KokkosAssembly.h.

Referenced by Moose::Kokkos::MaterialProperty< T, dimension >::allocate().

61  {
62  return _n_subdomain_qps[subdomain];
63  }
Array< dof_id_type > _n_subdomain_qps

◆ getNumQps() [2/2]

KOKKOS_FUNCTION unsigned int Moose::Kokkos::Assembly::getNumQps ( ElementInfo  info) const
inline

Get the number of quadrature points of an element.

Parameters
infoThe element information object
Returns
The number of quadrature points

Definition at line 69 of file KokkosAssembly.h.

69 { return _n_qps[info.id]; }
MPI_Info info
Array< unsigned int > _n_qps
Number of quadrature points.

◆ getPhi()

KOKKOS_FUNCTION const auto& Moose::Kokkos::Assembly::getPhi ( ContiguousSubdomainID  subdomain,
unsigned int  elem_type,
unsigned int  fe_type 
) const
inline

Get the shape functions of a FE type for an element type and subdomain.

Parameters
subdomainThe contiguous subdomain ID
elem_typeThe element type ID
fe_typeThe FE type ID
Returns
The shape functions at quadrature points

Definition at line 127 of file KokkosAssembly.h.

Referenced by Moose::Kokkos::VariablePhiValue::operator()(), and Moose::Kokkos::VariableTestValue::operator()().

128  {
129  return _phi(subdomain, elem_type, fe_type);
130  }
Array3D< Array2D< Real > > _phi
Shape functions for reference elements.

◆ getPhiFace()

KOKKOS_FUNCTION const auto& Moose::Kokkos::Assembly::getPhiFace ( ContiguousSubdomainID  subdomain,
unsigned int  elem_type,
unsigned int  fe_type 
) const
inline

Get the face shape functions of a FE type for an element type and subdomain.

Parameters
subdomainThe contiguous subdomain ID
elem_typeThe element type ID
fe_typeThe FE type ID
Returns
The shape functions of all sides at quadrature points

Definition at line 139 of file KokkosAssembly.h.

Referenced by Moose::Kokkos::System::getVectorQpValueFace(), Moose::Kokkos::VariablePhiValue::operator()(), and Moose::Kokkos::VariableTestValue::operator()().

140  {
141  return _phi_face(subdomain, elem_type, fe_type);
142  }
Array3D< Array< Array2D< Real > > > _phi_face

◆ getQpFaceOffset()

KOKKOS_FUNCTION dof_id_type Moose::Kokkos::Assembly::getQpFaceOffset ( ElementInfo  info,
unsigned int  side 
) const
inline

Get the starting offset of quadrature points of a side of an element into the global quadrature point index.

Parameters
infoThe element information object
sideThe side index
Returns
The starting offset

Definition at line 105 of file KokkosAssembly.h.

106  {
107  return _qp_offset_face(side, info.id);
108  }
MPI_Info info
Array2D< dof_id_type > _qp_offset_face

◆ getQpOffset()

KOKKOS_FUNCTION dof_id_type Moose::Kokkos::Assembly::getQpOffset ( ElementInfo  info) const
inline

Get the starting offset of quadrature points of an element into the global quadrature point index.

Parameters
infoThe element information object
Returns
The starting offset

Definition at line 97 of file KokkosAssembly.h.

Referenced by getJacobian(), getJxW(), and getQPoint().

97 { return _qp_offset[info.id]; }
MPI_Info info
Array< dof_id_type > _qp_offset
Starting offset into the global quadrature point index.

◆ getQPoint()

KOKKOS_FUNCTION Real3 Moose::Kokkos::Assembly::getQPoint ( ElementInfo  info,
unsigned int  qp 
) const
inline

Get the coordinate of an element quadrature point.

Parameters
infoThe element information object
qpThe local quadrature point index
Returns
The inverse of Jacobian matrix

Definition at line 194 of file KokkosAssembly.h.

Referenced by Moose::Kokkos::Datum::reinitTransform().

195  {
196  return _xyz[info.subdomain][getQpOffset(info) + qp];
197  }
KOKKOS_FUNCTION dof_id_type getQpOffset(ElementInfo info) const
Get the starting offset of quadrature points of an element into the global quadrature point index...
MPI_Info info
Array< Array< Real3 > > _xyz

◆ init()

void Moose::Kokkos::Assembly::init ( )

Initialize assembly.

◆ initQuadrature()

void Moose::Kokkos::Assembly::initQuadrature ( )
private

Initialize quadrature data.

◆ initShape()

void Moose::Kokkos::Assembly::initShape ( )
private

Initialize shape data.

◆ kokkosMesh()

KOKKOS_FUNCTION const Mesh& Moose::Kokkos::MeshHolder::kokkosMesh ( ) const
inlineinherited

Get the const reference of the Kokkos mesh.

Returns
The const reference of the Kokkos mesh depending on the architecture this function is being called on

Definition at line 360 of file KokkosMesh.h.

Referenced by computePhysicalMap().

361  {
362  KOKKOS_IF_ON_HOST(return _mesh_host;)
363 
364  return _mesh_device;
365  }
const Mesh _mesh_device
Device copy of the Kokkos mesh.
Definition: KokkosMesh.h:376
const Mesh & _mesh_host
Host reference of the Kokkos mesh.
Definition: KokkosMesh.h:372

◆ operator()()

KOKKOS_FUNCTION void Moose::Kokkos::Assembly::operator() ( const ThreadID  tid) const

Kokkos function for caching physical maps on element quadrature points.

Member Data Documentation

◆ _coord_type

Array<Moose::CoordinateSystemType> Moose::Kokkos::Assembly::_coord_type
private

Coordinate system type of each subdomain.

Definition at line 282 of file KokkosAssembly.h.

Referenced by coordTransformFactor().

◆ _dimension

const unsigned int Moose::Kokkos::Assembly::_dimension
private

Mesh dimension.

Definition at line 278 of file KokkosAssembly.h.

Referenced by computePhysicalMap().

◆ _fe_type_map

std::map<FEType, unsigned int> Moose::Kokkos::Assembly::_fe_type_map
private

FE type ID map.

Definition at line 273 of file KokkosAssembly.h.

Referenced by getFETypeID().

◆ _grad_phi

Array3D<Array2D<Real3> > Moose::Kokkos::Assembly::_grad_phi
private

Definition at line 326 of file KokkosAssembly.h.

Referenced by getGradPhi().

◆ _grad_phi_face

Array3D<Array<Array2D<Real3> > > Moose::Kokkos::Assembly::_grad_phi_face
private

Definition at line 327 of file KokkosAssembly.h.

Referenced by getGradPhiFace().

◆ _jacobian

Array<Array<Real33> > Moose::Kokkos::Assembly::_jacobian
private

Cached physical maps on element quadrature points.

Definition at line 345 of file KokkosAssembly.h.

Referenced by getJacobian().

◆ _jxw

Array<Array<Real> > Moose::Kokkos::Assembly::_jxw
private

Definition at line 346 of file KokkosAssembly.h.

Referenced by getJxW().

◆ _map_grad_phi

Array2D<Array2D<Real3> > Moose::Kokkos::Assembly::_map_grad_phi
private

Definition at line 337 of file KokkosAssembly.h.

Referenced by computePhysicalMap().

◆ _map_grad_phi_face

Array2D<Array<Array2D<Real3> > > Moose::Kokkos::Assembly::_map_grad_phi_face
private

Definition at line 338 of file KokkosAssembly.h.

Referenced by computePhysicalMap().

◆ _map_grad_psi_face

Array2D<Array<Array2D<Real3> > > Moose::Kokkos::Assembly::_map_grad_psi_face
private

Definition at line 339 of file KokkosAssembly.h.

Referenced by computePhysicalMap().

◆ _map_phi

Array2D<Array2D<Real> > Moose::Kokkos::Assembly::_map_phi
private

Shape functions for computing reference-to-physical maps.

Definition at line 334 of file KokkosAssembly.h.

Referenced by computePhysicalMap().

◆ _map_phi_face

Array2D<Array<Array2D<Real> > > Moose::Kokkos::Assembly::_map_phi_face
private

Definition at line 335 of file KokkosAssembly.h.

Referenced by computePhysicalMap().

◆ _map_psi_face

Array2D<Array<Array2D<Real> > > Moose::Kokkos::Assembly::_map_psi_face
private

Definition at line 336 of file KokkosAssembly.h.

◆ _material_boundaries

std::set<BoundaryID> Moose::Kokkos::Assembly::_material_boundaries
private

Boundaries to cache face material properties.

Definition at line 353 of file KokkosAssembly.h.

Referenced by getMaterialBoundaries().

◆ _max_qps_per_elem

unsigned int Moose::Kokkos::Assembly::_max_qps_per_elem = 0
private

Definition at line 306 of file KokkosAssembly.h.

Referenced by getMaxQpsPerElem().

◆ _mesh

MooseMesh& Moose::Kokkos::Assembly::_mesh
private

Reference of the MOOSE mesh.

Definition at line 269 of file KokkosAssembly.h.

◆ _n_dofs

Array2D<unsigned int> Moose::Kokkos::Assembly::_n_dofs
private

Definition at line 328 of file KokkosAssembly.h.

Referenced by getNumDofs().

◆ _n_qps

Array<unsigned int> Moose::Kokkos::Assembly::_n_qps
private

Number of quadrature points.

Definition at line 303 of file KokkosAssembly.h.

Referenced by getNumQps().

◆ _n_qps_face

Array2D<unsigned int> Moose::Kokkos::Assembly::_n_qps_face
private

Definition at line 304 of file KokkosAssembly.h.

Referenced by getNumFaceQps().

◆ _n_subdomain_qps

Array<dof_id_type> Moose::Kokkos::Assembly::_n_subdomain_qps
private

Definition at line 308 of file KokkosAssembly.h.

Referenced by getNumQps().

◆ _n_subdomain_qps_face

Array<dof_id_type> Moose::Kokkos::Assembly::_n_subdomain_qps_face
private

Definition at line 309 of file KokkosAssembly.h.

Referenced by getNumFaceQps().

◆ _phi

Array3D<Array2D<Real> > Moose::Kokkos::Assembly::_phi
private

Shape functions for reference elements.

Definition at line 324 of file KokkosAssembly.h.

Referenced by getPhi().

◆ _phi_face

Array3D<Array<Array2D<Real> > > Moose::Kokkos::Assembly::_phi_face
private

Definition at line 325 of file KokkosAssembly.h.

Referenced by getPhiFace().

◆ _problem

FEProblemBase& Moose::Kokkos::Assembly::_problem
private

Reference of the MOOSE problem.

Definition at line 265 of file KokkosAssembly.h.

◆ _q_points

Array2D<Array<Real3> > Moose::Kokkos::Assembly::_q_points
private

Quadrature points and weights for reference elements.

Definition at line 315 of file KokkosAssembly.h.

◆ _q_points_face

Array2D<Array<Array<Real3> > > Moose::Kokkos::Assembly::_q_points_face
private

Definition at line 316 of file KokkosAssembly.h.

◆ _qp_offset

Array<dof_id_type> Moose::Kokkos::Assembly::_qp_offset
private

Starting offset into the global quadrature point index.

Definition at line 296 of file KokkosAssembly.h.

Referenced by getQpOffset().

◆ _qp_offset_face

Array2D<dof_id_type> Moose::Kokkos::Assembly::_qp_offset_face
private

Definition at line 297 of file KokkosAssembly.h.

Referenced by getQpFaceOffset().

◆ _rz_axis

Array<Pair<Real3, Real3> > Moose::Kokkos::Assembly::_rz_axis
private

General axisymmetric axis of each subdomain in cylindrical coordinate system.

Definition at line 290 of file KokkosAssembly.h.

Referenced by coordTransformFactor().

◆ _rz_radial_coord

unsigned int Moose::Kokkos::Assembly::_rz_radial_coord = libMesh::invalid_uint
private

Radial coordinate index in cylindrical coordinate system.

Definition at line 286 of file KokkosAssembly.h.

Referenced by coordTransformFactor().

◆ _weights

Array2D<Array<Real> > Moose::Kokkos::Assembly::_weights
private

Definition at line 317 of file KokkosAssembly.h.

Referenced by computePhysicalMap().

◆ _weights_face

Array2D<Array<Array<Real> > > Moose::Kokkos::Assembly::_weights_face
private

Definition at line 318 of file KokkosAssembly.h.

Referenced by computePhysicalMap().

◆ _xyz

Array<Array<Real3> > Moose::Kokkos::Assembly::_xyz
private

Definition at line 347 of file KokkosAssembly.h.

Referenced by getQPoint().


The documentation for this class was generated from the following file: