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 /// @file 9 /// Internal header for CUDA strided element restriction kernels 10 #ifndef CEED_CUDA_REF_RESTRICTION_STRIDED_H 11 #define CEED_CUDA_REF_RESTRICTION_STRIDED_H 12 13 #include <ceed.h> 14 15 //------------------------------------------------------------------------------ 16 // L-vector -> E-vector, strided 17 //------------------------------------------------------------------------------ 18 extern "C" __global__ void StridedNoTranspose(const CeedScalar *__restrict__ u, CeedScalar *__restrict__ v) { 19 for (CeedInt node = blockIdx.x * blockDim.x + threadIdx.x; node < RSTR_NUM_ELEM * RSTR_ELEM_SIZE; node += blockDim.x * gridDim.x) { 20 const CeedInt loc_node = node % RSTR_ELEM_SIZE; 21 const CeedInt elem = node / RSTR_ELEM_SIZE; 22 23 for (CeedInt comp = 0; comp < RSTR_NUM_COMP; comp++) { 24 v[loc_node + comp * RSTR_ELEM_SIZE * RSTR_NUM_ELEM + elem * RSTR_ELEM_SIZE] = 25 u[loc_node * RSTR_STRIDE_NODES + comp * RSTR_STRIDE_COMP + elem * RSTR_STRIDE_ELEM]; 26 } 27 } 28 } 29 30 //------------------------------------------------------------------------------ 31 // E-vector -> L-vector, strided 32 //------------------------------------------------------------------------------ 33 extern "C" __global__ void StridedTranspose(const CeedScalar *__restrict__ u, CeedScalar *__restrict__ v) { 34 for (CeedInt node = blockIdx.x * blockDim.x + threadIdx.x; node < RSTR_NUM_ELEM * RSTR_ELEM_SIZE; node += blockDim.x * gridDim.x) { 35 const CeedInt loc_node = node % RSTR_ELEM_SIZE; 36 const CeedInt elem = node / RSTR_ELEM_SIZE; 37 38 for (CeedInt comp = 0; comp < RSTR_NUM_COMP; comp++) { 39 v[loc_node * RSTR_STRIDE_NODES + comp * RSTR_STRIDE_COMP + elem * RSTR_STRIDE_ELEM] += 40 u[loc_node + comp * RSTR_ELEM_SIZE * RSTR_NUM_ELEM + elem * RSTR_ELEM_SIZE]; 41 } 42 } 43 } 44 45 //------------------------------------------------------------------------------ 46 47 #endif // CEED_CUDA_REF_RESTRICTION_STRIDED_H 48