xref: /libCEED/backends/hip-ref/ceed-hip-ref-basis.c (revision 391f7d98b23119a3db61bce3e938b5dca7339952)
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 <string.h>
12 #include <hip/hip_runtime.h>
13 
14 #include "../hip/ceed-hip-common.h"
15 #include "../hip/ceed-hip-compile.h"
16 #include "ceed-hip-ref.h"
17 
18 //------------------------------------------------------------------------------
19 // Basis apply - tensor
20 //------------------------------------------------------------------------------
21 static int CeedBasisApplyCore_Hip(CeedBasis basis, bool apply_add, const CeedInt num_elem, CeedTransposeMode t_mode, CeedEvalMode eval_mode,
22                                   CeedVector u, CeedVector v) {
23   Ceed              ceed;
24   CeedInt           Q_1d, dim;
25   const CeedInt     is_transpose   = t_mode == CEED_TRANSPOSE;
26   const int         max_block_size = 64;
27   const CeedScalar *d_u;
28   CeedScalar       *d_v;
29   CeedBasis_Hip    *data;
30 
31   CeedCallBackend(CeedBasisGetCeed(basis, &ceed));
32   CeedCallBackend(CeedBasisGetData(basis, &data));
33 
34   // Get read/write access to u, v
35   if (u != CEED_VECTOR_NONE) CeedCallBackend(CeedVectorGetArrayRead(u, CEED_MEM_DEVICE, &d_u));
36   else CeedCheck(eval_mode == CEED_EVAL_WEIGHT, ceed, CEED_ERROR_BACKEND, "An input vector is required for this CeedEvalMode");
37   if (apply_add) CeedCallBackend(CeedVectorGetArray(v, CEED_MEM_DEVICE, &d_v));
38   else CeedCallBackend(CeedVectorGetArrayWrite(v, CEED_MEM_DEVICE, &d_v));
39 
40   // Clear v for transpose operation
41   if (is_transpose && !apply_add) {
42     CeedInt  num_comp, q_comp, num_nodes, num_qpts;
43     CeedSize length;
44 
45     CeedCallBackend(CeedBasisGetNumComponents(basis, &num_comp));
46     CeedCallBackend(CeedBasisGetNumQuadratureComponents(basis, eval_mode, &q_comp));
47     CeedCallBackend(CeedBasisGetNumNodes(basis, &num_nodes));
48     CeedCallBackend(CeedBasisGetNumQuadraturePoints(basis, &num_qpts));
49     length = (CeedSize)num_elem * (CeedSize)num_comp * (t_mode == CEED_TRANSPOSE ? (CeedSize)num_nodes : ((CeedSize)num_qpts * (CeedSize)q_comp));
50     CeedCallHip(ceed, hipMemset(d_v, 0, length * sizeof(CeedScalar)));
51   }
52   CeedCallBackend(CeedBasisGetNumQuadraturePoints1D(basis, &Q_1d));
53   CeedCallBackend(CeedBasisGetDimension(basis, &dim));
54 
55   // Basis action
56   switch (eval_mode) {
57     case CEED_EVAL_INTERP: {
58       void         *interp_args[] = {(void *)&num_elem, (void *)&is_transpose, &data->d_interp_1d, &d_u, &d_v};
59       const CeedInt block_size    = CeedIntMin(CeedIntPow(Q_1d, dim), max_block_size);
60 
61       CeedCallBackend(CeedRunKernel_Hip(ceed, data->Interp, num_elem, block_size, interp_args));
62     } break;
63     case CEED_EVAL_GRAD: {
64       void         *grad_args[] = {(void *)&num_elem, (void *)&is_transpose, &data->d_interp_1d, &data->d_grad_1d, &d_u, &d_v};
65       const CeedInt block_size  = max_block_size;
66 
67       CeedCallBackend(CeedRunKernel_Hip(ceed, data->Grad, num_elem, block_size, grad_args));
68     } break;
69     case CEED_EVAL_WEIGHT: {
70       CeedCheck(data->d_q_weight_1d, ceed, CEED_ERROR_BACKEND, "%s not supported; q_weights_1d not set", CeedEvalModes[eval_mode]);
71       void     *weight_args[] = {(void *)&num_elem, (void *)&data->d_q_weight_1d, &d_v};
72       const int block_size_x  = Q_1d;
73       const int block_size_y  = dim >= 2 ? Q_1d : 1;
74 
75       CeedCallBackend(CeedRunKernelDim_Hip(ceed, data->Weight, num_elem, block_size_x, block_size_y, 1, weight_args));
76     } break;
77     case CEED_EVAL_NONE: /* handled separately below */
78       break;
79     // LCOV_EXCL_START
80     case CEED_EVAL_DIV:
81     case CEED_EVAL_CURL:
82       return CeedError(ceed, CEED_ERROR_BACKEND, "%s not supported", CeedEvalModes[eval_mode]);
83       // LCOV_EXCL_STOP
84   }
85 
86   // Restore vectors, cover CEED_EVAL_NONE
87   CeedCallBackend(CeedVectorRestoreArray(v, &d_v));
88   if (eval_mode == CEED_EVAL_NONE) CeedCallBackend(CeedVectorSetArray(v, CEED_MEM_DEVICE, CEED_COPY_VALUES, (CeedScalar *)d_u));
89   if (eval_mode != CEED_EVAL_WEIGHT) CeedCallBackend(CeedVectorRestoreArrayRead(u, &d_u));
90   return CEED_ERROR_SUCCESS;
91 }
92 
93 static int CeedBasisApply_Hip(CeedBasis basis, const CeedInt num_elem, CeedTransposeMode t_mode, CeedEvalMode eval_mode, CeedVector u, CeedVector v) {
94   CeedCallBackend(CeedBasisApplyCore_Hip(basis, false, num_elem, t_mode, eval_mode, u, v));
95   return CEED_ERROR_SUCCESS;
96 }
97 
98 static int CeedBasisApplyAdd_Hip(CeedBasis basis, const CeedInt num_elem, CeedTransposeMode t_mode, CeedEvalMode eval_mode, CeedVector u,
99                                  CeedVector v) {
100   CeedCallBackend(CeedBasisApplyCore_Hip(basis, true, num_elem, t_mode, eval_mode, u, v));
101   return CEED_ERROR_SUCCESS;
102 }
103 
104 //------------------------------------------------------------------------------
105 // Basis apply - tensor AtPoints
106 //------------------------------------------------------------------------------
107 static int CeedBasisApplyAtPointsCore_Hip(CeedBasis basis, bool apply_add, const CeedInt num_elem, const CeedInt *num_points,
108                                           CeedTransposeMode t_mode, CeedEvalMode eval_mode, CeedVector x_ref, CeedVector u, CeedVector v) {
109   Ceed              ceed;
110   CeedInt           Q_1d, dim, max_num_points = num_points[0];
111   const CeedInt     is_transpose   = t_mode == CEED_TRANSPOSE;
112   const int         max_block_size = 32;
113   const CeedScalar *d_x, *d_u;
114   CeedScalar       *d_v;
115   CeedBasis_Hip    *data;
116 
117   CeedCallBackend(CeedBasisGetCeed(basis, &ceed));
118   CeedCallBackend(CeedBasisGetData(basis, &data));
119   CeedCallBackend(CeedBasisGetNumQuadraturePoints1D(basis, &Q_1d));
120   CeedCallBackend(CeedBasisGetDimension(basis, &dim));
121 
122   // Weight handled separately
123   if (eval_mode == CEED_EVAL_WEIGHT) {
124     CeedCallBackend(CeedVectorSetValue(v, 1.0));
125     return CEED_ERROR_SUCCESS;
126   }
127 
128   // Check padded to uniform number of points per elem
129   for (CeedInt i = 1; i < num_elem; i++) max_num_points = CeedIntMax(max_num_points, num_points[i]);
130   {
131     CeedInt  num_comp, q_comp;
132     CeedSize len, len_required;
133 
134     CeedCallBackend(CeedBasisGetNumComponents(basis, &num_comp));
135     CeedCallBackend(CeedBasisGetNumQuadratureComponents(basis, eval_mode, &q_comp));
136     CeedCallBackend(CeedVectorGetLength(is_transpose ? u : v, &len));
137     len_required = (CeedSize)num_comp * (CeedSize)q_comp * (CeedSize)num_elem * (CeedSize)max_num_points;
138     CeedCheck(len >= len_required, ceed, CEED_ERROR_BACKEND,
139               "Vector at points must be padded to the same number of points in each element for BasisApplyAtPoints on GPU backends."
140               " Found %" CeedSize_FMT ", Required %" CeedSize_FMT,
141               len, len_required);
142   }
143 
144   // Move num_points array to device
145   if (is_transpose) {
146     const CeedInt num_bytes = num_elem * sizeof(CeedInt);
147 
148     if (num_elem != data->num_elem_at_points) {
149       data->num_elem_at_points = num_elem;
150 
151       if (data->d_points_per_elem) CeedCallHip(ceed, hipFree(data->d_points_per_elem));
152       CeedCallHip(ceed, hipMalloc((void **)&data->d_points_per_elem, num_bytes));
153       CeedCallBackend(CeedFree(&data->h_points_per_elem));
154       CeedCallBackend(CeedCalloc(num_elem, &data->h_points_per_elem));
155     }
156     if (memcmp(data->h_points_per_elem, num_points, num_bytes)) {
157       memcpy(data->h_points_per_elem, num_points, num_bytes);
158       CeedCallHip(ceed, hipMemcpy(data->d_points_per_elem, num_points, num_bytes, hipMemcpyHostToDevice));
159     }
160   }
161 
162   // Build kernels if needed
163   if (data->num_points != max_num_points) {
164     CeedInt P_1d;
165 
166     CeedCallBackend(CeedBasisGetNumNodes1D(basis, &P_1d));
167     data->num_points = max_num_points;
168 
169     // -- Create interp matrix to Chebyshev coefficients
170     if (!data->d_chebyshev_interp_1d) {
171       CeedSize    interp_bytes;
172       CeedScalar *chebyshev_interp_1d;
173 
174       interp_bytes = P_1d * Q_1d * sizeof(CeedScalar);
175       CeedCallBackend(CeedCalloc(P_1d * Q_1d, &chebyshev_interp_1d));
176       CeedCallBackend(CeedBasisGetChebyshevInterp1D(basis, chebyshev_interp_1d));
177       CeedCallHip(ceed, hipMalloc((void **)&data->d_chebyshev_interp_1d, interp_bytes));
178       CeedCallHip(ceed, hipMemcpy(data->d_chebyshev_interp_1d, chebyshev_interp_1d, interp_bytes, hipMemcpyHostToDevice));
179       CeedCallBackend(CeedFree(&chebyshev_interp_1d));
180     }
181 
182     // -- Compile kernels
183     const char basis_kernel_source[] = "// AtPoints basis source\n#include <ceed/jit-source/hip/hip-ref-basis-tensor-at-points.h>\n";
184     CeedInt    num_comp;
185 
186     if (data->moduleAtPoints) CeedCallHip(ceed, hipModuleUnload(data->moduleAtPoints));
187     CeedCallBackend(CeedBasisGetNumComponents(basis, &num_comp));
188     CeedCallBackend(CeedCompile_Hip(ceed, basis_kernel_source, &data->moduleAtPoints, 9, "BASIS_Q_1D", Q_1d, "BASIS_P_1D", P_1d, "BASIS_BUF_LEN",
189                                     Q_1d * CeedIntPow(Q_1d > P_1d ? Q_1d : P_1d, dim - 1), "BASIS_DIM", dim, "BASIS_NUM_COMP", num_comp,
190                                     "BASIS_NUM_NODES", CeedIntPow(P_1d, dim), "BASIS_NUM_QPTS", CeedIntPow(Q_1d, dim), "BASIS_NUM_PTS",
191                                     max_num_points, "POINTS_BUFF_LEN", CeedIntPow(Q_1d, dim - 1)));
192     CeedCallBackend(CeedGetKernel_Hip(ceed, data->moduleAtPoints, "InterpAtPoints", &data->InterpAtPoints));
193     CeedCallBackend(CeedGetKernel_Hip(ceed, data->moduleAtPoints, "GradAtPoints", &data->GradAtPoints));
194   }
195 
196   // Get read/write access to u, v
197   CeedCallBackend(CeedVectorGetArrayRead(x_ref, CEED_MEM_DEVICE, &d_x));
198   if (u != CEED_VECTOR_NONE) CeedCallBackend(CeedVectorGetArrayRead(u, CEED_MEM_DEVICE, &d_u));
199   else CeedCheck(eval_mode == CEED_EVAL_WEIGHT, ceed, CEED_ERROR_BACKEND, "An input vector is required for this CeedEvalMode");
200   if (apply_add) CeedCallBackend(CeedVectorGetArray(v, CEED_MEM_DEVICE, &d_v));
201   else CeedCallBackend(CeedVectorGetArrayWrite(v, CEED_MEM_DEVICE, &d_v));
202 
203   // Clear v for transpose operation
204   if (is_transpose && !apply_add) {
205     CeedInt  num_comp, q_comp, num_nodes;
206     CeedSize length;
207 
208     CeedCallBackend(CeedBasisGetNumComponents(basis, &num_comp));
209     CeedCallBackend(CeedBasisGetNumQuadratureComponents(basis, eval_mode, &q_comp));
210     CeedCallBackend(CeedBasisGetNumNodes(basis, &num_nodes));
211     length =
212         (CeedSize)num_elem * (CeedSize)num_comp * (t_mode == CEED_TRANSPOSE ? (CeedSize)num_nodes : ((CeedSize)max_num_points * (CeedSize)q_comp));
213     CeedCallHip(ceed, hipMemset(d_v, 0, length * sizeof(CeedScalar)));
214   }
215 
216   // Basis action
217   switch (eval_mode) {
218     case CEED_EVAL_INTERP: {
219       void *interp_args[]      = {(void *)&num_elem, (void *)&is_transpose, &data->d_chebyshev_interp_1d, &data->d_points_per_elem, &d_x, &d_u, &d_v};
220       const CeedInt block_size = CeedIntMin(CeedIntPow(Q_1d, dim), max_block_size);
221 
222       CeedCallBackend(CeedRunKernel_Hip(ceed, data->InterpAtPoints, num_elem, block_size, interp_args));
223     } break;
224     case CEED_EVAL_GRAD: {
225       void *grad_args[]        = {(void *)&num_elem, (void *)&is_transpose, &data->d_chebyshev_interp_1d, &data->d_points_per_elem, &d_x, &d_u, &d_v};
226       const CeedInt block_size = CeedIntMin(CeedIntPow(Q_1d, dim), max_block_size);
227 
228       CeedCallBackend(CeedRunKernel_Hip(ceed, data->GradAtPoints, num_elem, block_size, grad_args));
229     } break;
230     case CEED_EVAL_WEIGHT:
231     case CEED_EVAL_NONE: /* handled separately below */
232       break;
233     // LCOV_EXCL_START
234     case CEED_EVAL_DIV:
235     case CEED_EVAL_CURL:
236       return CeedError(ceed, CEED_ERROR_BACKEND, "%s not supported", CeedEvalModes[eval_mode]);
237       // LCOV_EXCL_STOP
238   }
239 
240   // Restore vectors, cover CEED_EVAL_NONE
241   CeedCallBackend(CeedVectorRestoreArrayRead(x_ref, &d_x));
242   CeedCallBackend(CeedVectorRestoreArray(v, &d_v));
243   if (eval_mode == CEED_EVAL_NONE) CeedCallBackend(CeedVectorSetArray(v, CEED_MEM_DEVICE, CEED_COPY_VALUES, (CeedScalar *)d_u));
244   if (eval_mode != CEED_EVAL_WEIGHT) CeedCallBackend(CeedVectorRestoreArrayRead(u, &d_u));
245   return CEED_ERROR_SUCCESS;
246 }
247 
248 static int CeedBasisApplyAtPoints_Hip(CeedBasis basis, const CeedInt num_elem, const CeedInt *num_points, CeedTransposeMode t_mode,
249                                       CeedEvalMode eval_mode, CeedVector x_ref, CeedVector u, CeedVector v) {
250   CeedCallBackend(CeedBasisApplyAtPointsCore_Hip(basis, false, num_elem, num_points, t_mode, eval_mode, x_ref, u, v));
251   return CEED_ERROR_SUCCESS;
252 }
253 
254 static int CeedBasisApplyAddAtPoints_Hip(CeedBasis basis, const CeedInt num_elem, const CeedInt *num_points, CeedTransposeMode t_mode,
255                                          CeedEvalMode eval_mode, CeedVector x_ref, CeedVector u, CeedVector v) {
256   CeedCallBackend(CeedBasisApplyAtPointsCore_Hip(basis, true, num_elem, num_points, t_mode, eval_mode, x_ref, u, v));
257   return CEED_ERROR_SUCCESS;
258 }
259 
260 //------------------------------------------------------------------------------
261 // Basis apply - non-tensor
262 //------------------------------------------------------------------------------
263 static int CeedBasisApplyNonTensorCore_Hip(CeedBasis basis, bool apply_add, const CeedInt num_elem, CeedTransposeMode t_mode, CeedEvalMode eval_mode,
264                                            CeedVector u, CeedVector v) {
265   Ceed                    ceed;
266   CeedInt                 num_nodes, num_qpts;
267   const CeedInt           is_transpose    = t_mode == CEED_TRANSPOSE;
268   const int               elems_per_block = 1;
269   const int               grid            = CeedDivUpInt(num_elem, elems_per_block);
270   const CeedScalar       *d_u;
271   CeedScalar             *d_v;
272   CeedBasisNonTensor_Hip *data;
273 
274   CeedCallBackend(CeedBasisGetCeed(basis, &ceed));
275   CeedCallBackend(CeedBasisGetData(basis, &data));
276   CeedCallBackend(CeedBasisGetNumQuadraturePoints(basis, &num_qpts));
277   CeedCallBackend(CeedBasisGetNumNodes(basis, &num_nodes));
278 
279   // Get read/write access to u, v
280   if (u != CEED_VECTOR_NONE) CeedCallBackend(CeedVectorGetArrayRead(u, CEED_MEM_DEVICE, &d_u));
281   else CeedCheck(eval_mode == CEED_EVAL_WEIGHT, ceed, CEED_ERROR_BACKEND, "An input vector is required for this CeedEvalMode");
282   if (apply_add) CeedCallBackend(CeedVectorGetArray(v, CEED_MEM_DEVICE, &d_v));
283   else CeedCallBackend(CeedVectorGetArrayWrite(v, CEED_MEM_DEVICE, &d_v));
284 
285   // Clear v for transpose operation
286   if (is_transpose && !apply_add) {
287     CeedSize length;
288 
289     CeedCallBackend(CeedVectorGetLength(v, &length));
290     CeedCallHip(ceed, hipMemset(d_v, 0, length * sizeof(CeedScalar)));
291   }
292 
293   // Apply basis operation
294   switch (eval_mode) {
295     case CEED_EVAL_INTERP: {
296       void     *interp_args[] = {(void *)&num_elem, &data->d_interp, &d_u, &d_v};
297       const int block_size_x  = is_transpose ? num_nodes : num_qpts;
298 
299       if (is_transpose) {
300         CeedCallBackend(CeedRunKernelDim_Hip(ceed, data->InterpTranspose, grid, block_size_x, 1, elems_per_block, interp_args));
301       } else {
302         CeedCallBackend(CeedRunKernelDim_Hip(ceed, data->Interp, grid, block_size_x, 1, elems_per_block, interp_args));
303       }
304     } break;
305     case CEED_EVAL_GRAD: {
306       void     *grad_args[]  = {(void *)&num_elem, &data->d_grad, &d_u, &d_v};
307       const int block_size_x = is_transpose ? num_nodes : num_qpts;
308 
309       if (is_transpose) {
310         CeedCallBackend(CeedRunKernelDim_Hip(ceed, data->DerivTranspose, grid, block_size_x, 1, elems_per_block, grad_args));
311       } else {
312         CeedCallBackend(CeedRunKernelDim_Hip(ceed, data->Deriv, grid, block_size_x, 1, elems_per_block, grad_args));
313       }
314     } break;
315     case CEED_EVAL_DIV: {
316       void     *div_args[]   = {(void *)&num_elem, &data->d_div, &d_u, &d_v};
317       const int block_size_x = is_transpose ? num_nodes : num_qpts;
318 
319       if (is_transpose) {
320         CeedCallBackend(CeedRunKernelDim_Hip(ceed, data->DerivTranspose, grid, block_size_x, 1, elems_per_block, div_args));
321       } else {
322         CeedCallBackend(CeedRunKernelDim_Hip(ceed, data->Deriv, grid, block_size_x, 1, elems_per_block, div_args));
323       }
324     } break;
325     case CEED_EVAL_CURL: {
326       void     *curl_args[]  = {(void *)&num_elem, &data->d_curl, &d_u, &d_v};
327       const int block_size_x = is_transpose ? num_nodes : num_qpts;
328 
329       if (is_transpose) {
330         CeedCallBackend(CeedRunKernelDim_Hip(ceed, data->DerivTranspose, grid, block_size_x, 1, elems_per_block, curl_args));
331       } else {
332         CeedCallBackend(CeedRunKernelDim_Hip(ceed, data->Deriv, grid, block_size_x, 1, elems_per_block, curl_args));
333       }
334     } break;
335     case CEED_EVAL_WEIGHT: {
336       CeedCheck(data->d_q_weight, ceed, CEED_ERROR_BACKEND, "%s not supported; q_weights not set", CeedEvalModes[eval_mode]);
337       void *weight_args[] = {(void *)&num_elem, (void *)&data->d_q_weight, &d_v};
338 
339       CeedCallBackend(CeedRunKernelDim_Hip(ceed, data->Weight, grid, num_qpts, 1, elems_per_block, weight_args));
340     } break;
341     case CEED_EVAL_NONE: /* handled separately below */
342       break;
343   }
344 
345   // Restore vectors, cover CEED_EVAL_NONE
346   CeedCallBackend(CeedVectorRestoreArray(v, &d_v));
347   if (eval_mode == CEED_EVAL_NONE) CeedCallBackend(CeedVectorSetArray(v, CEED_MEM_DEVICE, CEED_COPY_VALUES, (CeedScalar *)d_u));
348   if (eval_mode != CEED_EVAL_WEIGHT) CeedCallBackend(CeedVectorRestoreArrayRead(u, &d_u));
349   return CEED_ERROR_SUCCESS;
350 }
351 
352 static int CeedBasisApplyNonTensor_Hip(CeedBasis basis, const CeedInt num_elem, CeedTransposeMode t_mode, CeedEvalMode eval_mode, CeedVector u,
353                                        CeedVector v) {
354   CeedCallBackend(CeedBasisApplyNonTensorCore_Hip(basis, false, num_elem, t_mode, eval_mode, u, v));
355   return CEED_ERROR_SUCCESS;
356 }
357 
358 static int CeedBasisApplyAddNonTensor_Hip(CeedBasis basis, const CeedInt num_elem, CeedTransposeMode t_mode, CeedEvalMode eval_mode, CeedVector u,
359                                           CeedVector v) {
360   CeedCallBackend(CeedBasisApplyNonTensorCore_Hip(basis, true, num_elem, t_mode, eval_mode, u, v));
361   return CEED_ERROR_SUCCESS;
362 }
363 
364 //------------------------------------------------------------------------------
365 // Destroy tensor basis
366 //------------------------------------------------------------------------------
367 static int CeedBasisDestroy_Hip(CeedBasis basis) {
368   Ceed           ceed;
369   CeedBasis_Hip *data;
370 
371   CeedCallBackend(CeedBasisGetCeed(basis, &ceed));
372   CeedCallBackend(CeedBasisGetData(basis, &data));
373   CeedCallHip(ceed, hipModuleUnload(data->module));
374   if (data->moduleAtPoints) CeedCallHip(ceed, hipModuleUnload(data->moduleAtPoints));
375   if (data->d_q_weight_1d) CeedCallHip(ceed, hipFree(data->d_q_weight_1d));
376   CeedCallBackend(CeedFree(&data->h_points_per_elem));
377   if (data->d_points_per_elem) CeedCallHip(ceed, hipFree(data->d_points_per_elem));
378   CeedCallHip(ceed, hipFree(data->d_interp_1d));
379   CeedCallHip(ceed, hipFree(data->d_grad_1d));
380   CeedCallHip(ceed, hipFree(data->d_chebyshev_interp_1d));
381   CeedCallBackend(CeedFree(&data));
382   return CEED_ERROR_SUCCESS;
383 }
384 
385 //------------------------------------------------------------------------------
386 // Destroy non-tensor basis
387 //------------------------------------------------------------------------------
388 static int CeedBasisDestroyNonTensor_Hip(CeedBasis basis) {
389   Ceed                    ceed;
390   CeedBasisNonTensor_Hip *data;
391 
392   CeedCallBackend(CeedBasisGetCeed(basis, &ceed));
393   CeedCallBackend(CeedBasisGetData(basis, &data));
394   CeedCallHip(ceed, hipModuleUnload(data->module));
395   if (data->d_q_weight) CeedCallHip(ceed, hipFree(data->d_q_weight));
396   CeedCallHip(ceed, hipFree(data->d_interp));
397   CeedCallHip(ceed, hipFree(data->d_grad));
398   CeedCallHip(ceed, hipFree(data->d_div));
399   CeedCallHip(ceed, hipFree(data->d_curl));
400   CeedCallBackend(CeedFree(&data));
401   return CEED_ERROR_SUCCESS;
402 }
403 
404 //------------------------------------------------------------------------------
405 // Create tensor
406 //------------------------------------------------------------------------------
407 int CeedBasisCreateTensorH1_Hip(CeedInt dim, CeedInt P_1d, CeedInt Q_1d, const CeedScalar *interp_1d, const CeedScalar *grad_1d,
408                                 const CeedScalar *q_ref_1d, const CeedScalar *q_weight_1d, CeedBasis basis) {
409   Ceed           ceed;
410   CeedInt        num_comp;
411   const CeedInt  q_bytes      = Q_1d * sizeof(CeedScalar);
412   const CeedInt  interp_bytes = q_bytes * P_1d;
413   CeedBasis_Hip *data;
414 
415   CeedCallBackend(CeedBasisGetCeed(basis, &ceed));
416   CeedCallBackend(CeedCalloc(1, &data));
417 
418   // Copy data to GPU
419   if (q_weight_1d) {
420     CeedCallHip(ceed, hipMalloc((void **)&data->d_q_weight_1d, q_bytes));
421     CeedCallHip(ceed, hipMemcpy(data->d_q_weight_1d, q_weight_1d, q_bytes, hipMemcpyHostToDevice));
422   }
423   CeedCallHip(ceed, hipMalloc((void **)&data->d_interp_1d, interp_bytes));
424   CeedCallHip(ceed, hipMemcpy(data->d_interp_1d, interp_1d, interp_bytes, hipMemcpyHostToDevice));
425   CeedCallHip(ceed, hipMalloc((void **)&data->d_grad_1d, interp_bytes));
426   CeedCallHip(ceed, hipMemcpy(data->d_grad_1d, grad_1d, interp_bytes, hipMemcpyHostToDevice));
427 
428   // Compile basis kernels
429   const char basis_kernel_source[] = "// Tensor basis source\n#include <ceed/jit-source/hip/hip-ref-basis-tensor.h>\n";
430 
431   CeedCallBackend(CeedBasisGetNumComponents(basis, &num_comp));
432   CeedCallBackend(CeedCompile_Hip(ceed, basis_kernel_source, &data->module, 7, "BASIS_Q_1D", Q_1d, "BASIS_P_1D", P_1d, "BASIS_BUF_LEN",
433                                   Q_1d * CeedIntPow(Q_1d > P_1d ? Q_1d : P_1d, dim - 1), "BASIS_DIM", dim, "BASIS_NUM_COMP", num_comp,
434                                   "BASIS_NUM_NODES", CeedIntPow(P_1d, dim), "BASIS_NUM_QPTS", CeedIntPow(Q_1d, dim)));
435   CeedCallBackend(CeedGetKernel_Hip(ceed, data->module, "Interp", &data->Interp));
436   CeedCallBackend(CeedGetKernel_Hip(ceed, data->module, "Grad", &data->Grad));
437   CeedCallBackend(CeedGetKernel_Hip(ceed, data->module, "Weight", &data->Weight));
438 
439   CeedCallBackend(CeedBasisSetData(basis, data));
440 
441   // Register backend functions
442   CeedCallBackend(CeedSetBackendFunction(ceed, "Basis", basis, "Apply", CeedBasisApply_Hip));
443   CeedCallBackend(CeedSetBackendFunction(ceed, "Basis", basis, "ApplyAdd", CeedBasisApplyAdd_Hip));
444   CeedCallBackend(CeedSetBackendFunction(ceed, "Basis", basis, "ApplyAtPoints", CeedBasisApplyAtPoints_Hip));
445   CeedCallBackend(CeedSetBackendFunction(ceed, "Basis", basis, "ApplyAddAtPoints", CeedBasisApplyAddAtPoints_Hip));
446   CeedCallBackend(CeedSetBackendFunction(ceed, "Basis", basis, "Destroy", CeedBasisDestroy_Hip));
447   return CEED_ERROR_SUCCESS;
448 }
449 
450 //------------------------------------------------------------------------------
451 // Create non-tensor H^1
452 //------------------------------------------------------------------------------
453 int CeedBasisCreateH1_Hip(CeedElemTopology topo, CeedInt dim, CeedInt num_nodes, CeedInt num_qpts, const CeedScalar *interp, const CeedScalar *grad,
454                           const CeedScalar *q_ref, const CeedScalar *q_weight, CeedBasis basis) {
455   Ceed                    ceed;
456   CeedInt                 num_comp, q_comp_interp, q_comp_grad;
457   const CeedInt           q_bytes = num_qpts * sizeof(CeedScalar);
458   CeedBasisNonTensor_Hip *data;
459 
460   CeedCallBackend(CeedBasisGetCeed(basis, &ceed));
461   CeedCallBackend(CeedCalloc(1, &data));
462 
463   // Copy basis data to GPU
464   CeedCallBackend(CeedBasisGetNumQuadratureComponents(basis, CEED_EVAL_INTERP, &q_comp_interp));
465   CeedCallBackend(CeedBasisGetNumQuadratureComponents(basis, CEED_EVAL_GRAD, &q_comp_grad));
466   if (q_weight) {
467     CeedCallHip(ceed, hipMalloc((void **)&data->d_q_weight, q_bytes));
468     CeedCallHip(ceed, hipMemcpy(data->d_q_weight, q_weight, q_bytes, hipMemcpyHostToDevice));
469   }
470   if (interp) {
471     const CeedInt interp_bytes = q_bytes * num_nodes * q_comp_interp;
472 
473     CeedCallHip(ceed, hipMalloc((void **)&data->d_interp, interp_bytes));
474     CeedCallHip(ceed, hipMemcpy(data->d_interp, interp, interp_bytes, hipMemcpyHostToDevice));
475   }
476   if (grad) {
477     const CeedInt grad_bytes = q_bytes * num_nodes * q_comp_grad;
478 
479     CeedCallHip(ceed, hipMalloc((void **)&data->d_grad, grad_bytes));
480     CeedCallHip(ceed, hipMemcpy(data->d_grad, grad, grad_bytes, hipMemcpyHostToDevice));
481   }
482 
483   // Compile basis kernels
484   const char basis_kernel_source[] = "// Nontensor basis source\n#include <ceed/jit-source/hip/hip-ref-basis-nontensor.h>\n";
485 
486   CeedCallBackend(CeedBasisGetNumComponents(basis, &num_comp));
487   CeedCallBackend(CeedCompile_Hip(ceed, basis_kernel_source, &data->module, 5, "BASIS_Q", num_qpts, "BASIS_P", num_nodes, "BASIS_Q_COMP_INTERP",
488                                   q_comp_interp, "BASIS_Q_COMP_DERIV", q_comp_grad, "BASIS_NUM_COMP", num_comp));
489   CeedCallBackend(CeedGetKernel_Hip(ceed, data->module, "Interp", &data->Interp));
490   CeedCallBackend(CeedGetKernel_Hip(ceed, data->module, "InterpTranspose", &data->InterpTranspose));
491   CeedCallBackend(CeedGetKernel_Hip(ceed, data->module, "Deriv", &data->Deriv));
492   CeedCallBackend(CeedGetKernel_Hip(ceed, data->module, "DerivTranspose", &data->DerivTranspose));
493   CeedCallBackend(CeedGetKernel_Hip(ceed, data->module, "Weight", &data->Weight));
494 
495   CeedCallBackend(CeedBasisSetData(basis, data));
496 
497   // Register backend functions
498   CeedCallBackend(CeedSetBackendFunction(ceed, "Basis", basis, "Apply", CeedBasisApplyNonTensor_Hip));
499   CeedCallBackend(CeedSetBackendFunction(ceed, "Basis", basis, "ApplyAdd", CeedBasisApplyAddNonTensor_Hip));
500   CeedCallBackend(CeedSetBackendFunction(ceed, "Basis", basis, "Destroy", CeedBasisDestroyNonTensor_Hip));
501   return CEED_ERROR_SUCCESS;
502 }
503 
504 //------------------------------------------------------------------------------
505 // Create non-tensor H(div)
506 //------------------------------------------------------------------------------
507 int CeedBasisCreateHdiv_Hip(CeedElemTopology topo, CeedInt dim, CeedInt num_nodes, CeedInt num_qpts, const CeedScalar *interp, const CeedScalar *div,
508                             const CeedScalar *q_ref, const CeedScalar *q_weight, CeedBasis basis) {
509   Ceed                    ceed;
510   CeedInt                 num_comp, q_comp_interp, q_comp_div;
511   const CeedInt           q_bytes = num_qpts * sizeof(CeedScalar);
512   CeedBasisNonTensor_Hip *data;
513 
514   CeedCallBackend(CeedBasisGetCeed(basis, &ceed));
515   CeedCallBackend(CeedCalloc(1, &data));
516 
517   // Copy basis data to GPU
518   CeedCallBackend(CeedBasisGetNumQuadratureComponents(basis, CEED_EVAL_INTERP, &q_comp_interp));
519   CeedCallBackend(CeedBasisGetNumQuadratureComponents(basis, CEED_EVAL_DIV, &q_comp_div));
520   if (q_weight) {
521     CeedCallHip(ceed, hipMalloc((void **)&data->d_q_weight, q_bytes));
522     CeedCallHip(ceed, hipMemcpy(data->d_q_weight, q_weight, q_bytes, hipMemcpyHostToDevice));
523   }
524   if (interp) {
525     const CeedInt interp_bytes = q_bytes * num_nodes * q_comp_interp;
526 
527     CeedCallHip(ceed, hipMalloc((void **)&data->d_interp, interp_bytes));
528     CeedCallHip(ceed, hipMemcpy(data->d_interp, interp, interp_bytes, hipMemcpyHostToDevice));
529   }
530   if (div) {
531     const CeedInt div_bytes = q_bytes * num_nodes * q_comp_div;
532 
533     CeedCallHip(ceed, hipMalloc((void **)&data->d_div, div_bytes));
534     CeedCallHip(ceed, hipMemcpy(data->d_div, div, div_bytes, hipMemcpyHostToDevice));
535   }
536 
537   // Compile basis kernels
538   const char basis_kernel_source[] = "// Nontensor basis source\n#include <ceed/jit-source/hip/hip-ref-basis-nontensor.h>\n";
539 
540   CeedCallBackend(CeedBasisGetNumComponents(basis, &num_comp));
541   CeedCallBackend(CeedCompile_Hip(ceed, basis_kernel_source, &data->module, 5, "BASIS_Q", num_qpts, "BASIS_P", num_nodes, "BASIS_Q_COMP_INTERP",
542                                   q_comp_interp, "BASIS_Q_COMP_DERIV", q_comp_div, "BASIS_NUM_COMP", num_comp));
543   CeedCallBackend(CeedGetKernel_Hip(ceed, data->module, "Interp", &data->Interp));
544   CeedCallBackend(CeedGetKernel_Hip(ceed, data->module, "InterpTranspose", &data->InterpTranspose));
545   CeedCallBackend(CeedGetKernel_Hip(ceed, data->module, "Deriv", &data->Deriv));
546   CeedCallBackend(CeedGetKernel_Hip(ceed, data->module, "DerivTranspose", &data->DerivTranspose));
547   CeedCallBackend(CeedGetKernel_Hip(ceed, data->module, "Weight", &data->Weight));
548 
549   CeedCallBackend(CeedBasisSetData(basis, data));
550 
551   // Register backend functions
552   CeedCallBackend(CeedSetBackendFunction(ceed, "Basis", basis, "Apply", CeedBasisApplyNonTensor_Hip));
553   CeedCallBackend(CeedSetBackendFunction(ceed, "Basis", basis, "ApplyAdd", CeedBasisApplyAddNonTensor_Hip));
554   CeedCallBackend(CeedSetBackendFunction(ceed, "Basis", basis, "Destroy", CeedBasisDestroyNonTensor_Hip));
555   return CEED_ERROR_SUCCESS;
556 }
557 
558 //------------------------------------------------------------------------------
559 // Create non-tensor H(curl)
560 //------------------------------------------------------------------------------
561 int CeedBasisCreateHcurl_Hip(CeedElemTopology topo, CeedInt dim, CeedInt num_nodes, CeedInt num_qpts, const CeedScalar *interp,
562                              const CeedScalar *curl, const CeedScalar *q_ref, const CeedScalar *q_weight, CeedBasis basis) {
563   Ceed                    ceed;
564   CeedInt                 num_comp, q_comp_interp, q_comp_curl;
565   const CeedInt           q_bytes = num_qpts * sizeof(CeedScalar);
566   CeedBasisNonTensor_Hip *data;
567 
568   CeedCallBackend(CeedBasisGetCeed(basis, &ceed));
569   CeedCallBackend(CeedCalloc(1, &data));
570 
571   // Copy basis data to GPU
572   CeedCallBackend(CeedBasisGetNumQuadratureComponents(basis, CEED_EVAL_INTERP, &q_comp_interp));
573   CeedCallBackend(CeedBasisGetNumQuadratureComponents(basis, CEED_EVAL_CURL, &q_comp_curl));
574   if (q_weight) {
575     CeedCallHip(ceed, hipMalloc((void **)&data->d_q_weight, q_bytes));
576     CeedCallHip(ceed, hipMemcpy(data->d_q_weight, q_weight, q_bytes, hipMemcpyHostToDevice));
577   }
578   if (interp) {
579     const CeedInt interp_bytes = q_bytes * num_nodes * q_comp_interp;
580 
581     CeedCallHip(ceed, hipMalloc((void **)&data->d_interp, interp_bytes));
582     CeedCallHip(ceed, hipMemcpy(data->d_interp, interp, interp_bytes, hipMemcpyHostToDevice));
583   }
584   if (curl) {
585     const CeedInt curl_bytes = q_bytes * num_nodes * q_comp_curl;
586 
587     CeedCallHip(ceed, hipMalloc((void **)&data->d_curl, curl_bytes));
588     CeedCallHip(ceed, hipMemcpy(data->d_curl, curl, curl_bytes, hipMemcpyHostToDevice));
589   }
590 
591   // Compile basis kernels
592   const char basis_kernel_source[] = "// Nontensor basis source\n#include <ceed/jit-source/hip/hip-ref-basis-nontensor.h>\n";
593 
594   CeedCallBackend(CeedBasisGetNumComponents(basis, &num_comp));
595   CeedCallBackend(CeedCompile_Hip(ceed, basis_kernel_source, &data->module, 5, "BASIS_Q", num_qpts, "BASIS_P", num_nodes, "BASIS_Q_COMP_INTERP",
596                                   q_comp_interp, "BASIS_Q_COMP_DERIV", q_comp_curl, "BASIS_NUM_COMP", num_comp));
597   CeedCallBackend(CeedGetKernel_Hip(ceed, data->module, "Interp", &data->Interp));
598   CeedCallBackend(CeedGetKernel_Hip(ceed, data->module, "InterpTranspose", &data->InterpTranspose));
599   CeedCallBackend(CeedGetKernel_Hip(ceed, data->module, "Deriv", &data->Deriv));
600   CeedCallBackend(CeedGetKernel_Hip(ceed, data->module, "DerivTranspose", &data->DerivTranspose));
601   CeedCallBackend(CeedGetKernel_Hip(ceed, data->module, "Weight", &data->Weight));
602 
603   CeedCallBackend(CeedBasisSetData(basis, data));
604 
605   // Register backend functions
606   CeedCallBackend(CeedSetBackendFunction(ceed, "Basis", basis, "Apply", CeedBasisApplyNonTensor_Hip));
607   CeedCallBackend(CeedSetBackendFunction(ceed, "Basis", basis, "ApplyAdd", CeedBasisApplyAddNonTensor_Hip));
608   CeedCallBackend(CeedSetBackendFunction(ceed, "Basis", basis, "Destroy", CeedBasisDestroyNonTensor_Hip));
609   return CEED_ERROR_SUCCESS;
610 }
611 
612 //------------------------------------------------------------------------------
613