nv-cuda-prll-prog-merge-0

Nvidia CUDA 并行编程笔记（一）

001：NVIDIA CUDA中的HelloWorld程序写法 🚀

在本节课中，我们将学习如何编写和运行CUDA程序。我们将从最简单的“Hello World”程序开始，逐步探索其不同变体，并了解如何执行多线程程序以及如何在GPU上分配内存。

概述

本节课是CUDA编程的入门。我们将首先了解一个基本的CUDA程序结构，然后分析其执行流程。接着，我们会探讨如何通过多线程来并行执行任务，并学习如何管理CPU与GPU之间的异步执行和数据传输。课程的核心目标是理解CUDA内核（Kernel）的启动机制以及线程层次结构。

CPU上的Hello World程序

首先，我们来看一个用C语言编写的标准“Hello World”程序。这个程序与CUDA无关，完全在CPU上运行。

#include <stdio.h>

int main() {
    printf("Hello World\n");
    return 0;
}

使用GCC编译器编译并运行此程序，会在屏幕上输出“Hello World”。

gcc hello_world.c -o hello_world
./hello_world

第一个CUDA Hello World程序

现在，我们来看一个在GPU上运行的“Hello World”程序。它与CPU版本有三个关键区别。

#include <stdio.h>
#include <cuda_runtime.h> // 区别1：包含CUDA头文件

// 区别2：使用 __global__ 关键字声明GPU内核函数
__global__ void d_kernel() {
    printf("Hello World\n");
}

![](https://github.com/OpenDocCN/dsai-notes-pt1-zh/raw/master/docs/nv-cuda-prll-prog/img/84d8c0d6900c565dc28757b2f2d2785e_5.png)

int main() {
    // 区别3：使用特殊语法启动内核，<<<1, 1>>> 表示用1个线程块，每个块1个线程
    d_kernel<<<1, 1>>>();
    cudaDeviceSynchronize(); // 等待GPU内核执行完成
    return 0;
}

关键概念解析

  1. 头文件：#include <cuda_runtime.h> 是编写CUDA代码所必需的，它包含了CUDA相关的函数和类型定义。

  2. 内核函数：__global__ 关键字告诉编译器，d_kernel 是一个GPU内核函数，它将在GPU上执行，而不是在CPU上。

  3. 内核启动：d_kernel<<<1, 1>>>() 是启动内核的语法。<<<1, 1>>> 中的两个参数定义了线程的配置：

    • 第一个 1 表示线程块（Block）的数量。

    • 第二个 1 表示每个线程块中的线程（Thread）数量。

    • 因此，<<<1, 1>>> 表示启动一个包含1个线程的线程块，即总共只有1个线程执行该内核。

注意：如果省略 cudaDeviceSynchronize()，程序可能不会输出任何内容。这是因为CPU在启动内核后不会等待GPU完成，而是继续执行并可能提前结束程序。cudaDeviceSynchronize() 的作用是让CPU等待所有GPU线程执行完毕。

多个内核的调用顺序

上一节我们介绍了单个内核的启动。本节中我们来看看当程序中有多个内核调用时，它们的执行顺序是怎样的。

__global__ void d_kernel() {
    printf("Hello World\n");
}

int main() {
    d_kernel<<<1, 1>>>();
    d_kernel<<<1, 1>>>();
    d_kernel<<<1, 1>>>();
    cudaDeviceSynchronize();
    printf("on CPU\n");
    return 0;
}

执行流程分析

  1. CPU按顺序将三个内核调用放入一个默认的队列（称为“流”，Stream）中。

  2. GPU从这个队列中按先进先出（FIFO）的顺序取出并执行内核。

  3. cudaDeviceSynchronize() 确保CPU等待队列中所有内核执行完毕。

  4. 最后，CPU执行 printf(“on CPU\n”)。

因此，输出结果是确定性的：三行“Hello World”，然后是一行“on CPU”。内核的执行是串行的（一个接一个），但这是GPU内部的串行，CPU在启动它们后就去等待了。

CPU与GPU的异步执行

我们调整一下代码顺序，观察CPU和GPU异步执行的影响。

int main() {
    d_kernel<<<1, 1>>>();
    d_kernel<<<1, 1>>>();
    printf("on CPU\n");
    cudaDeviceSynchronize(); // 注意：同步调用放到了打印语句之后
    return 0;
}

可能的输出结果

由于CPU和GPU是异步执行的，输出顺序变得不确定。以下是几种可能的情况：

  1. on CPU -> Hello World -> Hello World

  2. Hello World -> on CPU -> Hello World

  3. Hello World -> Hello World -> on CPU

关键点：

  • CPU启动内核后立即继续执行后面的 printf 语句。

  • 两个内核在默认流中仍按顺序执行（Hello World 总是先于第二个 Hello World 打印）。

  • cudaDeviceSynchronize() 保证了两个内核最终都会执行完毕，但无法控制它相对于CPU printf 的执行时机。

混合CPU与GPU打印任务

让我们看一个更复杂的例子，其中混合了CPU和GPU的打印任务。

int main() {
    d_kernel<<<1, 1>>>();
    printf("CPU1\n");
    d_kernel<<<1, 1>>>();
    printf("CPU2\n");
    d_kernel<<<1, 1>>>();
    printf("CPU3\n");
    cudaDeviceSynchronize();
    printf("on CPU\n");
    return 0;
}

并行性分析

在这个例子中，以下操作可以并行执行：

  • GPU执行三个 d_kernel 打印“Hello World”。

  • CPU执行三个 printf 打印“CPU1”、“CPU2”、“CPU3”。

因此，你可能会看到“CPU1/2/3”与“Hello World”交错打印的输出。然而，由于 cudaDeviceSynchronize() 的存在，最后的 ”on CPU” 一定会在所有GPU内核执行完毕后才打印。

一个常见的输出顺序是：CPU1 CPU2 CPU3 -> 三个 Hello World -> on CPU。这是因为CPU执行速度通常很快，在GPU完成大量工作前就可能完成了自己的打印任务。

使用多线程并行打印

前面的例子每个内核只使用一个线程。CUDA的强大之处在于大规模并行。以下是如何用多个线程打印“Hello World”。

__global__ void d_kernel() {
    printf("Hello World\n");
}

int main() {
    // 启动1个线程块，该块中包含32个线程
    d_kernel<<<1, 32>>>();
    cudaDeviceSynchronize();
    return 0;
}

编译并运行此程序，将会看到“Hello World”被打印了32次，每个线程执行一次 printf 语句。

线程配置公式：<<<NumBlocks, ThreadsPerBlock>>>

  • NumBlocks：线程网格（Grid）中线程块的数量。

  • ThreadsPerBlock：每个线程块中的线程数量。

  • 总线程数 = NumBlocks * ThreadsPerBlock

并行计算示例：打印平方数

上一节我们介绍了如何用多线程执行相同任务。本节中我们来看看如何让每个线程处理不同的数据。假设我们想并行计算并打印0到99的平方。

思路分析

目标是让线程i计算并打印 i * i。我们需要一种方法让每个线程知道自己独特的索引（ID）。

__global__ void square_kernel() {
    int tid = threadIdx.x; // 获取当前线程在线程块内的索引（0到ThreadsPerBlock-1）
    printf(“Square of %d is %d\n”, tid, tid * tid);
}

![](https://github.com/OpenDocCN/dsai-notes-pt1-zh/raw/master/docs/nv-cuda-prll-prog/img/84d8c0d6900c565dc28757b2f2d2785e_45.png)

int main() {
    int n = 100;
    // 启动1个线程块，包含100个线程
    square_kernel<<<1, n>>>();
    cudaDeviceSynchronize();
    printf(“on CPU\n”);
    return 0;
}

代码解释：

  • threadIdx.x 是一个内置变量，表示线程在其所属线程块中的x维度索引。

  • 由于我们只使用了一个线程块（<<<1, n>>>），threadIdx.x 的范围是0到n-1。

  • 每个线程根据自己独有的 tid 计算平方并打印。

注意：由于所有线程是并行执行的，打印输出的顺序可能是混乱的（即非从0到99的顺序）。这是并行编程中常见的现象，如果需要顺序输出，则需要额外的同步或排序操作。

核心概念回顾

在深入代码之前，我们先简要回顾一些核心的计算机系统概念，这有助于理解CUDA的线程模型。

  1. 进程（Process）：正在执行的程序实例。它拥有独立的地址空间，包含代码、数据、堆栈等资源。例如，打开的Word应用程序就是一个进程。

  2. 线程（Thread）：进程内的一个执行流。它是“轻量级进程”，与同进程的其他线程共享代码和全局数据，但拥有独立的栈和寄存器。例如，Word中的拼写检查功能可以是一个独立的线程。

  3. 操作系统（OS）：管理硬件资源和为用户程序提供服务的软件。例如，Linux， Windows。

  4. 硬件组件：

    • 缓存（Cache）：靠近处理器的小容量高速内存，用于加速数据访问。

    • 主存（RAM）：容量更大的内存，用于存储正在运行的程序和数据。

    • 核心（Core）：执行线程的物理处理单元。一个核心同一时刻只能执行一个线程。线程可能会在不同的核心之间调度（“跳跃”）。

在CUDA中，我们创建成千上万个线程，它们被组织成线程块和网格，在GPU的众多核心上并行执行。

总结

本节课我们一起学习了CUDA并行编程的基础。我们从最简单的“Hello World”程序出发，理解了CUDA内核函数（用 __global__ 修饰）的概念以及如何用 <<<…>>> 语法启动内核。我们探讨了CPU与GPU的异步执行特性，并学习了使用 cudaDeviceSynchronize() 进行同步。通过示例，我们看到了如何配置多线程（<<<NumBlocks, ThreadsPerBlock>>>）来并行执行任务，以及如何利用 threadIdx.x 让每个线程处理不同的数据。这些概念是构建更复杂、高性能CUDA应用程序的基石。

002：CUDA程序流程与cudaMemcpy操作 🚀

概述

在本节课中，我们将学习CUDA程序的基本流程，特别是CPU与GPU之间如何通过cudaMemcpy函数进行数据传输。我们将通过具体的代码示例，理解如何为GPU分配内存、启动内核以及将结果复制回CPU。

回顾与引入

上一节我们介绍了如何编写和启动一个简单的CUDA内核，并使用多个线程并行打印数字的平方。本节中，我们来看看如何将计算结果存储到数组中，并将其从GPU传回CPU。

我们将从一个简单的C语言串行代码开始，目标是将其并行化。

// 串行C代码示例
for (int i = 0; i < 100; i++) {
    a[i] = i * i;
}

以下是实现相同功能的CUDA并行代码。

#include <stdio.h>

![](https://github.com/OpenDocCN/dsai-notes-pt1-zh/raw/master/docs/nv-cuda-prll-prog/img/79f6aa2eef992bc8530785c087c51019_20.png)

__global__ void square(int *a) {
    int tid = threadIdx.x;
    a[tid] = tid * tid;
}

![](https://github.com/OpenDocCN/dsai-notes-pt1-zh/raw/master/docs/nv-cuda-prll-prog/img/79f6aa2eef992bc8530785c087c51019_22.png)

int main() {
    int n = 100;
    int a[n];
    int *d_a;
    cudaMalloc((void**)&d_a, n * sizeof(int));
    square<<<1, n>>>(d_a);
    cudaMemcpy(a, d_a, n * sizeof(int), cudaMemcpyDeviceToHost);
    for (int i = 0; i < n; i++) {
        printf("%d ", a[i]);
    }
    cudaFree(d_a);
    return 0;
}

代码解析

以下是上述代码的关键步骤解析：

  1. 主机（CPU）内存分配：在CPU上声明数组 a。

  2. 设备（GPU）内存分配：使用 cudaMalloc 在GPU上分配内存，d_a 是一个位于CPU上的指针，但它指向GPU上的内存地址。

  3. 内核启动：使用 <<<1, n>>> 语法启动 square 内核，其中 n 是线程数。每个线程计算其 threadIdx.x 的平方，并写入GPU数组 d_a 的对应位置。

  4. 数据回传：使用 cudaMemcpy 将计算结果从GPU内存 (d_a) 复制回CPU内存 (a)。第四个参数 cudaMemcpyDeviceToHost 指明了传输方向。

  5. 内存释放：使用 cudaFree 释放GPU上分配的内存。

核心概念与注意事项

在深入更多示例前，我们需要理解几个核心概念：

  • CPU与GPU内存分离：CPU和GPU拥有独立的内存空间（显存）。在离散GPU上，CPU无法直接访问GPU内存，反之亦然。因此，必须显式地进行内存分配和数据拷贝。

  • cudaMemcpy 是阻塞操作：cudaMemcpy 调用会阻塞CPU线程，直到所有先前启动的内核执行完毕且数据传输完成。因此，在这个例子中，我们不需要显式调用 cudaDeviceSynchronize()。

  • PCIe总线瓶颈：CPU与GPU之间的数据传输通过PCI Express总线进行，其带宽有限。频繁或大量数据传输可能成为程序性能瓶颈。因此，应确保在GPU上进行足够的计算，以抵消数据传输的开销。

CPU/GPU内存访问示例

为了加深对内存分离的理解，我们来看几个例子。

示例1：错误的全局变量访问

以下代码尝试在内核中访问一个在CPU上定义的全局变量字符串指针，这将导致编译错误。

char *msg = "Hello World"; // 在CPU内存中

![](https://github.com/OpenDocCN/dsai-notes-pt1-zh/raw/master/docs/nv-cuda-prll-prog/img/79f6aa2eef992bc8530785c087c51019_47.png)

![](https://github.com/OpenDocCN/dsai-notes-pt1-zh/raw/master/docs/nv-cuda-prll-prog/img/79f6aa2eef992bc8530785c087c51019_49.png)

__global__ void hello() {
    printf("%s\n", msg); // 错误！GPU内核无法直接访问CPU内存
}

![](https://github.com/OpenDocCN/dsai-notes-pt1-zh/raw/master/docs/nv-cuda-prll-prog/img/79f6aa2eef992bc8530785c087c51019_51.png)

int main() {
    hello<<<1, 32>>>();
    cudaDeviceSynchronize();
    return 0;
}

关键点：在离散GPU上，定义在主机（CPU）上的变量不能被设备（GPU）代码直接访问。

示例2：使用预处理器宏

如果使用预处理器宏 #define，代码则可以正常工作。

#define MSG "Hello World" // 预处理器指令，在编译前进行文本替换

![](https://github.com/OpenDocCN/dsai-notes-pt1-zh/raw/master/docs/nv-cuda-prll-prog/img/79f6aa2eef992bc8530785c087c51019_63.png)

![](https://github.com/OpenDocCN/dsai-notes-pt1-zh/raw/master/docs/nv-cuda-prll-prog/img/79f6aa2eef992bc8530785c087c51019_65.png)

__global__ void hello() {
    printf("%s\n", MSG); // 正确！编译前MSG已被替换为字符串字面量
}

int main() {
    hello<<<1, 32>>>();
    cudaDeviceSynchronize();
    return 0;
}

关键点：#define 是预处理器指令，它在编译之前将 MSG 直接替换为 "Hello World" 字符串字面量。对于内核来说，它只是在打印一个常量字符串，而非访问CPU内存地址。

典型的CUDA程序流程 📊

一个标准的CUDA程序通常遵循以下步骤，我们将其与所需的CUDA API对应起来：

  1. 从文件系统加载数据到CPU内存：使用标准C/C++文件操作（如fread）。

  2. 将数据从CPU内存传输到GPU内存：

    • 使用 cudaMalloc 在GPU上分配内存。

    • 使用 cudaMemcpy(..., cudaMemcpyHostToDevice) 将数据从主机复制到设备。

  3. 在GPU上执行内核：使用 kernel_name<<<grid_size, block_size>>>(parameters) 启动内核。

  4. 将结果从GPU内存传输回CPU内存：使用 cudaMemcpy(..., cudaMemcpyDeviceToHost)。

  5. 在CPU上使用或输出结果：使用 printf 打印或 fwrite 写回文件。

重要提示：程序员有责任维护CPU和GPU上数据的一致性。通常，建议使用不同的变量名来区分它们（例如，h_a 表示主机数组，d_a 表示设备数组），以避免错误传递指针。

实践练习：字符串操作与边界问题

让我们通过一个操作字符串的练习，来学习如何正确处理GPU内存和线程索引。

任务描述

我们有一个字符串 "GdkknVnqkc"，每个字符在字母表中都是前一个字符（例如，‘H’的前一个字符是‘G’）。我们想用GPU线程并行地将每个字符递增，最终得到 "HelloWorld"。

有问题的初始代码

以下代码存在一个常见的边界错误。

__global__ void decode(char *a, int len) {
    int tid = threadIdx.x;
    if (tid < len) {
        a[tid]++; // 所有线程都对字符进行加一操作
    }
}

int main() {
    char cpu_a[] = "GdkknVnqkc";
    int len = strlen(cpu_a);
    char *gpu_a;
    cudaMalloc((void**)&gpu_a, (len + 1) * sizeof(char)); // 分配len+1字节，包含'\0'
    cudaMemcpy(gpu_a, cpu_a, (len + 1) * sizeof(char), cudaMemcpyHostToDevice);
    decode<<<1, len + 1>>>(gpu_a, len + 1); // 启动了 len+1 个线程
    cudaMemcpy(cpu_a, gpu_a, (len + 1) * sizeof(char), cudaMemcpyDeviceToHost);
    printf("%s\n", cpu_a);
    cudaFree(gpu_a);
    return 0;
}

问题分析：字符串 "GdkknVnqkc" 长度为10，但C风格字符串以 \0 结尾，所以 cpu_a 实际占用11字节。代码分配了11字节，并启动了11个线程。第11个线程（tid=10）操作的是字符串终止符 \0，将其递增为 \1，导致字符串失去终止符。后续 printf 可能打印乱码或导致段错误。

正确的代码

正确的做法是只对有效字符进行操作，不改变终止符。

decode<<<1, len>>>(gpu_a, len); // 只启动 len 个线程，不处理终止符

关键点：在并行处理数组（尤其是字符串）时，必须仔细管理线程索引与数据边界，确保不会意外修改或访问越界。

线程组织与限制 ⚙️

最后，我们简要探讨一下CUDA的线程组织。目前我们只使用了 <<<1, n>>> 中的第二个参数，它定义了一个线程块（Thread Block）中的线程数量。

  • 线程块大小限制：每个线程块的最大线程数通常是1024（取决于GPU架构）。这意味着，使用 <<<1, n>>> 这种格式，n 不能超过1024。

  • 如何启动更多线程：为了启动超过1024个线程（例如，处理8000个元素的数组），我们需要使用第一个参数来指定网格（Grid）中线程块的数量。这是下节课的重点。

例如，以下调用试图启动8000个线程，但可能会失败，因为单个线程块无法容纳8000个线程。

kernel<<<1, 8000>>>(); // 可能失败，超出线程块最大限制

我们将学习如何使用网格和线程块来灵活地组织大量线程。

总结

本节课中我们一起学习了：

  1. CUDA程序基本流程：包括主机与设备内存分配、数据传输和内核执行。

  2. cudaMemcpy 函数：用于在CPU和GPU之间复制数据，方向由 cudaMemcpyHostToDevice 和 cudaMemcpyDeviceToHost 控制。

  3. CPU与GPU内存分离：这是CUDA编程的核心概念，必须为需在GPU上处理的数据在设备端创建副本。

  4. 数据边界管理：在线程化代码中，必须谨慎处理数组索引，防止越界访问。

  5. 线程组织初步：了解了线程块的概念及其线程数量限制，为学习网格化线程组织打下基础。

掌握这些基础是编写正确、高效CUDA程序的关键。下一节，我们将深入探讨CUDA的线程层次结构（网格和线程块），以支持大规模并行计算。

003：CUDA中的线程组织与线程块 🧵

概述

在本节课中，我们将深入学习CUDA中的线程组织模型。我们将探讨线程网格、线程块和线程的层次结构，理解如何通过多维索引来唯一标识每个线程，并通过实例学习如何为一维和二维数据结构编写并行初始化代码。

线程组织的层次结构

在上一讲中，我们学习了典型的CUDA程序流程。本节中，我们来看看CUDA如何组织海量的线程来执行并行计算。

当一个内核被启动时，它会以一个线程网格的形式执行。这个网格是一个三维数据结构，可以将其想象成一个魔方。整个魔方就是网格，而构成魔方的每一个小立方体就是一个线程块。

  • 网格：是内核启动的最高层次组织，由线程块构成。

  • 线程块：是网格的组成部分，内部包含线程。

  • 线程：是实际执行计算的最小单元，位于线程块内部。

网格和线程块都可以是一维、二维或三维的，这为处理不同维度的数据（如图像、矩阵、体积数据）提供了灵活性。

标识线程与线程块

为了在并行计算中准确定位每个线程，我们需要一套坐标系统。

网格维度与块索引

网格的尺寸（即每个维度上有多少个线程块）可以通过内置变量访问：

  • gridDim.x：网格在x方向的尺寸（线程块数量）。

  • gridDim.y：网格在y方向的尺寸。

  • gridDim.z：网格在z方向的尺寸。

网格中每个线程块的位置由块索引唯一确定：

  • blockIdx.x：当前线程块在网格x方向上的索引。

  • blockIdx.y：当前线程块在网格y方向上的索引。

  • blockIdx.z：当前线程块在网格z方向上的索引。

索引范围：如果gridDim.x = 10，那么blockIdx.x的取值范围是 0 到 9。

线程块维度与线程索引

每个线程块的尺寸（即每个维度上有多少个线程）可以通过内置变量访问：

  • blockDim.x：线程块在x方向的尺寸（线程数量）。

  • blockDim.y：线程块在y方向的尺寸。

  • blockDim.z：线程块在z方向的尺寸。

线程块中每个线程的位置由线程索引唯一确定：

  • threadIdx.x：当前线程在线程块x方向上的索引。

  • threadIdx.y：当前线程在线程块y方向上的索引。

  • threadIdx.z：当前线程在线程块z方向上的索引。

索引范围：如果blockDim.x = 5，那么threadIdx.x的取值范围是 0 到 4。

内核启动配置示例

以下代码展示了如何配置并启动一个多维的内核：

dim3 grid(2, 3, 4); // 定义一个2x3x4的网格（共24个线程块）
dim3 block(5, 6, 7); // 定义一个5x6x7的线程块（共210个线程）

![](https://github.com/OpenDocCN/dsai-notes-pt1-zh/raw/master/docs/nv-cuda-prll-prog/img/68e31b967d69cfea77cc142db62cac46_44.png)

// 启动内核，共创建 24 * 210 = 5040 个线程
myKernel<<<grid, block>>>();
cudaDeviceSynchronize();

在内核函数myKernel内部，我们可以使用上述变量。例如，以下条件只会在一个特定的线程（位于第一个线程块的第一个线程）中为真：

if (blockIdx.x == 0 && blockIdx.y == 0 && blockIdx.z == 0 &&
    threadIdx.x == 0 && threadIdx.y == 0 && threadIdx.z == 0) {
    // 只有这个线程会执行这里的代码
}

应用实例：计算全局唯一线程ID

在实际编程中，我们经常需要为每个线程计算一个在整个网格范围内唯一的ID，以便处理像数组这样的线性内存。

情况一：单线程块，二维线程组织

假设我们有一个大小为 N * M 的一维数组，我们使用一个二维线程块（N x M 个线程）来初始化它。

内核启动配置：

#define N 5
#define M 6
dim3 block(N, M); // 5x6的线程块
myKernel<<<1, block>>>(matrix); // 网格中只有1个线程块

在内核中，计算唯一ID的公式为：

int i = threadIdx.x * blockDim.y + threadIdx.y;
matrix[i] = i; // 为每个位置赋予唯一值

公式解释：threadIdx.x 是行索引，blockDim.y（固定为6）是总列数，threadIdx.y 是列索引。该公式将二维索引映射到一维数组。

情况二：多线程块，一维线程组织

同样初始化一个 N * M 的数组，现在我们使用多个一维线程块。每个线程块包含 M 个线程，总共启动 N 个线程块。

内核启动配置：

#define N 5
#define M 6
myKernel<<<N, M>>>(matrix); // N个块，每个块M个线程

在内核中，计算唯一ID的公式变为：

int i = blockIdx.x * blockDim.x + threadIdx.x;
matrix[i] = i;

公式解释：blockIdx.x 是线程块索引，blockDim.x（固定为6）是每个块的线程数，threadIdx.x 是块内线程索引。该公式结合了块索引和块内线程索引来生成全局唯一ID。

这两种方法都能成功地将数组初始化为 [0, 1, 2, ..., 29]，展示了组织线程的灵活性。

总结

本节课我们一起学习了CUDA中线程组织的核心概念：

  1. CUDA执行模型是分层的：网格(Grid) -> 线程块(Block) -> 线程(Thread)。

  2. 网格和线程块都可以配置为一维、二维或三维结构，以适应不同维度的数据。

  3. 使用 blockIdx、threadIdx、gridDim、blockDim 这些内置变量可以唯一标识每个线程并获取执行配置信息。

  4. 关键在于计算全局唯一的线程ID，其公式取决于网格和线程块的维度配置。我们通过两个初始化数组的实例演示了如何推导和应用这个公式。

理解线程组织是编写高效CUDA并行程序的基础，它直接决定了任务如何被分解并映射到GPU的数千个核心上执行。

004：GPU计算层次与线程发散

在本节课中，我们将学习如何为大规模数据启动CUDA内核，深入理解GPU的计算层次结构，并探讨线程发散问题及其对性能的影响。

概述

上一节我们介绍了线程在GPU上的组织方式，包括网格、线程块以及线程ID的计算。本节中，我们将看看如何为大规模数据启动内核，理解GPU的底层执行单元——线程束，并分析条件分支如何导致线程发散，从而影响性能。

为大规模数据启动内核

当数据规模巨大时，我们需要计算合适的线程块数量来启动内核。以下是一个示例，展示了如何根据用户指定的数据大小动态计算所需的线程块数量。

#include <stdio.h>
#include <cuda.h>

#define BLOCK_SIZE 1024

__global__ void initVectorKernel(unsigned int *vector, unsigned int vector_size) {
    unsigned int id = blockIdx.x * blockDim.x + threadIdx.x;
    if (id < vector_size) {
        vector[id] = id;
    }
}

int main(int argc, char *argv[]) {
    unsigned int n = atoi(argv[1]);
    unsigned int *d_vector, *h_vector;

    cudaMalloc(&d_vector, n * sizeof(unsigned int));
    h_vector = (unsigned int*)malloc(n * sizeof(unsigned int));

    unsigned int num_blocks = (n + BLOCK_SIZE - 1) / BLOCK_SIZE;
    printf("Number of blocks: %u\n", num_blocks);

    initVectorKernel<<<num_blocks, BLOCK_SIZE>>>(d_vector, n);
    cudaMemcpy(h_vector, d_vector, n * sizeof(unsigned int), cudaMemcpyDeviceToHost);

    for (int i = 0; i < n; i++) {
        printf("%u ", h_vector[i]);
    }
    printf("\n");

    free(h_vector);
    cudaFree(d_vector);
    return 0;
}

以下是代码中需要解决的两个关键问题及其解决方案：

  1. 整数除法问题：使用 (n + BLOCK_SIZE - 1) / BLOCK_SIZE 代替 n / BLOCK_SIZE，确保当 n 不能被 BLOCK_SIZE 整除时，能分配足够的线程块。

  2. 内存越界访问：在内核函数中添加条件判断 if (id < vector_size)，确保只有有效的线程ID才会执行内存写入操作，防止访问非法内存地址。

矩阵平方应用：从CPU到GPU的并行化

我们以一个矩阵平方运算为例，展示如何逐步将其从CPU版本并行化到GPU上，并分析性能变化。

CPU版本

CPU版本的矩阵平方使用三层嵌套循环，时间复杂度为 O(n³)。

void matrixSquareCPU(float *matrix, float *result, int n) {
    for (int i = 0; i < n; i++) {
        for (int j = 0; j < n; j++) {
            float sum = 0.0f;
            for (int k = 0; k < n; k++) {
                sum += matrix[i * n + k] * matrix[k * n + j];
            }
            result[i * n + j] = sum;
        }
    }
}

对于一个64x64的矩阵，CPU执行时间约为 1.52毫秒。

GPU版本1：仅并行化外层循环

在第一个GPU版本中，我们仅并行化了最外层的 i 循环。每个线程负责计算结果矩阵的一整行。

__global__ void matrixSquareGPUv1(float *matrix, float *result, int n) {
    int i = blockIdx.x * blockDim.x + threadIdx.x;
    if (i < n) {
        for (int j = 0; j < n; j++) {
            float sum = 0.0f;
            for (int k = 0; k < n; k++) {
                sum += matrix[i * n + k] * matrix[k * n + j];
            }
            result[i * n + j] = sum;
        }
    }
}
// 内核启动: <<<1, n>>>

这个版本存在两个问题：

  1. 每个线程内部仍有大量串行工作（j 和 k 循环）。

  2. 并行度仅为 n（例如64），远未充分利用GPU的数千个核心。

因此，其执行时间（6.39毫秒）甚至比CPU版本更慢。

GPU版本2：并行化两个外层循环

为了提升并行度，我们同时并行化 i 和 j 循环。每个线程负责计算结果矩阵中的一个单独元素。

__global__ void matrixSquareGPUv2(float *matrix, float *result, int n) {
    int id = blockIdx.x * blockDim.x + threadIdx.x;
    int i = id / n;
    int j = id % n;
    if (i < n && j < n) {
        float sum = 0.0f;
        for (int k = 0; k < n; k++) {
            sum += matrix[i * n + k] * matrix[k * n + j];
        }
        result[i * n + j] = sum;
    }
}
// 内核启动: <<<(n*n + BLOCK_SIZE -1)/BLOCK_SIZE, BLOCK_SIZE>>>

这个版本取得了显著改进：

  • 并行度：提升到 n * n（例如4096）。

  • 串行工作：每个线程仅剩一个 k 循环。

其执行时间大幅降低至 0.1毫秒，优于CPU版本。

注意：最内层的 k 循环涉及对同一内存位置（sum）的累加，如果直接并行化会导致数据竞争，需要同步机制，这将在后续课程中讨论。

GPU计算层次结构

理解GPU的硬件执行模型对编写高效CUDA程序至关重要。其计算层次从高到低如下：

  1. 线程：最小的执行单元，类似于一个“工人”。

  2. 线程束：由 32个线程 组成的基本执行单元。线程束内的所有线程以 SIMD 方式执行，即单指令，多数据。它们步调一致地执行相同的指令，但操作不同的数据。

  3. 线程块：由多个线程束组成。一个线程块最多可包含1024个线程。

    • 块内同步：线程块内的线程可以通过 __syncthreads() 函数进行同步。

    • 共享内存：线程块拥有一个高速的、块内线程可共享的片上内存（Shared Memory），能显著减少全局内存访问延迟。

  4. 流式多处理器：GPU上的一个处理单元，可同时管理多个线程块（通常1-8个）。一个SM包含大量核心（例如128个）。

  5. GPU设备：由多个SM组成，可同时运行数万甚至数十万个线程。

关键点：锁步执行仅发生在线程束内部。不同线程束之间的执行是异步的。

线程发散

线程发散是GPU编程中一个重要的性能考量因素。它发生在一个线程束内的线程需要执行不同的指令路径时。

发散的原因与影响

由于线程束以SIMD方式执行，当遇到条件分支（如 if-else）时，如果线程束内有些线程满足条件，有些不满足，GPU会先执行满足条件分支的线程，其他线程等待（执行空操作）；然后再执行另一分支的线程。这导致了串行化执行，降低了效率。

__global__ void divergentKernel(int *data) {
    int id = threadIdx.x;
    if (id % 2 == 0) {
        data[id] *= 2; // 偶数线程执行
    } else {
        data[id] += 1; // 奇数线程执行
    }
}

在上面的内核中，一个线程束内的奇偶线程将分两批执行，造成发散。

何时会发生发散？

导致发散的关键不是条件语句本身，而是条件在同一个线程束内是否对所有线程评估出相同的结果。

  • 不会发散的例子：if (threadIdx.x < 32)。对于一个线程束（0-31），条件对所有线程都为真，因此没有发散。

  • 会发散的例子：循环次数依赖于线程ID，for(int i=0; i<threadIdx.x; i++)。线程束内不同线程的循环次数不同，最后一个线程（迭代31次）将决定整个线程束的执行时间，导致其他线程等待，引发发散。

避免发散的一个技巧

在某些特定情况下，可以重写代码来消除条件分支。

原始代码（可能发散）：

if (x == y) {
    x = z;
} else if (x == z) {
    x = y;
}

优化后代码（无分支）：

x = (x == y) * z + (x == z) * y;
// 或者更通用的： x = z * (x == y) + y * (x == z);

优化后的代码通过逻辑运算和算术运算替代了条件判断，使得线程束内所有线程执行相同的指令流，从而避免了发散。但这种方法并非在所有情况下都适用。

总结

本节课中我们一起学习了：

  1. 如何为大规模数据计算并启动合适数量的线程块，并注意处理整数除法和内存越界问题。

  2. 通过矩阵平方的例子，实践了将CPU代码并行化到GPU的逐步优化过程，认识到提高并行度和减少线程内串行工作的重要性。

  3. 深入理解了GPU的计算层次结构，特别是线程束的SIMD执行模式。

  4. 分析了线程发散问题的成因、影响以及避免方法。核心在于意识到：线程束内执行路径的分化是性能杀手，编写代码时应尽量让同一个线程束内的线程执行相同的指令流。

掌握这些概念是编写高效CUDA程序的基础。下一节课，我们将探讨GPU的内存模型。

005：CUDA中的线程发散（续）与内存模型概述

在本节课中，我们将深入学习CUDA编程中的线程发散问题，并通过具体示例探讨如何优化代码以减少发散。随后，我们将从计算部分转向内存部分，初步了解GPU的内存层次结构。

回顾与概述

上一节我们介绍了GPU的计算层次结构、典型的矩阵平方应用并行化方法，以及处理大规模数据时内核启动的常见问题。我们还深入探讨了线程束（Warp）的概念、执行方式（锁步、SIMD），并初步了解了线程发散（Thread Divergence）问题。

本节我们将继续深入线程发散问题，通过更多示例学习如何识别和缓解发散。最后，我们将开始新的主题——GPU内存，了解其基本组织方式。

线程发散深入分析

示例一：多分支条件语句

首先，我们回顾上节课结束时的一个简单示例。以下内核启动了一个包含32个线程（即一个线程束）的核函数：

__global__ void D_kernel(float* vector, int vector_size) {
    int id = threadIdx.x + blockIdx.x * blockDim.x;
    if (id < vector_size) {
        if (id % 2 == 1) {
            // 语句 S1 (奇数ID线程执行)
        } else {
            // 语句 S2 (偶数ID线程执行)
        }
        // 语句 S3 (所有线程执行)
    }
}

关键观察点：对于这个 if-else 结构，无论线程ID是奇数还是偶数，线程束执行它都需要 2步。第一步，一部分线程执行 S1，另一部分闲置；第二步，另一部分线程执行 S2，之前执行的线程闲置。虽然存在发散，但执行步数是固定的。

练习：复杂条件分支

现在考虑一个更复杂的条件分支：

if (x == y) {
    a = z + w;
} else if (x == z) {
    a = y + w;
} else {
    a = y + z;
}

问题：

  1. 线程束执行此 if-else-if 结构最多需要多少步？

  2. 如何消除此处的线程发散？

分析与解答：

  1. 执行步数：在最坏情况下，线程束中的线程可能分别进入三个不同的条件分支。因此，线程束需要按顺序让进入每个分支的线程组依次执行，总共需要 3步。

  2. 消除发散：我们可以利用数学等价替换，将条件分支转换为一个无分支的表达式，从而让所有线程执行相同的操作：

    a = (x == y) * (z + w) + (x == z) * (y + w) + (x == w) * (y + z);

    但更优雅的解法是利用这个公式：a = (y * (x == z) + z * (x == y) + w * (x == w)) / ((x == y) + (x == z) + (x == w));

    当x等于y、z或w时，此公式能正确计算出对应的a值，且所有线程执行相同的计算，只需1步，彻底消除了发散。

示例二：包含循环的昂贵分支

考虑以下内核，它同样由32个线程（一个线程束）执行：

__global__ void T_kernel(float* vector, int vector_size) {
    int id = threadIdx.x + blockIdx.x * blockDim.x;
    if (id <= 0) { // 仅线程ID为0满足条件
        vector[0] = 0;
        for(int i = 1; i <= 100; i++) {
            vector[0] += i; // 计算1到100的和
        }
    } else {
        vector[id] = 1;
    }
}

问题：线程束执行这个 if-else 结构需要多少步？

分析与解答：

  • 只有 id=0 的线程会进入 if 分支，执行 vector[0] = 0;（1步）和100次循环（100步），共101步。

  • 其余31个线程 (id=1 到 31) 在 if 阶段闲置。

  • 然后，这31个线程进入 else 分支，执行 vector[id] = 1;（1步），而 id=0 的线程在此阶段闲置。

  • 因此，线程束总共需要 101 + 1 = 102步 来执行此结构。性能瓶颈在于 id=0 线程的昂贵循环。

优化：我们可以将循环计算（求和）替换为闭合公式（等差数列求和公式）：

if (id <= 0) { vector[0] = 100 * 101 / 2; }

这样，if 分支的耗时从101步降为1步。优化后，线程束执行整个 if-else 结构只需要 2步（if分支1步，else分支1步）。

进一步优化以消除发散：虽然步数减为2，但线程依然发散（一部分执行if，另一部分执行else）。我们可以利用整数运算技巧，将其合并为单一步骤：

vector[id] = 1 + ( ( -(id >> 31) ) & (100*101/2 - 1) );

这个表达式利用了整数右移和位操作，使得当 id=0 时，结果为 100*101/2；当 id>0 时，结果为 1。所有线程执行相同的指令序列，仅需1步，完全消除了发散。

示例三：Switch-Case 语句

考虑以下使用 switch-case 的内核：

__global__ void D_kernel(float* vector, int vector_size) {
    int id = threadIdx.x + blockIdx.x * blockDim.x;
    switch(id) {
        case 0: vector[0] = 0; break;
        case 1: vector[1] = vector[0]; break;
        case 2: vector[2] = vector[id-2]; break;
        case 3: vector[3] = vector[6]; break;
        case 4: vector[4] = vector[5] + 4; break;
        case 5: vector[5] = vector[4] - 4; break;
        case 6: vector[6] = vector[3]; break;
        case 7: vector[7] = vector[7] * 2; break;
        case 8: vector[8] = vector[4] + id; break;
        case 9: vector[9] = 9; break;
    }
}

问题：线程束执行此 switch 结构需要多少步？

分析与解答：

switch 语句在逻辑上等价于一系列的 if-else if 判断。线程束需要让进入每个 case 的线程组依次执行。因为有10个 case，所以在最坏情况下需要 10步。

优化：通过观察，我们可以合并一些执行相同或等价操作的 case，从而减少步数：

  • case 0 和 case 2：由于线程锁步执行，当 id=2 执行 vector[2] = vector[0]; 时，vector[0] 已被 id=0 的线程初始化为0。因此两者效果可合并（将vector对应位置赋0）。

  • case 3 和 case 6：都执行 vector[3] = vector[6]; 和 vector[6] = vector[3];（具体值取决于执行顺序，但操作对称）。

  • case 4 和 case 8：case 8 执行 vector[8] = vector[4] + 8;，这与 case 4 的 vector[4] = vector[5] + 4; 在数据流上相关，在某些条件下可考虑合并逻辑。

通过合并，可以将执行步数从10步减少到7步或更少。

示例四：屏障同步与死锁

考虑以下包含全局屏障的内核：

__global__ void example_kernel() {
    int id = threadIdx.x + blockIdx.x * blockDim.x;
    if (id < 16) {
        printf("Inside if\n");
        __syncthreads(); // 全局屏障
    } else {
        printf("Inside else\n");
    }
}

问题：这个程序的输出是什么？会发生什么？

分析与解答：

  • 假设内核启动了一个包含32个线程的线程束。

  • 前16个线程 (id=0 到 15) 进入 if 分支，打印 “Inside if”，然后到达 __syncthreads() 屏障并等待。

  • 后16个线程 (id=16 到 31) 进入 else 分支，打印 “Inside else”，但永远不会到达 __syncthreads() 屏障。

  • 因此，先到达屏障的16个线程将无限期等待另外16个线程到达屏障，而这是永远不会发生的。这就导致了 死锁（Deadlock）。

关键点：线程束内的线程如果因为条件分支而无法全部到达同步点，就会导致死锁。这是编写CUDA代码时需要特别注意的问题。

GPU内存模型概述

前面我们主要关注了GPU的计算方面。现在，让我们开始了解GPU的另一个核心方面——内存。

GPU拥有复杂的内存层次结构，设计目标是隐藏全局内存的高延迟并提高数据访问带宽。以下是其主要的层次：

线程 (Thread) -> 寄存器 (Registers) [最快，每个线程私有]
        |
        V
线程块 (Block) -> 共享内存 (Shared Memory) [较快，块内线程共享]
        |
        V
全局内存 (Global Memory) / 设备内存 (Device Memory) [较慢，所有线程可访问]
        |
        V
纹理内存 (Texture Memory) 和 常量内存 (Constant Memory) [特殊用途，只读]

各级内存详解

  1. 寄存器 (Registers)

    • 位置与速度：位于芯片上，是速度最快的内存。

    • 作用域：每个线程私有。编译器会尽可能将局部变量分配到寄存器中。

    • 容量：数量有限。如果寄存器不足，变量会“溢出”到更慢的本地或全局内存。

  1. 共享内存 (Shared Memory)

    • 位置与速度：位于芯片上，速度仅次于寄存器，但远快于全局内存。

    • 作用域：同一个线程块内的所有线程共享。是块内线程通信和协作的主要桥梁。

    • 用途：常用于存储需要被块内线程频繁访问的数据片段，例如矩阵乘法中的瓦片（Tile）。

  2. 全局内存 (Global Memory)

    • 位置与速度：即GPU的板载DRAM（设备内存）。容量大（数GB），但延迟非常高（数百个时钟周期）。

    • 作用域：对所有网格中的所有线程可见。是主机（CPU）与设备（GPU）之间传输数据的主要区域，也是不同线程块之间通信的媒介（尽管不直接同步）。

    • 带宽：现代GPU的全局内存带宽很高（可达数百GB/s），但高延迟仍是主要性能瓶颈。

  3. 纹理内存 (Texture Memory) 和 常量内存 (Constant Memory)

    • 特性：它们是全局内存中经过特殊优化、具有独立缓存的空间。

    • 只读性：内核通常只能读取它们的内容。

    • 纹理内存：针对具有空间局部性的访问模式（如图像处理中的邻近像素访问）进行了优化。

    • 常量内存：适合存储所有线程都需要读取的常量数据。当所有线程读取相同地址时，性能极高（广播机制）。

    • 容量：很小（常量内存通常64KB，纹理内存缓存也有限）。

性能启示

编程时应遵循以下原则以提升性能：

  • 优先使用寄存器：尽量使用局部变量。

  • 合理利用共享内存：将需要重复访问的数据从全局内存缓存到共享内存。

  • 优化全局内存访问：通过合并访问（Coalesced Access）来最大化内存带宽利用率（后续课程详述）。

  • 善用常量与纹理内存：将只读的常量或具有空间局部性的数据放置于此。

总结

本节课我们深入探讨了CUDA中的线程发散问题。通过多个示例，我们学习了：

  1. 如何分析条件分支、循环和switch语句导致的线程束执行步数。

  2. 如何运用数学等价变换和位操作等技巧，将发散代码重构为无分支代码，从而显著提升性能，甚至完全消除发散。

  3. 线程发散在同步点（如__syncthreads()）可能引发死锁，需要谨慎处理。

最后，我们开启了GPU内存主题的学习，概述了GPU的内存层次结构，包括寄存器、共享内存、全局内存、纹理内存和常量内存的特性与用途。理解这个层次结构对于编写高性能CUDA程序至关重要。

在接下来的课程中，我们将更深入地研究内存访问模式、共享内存的使用以及如何实现高效的内存访问。

006：GPU内存类型、延迟与局部性 🧠

概述

在本节课中，我们将深入学习GPU的内存层次结构，理解带宽与延迟的概念，并探讨如何利用局部性原理来优化CUDA程序的性能。

GPU内存层次结构回顾

上一节我们深入探讨了线程发散问题。本节中，我们来看看GPU的内存组织。

GPU的内存层次结构与计算层次结构非常相似，了解这种层次结构有助于我们优化应用程序性能，例如利用共享内存来降低延迟。

以下是GPU的典型内存层次结构：

  • 寄存器：每个线程私有，速度最快，容量有限。

  • 共享内存：每个线程块共享，可编程，速度快。

  • L1缓存/共享内存：每个流多处理器（SM）配置，部分可编程（作为共享内存）。

  • L2缓存：所有SM共享，用于实现快速原子操作等。

  • 全局内存：所有线程和网格共享，主机与设备通信的主要通道，容量大但延迟高。

  • 常量内存：只读，所有线程共享，容量小（64KB），带宽高。

  • 纹理内存：只读，所有线程共享，容量小（约2KB），针对2D空间局部性优化。

核心概念公式：

  • 全局内存延迟：400 - 800 个时钟周期

  • 共享内存延迟：20 - 30 个时钟周期

  • 全局内存带宽：~200 GB/s

  • 共享内存带宽：~1 TB/s

带宽与延迟

带宽

带宽衡量的是单位时间内能够传输的数据量，而不是数据传输的速度。在GPU中，高带宽意味着可以同时进行大量的数据读写操作。

类比：飞机运输少量人员速度快（低延迟），火车运输大量人员速度慢但运量大（高带宽）。GPU类似于火车，虽然单次数据访问可能较慢，但可以并行处理大量数据访问，从而获得高吞吐量。

以下是提高带宽的几种通用技术：

  • 数据重用：在设备上复用已存在的数据，减少主机与设备间的数据传输。

  • 数据压缩：在传输前压缩数据，减少传输量。条件是：压缩时间 + 传输压缩数据时间 < 传输原始数据时间。

  • 重计算：在某些情况下，与其将结果从设备传回主机，不如在主机上重新计算，以节省回传的带宽。

延迟

延迟是指完成一次操作（如内存读写）所需的时间。理想情况下延迟应为零，但内存访问延迟不可避免。

在CPU上，主要通过以下方式隐藏内存访问延迟：

  1. 缓存层次结构（L1, L2, L3, L4缓存）。

  2. 多线程/多进程：当一个线程等待IO时，调度其他线程执行。

  3. 写回缓存：先写入高速缓存，CPU无需等待数据写入慢速内存。

  4. 指令级并行：执行不依赖于正在等待IO的指令的下一条指令。

  5. 流水线：利用硬件不同单元并行工作。

在GPU上，由于线程数量极多，每个线程可用的缓存容量很小，因此通过缓存隐藏延迟的效果有限。GPU主要依靠大规模多线程来隐藏延迟。

延迟隐藏原理：当一个线程块因内存访问（IO）而暂停时，流多处理器（SM）会立即切换到另一个就绪的线程块执行计算。通过让大量线程块交替进行计算和IO，可以保持SM始终繁忙，从而隐藏IO延迟，甚至获得超线性加速。

超线性加速示例：

  • 顺序执行（无延迟隐藏）：4个线程块，每个需5个时间单位，共需 4 * 5 = 20 个单位时间。

  • 利用多线程隐藏延迟执行：4个线程块交错执行，仅需约8个单位时间即可全部完成。

  • 加速比：20 / 8 ≈ 2.5，超过了线程块数量（4）的线性比例（1），故为超线性加速。

局部性原理

局部性原理对CPU和GPU的性能都至关重要，缓存就是为利用局部性而设计的。GPU线程块内的线程会利用L1和L2缓存。

GPU的L1缓存（64KB）是可配置的，程序员可以指定一部分作为可编程的共享内存（Scratchpad Memory），另一部分作为硬件管理的L1缓存。程序员应积极利用共享内存来存储被频繁访问的数据片段，以利用局部性。

在GPU上，空间局部性尤为重要，它与我们后续将学到的内存合并访问技术紧密相关。

局部性的类型

在CPU架构中，主要利用两种局部性：

  1. 时间局部性：如果某个内存位置被访问，那么它很可能在不久的将来再次被访问。

    • 代码示例：

      a[i] = 10; // 第一次访问（写）
a[i] = a[i] + 5; // 再次访问（读+写）
b = a[i]; // 第三次访问（读）

  2. 空间局部性：如果某个内存位置被访问，那么它附近的内存位置很可能也会被访问。

    • 代码示例：

      for (int i = 0; i < N; i++) {
    sum += a[i]; // 访问 a[0], a[1], a[2]... 这些地址在内存中是连续的
}

缓存的设计就是为了利用这两种局部性。当访问某个字（word）时，其所在的整个内存块会被载入缓存。如果后续访问同一位置（时间局部性）或其附近位置（空间局部性），就可以快速从缓存中获取数据。

案例分析：矩阵乘法的循环顺序优化

以下是两个完成矩阵乘法 C = A * B 的代码片段，假设矩阵以行主序方式存储。

代码片段 1:

for (int i = 0; i < N; i++) {
    for (int j = 0; j < P; j++) {
        for (int k = 0; k < M; k++) {
            C[i][j] += A[i][k] * B[k][j]; // B 被以列主序访问
        }
    }
}

代码片段 2:

for (int i = 0; i < N; i++) {
    for (int k = 0; k < M; k++) {
        for (int j = 0; j < P; j++) {
            C[i][j] += A[i][k] * B[k][j]; // A, B, C 均被以行主序访问
        }
    }
}

问题：哪个代码片段性能更好？

分析与答案：

  • 在行主序存储下，按行顺序访问数组能充分利用空间局部性，因为相邻元素被一起载入缓存。

  • 代码片段1：C[i][j] 和 A[i][k] 是行主序访问，但 B[k][j] 是列主序访问（内层循环 k 变化最快）。这导致对 B 的访问无法利用缓存，性能较差。

  • 代码片段2：C[i][j]、A[i][k] 和 B[k][j] 全部是行主序访问（对于 B，外层 k 慢，内层 j 快，相当于访问 B 的行）。这能充分利用空间局部性，性能更好。

实测结果通常显示，代码片段2的执行时间显著短于代码片段1（例如，4.7秒 vs 9.5秒），因为它所有矩阵的访问模式都与内存布局一致，缓存命中率更高。

总结

本节课我们一起学习了GPU详细的内存层次结构及其特点，理解了带宽与延迟的区别以及GPU如何利用大规模多线程隐藏延迟。我们还回顾了局部性原理（时间局部性与空间局部性），并通过矩阵乘法的例子，看到了访问模式如何显著影响程序性能。理解这些概念是进行后续GPU内存优化（如使用共享内存、实现合并访问）的基础。

007：内存合并与AoS对SoA 🚀

概述

在本节课中，我们将要学习GPU内存访问中一个至关重要的性能优化概念——内存合并。我们将探讨什么是内存合并，以及如何通过选择合适的数据结构布局（数组结构体 与 结构体数组）来最大化内存合并的效率，从而显著提升GPU程序的性能。

上一节回顾

在上一节中，我们介绍了内存层次结构、带宽以及空间局部性和时间局部性的概念。我们通过CPU上的例子看到，遵循内存的存储顺序（如行主序）进行访问可以带来显著的性能提升。

本节中，我们将把焦点转移到GPU上，看看除了局部性之外，另一个对GPU内存性能影响巨大的因素——内存合并。

什么是内存合并？ 🤔

内存合并是GPU硬件层面的一项优化技术。其核心定义如下：

当一个线程束（32个线程）访问一个连续的内存块（例如，一个32字大小的块）时，这些零散的内存访问请求会被硬件“合并”成一个单一的、更宽的内存事务。

这意味着，原本可能需要32次独立内存访问的操作，现在可能只需要1次或几次合并后的访问即可完成。这极大地减少了内存延迟，提升了内存带宽的利用率。

合并访问示例

考虑一个经典的初始化数组的内核：

__global__ void initArray(int *A) {
    int tid = blockIdx.x * blockDim.x + threadIdx.x;
    A[tid] = 0; // 合并的内存访问
}

在这个例子中，线程0访问A[0]，线程1访问A[1]，...，线程31访问A[31]。由于这些元素在内存中是连续存储的，并且位于同一个内存块内，因此这32次访问会被硬件合并成一次高效的内存事务。

非合并访问的代价

如果没有内存合并，每个线程的加载/存储指令都需要独立的内存周期。在最坏情况下，一个线程束可能需要32个内存周期，这将导致严重的性能瓶颈，使得GPU的计算能力无法充分发挥。

内存合并的内部实现 🛠️

内存合并是由GPU硬件自动完成的，程序员无法直接控制其开关。其工作方式与线程束的执行模型紧密相关：

  • 线程束以锁步方式执行指令。

  • 当线程束执行一条内存指令（加载或存储）时，内存控制器会检查所有线程请求的内存地址。

  • 如果这些地址落在一个对齐的、连续的内存块（例如32字节、64字节或128字节，具体取决于GPU架构和访问类型）内，硬件就会将这些请求合并为一个内存事务。

关键要点：内存合并是硬件的固有特性。程序员的责任是通过组织数据和访问模式来促成合并访问的发生。

内存访问模式分析 📊

根据线程束内线程访问内存的方式，我们可以将访问模式分为三类：

以下是三种典型的内存访问模式图示：

  1. 完全合并访问：线程束中所有线程访问连续且位于同一内存块内的地址。这是最理想的情况。

  1. 部分合并访问：线程束中只有部分线程的访问地址是连续的，其他线程的访问地址分散在不同的内存块中。这会导致多个内存事务。

  2. 完全非合并访问：线程束中每个线程访问的地址都位于完全不同的、不连续的内存块中。这会导致最差的性能，需要最多（32次）的内存事务。

代码示例：步长访问与随机访问 💻

让我们通过具体代码来理解不同的访问模式。

1. 完全合并访问

// 经典示例：连续访问
A[tid] = value;

2. 完全非合并访问（通过大跨度步长实现）

// 假设 chunk_size = 32 (或更大)
int start = tid * chunk_size;
int end = start + chunk_size;
for (int i = start; i < end; i++) {
    A[i] = value; // 每个线程访问完全独立的内存块
}

当chunk_size为32时，线程0访问元素0-31，线程1访问元素32-63...，彼此毫无重叠，导致完全非合并。

3. 部分合并访问

// 假设 chunk_size = 2
int start = tid * chunk_size;
int end = start + chunk_size;
for (int i = start; i < end; i++) {
    A[i] = value;
}

当chunk_size为2时，前16个线程（tid 0-15）访问的元素（如0-1, 2-3, ..., 30-31）可能位于同一个或两个内存块内，从而产生部分合并。这通常需要2次内存事务。

4. 随机访问

// 索引来自另一个输入数组，访问模式不可预测
int index = input[tid];
A[index] = value;

这种模式在图形算法（如BFS、SSSP）中很常见，其合并效率完全取决于input数组中的数据分布，通常是低效的。

数据结构布局：AoS vs SoA 🏗️

数据在内存中的组织方式（数据结构布局）对内存合并有决定性影响。主要分为两种：

数组结构体

数组结构体 是CPU编程中常见的形式。它将一个对象的所有属性打包在一起，然后创建该对象的数组。

// AoS (Array of Structures)
struct Student {
    int roll_no;
    double marks;
    char name[50];
};
Student students[1000]; // 数组的每个元素都是一个完整的Student结构体

  • CPU优势：利于空间局部性。当CPU线程访问students[i].roll_no时，很可能紧接着访问students[i].marks，后者可能已被缓存。

  • GPU劣势：当GPU线程束访问不同学生的roll_no时，访问地址的跨度是sizeof(Student)，这会导致跨步访问，严重阻碍内存合并。

结构体数组

结构体数组 是GPU优化中推荐的形式。它将所有对象的同一属性分别存储在独立的数组中。

// SoA (Structure of Arrays)
struct StudentData {
    int roll_no[1000];
    double marks[1000];
    char name[1000][50];
};
StudentData students; // 一个结构体，其成员是多个数组

  • GPU优势：利于内存合并。当GPU线程束访问不同学生的roll_no时（即访问roll_no[0], roll_no[1], ...），这些地址是连续的，可以实现完美的合并访问。

  • CPU劣势：可能破坏空间局部性，因为访问一个学生的所有属性需要从三个不同的数组位置读取。

GPU上的AoS与SoA实现对比 ⚖️

以下是两种布局在GPU内核中访问方式的对比：

AoS 访问方式（性能较差）

// 内核函数
__global__ void processAoS(Student *d_students) {
    int tid = blockIdx.x * blockDim.x + threadIdx.x;
    // 跨步访问，不利于合并
    d_students[tid].roll_no = tid;
    d_students[tid].marks = tid * 10.0;
}
// 内存分配与拷贝
Student *d_students;
cudaMalloc(&d_students, N * sizeof(Student));
cudaMemcpy(d_students, h_students, N * sizeof(Student), cudaMemcpyHostToDevice);

SoA 访问方式（性能更优）

// 内核函数
__global__ void processSoA(int *roll_no, double *marks) {
    int tid = blockIdx.x * blockDim.x + threadIdx.x;
    // 连续访问，利于合并
    roll_no[tid] = tid;
    marks[tid] = tid * 10.0;
}
// 内存分配与拷贝（需要为每个数组成员单独分配）
int *d_roll_no;
double *d_marks;
cudaMalloc(&d_roll_no, N * sizeof(int));
cudaMalloc(&d_marks, N * sizeof(double));
cudaMemcpy(d_roll_no, h_roll_no, N * sizeof(int), cudaMemcpyHostToDevice);
cudaMemcpy(d_marks, h_marks, N * sizeof(double), cudaMemcpyHostToDevice);

性能对比：在实际测试中，使用SoA布局的内核执行时间可能仅为使用AoS布局内核的1/3，这巨大的性能提升主要归功于高效的内存合并。

总结

本节课中我们一起学习了GPU并行编程中两个核心的优化概念：

  1. 内存合并：GPU硬件将线程束内连续的、对齐的内存访问合并为单个宽内存事务的机制，是获得高内存带宽的关键。

  2. AoS与SoA：数据结构布局的选择直接影响内存合并的效率。

    • AoS 在CPU上利用局部性有优势，但在GPU上会导致跨步访问，破坏合并。

    • SoA 是GPU编程的推荐模式，它能确保对同一属性的访问是连续的，从而实现完美的内存合并，带来显著的性能提升。

作为程序员，在设计GPU程序时，应当时刻考虑数据的访问模式，并优先采用SoA布局来最大化内存合并的效益。

008：GPU中的共享内存 🧠

概述

在本节课中，我们将要学习CUDA编程中一个非常强大的特性：共享内存。我们将了解什么是共享内存、如何声明和使用共享内存变量，以及为什么它在某些场景下至关重要。我们还将通过示例来探讨使用共享内存时可能遇到的同步问题。

共享内存简介

上一节我们介绍了内存合并（Memory Coalescing）的概念，它对于提升GPU内存访问效率至关重要。本节中，我们来看看另一种可以显著提升性能的内存：共享内存。

共享内存是GPU上L1缓存的一部分，但它与CPU缓存有一个关键区别：共享内存是可由程序员显式控制的。在CPU上，数据如何进入L1、L2、L3缓存完全由硬件管理，程序员无法干预。而在GPU上，程序员可以指定哪些数据应该存放在共享内存中，这为我们优化程序提供了极大的灵活性。

共享内存主要有两个关键特性：

  1. 按线程块分配：每个线程块（Thread Block）都拥有自己独立的一块共享内存。即使多个线程块被调度到同一个流式多处理器（SM）上执行，它们各自的共享内存也是隔离的。

  2. 线程块内共享：在同一个线程块内的所有线程，可以看到并修改共享内存中的同一份数据副本。这使得共享内存成为线程块内线程间通信和协作的理想场所。

共享内存通常用于以下两种情况：

  • 数据重用：当线程块内需要反复访问一小块数据时，将其加载到共享内存中可以极大降低访问延迟。

  • 线程协作：当需要在线程块内的线程之间进行同步或协调时（例如，实现一个屏障或归约操作），共享内存可以作为共享的工作空间。

如何声明共享内存

在CUDA内核（Kernel）中声明共享内存变量非常简单，只需在变量声明前加上 __shared__ 限定符即可。

以下是声明共享内存数组和变量的语法：

// 在内核函数中声明一个共享内存数组
__shared__ float shared_array[1024];

![](https://github.com/OpenDocCN/dsai-notes-pt1-zh/raw/master/docs/nv-cuda-prll-prog/img/2a3df8b485a890f6aeea02b32796c9af_29.png)

![](https://github.com/OpenDocCN/dsai-notes-pt1-zh/raw/master/docs/nv-cuda-prll-prog/img/2a3df8b485a890f6aeea02b32796c9af_31.png)

![](https://github.com/OpenDocCN/dsai-notes-pt1-zh/raw/master/docs/nv-cuda-prll-prog/img/2a3df8b485a890f6aeea02b32796c9af_32.png)

![](https://github.com/OpenDocCN/dsai-notes-pt1-zh/raw/master/docs/nv-cuda-prll-prog/img/2a3df8b485a890f6aeea02b32796c9af_34.png)

![](https://github.com/OpenDocCN/dsai-notes-pt1-zh/raw/master/docs/nv-cuda-prll-prog/img/2a3df8b485a890f6aeea02b32796c9af_36.png)

![](https://github.com/OpenDocCN/dsai-notes-pt1-zh/raw/master/docs/nv-cuda-prll-prog/img/2a3df8b485a890f6aeea02b32796c9af_37.png)

![](https://github.com/OpenDocCN/dsai-notes-pt1-zh/raw/master/docs/nv-cuda-prll-prog/img/2a3df8b485a890f6aeea02b32796c9af_39.png)

![](https://github.com/OpenDocCN/dsai-notes-pt1-zh/raw/master/docs/nv-cuda-prll-prog/img/2a3df8b485a890f6aeea02b32796c9af_41.png)

// 在内核函数中声明一个共享内存变量
__shared__ unsigned int shared_variable;

重要说明：

  • 这些声明必须写在CUDA内核函数内部。

  • 如果声明时不加 __shared__，例如 float local_array[1024];，那么这个数组将是线程局部的（Thread-local），每个线程都有自己的副本，线程间无法直接访问彼此的数据。

  • 加上 __shared__ 后，该变量在线程块内是共享的。例如，一个包含1024个线程的线程块，所有线程都操作 shared_array 的同一份物理存储。

示例与分析：理解共享内存的并发访问

为了深入理解共享内存的共享特性以及可能引发的并发问题，我们来看一个具体的代码示例。

假设我们启动一个内核，只使用一个线程块，该线程块包含1024个线程。内核代码如下：

__global__ void myKernel() {
    __shared__ int s;
    int i = threadIdx.x;

    if (i == 0) {
        s = 0;          // 语句A：线程0初始化s为0
    }
    if (i == 1) {
        s += 1;         // 语句B：线程1给s加1
    }
    if (i == 100) {
        s += 2;         // 语句C：线程100给s加2
    }
    if (i == 0) {
        printf("%d\n", s); // 语句D：线程0打印s的值
    }
}

我们需要思考：这个程序可能的输出结果是什么？

关键概念回顾：线程束（Warp）与锁步执行

在分析之前，必须回忆GPU的执行模型：线程以32个为一组（称为一个Warp）进行调度。同一个Warp内的线程是“锁步”（Lock-step）执行的，即它们同时执行同一条指令（尽管可能因为分支而产生部分线程不活跃的情况）。

在我们的例子中：

  • 线程0和线程1属于同一个Warp（Warp 0）。

  • 线程100属于另一个Warp（Warp 3）。

并发执行与可能的输出

由于不同Warp之间的执行顺序是不确定的，语句A、B、C、D可能以多种方式交织执行，导致最终 s 的值不确定。以下是几种可能的执行顺序及其结果：

  1. 顺序执行 (A -> B -> C -> D)：s 最终为 0 + 1 + 2 = 3。输出 3。

  2. 线程100先执行 (C -> A -> B -> D)：线程100读取未初始化的 s（垃圾值），加2后写入。随后A将其覆盖为0，B加1。最终 s 为 1。输出 1。

  3. A和C并发执行产生竞争：线程0执行 s=0 和线程100执行 s+=2 同时发生。s+=2 操作在底层可能被分解为“读-改-写”三步。如果线程100读取了垃圾值，计算了垃圾值+2，但在写入前，线程0写入了0，线程1写入了1，最后线程100的写入覆盖了结果。最终 s 可能为 垃圾值+2。

  4. 另一种竞争情况：类似情况3，但线程1的 s+=1 在线程100最终写入后才执行，最终 s 可能为 垃圾值+3。

  5. B和C并发执行：线程1的 s+=1 和线程100的 s+=2 并发执行，且 s 已被初始化为0。如果线程100的写入后发生，s 最终为 2。输出 2。

结论：这个简单的程序由于存在数据竞争（Data Race），可能产生多种输出（如1, 2, 3, 垃圾值+2, 垃圾值+3），其行为是非确定性的。

这个例子清晰地展示了当多个线程不加协调地访问共享内存时会导致的问题。要获得确定性的、正确的结果，我们必须引入同步机制。

共享内存的典型使用模式

让我们回到一个更实际的例子，看看共享内存如何被应用，以及同步为何必不可少。

问题：有一个1024x1024的矩阵 M。我们想对每个元素 M[i][j] 进行如下操作：M[i][j] = M[i][j] + M[i][j+1]（对最后一列不操作）。我们计划将每一行分配给一个线程块（共1024个线程块），每个线程块内的1024个线程处理该行的1024个元素。

为了优化性能，我们决定将一行数据从全局内存加载到共享内存中进行计算。实现步骤可分为三步：

  1. 声明并加载数据：在线程块内声明一个共享内存数组，然后将全局内存中对应行的数据复制进来。

  2. 执行计算：所有线程使用共享内存中的数据执行 shared_data[j] = shared_data[j] + shared_data[j+1]。

  3. 写回结果：将计算后的共享内存数据写回全局内存的对应行。

这个流程中存在两个关键的同步问题：

问题一：步骤间的屏障

  • 在步骤1和步骤2之间，必须确保线程块内所有线程都已完成数据加载，才能开始计算。否则，某些线程可能还在加载数据，而另一些线程已经开始使用未完全加载的数据进行计算。

  • 同理，在步骤2和步骤3之间，必须确保所有线程都已完成计算，才能将结果写回全局内存。

问题二：计算步骤中的数据竞争

  • 在步骤2中，线程 j 需要读取 shared_data[j] 和 shared_data[j+1]。而线程 j+1 同时也在读取 shared_data[j+1] 和 shared_data[j+2]。如果操作是就地更新（in-place update），并且没有同步，就可能发生我们前面例子中类似的竞争条件，导致结果错误。

为了解决这些问题，CUDA提供了 __syncthreads() 函数。它是一个线程块级别的屏障，调用该函数会阻塞线程块内的所有线程，直到每个线程都到达这个调用点，然后所有线程才能继续执行后续指令。

一个正确的、使用了同步的伪代码结构如下：

__global__ void matrixKernel(float* M) {
    __shared__ float row_data[1024];
    int tid = threadIdx.x;
    int row_start = blockIdx.x * 1024;

    // 步骤1：协作加载数据
    row_data[tid] = M[row_start + tid];
    __syncthreads(); // 等待所有线程加载完毕

    // 步骤2：执行计算（注意处理边界）
    if (tid < 1023) { // 最后一列不计算
        float temp = row_data[tid] + row_data[tid + 1];
        __syncthreads(); // 确保所有读取操作完成后再写入
        row_data[tid] = temp;
    }
    __syncthreads(); // 等待所有计算完成

    // 步骤3：写回结果
    M[row_start + tid] = row_data[tid];
}

注意：实际处理此类相邻元素依赖的计算时，模式可能更复杂（例如使用双缓冲），但此示例展示了同步的基本用法。

总结

本节课中我们一起学习了CUDA中共享内存的核心知识：

posted @ 2026-10-17 13:17  绝不原创的飞龙  阅读(4)  评论(0)    收藏  举报