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 14 #include "../cuda/ceed-cuda-common.h" 15 #include "../cuda/ceed-cuda-compile.h" 16 #include "ceed-cuda-ref.h" 17 18 //------------------------------------------------------------------------------ 19 // Basis apply - tensor 20 //------------------------------------------------------------------------------ 21 int CeedBasisApply_Cuda(CeedBasis basis, const CeedInt num_elem, CeedTransposeMode t_mode, CeedEvalMode eval_mode, CeedVector u, CeedVector v) { 22 Ceed ceed; 23 CeedInt Q_1d, dim; 24 const CeedInt is_transpose = t_mode == CEED_TRANSPOSE; 25 const int max_block_size = 32; 26 const CeedScalar *d_u; 27 CeedScalar *d_v; 28 CeedBasis_Cuda *data; 29 30 CeedCallBackend(CeedBasisGetCeed(basis, &ceed)); 31 CeedCallBackend(CeedBasisGetData(basis, &data)); 32 33 // Get read/write access to u, v 34 if (u != CEED_VECTOR_NONE) CeedCallBackend(CeedVectorGetArrayRead(u, CEED_MEM_DEVICE, &d_u)); 35 else CeedCheck(eval_mode == CEED_EVAL_WEIGHT, ceed, CEED_ERROR_BACKEND, "An input vector is required for this CeedEvalMode"); 36 CeedCallBackend(CeedVectorGetArrayWrite(v, CEED_MEM_DEVICE, &d_v)); 37 38 // Clear v for transpose operation 39 if (is_transpose) { 40 CeedSize length; 41 42 CeedCallBackend(CeedVectorGetLength(v, &length)); 43 CeedCallCuda(ceed, cudaMemset(d_v, 0, length * sizeof(CeedScalar))); 44 } 45 CeedCallBackend(CeedBasisGetNumQuadraturePoints1D(basis, &Q_1d)); 46 CeedCallBackend(CeedBasisGetDimension(basis, &dim)); 47 48 // Basis action 49 switch (eval_mode) { 50 case CEED_EVAL_INTERP: { 51 void *interp_args[] = {(void *)&num_elem, (void *)&is_transpose, &data->d_interp_1d, &d_u, &d_v}; 52 const CeedInt block_size = CeedIntMin(CeedIntPow(Q_1d, dim), max_block_size); 53 54 CeedCallBackend(CeedRunKernel_Cuda(ceed, data->Interp, num_elem, block_size, interp_args)); 55 } break; 56 case CEED_EVAL_GRAD: { 57 void *grad_args[] = {(void *)&num_elem, (void *)&is_transpose, &data->d_interp_1d, &data->d_grad_1d, &d_u, &d_v}; 58 const CeedInt block_size = max_block_size; 59 60 CeedCallBackend(CeedRunKernel_Cuda(ceed, data->Grad, num_elem, block_size, grad_args)); 61 } break; 62 case CEED_EVAL_WEIGHT: { 63 CeedCheck(data->d_q_weight_1d, ceed, CEED_ERROR_BACKEND, "%s not supported; q_weights_1d not set", CeedEvalModes[eval_mode]); 64 void *weight_args[] = {(void *)&num_elem, (void *)&data->d_q_weight_1d, &d_v}; 65 const int block_size_x = Q_1d; 66 const int block_size_y = dim >= 2 ? Q_1d : 1; 67 68 CeedCallBackend(CeedRunKernelDim_Cuda(ceed, data->Weight, num_elem, block_size_x, block_size_y, 1, weight_args)); 69 } break; 70 case CEED_EVAL_NONE: /* handled separately below */ 71 break; 72 // LCOV_EXCL_START 73 case CEED_EVAL_DIV: 74 case CEED_EVAL_CURL: 75 return CeedError(ceed, CEED_ERROR_BACKEND, "%s not supported", CeedEvalModes[eval_mode]); 76 // LCOV_EXCL_STOP 77 } 78 79 // Restore vectors, cover CEED_EVAL_NONE 80 CeedCallBackend(CeedVectorRestoreArray(v, &d_v)); 81 if (eval_mode == CEED_EVAL_NONE) CeedCallBackend(CeedVectorSetArray(v, CEED_MEM_DEVICE, CEED_COPY_VALUES, (CeedScalar *)d_u)); 82 if (eval_mode != CEED_EVAL_WEIGHT) CeedCallBackend(CeedVectorRestoreArrayRead(u, &d_u)); 83 return CEED_ERROR_SUCCESS; 84 } 85 86 //------------------------------------------------------------------------------ 87 // Basis apply - non-tensor 88 //------------------------------------------------------------------------------ 89 int CeedBasisApplyNonTensor_Cuda(CeedBasis basis, const CeedInt num_elem, CeedTransposeMode t_mode, CeedEvalMode eval_mode, CeedVector u, 90 CeedVector v) { 91 Ceed ceed; 92 CeedInt num_nodes, num_qpts; 93 const CeedInt is_transpose = t_mode == CEED_TRANSPOSE; 94 const int elems_per_block = 1; 95 const int grid = CeedDivUpInt(num_elem, elems_per_block); 96 const CeedScalar *d_u; 97 CeedScalar *d_v; 98 CeedBasisNonTensor_Cuda *data; 99 100 CeedCallBackend(CeedBasisGetCeed(basis, &ceed)); 101 CeedCallBackend(CeedBasisGetData(basis, &data)); 102 CeedCallBackend(CeedBasisGetNumQuadraturePoints(basis, &num_qpts)); 103 CeedCallBackend(CeedBasisGetNumNodes(basis, &num_nodes)); 104 105 // Get read/write access to u, v 106 if (u != CEED_VECTOR_NONE) CeedCallBackend(CeedVectorGetArrayRead(u, CEED_MEM_DEVICE, &d_u)); 107 else CeedCheck(eval_mode == CEED_EVAL_WEIGHT, ceed, CEED_ERROR_BACKEND, "An input vector is required for this CeedEvalMode"); 108 CeedCallBackend(CeedVectorGetArrayWrite(v, CEED_MEM_DEVICE, &d_v)); 109 110 // Clear v for transpose operation 111 if (is_transpose) { 112 CeedSize length; 113 114 CeedCallBackend(CeedVectorGetLength(v, &length)); 115 CeedCallCuda(ceed, cudaMemset(d_v, 0, length * sizeof(CeedScalar))); 116 } 117 118 // Apply basis operation 119 switch (eval_mode) { 120 case CEED_EVAL_INTERP: { 121 void *interp_args[] = {(void *)&num_elem, &data->d_interp, &d_u, &d_v}; 122 const int block_size_x = is_transpose ? num_nodes : num_qpts; 123 124 if (is_transpose) { 125 CeedCallBackend(CeedRunKernelDim_Cuda(ceed, data->InterpTranspose, grid, block_size_x, 1, elems_per_block, interp_args)); 126 } else { 127 CeedCallBackend(CeedRunKernelDim_Cuda(ceed, data->Interp, grid, block_size_x, 1, elems_per_block, interp_args)); 128 } 129 } break; 130 case CEED_EVAL_GRAD: { 131 void *grad_args[] = {(void *)&num_elem, &data->d_grad, &d_u, &d_v}; 132 const int block_size_x = is_transpose ? num_nodes : num_qpts; 133 134 if (is_transpose) { 135 CeedCallBackend(CeedRunKernelDim_Cuda(ceed, data->DerivTranspose, grid, block_size_x, 1, elems_per_block, grad_args)); 136 } else { 137 CeedCallBackend(CeedRunKernelDim_Cuda(ceed, data->Deriv, grid, block_size_x, 1, elems_per_block, grad_args)); 138 } 139 } break; 140 case CEED_EVAL_DIV: { 141 void *div_args[] = {(void *)&num_elem, &data->d_div, &d_u, &d_v}; 142 const int block_size_x = is_transpose ? num_nodes : num_qpts; 143 144 if (is_transpose) { 145 CeedCallBackend(CeedRunKernelDim_Cuda(ceed, data->DerivTranspose, grid, block_size_x, 1, elems_per_block, div_args)); 146 } else { 147 CeedCallBackend(CeedRunKernelDim_Cuda(ceed, data->Deriv, grid, block_size_x, 1, elems_per_block, div_args)); 148 } 149 } break; 150 case CEED_EVAL_CURL: { 151 void *curl_args[] = {(void *)&num_elem, &data->d_curl, &d_u, &d_v}; 152 const int block_size_x = is_transpose ? num_nodes : num_qpts; 153 154 if (is_transpose) { 155 CeedCallBackend(CeedRunKernelDim_Cuda(ceed, data->DerivTranspose, grid, block_size_x, 1, elems_per_block, curl_args)); 156 } else { 157 CeedCallBackend(CeedRunKernelDim_Cuda(ceed, data->Deriv, grid, block_size_x, 1, elems_per_block, curl_args)); 158 } 159 } break; 160 case CEED_EVAL_WEIGHT: { 161 CeedCheck(data->d_q_weight, ceed, CEED_ERROR_BACKEND, "%s not supported; q_weights not set", CeedEvalModes[eval_mode]); 162 void *weight_args[] = {(void *)&num_elem, (void *)&data->d_q_weight, &d_v}; 163 164 CeedCallBackend(CeedRunKernelDim_Cuda(ceed, data->Weight, grid, num_qpts, 1, elems_per_block, weight_args)); 165 } break; 166 case CEED_EVAL_NONE: /* handled separately below */ 167 break; 168 } 169 170 // Restore vectors, cover CEED_EVAL_NONE 171 CeedCallBackend(CeedVectorRestoreArray(v, &d_v)); 172 if (eval_mode == CEED_EVAL_NONE) CeedCallBackend(CeedVectorSetArray(v, CEED_MEM_DEVICE, CEED_COPY_VALUES, (CeedScalar *)d_u)); 173 if (eval_mode != CEED_EVAL_WEIGHT) CeedCallBackend(CeedVectorRestoreArrayRead(u, &d_u)); 174 return CEED_ERROR_SUCCESS; 175 } 176 177 //------------------------------------------------------------------------------ 178 // Destroy tensor basis 179 //------------------------------------------------------------------------------ 180 static int CeedBasisDestroy_Cuda(CeedBasis basis) { 181 Ceed ceed; 182 CeedBasis_Cuda *data; 183 184 CeedCallBackend(CeedBasisGetCeed(basis, &ceed)); 185 CeedCallBackend(CeedBasisGetData(basis, &data)); 186 CeedCallCuda(ceed, cuModuleUnload(data->module)); 187 if (data->d_q_weight_1d) CeedCallCuda(ceed, cudaFree(data->d_q_weight_1d)); 188 CeedCallCuda(ceed, cudaFree(data->d_interp_1d)); 189 CeedCallCuda(ceed, cudaFree(data->d_grad_1d)); 190 CeedCallBackend(CeedFree(&data)); 191 return CEED_ERROR_SUCCESS; 192 } 193 194 //------------------------------------------------------------------------------ 195 // Destroy non-tensor basis 196 //------------------------------------------------------------------------------ 197 static int CeedBasisDestroyNonTensor_Cuda(CeedBasis basis) { 198 Ceed ceed; 199 CeedBasisNonTensor_Cuda *data; 200 201 CeedCallBackend(CeedBasisGetCeed(basis, &ceed)); 202 CeedCallBackend(CeedBasisGetData(basis, &data)); 203 CeedCallCuda(ceed, cuModuleUnload(data->module)); 204 if (data->d_q_weight) CeedCallCuda(ceed, cudaFree(data->d_q_weight)); 205 CeedCallCuda(ceed, cudaFree(data->d_interp)); 206 CeedCallCuda(ceed, cudaFree(data->d_grad)); 207 CeedCallCuda(ceed, cudaFree(data->d_div)); 208 CeedCallCuda(ceed, cudaFree(data->d_curl)); 209 CeedCallBackend(CeedFree(&data)); 210 return CEED_ERROR_SUCCESS; 211 } 212 213 //------------------------------------------------------------------------------ 214 // Create tensor 215 //------------------------------------------------------------------------------ 216 int CeedBasisCreateTensorH1_Cuda(CeedInt dim, CeedInt P_1d, CeedInt Q_1d, const CeedScalar *interp_1d, const CeedScalar *grad_1d, 217 const CeedScalar *q_ref_1d, const CeedScalar *q_weight_1d, CeedBasis basis) { 218 Ceed ceed; 219 char *basis_kernel_source; 220 const char *basis_kernel_path; 221 CeedInt num_comp; 222 const CeedInt q_bytes = Q_1d * sizeof(CeedScalar); 223 const CeedInt interp_bytes = q_bytes * P_1d; 224 CeedBasis_Cuda *data; 225 226 CeedCallBackend(CeedBasisGetCeed(basis, &ceed)); 227 CeedCallBackend(CeedCalloc(1, &data)); 228 229 // Copy data to GPU 230 if (q_weight_1d) { 231 CeedCallCuda(ceed, cudaMalloc((void **)&data->d_q_weight_1d, q_bytes)); 232 CeedCallCuda(ceed, cudaMemcpy(data->d_q_weight_1d, q_weight_1d, q_bytes, cudaMemcpyHostToDevice)); 233 } 234 CeedCallCuda(ceed, cudaMalloc((void **)&data->d_interp_1d, interp_bytes)); 235 CeedCallCuda(ceed, cudaMemcpy(data->d_interp_1d, interp_1d, interp_bytes, cudaMemcpyHostToDevice)); 236 CeedCallCuda(ceed, cudaMalloc((void **)&data->d_grad_1d, interp_bytes)); 237 CeedCallCuda(ceed, cudaMemcpy(data->d_grad_1d, grad_1d, interp_bytes, cudaMemcpyHostToDevice)); 238 239 // Compile basis kernels 240 CeedCallBackend(CeedBasisGetNumComponents(basis, &num_comp)); 241 CeedCallBackend(CeedGetJitAbsolutePath(ceed, "ceed/jit-source/cuda/cuda-ref-basis-tensor.h", &basis_kernel_path)); 242 CeedDebug256(ceed, CEED_DEBUG_COLOR_SUCCESS, "----- Loading Basis Kernel Source -----\n"); 243 CeedCallBackend(CeedLoadSourceToBuffer(ceed, basis_kernel_path, &basis_kernel_source)); 244 CeedDebug256(ceed, CEED_DEBUG_COLOR_SUCCESS, "----- Loading Basis Kernel Source Complete! -----\n"); 245 CeedCallBackend(CeedCompile_Cuda(ceed, basis_kernel_source, &data->module, 7, "BASIS_Q_1D", Q_1d, "BASIS_P_1D", P_1d, "BASIS_BUF_LEN", 246 num_comp * CeedIntPow(Q_1d > P_1d ? Q_1d : P_1d, dim), "BASIS_DIM", dim, "BASIS_NUM_COMP", num_comp, 247 "BASIS_NUM_NODES", CeedIntPow(P_1d, dim), "BASIS_NUM_QPTS", CeedIntPow(Q_1d, dim))); 248 CeedCallBackend(CeedGetKernel_Cuda(ceed, data->module, "Interp", &data->Interp)); 249 CeedCallBackend(CeedGetKernel_Cuda(ceed, data->module, "Grad", &data->Grad)); 250 CeedCallBackend(CeedGetKernel_Cuda(ceed, data->module, "Weight", &data->Weight)); 251 CeedCallBackend(CeedFree(&basis_kernel_path)); 252 CeedCallBackend(CeedFree(&basis_kernel_source)); 253 254 CeedCallBackend(CeedBasisSetData(basis, data)); 255 256 // Register backend functions 257 CeedCallBackend(CeedSetBackendFunction(ceed, "Basis", basis, "Apply", CeedBasisApply_Cuda)); 258 CeedCallBackend(CeedSetBackendFunction(ceed, "Basis", basis, "Destroy", CeedBasisDestroy_Cuda)); 259 return CEED_ERROR_SUCCESS; 260 } 261 262 //------------------------------------------------------------------------------ 263 // Create non-tensor H^1 264 //------------------------------------------------------------------------------ 265 int CeedBasisCreateH1_Cuda(CeedElemTopology topo, CeedInt dim, CeedInt num_nodes, CeedInt num_qpts, const CeedScalar *interp, const CeedScalar *grad, 266 const CeedScalar *q_ref, const CeedScalar *q_weight, CeedBasis basis) { 267 Ceed ceed; 268 char *basis_kernel_source; 269 const char *basis_kernel_path; 270 CeedInt num_comp, q_comp_interp, q_comp_grad; 271 const CeedInt q_bytes = num_qpts * sizeof(CeedScalar); 272 CeedBasisNonTensor_Cuda *data; 273 274 CeedCallBackend(CeedBasisGetCeed(basis, &ceed)); 275 CeedCallBackend(CeedCalloc(1, &data)); 276 277 // Copy basis data to GPU 278 CeedCallBackend(CeedBasisGetNumQuadratureComponents(basis, CEED_EVAL_INTERP, &q_comp_interp)); 279 CeedCallBackend(CeedBasisGetNumQuadratureComponents(basis, CEED_EVAL_GRAD, &q_comp_grad)); 280 if (q_weight) { 281 CeedCallCuda(ceed, cudaMalloc((void **)&data->d_q_weight, q_bytes)); 282 CeedCallCuda(ceed, cudaMemcpy(data->d_q_weight, q_weight, q_bytes, cudaMemcpyHostToDevice)); 283 } 284 if (interp) { 285 const CeedInt interp_bytes = q_bytes * num_nodes * q_comp_interp; 286 287 CeedCallCuda(ceed, cudaMalloc((void **)&data->d_interp, interp_bytes)); 288 CeedCallCuda(ceed, cudaMemcpy(data->d_interp, interp, interp_bytes, cudaMemcpyHostToDevice)); 289 } 290 if (grad) { 291 const CeedInt grad_bytes = q_bytes * num_nodes * q_comp_grad; 292 293 CeedCallCuda(ceed, cudaMalloc((void **)&data->d_grad, grad_bytes)); 294 CeedCallCuda(ceed, cudaMemcpy(data->d_grad, grad, grad_bytes, cudaMemcpyHostToDevice)); 295 } 296 297 // Compile basis kernels 298 CeedCallBackend(CeedBasisGetNumComponents(basis, &num_comp)); 299 CeedCallBackend(CeedGetJitAbsolutePath(ceed, "ceed/jit-source/cuda/cuda-ref-basis-nontensor.h", &basis_kernel_path)); 300 CeedDebug256(ceed, CEED_DEBUG_COLOR_SUCCESS, "----- Loading Basis Kernel Source -----\n"); 301 CeedCallBackend(CeedLoadSourceToBuffer(ceed, basis_kernel_path, &basis_kernel_source)); 302 CeedDebug256(ceed, CEED_DEBUG_COLOR_SUCCESS, "----- Loading Basis Kernel Source Complete! -----\n"); 303 CeedCallBackend(CeedCompile_Cuda(ceed, basis_kernel_source, &data->module, 5, "BASIS_Q", num_qpts, "BASIS_P", num_nodes, "BASIS_Q_COMP_INTERP", 304 q_comp_interp, "BASIS_Q_COMP_DERIV", q_comp_grad, "BASIS_NUM_COMP", num_comp)); 305 CeedCallBackend(CeedGetKernel_Cuda(ceed, data->module, "Interp", &data->Interp)); 306 CeedCallBackend(CeedGetKernel_Cuda(ceed, data->module, "InterpTranspose", &data->InterpTranspose)); 307 CeedCallBackend(CeedGetKernel_Cuda(ceed, data->module, "Deriv", &data->Deriv)); 308 CeedCallBackend(CeedGetKernel_Cuda(ceed, data->module, "DerivTranspose", &data->DerivTranspose)); 309 CeedCallBackend(CeedGetKernel_Cuda(ceed, data->module, "Weight", &data->Weight)); 310 CeedCallBackend(CeedFree(&basis_kernel_path)); 311 CeedCallBackend(CeedFree(&basis_kernel_source)); 312 313 CeedCallBackend(CeedBasisSetData(basis, data)); 314 315 // Register backend functions 316 CeedCallBackend(CeedSetBackendFunction(ceed, "Basis", basis, "Apply", CeedBasisApplyNonTensor_Cuda)); 317 CeedCallBackend(CeedSetBackendFunction(ceed, "Basis", basis, "Destroy", CeedBasisDestroyNonTensor_Cuda)); 318 return CEED_ERROR_SUCCESS; 319 } 320 321 //------------------------------------------------------------------------------ 322 // Create non-tensor H(div) 323 //------------------------------------------------------------------------------ 324 int CeedBasisCreateHdiv_Cuda(CeedElemTopology topo, CeedInt dim, CeedInt num_nodes, CeedInt num_qpts, const CeedScalar *interp, const CeedScalar *div, 325 const CeedScalar *q_ref, const CeedScalar *q_weight, CeedBasis basis) { 326 Ceed ceed; 327 char *basis_kernel_source; 328 const char *basis_kernel_path; 329 CeedInt num_comp, q_comp_interp, q_comp_div; 330 const CeedInt q_bytes = num_qpts * sizeof(CeedScalar); 331 CeedBasisNonTensor_Cuda *data; 332 333 CeedCallBackend(CeedBasisGetCeed(basis, &ceed)); 334 CeedCallBackend(CeedCalloc(1, &data)); 335 336 // Copy basis data to GPU 337 CeedCallBackend(CeedBasisGetNumQuadratureComponents(basis, CEED_EVAL_INTERP, &q_comp_interp)); 338 CeedCallBackend(CeedBasisGetNumQuadratureComponents(basis, CEED_EVAL_DIV, &q_comp_div)); 339 if (q_weight) { 340 CeedCallCuda(ceed, cudaMalloc((void **)&data->d_q_weight, q_bytes)); 341 CeedCallCuda(ceed, cudaMemcpy(data->d_q_weight, q_weight, q_bytes, cudaMemcpyHostToDevice)); 342 } 343 if (interp) { 344 const CeedInt interp_bytes = q_bytes * num_nodes * q_comp_interp; 345 346 CeedCallCuda(ceed, cudaMalloc((void **)&data->d_interp, interp_bytes)); 347 CeedCallCuda(ceed, cudaMemcpy(data->d_interp, interp, interp_bytes, cudaMemcpyHostToDevice)); 348 } 349 if (div) { 350 const CeedInt div_bytes = q_bytes * num_nodes * q_comp_div; 351 352 CeedCallCuda(ceed, cudaMalloc((void **)&data->d_div, div_bytes)); 353 CeedCallCuda(ceed, cudaMemcpy(data->d_div, div, div_bytes, cudaMemcpyHostToDevice)); 354 } 355 356 // Compile basis kernels 357 CeedCallBackend(CeedBasisGetNumComponents(basis, &num_comp)); 358 CeedCallBackend(CeedGetJitAbsolutePath(ceed, "ceed/jit-source/cuda/cuda-ref-basis-nontensor.h", &basis_kernel_path)); 359 CeedDebug256(ceed, CEED_DEBUG_COLOR_SUCCESS, "----- Loading Basis Kernel Source -----\n"); 360 CeedCallBackend(CeedLoadSourceToBuffer(ceed, basis_kernel_path, &basis_kernel_source)); 361 CeedDebug256(ceed, CEED_DEBUG_COLOR_SUCCESS, "----- Loading Basis Kernel Source Complete! -----\n"); 362 CeedCallBackend(CeedCompile_Cuda(ceed, basis_kernel_source, &data->module, 5, "BASIS_Q", num_qpts, "BASIS_P", num_nodes, "BASIS_Q_COMP_INTERP", 363 q_comp_interp, "BASIS_Q_COMP_DERIV", q_comp_div, "BASIS_NUM_COMP", num_comp)); 364 CeedCallBackend(CeedGetKernel_Cuda(ceed, data->module, "Interp", &data->Interp)); 365 CeedCallBackend(CeedGetKernel_Cuda(ceed, data->module, "InterpTranspose", &data->InterpTranspose)); 366 CeedCallBackend(CeedGetKernel_Cuda(ceed, data->module, "Deriv", &data->Deriv)); 367 CeedCallBackend(CeedGetKernel_Cuda(ceed, data->module, "DerivTranspose", &data->DerivTranspose)); 368 CeedCallBackend(CeedGetKernel_Cuda(ceed, data->module, "Weight", &data->Weight)); 369 CeedCallBackend(CeedFree(&basis_kernel_path)); 370 CeedCallBackend(CeedFree(&basis_kernel_source)); 371 372 CeedCallBackend(CeedBasisSetData(basis, data)); 373 374 // Register backend functions 375 CeedCallBackend(CeedSetBackendFunction(ceed, "Basis", basis, "Apply", CeedBasisApplyNonTensor_Cuda)); 376 CeedCallBackend(CeedSetBackendFunction(ceed, "Basis", basis, "Destroy", CeedBasisDestroyNonTensor_Cuda)); 377 return CEED_ERROR_SUCCESS; 378 } 379 380 //------------------------------------------------------------------------------ 381 // Create non-tensor H(curl) 382 //------------------------------------------------------------------------------ 383 int CeedBasisCreateHcurl_Cuda(CeedElemTopology topo, CeedInt dim, CeedInt num_nodes, CeedInt num_qpts, const CeedScalar *interp, 384 const CeedScalar *curl, const CeedScalar *q_ref, const CeedScalar *q_weight, CeedBasis basis) { 385 Ceed ceed; 386 char *basis_kernel_source; 387 const char *basis_kernel_path; 388 CeedInt num_comp, q_comp_interp, q_comp_curl; 389 const CeedInt q_bytes = num_qpts * sizeof(CeedScalar); 390 CeedBasisNonTensor_Cuda *data; 391 392 CeedCallBackend(CeedBasisGetCeed(basis, &ceed)); 393 CeedCallBackend(CeedCalloc(1, &data)); 394 395 // Copy basis data to GPU 396 CeedCallBackend(CeedBasisGetNumQuadratureComponents(basis, CEED_EVAL_INTERP, &q_comp_interp)); 397 CeedCallBackend(CeedBasisGetNumQuadratureComponents(basis, CEED_EVAL_CURL, &q_comp_curl)); 398 if (q_weight) { 399 CeedCallCuda(ceed, cudaMalloc((void **)&data->d_q_weight, q_bytes)); 400 CeedCallCuda(ceed, cudaMemcpy(data->d_q_weight, q_weight, q_bytes, cudaMemcpyHostToDevice)); 401 } 402 if (interp) { 403 const CeedInt interp_bytes = q_bytes * num_nodes * q_comp_interp; 404 405 CeedCallCuda(ceed, cudaMalloc((void **)&data->d_interp, interp_bytes)); 406 CeedCallCuda(ceed, cudaMemcpy(data->d_interp, interp, interp_bytes, cudaMemcpyHostToDevice)); 407 } 408 if (curl) { 409 const CeedInt curl_bytes = q_bytes * num_nodes * q_comp_curl; 410 411 CeedCallCuda(ceed, cudaMalloc((void **)&data->d_curl, curl_bytes)); 412 CeedCallCuda(ceed, cudaMemcpy(data->d_curl, curl, curl_bytes, cudaMemcpyHostToDevice)); 413 } 414 415 // Compile basis kernels 416 CeedCallBackend(CeedBasisGetNumComponents(basis, &num_comp)); 417 CeedCallBackend(CeedGetJitAbsolutePath(ceed, "ceed/jit-source/cuda/cuda-ref-basis-nontensor.h", &basis_kernel_path)); 418 CeedDebug256(ceed, CEED_DEBUG_COLOR_SUCCESS, "----- Loading Basis Kernel Source -----\n"); 419 CeedCallBackend(CeedLoadSourceToBuffer(ceed, basis_kernel_path, &basis_kernel_source)); 420 CeedDebug256(ceed, CEED_DEBUG_COLOR_SUCCESS, "----- Loading Basis Kernel Source Complete! -----\n"); 421 CeedCallBackend(CeedCompile_Cuda(ceed, basis_kernel_source, &data->module, 5, "BASIS_Q", num_qpts, "BASIS_P", num_nodes, "BASIS_Q_COMP_INTERP", 422 q_comp_interp, "BASIS_Q_COMP_DERIV", q_comp_curl, "BASIS_NUM_COMP", num_comp)); 423 CeedCallBackend(CeedGetKernel_Cuda(ceed, data->module, "Interp", &data->Interp)); 424 CeedCallBackend(CeedGetKernel_Cuda(ceed, data->module, "InterpTranspose", &data->InterpTranspose)); 425 CeedCallBackend(CeedGetKernel_Cuda(ceed, data->module, "Deriv", &data->Deriv)); 426 CeedCallBackend(CeedGetKernel_Cuda(ceed, data->module, "DerivTranspose", &data->DerivTranspose)); 427 CeedCallBackend(CeedGetKernel_Cuda(ceed, data->module, "Weight", &data->Weight)); 428 CeedCallBackend(CeedFree(&basis_kernel_path)); 429 CeedCallBackend(CeedFree(&basis_kernel_source)); 430 431 CeedCallBackend(CeedBasisSetData(basis, data)); 432 433 // Register backend functions 434 CeedCallBackend(CeedSetBackendFunction(ceed, "Basis", basis, "Apply", CeedBasisApplyNonTensor_Cuda)); 435 CeedCallBackend(CeedSetBackendFunction(ceed, "Basis", basis, "Destroy", CeedBasisDestroyNonTensor_Cuda)); 436 return CEED_ERROR_SUCCESS; 437 } 438 439 //------------------------------------------------------------------------------ 440