c++ - CUDA-优化表面检测内核

标签 c++ cuda

我正在尝试优化我的表面检测内核；给定输入二进制512w x 1024h图像，我想找到图像中的第一个亮表面。我编写的代码声明了512个线程，并在3x3附近搜索了第一个亮像素。该代码可以正常工作，但是~9.46 ms有点慢，我想使其运行得更快。

编辑1:，性能提高不到我的原始内核运行时间的一半。 Robert的内核在我的Quadro K6000上以4.032 ms运行。

编辑2:设法通过将线程数减少一半来进一步提高性能。现在，我的(Robert修改过的)内核在Quadro K6000上以2.125 ms运行。

使用以下命令调用内核:
firstSurfaceDetection <<< 1, 512 >>> (threshImg, firstSurfaceImg, actualImHeight, actualImWidth);
我想使用共享内存来改善内存获取；关于如何优化此代码补丁的任何想法？

__global__ void firstSurfaceDetection (float *threshImg, float *firstSurfaceImg, int height, int width) {

int col = threadIdx.x + (blockDim.x*blockIdx.x); 
int rows2skip = 10; 
float thresh = 1.0f;

 //thread Index: (0 -> 511)

if (col < width) {

    if( col == 0 ) { // first col - 0
        for (int row = 0 + rows2skip; row < height - 2; row++) { // skip first 30 rows
            int cnt = 0;
             float neibs[6]; // not shared mem as it reduces speed  

            // get six neighbours - three in same col, and three to the right 
            neibs[0] = threshImg[((row)*width) +(col)];             if(neibs[0] == thresh) { cnt++; }   // current position
            neibs[1] = threshImg[((row)*width) +(col+1)];           if(neibs[1] == thresh) { cnt++; }   // right
            neibs[2] = threshImg[((row+1)*width) +(col)];           if(neibs[2] == thresh) { cnt++; }   // bottom
            neibs[3] = threshImg[((row+1)*width) +(col+1)];         if(neibs[3] == thresh) { cnt++; }   // bottom right
            neibs[4] = threshImg[((row+2)*width) +(col)];           if(neibs[4] == thresh) { cnt++; }   // curr offset by 2 - bottom
            neibs[5] = threshImg[((row+2)*width) +(col+1)];         if(neibs[5] == thresh) { cnt++; }   // curr offset by 2 - bottom right

            if(cnt == 6) { // if all neighbours are bright, we are at the edge boundary
                firstSurfaceImg[(row)*width + col] = 1.0f;
                row = height;
            }
        }
    }

    else if ( col == (width-1) ) { // last col 
        for (int row = 0 + rows2skip; row < height -2; row++) { 
            int cnt = 0;
             float neibs[6]; // not shared mem as it reduces speed  

            // get six neighbours - three in same col, and three to the left
            neibs[0] = threshImg[((row)*width) +(col)];             if(neibs[0] == thresh) { cnt++; }   // current position
            neibs[1] = threshImg[((row)*width) +(col-1)];           if(neibs[1] == thresh) { cnt++; }   // left
            neibs[2] = threshImg[((row+1)*width) +(col)];           if(neibs[2] == thresh) { cnt++; }   // bottom
            neibs[3] = threshImg[((row+1)*width) +(col-1)];         if(neibs[3] == thresh) { cnt++; }   // bottom left
            neibs[4] = threshImg[((row+2)*width) +(col)];           if(neibs[4] == thresh) { cnt++; }   // curr offset by 2 - bottom
            neibs[5] = threshImg[((row+2)*width) +(col-1)];         if(neibs[5] == thresh) { cnt++; }   // curr offset by 2 - bottom left

            if(cnt == 6) { // if all neighbours are bright, we are at the edge boundary
                firstSurfaceImg[(row)*width + col] = 1.0f;
                row = height;
            }
        }       
    }

    // remaining threads are: (1 -> 510) 

    else { // any col other than first or last column
        for (int row = 0 + rows2skip; row < height - 2; row++) { 

            int cnt = 0;
            float neibs[9]; // not shared mem as it reduces speed   

            // for threads < width/4, get the neighbors
            // get nine neighbours - three in curr col, three each to left and right
            neibs[0] = threshImg[((row)*width) +(col-1)];           if(neibs[0] == thresh) { cnt++; } 
            neibs[1] = threshImg[((row)*width) +(col)];             if(neibs[1] == thresh) { cnt++; } 
            neibs[2] = threshImg[((row)*width) +(col+1)];           if(neibs[2] == thresh) { cnt++; }           
            neibs[3] = threshImg[((row+1)*width) +(col-1)];         if(neibs[3] == thresh) { cnt++; }           
            neibs[4] = threshImg[((row+1)*width) +(col)];           if(neibs[4] == thresh) { cnt++; }           
            neibs[5] = threshImg[((row+1)*width) +(col+1)];         if(neibs[5] == thresh) { cnt++; }           
            neibs[6] = threshImg[((row+2)*width) +(col-1)];         if(neibs[6] == thresh) { cnt++; }           
            neibs[7] = threshImg[((row+2)*width) +(col)];           if(neibs[7] == thresh) { cnt++; }           
            neibs[8] = threshImg[((row+2)*width) +(col+1)];         if(neibs[8] == thresh) { cnt++; }

            if(cnt == 9) { // if all neighbours are bright, we are at the edge boundary

                firstSurfaceImg[(row)*width + col] = 1.0f;
                row = height;
                }
            }
        }       
    }           

__syncthreads();
}

最佳答案

这是一个工作示例，演示了注释中讨论的3个概念中的2个:

  • 首先要考虑的优化是512个线程不足以使任何GPU繁忙。我们希望定位10000个或更多线程。 GPU是一个隐藏延迟的机器，当您的线程不足以帮助GPU隐藏延迟时，内核就会受到延迟的限制，这是一种内存受限的问题。最简单的方法是让每个线程在图像中处理一个像素(总共允许512 * 1024个线程)，而不是一列(仅允许总共512个线程)。但是，由于这似乎“破坏”了我们的“第一表面检测”算法，因此我们也必须进行另一种修改(2)。
  • 一旦我们并行处理了所有像素，则上面第1项的直接修改意味着我们不再知道哪个表面是“第一个”，即哪个“明亮”表面(每列)最接近第0行。该算法将问题从简单的转换更改为缩小(实际上，图像的每列减少一遍。)我们将通过为每个像素分配1个线程来并行处理每一列，但我们将选择结果满足最接近零行的亮度测试的像素。相对简单的方法是在发现合适的明亮像素邻域的最小行(每列)的每列一个数组中使用atomicMin。

    • 以下代码演示了上述2种更改(并且不包括共享内存的任何使用)，并演示了(对我而言，在Tesla K40上)与PC原始内核相比，加速范围是1x-20x。加速范围是由于算法的工作取决于图像结构而变化的事实。两种算法都有提早退出策略。由于for循环上的早期退出结构，原始内核可以完成或多或少的工作，具体取决于每列中发现“明亮”像素邻域的位置(如果有)。因此，如果所有列在第0行附近都具有明亮的邻域，我会看到大约1倍的“改进”(即，我的内核以与原始内核相同的速度运行)。如果所有列(仅)在图像的另一“端”附近都具有明亮的邻域，则可以看到大约提高了20倍。这可能会因GPU而异，因为开普勒GPU已提高了我正在使用的全局原子吞吐量。 编辑:由于这种可变的工作方式，我添加了一个粗略的“提前退出”策略作为对我的代码的琐碎修改。这使最短的执行时间接近两个内核之间的奇偶校验(即大约1倍)。

    其余的优化可能包括:
  • 使用共享内存。这应该是对例如矩阵乘法所使用的相同的基于图块的共享内存方法的简单改编。如果您使用方形方块，那么您将需要调整内核块/网格的尺寸以使其成为“方形方块”。
  • 改进的还原技术。对于某些图像结构，原子方法可能会比较慢。可以通过切换为每列适当的并行缩减来改善这一点。通过将测试镜像设置为各处的阈值，可以对我的内核进行“最坏情况”测试。这应该导致原始内核运行最快，而我的内核运行最慢，但是在这种情况下，我没有观察到我的内核有任何明显的减速。我的内核的执行时间非常稳定。同样，这可能取决于GPU。

    • 例:
    #include <stdlib.h>
#include <stdio.h>

#define SKIP_ROWS 10
#define THRESH 1.0f

#include <time.h>
#include <sys/time.h>
#define USECPSEC 1000000ULL

unsigned long long dtime_usec(unsigned long long start){

  timeval tv;
  gettimeofday(&tv, 0);
  return ((tv.tv_sec*USECPSEC)+tv.tv_usec)-start;
}

__global__ void firstSurfaceDetection (float *threshImg, float *firstSurfaceImg, int height, int width) {

int col = threadIdx.x + (blockDim.x*blockIdx.x); 
int rows2skip = SKIP_ROWS; 
float thresh = THRESH;

 //thread Index: (0 -> 511)

if (col < width) {

    if( col == 0 ) { // first col - 0
        for (int row = 0 + rows2skip; row < height; row++) { // skip first 30 rows
            int cnt = 0;
             float neibs[6]; // not shared mem as it reduces speed  

            // get six neighbours - three in same col, and three to the right 
            neibs[0] = threshImg[((row)*width) +(col)];             if(neibs[0] == thresh) { cnt++; }   // current position
            neibs[1] = threshImg[((row)*width) +(col+1)];           if(neibs[1] == thresh) { cnt++; }   // right
            neibs[2] = threshImg[((row+1)*width) +(col)];           if(neibs[2] == thresh) { cnt++; }   // bottom
            neibs[3] = threshImg[((row+1)*width) +(col+1)];         if(neibs[3] == thresh) { cnt++; }   // bottom right
            neibs[4] = threshImg[((row+2)*width) +(col)];           if(neibs[4] == thresh) { cnt++; }   // curr offset by 2 - bottom
            neibs[5] = threshImg[((row+2)*width) +(col+1)];         if(neibs[5] == thresh) { cnt++; }   // curr offset by 2 - bottom right

            if(cnt == 6) { // if all neighbours are bright, we are at the edge boundary
                firstSurfaceImg[(row)*width + col] = 1.0f;
                row = height;
            }
        }
    }

    else if ( col == (width-1) ) { // last col 
        for (int row = 0 + rows2skip; row < height; row++) { 
            int cnt = 0;
             float neibs[6]; // not shared mem as it reduces speed  

            // get six neighbours - three in same col, and three to the left
            neibs[0] = threshImg[((row)*width) +(col)];             if(neibs[0] == thresh) { cnt++; }   // current position
            neibs[1] = threshImg[((row)*width) +(col-1)];           if(neibs[1] == thresh) { cnt++; }   // left
            neibs[2] = threshImg[((row+1)*width) +(col)];           if(neibs[2] == thresh) { cnt++; }   // bottom
            neibs[3] = threshImg[((row+1)*width) +(col-1)];         if(neibs[3] == thresh) { cnt++; }   // bottom left
            neibs[4] = threshImg[((row+2)*width) +(col)];           if(neibs[4] == thresh) { cnt++; }   // curr offset by 2 - bottom
            neibs[5] = threshImg[((row+2)*width) +(col-1)];         if(neibs[5] == thresh) { cnt++; }   // curr offset by 2 - bottom left

            if(cnt == 6) { // if all neighbours are bright, we are at the edge boundary
                firstSurfaceImg[(row)*width + col] = 1.0f;
                row = height;
            }
        }       
    }

    // remaining threads are: (1 -> 510) 

    else { // any col other than first or last column
        for (int row = 0 + rows2skip; row < height; row++) { 

            int cnt = 0;
            float neibs[9]; // not shared mem as it reduces speed   

            // for threads < width/4, get the neighbors
            // get nine neighbours - three in curr col, three each to left and right
            neibs[0] = threshImg[((row)*width) +(col-1)];           if(neibs[0] == thresh) { cnt++; } 
            neibs[1] = threshImg[((row)*width) +(col)];             if(neibs[1] == thresh) { cnt++; } 
            neibs[2] = threshImg[((row)*width) +(col+1)];           if(neibs[2] == thresh) { cnt++; }           
            neibs[3] = threshImg[((row+1)*width) +(col-1)];         if(neibs[3] == thresh) { cnt++; }           
            neibs[4] = threshImg[((row+1)*width) +(col)];           if(neibs[4] == thresh) { cnt++; }           
            neibs[5] = threshImg[((row+1)*width) +(col+1)];         if(neibs[5] == thresh) { cnt++; }           
            neibs[6] = threshImg[((row+2)*width) +(col-1)];         if(neibs[6] == thresh) { cnt++; }           
            neibs[7] = threshImg[((row+2)*width) +(col)];           if(neibs[7] == thresh) { cnt++; }           
            neibs[8] = threshImg[((row+2)*width) +(col+1)];         if(neibs[8] == thresh) { cnt++; }

            if(cnt == 9) { // if all neighbours are bright, we are at the edge boundary

                firstSurfaceImg[(row)*width + col] = 1.0f;
                row = height;
                }
            }
        }       
    }           

__syncthreads();
}

__global__ void firstSurfaceDetection_opt (const float * __restrict__ threshImg, int *firstSurfaceImgRow, int height, int width) {

  int col = threadIdx.x + (blockDim.x*blockIdx.x); 
  int row = threadIdx.y + (blockDim.y*blockIdx.y);

  int rows2skip = SKIP_ROWS; 
  float thresh = THRESH;

  if ((row >= rows2skip) && (row < height-2) && (col < width) && (row < firstSurfaceImgRow[col])) {

    int cnt = 0;
    int inc = 0;
    if (col == 0) inc = +1;
    if (col == (width-1)) inc = -1;
    if (inc){
            cnt = 3;
            if (threshImg[((row)*width)   +(col)]     == thresh) cnt++;
            if (threshImg[((row)*width)   +(col+inc)] == thresh) cnt++;
            if (threshImg[((row+1)*width) +(col)]     == thresh) cnt++;   
            if (threshImg[((row+1)*width) +(col+inc)] == thresh) cnt++;      
            if (threshImg[((row+2)*width) +(col)]     == thresh) cnt++;     
            if (threshImg[((row+2)*width) +(col+inc)] == thresh) cnt++;
            }
    else {
            // get nine neighbours - three in curr col, three each to left and right
            if (threshImg[((row)*width)   +(col-1)] == thresh) cnt++;
            if (threshImg[((row)*width)   +(col)]   == thresh) cnt++;
            if (threshImg[((row)*width)   +(col+1)] == thresh) cnt++;
            if (threshImg[((row+1)*width) +(col-1)] == thresh) cnt++;
            if (threshImg[((row+1)*width) +(col)]   == thresh) cnt++;   
            if (threshImg[((row+1)*width) +(col+1)] == thresh) cnt++;      
            if (threshImg[((row+2)*width) +(col-1)] == thresh) cnt++;
            if (threshImg[((row+2)*width) +(col)]   == thresh) cnt++;     
            if (threshImg[((row+2)*width) +(col+1)] == thresh) cnt++;
            }
    if(cnt == 9) { // if all neighbours are bright, we are at the edge boundary
            atomicMin(firstSurfaceImgRow + col, row);
            }
    }
}


int main(int argc, char *argv[]){

  float *threshImg, *h_threshImg, *firstSurfaceImg, *h_firstSurfaceImg;
  int *firstSurfaceImgRow, *h_firstSurfaceImgRow;
  int actualImHeight = 1024;
  int actualImWidth = 512;
  int row_set = 512;
  if (argc > 1){
    int my_val = atoi(argv[1]);
    if ((my_val > SKIP_ROWS) && (my_val < actualImHeight - 3)) row_set = my_val;
    }

  h_firstSurfaceImg = (float *)malloc(actualImHeight*actualImWidth*sizeof(float));
  h_threshImg = (float *)malloc(actualImHeight*actualImWidth*sizeof(float));
  h_firstSurfaceImgRow = (int *)malloc(actualImWidth*sizeof(int));
  cudaMalloc(&threshImg, actualImHeight*actualImWidth*sizeof(float));
  cudaMalloc(&firstSurfaceImg, actualImHeight*actualImWidth*sizeof(float));
  cudaMalloc(&firstSurfaceImgRow, actualImWidth*sizeof(int));
  cudaMemset(firstSurfaceImgRow, 1, actualImWidth*sizeof(int));
  cudaMemset(firstSurfaceImg, 0, actualImHeight*actualImWidth*sizeof(float));

  for (int i = 0; i < actualImHeight*actualImWidth; i++) h_threshImg[i] = 0.0f;
  // insert "bright row" here
  for (int i = (row_set*actualImWidth); i < ((row_set+3)*actualImWidth); i++) h_threshImg[i] = THRESH;

  cudaMemcpy(threshImg, h_threshImg, actualImHeight*actualImWidth*sizeof(float), cudaMemcpyHostToDevice);


  dim3 grid(1,1024);
  //warm-up run
  firstSurfaceDetection_opt <<< grid, 512 >>> (threshImg, firstSurfaceImgRow, actualImHeight, actualImWidth);
  cudaDeviceSynchronize();
  cudaMemset(firstSurfaceImgRow, 1, actualImWidth*sizeof(int));
  cudaDeviceSynchronize();
  unsigned long long t2 = dtime_usec(0);
  firstSurfaceDetection_opt <<< grid, 512 >>> (threshImg, firstSurfaceImgRow, actualImHeight, actualImWidth);
  cudaDeviceSynchronize();
  t2 = dtime_usec(t2);
  cudaMemcpy(h_firstSurfaceImgRow, firstSurfaceImgRow, actualImWidth*sizeof(float), cudaMemcpyDeviceToHost);
  unsigned long long t1 = dtime_usec(0);
  firstSurfaceDetection <<< 1, 512 >>> (threshImg, firstSurfaceImg, actualImHeight, actualImWidth);
  cudaDeviceSynchronize();
  t1 = dtime_usec(t1);
  cudaMemcpy(h_firstSurfaceImg, firstSurfaceImg, actualImWidth*actualImHeight*sizeof(float), cudaMemcpyDeviceToHost); 

  printf("t1 = %fs, t2 = %fs\n", t1/(float)USECPSEC, t2/(float)USECPSEC);
  // validate results
  for (int i = 0; i < actualImWidth; i++) 
    if (h_firstSurfaceImgRow[i] < actualImHeight) 
      if (h_firstSurfaceImg[(h_firstSurfaceImgRow[i]*actualImWidth)+i] != THRESH)
        {printf("mismatch at %d, was %f, should be %d\n", i, h_firstSurfaceImg[(h_firstSurfaceImgRow[i]*actualImWidth)+i], THRESH); return 1;}
  return 0;
}
$ nvcc -arch=sm_35 -o t667 t667.cu
$ ./t667
t1 = 0.000978s, t2 = 0.000050s
$

    笔记:

    上面的示例
  • 在行= 512的整个图像上插入了一个“明亮的邻域”，因此在我的情况下(K40c)给出的道路中间加速因子几乎是20倍。
  • 为了简洁起见，我省略了proper cuda error checking。我推荐它。
  • 每个内核的执行性能在很大程度上取决于它是否首次运行。这可能与缓存和一般的预热效果有关。因此，为了给出合理的结果，我首先运行了内核，以进行额外的不定时预热。
  • 我没有进行共享内存优化的原因之一是这个问题已经很小了，至少对于像K40这样的大型开普勒GPU，它几乎完全适合L2缓存(尤其是我的内核，因为它使用较小的输出数据结构。)鉴于此，共享内存可能无法提供很多性能提升。

    • 编辑:我再次修改了代码，以便可以将插入明亮边界的测试图像中的行(行)作为命令行参数传递，并且我已经在3种不同的设备上测试了该代码，对亮行使用3种不同的设置:
    execution time on:     K40    C2075    Quadro NVS 310
bright row =   15:   31/33    23/45       29/314
bright row =  512:  978/50  929/112     1232/805
bright row = 1000: 2031/71 2033/149    2418/1285

all times are microseconds (original kernel/optimized kernel)
CUDA 6.5, CentOS 6.2

    关于c++ - CUDA-优化表面检测内核，我们在Stack Overflow上找到一个类似的问题： https://stackoverflow.com/questions/28509499/

    上一篇：C++ - 具有任意数量参数的函数指针

    下一篇：c++ - 调试断言失败错误

    相关文章：

    在 CUDA 中对许多小数组进行排序

    c++ - 不同应用的CUDA设备堆内存

    c++ - undefined reference CUDA 和 CMAKE

    c++ - 用于创建带有 cuda 函数的 C++ 共享库的 Makefile

    c++ - 将 protobuf 与 boost asio 结合使用

    c++ - 如何在 openGL 中乘以矩阵

    c++ - 在代码中将 C++ 应用程序设置为使用最大 CPU 使用率

    c++ - 在 Factory : Pros and Cons? 中使用静态方法

    c++ - 将数组传递给 Cuda

    具有可变键数的 C++ 键值容器

    ©2024 IT工具网  联系我们