Lines Matching +full:- +full:- +full:ceed
1 // Copyright (c) 2017-2026, Lawrence Livermore National Security, LLC and other CEED contributors.
2 // All Rights Reserved. See the top-level LICENSE and NOTICE files for details.
4 // SPDX-License-Identifier: BSD-2-Clause
6 // This file is part of CEED: http://github.com/ceed
8 #include <ceed.h>
9 #include <ceed/backend.h>
10 #include <ceed/jit-tools.h>
15 #include "../cuda/ceed-cuda-common.h"
16 #include "../cuda/ceed-cuda-compile.h"
17 #include "ceed-cuda-ref.h"
19 //------------------------------------------------------------------------------
20 // Basis apply - tensor
21 //------------------------------------------------------------------------------
24 Ceed ceed; in CeedBasisApplyCore_Cuda() local
32 CeedCallBackend(CeedBasisGetCeed(basis, &ceed)); in CeedBasisApplyCore_Cuda()
37 …else CeedCheck(eval_mode == CEED_EVAL_WEIGHT, ceed, CEED_ERROR_BACKEND, "An input vector is requir… in CeedBasisApplyCore_Cuda()
51 …void *interp_args[] = {(void *)&num_elem, (void *)&is_transpose, &data->d_interp_1d, &d_u,… in CeedBasisApplyCore_Cuda()
54 CeedCallBackend(CeedRunKernel_Cuda(ceed, data->Interp, num_elem, block_size, interp_args)); in CeedBasisApplyCore_Cuda()
57 … *grad_args[] = {(void *)&num_elem, (void *)&is_transpose, &data->d_interp_1d, &data->d_grad_1d,… in CeedBasisApplyCore_Cuda()
60 CeedCallBackend(CeedRunKernel_Cuda(ceed, data->Grad, num_elem, block_size, grad_args)); in CeedBasisApplyCore_Cuda()
63 …CeedCheck(data->d_q_weight_1d, ceed, CEED_ERROR_BACKEND, "%s not supported; q_weights_1d not set",… in CeedBasisApplyCore_Cuda()
64 void *weight_args[] = {(void *)&num_elem, (void *)&data->d_q_weight_1d, &d_v}; in CeedBasisApplyCore_Cuda()
68 …CeedCallBackend(CeedRunKernelDim_Cuda(ceed, data->Weight, num_elem, block_size_x, block_size_y, 1,… in CeedBasisApplyCore_Cuda()
75 return CeedError(ceed, CEED_ERROR_BACKEND, "%s not supported", CeedEvalModes[eval_mode]); in CeedBasisApplyCore_Cuda()
83 CeedCallBackend(CeedDestroy(&ceed)); in CeedBasisApplyCore_Cuda()
99 //------------------------------------------------------------------------------
100 // Basis apply - tensor AtPoints
101 //------------------------------------------------------------------------------
104 Ceed ceed; in CeedBasisApplyAtPointsCore_Cuda() local
122 CeedCallBackend(CeedBasisGetCeed(basis, &ceed)); in CeedBasisApplyAtPointsCore_Cuda()
134 CeedCheck(len >= len_required, ceed, CEED_ERROR_BACKEND, in CeedBasisApplyAtPointsCore_Cuda()
144 if (num_elem != data->num_elem_at_points) { in CeedBasisApplyAtPointsCore_Cuda()
145 data->num_elem_at_points = num_elem; in CeedBasisApplyAtPointsCore_Cuda()
147 if (data->d_points_per_elem) CeedCallCuda(ceed, cudaFree(data->d_points_per_elem)); in CeedBasisApplyAtPointsCore_Cuda()
148 CeedCallCuda(ceed, cudaMalloc((void **)&data->d_points_per_elem, num_bytes)); in CeedBasisApplyAtPointsCore_Cuda()
149 CeedCallBackend(CeedFree(&data->h_points_per_elem)); in CeedBasisApplyAtPointsCore_Cuda()
150 CeedCallBackend(CeedCalloc(num_elem, &data->h_points_per_elem)); in CeedBasisApplyAtPointsCore_Cuda()
152 if (memcmp(data->h_points_per_elem, num_points, num_bytes)) { in CeedBasisApplyAtPointsCore_Cuda()
153 memcpy(data->h_points_per_elem, num_points, num_bytes); in CeedBasisApplyAtPointsCore_Cuda()
154 …CeedCallCuda(ceed, cudaMemcpy(data->d_points_per_elem, num_points, num_bytes, cudaMemcpyHostToDevi… in CeedBasisApplyAtPointsCore_Cuda()
159 if (data->num_points != max_num_points) { in CeedBasisApplyAtPointsCore_Cuda()
163 data->num_points = max_num_points; in CeedBasisApplyAtPointsCore_Cuda()
165 // -- Create interp matrix to Chebyshev coefficients in CeedBasisApplyAtPointsCore_Cuda()
166 if (!data->d_chebyshev_interp_1d) { in CeedBasisApplyAtPointsCore_Cuda()
173 CeedCallCuda(ceed, cudaMalloc((void **)&data->d_chebyshev_interp_1d, interp_bytes)); in CeedBasisApplyAtPointsCore_Cuda()
174 …CeedCallCuda(ceed, cudaMemcpy(data->d_chebyshev_interp_1d, chebyshev_interp_1d, interp_bytes, cuda… in CeedBasisApplyAtPointsCore_Cuda()
178 // -- Compile kernels in CeedBasisApplyAtPointsCore_Cuda()
179 …rnel_source[] = "// AtPoints basis source\n#include <ceed/jit-source/cuda/cuda-ref-basis-tensor-at… in CeedBasisApplyAtPointsCore_Cuda()
182 if (data->moduleAtPoints) CeedCallCuda(ceed, cuModuleUnload(data->moduleAtPoints)); in CeedBasisApplyAtPointsCore_Cuda()
184 …CeedCallBackend(CeedCompile_Cuda(ceed, basis_kernel_source, &data->moduleAtPoints, 9, "BASIS_Q_1D"… in CeedBasisApplyAtPointsCore_Cuda()
185 …Q_1d * CeedIntPow(Q_1d > P_1d ? Q_1d : P_1d, dim - 1), "BASIS_DIM", dim, "BASIS_NUM_COMP", num_com… in CeedBasisApplyAtPointsCore_Cuda()
187 max_num_points, "POINTS_BUFF_LEN", CeedIntPow(Q_1d, dim - 1))); in CeedBasisApplyAtPointsCore_Cuda()
188 …CeedCallBackend(CeedGetKernel_Cuda(ceed, data->moduleAtPoints, "InterpAtPoints", &data->InterpAtPo… in CeedBasisApplyAtPointsCore_Cuda()
189 …CeedCallBackend(CeedGetKernel_Cuda(ceed, data->moduleAtPoints, "InterpTransposeAtPoints", &data->I… in CeedBasisApplyAtPointsCore_Cuda()
190 …CeedCallBackend(CeedGetKernel_Cuda(ceed, data->moduleAtPoints, "GradAtPoints", &data->GradAtPoints… in CeedBasisApplyAtPointsCore_Cuda()
191 …CeedCallBackend(CeedGetKernel_Cuda(ceed, data->moduleAtPoints, "GradTransposeAtPoints", &data->Gra… in CeedBasisApplyAtPointsCore_Cuda()
197 …else CeedCheck(eval_mode == CEED_EVAL_WEIGHT, ceed, CEED_ERROR_BACKEND, "An input vector is requir… in CeedBasisApplyAtPointsCore_Cuda()
209 …void *interp_args[] = {(void *)&num_elem, &data->d_chebyshev_interp_1d, &data->d_points_pe… in CeedBasisApplyAtPointsCore_Cuda()
212 …CeedCallBackend(CeedRunKernel_Cuda(ceed, is_transpose ? data->InterpTransposeAtPoints : data->Inte… in CeedBasisApplyAtPointsCore_Cuda()
216 …void *grad_args[] = {(void *)&num_elem, &data->d_chebyshev_interp_1d, &data->d_points_per_… in CeedBasisApplyAtPointsCore_Cuda()
219 …CeedCallBackend(CeedRunKernel_Cuda(ceed, is_transpose ? data->GradTransposeAtPoints : data->GradAt… in CeedBasisApplyAtPointsCore_Cuda()
227 return CeedError(ceed, CEED_ERROR_BACKEND, "%s not supported", CeedEvalModes[eval_mode]); in CeedBasisApplyAtPointsCore_Cuda()
236 CeedCallBackend(CeedDestroy(&ceed)); in CeedBasisApplyAtPointsCore_Cuda()
252 //------------------------------------------------------------------------------
253 // Basis apply - non-tensor
254 //------------------------------------------------------------------------------
257 Ceed ceed; in CeedBasisApplyNonTensorCore_Cuda() local
266 CeedCallBackend(CeedBasisGetCeed(basis, &ceed)); in CeedBasisApplyNonTensorCore_Cuda()
273 …else CeedCheck(eval_mode == CEED_EVAL_WEIGHT, ceed, CEED_ERROR_BACKEND, "An input vector is requir… in CeedBasisApplyNonTensorCore_Cuda()
285 void *interp_args[] = {(void *)&num_elem, &data->d_interp, &d_u, &d_v}; in CeedBasisApplyNonTensorCore_Cuda()
289 …CeedCallBackend(CeedRunKernelDim_Cuda(ceed, data->InterpTranspose, grid, block_size_x, 1, elems_pe… in CeedBasisApplyNonTensorCore_Cuda()
291 …CeedCallBackend(CeedRunKernelDim_Cuda(ceed, data->Interp, grid, block_size_x, 1, elems_per_block, … in CeedBasisApplyNonTensorCore_Cuda()
295 void *grad_args[] = {(void *)&num_elem, &data->d_grad, &d_u, &d_v}; in CeedBasisApplyNonTensorCore_Cuda()
299 …CeedCallBackend(CeedRunKernelDim_Cuda(ceed, data->DerivTranspose, grid, block_size_x, 1, elems_per… in CeedBasisApplyNonTensorCore_Cuda()
301 …CeedCallBackend(CeedRunKernelDim_Cuda(ceed, data->Deriv, grid, block_size_x, 1, elems_per_block, g… in CeedBasisApplyNonTensorCore_Cuda()
305 void *div_args[] = {(void *)&num_elem, &data->d_div, &d_u, &d_v}; in CeedBasisApplyNonTensorCore_Cuda()
309 …CeedCallBackend(CeedRunKernelDim_Cuda(ceed, data->DerivTranspose, grid, block_size_x, 1, elems_per… in CeedBasisApplyNonTensorCore_Cuda()
311 …CeedCallBackend(CeedRunKernelDim_Cuda(ceed, data->Deriv, grid, block_size_x, 1, elems_per_block, d… in CeedBasisApplyNonTensorCore_Cuda()
315 void *curl_args[] = {(void *)&num_elem, &data->d_curl, &d_u, &d_v}; in CeedBasisApplyNonTensorCore_Cuda()
319 …CeedCallBackend(CeedRunKernelDim_Cuda(ceed, data->DerivTranspose, grid, block_size_x, 1, elems_per… in CeedBasisApplyNonTensorCore_Cuda()
321 …CeedCallBackend(CeedRunKernelDim_Cuda(ceed, data->Deriv, grid, block_size_x, 1, elems_per_block, c… in CeedBasisApplyNonTensorCore_Cuda()
325 …CeedCheck(data->d_q_weight, ceed, CEED_ERROR_BACKEND, "%s not supported; q_weights not set", CeedE… in CeedBasisApplyNonTensorCore_Cuda()
326 void *weight_args[] = {(void *)&num_elem, (void *)&data->d_q_weight, &d_v}; in CeedBasisApplyNonTensorCore_Cuda()
328 …CeedCallBackend(CeedRunKernelDim_Cuda(ceed, data->Weight, grid, num_qpts, 1, elems_per_block, weig… in CeedBasisApplyNonTensorCore_Cuda()
338 CeedCallBackend(CeedDestroy(&ceed)); in CeedBasisApplyNonTensorCore_Cuda()
354 //------------------------------------------------------------------------------
356 //------------------------------------------------------------------------------
358 Ceed ceed; in CeedBasisDestroy_Cuda() local
361 CeedCallBackend(CeedBasisGetCeed(basis, &ceed)); in CeedBasisDestroy_Cuda()
363 CeedCallCuda(ceed, cuModuleUnload(data->module)); in CeedBasisDestroy_Cuda()
364 if (data->moduleAtPoints) CeedCallCuda(ceed, cuModuleUnload(data->moduleAtPoints)); in CeedBasisDestroy_Cuda()
365 if (data->d_q_weight_1d) CeedCallCuda(ceed, cudaFree(data->d_q_weight_1d)); in CeedBasisDestroy_Cuda()
366 CeedCallBackend(CeedFree(&data->h_points_per_elem)); in CeedBasisDestroy_Cuda()
367 if (data->d_points_per_elem) CeedCallCuda(ceed, cudaFree(data->d_points_per_elem)); in CeedBasisDestroy_Cuda()
368 CeedCallCuda(ceed, cudaFree(data->d_interp_1d)); in CeedBasisDestroy_Cuda()
369 CeedCallCuda(ceed, cudaFree(data->d_grad_1d)); in CeedBasisDestroy_Cuda()
370 CeedCallCuda(ceed, cudaFree(data->d_chebyshev_interp_1d)); in CeedBasisDestroy_Cuda()
372 CeedCallBackend(CeedDestroy(&ceed)); in CeedBasisDestroy_Cuda()
376 //------------------------------------------------------------------------------
377 // Destroy non-tensor basis
378 //------------------------------------------------------------------------------
380 Ceed ceed; in CeedBasisDestroyNonTensor_Cuda() local
383 CeedCallBackend(CeedBasisGetCeed(basis, &ceed)); in CeedBasisDestroyNonTensor_Cuda()
385 CeedCallCuda(ceed, cuModuleUnload(data->module)); in CeedBasisDestroyNonTensor_Cuda()
386 if (data->d_q_weight) CeedCallCuda(ceed, cudaFree(data->d_q_weight)); in CeedBasisDestroyNonTensor_Cuda()
387 CeedCallCuda(ceed, cudaFree(data->d_interp)); in CeedBasisDestroyNonTensor_Cuda()
388 CeedCallCuda(ceed, cudaFree(data->d_grad)); in CeedBasisDestroyNonTensor_Cuda()
389 CeedCallCuda(ceed, cudaFree(data->d_div)); in CeedBasisDestroyNonTensor_Cuda()
390 CeedCallCuda(ceed, cudaFree(data->d_curl)); in CeedBasisDestroyNonTensor_Cuda()
392 CeedCallBackend(CeedDestroy(&ceed)); in CeedBasisDestroyNonTensor_Cuda()
396 //------------------------------------------------------------------------------
398 //------------------------------------------------------------------------------
401 Ceed ceed; in CeedBasisCreateTensorH1_Cuda() local
407 CeedCallBackend(CeedBasisGetCeed(basis, &ceed)); in CeedBasisCreateTensorH1_Cuda()
412 CeedCallCuda(ceed, cudaMalloc((void **)&data->d_q_weight_1d, q_bytes)); in CeedBasisCreateTensorH1_Cuda()
413 … CeedCallCuda(ceed, cudaMemcpy(data->d_q_weight_1d, q_weight_1d, q_bytes, cudaMemcpyHostToDevice)); in CeedBasisCreateTensorH1_Cuda()
415 CeedCallCuda(ceed, cudaMalloc((void **)&data->d_interp_1d, interp_bytes)); in CeedBasisCreateTensorH1_Cuda()
416 …CeedCallCuda(ceed, cudaMemcpy(data->d_interp_1d, interp_1d, interp_bytes, cudaMemcpyHostToDevice)); in CeedBasisCreateTensorH1_Cuda()
417 CeedCallCuda(ceed, cudaMalloc((void **)&data->d_grad_1d, interp_bytes)); in CeedBasisCreateTensorH1_Cuda()
418 CeedCallCuda(ceed, cudaMemcpy(data->d_grad_1d, grad_1d, interp_bytes, cudaMemcpyHostToDevice)); in CeedBasisCreateTensorH1_Cuda()
421 …asis_kernel_source[] = "// Tensor basis source\n#include <ceed/jit-source/cuda/cuda-ref-basis-tens… in CeedBasisCreateTensorH1_Cuda()
424 …CeedCallBackend(CeedCompile_Cuda(ceed, basis_kernel_source, &data->module, 7, "BASIS_Q_1D", Q_1d, … in CeedBasisCreateTensorH1_Cuda()
425 …Q_1d * CeedIntPow(Q_1d > P_1d ? Q_1d : P_1d, dim - 1), "BASIS_DIM", dim, "BASIS_NUM_COMP", num_com… in CeedBasisCreateTensorH1_Cuda()
427 CeedCallBackend(CeedGetKernel_Cuda(ceed, data->module, "Interp", &data->Interp)); in CeedBasisCreateTensorH1_Cuda()
428 CeedCallBackend(CeedGetKernel_Cuda(ceed, data->module, "Grad", &data->Grad)); in CeedBasisCreateTensorH1_Cuda()
429 CeedCallBackend(CeedGetKernel_Cuda(ceed, data->module, "Weight", &data->Weight)); in CeedBasisCreateTensorH1_Cuda()
434 CeedCallBackend(CeedSetBackendFunction(ceed, "Basis", basis, "Apply", CeedBasisApply_Cuda)); in CeedBasisCreateTensorH1_Cuda()
435 CeedCallBackend(CeedSetBackendFunction(ceed, "Basis", basis, "ApplyAdd", CeedBasisApplyAdd_Cuda)); in CeedBasisCreateTensorH1_Cuda()
436 …CeedCallBackend(CeedSetBackendFunction(ceed, "Basis", basis, "ApplyAtPoints", CeedBasisApplyAtPoin… in CeedBasisCreateTensorH1_Cuda()
437 …CeedCallBackend(CeedSetBackendFunction(ceed, "Basis", basis, "ApplyAddAtPoints", CeedBasisApplyAdd… in CeedBasisCreateTensorH1_Cuda()
438 CeedCallBackend(CeedSetBackendFunction(ceed, "Basis", basis, "Destroy", CeedBasisDestroy_Cuda)); in CeedBasisCreateTensorH1_Cuda()
439 CeedCallBackend(CeedDestroy(&ceed)); in CeedBasisCreateTensorH1_Cuda()
443 //------------------------------------------------------------------------------
444 // Create non-tensor H^1
445 //------------------------------------------------------------------------------
448 Ceed ceed; in CeedBasisCreateH1_Cuda() local
453 CeedCallBackend(CeedBasisGetCeed(basis, &ceed)); in CeedBasisCreateH1_Cuda()
460 CeedCallCuda(ceed, cudaMalloc((void **)&data->d_q_weight, q_bytes)); in CeedBasisCreateH1_Cuda()
461 CeedCallCuda(ceed, cudaMemcpy(data->d_q_weight, q_weight, q_bytes, cudaMemcpyHostToDevice)); in CeedBasisCreateH1_Cuda()
466 CeedCallCuda(ceed, cudaMalloc((void **)&data->d_interp, interp_bytes)); in CeedBasisCreateH1_Cuda()
467 CeedCallCuda(ceed, cudaMemcpy(data->d_interp, interp, interp_bytes, cudaMemcpyHostToDevice)); in CeedBasisCreateH1_Cuda()
472 CeedCallCuda(ceed, cudaMalloc((void **)&data->d_grad, grad_bytes)); in CeedBasisCreateH1_Cuda()
473 CeedCallCuda(ceed, cudaMemcpy(data->d_grad, grad, grad_bytes, cudaMemcpyHostToDevice)); in CeedBasisCreateH1_Cuda()
477 …s_kernel_source[] = "// Nontensor basis source\n#include <ceed/jit-source/cuda/cuda-ref-basis-nont… in CeedBasisCreateH1_Cuda()
480 …CeedCallBackend(CeedCompile_Cuda(ceed, basis_kernel_source, &data->module, 5, "BASIS_Q", num_qpts,… in CeedBasisCreateH1_Cuda()
482 CeedCallBackend(CeedGetKernel_Cuda(ceed, data->module, "Interp", &data->Interp)); in CeedBasisCreateH1_Cuda()
483 …CeedCallBackend(CeedGetKernel_Cuda(ceed, data->module, "InterpTranspose", &data->InterpTranspose)); in CeedBasisCreateH1_Cuda()
484 CeedCallBackend(CeedGetKernel_Cuda(ceed, data->module, "Deriv", &data->Deriv)); in CeedBasisCreateH1_Cuda()
485 CeedCallBackend(CeedGetKernel_Cuda(ceed, data->module, "DerivTranspose", &data->DerivTranspose)); in CeedBasisCreateH1_Cuda()
486 CeedCallBackend(CeedGetKernel_Cuda(ceed, data->module, "Weight", &data->Weight)); in CeedBasisCreateH1_Cuda()
491 …CeedCallBackend(CeedSetBackendFunction(ceed, "Basis", basis, "Apply", CeedBasisApplyNonTensor_Cuda… in CeedBasisCreateH1_Cuda()
492 …CeedCallBackend(CeedSetBackendFunction(ceed, "Basis", basis, "ApplyAdd", CeedBasisApplyAddNonTenso… in CeedBasisCreateH1_Cuda()
493 …CeedCallBackend(CeedSetBackendFunction(ceed, "Basis", basis, "Destroy", CeedBasisDestroyNonTensor_… in CeedBasisCreateH1_Cuda()
494 CeedCallBackend(CeedDestroy(&ceed)); in CeedBasisCreateH1_Cuda()
498 //------------------------------------------------------------------------------
499 // Create non-tensor H(div)
500 //------------------------------------------------------------------------------
503 Ceed ceed; in CeedBasisCreateHdiv_Cuda() local
508 CeedCallBackend(CeedBasisGetCeed(basis, &ceed)); in CeedBasisCreateHdiv_Cuda()
515 CeedCallCuda(ceed, cudaMalloc((void **)&data->d_q_weight, q_bytes)); in CeedBasisCreateHdiv_Cuda()
516 CeedCallCuda(ceed, cudaMemcpy(data->d_q_weight, q_weight, q_bytes, cudaMemcpyHostToDevice)); in CeedBasisCreateHdiv_Cuda()
521 CeedCallCuda(ceed, cudaMalloc((void **)&data->d_interp, interp_bytes)); in CeedBasisCreateHdiv_Cuda()
522 CeedCallCuda(ceed, cudaMemcpy(data->d_interp, interp, interp_bytes, cudaMemcpyHostToDevice)); in CeedBasisCreateHdiv_Cuda()
527 CeedCallCuda(ceed, cudaMalloc((void **)&data->d_div, div_bytes)); in CeedBasisCreateHdiv_Cuda()
528 CeedCallCuda(ceed, cudaMemcpy(data->d_div, div, div_bytes, cudaMemcpyHostToDevice)); in CeedBasisCreateHdiv_Cuda()
532 …s_kernel_source[] = "// Nontensor basis source\n#include <ceed/jit-source/cuda/cuda-ref-basis-nont… in CeedBasisCreateHdiv_Cuda()
535 …CeedCallBackend(CeedCompile_Cuda(ceed, basis_kernel_source, &data->module, 5, "BASIS_Q", num_qpts,… in CeedBasisCreateHdiv_Cuda()
537 CeedCallBackend(CeedGetKernel_Cuda(ceed, data->module, "Interp", &data->Interp)); in CeedBasisCreateHdiv_Cuda()
538 …CeedCallBackend(CeedGetKernel_Cuda(ceed, data->module, "InterpTranspose", &data->InterpTranspose)); in CeedBasisCreateHdiv_Cuda()
539 CeedCallBackend(CeedGetKernel_Cuda(ceed, data->module, "Deriv", &data->Deriv)); in CeedBasisCreateHdiv_Cuda()
540 CeedCallBackend(CeedGetKernel_Cuda(ceed, data->module, "DerivTranspose", &data->DerivTranspose)); in CeedBasisCreateHdiv_Cuda()
541 CeedCallBackend(CeedGetKernel_Cuda(ceed, data->module, "Weight", &data->Weight)); in CeedBasisCreateHdiv_Cuda()
546 …CeedCallBackend(CeedSetBackendFunction(ceed, "Basis", basis, "Apply", CeedBasisApplyNonTensor_Cuda… in CeedBasisCreateHdiv_Cuda()
547 …CeedCallBackend(CeedSetBackendFunction(ceed, "Basis", basis, "ApplyAdd", CeedBasisApplyAddNonTenso… in CeedBasisCreateHdiv_Cuda()
548 …CeedCallBackend(CeedSetBackendFunction(ceed, "Basis", basis, "Destroy", CeedBasisDestroyNonTensor_… in CeedBasisCreateHdiv_Cuda()
549 CeedCallBackend(CeedDestroy(&ceed)); in CeedBasisCreateHdiv_Cuda()
553 //------------------------------------------------------------------------------
554 // Create non-tensor H(curl)
555 //------------------------------------------------------------------------------
558 Ceed ceed; in CeedBasisCreateHcurl_Cuda() local
563 CeedCallBackend(CeedBasisGetCeed(basis, &ceed)); in CeedBasisCreateHcurl_Cuda()
570 CeedCallCuda(ceed, cudaMalloc((void **)&data->d_q_weight, q_bytes)); in CeedBasisCreateHcurl_Cuda()
571 CeedCallCuda(ceed, cudaMemcpy(data->d_q_weight, q_weight, q_bytes, cudaMemcpyHostToDevice)); in CeedBasisCreateHcurl_Cuda()
576 CeedCallCuda(ceed, cudaMalloc((void **)&data->d_interp, interp_bytes)); in CeedBasisCreateHcurl_Cuda()
577 CeedCallCuda(ceed, cudaMemcpy(data->d_interp, interp, interp_bytes, cudaMemcpyHostToDevice)); in CeedBasisCreateHcurl_Cuda()
582 CeedCallCuda(ceed, cudaMalloc((void **)&data->d_curl, curl_bytes)); in CeedBasisCreateHcurl_Cuda()
583 CeedCallCuda(ceed, cudaMemcpy(data->d_curl, curl, curl_bytes, cudaMemcpyHostToDevice)); in CeedBasisCreateHcurl_Cuda()
587 …s_kernel_source[] = "// Nontensor basis source\n#include <ceed/jit-source/cuda/cuda-ref-basis-nont… in CeedBasisCreateHcurl_Cuda()
590 …CeedCallBackend(CeedCompile_Cuda(ceed, basis_kernel_source, &data->module, 5, "BASIS_Q", num_qpts,… in CeedBasisCreateHcurl_Cuda()
592 CeedCallBackend(CeedGetKernel_Cuda(ceed, data->module, "Interp", &data->Interp)); in CeedBasisCreateHcurl_Cuda()
593 …CeedCallBackend(CeedGetKernel_Cuda(ceed, data->module, "InterpTranspose", &data->InterpTranspose)); in CeedBasisCreateHcurl_Cuda()
594 CeedCallBackend(CeedGetKernel_Cuda(ceed, data->module, "Deriv", &data->Deriv)); in CeedBasisCreateHcurl_Cuda()
595 CeedCallBackend(CeedGetKernel_Cuda(ceed, data->module, "DerivTranspose", &data->DerivTranspose)); in CeedBasisCreateHcurl_Cuda()
596 CeedCallBackend(CeedGetKernel_Cuda(ceed, data->module, "Weight", &data->Weight)); in CeedBasisCreateHcurl_Cuda()
601 …CeedCallBackend(CeedSetBackendFunction(ceed, "Basis", basis, "Apply", CeedBasisApplyNonTensor_Cuda… in CeedBasisCreateHcurl_Cuda()
602 …CeedCallBackend(CeedSetBackendFunction(ceed, "Basis", basis, "ApplyAdd", CeedBasisApplyAddNonTenso… in CeedBasisCreateHcurl_Cuda()
603 …CeedCallBackend(CeedSetBackendFunction(ceed, "Basis", basis, "Destroy", CeedBasisDestroyNonTensor_… in CeedBasisCreateHcurl_Cuda()
604 CeedCallBackend(CeedDestroy(&ceed)); in CeedBasisCreateHcurl_Cuda()
608 //------------------------------------------------------------------------------