这次回顾ECE408 Lecture 3，这次的介绍了CUDA并行执行模型。

课程主页：

搬运视频：

本讲内容

  • 了解更多关于CUDA线程的多维逻辑组织；
  • 学习使用控制结构，例如kernel中的循环；
  • 学习线程调度、延迟容限和硬件占用的概念；

内容回顾

这里回顾向量加法。

基本情形

矩阵加法，假设$n=1000$，block size为256：

vecAdd<<<ceil(N/256.0), 256>>>(...)
    
i = blockIdx.x * blockDim.x + threadIdx.x;
if (i<N) C[i] = A[i] + B[i];

图示：

红色部分解释了为什么要增加i<N的判断：为了在数组大小不整除block size时结果依然正确。

该代码对应一个线程处理1个元素，那么处理两个元素时该如何处理呢。

一个线程处理两个元素

代码：

vecAdd<<<ceil(N/(2*256.0)), 256>>>(...)
    
i = blockIdx.x * (2*blockDim.x) + threadIdx.x;
if (i<N) C[i] = A[i] + B[i];
i = i + blockDim.x
if (i<N) C[i] = A[i] + B[i];

说明：

  • 因为一个线程处理两个元素，所以block数量为N/(2*256.0)；
  • 后续代码可以理解为处理$2i, 2i+1$位置的元素；

Kernel, Grid, Block回顾

案例1：将彩色图像转换为灰度图像

案例1是将彩色图转换为灰度图：

处理逻辑

假设我们用2D Grid处理图片：

那么图像可以用许多个block覆盖：

因为图片长宽不一定整除block size，所以多余的部分需要做判断：

C/C++中二维数组的行优先布局

在C/C++中，二维数组是按行优先的顺序，以一维数组的形式存储在内存中，对应的二维坐标可以转换为一维坐标：

代码

// we have 3 channels corresponding to RGB
// The input image is encoded as unsigned characters [0, 255]
__global__
void RGBToGrayscale(unsigned char * grayImage, unsigned char * rgbImage, int width, int height)
{
    int Col = threadIdx.x + blockIdx.x * blockDim.x;
    int Row = threadIdx.y + blockIdx.y * blockDim.y;
    if (Col < width && Row < height) {
        // get 1D coordinate for the grayscale image
        int grayOffset = Row*width + Col;
        // one can think of the RGB image having
        // CHANNEL times columns of the gray scale image
        int rgbOffset = grayOffset*CHANNELS;
        unsigned char r = rgbImage[rgbOffset    ]; // red value for pixel
        unsigned char g = rgbImage[rgbOffset + 1]; // green value for pixel
        unsigned char b = rgbImage[rgbOffset + 2]; // blue value for pixel
        // perform the rescaling and store it
        // We multiply by floating point constants
        grayImage[grayOffset] = 0.21f*r + 0.71f*g + 0.07f*b;
    }
}

案例2：图像模糊

图像模糊将图片转换为另一幅图片：

每个输出像素是（原图）周围像素的平均值：

注意一些边界情况需要特殊处理：

代码

__global__
void blurKernel(unsigned char * in, unsigned char * out, int w, int h) {
    int Col = blockIdx.x * blockDim.x + threadIdx.x;
    int Row = blockIdx.y * blockDim.y + threadIdx.y;
    if (Col < w && Row < h) {
        int pixVal = 0;
        int pixels = 0;
        // Get the average of the surrounding BLUR_SIZE x BLUR_SIZE box
        for(int blurRow = -BLUR_SIZE; blurRow <= BLUR_SIZE; ++blurRow) {
            for(int blurCol = -BLUR_SIZE; blurCol <= BLUR_SIZE; ++blurCol) {
                int curRow = Row + blurRow;
                int curCol = Col + blurCol;
                // Verify we have a valid image pixel
            	if(curRow > -1 && curRow < h && curCol > -1 && curCol < w) {
                    pixVal += in[curRow * w + curCol];
                    pixels++; // Keep track of number of pixels in the avg
                }
            }
        }
        // Write our new pixel value out
        out[Row * w + Col] = (unsigned char)(pixVal / pixels);
    }
}

CUDA执行模型：线程块

  • 块中的所有线程执行相同的内核程序(SPMD)；
  • 同一块中的线程共享数据并在完成工作时进行同步；
  • 不同块中的线程不能协作；
  • 块以任意顺序执行；
  • 同一块内的线程按warp顺序执行；

执行线程块

  • 线程以块粒度分配给流式多处理器（SM）
    • 每个SM最多32个块（Maxwell 的资源限制）；
    • Maxwell/Pascal/Turing SM最多可占用2048个线程；
  • 线程以wrap的形式并发运行
    • SM维护线程/块ID；
    • SM管理/调度线程执行；

图示：

线程调度

  • 每个块作为32线程线程束（warp）执行
    • 是实现决定的，不是CUDA编程模型的一部分
    • warp根据其线性化线程索引进行划分
      • 线程0-31：warp 0；
      • 线程32-63：warp 1；
      • 按照X, Y, Z维度的顺序；
    • warps是SM中的调度单位
  • 如果将3个block分配给一个SM，每个块有256个线程，那么一个SM中有多少个 warp？
    • 每个块被分成256/32 = 8 warps；
    • 8 warps/block * 3 block = 24 warps；
  • SM实现零开销warp调度
    • 其下一条指令的操作数在执行前已经准备好；
    • 选择符合条件的warp以根据优先调度策略执行；
    • warp中的所有线程在被选中时执行相同的指令；

分支（Branch Divergence）

  • 分支分歧
    • warp中的线程在程序中采用不同路径；
    • 控制流的主要性能问题；
  • GPU使用预测执行
    • 每个线程计算是否走某个路径；
    • warp中线程采用的多条路径是串行执行的；

Branch Divergence示例

使用thread ID作为分支条件：

if (threadIdx.x > 2) {
	// THEN path (lots of lines)
} else {
	// ELSE path (lots more lines)
}
  • 所有在warp中的进程，有两条控制路径（THEN/ELSE）；
  • 所有线程都执行两条路径（仅当线程的谓词为真时才保留结果）；
  • 这样效率很低。

尝试使分支粒度成为warp大小的倍数（记住，它可能并不总是32！）：

if ((threadIdx.x / WARP_SIZE) > 2) {
	// THEN path (lots of lines)
} else {
	// ELSE path (lots of lines)
}
  • 仍然有两个控制路径；
  • 但是任何warp中的所有线程都只遵循一条路径；
  • 效率更高；
    • 以warp为单位划分分支；

块粒度注意事项

  • 对于RGBToGrayscaleConversion，我们应该使用8×8、16×16还是32×32块？假设GPU可以有1,536个线程和每个 SM最多8个块。
    • 对于8×8，我们每个块有64个线程。每个SM最多可以占用1,536个线程，也就是1,536/64=24个块。但是每个 SM最多只能占用8个块！
    • 对于16×16，我们每个块有256个线程。每个SM最多可以占用1,536个线程（1536/32=48个warp），即6个块（在8个块限制内）。因此，我们使用SM的全部线程容量。
    • 对于32×32，每个块有1,024个线程。一个SM只能容纳一个块，仅使用SM线程容量的2/3。