1 // Copyright (c) 2017-2024, 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 #include <ceed.h> 9 #include <ceed/backend.h> 10 #include <ceed/jit-tools.h> 11 #include <cuda.h> 12 #include <cuda_runtime.h> 13 #include <stdbool.h> 14 #include <stddef.h> 15 16 #include "../cuda/ceed-cuda-common.h" 17 #include "../cuda/ceed-cuda-compile.h" 18 #include "ceed-cuda-shared.h" 19 20 //------------------------------------------------------------------------------ 21 // Device initalization 22 //------------------------------------------------------------------------------ 23 int CeedInit_CudaInterp(CeedScalar *d_B, CeedInt P_1d, CeedInt Q_1d, CeedScalar **c_B); 24 int CeedInit_CudaGrad(CeedScalar *d_B, CeedScalar *d_G, CeedInt P_1d, CeedInt Q_1d, CeedScalar **c_B_ptr, CeedScalar **c_G_ptr); 25 int CeedInit_CudaCollocatedGrad(CeedScalar *d_B, CeedScalar *d_G, CeedInt P_1d, CeedInt Q_1d, CeedScalar **c_B_ptr, CeedScalar **c_G_ptr); 26 27 //------------------------------------------------------------------------------ 28 // Apply basis 29 //------------------------------------------------------------------------------ 30 int CeedBasisApplyTensor_Cuda_shared(CeedBasis basis, const CeedInt num_elem, CeedTransposeMode t_mode, CeedEvalMode eval_mode, CeedVector u, 31 CeedVector v) { 32 Ceed ceed; 33 Ceed_Cuda *ceed_Cuda; 34 CeedInt dim, num_comp; 35 const CeedScalar *d_u; 36 CeedScalar *d_v; 37 CeedBasis_Cuda_shared *data; 38 39 CeedCallBackend(CeedBasisGetCeed(basis, &ceed)); 40 CeedCallBackend(CeedGetData(ceed, &ceed_Cuda)); 41 CeedCallBackend(CeedBasisGetData(basis, &data)); 42 CeedCallBackend(CeedBasisGetDimension(basis, &dim)); 43 CeedCallBackend(CeedBasisGetNumComponents(basis, &num_comp)); 44 45 // Get read/write access to u, v 46 if (u != CEED_VECTOR_NONE) CeedCallBackend(CeedVectorGetArrayRead(u, CEED_MEM_DEVICE, &d_u)); 47 else CeedCheck(eval_mode == CEED_EVAL_WEIGHT, ceed, CEED_ERROR_BACKEND, "An input vector is required for this CeedEvalMode"); 48 CeedCallBackend(CeedVectorGetArrayWrite(v, CEED_MEM_DEVICE, &d_v)); 49 50 // Apply basis operation 51 switch (eval_mode) { 52 case CEED_EVAL_INTERP: { 53 CeedInt P_1d, Q_1d; 54 55 CeedCallBackend(CeedBasisGetNumNodes1D(basis, &P_1d)); 56 CeedCallBackend(CeedBasisGetNumQuadraturePoints1D(basis, &Q_1d)); 57 CeedInt thread_1d = CeedIntMax(Q_1d, P_1d); 58 59 CeedCallBackend(CeedInit_CudaInterp(data->d_interp_1d, P_1d, Q_1d, &data->c_B)); 60 void *interp_args[] = {(void *)&num_elem, &data->c_B, &d_u, &d_v}; 61 62 if (dim == 1) { 63 CeedInt elems_per_block = CeedIntMin(ceed_Cuda->device_prop.maxThreadsDim[2], CeedIntMax(512 / thread_1d, 64 1)); // avoid >512 total threads 65 CeedInt grid = num_elem / elems_per_block + ((num_elem / elems_per_block * elems_per_block < num_elem) ? 1 : 0); 66 CeedInt shared_mem = elems_per_block * thread_1d * sizeof(CeedScalar); 67 68 if (t_mode == CEED_TRANSPOSE) { 69 CeedCallBackend(CeedRunKernelDimShared_Cuda(ceed, data->InterpTranspose, grid, thread_1d, 1, elems_per_block, shared_mem, interp_args)); 70 } else { 71 CeedCallBackend(CeedRunKernelDimShared_Cuda(ceed, data->Interp, grid, thread_1d, 1, elems_per_block, shared_mem, interp_args)); 72 } 73 } else if (dim == 2) { 74 const CeedInt opt_elems[7] = {0, 32, 8, 6, 4, 2, 8}; 75 // elems_per_block must be at least 1 76 CeedInt elems_per_block = CeedIntMax(thread_1d < 7 ? opt_elems[thread_1d] / num_comp : 1, 1); 77 CeedInt grid = num_elem / elems_per_block + ((num_elem / elems_per_block * elems_per_block < num_elem) ? 1 : 0); 78 CeedInt shared_mem = elems_per_block * thread_1d * thread_1d * sizeof(CeedScalar); 79 80 if (t_mode == CEED_TRANSPOSE) { 81 CeedCallBackend( 82 CeedRunKernelDimShared_Cuda(ceed, data->InterpTranspose, grid, thread_1d, thread_1d, elems_per_block, shared_mem, interp_args)); 83 } else { 84 CeedCallBackend(CeedRunKernelDimShared_Cuda(ceed, data->Interp, grid, thread_1d, thread_1d, elems_per_block, shared_mem, interp_args)); 85 } 86 } else if (dim == 3) { 87 CeedInt elems_per_block = 1; 88 CeedInt grid = num_elem / elems_per_block + ((num_elem / elems_per_block * elems_per_block < num_elem) ? 1 : 0); 89 CeedInt shared_mem = elems_per_block * thread_1d * thread_1d * sizeof(CeedScalar); 90 91 if (t_mode == CEED_TRANSPOSE) { 92 CeedCallBackend( 93 CeedRunKernelDimShared_Cuda(ceed, data->InterpTranspose, grid, thread_1d, thread_1d, elems_per_block, shared_mem, interp_args)); 94 } else { 95 CeedCallBackend(CeedRunKernelDimShared_Cuda(ceed, data->Interp, grid, thread_1d, thread_1d, elems_per_block, shared_mem, interp_args)); 96 } 97 } 98 } break; 99 case CEED_EVAL_GRAD: { 100 CeedInt P_1d, Q_1d; 101 102 CeedCallBackend(CeedBasisGetNumNodes1D(basis, &P_1d)); 103 CeedCallBackend(CeedBasisGetNumQuadraturePoints1D(basis, &Q_1d)); 104 CeedInt thread_1d = CeedIntMax(Q_1d, P_1d); 105 106 if (data->d_collo_grad_1d) { 107 CeedCallBackend(CeedInit_CudaCollocatedGrad(data->d_interp_1d, data->d_collo_grad_1d, P_1d, Q_1d, &data->c_B, &data->c_G)); 108 } else { 109 CeedCallBackend(CeedInit_CudaGrad(data->d_interp_1d, data->d_grad_1d, P_1d, Q_1d, &data->c_B, &data->c_G)); 110 } 111 void *grad_args[] = {(void *)&num_elem, &data->c_B, &data->c_G, &d_u, &d_v}; 112 if (dim == 1) { 113 CeedInt elems_per_block = CeedIntMin(ceed_Cuda->device_prop.maxThreadsDim[2], CeedIntMax(512 / thread_1d, 114 1)); // avoid >512 total threads 115 CeedInt grid = num_elem / elems_per_block + ((num_elem / elems_per_block * elems_per_block < num_elem) ? 1 : 0); 116 CeedInt shared_mem = elems_per_block * thread_1d * sizeof(CeedScalar); 117 118 if (t_mode == CEED_TRANSPOSE) { 119 CeedCallBackend(CeedRunKernelDimShared_Cuda(ceed, data->GradTranspose, grid, thread_1d, 1, elems_per_block, shared_mem, grad_args)); 120 } else { 121 CeedCallBackend(CeedRunKernelDimShared_Cuda(ceed, data->Grad, grid, thread_1d, 1, elems_per_block, shared_mem, grad_args)); 122 } 123 } else if (dim == 2) { 124 const CeedInt opt_elems[7] = {0, 32, 8, 6, 4, 2, 8}; 125 // elems_per_block must be at least 1 126 CeedInt elems_per_block = CeedIntMax(thread_1d < 7 ? opt_elems[thread_1d] / num_comp : 1, 1); 127 CeedInt grid = num_elem / elems_per_block + ((num_elem / elems_per_block * elems_per_block < num_elem) ? 1 : 0); 128 CeedInt shared_mem = elems_per_block * thread_1d * thread_1d * sizeof(CeedScalar); 129 130 if (t_mode == CEED_TRANSPOSE) { 131 CeedCallBackend(CeedRunKernelDimShared_Cuda(ceed, data->GradTranspose, grid, thread_1d, thread_1d, elems_per_block, shared_mem, grad_args)); 132 } else { 133 CeedCallBackend(CeedRunKernelDimShared_Cuda(ceed, data->Grad, grid, thread_1d, thread_1d, elems_per_block, shared_mem, grad_args)); 134 } 135 } else if (dim == 3) { 136 CeedInt elems_per_block = 1; 137 CeedInt grid = num_elem / elems_per_block + ((num_elem / elems_per_block * elems_per_block < num_elem) ? 1 : 0); 138 CeedInt shared_mem = elems_per_block * thread_1d * thread_1d * sizeof(CeedScalar); 139 140 if (t_mode == CEED_TRANSPOSE) { 141 CeedCallBackend(CeedRunKernelDimShared_Cuda(ceed, data->GradTranspose, grid, thread_1d, thread_1d, elems_per_block, shared_mem, grad_args)); 142 } else { 143 CeedCallBackend(CeedRunKernelDimShared_Cuda(ceed, data->Grad, grid, thread_1d, thread_1d, elems_per_block, shared_mem, grad_args)); 144 } 145 } 146 } break; 147 case CEED_EVAL_WEIGHT: { 148 CeedInt Q_1d; 149 CeedInt block_size = 32; 150 151 CeedCheck(data->d_q_weight_1d, ceed, CEED_ERROR_BACKEND, "%s not supported; q_weights_1d not set", CeedEvalModes[eval_mode]); 152 CeedCallBackend(CeedBasisGetNumQuadraturePoints1D(basis, &Q_1d)); 153 void *weight_args[] = {(void *)&num_elem, (void *)&data->d_q_weight_1d, &d_v}; 154 if (dim == 1) { 155 const CeedInt elems_per_block = block_size / Q_1d; 156 const CeedInt grid_size = num_elem / elems_per_block + ((num_elem / elems_per_block * elems_per_block < num_elem) ? 1 : 0); 157 158 CeedCallBackend(CeedRunKernelDim_Cuda(ceed, data->Weight, grid_size, Q_1d, elems_per_block, 1, weight_args)); 159 } else if (dim == 2) { 160 const CeedInt opt_elems = block_size / (Q_1d * Q_1d); 161 const CeedInt elems_per_block = opt_elems > 0 ? opt_elems : 1; 162 const CeedInt grid_size = num_elem / elems_per_block + ((num_elem / elems_per_block * elems_per_block < num_elem) ? 1 : 0); 163 164 CeedCallBackend(CeedRunKernelDim_Cuda(ceed, data->Weight, grid_size, Q_1d, Q_1d, elems_per_block, weight_args)); 165 } else if (dim == 3) { 166 const CeedInt opt_elems = block_size / (Q_1d * Q_1d); 167 const CeedInt elems_per_block = opt_elems > 0 ? opt_elems : 1; 168 const CeedInt grid_size = num_elem / elems_per_block + ((num_elem / elems_per_block * elems_per_block < num_elem) ? 1 : 0); 169 170 CeedCallBackend(CeedRunKernelDim_Cuda(ceed, data->Weight, grid_size, Q_1d, Q_1d, elems_per_block, weight_args)); 171 } 172 } break; 173 case CEED_EVAL_NONE: /* handled separately below */ 174 break; 175 // LCOV_EXCL_START 176 case CEED_EVAL_DIV: 177 case CEED_EVAL_CURL: 178 return CeedError(ceed, CEED_ERROR_BACKEND, "%s not supported", CeedEvalModes[eval_mode]); 179 // LCOV_EXCL_STOP 180 } 181 182 // Restore vectors, cover CEED_EVAL_NONE 183 CeedCallBackend(CeedVectorRestoreArray(v, &d_v)); 184 if (eval_mode == CEED_EVAL_NONE) CeedCallBackend(CeedVectorSetArray(v, CEED_MEM_DEVICE, CEED_COPY_VALUES, (CeedScalar *)d_u)); 185 if (eval_mode != CEED_EVAL_WEIGHT) CeedCallBackend(CeedVectorRestoreArrayRead(u, &d_u)); 186 return CEED_ERROR_SUCCESS; 187 } 188 189 //------------------------------------------------------------------------------ 190 // Basis apply - tensor AtPoints 191 //------------------------------------------------------------------------------ 192 int CeedBasisApplyAtPoints_Cuda_shared(CeedBasis basis, const CeedInt num_elem, const CeedInt *num_points, CeedTransposeMode t_mode, 193 CeedEvalMode eval_mode, CeedVector x_ref, CeedVector u, CeedVector v) { 194 Ceed ceed; 195 CeedInt Q_1d, dim, max_num_points = num_points[0]; 196 const CeedInt is_transpose = t_mode == CEED_TRANSPOSE; 197 const int max_block_size = 32; 198 const CeedScalar *d_x, *d_u; 199 CeedScalar *d_v; 200 CeedBasis_Cuda_shared *data; 201 202 CeedCallBackend(CeedBasisGetCeed(basis, &ceed)); 203 CeedCallBackend(CeedBasisGetData(basis, &data)); 204 CeedCallBackend(CeedBasisGetNumQuadraturePoints1D(basis, &Q_1d)); 205 CeedCallBackend(CeedBasisGetDimension(basis, &dim)); 206 207 // Check uniform number of points per elem 208 for (CeedInt i = 1; i < num_elem; i++) { 209 CeedCheck(max_num_points == num_points[i], ceed, CEED_ERROR_BACKEND, 210 "BasisApplyAtPoints only supported for the same number of points in each element"); 211 } 212 213 // Weight handled separately 214 if (eval_mode == CEED_EVAL_WEIGHT) { 215 CeedCall(CeedVectorSetValue(v, 1.0)); 216 return CEED_ERROR_SUCCESS; 217 } 218 219 // Build kernels if needed 220 if (data->num_points != max_num_points) { 221 CeedInt P_1d; 222 223 CeedCallBackend(CeedBasisGetNumNodes1D(basis, &P_1d)); 224 data->num_points = max_num_points; 225 226 // -- Create interp matrix to Chebyshev coefficients 227 if (!data->d_chebyshev_interp_1d) { 228 CeedSize interp_bytes; 229 CeedScalar *chebyshev_interp_1d; 230 231 interp_bytes = P_1d * Q_1d * sizeof(CeedScalar); 232 CeedCallBackend(CeedCalloc(P_1d * Q_1d, &chebyshev_interp_1d)); 233 CeedCall(CeedBasisGetChebyshevInterp1D(basis, chebyshev_interp_1d)); 234 CeedCallCuda(ceed, cudaMalloc((void **)&data->d_chebyshev_interp_1d, interp_bytes)); 235 CeedCallCuda(ceed, cudaMemcpy(data->d_chebyshev_interp_1d, chebyshev_interp_1d, interp_bytes, cudaMemcpyHostToDevice)); 236 CeedCallBackend(CeedFree(&chebyshev_interp_1d)); 237 } 238 239 // -- Compile kernels 240 char *basis_kernel_source; 241 const char *basis_kernel_path; 242 CeedInt num_comp; 243 244 if (data->moduleAtPoints) CeedCallCuda(ceed, cuModuleUnload(data->moduleAtPoints)); 245 CeedCallBackend(CeedBasisGetNumComponents(basis, &num_comp)); 246 CeedCallBackend(CeedGetJitAbsolutePath(ceed, "ceed/jit-source/cuda/cuda-ref-basis-tensor-at-points.h", &basis_kernel_path)); 247 CeedDebug256(ceed, CEED_DEBUG_COLOR_SUCCESS, "----- Loading Basis Kernel Source -----\n"); 248 CeedCallBackend(CeedLoadSourceToBuffer(ceed, basis_kernel_path, &basis_kernel_source)); 249 CeedDebug256(ceed, CEED_DEBUG_COLOR_SUCCESS, "----- Loading Basis Kernel Source Complete! -----\n"); 250 CeedCallBackend(CeedCompile_Cuda(ceed, basis_kernel_source, &data->moduleAtPoints, 9, "BASIS_Q_1D", Q_1d, "BASIS_P_1D", P_1d, "BASIS_BUF_LEN", 251 Q_1d * CeedIntPow(Q_1d > P_1d ? Q_1d : P_1d, dim - 1), "BASIS_DIM", dim, "BASIS_NUM_COMP", num_comp, 252 "BASIS_NUM_NODES", CeedIntPow(P_1d, dim), "BASIS_NUM_QPTS", CeedIntPow(Q_1d, dim), "BASIS_NUM_PTS", 253 max_num_points, "POINTS_BUFF_LEN", CeedIntPow(Q_1d, dim - 1))); 254 CeedCallBackend(CeedGetKernel_Cuda(ceed, data->moduleAtPoints, "InterpAtPoints", &data->InterpAtPoints)); 255 CeedCallBackend(CeedGetKernel_Cuda(ceed, data->moduleAtPoints, "GradAtPoints", &data->GradAtPoints)); 256 CeedCallBackend(CeedFree(&basis_kernel_path)); 257 CeedCallBackend(CeedFree(&basis_kernel_source)); 258 } 259 260 // Get read/write access to u, v 261 CeedCallBackend(CeedVectorGetArrayRead(x_ref, CEED_MEM_DEVICE, &d_x)); 262 if (u != CEED_VECTOR_NONE) CeedCallBackend(CeedVectorGetArrayRead(u, CEED_MEM_DEVICE, &d_u)); 263 else CeedCheck(eval_mode == CEED_EVAL_WEIGHT, ceed, CEED_ERROR_BACKEND, "An input vector is required for this CeedEvalMode"); 264 CeedCallBackend(CeedVectorGetArrayWrite(v, CEED_MEM_DEVICE, &d_v)); 265 266 // Clear v for transpose operation 267 if (is_transpose) { 268 CeedSize length; 269 270 CeedCallBackend(CeedVectorGetLength(v, &length)); 271 CeedCallCuda(ceed, cudaMemset(d_v, 0, length * sizeof(CeedScalar))); 272 } 273 274 // Basis action 275 switch (eval_mode) { 276 case CEED_EVAL_INTERP: { 277 void *interp_args[] = {(void *)&num_elem, (void *)&is_transpose, &data->d_chebyshev_interp_1d, &d_x, &d_u, &d_v}; 278 const CeedInt block_size = CeedIntMin(CeedIntPow(Q_1d, dim), max_block_size); 279 280 CeedCallBackend(CeedRunKernel_Cuda(ceed, data->InterpAtPoints, num_elem, block_size, interp_args)); 281 } break; 282 case CEED_EVAL_GRAD: { 283 void *grad_args[] = {(void *)&num_elem, (void *)&is_transpose, &data->d_chebyshev_interp_1d, &d_x, &d_u, &d_v}; 284 const CeedInt block_size = CeedIntMin(CeedIntPow(Q_1d, dim), max_block_size); 285 286 CeedCallBackend(CeedRunKernel_Cuda(ceed, data->GradAtPoints, num_elem, block_size, grad_args)); 287 } break; 288 case CEED_EVAL_WEIGHT: 289 case CEED_EVAL_NONE: /* handled separately below */ 290 break; 291 // LCOV_EXCL_START 292 case CEED_EVAL_DIV: 293 case CEED_EVAL_CURL: 294 return CeedError(ceed, CEED_ERROR_BACKEND, "%s not supported", CeedEvalModes[eval_mode]); 295 // LCOV_EXCL_STOP 296 } 297 298 // Restore vectors, cover CEED_EVAL_NONE 299 CeedCallBackend(CeedVectorRestoreArrayRead(x_ref, &d_x)); 300 CeedCallBackend(CeedVectorRestoreArray(v, &d_v)); 301 if (eval_mode == CEED_EVAL_NONE) CeedCallBackend(CeedVectorSetArray(v, CEED_MEM_DEVICE, CEED_COPY_VALUES, (CeedScalar *)d_u)); 302 if (eval_mode != CEED_EVAL_WEIGHT) CeedCallBackend(CeedVectorRestoreArrayRead(u, &d_u)); 303 return CEED_ERROR_SUCCESS; 304 } 305 306 //------------------------------------------------------------------------------ 307 // Destroy basis 308 //------------------------------------------------------------------------------ 309 static int CeedBasisDestroy_Cuda_shared(CeedBasis basis) { 310 Ceed ceed; 311 CeedBasis_Cuda_shared *data; 312 313 CeedCallBackend(CeedBasisGetCeed(basis, &ceed)); 314 CeedCallBackend(CeedBasisGetData(basis, &data)); 315 CeedCallCuda(ceed, cuModuleUnload(data->module)); 316 if (data->moduleAtPoints) CeedCallCuda(ceed, cuModuleUnload(data->moduleAtPoints)); 317 if (data->d_q_weight_1d) CeedCallCuda(ceed, cudaFree(data->d_q_weight_1d)); 318 CeedCallCuda(ceed, cudaFree(data->d_interp_1d)); 319 CeedCallCuda(ceed, cudaFree(data->d_grad_1d)); 320 CeedCallCuda(ceed, cudaFree(data->d_collo_grad_1d)); 321 CeedCallCuda(ceed, cudaFree(data->d_chebyshev_interp_1d)); 322 CeedCallBackend(CeedFree(&data)); 323 return CEED_ERROR_SUCCESS; 324 } 325 326 //------------------------------------------------------------------------------ 327 // Create tensor basis 328 //------------------------------------------------------------------------------ 329 int CeedBasisCreateTensorH1_Cuda_shared(CeedInt dim, CeedInt P_1d, CeedInt Q_1d, const CeedScalar *interp_1d, const CeedScalar *grad_1d, 330 const CeedScalar *q_ref_1d, const CeedScalar *q_weight_1d, CeedBasis basis) { 331 Ceed ceed; 332 char *basis_kernel_source; 333 const char *basis_kernel_path; 334 CeedInt num_comp; 335 const CeedInt q_bytes = Q_1d * sizeof(CeedScalar); 336 const CeedInt interp_bytes = q_bytes * P_1d; 337 CeedBasis_Cuda_shared *data; 338 339 CeedCallBackend(CeedBasisGetCeed(basis, &ceed)); 340 CeedCallBackend(CeedCalloc(1, &data)); 341 342 // Copy basis data to GPU 343 if (q_weight_1d) { 344 CeedCallCuda(ceed, cudaMalloc((void **)&data->d_q_weight_1d, q_bytes)); 345 CeedCallCuda(ceed, cudaMemcpy(data->d_q_weight_1d, q_weight_1d, q_bytes, cudaMemcpyHostToDevice)); 346 } 347 CeedCallCuda(ceed, cudaMalloc((void **)&data->d_interp_1d, interp_bytes)); 348 CeedCallCuda(ceed, cudaMemcpy(data->d_interp_1d, interp_1d, interp_bytes, cudaMemcpyHostToDevice)); 349 CeedCallCuda(ceed, cudaMalloc((void **)&data->d_grad_1d, interp_bytes)); 350 CeedCallCuda(ceed, cudaMemcpy(data->d_grad_1d, grad_1d, interp_bytes, cudaMemcpyHostToDevice)); 351 352 // Compute collocated gradient and copy to GPU 353 data->d_collo_grad_1d = NULL; 354 bool has_collocated_grad = dim == 3 && Q_1d >= P_1d; 355 356 if (has_collocated_grad) { 357 CeedScalar *collo_grad_1d; 358 359 CeedCallBackend(CeedMalloc(Q_1d * Q_1d, &collo_grad_1d)); 360 CeedCallBackend(CeedBasisGetCollocatedGrad(basis, collo_grad_1d)); 361 CeedCallCuda(ceed, cudaMalloc((void **)&data->d_collo_grad_1d, q_bytes * Q_1d)); 362 CeedCallCuda(ceed, cudaMemcpy(data->d_collo_grad_1d, collo_grad_1d, q_bytes * Q_1d, cudaMemcpyHostToDevice)); 363 CeedCallBackend(CeedFree(&collo_grad_1d)); 364 } 365 366 // Compile basis kernels 367 CeedCallBackend(CeedBasisGetNumComponents(basis, &num_comp)); 368 CeedCallBackend(CeedGetJitAbsolutePath(ceed, "ceed/jit-source/cuda/cuda-shared-basis-tensor.h", &basis_kernel_path)); 369 CeedDebug256(ceed, CEED_DEBUG_COLOR_SUCCESS, "----- Loading Basis Kernel Source -----\n"); 370 CeedCallBackend(CeedLoadSourceToBuffer(ceed, basis_kernel_path, &basis_kernel_source)); 371 CeedDebug256(ceed, CEED_DEBUG_COLOR_SUCCESS, "----- Loading Basis Kernel Source Complete -----\n"); 372 CeedCallBackend(CeedCompile_Cuda(ceed, basis_kernel_source, &data->module, 8, "BASIS_Q_1D", Q_1d, "BASIS_P_1D", P_1d, "T_1D", 373 CeedIntMax(Q_1d, P_1d), "BASIS_DIM", dim, "BASIS_NUM_COMP", num_comp, "BASIS_NUM_NODES", CeedIntPow(P_1d, dim), 374 "BASIS_NUM_QPTS", CeedIntPow(Q_1d, dim), "BASIS_HAS_COLLOCATED_GRAD", has_collocated_grad)); 375 CeedCallBackend(CeedGetKernel_Cuda(ceed, data->module, "Interp", &data->Interp)); 376 CeedCallBackend(CeedGetKernel_Cuda(ceed, data->module, "InterpTranspose", &data->InterpTranspose)); 377 CeedCallBackend(CeedGetKernel_Cuda(ceed, data->module, "Grad", &data->Grad)); 378 CeedCallBackend(CeedGetKernel_Cuda(ceed, data->module, "GradTranspose", &data->GradTranspose)); 379 CeedCallBackend(CeedGetKernel_Cuda(ceed, data->module, "Weight", &data->Weight)); 380 CeedCallBackend(CeedFree(&basis_kernel_path)); 381 CeedCallBackend(CeedFree(&basis_kernel_source)); 382 383 CeedCallBackend(CeedBasisSetData(basis, data)); 384 385 // Register backend functions 386 CeedCallBackend(CeedSetBackendFunction(ceed, "Basis", basis, "Apply", CeedBasisApplyTensor_Cuda_shared)); 387 CeedCallBackend(CeedSetBackendFunction(ceed, "Basis", basis, "ApplyAtPoints", CeedBasisApplyAtPoints_Cuda_shared)); 388 CeedCallBackend(CeedSetBackendFunction(ceed, "Basis", basis, "Destroy", CeedBasisDestroy_Cuda_shared)); 389 return CEED_ERROR_SUCCESS; 390 } 391 392 //------------------------------------------------------------------------------ 393