-1

我正在尝试实现三相并行扫描，如 Programming Massively Parallel Processors 3rd edition 的第 8 章所述（有任何代码行，但只有指令）。该算法只允许使用 1 个块中的最大线程数，并且它受到共享内存大小的限制，因为所有元素都必须适合共享内存

经过一些调试，当我使用大量元素（例如 8192 和 1 个以上的线程）时，我在阶段 3 的求和过程中遇到了一个问题（参见下面的代码）。

算法的图形概念如下： 在此处输入图像描述

在下面你可以看到内核代码：

__global__ 
void efficient_Kogge_Stone_scan_kernel(float *X, float *Y, int InputSize) {
    __shared__ float XY[SECTION_SIZE];
    __shared__ float AUS[BLOCK_DIM];
    //int i = blockIdx.x * blockDim.x + threadIdx.x;

    // Keep mind: Partition the input into blockDim.x subsections: i.e. for 8 threads --> 8 subsections

    // collaborative load in a coalesced manner
    for (int j = 0; j < SECTION_SIZE; j += blockDim.x) {
        XY[threadIdx.x + j] = X[threadIdx.x + j];
    }
    __syncthreads();


    // PHASE 1: scan inner own subsection
    // At the end of this phase the last element of each subsection contains the sum of all alements in own subsection
    for (int j = 1; j < SUBSECTION_SIZE; j++) {
        XY[threadIdx.x * (SUBSECTION_SIZE)+j] += XY[threadIdx.x * (SUBSECTION_SIZE)+j - 1];
    }
    __syncthreads();


    // PHASE 2: perform iterative kogge_stone_scan of the last elements of each subsections of XY loaded first in AUS
    AUS[threadIdx.x] = XY[threadIdx.x * (SUBSECTION_SIZE)+(SUBSECTION_SIZE)-1];
    float in;
    for (unsigned int stride = 1; stride < blockDim.x; stride *= 2) {
        __syncthreads();
        if (threadIdx.x >= stride) {
            in = AUS[threadIdx.x - stride];
        }
        __syncthreads();
        if (threadIdx.x >= stride) {
            AUS[threadIdx.x] += in;
        }
    }
    __syncthreads();


    // PHASE 3: each thread adds to its elements the new value of the last element of its predecessor's section
    if (threadIdx.x > 0) {
        for (unsigned int stride = 0; stride < (SUBSECTION_SIZE); stride++) {
            XY[threadIdx.x * (SUBSECTION_SIZE)+stride] += AUS[threadIdx.x - 1];  // <-- 
        }
    }
    __syncthreads();


    // store the result into output vector
    for (int j = 0; j < SECTION_SIZE; j += blockDim.x) {
        Y[threadIdx.x + j] = XY[threadIdx.x + j];
    }
}

如果我在块中使用一个线程和 8192 个元素，它可以完美运行，但如果我使用超过一个线程，则 XY[5793] （或 X[5793] 完成并存储结果时）的结果是错误的。它有 4096 个元素和一个或多个线程，最多 1024 个线程。如果我使用 int 而不是浮点数，它甚至适用于具有一个或多个线程的 8192 个元素。

我也尝试在 MATLAB 中进行验证，这些是输出比较：

  • X[5973] = 1678811 5 ---- MATLAB
  • X[5973] = 1678811 4 ---- CPU
  • X[5973] = 1678811 6 ---- GPU

正如我们所看到的，CPU 结果也与 MATLAB 不同，所以在这些结果之后，我认为问题出在浮点加法上，但我告诉你，我用有序的“x.00”浮点数填充了输入数组（例如 {1.00, 2.00, 3.00, 4.00 ..... 8192.00}）。

另一个问题是关于性能的，主机代码总是比内核代码快，有这些配置参数和这些输入，是否正常？

如果您需要完整的源代码，可以在这里找到

4

2 回答 2

2

8192 是 2^13

sum(1..8192) 接近 8192^2/2: 8192*8193/2，比 2^25 多一点。

因此，您需要 26 位来表示它（参见下面的注释）。

单精度 IEEE 754 浮点数只有 24 位有效位，因此，取决于求和的执行方式（按顺序），以及最终的舍入方向（通常是默认舍入到最接近，平到偶数），结果可能会有所不同。

请注意，严格来说，精确的和可以用浮点数表示而无需四舍五入，因为最后 12 位为零，因此有效位仅跨越 14 位。但部分总和并非如此。

于 2017-04-28T08:17:28.680 回答
1

第一次扫描可能有问题：

XY[threadIdx.x * (SUBSECTION_SIZE)+j] += XY[threadIdx.x * (SUBSECTION_SIZE)+j - 1];

这可能导致共享内存中元素的读取不一致。当您阅读前一个元素时，不能保证任何其他线程都没有更新该值。

尝试通过将值存储在寄存器中来将此部分分成几部分。例子：

int t =  XY[threadIdx.x * (SUBSECTION_SIZE)+j - 1];
 __syncthreads();
 XY[threadIdx.x * (SUBSECTION_SIZE)+j] += t;
于 2017-04-27T16:03:50.780 回答