https://mooseframework.inl.gov
KokkosKernelValue.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 
40 class KernelValue : public Kernel
41 {
42 public:
44 
49 
54 
63  KOKKOS_FUNCTION Real computeQpJacobian(const unsigned int /* j */,
64  const unsigned int /* qp */,
65  ResidualDatum & /* datum */) const
66  {
67  return 0;
68  }
75 
80  template <typename Derived>
82  KOKKOS_FUNCTION void computeResidualInternal(const Derived & kernel, ResidualDatum & datum) const;
83  template <typename Derived>
84  KOKKOS_FUNCTION void computeJacobianInternal(const Derived & kernel, ResidualDatum & datum) const;
86 };
87 
88 template <typename Derived>
89 KOKKOS_FUNCTION void
90 KernelValue::computeResidualInternal(const Derived & kernel, ResidualDatum & datum) const
91 {
93  datum,
94  [&](Real * local_re, const unsigned int ib, const unsigned int ie)
95  {
96  for (unsigned int qp = 0; qp < datum.n_qps(); ++qp)
97  {
98  datum.reinit();
99 
100  Real value = datum.JxW(qp) * kernel.computeQpResidual(qp, datum);
101 
102  for (unsigned int i = ib; i < ie; ++i)
103  local_re[i] += value * _test(datum, i, qp);
104  }
105  });
106 }
107 
108 template <typename Derived>
109 KOKKOS_FUNCTION void
110 KernelValue::computeJacobianInternal(const Derived & kernel, ResidualDatum & datum) const
111 {
112  Real value = 0;
113 
115  datum,
116  [&](Real * local_ke, const unsigned int ijb, const unsigned int ije)
117  {
118  for (unsigned int qp = 0; qp < datum.n_qps(); ++qp)
119  {
120  datum.reinit();
121 
122  unsigned int j_old = libMesh::invalid_uint;
123 
124  for (unsigned int ij = ijb; ij < ije; ++ij)
125  {
126  unsigned int i = ij % datum.n_jdofs();
127  unsigned int j = ij / datum.n_jdofs();
128 
129  if (j != j_old)
130  {
131  value = datum.JxW(qp) * kernel.computeQpJacobian(j, qp, datum);
132  j_old = j;
133  }
134 
135  local_ke[ij] += value * _test(datum, i, qp);
136  }
137  }
138  });
139 }
140 
141 } // namespace Kokkos
142 } // namespace Moose
static auto defaultJacobian()
Get the function pointer of the default computeQpJacobian()
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
const VariableTestValue _test
Current test function.
Definition: KokkosKernel.h:155
The main MOOSE class responsible for handling user-defined parameters in almost every MOOSE system...
The base class for a user to derive their own Kokkos kernels where the residual is of the form...
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 te...
KOKKOS_FUNCTION void computeJacobianInternal(const Derived &kernel, ResidualDatum &datum) const
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
KOKKOS_FUNCTION Real 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...
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
KernelValue(const InputParameters &parameters)
Constructor.
KOKKOS_FUNCTION void computeResidualInternal(ResidualDatum &datum, function body) const
The common loop structure template for computing elemental residual.
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 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
static InputParameters validParams()