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 <stdbool.h> 12 #include <stddef.h> 13 #include <hip/hip_runtime.h> 14 15 #include "../hip/ceed-hip-common.h" 16 #include "../hip/ceed-hip-compile.h" 17 #include "ceed-hip-shared.h" 18 19 //------------------------------------------------------------------------------ 20 // Compute a block size based on required minimum threads 21 //------------------------------------------------------------------------------ 22 static CeedInt ComputeBlockSizeFromRequirement(const CeedInt required) { 23 CeedInt maxSize = 1024; // Max total threads per block 24 CeedInt currentSize = 64; // Start with one group 25 26 while (currentSize < maxSize) { 27 if (currentSize > required) break; 28 else currentSize = currentSize * 2; 29 } 30 return currentSize; 31 } 32 33 //------------------------------------------------------------------------------ 34 // Compute required thread block sizes for basis kernels given P, Q, dim, and 35 // num_comp (num_comp not currently used, but may be again in other basis 36 // parallelization options) 37 //------------------------------------------------------------------------------ 38 static int ComputeBasisThreadBlockSizes(const CeedInt dim, const CeedInt P_1d, const CeedInt Q_1d, const CeedInt num_comp, CeedInt *block_sizes) { 39 // Note that this will use the same block sizes for all dimensions when compiling, 40 // but as each basis object is defined for a particular dimension, we will never 41 // call any kernels except the ones for the dimension for which we have computed the 42 // block sizes. 43 const CeedInt thread_1d = CeedIntMax(P_1d, Q_1d); 44 45 switch (dim) { 46 case 1: { 47 // Interp kernels: 48 block_sizes[0] = 256; 49 50 // Grad kernels: 51 block_sizes[1] = 256; 52 53 // Weight kernels: 54 block_sizes[2] = 256; 55 } break; 56 case 2: { 57 // Interp kernels: 58 CeedInt required = thread_1d * thread_1d; 59 60 block_sizes[0] = CeedIntMax(256, ComputeBlockSizeFromRequirement(required)); 61 62 // Grad kernels: currently use same required minimum threads 63 block_sizes[1] = CeedIntMax(256, ComputeBlockSizeFromRequirement(required)); 64 65 // Weight kernels: 66 required = CeedIntMax(64, Q_1d * Q_1d); 67 block_sizes[2] = CeedIntMax(256, ComputeBlockSizeFromRequirement(required)); 68 69 } break; 70 case 3: { 71 // Interp kernels: 72 CeedInt required = thread_1d * thread_1d; 73 74 block_sizes[0] = CeedIntMax(256, ComputeBlockSizeFromRequirement(required)); 75 76 // Grad kernels: currently use same required minimum threads 77 block_sizes[1] = CeedIntMax(256, ComputeBlockSizeFromRequirement(required)); 78 79 // Weight kernels: 80 required = Q_1d * Q_1d * Q_1d; 81 block_sizes[2] = CeedIntMax(256, ComputeBlockSizeFromRequirement(required)); 82 } 83 } 84 return CEED_ERROR_SUCCESS; 85 } 86 87 //------------------------------------------------------------------------------ 88 // Apply basis 89 //------------------------------------------------------------------------------ 90 int CeedBasisApplyTensor_Hip_shared(CeedBasis basis, const CeedInt num_elem, CeedTransposeMode t_mode, CeedEvalMode eval_mode, CeedVector u, 91 CeedVector v) { 92 Ceed ceed; 93 Ceed_Hip *ceed_Hip; 94 CeedInt dim, num_comp; 95 const CeedScalar *d_u; 96 CeedScalar *d_v; 97 CeedBasis_Hip_shared *data; 98 99 CeedCallBackend(CeedBasisGetCeed(basis, &ceed)); 100 CeedCallBackend(CeedGetData(ceed, &ceed_Hip)); 101 CeedCallBackend(CeedBasisGetData(basis, &data)); 102 CeedCallBackend(CeedBasisGetDimension(basis, &dim)); 103 CeedCallBackend(CeedBasisGetNumComponents(basis, &num_comp)); 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 // Apply basis operation 111 switch (eval_mode) { 112 case CEED_EVAL_INTERP: { 113 CeedInt P_1d, Q_1d; 114 CeedInt block_size = data->block_sizes[0]; 115 116 CeedCallBackend(CeedBasisGetNumNodes1D(basis, &P_1d)); 117 CeedCallBackend(CeedBasisGetNumQuadraturePoints1D(basis, &Q_1d)); 118 CeedInt thread_1d = CeedIntMax(Q_1d, P_1d); 119 void *interp_args[] = {(void *)&num_elem, &data->d_interp_1d, &d_u, &d_v}; 120 121 if (dim == 1) { 122 CeedInt elems_per_block = 64 * thread_1d > 256 ? 256 / thread_1d : 64; 123 elems_per_block = elems_per_block > 0 ? elems_per_block : 1; 124 CeedInt grid = num_elem / elems_per_block + ((num_elem / elems_per_block * elems_per_block < num_elem) ? 1 : 0); 125 CeedInt shared_mem = elems_per_block * thread_1d * sizeof(CeedScalar); 126 127 if (t_mode == CEED_TRANSPOSE) { 128 CeedCallBackend(CeedRunKernelDimShared_Hip(ceed, data->InterpTranspose, grid, thread_1d, 1, elems_per_block, shared_mem, interp_args)); 129 } else { 130 CeedCallBackend(CeedRunKernelDimShared_Hip(ceed, data->Interp, grid, thread_1d, 1, elems_per_block, shared_mem, interp_args)); 131 } 132 } else if (dim == 2) { 133 // Check if required threads is small enough to do multiple elems 134 const CeedInt elems_per_block = CeedIntMax(block_size / (thread_1d * thread_1d), 1); 135 CeedInt grid = num_elem / elems_per_block + ((num_elem / elems_per_block * elems_per_block < num_elem) ? 1 : 0); 136 CeedInt shared_mem = elems_per_block * thread_1d * thread_1d * sizeof(CeedScalar); 137 138 if (t_mode == CEED_TRANSPOSE) { 139 CeedCallBackend( 140 CeedRunKernelDimShared_Hip(ceed, data->InterpTranspose, grid, thread_1d, thread_1d, elems_per_block, shared_mem, interp_args)); 141 } else { 142 CeedCallBackend(CeedRunKernelDimShared_Hip(ceed, data->Interp, grid, thread_1d, thread_1d, elems_per_block, shared_mem, interp_args)); 143 } 144 } else if (dim == 3) { 145 const CeedInt elems_per_block = CeedIntMax(block_size / (thread_1d * thread_1d), 1); 146 CeedInt grid = num_elem / elems_per_block + ((num_elem / elems_per_block * elems_per_block < num_elem) ? 1 : 0); 147 CeedInt shared_mem = elems_per_block * thread_1d * thread_1d * sizeof(CeedScalar); 148 149 if (t_mode == CEED_TRANSPOSE) { 150 CeedCallBackend( 151 CeedRunKernelDimShared_Hip(ceed, data->InterpTranspose, grid, thread_1d, thread_1d, elems_per_block, shared_mem, interp_args)); 152 } else { 153 CeedCallBackend(CeedRunKernelDimShared_Hip(ceed, data->Interp, grid, thread_1d, thread_1d, elems_per_block, shared_mem, interp_args)); 154 } 155 } 156 } break; 157 case CEED_EVAL_GRAD: { 158 CeedInt P_1d, Q_1d; 159 CeedInt block_size = data->block_sizes[1]; 160 161 CeedCallBackend(CeedBasisGetNumNodes1D(basis, &P_1d)); 162 CeedCallBackend(CeedBasisGetNumQuadraturePoints1D(basis, &Q_1d)); 163 CeedInt thread_1d = CeedIntMax(Q_1d, P_1d); 164 CeedScalar *d_grad_1d = data->d_grad_1d; 165 166 if (data->d_collo_grad_1d) { 167 d_grad_1d = data->d_collo_grad_1d; 168 } 169 void *grad_args[] = {(void *)&num_elem, &data->d_interp_1d, &d_grad_1d, &d_u, &d_v}; 170 if (dim == 1) { 171 CeedInt elems_per_block = 64 * thread_1d > 256 ? 256 / thread_1d : 64; 172 elems_per_block = elems_per_block > 0 ? elems_per_block : 1; 173 CeedInt grid = num_elem / elems_per_block + ((num_elem / elems_per_block * elems_per_block < num_elem) ? 1 : 0); 174 CeedInt shared_mem = elems_per_block * thread_1d * sizeof(CeedScalar); 175 176 if (t_mode == CEED_TRANSPOSE) { 177 CeedCallBackend(CeedRunKernelDimShared_Hip(ceed, data->GradTranspose, grid, thread_1d, 1, elems_per_block, shared_mem, grad_args)); 178 } else { 179 CeedCallBackend(CeedRunKernelDimShared_Hip(ceed, data->Grad, grid, thread_1d, 1, elems_per_block, shared_mem, grad_args)); 180 } 181 } else if (dim == 2) { 182 // Check if required threads is small enough to do multiple elems 183 const CeedInt elems_per_block = CeedIntMax(block_size / (thread_1d * thread_1d), 1); 184 CeedInt grid = num_elem / elems_per_block + ((num_elem / elems_per_block * elems_per_block < num_elem) ? 1 : 0); 185 CeedInt shared_mem = elems_per_block * thread_1d * thread_1d * sizeof(CeedScalar); 186 187 if (t_mode == CEED_TRANSPOSE) { 188 CeedCallBackend(CeedRunKernelDimShared_Hip(ceed, data->GradTranspose, grid, thread_1d, thread_1d, elems_per_block, shared_mem, grad_args)); 189 } else { 190 CeedCallBackend(CeedRunKernelDimShared_Hip(ceed, data->Grad, grid, thread_1d, thread_1d, elems_per_block, shared_mem, grad_args)); 191 } 192 } else if (dim == 3) { 193 const CeedInt elems_per_block = CeedIntMax(block_size / (thread_1d * thread_1d), 1); 194 CeedInt grid = num_elem / elems_per_block + ((num_elem / elems_per_block * elems_per_block < num_elem) ? 1 : 0); 195 CeedInt shared_mem = elems_per_block * thread_1d * thread_1d * sizeof(CeedScalar); 196 197 if (t_mode == CEED_TRANSPOSE) { 198 CeedCallBackend(CeedRunKernelDimShared_Hip(ceed, data->GradTranspose, grid, thread_1d, thread_1d, elems_per_block, shared_mem, grad_args)); 199 } else { 200 CeedCallBackend(CeedRunKernelDimShared_Hip(ceed, data->Grad, grid, thread_1d, thread_1d, elems_per_block, shared_mem, grad_args)); 201 } 202 } 203 } break; 204 case CEED_EVAL_WEIGHT: { 205 CeedInt Q_1d; 206 CeedInt block_size = data->block_sizes[2]; 207 208 CeedCheck(data->d_q_weight_1d, ceed, CEED_ERROR_BACKEND, "%s not supported; q_weights_1d not set", CeedEvalModes[eval_mode]); 209 CeedCallBackend(CeedBasisGetNumQuadraturePoints1D(basis, &Q_1d)); 210 void *weight_args[] = {(void *)&num_elem, (void *)&data->d_q_weight_1d, &d_v}; 211 212 if (dim == 1) { 213 const CeedInt opt_elems = block_size / Q_1d; 214 const CeedInt elems_per_block = opt_elems > 0 ? opt_elems : 1; 215 const CeedInt grid_size = num_elem / elems_per_block + ((num_elem / elems_per_block * elems_per_block < num_elem) ? 1 : 0); 216 217 CeedCallBackend(CeedRunKernelDim_Hip(ceed, data->Weight, grid_size, Q_1d, elems_per_block, 1, weight_args)); 218 } else if (dim == 2) { 219 const CeedInt opt_elems = block_size / (Q_1d * Q_1d); 220 const CeedInt elems_per_block = opt_elems > 0 ? opt_elems : 1; 221 const CeedInt grid_size = num_elem / elems_per_block + ((num_elem / elems_per_block * elems_per_block < num_elem) ? 1 : 0); 222 223 CeedCallBackend(CeedRunKernelDim_Hip(ceed, data->Weight, grid_size, Q_1d, Q_1d, elems_per_block, weight_args)); 224 } else if (dim == 3) { 225 const CeedInt opt_elems = block_size / (Q_1d * Q_1d); 226 const CeedInt elems_per_block = opt_elems > 0 ? opt_elems : 1; 227 const CeedInt grid_size = num_elem / elems_per_block + ((num_elem / elems_per_block * elems_per_block < num_elem) ? 1 : 0); 228 229 CeedCallBackend(CeedRunKernelDim_Hip(ceed, data->Weight, grid_size, Q_1d, Q_1d, elems_per_block, weight_args)); 230 } 231 } break; 232 case CEED_EVAL_NONE: /* handled separately below */ 233 break; 234 // LCOV_EXCL_START 235 case CEED_EVAL_DIV: 236 case CEED_EVAL_CURL: 237 return CeedError(ceed, CEED_ERROR_BACKEND, "%s not supported", CeedEvalModes[eval_mode]); 238 // LCOV_EXCL_STOP 239 } 240 241 // Restore vectors, cover CEED_EVAL_NONE 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 //------------------------------------------------------------------------------ 249 // Destroy basis 250 //------------------------------------------------------------------------------ 251 static int CeedBasisDestroy_Hip_shared(CeedBasis basis) { 252 Ceed ceed; 253 CeedBasis_Hip_shared *data; 254 255 CeedCallBackend(CeedBasisGetCeed(basis, &ceed)); 256 CeedCallBackend(CeedBasisGetData(basis, &data)); 257 CeedCallHip(ceed, hipModuleUnload(data->module)); 258 if (data->d_q_weight_1d) CeedCallHip(ceed, hipFree(data->d_q_weight_1d)); 259 CeedCallHip(ceed, hipFree(data->d_interp_1d)); 260 CeedCallHip(ceed, hipFree(data->d_grad_1d)); 261 CeedCallHip(ceed, hipFree(data->d_collo_grad_1d)); 262 CeedCallBackend(CeedFree(&data)); 263 return CEED_ERROR_SUCCESS; 264 } 265 266 //------------------------------------------------------------------------------ 267 // Create tensor basis 268 //------------------------------------------------------------------------------ 269 int CeedBasisCreateTensorH1_Hip_shared(CeedInt dim, CeedInt P_1d, CeedInt Q_1d, const CeedScalar *interp_1d, const CeedScalar *grad_1d, 270 const CeedScalar *q_ref_1d, const CeedScalar *q_weight_1d, CeedBasis basis) { 271 Ceed ceed; 272 char *basis_kernel_source; 273 const char *basis_kernel_path; 274 CeedInt num_comp; 275 const CeedInt q_bytes = Q_1d * sizeof(CeedScalar); 276 const CeedInt interp_bytes = q_bytes * P_1d; 277 CeedBasis_Hip_shared *data; 278 279 CeedCallBackend(CeedBasisGetCeed(basis, &ceed)); 280 CeedCallBackend(CeedCalloc(1, &data)); 281 282 // Copy basis data to GPU 283 if (q_weight_1d) { 284 CeedCallHip(ceed, hipMalloc((void **)&data->d_q_weight_1d, q_bytes)); 285 CeedCallHip(ceed, hipMemcpy(data->d_q_weight_1d, q_weight_1d, q_bytes, hipMemcpyHostToDevice)); 286 } 287 CeedCallHip(ceed, hipMalloc((void **)&data->d_interp_1d, interp_bytes)); 288 CeedCallHip(ceed, hipMemcpy(data->d_interp_1d, interp_1d, interp_bytes, hipMemcpyHostToDevice)); 289 CeedCallHip(ceed, hipMalloc((void **)&data->d_grad_1d, interp_bytes)); 290 CeedCallHip(ceed, hipMemcpy(data->d_grad_1d, grad_1d, interp_bytes, hipMemcpyHostToDevice)); 291 292 // Compute collocated gradient and copy to GPU 293 data->d_collo_grad_1d = NULL; 294 bool has_collocated_grad = dim == 3 && Q_1d >= P_1d; 295 296 if (has_collocated_grad) { 297 CeedScalar *collo_grad_1d; 298 299 CeedCallBackend(CeedMalloc(Q_1d * Q_1d, &collo_grad_1d)); 300 CeedCallBackend(CeedBasisGetCollocatedGrad(basis, collo_grad_1d)); 301 CeedCallHip(ceed, hipMalloc((void **)&data->d_collo_grad_1d, q_bytes * Q_1d)); 302 CeedCallHip(ceed, hipMemcpy(data->d_collo_grad_1d, collo_grad_1d, q_bytes * Q_1d, hipMemcpyHostToDevice)); 303 CeedCallBackend(CeedFree(&collo_grad_1d)); 304 } 305 306 // Set number of threads per block for basis kernels 307 CeedCallBackend(CeedBasisGetNumComponents(basis, &num_comp)); 308 CeedCallBackend(ComputeBasisThreadBlockSizes(dim, P_1d, Q_1d, num_comp, data->block_sizes)); 309 310 // Compile basis kernels 311 CeedCallBackend(CeedGetJitAbsolutePath(ceed, "ceed/jit-source/hip/hip-shared-basis-tensor.h", &basis_kernel_path)); 312 CeedDebug256(ceed, CEED_DEBUG_COLOR_SUCCESS, "----- Loading Basis Kernel Source -----\n"); 313 CeedCallBackend(CeedLoadSourceToBuffer(ceed, basis_kernel_path, &basis_kernel_source)); 314 CeedDebug256(ceed, CEED_DEBUG_COLOR_SUCCESS, "----- Loading Basis Kernel Source Complete! -----\n"); 315 CeedCallBackend(CeedCompile_Hip(ceed, basis_kernel_source, &data->module, 11, "BASIS_Q_1D", Q_1d, "BASIS_P_1D", P_1d, "T_1D", 316 CeedIntMax(Q_1d, P_1d), "BASIS_DIM", dim, "BASIS_NUM_COMP", num_comp, "BASIS_NUM_NODES", CeedIntPow(P_1d, dim), 317 "BASIS_NUM_QPTS", CeedIntPow(Q_1d, dim), "BASIS_INTERP_BLOCK_SIZE", data->block_sizes[0], "BASIS_GRAD_BLOCK_SIZE", 318 data->block_sizes[1], "BASIS_WEIGHT_BLOCK_SIZE", data->block_sizes[2], "BASIS_HAS_COLLOCATED_GRAD", 319 has_collocated_grad)); 320 CeedCallBackend(CeedGetKernel_Hip(ceed, data->module, "Interp", &data->Interp)); 321 CeedCallBackend(CeedGetKernel_Hip(ceed, data->module, "InterpTranspose", &data->InterpTranspose)); 322 CeedCallBackend(CeedGetKernel_Hip(ceed, data->module, "Grad", &data->Grad)); 323 CeedCallBackend(CeedGetKernel_Hip(ceed, data->module, "GradTranspose", &data->GradTranspose)); 324 CeedCallBackend(CeedGetKernel_Hip(ceed, data->module, "Weight", &data->Weight)); 325 CeedCallBackend(CeedFree(&basis_kernel_path)); 326 CeedCallBackend(CeedFree(&basis_kernel_source)); 327 328 CeedCallBackend(CeedBasisSetData(basis, data)); 329 330 // Register backend functions 331 CeedCallBackend(CeedSetBackendFunction(ceed, "Basis", basis, "Apply", CeedBasisApplyTensor_Hip_shared)); 332 CeedCallBackend(CeedSetBackendFunction(ceed, "Basis", basis, "Destroy", CeedBasisDestroy_Hip_shared)); 333 return CEED_ERROR_SUCCESS; 334 } 335 336 //------------------------------------------------------------------------------ 337