各位同事,想请教一个关于CUDA归约的问题

  • 7 replies
  • 548 views
各位同事,想请教一个关于CUDA归约的问题
« 于: 六月 21, 2021, 03:40:41 pm »
在学习cuda的过程中遇到了一个有关于并行归约的困难,希望有朋友能解答一下,不胜感激。我遇到的问题如下:在展开归约时(展开因子为8),末尾的几个数据总是丢失。导致我的归约结果比正确值小一些。代码如下:主机端:
               cudaMemcpy(dev_p, p, N * sizeof(float), cudaMemcpyHostToDevice);
               cudaMemcpy(dev_re12, re12, N * sizeof(float), cudaMemcpyHostToDevice);
               RN << <blocksPerGrid/8, threadsPerBlock >> >(dev_p, dev_re12, dev_partial_rn);
               cudaMemcpy(partial_rn, dev_partial_rn, (blocksPerGrid/8) * sizeof(float), cudaMemcpyDeviceToHost);
               rn = 0;
               for (int i = 0; i<blocksPerGrid/8; i++)
                {
                    rn += partial_rn;
                }
其中threadsPerBlock=1024,blocksPerGrid=128, N=34329.
设备端:
__global__ void RN(float *p, float *re12, float *rn)
{
        __shared__ float cache[threadsPerBlock];
        int tid = threadIdx.x + blockIdx.x * blockDim.x * 8;
        int cacheIndex = threadIdx.x;
        float temp = 0;
        if(tid + 7 * blockDim.x < N)
        {       
                float a1 = p[tid] * re12[tid];
                float a2 = p[tid + blockDim.x] * re12[tid + blockDim.x];
                float a3 = p[tid + 2 * blockDim.x] * re12[tid + 2 * blockDim.x];
                float a4 = p[tid + 3 * blockDim.x] * re12[tid + 3 * blockDim.x];
                float a5 = p[tid + 4 * blockDim.x] * re12[tid + 4 * blockDim.x];
                float a6 = p[tid + 5 * blockDim.x] * re12[tid + 5 * blockDim.x];
                float a7 = p[tid + 6 * blockDim.x] * re12[tid + 6 * blockDim.x];
                float a8 = p[tid + 7 * blockDim.x] * re12[tid + 7 * blockDim.x];
                temp = a1 + a2 + a3 + a4 + a5 + a6 + a7 + a8;
        }
        cache[cacheIndex] = temp;
        // synchronize threads in this block
        __syncthreads();
        //unrolling warp
        if (blockDim.x >= 1024 && cacheIndex < 512)
        {
                cache[cacheIndex] += cache[cacheIndex + 512];
        }
        __syncthreads();
        if (blockDim.x >= 512 && cacheIndex < 256)
        {
                cache[cacheIndex] += cache[cacheIndex + 256];
        }
        __syncthreads();
        if (blockDim.x >= 256 && cacheIndex < 128)
        {
                cache[cacheIndex] += cache[cacheIndex + 128];
        }
        __syncthreads();
        if (blockDim.x >= 128 && cacheIndex < 64)
        {
                cache[cacheIndex] += cache[cacheIndex + 64];
        }
        __syncthreads();
        if (cacheIndex < 32)
        {
                volatile float *vcache = cache;
                vcache[cacheIndex] += vcache[cacheIndex + 32];
                vcache[cacheIndex] += vcache[cacheIndex + 16];
                vcache[cacheIndex] += vcache[cacheIndex + 8];
                vcache[cacheIndex] += vcache[cacheIndex + 4];
                vcache[cacheIndex] += vcache[cacheIndex + 2];
                vcache[cacheIndex] += vcache[cacheIndex + 1];
        }
        if (cacheIndex == 0)
                rn[blockIdx.x] = cache[0];
}

被这个问题困扰很久了,课题止步不前,请有能力的前辈不吝赐教,再次感谢~!

Re: 各位同事,想请教一个关于CUDA归约的问题
« 回复 #1 于: 六月 22, 2021, 07:20:24 pm »
平台版本是cuda10,硬件是2060。请有经验的前辈指点迷津。感谢!

Re: 各位同事,想请教一个关于CUDA归约的问题
« 回复 #2 于: 六月 22, 2021, 07:36:03 pm »
发现一个问题:在注释掉以下程序后,结果仍不变。注释掉的程序如下:
 if (blockDim.x >= 1024 && cacheIndex < 512)
        {
                cache[cacheIndex] += cache[cacheIndex + 512];
        }
        __syncthreads();
晚辈我分段输出后,认为这段程序也正好对应少了的那段数据。这是为什么呢?

Re: 各位同事,想请教一个关于CUDA归约的问题
« 回复 #3 于: 六月 23, 2021, 09:13:20 pm »
大概找到问题出在哪里了,
   if(tid + 7 * blockDim.x <= N)
   {   
      float a1 = p[tid] * re12[tid];
      float a2 = p[tid + blockDim.x] * re12[tid + blockDim.x];
      float a3 = p[tid + 2 * blockDim.x] * re12[tid + 2 * blockDim.x];
      float a4 = p[tid + 3 * blockDim.x] * re12[tid + 3 * blockDim.x];
      float a5 = p[tid + 4 * blockDim.x] * re12[tid + 4 * blockDim.x];
      float a6 = p[tid + 5 * blockDim.x] * re12[tid + 5 * blockDim.x];
      float a7 = p[tid + 6 * blockDim.x] * re12[tid + 6 * blockDim.x];
      float a8 = p[tid + 7 * blockDim.x] * re12[tid + 7 * blockDim.x];
      //temp = a1 + a2 + a3 + a4;
      temp = a1 + a2 + a3 + a4 + a5 + a6 + a7 + a8;
        }
这段里面,N原本是34329,但是到了32768之后数据戛然而止,没有了。我已经查过输入到cuda中的数组,所有数据都传过来了没有问题。这是为何呢?

Re: 各位同事,想请教一个关于CUDA归约的问题
« 回复 #4 于: 六月 23, 2021, 09:26:16 pm »
32768=512*64,难道这是cuda的极限?应该不会这么弱吧,

Re: 各位同事,想请教一个关于CUDA归约的问题
« 回复 #5 于: 六月 25, 2021, 12:56:19 pm »
在学习cuda的过程中遇到了一个有关于并行归约的困难,希望有朋友能解答一下,不胜感激。我遇到的问题如下:在展开归约时(展开因子为8),末尾的几个数据总是丢失。导致我的归约结果比正确值小一些。代码如下:主机端:
               cudaMemcpy(dev_p, p, N * sizeof(float), cudaMemcpyHostToDevice);
               cudaMemcpy(dev_re12, re12, N * sizeof(float), cudaMemcpyHostToDevice);
               RN << <blocksPerGrid/8, threadsPerBlock >> >(dev_p, dev_re12, dev_partial_rn);
               cudaMemcpy(partial_rn, dev_partial_rn, (blocksPerGrid/8) * sizeof(float), cudaMemcpyDeviceToHost);
               rn = 0;
               for (int i = 0; i<blocksPerGrid/8; i++)
                {
                    rn += partial_rn;
                }
其中threadsPerBlock=1024,blocksPerGrid=128, N=34329.
设备端:
__global__ void RN(float *p, float *re12, float *rn)
{
        __shared__ float cache[threadsPerBlock];
        int tid = threadIdx.x + blockIdx.x * blockDim.x * 8;
        int cacheIndex = threadIdx.x;
        float temp = 0;
        if(tid + 7 * blockDim.x < N)
        {       
                float a1 = p[tid] * re12[tid];
                float a2 = p[tid + blockDim.x] * re12[tid + blockDim.x];
                float a3 = p[tid + 2 * blockDim.x] * re12[tid + 2 * blockDim.x];
                float a4 = p[tid + 3 * blockDim.x] * re12[tid + 3 * blockDim.x];
                float a5 = p[tid + 4 * blockDim.x] * re12[tid + 4 * blockDim.x];
                float a6 = p[tid + 5 * blockDim.x] * re12[tid + 5 * blockDim.x];
                float a7 = p[tid + 6 * blockDim.x] * re12[tid + 6 * blockDim.x];
                float a8 = p[tid + 7 * blockDim.x] * re12[tid + 7 * blockDim.x];
                temp = a1 + a2 + a3 + a4 + a5 + a6 + a7 + a8;
        }
        cache[cacheIndex] = temp;
        // synchronize threads in this block
        __syncthreads();
        //unrolling warp
        if (blockDim.x >= 1024 && cacheIndex < 512)
        {
                cache[cacheIndex] += cache[cacheIndex + 512];
        }
        __syncthreads();
        if (blockDim.x >= 512 && cacheIndex < 256)
        {
                cache[cacheIndex] += cache[cacheIndex + 256];
        }
        __syncthreads();
        if (blockDim.x >= 256 && cacheIndex < 128)
        {
                cache[cacheIndex] += cache[cacheIndex + 128];
        }
        __syncthreads();
        if (blockDim.x >= 128 && cacheIndex < 64)
        {
                cache[cacheIndex] += cache[cacheIndex + 64];
        }
        __syncthreads();
        if (cacheIndex < 32)
        {
                volatile float *vcache = cache;
                vcache[cacheIndex] += vcache[cacheIndex + 32];
                vcache[cacheIndex] += vcache[cacheIndex + 16];
                vcache[cacheIndex] += vcache[cacheIndex + 8];
                vcache[cacheIndex] += vcache[cacheIndex + 4];
                vcache[cacheIndex] += vcache[cacheIndex + 2];
                vcache[cacheIndex] += vcache[cacheIndex + 1];
        }
        if (cacheIndex == 0)
                rn[blockIdx.x] = cache[0];
}

被这个问题困扰很久了,课题止步不前,请有能力的前辈不吝赐教,再次感谢~!

你要考虑到数据和线程/blocks数量之间能否在特定的对应下,整除的问题。
例如最后多出来1个block用来处理后续的残余数据。例如可能会多出来一些线程/少出来一些数据的处理。

至于只能处理32K个,不存在这种极限。

以及,直接抄手册上的代码更快。网上这种经典代码能找出多种实现(你的两步样式的,第二步上单一block或者CPU扫尾;用原子操作一步到位的;不用原子操作单步的)。直接抄?
« 最后编辑时间: 六月 25, 2021, 12:57:27 pm 作者 屠戮人神 »

Re: 各位同事,想请教一个关于CUDA归约的问题
« 回复 #6 于: 六月 25, 2021, 01:01:37 pm »
你要考虑到数据和线程/blocks数量之间能否在特定的对应下,整除的问题。
例如最后多出来1个block用来处理后续的残余数据。例如可能会多出来一些线程/少出来一些数据的处理。

至于只能处理32K个,不存在这种极限。

以及,直接抄手册上的代码更快。网上这种经典代码能找出多种实现(你的两步样式的,第二步上单一block或者CPU扫尾;用原子操作一步到位的;不用原子操作单步的)。直接抄?

我举个例子说,最后最后一个对应的block:
很可能不能满足:  if(tid + 7 * blockDim.x < N)的条件(整除方面的原因)
但是会满足:
float a1 = p[tid] * re12[tid]; <--你假定的完全能进入的body
float a2 = p[tid + blockDim.x] * re12[tid + blockDim.x];
float a3 = p[tid + 2 * blockDim.x] * re12[tid + 2 * blockDim.x];
float a4 = p[tid + 3 * blockDim.x] * re12[tid + 3 * blockDim.x];
float a5 = p[tid + 4 * blockDim.x] * re12[tid + 4 * blockDim.x];
<---->实际可能a1 - a5都可以,但是条件是a8, 就会漏掉
float a6 = p[tid + 5 * blockDim.x] * re12[tid + 5 * blockDim.x];
float a7 = p[tid + 6 * blockDim.x] * re12[tid + 6 * blockDim.x];
float a8 = p[tid + 7 * blockDim.x] * re12[tid + 7 * blockDim.x];
(你假定以a8的索引条件判断)

Re: 各位同事,想请教一个关于CUDA归约的问题
« 回复 #7 于: 六月 25, 2021, 04:05:24 pm »
非常感谢屠戮人神前辈的帮助,您的建议非常有用!感谢!