您在已删除答案中的主要结论是正确的:您发布的内核不理解这样一个事实,即在内核执行结束时,您已经完成了大量的整体减少,但结果并不完全。必须(以某种方式)组合每个块的结果。正如评论中所指出的,您的代码还有一些其他问题。让我们看一下它的修改版本:
__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]);
}
- 正如 Pavan 所指出的,您需要初始化共享内存数组。
gridDim.x*blockDim.x
如果大于,则启动的最后一个块可能不是“完整”块elements
。
- 请注意,在这一行中,即使我们正在检查线程操作 (
gid
) 是否小于elements
,当我们添加s
用于gid
索引到共享内存中时,我们仍然可以在最后一个块中复制到共享内存中的合法值之外进行索引。因此,我们需要注释 1 中指出的共享内存初始化。
- 正如您已经发现的那样,您的最后一行不正确。每个块产生它自己的结果,我们必须以某种方式组合它们。如果启动的块数量很少(稍后会详细介绍),您可能会考虑的一种方法是使用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 示例一起提供。