网格跨步循环CUDA核：如何选取最优网格维度

假设有大量独立计算任务，不想为每个工作项分配一个线程，而是让CUDA线程通过循环处理任务——可能因为Z维度超出CUDA设备限制，或是为了实现计算复用。现有如下核函数：

__global__ void kernel(dim3 size, float* arr, float div)
{
    int xoff = blockIdx.x * blockDim.x + threadIdx.x;
    int yoff = blockIdx.y * blockDim.y + threadIdx.y;
    int zoff = blockIdx.z * blockDim.z + threadIdx.z;
    int xstride = blockDim.x * gridDim.x;
    int ystride = blockDim.y * gridDim.y;
    int zstride = blockDim.z * gridDim.z;
    for(int z = zoff; z < size.z; z += zstride)
        for(int y = yoff; y < size.y; y += ystride)
            for(int x = xoff; x < size.x; x += xstride)
                arr[(z * size.y + y) * size.x + x] /= div;
}

请问如何设置核启动的网格维度，以避免工作分配不均的问题？当前解决方案如下：

void call(dim3 size, float* arr, float div)
{
    dim3 blocksize(16, 16);
    int blocksize1d = blocksize.x * blocksize.y * blocksize.z;
    int device, sm_count, sm_blocks;
    cudaGetDevice(&device);
    cudaDeviceGetAttribute(&sm_count, cudaDevAttrMultiProcessorCount, device);
    cudaOccupancyMaxActiveBlocksPerMultiprocessor(
            &sm_blocks, kernel, blocksize1d, 0 /*shared memory*/);
    int active_blocks = sm_count * sm_blocks;
    dim3 blocks;
    blocks.x = min(active_blocks, (size.x + blocksize.x - 1) / blocksize.x);
    active_blocks /= blocks.x;
    blocks.y = min(active_blocks, (size.y + blocksize.y - 1) / blocksize.y);
    active_blocks /= blocks.y;
    blocks.z = min(active_blocks, (size.z + blocksize.z - 1) / blocksize.z);
    kernel<<<blocks, blocksize>>>(size, arr, div);
}

更优/标准的核启动方法 #

当前方案的核心问题是按X→Y→Z的顺序分配活动块，容易导致维度分配失衡：若X维度所需块数占满大部分活动块，Y、Z维度的块数会被过度压缩，进而让不同线程的循环迭代次数差异过大，引发负载不均。以下是几种更合理的优化思路：

1. 基于任务量均匀分配的网格规划 #

核心目标是让每个线程（或线程块）处理的任务量尽可能接近：

• 先计算总任务数：total_tasks = size.x * size.y * size.z
• 结合线程块总线程数block_size = blocksize.x * blocksize.y * blocksize.z，得到理论所需总块数：total_blocks = (total_tasks + block_size - 1) / block_size
• 取设备支持的最大活动块数active_blocks与total_blocks的较小值作为实际启动的总块数launch_blocks
• 将launch_blocks按维度任务占比分配到三维网格，或优先保证高并行性维度的块数，避免单一维度占用过多资源

2. 优化三维网格的分配逻辑 #

基于现有occupancy计算逻辑，调整维度分配顺序或比例：

• 先计算各维度单独所需的块数：

  int x_blocks = (size.x + blocksize.x - 1) / blocksize.x;
  int y_blocks = (size.y + blocksize.y - 1) / blocksize.y;
  int z_blocks = (size.z + blocksize.z - 1) / blocksize.z;
  

• 若x_blocks * y_blocks * z_blocks <= active_blocks，直接使用该三维网格；若超过，则按比例缩减各维度块数，优先保证X、Y维度的并行度（CUDA调度更偏好XY维度的线程块并行），剩余任务由线程循环处理

3. 扁平化任务维度（推荐方案） #

将三维任务映射为一维索引，简化核函数逻辑的同时，从根源上避免负载不均：

• 修改核函数为一维任务遍历：

  __global__ void kernel(dim3 size, float* arr, float div)
  {
      int total_tasks = size.x * size.y * size.z;
      int global_idx = blockIdx.x * blockDim.x + threadIdx.x;
      int stride = gridDim.x * blockDim.x;
      for(int idx = global_idx; idx < total_tasks; idx += stride)
      {
          int z = idx / (size.x * size.y);
          int rem = idx % (size.x * size.y);
          int y = rem / size.x;
          int x = rem % size.x;
          arr[(z * size.y + y) * size.x + x] /= div;
      }
  }
  

• 启动时设置2D线程块（如16×16=256线程，32的倍数适配warp调度），计算网格维度：

  void call(dim3 size, float* arr, float div)
  {
      dim3 blocksize(16, 16);
      int block_size = blocksize.x * blocksize.y;
      int total_tasks = size.x * size.y * size.z;
      int device, sm_count, sm_blocks;
      cudaGetDevice(&device);
      cudaDeviceGetAttribute(&sm_count, cudaDevAttrMultiProcessorCount, device);
      cudaOccupancyMaxActiveBlocksPerMultiprocessor(
              &sm_blocks, kernel, block_size, 0);
      int max_active_blocks = sm_count * sm_blocks;
      int total_blocks = (total_tasks + block_size - 1) / block_size;
      dim3 blocks(min(max_active_blocks, total_blocks), 1, 1);
      kernel<<<blocks, blocksize>>>(size, arr, div);
  }
  

  这种方式下，每个线程处理的任务数差异最多为1，完美解决负载不均问题。

4. 适配硬件的线程块调整 #

线程块大小优先选择32的倍数（如256、512），保证warp内线程利用率；同时可根据硬件缓存特性调整线程块的XY比例，比如32×8的线程块可能比16×16更适合某些架构的缓存布局。

内容的提问来源于stack exchange，提问作者Homer512