https://mooseframework.inl.gov
KokkosIntegratedBC.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 "KokkosIntegratedBCBase.h"
13 
14 namespace Moose
15 {
16 namespace Kokkos
17 {
18 
40 {
41 public:
43 
48 
52  virtual void computeResidual() override;
56  virtual void computeJacobian() override;
57 
62 
71  KOKKOS_FUNCTION Real computeQpJacobian(const unsigned int /* i */,
72  const unsigned int /* j */,
73  const unsigned int /* qp */,
74  ResidualDatum & /* datum */) const
75  {
76  return 0;
77  }
87  KOKKOS_FUNCTION Real computeQpOffDiagJacobian(const unsigned int /* i */,
88  const unsigned int /* j */,
89  const unsigned int /* jvar */,
90  const unsigned int /* qp */,
91  ResidualDatum & /* datum */) const
92  {
93  return 0;
94  }
106 
110  template <typename Derived>
112  KOKKOS_FUNCTION void operator()(ResidualLoop, const ThreadID tid, const Derived & bc) const;
113  template <typename Derived>
114  KOKKOS_FUNCTION void operator()(JacobianLoop, const ThreadID tid, const Derived & bc) const;
115  template <typename Derived>
116  KOKKOS_FUNCTION void
117  operator()(OffDiagJacobianLoop, const ThreadID tid, const Derived & bc) const;
119 
125 
131  template <typename Derived>
132  KOKKOS_FUNCTION void computeResidualInternal(const Derived & bc, ResidualDatum & datum) const;
138  template <typename Derived>
139  KOKKOS_FUNCTION void computeJacobianInternal(const Derived & bc, ResidualDatum & datum) const;
145  template <typename Derived>
146  KOKKOS_FUNCTION void computeOffDiagJacobianInternal(const Derived & bc,
147  ResidualDatum & datum) const;
149 
150 protected:
175 };
176 
177 template <typename Derived>
178 KOKKOS_FUNCTION void
179 IntegratedBC::operator()(ResidualLoop, const ThreadID tid, const Derived & bc) const
180 {
181  auto [elem, side] = kokkosBoundaryElementSideID(tid);
182 
183  ResidualDatum datum(
185 
186  bc.computeResidualInternal(bc, datum);
187 }
188 
189 template <typename Derived>
190 KOKKOS_FUNCTION void
191 IntegratedBC::operator()(JacobianLoop, const ThreadID tid, const Derived & bc) const
192 {
193  auto [elem, side] = kokkosBoundaryElementSideID(tid);
194 
195  ResidualDatum datum(
197 
198  bc.computeJacobianInternal(bc, datum);
199 }
200 
201 template <typename Derived>
202 KOKKOS_FUNCTION void
203 IntegratedBC::operator()(OffDiagJacobianLoop, const ThreadID tid, const Derived & bc) const
204 {
205  auto [elem, side] = kokkosBoundaryElementSideID(_thread(tid, 1));
206 
207  auto & sys = kokkosSystem(_kokkos_var.sys());
208  auto jvar = sys.getCoupling(_kokkos_var.var())[_thread(tid, 0)];
209 
210  if (!sys.isVariableActive(jvar, kokkosMesh().getElementInfo(elem).subdomain))
211  return;
212 
213  ResidualDatum datum(elem, side, kokkosAssembly(), kokkosSystems(), _kokkos_var, jvar);
214 
215  bc.computeOffDiagJacobianInternal(bc, datum);
216 }
217 
218 template <typename Derived>
219 KOKKOS_FUNCTION void
220 IntegratedBC::computeResidualInternal(const Derived & bc, ResidualDatum & datum) const
221 {
223  datum,
224  [&](Real * local_re, const unsigned int ib, const unsigned int ie)
225  {
226  for (unsigned int qp = 0; qp < datum.n_qps(); ++qp)
227  {
228  datum.reinit();
229 
230  for (unsigned int i = ib; i < ie; ++i)
231  local_re[i] += datum.JxW(qp) * bc.computeQpResidual(i, qp, datum);
232  }
233  });
234 }
235 
236 template <typename Derived>
237 KOKKOS_FUNCTION void
238 IntegratedBC::computeJacobianInternal(const Derived & bc, ResidualDatum & datum) const
239 {
241  datum,
242  [&](Real * local_ke, const unsigned int ijb, const unsigned int ije)
243  {
244  for (unsigned int qp = 0; qp < datum.n_qps(); ++qp)
245  {
246  datum.reinit();
247 
248  for (unsigned int ij = ijb; ij < ije; ++ij)
249  {
250  unsigned int i = ij % datum.n_jdofs();
251  unsigned int j = ij / datum.n_jdofs();
252 
253  local_ke[ij] += datum.JxW(qp) * bc.computeQpJacobian(i, j, qp, datum);
254  }
255  }
256  });
257 }
258 
259 template <typename Derived>
260 KOKKOS_FUNCTION void
262 {
264  datum,
265  [&](Real * local_ke, const unsigned int ijb, const unsigned int ije)
266  {
267  for (unsigned int qp = 0; qp < datum.n_qps(); ++qp)
268  {
269  datum.reinit();
270 
271  for (unsigned int ij = ijb; ij < ije; ++ij)
272  {
273  unsigned int i = ij % datum.n_jdofs();
274  unsigned int j = ij / datum.n_jdofs();
275 
276  local_ke[ij] +=
277  datum.JxW(qp) * bc.computeQpOffDiagJacobian(i, j, datum.jvar(), qp, datum);
278  }
279  }
280  });
281 }
282 
283 } // namespace Kokkos
284 } // namespace Moose
KOKKOS_FUNCTION unsigned int n_jdofs() const
Get the number of local DOFs for the coupled variable.
Definition: KokkosDatum.h:315
KOKKOS_FUNCTION void computeJacobianInternal(const Derived &bc, ResidualDatum &datum) const
Compute diagonal Jacobian.
KOKKOS_FUNCTION unsigned int sys(unsigned int comp=0) const
Get the system number of a component.
KOKKOS_FUNCTION const Assembly & kokkosAssembly() const
Get the const reference of the Kokkos assembly.
static auto defaultJacobian()
Get the function pointer of the default computeQpJacobian()
dof_id_type ThreadID
Definition: KokkosThread.h:18
KOKKOS_FUNCTION void computeJacobianInternal(ResidualDatum &datum, function body) const
The common loop structure template for computing elemental Jacobian.
The base class for a user to derive their own Kokkos integrated boundary conditions.
const InputParameters & parameters() const
Get the parameters of the object.
Definition: MooseBase.h:131
The main MOOSE class responsible for handling user-defined parameters in almost every MOOSE system...
Thread _thread
Kokkos thread object.
The base class for Kokkos integrated boundary conditions.
const VariableValue _u
Current solution at quadrature points.
The Kokkos wrapper classes for MOOSE-like shape function access.
const VariableGradient _grad_u
Gradient of the current solution at quadrature points.
const VariablePhiGradient _grad_phi
Gradient of the current shape function.
KOKKOS_FUNCTION unsigned int n_qps() const
Get the number of local quadrature points.
Definition: KokkosDatum.h:100
KOKKOS_FUNCTION void operator()(ResidualLoop, const ThreadID tid, const Derived &bc) const
The parallel computation entry functions called by Kokkos.
KOKKOS_FUNCTION const System & kokkosSystem(unsigned int sys) const
Get the const reference of a Kokkos system.
Definition: KokkosSystem.h:615
const VariableTestGradient _grad_test
Gradient of the current test function.
static auto defaultOffDiagJacobian()
Get the function pointer of the default computeQpOffDiagJacobian()
const VariableTestValue _test
Current test function.
KOKKOS_FUNCTION const Mesh & kokkosMesh() const
Get the const reference of the Kokkos mesh.
Definition: KokkosMesh.h:360
KOKKOS_FUNCTION void computeResidualInternal(const Derived &bc, ResidualDatum &datum) const
The parallel computation bodies that can be customized in the derived class by defining them in the d...
KOKKOS_FUNCTION Real computeQpOffDiagJacobian(const unsigned int, const unsigned int, const unsigned int, const unsigned int, ResidualDatum &) const
Compute off-diagonal Jacobian contribution on a quadrature point.
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
KOKKOS_FUNCTION auto kokkosBoundaryElementSideID(ThreadID tid) const
Get the contiguous element ID - side index pair this Kokkos thread is operating on.
KOKKOS_FUNCTION void computeResidualInternal(ResidualDatum &datum, function body) const
The common loop structure template for computing elemental residual.
virtual void computeJacobian() override
Dispatch diagonal and off-diagonal Jacobian calculation.
const VariablePhiValue _phi
Current shape function.
virtual void computeResidual() override
Dispatch residual calculation.
KOKKOS_FUNCTION Real JxW(const unsigned int qp)
Get the transformed Jacobian weight.
Definition: KokkosDatum.h:139
Variable _kokkos_var
Kokkos variable.
static InputParameters validParams()
MOOSE now contains C++17 code, so give a reasonable error message stating what the user can do to add...
KOKKOS_FUNCTION Real computeQpJacobian(const unsigned int, const unsigned int, const unsigned int, ResidualDatum &) const
Default methods to prevent compile errors even when these methods were not defined in the derived cla...
The Kokkos wrapper classes for MOOSE-like variable value access.
KOKKOS_FUNCTION unsigned int var(unsigned int comp=0) const
Get the variable number of a component.
KOKKOS_FUNCTION void computeOffDiagJacobianInternal(const Derived &bc, ResidualDatum &datum) const
Compute off-diagonal Jacobian.
KOKKOS_FUNCTION const auto & getCoupling(unsigned int var) const
Get the list of off-diagonal coupled variable numbers of a variable.
Definition: KokkosSystem.h:140
KOKKOS_FUNCTION unsigned int jvar() const
Get the coupled variable number.
Definition: KokkosDatum.h:330
The Kokkos object that holds thread-private data in the parallel operations of Kokkos residual object...
Definition: KokkosDatum.h:245
KOKKOS_FUNCTION void reinit()
Reset the reinit flag.
Definition: KokkosDatum.h:166
IntegratedBC(const InputParameters &parameters)
Constructor.
KOKKOS_FUNCTION const Array< System > & kokkosSystems() const
Get the const reference of the Kokkos systems.
Definition: KokkosSystem.h:598