L
O
A
D
I
N
G

cuda 入门

cuda
笔记
发布日期:   2024-01-23
更新日期:   2026-06-08
文章字数:   7.5k
阅读时长:   29 分
阅读次数:  

CUDA入门

前言

06年，NVIDIA公司发布了CUDA，它是建立在NVIDIA的GPU上的一个通用并行计算平台和编程模型。CUDA编程使得可以更高效地利用GPU的并行计算引擎来解决复杂的计算问题。近年来，GPU在深度学习领域取得了巨大成功，基于GPU的并行计算已成为训练深度学习模型的标准。

需要强调的是，GPU并非独立运行的计算平台，而是需要与CPU协同工作的。它可以视作CPU的协处理器，在提及GPU并行计算时，实际上指的是基于CPU+GPU的异构计算架构。在这种异构计算架构中，GPU与CPU通过PCIe总线连接在一起，协同完成计算任务。CPU所在位置被称为主机端（host），而GPU所在位置被称为设备端（device），如下图所示：

image-20231219141540082

当比较GPU和CPU时，GPU拥有更多的计算核心，使其非常适合处理数据并行的计算密集型任务，例如大型矩阵运算。相比之下，CPU的核心数量较少，但它们擅长执行复杂的逻辑运算，适用于控制密集型任务。此外，CPU上的线程通常是重型的，导致较大的上下文切换开销，而GPU拥有许多轻量级线程。

因此，基于CPU+GPU的异构计算平台能够充分发挥优势。CPU负责处理串行程序中的逻辑复杂部分，而GPU则专注于处理数据密集型、并行计算的部分，这样能够实现最佳性能。

image-20231219143336313

CUDA是NVIDIA公司开发的GPU编程模型，它提供了简单易用的接口，可用于构建基于GPU计算的应用程序。CUDA支持多种编程语言，如C/C++、Python和Fortran等。下面我们将使用CUDA C/C++接口来讲解CUDA编程。Windows系统下的CUDA安装教程可以参考给出的链接进行安装：CUDA安装及环境配置——最新详细版-CSDN博客

入门基础

在CUDA编程模型中，涉及两个核心概念：host和device。

host代表着CPU及其相关内存，而device代表着GPU及其相关内存。CUDA程序中同时存在host程序和device程序，分别在CPU和GPU上运行。此外，host和device之间可以进行通信，允许数据在它们之间进行传输。

一般而言，典型的CUDA程序执行流程如下：

  1. 分配host内存并对数据进行初始化
  2. 分配device内存，并将数据从host拷贝到device上
  3. 调用CUDA的核函数（kernel）在device上完成指定的并行运算
  4. 将device上计算得到的结果拷贝回host
  5. 释放device和host上分配的内存空间

核函数（kernel）是CUDA中的关键概念，它是在device上并行执行的函数。核函数以 __global__ 修饰符声明，在调用时需要使用 <<<grid, block>>> 来指定kernel要执行的线程数量。在CUDA中，每个线程都会执行核函数，并且每个线程都会有一个唯一的线程ID，可以通过核函数的内置变量 threadIdx 来获取。

由于GPU是异构计算模型，因此需要区分host和device上的代码。在CUDA中，可以通过函数类型限定词来区分host和device上的函数。主要的三种函数类型限定词如下：

  • __global__：在device上执行，可以从host调用（一些特定的GPU也能从device上调用）。返回类型必须是 void，不支持可变参数，也不能成为类的成员函数。需要注意的是，用 __global__ 定义的kernel是异步执行的，这意味着host不会等待kernel执行完毕再执行下一步。
  • __device__：在device上执行，仅可以从device中调用，不能与 __global__ 同时使用。
  • __host__：在host上执行，仅可以从host中调用。一般情况下，此修饰符可以省略不写。它不能与 __global__ 同时使用，但可以与 __device__ 结合使用，此时函数可在device和host上都编译。

要深入理解kernel，需要对kernel的线程层次结构有清晰的认识。

  • 线程结构的第一层次：GPU上存在大量并行化的轻量级线程。当kernel在device上执行时，实际上是启动了大量线程，一个kernel启动的所有线程构成一个网格（grid）
  • 第二个层次：在同一个网格内，一些线程共享相同的全局内存空间，网格被分为多个部分，每一个都是一个线程块（block）一个线程块包含多个线程。
  • 第三个层次：在线程块内，线程被划分为更小的单元，称为线程（thread）。这是GPU上的最细粒度的并行层次。线程是执行计算的基本单元，它们共享线程块内的资源，包括寄存器文件和共享内存。通过合理的组织和协同工作，线程能够高效地执行各自的任务，实现更细致的并行化。

线程的层次结构如下图所示，以二维的grid和block为例。grid和block都被定义为 dim3 类型的变量，dim3 可以看作包含三个无符号整数（x，y，z）成员的结构体变量，在定义时，缺省值默认为1。因此，grid和block可以灵活定义为1维、2维或3维结构。在CUDA编程中，这些概念和规则是关键，对于理解和利用GPU并行计算极为重要。

要调用kernel，必须通过执行配置 <<<grid, block>>> 来指定kernel所使用的线程数及结构。

1
2
3
dim3 grid(3, 2);
dim3 block(5, 3);
kernel_fun<<< grid, block >>>(prams...);
Kernel的两层线程组织结构（2-dim）

所以，一个线程需要两个内置的坐标变量（blockIdx，threadIdx）来唯一标识，它们都是dim3类型变量，其中blockIdx指明线程所在grid中的位置，而threaIdx指明线程所在block中的位置，如图中的Thread (1,1)满足：

1
2
3
4
threadIdx.x = 1
threadIdx.y = 1
blockIdx.x = 1
blockIdx.y = 1

在CUDA编程中，线程块（block）中的线程被分配到同一个流式多处理器（SM）上，但是每个SM的资源是有限的，因此线程块中的线程数量受到限制。现代GPU支持的线程块中的线程数量通常可达1024个。有时候，需要知道一个线程在block中的全局ID，这就需要了解block的组织结构，可以通过线程的内置变量 blockDim 来获取。它获取线程块在各个维度上的大小。对于一个2维的block( Dx, Dy)，线程(x, y)的ID值为 (x + y * Dx )，如果是3维的block( Dx, Dy, Dz)，线程(x, y, z)的ID值为 (x + y * Dx + z * Dx * Dy)。此外，线程还具有内置变量 gridDim，用于获取网格块在各个维度上的大小。

kernel的这种线程组织结构天然适合vector、matrix等运算，如我们将利用上图2-dim结构实现两个矩阵的加法，每个线程负责处理每个位置的两个元素相加，代码如下所示。线程块大小为(16, 16)，然后将N * N大小的矩阵均分为不同的线程块来执行加法运算。

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
__global__ void matrixAddition(float *a, float *b, float *c, int N) {
    int row = blockIdx.y * blockDim.y + threadIdx.y;
    int col = blockIdx.x * blockDim.x + threadIdx.x;

    if (row < N && col < N) {
        int index = row * N + col;
        c[index] = a[index] + b[index];
    }
}

int main() {
    // 其他代码...

    dim3 blockSize(16, 16);
    dim3 gridSize((N + blockSize.x - 1) / blockSize.x, (N + blockSize.y - 1) / blockSize.y);

    matrixAddition<<<gridSize, blockSize>>>(device_a, device_b, device_c, N);

    // 其他代码...
}

以上代码中，kernel函数 matrixAddition 实现了矩阵加法操作。通过计算每个线程的行列索引，确保每个线程仅处理有效的矩阵元素。主函数中，根据矩阵大小 N 和设定的线程块大小，计算出适当的 grid 和 block 尺寸，并调用 kernel 函数执行矩阵加法。

内存模型

这里简单介绍一下CUDA的内存模型，如下图所示：

CUDA内存模型

可以看到，每个线程有自己的私有本地内存（Local Memory），而每个线程块有包含共享内存（Shared Memory）,可以被线程块中所有线程共享，其生命周期与线程块一致。此外，所有的线程都可以访问全局内存（Global Memory）。还可以访问一些只读内存块：常量内存（Constant Memory）和纹理内存（Texture Memory）

这些内存类型的特点如下：

  1. 私有本地内存（Local Memory）：每个线程都有自己的私有本地内存，用于存储线程私有的变量和数据。这些变量通常是在寄存器不足时或者动态分配内存时使用。访问本地内存的速度较慢，因此需要尽量减少对本地内存的访问，以提高性能。

  2. 共享内存（Shared Memory）：每个线程块都拥有共享内存，它可以被该线程块中的所有线程共享。共享内存的生命周期与线程块一致，对于需要在线程块内多个线程之间共享的数据，可以存储在共享内存中。共享内存的访问速度比全局内存快得多，因此可以用于提高线程块内的数据访问效率和并行计算速度。

  3. 全局内存（Global Memory）：所有的线程都可以访问全局内存，它是GPU中容量最大的内存类型。全局内存通常用于存储整个程序的数据，但其访问速度相对较慢。因此，在CUDA编程中，需要谨慎优化全局内存的访问，以最大化性能。

  4. 常量内存（Constant Memory）：常量内存用于存储对所有线程只读的数据，其特点是具有缓存，提供对只读内存的快速访问，适用于线程内的数据共享。

  5. 纹理内存（Texture Memory）：纹理内存是一种特殊的只读内存，适用于图形处理领域，提供对数据的高效访问和插值。它通常用于图像处理等应用场景。

CUDA内存模型中不同类型的内存具有不同的访问特性和速度，合理地利用和管理这些内存类型可以对CUDA程序的性能产生显著影响。优化内存访问是CUDA编程中的一个重要方面，能够有效地提高程序的执行效率。

硬件认识

还有重要一点，你需要对GPU的硬件实现有一个基本的认识。上面说到了kernel的线程组织层次，那么一个kernel实际上会启动很多线程，这些线程是逻辑上并行的，但是在物理层却并不一定。

在GPU硬件中，一个关键组件是SM（Streaming Multiprocessor，流式多处理器）。SM包含CUDA核心、共享内存、寄存器等。SM能够同时执行数百个线程，其并发能力取决于SM拥有的资源数量。当一个kernel被执行时，其线程块被分配到SM上。一个线程块只能在一个SM上被调度，而一个SM通常能够调度多个线程块，这取决于SM本身的能力。因此，一个kernel的线程块可能被分配到多个SM上，所以grid只是逻辑层，而SM才是执行的物理层。

SM采用的是SIMT (Single-Instruction, Multiple-Thread，单指令多线程)架构，基本的执行单元是线程束（warps)，每个线程束包含32个线程。这些线程同时执行相同的指令，但每个线程都有自己的指令地址计数器和寄存器状态，拥有独立的执行路径。因此，尽管线程束中的线程同时执行相同的指令，但它们可能具有不同的行为，特别是在遇到分支结构时，可能导致一些线程进入分支，而其他线程不执行，这些非执行的线程只能等待，因为GPU规定线程束中所有线程在同一周期执行相同的指令，线程束分化会导致性能下降。

当一个线程块被划分到SM上时，它会进一步划分为多个线程束，因为这才是SM的基本执行单元。然而，一个SM同时并发的线程束数量是有限的，这受限于资源分配，因为SM要为每个线程块分配共享内存，并为每个线程束中的线程分配独立的寄存器。因此，SM的配置会影响其所支持的线程块和线程束的并发数量。总的来说，虽然kernel的grid和block只是逻辑上的划分，但在物理层面，一个kernel的所有线程并非一定同时并发。因此，需要注意调整grid和block的配置，以充分发挥GPU的性能。此外，由于线程束的基本大小是32，因此block大小最好设置为32的倍数。

逻辑层和物理层

个人配置

在进行CUDA编程前，可以先检查一下自己的GPU硬件配置，这样才能有的放矢，通过以下代码获得GPU的配置属性：

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
#include <stdio.h>
#include <cuda_runtime.h>

int main()
{
    cudaDeviceProp deviceProp;
    cudaGetDeviceProperties(&deviceProp, 0);

    printf("Device 0 information:\n");
    printf("设备名称与型号: %s\n", deviceProp.name);
    printf("显存大小: %d MB\n", (int)(deviceProp.totalGlobalMem / 1024 / 1024));
    printf("计算能力: %d.%d\n", deviceProp.major, deviceProp.minor);
    printf("含有的SM数量: %d\n", deviceProp.multiProcessorCount);
    printf("CUDA CORE数量: %d\n", deviceProp.multiProcessorCount * 192);
    printf("每个线程块的共享内存大小: %.2f KB\n", deviceProp.sharedMemPerBlock / 1024.0);
    printf("每个线程块的最大线程数: %d\n", deviceProp.maxThreadsPerBlock);
    printf("每个EM的最大线程数: %d\n", deviceProp.maxThreadsPerMultiProcessor);
    printf("每个SM的最大线程束数: %d\n", deviceProp.maxThreadsPerMultiProcessor / 32);

    return 0;
}

向量加法实例

知道了CUDA编程基础，我们就来个简单的实战，利用CUDA编程实现两个向量的加法，在实现之前，先简单介绍一下CUDA编程中内存管理API。

在CUDA编程中，内存管理是至关重要的。以下是常用的内存管理函数的简要介绍：

  1. cudaMalloc() 函数：

    1
    cudaError_t cudaMalloc(void** devPtr, size_t size);

    cudaMalloc 用于在设备（GPU）上分配内存。它类似于C语言中的 malloc() 函数，但是在GPU的全局内存中分配内存。devPtr 是指向分配内存起始位置的指针，size 表示要分配的内存字节数。分配成功时返回 cudaSuccess，否则返回相应的错误码。需要注意，使用完内存后，必须使用 cudaFree() 函数释放相应的设备内存，以避免内存泄漏。

    这个函数和C语言中的malloc类似，但是在device上申请一定字节大小的显存，其中devPtr是指向所分配内存的指针。

    同时要释放分配的内存使用cudaFree函数，这和C语言中的free函数对应。

  2. cudaFree() 函数：

    1
    cudaError_t cudaFree(void* devPtr);

    cudaFree 用于释放之前由 cudaMalloc 分配的设备内存。传入的参数 devPtr 是之前分配的设备内存的指针。释放成功时返回 cudaSuccess，否则返回相应的错误码。

  3. cudaMemcpy() 函数：

    1
    cudaError_t cudaMemcpy(void* dst, const void* src, size_t count, cudaMemcpyKind kind);

    cudaMemcpy 用于在主机（host）和设备（device）之间进行数据传输。它可以控制数据的复制方向，例如从主机到设备、从设备到主机、设备到设备等。dst 和 src 分别是目标和源的指针，count表示要复制的字节数，kind指定了数据传输的方向。常用的传输方向有：cudaMemcpyHostToHost、cudaMemcpyHostToDevice、cudaMemcpyDeviceToHost、cudaMemcpyDeviceToDevice，根据英文即可知道传输方向。

在CUDA程序中，合理使用这些函数进行内存分配和数据传输是保证程序正确性和性能的重要一步。

现在我们来实现一个向量加法的实例，这里grid和block都设计为1-dim，首先定义kernel如下：

1
2
3
4
5
6
7
8
9
10
11
12
// 两个向量加法kernel，grid和block均为一维
__global__ void add(float* x, float * y, float* z, int n)
{
    // 获取全局索引
    int index = threadIdx.x + blockIdx.x * blockDim.x;
    // 步长
    int stride = blockDim.x * gridDim.x;
    for (int i = index; i < n; i += stride)
    {
        z[i] = x[i] + y[i];
    }
}

其中stride是整个grid的线程数，有时候向量的元素数很多，这时候可以将在每个线程实现多个元素（元素总数/线程总数）的加法，相当于使用了多个grid来处理，这是一种grid-stride loop方式，不过下面的例子一个线程只处理一个元素，所以kernel里面的循环是不执行的。

下面是一个简单的向量加法示例，使用CUDA编程中的基本概念和函数进行向量相加操作：

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
#include <stdio.h>
#include <stdlib.h>
#include <math.h>
#include <cuda_runtime.h>

// 两个向量加法kernel，grid和block均为一维
__global__ void add(float* x, float* y, float* z, int n) {
    // 获取全局索引
    int index = threadIdx.x + blockIdx.x * blockDim.x;
    // 步长
    int stride = blockDim.x * gridDim.x;

    for (int i = index; i < n; i += stride) {
        z[i] = x[i] + y[i];
    }
}

int main() {
    int n = 1 << 20; // 向量大小
    int nBytes = n * sizeof(float);               
    float *x, *y, *z; // 定义主机上的向量指针
    float *dev_x, *dev_y, *dev_z; // 定义设备上的向量指针

    // 在主机上分配内存
    x = (float *)malloc(nBytes);
    y = (float *)malloc(nBytes);
    z = (float *)malloc(nBytes);

    // 在设备上分配内存
    cudaMalloc((void **)&dev_x, nBytes);
    cudaMalloc((void **)&dev_y, nBytes);
    cudaMalloc((void **)&dev_z, nBytes);

    // 初始化主机上的向量数据
    for (int i = 0; i < n; i++) {
        x[i] = 1.0f;
        y[i] = 2.0f;
    }

    // 将主机上的向量数据复制到设备
    cudaMemcpy(dev_x, x, nBytes, cudaMemcpyHostToDevice);
    cudaMemcpy(dev_y, y, nBytes, cudaMemcpyHostToDevice);

    // 定义 grid 和 block 的大小
    int blockSize = 256;
    int numBlocks = (n + blockSize - 1) / blockSize;

    // 调用 kernel 函数
    add<<<numBlocks, blockSize>>>(dev_x, dev_y, dev_z, n);

    // 将设备上的结果复制回主机
    cudaMemcpy(z, dev_z, nBytes, cudaMemcpyDeviceToHost);

    // 验证结果
    float maxError = 0.0;
    for (int i = 0; i < n; i++) {
        maxError = fmaxf(maxError, fabsf(z[i] - 3.0f));
    }
    printf("最大误差: %f\n", maxError);

    // 释放设备和主机上的内存
    cudaFree(dev_x);
    cudaFree(dev_y);
    cudaFree(dev_z);
    free(x);
    free(y);
    free(z);

    return 0;
}

这段代码实现了向量加法操作，它使用了CUDA的核函数 add，将两个向量 x 和 y 相加得到结果向量 z。在主机和设备之间进行内存分配、数据传输和核函数调用，最后将结果从设备复制回主机进行验证。这里的向量大小为1<<20，而block大小为256，grid大小是4096，kernel的线程层级结构如下图所示：

kernel的线程层次结构

nvprof

NVIDIA nvprof / nvvp工具是英伟达N卡GPU编程中用于观察的利器。全称是NVIDIA Visual Profiler，是由2008年起开始支持的性能分析器。交互性好，利于使用。其中记录运行日志时使用命令nvprof，可视化显示日志时使用命令nvvp。

不过在最近几年，英伟达官方推出了新的性能分析工具NSight，官方更加建议使用新的工具，给出的原因是NSight运行时消耗的资源更少，统计的数据更加贴近实际运行情况的数据。相比之下使用nvprof/nvvp方式运行时消耗资源较多，数据统计容易不准确。

使用

nvprof 是 NVIDIA 提供的一个用于分析 CUDA 应用程序性能的命令行工具。它可以用来收集和分析 CUDA 应用程序的性能数据，包括核函数执行时间、内存使用情况等。

使用 nvprof 工具来分析 CUDA 核函数运行情况的基本步骤：

  1. 编译代码：首先，确保代码已经被编译成可执行文件。

  2. 运行 nvprof：在命令行中运行 nvprof 并指定要分析的可执行文件。

    1
    nvprof ./your_executable

    这将启动您的 CUDA 可执行文件并开始收集性能数据。

  3. 查看分析结果：nvprof 会输出核函数的执行时间、内存使用情况等相关信息。通常，您可以在输出中找到有关核函数执行时间、内存传输时间、资源使用情况等方面的详细信息。该输出可以帮助您确定哪些部分的代码可能需要优化或哪些操作占用了较多的时间和资源。

  4. 选择特定的分析选项：nvprof 还提供了许多选项，可用于指定收集的数据类型、输出格式、分析范围等。您可以根据需要使用不同的选项进行详细的分析。

例如，您可以通过以下方式使用 nvprof 对示例程序进行分析：

1
nvprof ./your_executable

这会收集并显示有关您的 CUDA 可执行文件的基本性能数据。您还可以使用其他选项来进行更详细的分析，例如收集内存信息、指定要分析的核函数等等。可以通过 nvprof --help 查看更多可用选项的详细信息。

请注意，nvprof 是一个强大的工具，可以帮助您深入了解 CUDA 应用程序的性能瓶颈和优化机会。

使用nvprof工具可以分析kernel运行情况，结果如下所示，可以看到kernel函数费时约1.5ms。

1
2
3
4
5
6
7
8
nvprof cuda9.exe
==7244== NVPROF is profiling process 7244, command: cuda9.exe
最大误差: 4.31602e+008
==7244== Profiling application: cuda9.exe
==7244== Profiling result:
            Type  Time(%)      Time     Calls       Avg       Min       Max  Name
 GPU activities:   67.57%  3.2256ms         2  1.6128ms  1.6017ms  1.6239ms  [CUDA memcpy HtoD]
                   32.43%  1.5478ms         1  1.5478ms  1.5478ms  1.5478ms  add(float*, float*, float*, int)

你调整block的大小，对比不同配置下的kernel运行情况，我这里测试的是当block为128时，kernel费时约1.6ms，而block为512时kernel费时约1.7ms，当block为64时，kernel费时约2.3ms。看来不是block越大越好，而要适当选择。

统一内存管理

在前面的实现中，我们需要单独在host和device上进行内存分配，并且要进行数据拷贝，这是很容易出错的。好在CUDA 6.0引入统一内存（Unified Memory）来避免这种麻烦，简单来说就是统一内存使用一个托管内存来共同管理host和device中的内存，并且自动在host和device中进行数据传输。

cudaMallocManaged是CUDA编程中用于动态分配内存的一个函数，它为设备和主机提供统一的、可共享的内存空间。这个函数属于CUDA运行时API的一部分。

函数：

1
cudaError_t cudaMallocManaged(void **devPtr, size_t size, unsigned int flags=0);

参数解析：

  1. void **devPtr：输出参数，指向设备内存的指针。调用成功后，devPtr将被设置为新分配内存区域的首地址，该内存区域可供GPU和CPU访问。

  2. size_t size：输入参数，要分配的内存大小（以字节为单位）。

  3. unsigned int flags：可选输入参数，默认值为0。可以使用此参数指定内存管理行为，例如：

    • cudaMemAttachHost: 表示内存应该同时对主机和设备可见，并且在需要时自动同步。
    • cudaMemAttachGlobal: 与cudaMemAttachHost类似，但不保证内存始终驻留在主机上。

返回值：

  • cudaSuccess：如果成功分配了内存，则返回此状态。
  • 其他错误代码（如cudaErrorMemoryAllocation等）：如果发生错误，如内存不足或其他CUDA运行时错误，则返回相应的错误代码。

通过cudaMallocManaged分配的内存会遵循统一内存（Unified Memory）模型，使得开发者无需手动管理数据在主机和设备之间的迁移，大大简化了并行计算中的内存管理复杂度。然而，需要注意的是，尽管统一内存提供了便捷性，但在某些情况下可能会影响性能，因为它依赖于CUDA驱动程序和硬件进行透明的数据迁移。

利用统一内存，可以将前面的程序简化如下：

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
#include <stdio.h>
#include <math.h>
#include <cuda_runtime.h>

__global__ void add(float* x, float* y, float* z, int n) {
    int index = threadIdx.x + blockIdx.x * blockDim.x;
    int stride = blockDim.x * gridDim.x;

    for (int i = index; i < n; i += stride) {
        z[i] = x[i] + y[i];
    }
}

int main() {
    int n = 1 << 20; // 向量大小
    int nBytes = n * sizeof(float);

    float *z; // 定义统一内存中的向量指针

    // 在统一内存中分配内存
    cudaMallocManaged((void**)&z, nBytes * 3); // 分配统一内存大小为三倍，包括 x, y, z

    float *x = z; // x指针指向统一内存
    float *y = x + n; // y指针指向统一内存中x的末尾
    float *dev_z = y + n; // dev_z指针指向统一内存中y的末尾

    // 初始化主机上的向量数据
    for (int i = 0; i < n; i++) {
        x[i] = 1.0f;
        y[i] = 2.0f;
    }

    // 定义 grid 和 block 的大小
    int blockSize = 256;
    int numBlocks = (n + blockSize - 1) / blockSize;

    // 调用 kernel 函数
    add<<<numBlocks, blockSize>>>(x, y, dev_z, n);

    // 同步device 保证结果能正确访问
    cudaDeviceSynchronize();

    // 验证结果
    float maxError = 0.0;
    for (int i = 0; i < n; i++) {
        maxError = fmaxf(maxError, fabsf(dev_z[i] - 3.0f));
    }
    printf("最大误差: %f\n", maxError);

    // 释放统一内存
    cudaFree(z);

    return 0;
}

相比之前的代码，使用统一内存更简洁了，值得注意的是kernel执行是与host异步的，由于托管内存自动进行数据传输，这里要用cudaDeviceSynchronize()函数保证device和host同步，这样后面才可以正确访问kernel计算的结果。

矩阵乘法实例

最后再实现一个稍微复杂一些的例子，就是两个矩阵的乘法，设输入矩阵为 A 和 B ，要得到 C=A×B 。实现思路是每个线程计算 C的一个元素值 Ci,j ，对于矩阵运算，应该选用grid和block为2-D的。首先定义矩阵的结构体：

1
2
3
4
5
6
7
// 矩阵类型，行优先，M(row, col) = *(M.elements + row * M.width + col)
struct Matrix
{
    int width;
    int height;
    float *elements;
};
矩阵乘法

然后实现矩阵乘法的核函数，这里定义了两个辅助的__device__函数分别用于获取矩阵的元素值和为矩阵元素赋值，具体代码如下：

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
// 获取矩阵A的(row, col)元素
__device__ float getElement(const Matrix* A, int row, int col)
{
    return A->elements[row * A->width + col];
}

// 为矩阵A的(row, col)元素赋值
__device__ void setElement(Matrix* A, int row, int col, float value)
{
    A->elements[row * A->width + col] = value;
}

// 矩阵相乘kernel，2-D，每个线程计算一个元素
__global__ void matMulKernel(const Matrix* A, const Matrix* B, Matrix* C)
{
    float Cvalue = 0.0;
    int row = threadIdx.y + blockIdx.y * blockDim.y;
    int col = threadIdx.x + blockIdx.x * blockDim.x;
    for (int i = 0; i < A->width; ++i)
    {
        Cvalue += getElement(A, row, i) * getElement(B, i, col);
    }
    setElement(C, row, col, Cvalue);
}

最后采用统一内存编写矩阵相乘的测试实例：

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
int main()
{
    int width = 1 << 10;
    int height = 1 << 10;
    Matrix *A, *B, *C;
    // 申请托管内存
    cudaMallocManaged((void**)&A, sizeof(Matrix));
    cudaMallocManaged((void**)&B, sizeof(Matrix));
    cudaMallocManaged((void**)&C, sizeof(Matrix));
    int nBytes = width * height * sizeof(float);
    cudaMallocManaged((void**)&A->elements, nBytes);
    cudaMallocManaged((void**)&B->elements, nBytes);
    cudaMallocManaged((void**)&C->elements, nBytes);

    // 初始化数据
    A->height = height;
    A->width = width;
    B->height = height;
    B->width = width;
    C->height = height;
    C->width = width;
    for (int i = 0; i < width * height; ++i)
    {
        A->elements[i] = 1.0;
        B->elements[i] = 2.0;
    }

    // 定义kernel的执行配置
    dim3 blockSize(32, 32);
    dim3 gridSize((width + blockSize.x - 1) / blockSize.x, 
        (height + blockSize.y - 1) / blockSize.y);
    // 执行kernel
    matMulKernel << < gridSize, blockSize >> >(A, B, C);


    // 同步device 保证结果能正确访问
    cudaDeviceSynchronize();
    // 检查执行结果
    float maxError = 0.0;
    for (int i = 0; i < width * height; ++i)
        maxError = fmax(maxError, fabs(C->elements[i] - 2 * width));
    std::cout << "最大误差: " << maxError << std::endl;

    return 0;
}

这里矩阵大小为1024，设计的线程的block大小为(32, 32)，那么grid大小为(32, 32)，最终测试结果如下：

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
nvprof cuda9.exe
==16304== NVPROF is profiling process 16304, command: cuda9.exe
最大误差: 0
==16304== Profiling application: cuda9.exe
==16304== Profiling result:
            Type  Time(%)      Time     Calls       Avg       Min       Max  Name
 GPU activities:  100.00%  1.32752s         1  1.32752s  1.32752s  1.32752s  matMulKernel(Matrix*, Matrix*, Matrix*)
      API calls:   83.11%  1.32762s         1  1.32762s  1.32762s  1.32762s  cudaDeviceSynchronize
                   13.99%  223.40ms         6  37.233ms  37.341us  217.66ms  cudaMallocManaged
                    2.81%  44.810ms         1  44.810ms  44.810ms  44.810ms  cudaLaunch
                    0.08%  1.3300ms        94  14.149us       0ns  884.64us  cuDeviceGetAttribute
                    0.01%  199.03us         1  199.03us  199.03us  199.03us  cuDeviceGetName
                    0.00%  10.009us         1  10.009us  10.009us  10.009us  cuDeviceTotalMem
                    0.00%  6.5440us         1  6.5440us  6.5440us  6.5440us  cudaConfigureCall
                    0.00%  3.0800us         3  1.0260us     385ns  1.5400us  cudaSetupArgument
                    0.00%  2.6940us         3     898ns     385ns  1.5390us  cuDeviceGetCount
                    0.00%  1.9250us         2     962ns     385ns  1.5400us  cuDeviceGet

==16304== Unified Memory profiling result:
Device "GeForce GT 730 (0)"
   Count  Avg Size  Min Size  Max Size  Total Size  Total Time  Name
    2051  4.0000KB  4.0000KB  4.0000KB  8.011719MB  21.20721ms  Host To Device
     270  45.570KB  4.0000KB  1.0000MB  12.01563MB  7.032508ms  Device To Host

小结

CUDA编程入门相对简单，通过掌握基本的线程管理和内存操作，可以快速实现并行计算的基本应用。然而，随着问题复杂度的提升，深入CUDA领域则需面对更多挑战：例如动态并行性管理、流式多处理器的细粒度同步控制、全局内存与共享内存的有效利用（包括bank冲突处理）、纹理和常量内存优化、以及先进的并行算法设计如归约、扫瞄等。此外，针对现代GPU架构的特性进行性能调优，理解并解决内存层次结构带来的瓶颈，亦是深入CUDA技术栈的重要环节。只有不断钻研这些高级技术和最佳实践，才能在大规模并行计算中发挥出CUDA的真正威力，实现从“入门”到“精通”的跨越。

参考文章

CUDA编程入门极简教程 - 知乎 (zhihu.com)

https://blog.csdn.net/TracelessLe/article/details/110880135

文章作者: loyeh
文章链接: https://Lo-Y-eH.github.io/posts/c7d8/
版权声明: 本博客所有文章除特別声明外，均采用 CC BY 4.0 许可协议。转载请注明来源 loyeh !
cuda
赏
评论
  目录