9

我一直在学习 Cuda,但我仍在掌握并行性。我目前遇到的问题是对一组值实现最大减少。这是我的内核

__global__ void max_reduce(const float* const d_array,
                     float* d_max,
                     const size_t elements)
{
    extern __shared__ float shared[];

    int tid = threadIdx.x;
    int gid = (blockDim.x * blockIdx.x) + tid;

    if (gid < elements)
        shared[tid] = d_array[gid];
    __syncthreads();

    for (unsigned int s=blockDim.x/2; s>0; s>>=1) 
    {
        if (tid < s && gid < elements)
            shared[tid] = max(shared[tid], shared[tid + s]);
        __syncthreads();
    }

    if (gid == 0)
        *d_max = shared[tid];
}

我使用相同的方法(用 min 替换 max 函数)实现了 min reduce,效果很好。

为了测试内核,我使用串行 for 循环找到了最小值和最大值。内核中的最小值和最大值总是相同的,但只有最小值减少匹配。

有什么明显的我错过/做错了吗?

4

2 回答 2

19

您在已删除答案中的主要结论是正确的:您发布的内核不理解这样一个事实,即在内核执行结束时,您已经完成了大量的整体减少,但结果并不完全。必须(以某种方式)组合每个块的结果。正如评论中所指出的,您的代码还有一些其他问题。让我们看一下它的修改版本:

__device__ float atomicMaxf(float* address, float val)
{
    int *address_as_int =(int*)address;
    int old = *address_as_int, assumed;
    while (val > __int_as_float(old)) {
        assumed = old;
        old = atomicCAS(address_as_int, assumed,
                        __float_as_int(val));
        }
    return __int_as_float(old);
}


__global__ void max_reduce(const float* const d_array, float* d_max, 
                                              const size_t elements)
{
    extern __shared__ float shared[];

    int tid = threadIdx.x;
    int gid = (blockDim.x * blockIdx.x) + tid;
    shared[tid] = -FLOAT_MAX;  // 1

    if (gid < elements)
        shared[tid] = d_array[gid];
    __syncthreads();

    for (unsigned int s=blockDim.x/2; s>0; s>>=1) 
    {
        if (tid < s && gid < elements)
            shared[tid] = max(shared[tid], shared[tid + s]);  // 2
        __syncthreads();
    }
    // what to do now?
    // option 1: save block result and launch another kernel
    if (tid == 0)        
        d_max[blockIdx.x] = shared[tid]; // 3
    // option 2: use atomics
    if (tid == 0)
      atomicMaxf(d_max, shared[0]);
}
  1. 正如 Pavan 所指出的,您需要初始化共享内存数组。gridDim.x*blockDim.x如果大于,则启动的最后一个块可能不是“完整”块elements
  2. 请注意,在这一行中,即使我们正在检查线程操作 ( gid) 是否小于elements当我们添加s用于gid索引到共享内存中时,我们仍然可以在最后一个块中复制到共享内存中的合法值之外进行索引。因此,我们需要注释 1 中指出的共享内存初始化。
  3. 正如您已经发现的那样,您的最后一行不正确。每个块产生它自己的结果,我们必须以某种方式组合它们。如果启动的块数量很少(稍后会详细介绍),您可能会考虑的一种方法是使用atomics。通常我们会引导人们远离使用原子,因为它们在执行时间方面是“昂贵的”。然而,我们面临的另一个选择是将块结果保存在全局内存中,完成内核,然后可能启动另一个内核来组合单个块结果。如果我最初启动了大量块(比如超过 1024 个),那么如果我遵循这种方法,我最终可能会启动两个额外的内核。因此考虑原子。如前所述,没有原生atomicMax用于浮点数的函数,但如文档中所述,您可以使用它atomicCAS来生成任意原子函数,我提供了一个示例,atomicMaxf其中为float.

但是运行 1024 个或更多原子函数(每个块一个)是最好的方法吗?可能不是。

当启动线程块内核时,我们真的只需要启动足够多的线程块来保持机器忙碌。根据经验,我们希望每个 SM 至少运行 4-8 个经线,多一些可能是个好主意。但是从机器利用率的角度来看,最初启动数千个线程块并没有特别的好处。如果我们选择一个数字,例如每个 SM 8 个线程块,并且我们的 GPU 中最多有 14-16 个 SM,这给我们提供了相对较少的 8*14 = 112 个线程块。让我们选择 128 (8*16) 作为一个不错的整数。这并没有什么神奇之处,它足以让 GPU 保持忙碌。如果我们让这 128 个线程块中的每一个都做额外的工作来解决整个问题问题,然后我们可以利用我们对原子的使用而无需(也许)为此付出太多代价,并避免多次内核启动。那么这看起来如何?:

__device__ float atomicMaxf(float* address, float val)
{
    int *address_as_int =(int*)address;
    int old = *address_as_int, assumed;
    while (val > __int_as_float(old)) {
        assumed = old;
        old = atomicCAS(address_as_int, assumed,
                        __float_as_int(val));
        }
    return __int_as_float(old);
}


__global__ void max_reduce(const float* const d_array, float* d_max, 
                                              const size_t elements)
{
    extern __shared__ float shared[];

    int tid = threadIdx.x;
    int gid = (blockDim.x * blockIdx.x) + tid;
    shared[tid] = -FLOAT_MAX; 

    while (gid < elements) {
        shared[tid] = max(shared[tid], d_array[gid]);
        gid += gridDim.x*blockDim.x;
        }
    __syncthreads();
    gid = (blockDim.x * blockIdx.x) + tid;  // 1
    for (unsigned int s=blockDim.x/2; s>0; s>>=1) 
    {
        if (tid < s && gid < elements)
            shared[tid] = max(shared[tid], shared[tid + s]);
        __syncthreads();
    }

    if (tid == 0)
      atomicMaxf(d_max, shared[0]);
}

elements使用这个修改过的内核,在创建内核启动时,我们不会根据整体数据大小( )来决定要启动多少线程块。相反,我们正在启动固定数量的块(例如,128,您可以修改这个数字以找出运行最快的块),并让每个线程块(以及整个网格)循环通过内存,计算每个元素的部分最大操作共享内存。然后,在注释 1 标记的行中,我们必须将变量重新设置gid为其初始值。这实际上是不必要的,如果我们保证网格 ( gridDim.x*blockDim.x) 的大小小于elements,则可以进一步简化块缩减循环代码,这在内核启动时并不难做到。

请注意,使用此原子方法时,有必要将结果(*d_max在这种情况下)初始化为适当的值,例如-FLOAT_MAX.

同样,我们通常会引导人们远离原子使用,但在这种情况下,如果我们仔细管理它是值得考虑的,它可以让我们节省额外内核启动的开销。

有关如何进行快速并行缩减的忍者级分析,请查看 Mark Harris 的优秀白皮书,该白皮书与相关CUDA 示例一起提供。

于 2013-06-29T15:40:36.107 回答
-1

这是一个看似幼稚但并非如此的问题。这不会推广到其他函数,例如sum(),但它适用于min()max()

__device__ const float float_min = -3.402e+38;

__global__ void maxKernel(float* d_data)
{ 
    // compute max over all threads, store max in d_data[0]
    int i = threadIdx.x;
    __shared__ float max_value;

    if (i == 0) max_value = float_min;
    float v = d_data[i];
    __syncthreads();

    while (max_value < v) max_value = v;

    __syncthreads();
    if (i == 0) d_data[0] = max_value;
}

是的,没错,只在初始化后同步一次,在写入结果前同步一次。该死的比赛条件!全速前进!

在你告诉我它不会工作之前,请先试一试。我已经彻底测试过,它每次都适用于各种任意内核大小。事实证明,在这种情况下,竞争条件无关紧要,因为 while 循环解决了它。

它的工作速度明显快于传统的减少。另一个令人惊讶的是,内核大小为 32 的平均通过次数为 4。是的,这是 (log(n)-1),这似乎违反直觉。这是因为比赛条件提供了好运的机会。除了消除传统减少的开销之外,还有这个好处。

对于较大的 n,无法避免每个 warp 至少进行一次迭代,但该迭代仅涉及一个比较操作,当 max_value 位于分布的高端时,该操作通常在整个 warp 中立即为假。您可以修改它以使用多个 SM,但这会大大增加总工作量并增加通信成本,因此不太可能有帮助。

为简洁起见,我省略了大小和输出参数。大小只是线程数(可以是 137 或任何你喜欢的)。输出以d_data[0].

我在这里上传了工作文件:https ://github.com/kenseehart/YAMR

于 2021-02-28T13:43:57.337 回答