0. 简介

最近作者希望系统性的去学习一下CUDA加速的相关知识，正好看到深蓝学院有这一门课程。所以这里作者以此课程来作为主线来进行记录分享，方便能给CUDA网络加速学习的萌新们去提供一定的帮助。

1. 基础矩阵乘法

下图是矩阵乘法的示意图，下面我们来看一下在CPU和GPU上是怎么表达的。
在这里插入图片描述
CPU代码示意流程：

// Matrix multiplication on the (CPU) host
void main(){
    define A, B, C
    for i= 0 to M-1 do
        for j = 0 to N-1 do
        /* compute element C(i,j) */
            for k = 0 to K-1 do
                C(i,j) <= C(i,j) + A(i,k) * B(k,j)
            end
        end
    end
}

GPU代码示意流程：

void main(){
    define A_cpu, B_cpu, C_cpuin the CPU memory
    define A_gpu, B_gpu, C_gpuin the GPU memory
    memcopyA_cputo A_gpu
    memcopyB_cputo B_gpu
    dim3 dimBlock(16, 16)
    dim3 dimGrid(N/dimBlock.x, M/dimBlock.y)
    matrixMul<<<dimGrid, dimBlock>>>(A_gpu,B_gpu,C_gpu,K)
    memcopyC_gputo C_cpu
}

__global__ void matrixMul(A_gpu,B_gpu,C_gpu,K){
    temp <= 0
    i<= blockIdx.y* blockDim.y+ threadIdx.y// Row iof matrix C
    j <= blockIdx.x* blockDim.x+ threadIdx.x// Column j of matrix C
    for k = 0 to K-1 do
    accu<= accu+ A_gpu(i,k) * B_gpu(k,j)
    end
    C_gpu(i,j) <= accu
}

在GPU中，一个线程负责计算C中的一个元素，其中A中的每一行从全局内存中载入N次，B中的每一列从全局内存中载入M次。总共的次数为2mnk的读取次数（因为在每一个k处都要从A和B中读取一次，所以要乘二），值得注意的是，可以将多次访问的数据放到共享内存中，减少重复读取的次数，并充分利用共享内存的延迟低的优势
在这里插入图片描述
所以我们可以看到共享内存可以被然被内存访问指令访问，同时拥有更高的速度访问(延迟&吞吐)。在共享内存中存在有两种申请空间方式，静态申请和动态申请，但是共享内存的大小只有几十K，过度使用共享内存会降低程序的并行性，共享内存在使用时候需要注意的是：

  • 使用shared关键字，同时有静态申请和动态申请两种方式；
  • 将每个线程从全局索引位置读取元素，将它存储到共享内存之中；
  • 注意数据存在着交叉，应该将边界上的数据拷贝进来；
  • 块内线程同步：__syncthreads()

2. 线程同步函数

上面讲到了syncthreads()函数，该函数是cuda的内建函数，用于块内线程通信。syncthreads()函数是对线程进行同步，需要保证的是需要对所有的共享内存需要同步的线程都被同步到，下面是两个示例。
在这里插入图片描述
申请共享内存会存在有两种申请方式：
1.静态申请

__global__ void staticReverse(int *d, int n) {
    __shared__ int s[64];
    int t = threadIdx.x;
    int tr = n-t-1;
    s[t] = d[t];
    __syncthreads();
    d[t] = s[tr];
}
staticReverse<<<1,n>>>(d_d, n);

2.动态申请

相较于静态申请来看在s后面是没有明确大小的，其次加了extern来强调了动态申请，最后在最底下除了限定维度外，还需要设置sizeof(int)来设置空间。

__global__ void dynamicReverse(int *d, int n) {
    extern __shared__ int s[];
    int t = threadIdx.x;
    int tr = n-t-1;
    s[t] = d[t];
    __syncthreads();
    d[t] = s[tr];
}
dynamicReverse<<<1,n,n*sizeof(int)>>>(d_d, n);

3. 进阶版矩阵乘法

之前第一章提到的每个输入元素被WIDTH个线程读取，这就导致被重复多次加载。所以可以将每个元素加载到共享内存中并让多个线程使用本地版本以减少内存带宽；从而让每次计算都会计算成一个矩阵块。从而解决了基础版矩阵乘法的两个点：

  1. 基础版矩阵乘法里面的参数被频繁读取，所以放到共享内存中解决这样的问题
  2. 基础版矩阵乘法里面的C每一次计算只能计算出一个，所以通过多线程让一次计算计算出多个C的元素。，使每个阶段的数据访问集中在一个子集上(tile) of Md and Nd
    在这里插入图片描述
    上图就是使用block完成多线程一次计算出多个C元素的操作，具体大小为block_size第一个平铺矩阵元素M[Row][tx]、N[ty][Col]；下一个平铺矩阵元素就是M[Row][1*TILE_WIDTH+tx]、N[1*TILE_WIDTH+ty][Col] 。具体代码如下：
/* Accumulate C tile by tile. */
for tileIdx= 0 to (K/blockDim.x-1) do
    /* Load one tile of A and one tile of B into shared mem */
    i<= blockIdx.y* blockDim.y+ threadIdx.y// Row iof matrix A
    j <= tileIdx* blockDim.x+ threadIdx.x// Column j of matrix A
    A_tile(threadIdx.y, threadIdx.x) <= A_gpu(i,j) // Load A(i,j) to shared mem
    B_tile(threadIdx.x, threadIdx.y) <= B_gpu(j,i) // Load B(j,i) to shared mem
    __sync() // Synchronize before computation
    /* Accumulate one tile of C from tiles of A and B in shared mem */
    for k = 0 to threadDim.xdo
        accu<= accu+ A_tile(threadIdx.y,k) * B_tile(k,threadIdx.x)
    end
    __sync()
end

我们可以看到原始矩阵乘法需要取2mnk次，但是平铺矩阵算法则可以简化为2*Block_size*Block_size*[k/Block_size]*[m/Block_size]*[n/Block_size]=2mnk/Block_size。

本系列接下来的内容，将会从CUDA的Stream和Event两个部分来完成对CUDA的进一步深入。