我有一个矩阵，我想使用CUDA并以最快的方式计算列方式的均值（简化为简单的总和），即返回包含该矩阵中每列的平均值的行向量 . 用于计算单列向量之和的总和减少实现如下所示：
template<typename T>
__global__ void kernelSum(const T* __restrict__ input, T* __restrict__ per_block_results, const size_t n) {
extern __shared__ T sdata[];
size_t tid = blockIdx.x * blockDim.x + threadIdx.x;
// load input into __shared__ memory
T x = 0.0;
if (tid < n) {
x = input[tid];
}
sdata[threadIdx.x] = x;
__syncthreads();
// contiguous range pattern
for(int offset = blockDim.x / 2; offset > 0; offset >>= 1) {
if(threadIdx.x < offset) {
// add a partial sum upstream to our own
sdata[threadIdx.x] += sdata[threadIdx.x + offset];
}
// wait until all threads in the block have
// updated their partial sums
__syncthreads();
}
// thread 0 writes the final result
if(threadIdx.x == 0) {
per_block_results[blockIdx.x] = sdata[0];
}
}
这被调用为：
int n = ... // vector size
const int BLOCK_SIZE = 1024;
int number_of_blocks = (n + BLOCK_SIZE - 1) / BLOCK_SIZE;
double* per_block_results = NULL;
cudaMalloc((void**) &per_block_results, sizeof(double)*(number_of_blocks + 1));
// launch one kernel to compute, per-block, a partial sum
kernelSum<double> <<<number_of_blocks, BLOCK_SIZE, BLOCK_SIZE*sizeof(double)>>>(a, per_block_results, n);
// launch a single block to compute the sum of the partial sums
kernelSum<double> <<<1, number_of_blocks, number_of_blocks*sizeof(double)>>>(per_block_results, per_block_results + number_of_blocks, number_of_blocks);
我可以将这个内核推广到任意数量的列的矩阵，但我受共享内存的限制 . 我的GPU具有计算能力 3.5 ，因此它具有 48KB 共享内存，最大块大小为 1024 ，即每个块的线程数 . 由于我对双精度感兴趣，我有 48*1024/8= 6144 共享内存的最大双精度数 . 由于每个块都进行了减少，因此我可以最多使用 6144 (doubles in shared memory) / 1024 (block size) = 6 列，我可以同时计算减少的总和 . 然后减小块大小将允许同时计算更多列，例如 6144 (doubles in shared memory) / 512 (block size) = 12 .
这种更复杂的策略是否会超过矩阵每列的简单CPU循环并调用总和减少量 . 还有另一种更好的方法吗？
2 回答
什么阻止你做这样的事情：
[标准免责声明：用浏览器编写，从未编译或测试，使用风险自负]
即 . 对于共享内存缩减，块中的每个线程只需要
sdata中的一个条目 . 每个线程总和所需的值以覆盖整列输入 . 然后没有共享内存限制，您可以使用相同的块大小对任何大小的列求和 .编辑：显然使用每个线程的部分总和的想法对你来说是新的，所以这里有一个完整的例子来研究：
您应该尝试相对于矩阵大小的块大小以获得最佳性能，但通常内核每个线程的工作量越多，整体性能就越好（因为共享内存减少非常昂贵） . 您可以在this answer中看到一种阻止和网格大小启发式的方法，用于类似的内存带宽限制问题 .
作为talonmies已经提供的答案的替代方案，我在这里报告
4其他减少列的方法，3基于使用CUDA Thrust和1基于使用cublas<t>gemv()和1列，如我在上面的评论中所建议的 .CUDA推力方法类似于Reduce matrix rows with CUDA，通过获得隐式转置
这是完整的代码：