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 -- 1D 10 template <typename T, int Q_> 11 __device__ __inline__ void magma_weight_1d_device(const T *sTweight, T *sV, const int tx) { 12 // Assumptions 13 // 1. 1D thread configuration of size Q_ 14 // 2. The output sV is in shared memory -- size 1xQ_ 15 if (tx < Q_) { 16 sV[tx] = sTweight[tx]; 17 } 18 } 19 20 ////////////////////////////////////////////////////////////////////////////////////////// 21 extern "C" __launch_bounds__(MAGMA_BASIS_BOUNDS(Q, MAGMA_MAXTHREADS_1D)) __global__ 22 void magma_weight_1d_kernel(const CeedScalar *dqweight1d, CeedScalar *dV, const int v_stride, const int nelem) { 23 MAGMA_DEVICE_SHARED(CeedScalar, shared_data) 24 25 const int tx = threadIdx.x; 26 const int ty = threadIdx.y; 27 const int elem_id = (blockIdx.x * blockDim.y) + ty; 28 29 if (elem_id >= nelem) return; 30 31 // global memory pointers 32 dV += elem_id * v_stride; 33 34 // shared memory pointers 35 CeedScalar *sTweight = (CeedScalar *)shared_data; 36 CeedScalar *sV = sTweight + Q; 37 sV += ty * Q; 38 39 // read dqweight_1d 40 if (ty == 0 && tx < Q) { 41 sTweight[tx] = dqweight1d[tx]; 42 } 43 44 __syncthreads(); 45 magma_weight_1d_device<CeedScalar, Q>(sTweight, sV, tx); 46 __syncthreads(); 47 48 // write V 49 dV[tx] = sV[tx]; 50 } 51