1 // Copyright (c) 2017-2022, Lawrence Livermore National Security, LLC and other CEED contributors. 2 // All Rights Reserved. See the top-level LICENSE and NOTICE files for details. 3 // 4 // SPDX-License-Identifier: BSD-2-Clause 5 // 6 // This file is part of CEED: http://github.com/ceed 7 8 ////////////////////////////////////////////////////////////////////////////////////////// 9 // weight basis action -- 3D 10 template <typename T, int DIM_, int NCOMP_, int Q_, int iDIM, int iCOMP> 11 __device__ __inline__ void magma_weight_3d_device(const T *sTweight, T rV[DIM_][NCOMP_][Q_], const int tx) { 12 // Assumptions 13 // 1. 1D thread configuration of size Q_^2 14 // 2. rV[][][] matches the storage used in other actions (interp, grad, ... etc) 15 // 3. iDIM and iCOMP specify which indexes to use in rV, 16 // since the output per thread is a register array of size Q_ 17 // 4. Sync is recommended after the call (to make sure sTweight can be overwritten) 18 19 if (tx < (Q_ * Q_)) { 20 // x sTweight[j] for first update 21 // x sTweight[tx%Q_] for second update 22 // x sTweight[tx/Q_] for third update 23 for (int j = 0; j < Q_; j++) { 24 rV[iDIM][iCOMP][j] = sTweight[j] * sTweight[tx % Q_] * sTweight[tx / Q_]; 25 } 26 } 27 } 28 29 ////////////////////////////////////////////////////////////////////////////////////////// 30 extern "C" __launch_bounds__(MAGMA_BASIS_BOUNDS(Q *Q, MAGMA_MAXTHREADS_3D)) __global__ 31 void magma_weight_3d_kernel(const CeedScalar *dqweight1d, CeedScalar *dV, const int v_stride, const int nelem) { 32 MAGMA_DEVICE_SHARED(CeedScalar, shared_data) 33 34 const int tx = threadIdx.x; 35 const int ty = threadIdx.y; 36 const int elem_id = (blockIdx.x * blockDim.y) + ty; 37 38 if (elem_id >= nelem) return; 39 40 CeedScalar rV[1][1][Q]; // allocate with DIM=NCOMP=1, but sizes may differ for a fused operator 41 // global memory pointers 42 dV += elem_id * v_stride; 43 44 // shared memory pointers 45 CeedScalar *sTweight = (CeedScalar *)shared_data; 46 47 // read dqweight_1d 48 if (tx < Q) { 49 sTweight[tx] = dqweight1d[tx]; 50 } 51 __syncthreads(); 52 53 magma_weight_3d_device<CeedScalar, 1, 1, Q, 0, 0>(sTweight, rV, tx); 54 55 // write V 56 if (tx < (Q * Q)) { 57 for (int j = 0; j < Q; j++) { 58 dV[j * (Q * Q) + tx] = rV[0][0][j]; 59 } 60 } 61 } 62