https://mooseframework.inl.gov
KokkosKernelGrad.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 "KokkosKernel.h"
13 
14 namespace Moose
15 {
16 namespace Kokkos
17 {
18 
41 class KernelGrad : public Kernel
42 {
43 public:
45 
50 
55 
64  KOKKOS_FUNCTION Real3 computeQpJacobian(const unsigned int /* j */,
65  const unsigned int /* qp */,
66  ResidualDatum & /* datum */) const
67  {
68  return Real3(0);
69  }
76 
81  template <typename Derived>
83  KOKKOS_FUNCTION void computeResidualInternal(const Derived & kernel, ResidualDatum & datum) const;
84  template <typename Derived>
85  KOKKOS_FUNCTION void computeJacobianInternal(const Derived & kernel, ResidualDatum & datum) const;
87 };
88 
89 template <typename Derived>
90 KOKKOS_FUNCTION void
91 KernelGrad::computeResidualInternal(const Derived & kernel, ResidualDatum & datum) const
92 {
94  datum,
95  [&](Real * local_re, const unsigned int ib, const unsigned int ie)
96  {
97  for (unsigned int qp = 0; qp < datum.n_qps(); ++qp)
98  {
99  datum.reinit();
100 
101  Real3 value = datum.JxW(qp) * kernel.computeQpResidual(qp, datum);
102 
103  for (unsigned int i = ib; i < ie; ++i)
104  local_re[i] += value * _grad_test(datum, i, qp);
105  }
106  });
107 }
108 
109 template <typename Derived>
110 KOKKOS_FUNCTION void
111 KernelGrad::computeJacobianInternal(const Derived & kernel, ResidualDatum & datum) const
112 {
113  Real3 value;
114 
116  datum,
117  [&](Real * local_ke, const unsigned int ijb, const unsigned int ije)
118  {
119  for (unsigned int qp = 0; qp < datum.n_qps(); ++qp)
120  {
121  datum.reinit();
122 
123  unsigned int j_old = libMesh::invalid_uint;
124 
125  for (unsigned int ij = ijb; ij < ije; ++ij)
126  {
127  unsigned int i = ij % datum.n_jdofs();
128  unsigned int j = ij / datum.n_jdofs();
129 
130  if (j != j_old)
131  {
132  value = datum.JxW(qp) * kernel.computeQpJacobian(j, qp, datum);
133  j_old = j;
134  }
135 
136  local_ke[ij] += value * _grad_test(datum, i, qp);
137  }
138  }
139  });
140 }
141 
142 } // namespace Kokkos
143 } // namespace Moose
KOKKOS_FUNCTION unsigned int n_jdofs() const
Get the number of local DOFs for the coupled variable.
Definition: KokkosDatum.h:315
const unsigned int invalid_uint
KOKKOS_FUNCTION void computeJacobianInternal(ResidualDatum &datum, function body) const
The common loop structure template for computing elemental Jacobian.
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...
KOKKOS_FUNCTION void computeResidualInternal(const Derived &kernel, ResidualDatum &datum) const
The parallel computation bodies that hide the base class methods to optimize for factoring out the gr...
KernelGrad(const InputParameters &parameters)
Constructor.
static auto defaultJacobian()
Get the function pointer of the default computeQpJacobian()
Real value(unsigned n, unsigned alpha, unsigned beta, Real x)
KOKKOS_FUNCTION unsigned int n_qps() const
Get the number of local quadrature points.
Definition: KokkosDatum.h:100
static InputParameters validParams()
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
const VariableTestGradient _grad_test
Gradient of the current test function.
Definition: KokkosKernel.h:159
KOKKOS_FUNCTION void computeResidualInternal(ResidualDatum &datum, function body) const
The common loop structure template for computing elemental residual.
KOKKOS_FUNCTION Real3 computeQpJacobian(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...
KOKKOS_FUNCTION Real JxW(const unsigned int qp)
Get the transformed Jacobian weight.
Definition: KokkosDatum.h:139
MOOSE now contains C++17 code, so give a reasonable error message stating what the user can do to add...
The base class for a user to derive their own Kokkos kernels.
Definition: KokkosKernel.h:40
The base class for a user to derive their own Kokkos kernels where the residual is of the form...
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
KOKKOS_FUNCTION void computeJacobianInternal(const Derived &kernel, ResidualDatum &datum) const