Line data Source code
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 : 19 : /** 20 : * The base class for a user to derive their own Kokkos kernels where the residual is of the form 21 : * 22 : * $(\dots, \psi_i)$ 23 : * 24 : * i.e. the test function $(\psi_i)$ can be factored out for optimization. 25 : * 26 : * The user should still define computeQpResidual() and computeQpJacobian(), but their signatures 27 : * are different from the base class. The signature of computeQpResidual() expected to be defined in 28 : * the derived class is as follows: 29 : * 30 : * @param qp The local quadrature point index 31 : * @param datum The ResidualDatum object of the current thread 32 : * @returns The component of the residual contribution that will be multiplied by the test function 33 : * 34 : * KOKKOS_FUNCTION Real computeQpResidual(const unsigned int qp, 35 : * ResidualDatum & datum) const; 36 : * 37 : * The signature of computeQpJacobian() can be found in the code below. The definition of 38 : * computeQpOffDiagJacobian() is still the same with the original Kokkos kernel. 39 : */ 40 : class KernelValue : public Kernel 41 : { 42 : public: 43 : static InputParameters validParams(); 44 : 45 : /** 46 : * Constructor 47 : */ 48 : KernelValue(const InputParameters & parameters); 49 : 50 : /** 51 : * Default methods to prevent compile errors even when these methods were not defined in the 52 : * derived class 53 : */ 54 : ///@{ 55 : /** 56 : * Compute diagonal Jacobian contribution on a quadrature point 57 : * @param j The trial function DOF index 58 : * @param qp The local quadrature point index 59 : * @param datum The ResidualDatum object of the current thread 60 : * @returns The component of the Jacobian contribution that will be multiplied by the test 61 : * function 62 : */ 63 0 : KOKKOS_FUNCTION Real computeQpJacobian(const unsigned int /* j */, 64 : const unsigned int /* qp */, 65 : ResidualDatum & /* datum */) const 66 : { 67 0 : return 0; 68 : } 69 : /** 70 : * Get the function pointer of the default computeQpJacobian() 71 : * @returns The function pointer 72 : */ 73 37361 : static auto defaultJacobian() { return &KernelValue::computeQpJacobian; } 74 : ///@} 75 : 76 : /** 77 : * The parallel computation bodies that hide the base class methods to optimize for factoring 78 : * out the test function 79 : */ 80 : ///@{ 81 : 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; 85 : ///@} 86 : }; 87 : 88 : template <typename Derived> 89 : KOKKOS_FUNCTION void 90 96 : KernelValue::computeResidualInternal(const Derived & kernel, ResidualDatum & datum) const 91 : { 92 96 : ResidualObject::computeResidualInternal( 93 : datum, 94 192 : [&](Real * local_re, const unsigned int ib, const unsigned int ie) 95 : { 96 480 : for (unsigned int qp = 0; qp < datum.n_qps(); ++qp) 97 : { 98 384 : datum.reinit(); 99 : 100 384 : Real value = datum.JxW(qp) * kernel.computeQpResidual(qp, datum); 101 : 102 1920 : for (unsigned int i = ib; i < ie; ++i) 103 1536 : local_re[i] += value * _test(datum, i, qp); 104 : } 105 : }); 106 96 : } 107 : 108 : template <typename Derived> 109 : KOKKOS_FUNCTION void 110 24 : KernelValue::computeJacobianInternal(const Derived & kernel, ResidualDatum & datum) const 111 : { 112 24 : Real value = 0; 113 : 114 24 : ResidualObject::computeJacobianInternal( 115 : datum, 116 48 : [&](Real * local_ke, const unsigned int ijb, const unsigned int ije) 117 : { 118 120 : for (unsigned int qp = 0; qp < datum.n_qps(); ++qp) 119 : { 120 96 : datum.reinit(); 121 : 122 96 : unsigned int j_old = libMesh::invalid_uint; 123 : 124 1632 : for (unsigned int ij = ijb; ij < ije; ++ij) 125 : { 126 1536 : unsigned int i = ij % datum.n_jdofs(); 127 1536 : unsigned int j = ij / datum.n_jdofs(); 128 : 129 1536 : if (j != j_old) 130 : { 131 384 : value = datum.JxW(qp) * kernel.computeQpJacobian(j, qp, datum); 132 384 : j_old = j; 133 : } 134 : 135 1536 : local_ke[ij] += value * _test(datum, i, qp); 136 : } 137 : } 138 : }); 139 24 : } 140 : 141 : } // namespace Kokkos 142 : } // namespace Moose