https://mooseframework.inl.gov
KokkosVariableValue.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 "KokkosDatum.h"
13 
14 #include "MooseVariableBase.h"
15 
16 namespace Moose
17 {
18 namespace Kokkos
19 {
20 
24 class VariablePhiValue
26 {
27 public:
35  KOKKOS_FUNCTION Real operator()(AssemblyDatum & datum, unsigned int i, unsigned int qp) const
36  {
37  auto & elem = datum.elem();
38  auto side = datum.side();
39  auto fe = datum.jfe();
40 
41  return side == libMesh::invalid_uint
42  ? datum.assembly().getPhi(elem.subdomain, elem.type, fe)(i, qp)
43  : datum.assembly().getPhiFace(elem.subdomain, elem.type, fe)(side)(i, qp);
44  }
45 };
46 
48 {
49 public:
57  KOKKOS_FUNCTION Real3 operator()(AssemblyDatum & datum, unsigned int i, unsigned int qp) const
58  {
59  auto & elem = datum.elem();
60  auto side = datum.side();
61  auto fe = datum.jfe();
62 
63  return datum.J(qp) *
64  (side == libMesh::invalid_uint
65  ? datum.assembly().getGradPhi(elem.subdomain, elem.type, fe)(i, qp)
66  : datum.assembly().getGradPhiFace(elem.subdomain, elem.type, fe)(side)(i, qp));
67  }
68 };
69 
71 {
72 public:
80  KOKKOS_FUNCTION Real operator()(AssemblyDatum & datum, unsigned int i, unsigned int qp) const
81  {
82  auto & elem = datum.elem();
83  auto side = datum.side();
84  auto fe = datum.ife();
85 
86  return side == libMesh::invalid_uint
87  ? datum.assembly().getPhi(elem.subdomain, elem.type, fe)(i, qp)
88  : datum.assembly().getPhiFace(elem.subdomain, elem.type, fe)(side)(i, qp);
89  }
90 };
91 
93 {
94 public:
102  KOKKOS_FUNCTION Real3 operator()(AssemblyDatum & datum, unsigned int i, unsigned int qp) const
103  {
104  auto & elem = datum.elem();
105  auto side = datum.side();
106  auto fe = datum.ife();
107 
108  return datum.J(qp) *
109  (side == libMesh::invalid_uint
110  ? datum.assembly().getGradPhi(elem.subdomain, elem.type, fe)(i, qp)
111  : datum.assembly().getGradPhiFace(elem.subdomain, elem.type, fe)(side)(i, qp));
112  }
113 };
115 
119 class VariableValue
121 {
122 public:
126  VariableValue() = default;
132  VariableValue(Variable var, bool dof = false) : _var(var), _dof(dof) {}
140  const TagName & tag = Moose::SOLUTION_TAG,
141  bool dof = false)
142  : _var(var, tag), _dof(dof)
143  {
144  }
145 
150  KOKKOS_FUNCTION operator bool() const { return _var.coupled(); }
151 
159  KOKKOS_FUNCTION Real operator()(Datum & datum, unsigned int qp, unsigned int comp = 0) const;
160 
161 private:
169  bool _dof = false;
170 };
171 
173 {
174 public:
178  VariableGradient() = default;
189  VariableGradient(const MooseVariableBase & var, const TagName & tag = Moose::SOLUTION_TAG)
190  : _var(var, tag)
191  {
192  }
193 
198  KOKKOS_FUNCTION operator bool() const { return _var.coupled(); }
199 
207  KOKKOS_FUNCTION Real3 operator()(Datum & datum, unsigned int idx, unsigned int comp = 0) const;
208 
209 private:
214 };
216 
217 KOKKOS_FUNCTION inline Real
218 VariableValue::operator()(Datum & datum, unsigned int idx, unsigned int comp) const
219 {
220  KOKKOS_ASSERT(_var.initialized());
221 
222  if (_var.coupled())
223  {
224  auto & sys = datum.system(_var.sys(comp));
225  auto var = _var.var(comp);
226  auto tag = _var.tag();
227 
228  if (_dof)
229  {
230  unsigned int dof;
231 
232  if (datum.isNodal())
233  {
234  auto node = datum.node();
235  dof = sys.getNodeLocalDofIndex(node, 0, var);
236  }
237  else
238  {
239  auto elem = datum.elem().id;
240  dof = sys.getElemLocalDofIndex(elem, idx, var);
241  }
242 
243  return sys.getVectorDofValue(dof, tag);
244  }
245  else
246  {
247  auto & elem = datum.elem();
248  auto side = datum.side();
249  auto offset = datum.qpOffset();
250 
251  return side == libMesh::invalid_uint ? sys.getVectorQpValue(elem, offset + idx, var, tag)
252  : sys.getVectorQpValueFace(elem, side, idx, var, tag);
253  }
254  }
255  else
256  return _var.value(comp);
257 }
258 
259 KOKKOS_FUNCTION inline Real3
260 VariableGradient::operator()(Datum & datum, unsigned int qp, unsigned int comp) const
261 {
262  KOKKOS_ASSERT(_var.initialized());
263 
264  if (_var.coupled())
265  {
266  KOKKOS_ASSERT(!datum.isNodal());
267 
268  auto & elem = datum.elem();
269  auto side = datum.side();
270  auto offset = datum.qpOffset();
271 
272  return side == libMesh::invalid_uint
273  ? datum.system(_var.sys(comp))
274  .getVectorQpGrad(elem, offset + qp, _var.var(comp), _var.tag())
275  : datum.system(_var.sys(comp))
276  .getVectorQpGradFace(elem, side, datum.J(qp), qp, _var.var(comp), _var.tag());
277  }
278  else
279  return Real3(0);
280 }
281 
282 } // namespace Kokkos
283 } // namespace Moose
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.
KOKKOS_FUNCTION TagID tag() const
Get the vector tag ID.
Variable _var
Coupled Kokkos variable.
KOKKOS_FUNCTION unsigned int sys(unsigned int comp=0) const
Get the system number of a component.
The Kokkos object that holds thread-private data in the parallel operations of any Kokkos object...
Definition: KokkosDatum.h:25
KOKKOS_FUNCTION Real operator()(Datum &datum, unsigned int qp, unsigned int comp=0) const
Get the current variable value.
const unsigned int invalid_uint
KOKKOS_FUNCTION Real3 operator()(Datum &datum, unsigned int idx, unsigned int comp=0) const
Get the current variable gradient.
OutputTools< Real >::VariablePhiValue VariablePhiValue
Definition: MooseTypes.h:320
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.
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.
KOKKOS_FUNCTION bool initialized() const
Get whether the variable is initialized.
KOKKOS_FUNCTION Real3 operator()(AssemblyDatum &datum, unsigned int i, unsigned int qp) const
Get the gradient of the current shape function.
VariableValue()=default
Default constructor.
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.
KOKKOS_FUNCTION bool coupled() const
Get whether the variable is coupled.
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 side() const
Get the side index.
Definition: KokkosDatum.h:90
KOKKOS_FUNCTION unsigned int jfe() const
Get the coupled variable FE type ID.
Definition: KokkosDatum.h:344
VariableValue(Variable var, bool dof=false)
Constructor.
KOKKOS_FUNCTION const Assembly & assembly() const
Get the Kokkos assembly.
Definition: KokkosDatum.h:68
VariableGradient(const MooseVariableBase &var, const TagName &tag=Moose::SOLUTION_TAG)
Constructor.
KOKKOS_FUNCTION Real operator()(AssemblyDatum &datum, unsigned int i, unsigned int qp) const
Get the current shape function.
KOKKOS_FUNCTION Real3 operator()(AssemblyDatum &datum, unsigned int i, unsigned int qp) const
Get the gradient of the current test function.
VariableValue(const MooseVariableBase &var, const TagName &tag=Moose::SOLUTION_TAG, bool dof=false)
Constructor.
OutputTools< Real >::VariableValue VariableValue
Definition: MooseTypes.h:315
KOKKOS_FUNCTION ContiguousNodeID node() const
Get the contiguous node ID.
Definition: KokkosDatum.h:95
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
KOKKOS_FUNCTION Real value(unsigned int comp=0) const
Get the default value of a component.
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
const TagName SOLUTION_TAG
Definition: MooseTypes.C:25
Variable _var
Coupled Kokkos variable.
MOOSE now contains C++17 code, so give a reasonable error message stating what the user can do to add...
KOKKOS_FUNCTION unsigned int var(unsigned int comp=0) const
Get the variable number of a component.
VariableGradient(Variable var)
Constructor.
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
KOKKOS_FUNCTION Real operator()(AssemblyDatum &datum, unsigned int i, unsigned int qp) const
Get the current test function.
VariableGradient()=default
Default constructor.
bool _dof
Flag whether DOF values are requested.
ContiguousElementID id
Contiguous element ID.
Definition: KokkosMesh.h:44
Base variable class.
unsigned int idx(const ElemType type, const unsigned int nx, const unsigned int i, const unsigned int j)