前言

自动驾驶之心推出的 《CUDA与TensorRT部署实战课程》，链接。记录下个人学习笔记，仅供自己参考

Note：关于 CUDA 中错误的处理杜老师之前也讲过，感兴趣的可以看看 3.9.cuda运行时API-错误处理的理解以及错误的传播特性

本次课程我们来学习课程第二章————CUDA 编程入门，一起来学习使用 CUDA 进行矩阵乘法的加速

课程大纲可以看下面的思维导图

0. 简述

这节课程我们来讲第二章第 2 小节，使用 CUDA 进行矩阵乘法的加速，这个小节我们主要分为以下三个部分讲解：

• 初步计算 MatMul
• CUDA 中的 error handler
• 获取 GPU 的硬件信息

1. 初步计算 MatMul

本小节目标：理解使用 cuda 进行矩阵乘法的加速方法，tile 的用意

这个部分主要讲解初步计算 Matmul 矩阵乘法，这个部分学完之后希望大家能够去理解 CUDA 编程中矩阵乘法的一些基本操作，加速方法，还有矩阵乘法中的 tile 是什么如何去使用，GPU 的并行运算里面它是如何体现 tile 的。同时会教大家一些基本的 CUDA Runtime API 的使用，一些参数的配置等待

1.1 执行一下我们的第三个CUDA程序

源代码获取地址：https://github.com/kalfazed/tensorrt_starter

我们先来看这节课的案例代码 2.3-matmul-basic，如下图所示：

2.3-matmul-basic

src 源文件里面东西比较多，首先这边实现了一个 timer 计时器（timer.hpp）以及 CPU 和 GPU 上矩阵乘法的实现，main 函数里面有对两种不同方式的矩阵乘法进行性能测试，具体的代码分析我们下面会给大家进行讲解，这里大家简单了解下就行。

Makefile

Makefile 文件我们也有做简单修改，因为这里面需要编译的源文件越来越多了，所以我们希望每次在编译的过程中去建立依赖关系，这样我们就可以每次在修改一些程序的时候，我们只要去编译修改的那个程序，没有修改的部分我们就不要让它去编译了，这个其实可以通过你的 g++ 和 nvcc 里面的 -M、-MF 这些选项可以建立起依赖关系，这个写法主要参考了 tensorRT_Pro，大家感兴趣的可以看下

utils.cpp

这里还给大家实现了一个矩阵的初始化，矩阵比较的一些实现

timer.hpp

刚才说的 timer 计时器的也有实现，之后我们进入代码分析再去讲解

1.2 host端与device端的数据传输

这个案例我们需要做的是矩阵乘法，所以也就意味着我们需要在 host 端和 device 端进行一些内存空间的分配和内存的一些拷贝，这样我们在 CPU(host) 端其实需要做的东西还是比较多，主要包括以下几个部分：

• 分配 host 与 device 端的内存空间
• 将数据传送到 GPU
• 配置核函数的参数
  • grid dim
  • block dim
  • shared memory size
  • stream
• 启动核函数
• 将数据从 GPU 传入回来

首先你需要分配你的 host 端和 device 端一个内存空间，这里 CUDA 提供了一些 API 可以去使用：

• cudaMalloc：在 device 端进行空间分配，是一种 cuda runtime api
• cudaMallocHost：在 host 端的 pinned memory 上分配空间

这两个 API 都属于 CUDA Runtime API，大家如果第一次听到 CUDA Runtime API 可能不是很清楚，一般来说大家在 CUDA 编程中会使用到一些 API 函数，如果是以 cudaxxx 开头的一般就是指 CUDA Runtime API。

另外一个大家可能接触的并不是很多，另外一种其实是以 CUxxx 开头的，这种 API 我们会叫 CUDA Driver API，这两种 API 其实还是有一点不同的，它们的层级关系可以看下图：

上图中你 CPU 有一个 Application 应用程序，你从你的程序端其实可以调用的东西大概有三种，一个是你的 cuda library，这个的话就是比如说 cuDNN 还有一些其它的 library 库我们都可以从这里面去访问，另外一个就是大家用得比较多的 CUDA Runtime API，比如 cudaMalloc、cudaMallocHost、cudaMemcpy、cudaMemcpyAsync 这些都是可以在这里面使用的

另外一个就是偏底层的 CUDA Driver API，这个是驱动级的 API，是和硬件沟通的偏底层的东西，我们写程序时一般不会用到，因为它内部的结构比较复杂，并且很难去进行调试，因为它是直接控制你的 GPU 端的，而 Runtime API 其实就是在 Driver API 上面进行的一个比较好的封装，这样方便我们对底层的操作，包括 implicit initialization（隐式初始化）、context management（上下文管理）、module management（模块管理）都是在 Runtime API 里面做好了的，大家记住一下这个区别就好了。

OK，回到正题，我们把数据分配完空间之后就可以进行数据传输，将数据从 CPU 传输到 GPU，数据传输一般有两个方式，一个是 cudaMemcpy，一个是 cudaMemcpyAsync，两个 API 的不同点在于是同步还是异步方式进行数据拷贝传输，数据拷贝的方向可以分为把数据从 Host 端拷贝到 Device 端，或者是 Device 端之间进行传输，亦或是将 Device 端数据拷贝到 Host 端，通过参数指定即可

传输完数据之后，我们就可以开始配置核函数的参数了，包括 (Grid, Block, shared_memory, stream) 这四个参数其实是需要配置的，但一般来说前面的 Grid, Block 大小是必须要配置的，后面的 shared_memory, stream 这两个参数并不是一定要配置的，因为它们都有一些默认值，具体这两个参数怎么用呢？我们之后的课程中再去讲解，这里大家不用去考虑

配置完参数之后我们就需要从你的 Host 端去启动 Device 端的核函数，记住这个核函数启动的时候一般都是异步的，所以我们启动完核函数之后就要进行一个同步，核函数执行完成之后我们需要将结果从你的 Deivce 端传入回你的 Host 端，这个也是可以用 cudaMemcpy 这个 API 去实现的

OK，这是 Host 端需要做的事情，其次就是你的 Device 端，Deivce 端的话它其实做的东西比较简单，就是根据你 Host 端配置的这些参数来执行一个核函数，利用多个 thread 来并行计算

以上就是一个 CUDA 程序中数据传输的基本流程，大家了解下即可

1.3 CUDA Core的矩阵乘法计算

我们下面来给大家简单介绍一下你的 CUDA Core 它的矩阵乘法计算是怎么做的，这里给大家举个 A * B = C 例子，如下图所示：

A 是 4x8 的矩阵，B 是 8x4 的矩阵，C 是 4x4 的矩阵，我们以 block size = 1 为例来讲解，首先我们先要计算红色小框框，它涉及到 A 的第一行和 B 的第一列，把它们用公式展开的话就是 c(0,0) = a(0,0) * b(0,0) + … + a(0,7)*b(7.0)

这里面我们需要八个乘法和八个加法，这个在 GPU 端我们会给它叫做 FMA（Fused Multiply-Add），就是一个乘加法指令。这里稍微扩展一下，一般来说我们讲一个 clk 在做一个浮点运算，指的是一个乘法的浮点运算或者加法的浮点运算，但是如果说我们的指令是 FMA，也就是你的乘法和加法混合在一起进行计算的话，那么我们可以说一个 clk 计算两个浮点数的计算，这个的话是可以提高吞吐量，所以说这个例子如果按照 FMA 来计算，如果想要计算 c(0,0) 我们需要 8 个 clk 来完成。

如果说我们用一个 thread 要完成 4x8 与 8x4 的矩阵乘法计算，需要 8 * 16 = 128个 clk 才可以完成，这个效率就相当低下，那么我们怎么办呢？

这个其实我们在并行处理中讲过，我们 c(0,0) 的计算和 c(0,1)、c(0,2) 等等这些计算其实是没有任何依赖关系的，因此我们可以把这些东西一起计算，我们刚才不是说分配一个 thread 来处理吗，那么我们现在可以直接分配 16 个 thread，每个 thread 负责 C 其中一个元素的计算，那么 16 个 thread 一起计算它所需要的时间其实跟你计算一个小元素所需要的时间是一样的，8 个 clk 就完成了，因为我们是并行处理的，这个就是利用 CUDA Core 完成的计算

值得注意的是 CUDA Core 和 Tensor Core 进行矩阵乘法的方法还是有一点不同的，Tensor Core 会比 CUDA Core 更快一点，这个我们在之后的课程中再跟大家去讲解

下面我们再扩展一下，现在有一个 8x16 的 A 矩阵和一个 16x8 的 B 矩阵相乘得到一个 8x8 的 C 矩阵，如下图所示：

现在需要 CUDA Core 完成这个矩阵的计算，那现在我们就可以有很多策略了，比如说我们可以把 block 分成 4x4，那么相应的我们的 grid 就是 2x2，也就是 blockDim(4,4)，gridDim(2,2)，计算过程如下图所示：

红色小块就是一个 block 的计算，是 16 个 thread 一起进行的，它所涉及到的数据就是矩阵 A 中绿色的四行和矩阵 B 中紫色的四列，由于我们设置的 grid 是 2x2，意味着我们还有其它三个 block 也是做类似的计算，四个 block 共同完成这个矩阵乘法的计算。

当然我们对 blockDim 和 gridDim 也可以有不同的划分，如下图所示：

在上图中我们的 blockDim 被设置为 8x8，gridDim 被设置为 1x1，这意味着我们只通过一个 block 中的 84 个 thread 就可以完成矩阵 C 的计算，每个 thread 负责 C 中的一个元素即 A 中的一行数据和 B 中的一列数据的乘加运算

大家可能会想既然可以给一个 block 分配 64 个 thread 完成计算，那是不是说我们也可以分配 128、256、512 甚至 1024 个 thread 完成相应的计算呢？这个是完全可以的，但是有一个点需要大家注意，就是一个 block 可分配的 thread 数量是有限的

CUDA 编程中有个规定，就是一个 block 中可以分配的 thread 的数量最大是 1024 个线程，如果大于 1024 会显示配置错误。

1.4 代码分析

我们讲了这么多，我们开始进入到代码里面看下它具体是怎么实现的，我们先 make 执行下 2.3 小节的案例，输出结果如下：

可以看到有三个结果，一个是 CPU 端的矩阵乘法耗时，大概是 2186.26 ms，一个是 GPU warmup 的耗时是 290.353 ms，还有 GPU 端的矩阵乘法耗时 5.31113 ms，最后是二者的精度比较，如果误差允许范围是 10-4 那么二者精度是一样的

我们先来看计时器代码 timer.hpp，代码如下：

#include <chrono>
#include <cstdio>
#include <ratio>
#include <string>
#include <iostream>


class Timer {
public:
    using s  = std::ratio<1, 1>;
    using ms = std::ratio<1, 1000>;
    using us = std::ratio<1, 1000000>;
    using ns = std::ratio<1, 1000000000>;

public:
    Timer(){};

public:
    void start() {mStart = std::chrono::high_resolution_clock::now();}
    void stop()  {mStop  = std::chrono::high_resolution_clock::now();}

    template <typename span>
    void duration(std::string msg);

private:
    std::chrono::time_point<std::chrono::high_resolution_clock> mStart;
    std::chrono::time_point<std::chrono::high_resolution_clock> mStop;
};

/* 
 * 注意：这个实现是不能够非常精准的获取kernel函数的执行时间 
 * 要如果想要精准的获取kernel实现需要通过cuda event来进行测量，这个在后面的案例中会讲
 */
template <typename span>
void Timer::duration(std::string msg){
    std::string str;
    char fMsg[100];
    std::sprintf(fMsg, "%-30s", msg.c_str());

    if(std::is_same<span, s>::value) { str = " s"; }
    else if(std::is_same<span, ms>::value) { str = " ms"; }
    else if(std::is_same<span, us>::value) { str = " us"; }
    else if(std::is_same<span, ns>::value) { str = " ns"; }

    std::chrono::duration<double, span> time = mStop - mStart;
    std::cout << fMsg << " uses " << time.count() << str << std::endl;
}

timer.hpp 是一个用于测量代码执行时间的 C++ 头文件，它基于 C++11 的 <chrono> 库构建。该头文件定义了一个 Timer 类，允许用户通过简单的接口来记录和计算代码执行的开始和结束时间点。Timer 类提供了 start() 和 stop() 成员函数用于分别记录时间点，以及一个模板成员函数 duration() 用于输出执行时间。该类支持秒、毫秒、微秒和纳秒四种时间单位，通过模板参数指定。(from chatGPT)

我们再来看工具代码 utils.cpp，代码如下：

#include "utils.hpp"
#include <math.h>
#include <random>


void initMatrix(float* data, int size, int min, int max, int seed) {
    srand(seed);
    for (int i = 0; i < size; i ++) {
        data[i] = float(rand()) * float(max - min) / RAND_MAX;
    }
}

void printMat(float* data, int size) {
    for (int i = 0; i < size; i ++) {
        printf("%.8lf", data[i]);
        if (i != size - 1) {
            printf(", ");
        } else {
            printf("\n");
        }
    }
}

void compareMat(float* h_data, float* d_data, int size) {
    double precision = 1.0E-4;
    bool error = false;
    /* 
     * 这里注意，浮点数运算时CPU和GPU之间的计算结果是有误差的
     * 一般来说误差保持在1.0E-4之内是可以接受的
    */
    for (int i = 0; i < size; i ++) {
        if (abs(h_data[i] - d_data[i]) > precision) {
            error = true;
            printf("res is different in %d, cpu: %.8lf, gpu: %.8lf\n",i, h_data[i], d_data[i]);
            break;
        }
    }
    if (error) 
      printf("Matmul result is different\n");
    else
      printf("Matmul result is same, precision is 1.0E-4\n");
}

utils.cpp 是 CUDA 矩阵乘法示例程序中的一个工具文件，包含了初始化矩阵、打印矩阵以及比较矩阵的功能。该文件中定义了三个函数：(from chatGPT)

1. initMatrix：用于初始化一个浮点数矩阵。它接收一个浮点数组指针、矩阵大小、值的最小和最大范围以及一个种子用于随机数生成，确保矩阵中的值随机分布在给定的范围内。
2. printMat：打印矩阵的函数。给定一个浮点数数组和大小，它将遍历数组并打印每个元素，元素之间用逗号分隔，以方便查看矩阵内容。
3. compareMat：比较两个矩阵的函数，用于验证 CUDA 计算结果的正确性。它接收两个浮点数数组（一个从 CPU 计算得到，另一个从 GPU 计算得到）和数组大小，然后逐元素比较两个矩阵。如果差异超出预设的精度阈值（1.0E-4），则报告不一致，并指出首个不一致元素的位置和值。

下面我们来看下 CPU 端矩阵乘法的实现代码 matmul_cpu.cpp，如下所示：

#include "matmul.hpp"

void MatmulOnHost(float *M, float *N, float *P, int width){
    for (int i = 0; i < width; i ++)
        for (int j = 0; j < width; j ++){
            float sum = 0;
            for (int k = 0; k < width; k++){
                float a = M[i * width + k];
                float b = N[k * width + j];
                sum += a * b;
            }
            P[i * width + j] = sum;
        }
}

这个就是一个正常的 A * B = C 的矩阵乘法实现，三次 for 循环串行执行的 CPU 实现

我们重点来看下 GPU 端矩阵乘法的实现代码 matmul_gpu_basic.cpp，如下所示：

#include "cuda_runtime.h"
#include "cuda.h"
#include "stdio.h"

/* matmul的函数实现*/
__global__ void MatmulKernel(float *M_device, float *N_device, float *P_device, int width){
    /* 
        我们设定每一个thread负责P中的一个坐标的matmul
        所以一共有width * width个thread并行处理P的计算
    */
    int y = blockIdx.y * blockDim.y + threadIdx.y;
    int x = blockIdx.x * blockDim.x + threadIdx.x;

    float P_element = 0;

    /* 对于每一个P的元素，我们只需要循环遍历width次M和N中的元素就可以了*/
    for (int k = 0; k < width; k ++){
        float M_element = M_device[y * width + k];
        float N_element = N_device[k * width + x];
        P_element += M_element * N_element;
    }

    P_device[y * width + x] = P_element;
}

/*
    CUDA中使用block对矩阵中某一片区域进行集中计算。这个类似于loop中的tile
    感兴趣的同学可以试着改一下blockSize，也就是tileSize，看看速度会发生什么样子的变化
    当blockSize达到一个数量的时候，这个程序会出错。下一个案例中我们会分析
*/
void MatmulOnDevice(float *M_host, float *N_host, float* P_host, int width, int blockSize){
    /* 设置矩阵大小 */
    int size = width * width * sizeof(float);

    /* 分配M, N在GPU上的空间*/
    float *M_device;
    float *N_device;

    cudaMalloc(&M_device, size);
    cudaMalloc(&N_device, size);

    /* 分配M, N拷贝到GPU上*/
    cudaMemcpy(M_device, M_host, size, cudaMemcpyHostToDevice);
    cudaMemcpy(N_device, N_host, size, cudaMemcpyHostToDevice);

    /* 分配P在GPU上的空间*/
    float *P_device;
    cudaMalloc(&P_device, size);

    /* 调用kernel来进行matmul计算, 在这个例子中我们用的方案是：将一个矩阵切分成多个blockSize * blockSize的大小 */
    dim3 dimBlock(blockSize, blockSize);
    dim3 dimGrid(width / blockSize, width / blockSize);
    MatmulKernel <<<dimGrid, dimBlock>>> (M_device, N_device, P_device, width);

    /* 将结果从device拷贝回host*/
    cudaMemcpy(P_host, P_device, size, cudaMemcpyDeviceToHost);
    cudaDeviceSynchronize();

    /* Free */
    cudaFree(P_device);
    cudaFree(N_device);
    cudaFree(M_device);
}

这段代码展示了如何使用 CUDA 编程模型在 GPU 上实现矩阵乘法。矩阵乘法是并行计算中一个常见的案例，适合于展示 GPU 的高并行处理能力。下面是对这段代码的详细分析：(from chatGPT)

核函数 MatmulKernel：

• 功能：核函数 MatmulKernel 负责计算矩阵乘法的一个元素。每个线程计算结果矩阵 P 中的一个元素，通过对矩阵 M 和 N 进行逐元素相乘和累加来实现。

• 索引计算：线程的全局索引 (x, y) 通过 blockIdx, blockDim, 和 threadIdx 计算得到，这些变量是 CUDA 提供的内建变量，用于定位当前线程在整个执行网格中的位置。x 和 y 表示当前线程负责计算 P 矩阵中位置为 (y, x) 的元素。

• 计算逻辑：对于 P 矩阵中的每个元素，核函数内部通过一个循环遍历 width 次，每次循环计算 M 的一行与 N 的一列对应元素的乘积并累加到 P_element 上。最后，计算结果 P_element 被写回到 P_device 对应的位置。

主机函数 MatmulOnDevice：

• 功能：这个函数封装了矩阵乘法在 GPU 上的整个过程，包括内存分配、数据传输、核函数调用和清理资源。

• 内存操作：首先，函数为输入矩阵 M 和 N，以及输出矩阵 P 在 GPU 上分配内存。使用 cudaMalloc 分配，cudaMemcpy 用于将矩阵从主机（CPU）内存复制到设备（GPU）内存。

• 核函数调用：dim3 dimBlock(blockSize, blockSize) 和 dim3 dimGrid(width / blockSize, width / blockSize) 定义了核函数的执行配置，包括每个块的维度和网格的维度。这里 blockSize 是一个参数，决定了每个块中线程的布局（blockSize × blockSize），而网格的维度根据矩阵大小和块大小计算得出。MatmulKernel 核函数随后被并行调用，每个线程负责计算输出矩阵 P 的一个元素。

• 结果回传：矩阵乘法完成后，P 矩阵从 GPU 内存复制回主机内存，使用 cudaMemcpyDeviceToHost。

• 资源释放：最后，使用 cudaFree 释放 GPU 上分配的内存。

注意事项：

• 性能考虑：调整 blockSize（即每个块中的线程数）可以影响性能。较大的块可以减少启动核函数的次数，但太大可能会超出 GPU 对每块线程数的限制，导致程序错误。
• 同步操作：cudaDeviceSynchronize 调用确保了所有的 GPU 操作完成后，程序才会继续执行。这对于在数据被复制回主机之前确保所有的 GPU 计算完成是必要的。
• 精度和误差：GPU 上的浮点计算可能与 CPU 上的浮点计算有细微差别，这是由于不同的硬件架构和优化引起的。
• 错误处理：示例代码中未包含 CUDA 调用的错误处理。在实践中，应检查几乎所有 CUDA 调用的返回值，确保正确地处理潜在的错误情况。错误处理对于调试和保证程序稳定运行非常关键。

小结：

在 matmul_gpu_basic.cpp 文件中，通过 CUDA 编程模型实现了矩阵乘法，展示了如何利用 GPU 的高并行能力进行大规模数据计算。该示例涵盖了从数据在主机与设备之间的传输、内存管理到并行计算的核心步骤，为理解和使用 CUDA 进行高性能计算提供了一个基础的框架。此外，通过修改 blockSize 参数，开发者可以实验不同的配置，优化计算性能。

重要的是，要认识到 GPU 编程涉及许多细节，如线程布局的选择、内存访问模式的优化、计算精度的考量以及错误的管理，这些都对最终的计算性能和结果的正确性有着重要影响。因此，开发高效且正确的 GPU 程序不仅需要深入理解 CUDA 编程模型，还需要对所解决问题的计算特性有深入的理解。

最后我们再来看下 main.cpp：

#include <stdio.h>
#include <cuda_runtime.h>

#include "utils.hpp"
#include "timer.hpp"
#include "matmul.hpp"


int seed;
int main(){
    Timer timer;
    int width     = 1<<10; // 1,024
    int min       = 0;
    int max       = 1;
    int size      = width * width;
    int blockSize = 1;

    float* h_matM = (float*)malloc(size * sizeof(float));
    float* h_matN = (float*)malloc(size * sizeof(float));
    float* h_matP = (float*)malloc(size * sizeof(float));
    float* d_matP = (float*)malloc(size * sizeof(float));
    
    seed = 1;
    initMatrix(h_matM, size, min, max, seed);
    seed += 1;
    initMatrix(h_matN, size, min, max, seed);
    
    /* CPU */
    timer.start();
    MatmulOnHost(h_matM, h_matN, h_matP, width);
    timer.stop();
    timer.duration<Timer::ms>("matmul in cpu");

    /* GPU warmup */
    timer.start();
    MatmulOnDevice(h_matM, h_matN, d_matP, width, blockSize);
    timer.stop();
    timer.duration<Timer::ms>("matmul in gpu(warmup)");

    /* GPU general implementation, bs = 16*/
    blockSize = 16;
    timer.start();
    MatmulOnDevice(h_matM, h_matN, d_matP, width, blockSize);
    timer.stop();
    timer.duration<Timer::ms>("matmul in gpu(bs = 16)");
    compareMat(h_matP, d_matP, size);

    /* GPU general implementation, bs = 1*/
    blockSize = 1;
    timer.start();
    MatmulOnDevice(h_matM, h_matN, d_matP, width, blockSize);
    timer.stop();
    timer.duration<Timer::ms>("matmul in gpu(bs = 1)");
    compareMat(h_matP, d_matP, size);

    /* GPU general implementation, bs = 32*/
    blockSize = 32;
    timer.start();
    MatmulOnDevice(h_matM, h_matN, d_matP, width, blockSize);
    timer.stop();
    timer.duration<Timer::ms>("matmul in gpu(bs = 32)");
    compareMat(h_matP, d_matP, size);
    return 0;
}

main.cpp 是 CUDA 矩阵乘法示例程序的入口文件，其中包括了矩阵乘法在 CPU 和 GPU 上的实现和性能比较。程序定义了几个关键变量，如矩阵大小 width、值的范围 min 和 max、以及 blockSize，后者用于在 GPU 上的矩阵乘法中定义每个块的线程布局。接着，程序分配了 CPU 上用于存储两个输入矩阵和一个输出矩阵的内存。(from chatGPT)

通过调用 initMatrix 函数初始化两个输入矩阵，程序首先在 CPU 上执行矩阵乘法，记录并打印出执行时间。随后，程序执行了几轮 GPU 上的矩阵乘法以测试不同的 blockSize 参数对性能的影响。每次 GPU 计算之前，都有一个“预热”运行（blockSize 为 1），以确保 GPU 处于激活状态，随后分别以 blockSize 为 16、1 和 32 进行实际的计算，每次都测量并输出执行时间。最后，使用 compareMat 函数比较 CPU 和 GPU 计算结果的正确性。

值得注意的是我们在测试核函数的执行时间时刻会有一个 warmup 的操作，在 GPU 编程中 warmup 操作是一个常见的实践，目的是为了将 GPU 从低功耗模式唤醒到高性能模式。当 GPU 首次被用于计算时，它可能需要一段时间来从较低的功率状态过渡到完全激活的状态，这个过程中可能会有较高的延迟和较低的性能。通过执行一个轻量级的计算任务作为 warmup，可以确保 GPU 在随后的计算任务中能够以最佳状态运行。

还有这里的时间计算实现是不能够非常精准的获取 kernel 函数的执行时间的，如果想要精准的获取 kernel 实现需要通过 cuda event 来进行测量，这个在后面的案例中会讲

这里我们使用的 blockSize = 16，我们可以看看 blockSize 变化时的结果，如下图所示：

blockSize=1

可以看到 blockSize=1 时 GPU 上核函数的执行时间明显变长了是 113.298ms，效率不是很高，大家也可以试下其它的 blockSize 看下输出时间结果的变化

2. CUDA中的error handler

本小节目标：养成使用 cuda 的 error handler 进行良好的编程习惯

这个部分给大家去讲 CUDA 中的 error handler 错误排查，通过这部分的内容讲解希望大家能养成一个良好的 CUDA 编程习惯，也就是如何利用 error handler 定位错误以及错误产生的原因，然后跟大家分享一下 CUDA 中的一些错误的表现形式和传播形式

这个小节的案例是 2.4-error-hanlder，如下图所示：

这个小节跟 2.3 小节内容一样，只不过在上面做了一些 error handler 上的扩展，大家在矩阵乘法的案例可能会做很多其它的测试，比如将 block size 设置为 2x2，4x4，8x8 等等，测试后我们会发现当设置为 64x64 时会报错，如果我们不使用 error handler 它的错误其实很隐藏，我们发现不出来，如果使用了 error handler 我们可以得到如下的信息：

在终端上我们会打印错误发生的地方，错误发生的原因等等，可以帮助我们进行调试，

2.4 小节案例中新添加的东西主要有两部分：

utils.hpp

首先在 utils.hpp 中添加了两个宏定义，一个是对 cuda runtime api 使用的 error handling，一个是对核函数的 error handling

matmul_gpu_basic.cu

然后在 matmul_gpu_basic.cu 中对于每一个 cuda_runtime_api 的使用都进行了错误检查，以及 kernel 执行结束后对 kernel 进行错误排查

2.1 为什么需要有error handler

我们为什么需要有 error handler 呢？它对我们有什么帮助呢？下面我们通过一个例子来简单讲下：

我们在代码中人为设置了一个错误，将 block size 的大小设置为 64x64，我们知道一个 block 中的 thread 数量是有上限的最多设置为 1024，如果设置为 64x64=4096 个线程的话它在 CUDA 核函数的配置上是有错误的，这个配置是无法启动你的核函数的

下面我们来对比看下没有 error handler 以及有 error handler 时程序执行的情况

没有error handling的时候

有error handling的时候

我们对比后可以发现，如果程序没有做任何 error handler 的处理，我们 GPU 上的执行结果是 0 也就是 cudaMalloc 初始化之后的一个结果，它的值没有发生任何改变，错误发生的坐标也是 [0,0]，也就是一开始就出错了，核函数压根没有执行。那如果我们程序加入了 error handler，我们可以发现它在 src/matmul_gpu_basic.cu 的 63 行发现了错误，错误类型是配置错误，错误的原因的无效的配置参数导致的，通过 error handler 我们就能够精准的获取错误信息，这对我们编程调试都是非常有帮助的，也是推荐大家去这么做的。

关于 error handler 的实现，这里还有几点需要大家注意：

一般来说，我们习惯把 cuda 的 error handler 定义成宏，如上图所示，因为这样可以避免在执行的时候发生调用 error handler 而引起的 overhead，此外代码中的 __FILE__ 和 __LINE__ 的含义如下：

• __FILE__：编译器内部定义的一个宏，表示的是当前文件的文件名
• __LINE__：编译器内部定义的一个宏，表示的是当前文件的行

一个良好的 CUDA 编程习惯就是在调用任何一个 cuda runtime api 的时候使用 error handler 进行包装，这样可以方便我们排查错误的来源

上图是 NVIDIA CUDA 官方文档中关于两种错误的一个解释，cudaGetLastError 它会返回一个最近的错误，同时它还会 reset 系统的错误状态，帮你 reset 到 cuda success。cudaPeekAtLastError 它也会返回一个错误，但它不会把你整个系统的 error 状态给返回到 cuda success

二者的差别在于错误是否传播，对于不可恢复的错误，我们需要使用 cudaGetLastError，因为如果发生了这种类型的错误并且不把系统的状态进行 reset 的话，错误会一直传播下去，导致后面的即便正确的 api 使用也会产生同样的错误

那大家可能会想什么叫做不可恢复的错误呢？其实不可恢复的（non-recoverable/sticky）一般指的是核函数内部的执行错误，比较典型的例子就是内存访问越界。相比之下，比如像 block size 以及 shared memory size 这种配置错误就属于可恢复的（recoverable/non-sticky）错误

我们除了可以将错误分为不可恢复/可恢复错误（sticky/non-sticky error）之外，还可以将它们分为同步/异步错误（synchronous/asynchronous error），我们这里稍微扩展一下，我们知道核函数是异步执行的，如果你的 CPU 端调用核函数后不做同步的话，我们 CPU 就会往下执行其它相关指令了，所以我们错误也有同步和异步之分

比如我们在启动核函数的时候发生的一些 block size、shared memory size 等的配置错误，它属于启动时就会发生的错误，不需要我们等待就可以直接找到错误，这个就是异步错误。但是也有一些情况就是我们需要等待你的核函数全部执行完我们才知道到底有没有出现诸如内存访问越界的错误，这种我们可以认为是同步错误

以上就是关于 CUDA 中错误处理的方式，大家可以了解下

3. 获取GPU信息

本小节目标：学习如何使用 cuda runtime api 显示 GPU 硬件信息，以及理解 GPU 硬件信息的重要性

这部分教大家如何去获取 GPU 的硬件信息，当然大家根据自己 GPU 的架构可以从 NVIDIA 官网或者网络上搜索出很多相关的信息，但是我们也需要知道可以使用 CUDA Runtime API 提供的一些接口去获取你的 GPU 信息，包括你的每个 block 中 thread 的最大数量、SM 数量、shared memory 大小等等，这些信息可能是我们在做模型部署优化的时候需要考虑的一些东西，它们的信息都可以通过 CUDA Runtime API 获取到

3.1 执行一下我们的第五个CUDA程序

本小节对应的案例是 2.5-device-info，其目录结构如下图所示：

我们先 make run 执行看下输出结果，如下图所示：

我们可以看到终端上打印了当前 GPU 的各种硬件信息，开源看到博主使用的 GPU 是 RTX3060，它的各种信息都被打印出来，相关的核心代码如下所示：

// main.cpp
#include <stdio.h>
#include <cuda_runtime.h>
#include <string>

#include "utils.hpp"

int main(){
    int count;
    int index = 0;
    cudaGetDeviceCount(&count);
    while (index < count) {
        cudaSetDevice(index);
        cudaDeviceProp prop;
        cudaGetDeviceProperties(&prop, index);
        LOG("%-40s",             "*********************Architecture related**********************");
        LOG("%-40s%d%s",         "Device id: ",                   index, "");
        LOG("%-40s%s%s",         "Device name: ",                 prop.name, "");
        LOG("%-40s%.1f%s",       "Device compute capability: ",   prop.major + (float)prop.minor / 10, "");
        LOG("%-40s%.2f%s",       "GPU global meory size: ",       (float)prop.totalGlobalMem / (1<<30), "GB");
        LOG("%-40s%.2f%s",       "L2 cache size: ",               (float)prop.l2CacheSize / (1<<20), "MB");
        LOG("%-40s%.2f%s",       "Shared memory per block: ",     (float)prop.sharedMemPerBlock / (1<<10), "KB");
        LOG("%-40s%.2f%s",       "Shared memory per SM: ",        (float)prop.sharedMemPerMultiprocessor / (1<<10), "KB");
        LOG("%-40s%.2f%s",       "Device clock rate: ",           prop.clockRate*1E-6, "GHz");
        LOG("%-40s%.2f%s",       "Device memory clock rate: ",    prop.memoryClockRate*1E-6, "Ghz");
        LOG("%-40s%d%s",         "Number of SM: ",                prop.multiProcessorCount, "");
        LOG("%-40s%d%s",         "Warp size: ",                   prop.warpSize, "");

        LOG("%-40s",             "*********************Parameter related************************");
        LOG("%-40s%d%s",         "Max block numbers: ",           prop.maxBlocksPerMultiProcessor, "");
        LOG("%-40s%d%s",         "Max threads per block: ",       prop.maxThreadsPerBlock, "");
        LOG("%-40s%d:%d:%d%s",   "Max block dimension size:",     prop.maxThreadsDim[0], prop.maxThreadsDim[1], prop.maxThreadsDim[2], "");
        LOG("%-40s%d:%d:%d%s",   "Max grid dimension size: ",     prop.maxGridSize[0], prop.maxGridSize[1], prop.maxGridSize[2], "");
        index ++;
        printf("\n");
    }
    return 0;
}

上述代码通过 CUDA Runtime API 来查询和展示系统中所有 CUDA 兼容 GPU 设备信息的程序。首先，它调用 cudaGetDeviceCount 以获取系统中 GPU 的数量。接着，通过一个循环，针对每个 GPU 设备执行以下步骤：设置当前设备上下文、获取该设备的属性（如设备名称、计算能力、全局内存大小、L2 缓存大小等），并使用一个自定义的 LOG 宏打印出这些信息。此外，还包括了设备的性能参数，比如每个多处理器的最大块数、每个块的最大线程数、最大块维度和网格维度大小等。

代码中展示了一些比较重要的 device info，比如 shared memory 大小，SM 数量，thread per block 的最大值等等，还要很多其它的没有展示，大家可以自行调查和使用一下。大家可能会发现我们在打印的时候并没有使用 printf 或者 cout 这种，而是使用的 LOG 宏，这个宏定义在 util.hpp 文件中，其内容如下：

#define LOG(...)                     __log_info(__VA_ARGS__)

// 使用变参进行LOG的打印。比较推荐的打印log的写法
static void __log_info(const char* format, ...) {
    char msg[1000];
    va_list args;
    va_start(args, format);

    vsnprintf(msg, sizeof(msg), format, args);

    fprintf(stdout, "%s\n", msg);
    va_end(args);
}

这段代码定义了一个用于日志打印的宏 LOG 和一个支持变参的函数 __log_info，旨在以格式化的方式输出日志信息到标准输出。大家如果看过一些开源项目会发现很多时候都会把打印日志的方法定义成一个宏封装起来方便使用，比如说 tensorRT_Pro 项目中的 INFO、INFOE 等等，那大家可能对 __log_info 函数中的 __VA_ARGS_ 参数有所困惑，它其实是编译器内部定义的一个宏，表示的是变参，配合 vsnprintf 函数开源将 LOG 中的变参信息存入到 msg 的这个 buffer 中，最终在一起打印出来。

3.2 为什么要注意硬件信息

最后我们给大家讲一下为什么要去关注硬件信息

上图是后面课程中一个案例的 Makefile 文件，它需要根据 GPU 来指定 compute capability 参数，供 nvcc 编译 CUDA 程序时使用，nvcc 会根据你的 archer 信息来选择不同的编译策略，比如说我们这里的 3060 的话，它在 computer capability 的 architecture 是 86，而像 Tesla A100 或者 Tesla V100 显卡它的 architecture 就是 compute_80、compute_70，你的 code 就是 SM80、SM70

那我们在 nvcc 编译的时候其实需要把这个给添加进去，这个信息的话就可以通过我们上面刚讲的 CUDA Runtime API 提供的接口去获取，除了编译时需要知道 GPU 硬件的 compute capability 信息外，我们在编程时还需要知道以下信息

首先我们需要知道我们在启动核函数的时候，配置信息的规定都有哪些，比如 block 的最大数量，每个 block 中 thread 的最大数量，这样我们在配置 block size 和 grid size 的时候就可以根据获取的这些信息进行相应的设置

*********************Parameter related************************
Max block numbers:                      16
Max threads per block:                  1024
Max block dimension size:               1024:1024:64
Max grid dimension size:                2147483647:65535:65535

此外，我们都知道 shared memory 的使用对 CUDA 程序的加速很重要，那当我们在使用 shared memory 的时候是不是需要知道它的大小上限是多少呢

Shared memory per block:                48.00KB
Shared memory per SM:                   100.00KB

还有你的一个 warp 线程束是由多少个thread 组成的等等这些信息也是非常重要的

Warp size:                              32

最后当我们在进行性能调优的时候，内存的大小和内存的带宽是我们需要考虑的一个很重要的因素，结合 roofline model（后续会讲），我们需要寻找想要隐藏 memory 的数据传输所造成的 overhead，需要多少的计算量和计算效率，那我们需要考虑的东西就有 memory bandwidth 内存带宽、clock rate 时钟频率、memory clock rate 内存时钟频率，SM 的数量等等，我们根据当前 GPU 架构的这些信息去计算我们的 roofline model 是什么样子的，如下图所示：

Device clock rate:                      1.84GHz
Device memory clock rate:               7.50Ghz
Number of SM:                           28

以上就是使用 cuda runtime API 来获取 GPU 信息的全部内容了，大家可以在自己的 GPU 上进行相关测试

总结

本次课程我们学习了使用 CUDA 进行矩阵乘法的加速，首先我们分析了如何利用 CUDA Core 去进行矩阵乘法的加速，然后我们分析了相关代码，对比了在 CPU 和 GPU 端进行矩阵乘法的耗时，结果表明 GPU 上的并行处理对矩阵乘法的加速非常明显。然后我们聊了一下 CUDA 中的错误处理，对于可恢复错误和不可恢复错误进行了讲解，一个良好的编程习惯是需要有 error handler 的，这可以帮助我们调试程序。最后我们学习了如何利用 cuda runtime api 来获取我们 GPU 的各自信息，这对于我们后续优化核函数非常有帮助

OK，以上就是第 2 小节有关使用 CUDA 进行矩阵乘法加速的全部内容了，下节我们来学习共享内存以及 Bank Conflict，敬请期待😄

参考

• 3.9.cuda运行时API-错误处理的理解以及错误的传播特性

• CUDA software stack (NVIDIA, 2007)

• NVIDIA, cuda toolkit documentation(cuda v12.2.0)

• Roofline Performance Model