CUDA中的热方程矩阵 - 非法地址错误

Question

按照this question参考官方指南中的shared memory example，我正在尝试构建热方程矩阵，就像我制作的这张画得不好的图片一样

这是我到目前为止所做的，最小的例子

#define N 32
#define BLOCK_SIZE 16
#define NUM_BLOCKS ((N + BLOCK_SIZE - 1)/ BLOCK_SIZE)

__global__ void heat_matrix(int* A)
{
    const unsigned int tid = threadIdx.x + blockIdx.x * blockDim.x; 
    __shared__ int temp_sm_A[N*N];
    int* temp_A = &temp_sm_A[0]; memset(temp_A, 0, N*N*sizeof(int));

    if (tid < N) //(*)
    {
        #pragma unroll
        for (unsigned int m = 0; m < NUM_BLOCKS; ++m)
        {           
            #pragma unroll
            for (unsigned int e = 0; e < BLOCK_SIZE ; ++e) 
            {
                if ( (tid == 0 && e == 0) || (tid == (N-1) && e == (BLOCK_SIZE-1) ) )
                {
                    temp_A[tid + (e + BLOCK_SIZE * m) * N] = -2;
                    temp_A[tid + (e + BLOCK_SIZE * m) * N + ( tid==0 ? 1 : -1 )] = 1;
                }
                if ( tid == e )
                {
                    temp_A[tid + (e + BLOCK_SIZE * m) * N - 1] = 1;
                    //printf("temp_A[%d] = 1;\n", (tid + (e + BLOCK_SIZE * m) * N -1));
                    temp_A[tid + (e + BLOCK_SIZE * m) * N] = -2;
                    //printf("temp_A[%d] = -2;\n", (tid + (e + BLOCK_SIZE * m) * N));
                    temp_A[tid + (e + BLOCK_SIZE * m) * N + 1] = 1;
                    //printf("temp_A[%d] = 1;\n", (tid + (e + BLOCK_SIZE * m) * N +1));
                }
            }
        }
        __syncthreads(); //(**)
        memcpy(A, temp_A, N*N*sizeof(int));
    }
}
int main(){
    int* h_A = (int*)malloc(N*N*sizeof(int)); memset(h_A, 0, N*N*sizeof(int));
    int* d_A; 
    checkCudaErrors(cudaMalloc((void**)&d_A, N*N*sizeof(int)));
    checkCudaErrors(cudaMemcpy(d_A, h_A, N*N*sizeof(int), cudaMemcpyHostToDevice));
    dim3 dim_grid((N/2 + BLOCK_SIZE -1)/ BLOCK_SIZE);
    dim3 dim_block(BLOCK_SIZE);

    heat_matrix <<< dim_grid, dim_block >>> (d_A);
    checkCudaErrors(cudaMemcpy(h_A, d_A, N*N*sizeof(int), cudaMemcpyDeviceToHost));
...
}

代码的排列适合较大的 N（大于 32）。我利用了块划分。执行nvcc 时产生

CUDA error at matrix.cu:102 code=77(cudaErrorIllegalAddress) "cudaMemcpy(h_A, d_A, N*N*sizeof(int), cudaMemcpyDeviceToHost)"

而cuda-memcheck只提供了一个错误（其实还有一个，不过来自cudasuccess=checkCudaErrors(cudaDeviceReset()); ...）

========= CUDA-MEMCHECK
========= Invalid __shared__ write of size 4
=========     at 0x00000cd0 in heat_matrix(int*)
=========     by thread (0,0,0) in block (0,0,0)
=========     Address 0xfffffffc is out of bounds
...

我看不出我在代码中哪里做错了。第一个块中的线程0如何引发非法访问？甚至还有具体的if案例来处理，也没有报错的代码行。

此外，我的代码有比处理所有ifs 更有效的方法吗？当然有，但我找不到更好的并行表达式来将案例拆分为第二个for。

---

在旁注中，(*) 对我来说似乎没有必要；如果我想跟随其他 GPU 函数调用，则需要 (**)。我说的对吗？

Answer 1

在您的内核中，temp_A 是指向共享内存数组开头的本地指针。考虑：

N = 32;

BLOCK_SIZE = 16;

m (0,1);

e (0,BLOCK_SIZE)

temp_A[tid + (e + BLOCK_SIZE * m) * N] 之类的访问很容易超出 1024 个元素的长数组的范围。

Answer 2

1. 看看这一行：

               temp_A[tid + (e + BLOCK_SIZE * m) * N - 1] = 1;
   

   对于第一次迭代期间tid 等于零的线程，tid + (e + BLOCK_SIZE * m) * N - 1 的计算结果为索引 -1。这正是 cuda-memcheck 输出所抱怨的（地址由于下溢而回绕）。

2. 稍后将对该行进行类似的越界访问

               temp_A[tid + (e + BLOCK_SIZE * m) * N + 1] = 1;
   

   当tid、e 和m 都取最大值时。

3. 您有多个线程写入同一内​​存位置。每个线程应该在每次内部循环迭代中写入一个数组元素。不需要写出相邻的元素，因为它们已经被自己的线程覆盖了。

4. 在初始化 memset() 和主循环内的存储之间存在竞争条件。在memset() 后面加上syncthreads()。

5. 对memset() 和memcpy() 的调用将导致每个线程执行完整的设置/复制，执行N 次操作而不是一次。
   处理这种情况的常用方法是显式写出操作，在块的线程之间划分工作。
   不过……

6. 首先在共享内存中创建矩阵，然后再将其复制到全局内存中没有任何好处。直接写入全局内存中的A，完全不需要memset()、memcpy() 和syncthreads()。

7. 使用只有 16 个线程的块大小会留下一半的资源未使用，因为线程块是以 32 个线程为单位分配的（一个扭曲）。

您可能需要重新阅读 CUDA C 编程指南中有关 Thread Hierarchy 的部分。