0. 从C++到CUDA


从 C++ 到 CUDA:完整入门笔记

第一部分:核心思维转变

1.1 设计哲学的差异

CPU GPU
核心数量 少(8~32) 多(数千)
单核能力 强(缓存、分支预测、乱序执行) 弱(砍掉了大部分控制逻辑)
擅长 复杂的串行逻辑 对海量数据做同样的简单操作
设计目标 低延迟 高吞吐

1.2 四个关键思维转变

从”循环”到”线程网格”

1
2
3
4
5
6
// CPU 视角:我循环 N 次
for (int i = 0; i < N; ++i) C[i] = A[i] + B[i];

// GPU 视角:启动 N 个线程,每个线程做一件事
int idx = blockIdx.x * blockDim.x + threadIdx.x;
if (idx < N) C[idx] = A[idx] + B[idx];

每个线程相当于循环的一次迭代被”展开”。

从”我”到”我们”

写 kernel 时,代码会被几万个线程同时执行。要时刻想”这一万个’我’各自在做什么,会不会撞车”。读写同一地址要用原子操作,block 内同步用 __syncthreads()

从”算得快”到”喂得饱”

GPU 编程中算力几乎从来不是瓶颈,瓶颈是内存带宽——能不能让那几千个核心都有活干,不是在等数据。这是和 ICPC/CPU 优化最大的不同。

从”按需分配”到”批量并行”

GPU 不擅长动态内存分配、复杂分支、不规则访问。算法常常被重新设计成规则的、数组化的形式。


第二部分:CUDA 的两套内存空间

2.1 物理上分离的两套内存

1
2
3
4
5
CPU 主内存 (Host)              GPU 显存 (Device)
┌──────────────┐ ┌──────────────────┐
│ malloc 分配 │ PCIe / NVLink │ cudaMalloc 分配 │
│ h_ 前缀指针 │ ◄────────────► │ d_ 前缀指针 │
└──────────────┘ cudaMemcpy └──────────────────┘

关键约束

  • Host 指针不能在 GPU kernel 里解引用
  • Device 指针不能在 CPU 代码里解引用
  • 编译器不会帮你检查,错用会段错误或更怪异的行为
  • 数据必须用 cudaMemcpy 显式搬运

2.2 PCIe 是性能瓶颈

  • PCIe 4.0 x16:~32 GB/s
  • GPU 显存带宽:~1000+ GB/s(差几十倍)

反模式:拷数据上去 → 跑简单 kernel → 拷回来——传输时间可能比计算时间还长。
最佳实践:数据一次上传后,留在显存上多次复用。

LeetGPU 的题目数据已经在 device 上,不涉及这层。

2.3 指针变量 vs 指向的内存

1
2
float *d_A;                  // (1) 栈上的指针变量,初值随机(脏数据)
cudaMalloc(&d_A, bytes); // (2) 把 GPU 显存地址写进 d_A,覆盖随机值
1
2
3
4
5
CPU 栈                          GPU 显存
┌─────────────┐ ┌──────────────────────┐
│ d_A: 0x7f.. │ ─────────────▶ │ 100MB 的连续区域 │
└─────────────┘ └──────────────────────┘
8 字节 cudaMalloc 保证虚拟连续

两件事独立

  • 指针变量本身的初值(随机)→ 被 cudaMalloc 覆盖
  • 分配的空间连续性 → cudaMalloc 的契约保证(虚拟地址连续,物理上由 GPU MMU 处理)

第三部分: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
#include <cuda_runtime.h>

__global__ void kernel(const float* A, float* C, int N) {
int idx = blockIdx.x * blockDim.x + threadIdx.x;
if (idx < N) C[idx] = A[idx] * 2;
}

int main() {
const int N = 1 << 20;
const size_t bytes = N * sizeof(float);

// 1. Host 准备数据
float* h_A = (float*)malloc(bytes);
float* h_C = (float*)malloc(bytes);
for (int i = 0; i < N; ++i) h_A[i] = 1.0f * i;

// 2. cudaMalloc 分配显存
float *d_A, *d_C;
cudaMalloc(&d_A, bytes);
cudaMalloc(&d_C, bytes);

// 3. cudaMemcpy 上传数据
cudaMemcpy(d_A, h_A, bytes, cudaMemcpyHostToDevice);

// 4. 启动 kernel
int threadsPerBlock = 256;
int blocksPerGrid = (N + threadsPerBlock - 1) / threadsPerBlock;
kernel<<<blocksPerGrid, threadsPerBlock>>>(d_A, d_C, N);

// 5. cudaMemcpy 下载结果(隐式同步)
cudaMemcpy(h_C, d_C, bytes, cudaMemcpyDeviceToHost);

// 6. 使用结果
printf("h_C[100] = %f\n", h_C[100]);

// 7. 释放资源
cudaFree(d_A); cudaFree(d_C);
free(h_A); free(h_C);
return 0;
}

LeetGPU 帮你做了 1、2、3、5、6、7,你只写第 4 步(kernel + 启动)。

3.1 cudaMemcpy 方向参数

  • cudaMemcpyHostToDevice:H → D
  • cudaMemcpyDeviceToHost:D → H
  • cudaMemcpyDeviceToDevice:D 内部
  • cudaMemcpyHostToHost:基本用不到

cudaMemcpy(同步版)会隐式同步——会等队列里所有 kernel 跑完。这就是为什么很多简单程序不需要显式 cudaDeviceSynchronize()


第四部分:CUDA 语法层面的差异

4.1 编译

C++ CUDA
文件后缀 .cpp .cu
编译器 g++ / clang++ nvcc

nvcc 是驱动程序:把 .cu 拆成 host 部分(丢给 g++)和 device 部分(自己编译成 PTX),最后链接。可以在 .cu 里自由混写 C++ 和 CUDA。

4.2 函数修饰符

修饰符 在哪执行 从哪调用
__host__(默认) CPU CPU
__global__ GPU CPU
__device__ GPU GPU

__host__ __device__ 可叠加,让函数两边都能编译:

1
__host__ __device__ float square(float x) { return x * x; }

4.3 kernel 调用语法

1
2
kernel<<<grid, block>>>(args);
kernel<<<grid, block, sharedMemBytes, stream>>>(args); // 完整形式

<<<...>>> 是 nvcc 特有的语法扩展。

4.4 GPU 代码的限制

  • 不能 throw / try-catch
  • 不能用 std::vectorstd::string
  • 标准库基本不能用(有 GPU 版的 Thrust 库)
  • 错误检查靠返回值

推荐宏

1
2
3
4
5
6
7
8
#define CUDA_CHECK(call) do {                                  \
cudaError_t err = (call); \
if (err != cudaSuccess) { \
fprintf(stderr, "CUDA error %s:%d: %s\n", \
__FILE__, __LINE__, cudaGetErrorString(err)); \
exit(1); \
} \
} while (0)

CUDA 错误经常是悄无声息的,不查返回值会让你 debug 时怀疑人生。

4.5 异步执行模型

1
2
3
kernel<<<...>>>(...);          // 不阻塞,CPU 立即继续
cpu_do_something(); // 可能在 kernel 还没跑完时就执行
cudaDeviceSynchronize(); // 显式等 GPU 跑完

陷阱:错误也是异步报告的。kernel 内越界访问,<<<...>>> 调用本身不会报错;要等下一次 CUDA API 调用才能拿到错误。

debug 时的标准做法:

1
2
3
kernel<<<...>>>(...);
CUDA_CHECK(cudaGetLastError()); // 检查启动配置错误
CUDA_CHECK(cudaDeviceSynchronize()); // 等 kernel 跑完,捕获运行时错误

第五部分:执行模型——grid / block / thread / warp

5.1 三层结构

1
2
3
4
Grid (网格)
└── Block (线程块) × blocksPerGrid 个
└── Thread (线程) × threadsPerBlock 个
└── 32 个组成一个 Warp(硬件层面)

类比:

  • thread = 一个工人
  • block = 一个工程队(共用工棚 = shared memory,可互相喊话 = __syncthreads()
  • grid = 整个工地

5.2 Block 和 SM 的关系

GPU 硬件层面的执行单元叫 SM (Streaming Multiprocessor)

  • 一个 block 整体调度到某个 SM 上执行
  • block 内线程共享该 SM 上的 shared memory
  • block 之间可能在不同 SM 上,互相隔离不能直接通信

5.3 Warp:硬件真实调度单位

block 里的线程在硬件上32 个一组(一个 warp)调度,整组一起执行同一条指令——SIMT (Single Instruction Multiple Threads)。

推论

  • threadsPerBlock 必须是 32 的倍数,否则有”死线程”占位
  • 常用值:128、256、512
  • 256 是不错的默认值

5.4 Warp Divergence

warp 里所有线程共用一个 PC(程序计数器),必须执行同一条指令。遇到 if-else:

1
2
if (cond) { A; }
else { B; }

执行流程:

  1. mask = “满足 cond 的线程”,所有线程一起走 A,不满足的结果被丢弃
  2. mask 反过来,再走一遍 B
  3. A 和 B 串行执行,时间是两者之和

这就是 warp divergence 的代价。GPU 不擅长复杂分支逻辑的根本原因。

Volta(2017)后引入 Independent Thread Scheduling,每线程独立 PC,但代价模型基本不变。


第六部分:内置变量——CUDA 最重要的疑问

6.1 kernel 调用的两组信息

1
2
3
4
5
vector_add  <<<grid, block>>>  (A, B, C, N);
^^^^^^^^^^^^^^^^^ ^^^^^^^^^^^^
执行配置 函数参数
传给硬件调度器 和普通 C++ 一样
告诉启动多少线程 所有线程收到同样的参数

6.2 内置变量从哪来

blockIdxthreadIdx不是函数参数,是硬件给每个线程提供的”内置变量”——本质是特殊寄存器。

内置变量 含义 类型
threadIdx 我在 block 内的坐标 dim3
blockIdx 我在 grid 内的坐标 dim3
blockDim block 的尺寸 dim3
gridDim grid 的尺寸 dim3
warpSize warp 大小(永远 32) int

PTX 层面:

1
2
3
4
mov.u32  %r1, %ctaid.x;     // blockIdx.x
mov.u32 %r2, %ntid.x; // blockDim.x
mov.u32 %r3, %tid.x; // threadIdx.x
mad.lo.s32 %r4, %r1, %r2, %r3;

%ctaid%ntid%tid 是 PTX 的特殊寄存器——硬件保证读它们时返回当前线程的身份。

类比 CPU 的 TLS、ARM 的 MPIDR——硬件给执行单元注入身份信息。

6.3 为什么这种设计

如果让 CUDA 把 idx 作为参数传进去,就要预先生成 25M 组参数、启动 25M 次函数调用——CPU 多线程模型,开销巨大。

GPU 的做法:你描述启动多少线程,硬件自己给每个线程编号,每个线程自己查自己的编号。

  • 启动开销极低
  • 编号天然分布式
  • 代码可以写成”我是谁,我做什么”的视角

6.4 心智模型:CUDA 的”两个世界”

逻辑世界(业务):

  • 数据:A、B、C 是什么
  • 参数:N 是多少
  • 通过 (A, B, C, N) 函数参数传

硬件世界(GPU 调度):

  • 启动多少线程:<<<grid, block>>>
  • 我是谁:threadIdxblockIdx(硬件注入)
  • 总共多少人:blockDimgridDim(硬件注入)

这两个世界通过 idx = blockIdx.x * blockDim.x + threadIdx.x 桥接——把”我是哪个线程”翻译成”我处理哪份数据”。


第七部分:线程到数据的映射

7.1 重要澄清

CUDA 不会自动把数组元素分配给线程。它做的事:

  1. 把参数 (d_A, d_B, d_C, N) 复制给每个线程——所有线程看到相同的参数
  2. 启动 blocksPerGrid * threadsPerBlock 个线程
  3. 给每个线程分配唯一的 (blockIdx, threadIdx)
  4. 让每个线程跑同一份 kernel 代码

“线程 i 处理数据 i” 这种映射是你在 kernel 里写的

1
2
int idx = blockIdx.x * blockDim.x + threadIdx.x;
if (idx < N) C[idx] = A[idx] + B[idx];

如果你写 int idx = 0,所有 25M 线程会争抢写 C[0]——CUDA 不报错,结果就是错的。

7.2 各种映射模式

不同问题有不同的映射规则。这是 CUDA 编程的核心设计决策。

1D 一对一(向量加法):

1
int idx = blockIdx.x * blockDim.x + threadIdx.x;

Grid-stride loop(一个线程处理多个元素):

1
2
3
4
int stride = gridDim.x * blockDim.x;
for (int idx = blockIdx.x * blockDim.x + threadIdx.x;
idx < N;
idx += stride) { ... }

关键点:用 idx += stride 而非 idx += 1,是为了保证同一时刻 warp 里 32 个线程访问连续地址(合并访存)。

2D 网格(图像、矩阵):

1
2
3
int col = blockIdx.x * blockDim.x + threadIdx.x;
int row = blockIdx.y * blockDim.y + threadIdx.y;
int idx = row * W + col;

关键点:让 threadIdx.x 对应 col(连续维度),保证合并访存。

Block 协作(reduction、矩阵乘法):

每个 block 负责一段数据,block 内线程通过 shared memory 协作。

7.3 设计选择遵循的原则

不只是”逻辑上对”,还要:

  • 合并访存(warp 内 32 线程访问连续地址)
  • 避免 divergence(warp 内尽量走同一分支)
  • 充分并行(避免大量线程闲置)
  • 数据复用(用 shared memory 减少 global memory 访问)

第八部分:合并访存(Coalesced Access)—— GPU 优化第一原则

8.1 GPU 显存的特性

GPU 显存(GDDR/HBM)按 burst 读取,一次最少读一行(32~128 字节)。

CPU 用 cache line 处理(64 字节自动缓存附近内容)。

GPU 用 coalescing 处理:硬件检查 warp 内 32 个线程的访存地址,拼起来一次读

8.2 合并访存的规则

  • 理想:32 线程访问连续的、对齐的 128 字节 → 1 次内存事务
  • 最坏:32 线程访问 32 个不同的 128 字节区域 → 32 次事务,带宽利用率 1/32

8.3 例子

1
2
3
4
5
6
7
// 好:warp 内线程的 idx 是 0,1,2,...,31,地址连续
int idx = blockIdx.x * blockDim.x + threadIdx.x;
C[idx] = A[idx] + B[idx];

// 坏:warp 内线程访问跨度大
int idx = threadIdx.x * gridDim.x + blockIdx.x;
C[idx] = A[idx] + B[idx];

8.4 二维数据的陷阱

1
2
3
4
5
6
7
8
9
// 错(按 C 风格行主序存储时):threadIdx.x 跨行访问
int row = threadIdx.x;
int col = threadIdx.y;
M[row * W + col] // warp 内线程访问 M[0*W], M[1*W], M[2*W], ... 跨度巨大

// 对:threadIdx.x 沿行内访问
int col = threadIdx.x;
int row = threadIdx.y;
M[row * W + col] // warp 内线程访问连续地址

8.5 对齐要求

cudaMalloc 返回的指针对齐到 256 字节,从开头访问肯定对齐。但 A + 1 这种偏移可能跨两个 128 字节块,需要 2 次事务。


第九部分:性能模型——Roofline

9.1 算术强度(Arithmetic Intensity)

每读 1 字节能做多少次浮点运算。这是程序本身的属性,与硬件无关

向量加法:$\text{AI} = \frac{N}{12N} = \frac{1}{12} \approx 0.083$ FLOP/byte(极低)

矩阵乘法:$\text{AI} = O(N)$(很高)

9.2 Roofline 模型

每个硬件有两个指标:

  • 峰值算力 $P_{\max}$(FLOPs/s)
  • 峰值带宽 $B_{\max}$(bytes/s)

可达最大性能:

  • $\text{AI} < \text{AI}_{\text{crit}}$ → memory-bound,性能上限由带宽决定
  • $\text{AI} > \text{AI}_{\text{crit}}$ → compute-bound,性能上限由算力决定

9.3 A100 实例

  • FP32 峰值算力 ≈ 19.5 TFLOPs/s
  • HBM2e 带宽 ≈ 1555 GB/s
  • 临界 AI ≈ 12.5 FLOP/byte

向量加法 AI = 0.083 → memory-bound,性能上限 ≈ 130 GFLOPs/s(峰值算力的 0.7%)。

9.4 优化导向

对 memory-bound 任务:

  • ✅ 合并访存
  • ✅ 减少数据搬运(kernel fusion)
  • ❌ shared memory(每元素只读一次,无复用)
  • ❌ 复杂算法(计算已最少)

第十部分:内存管理细节

10.1 malloc 的底层

1
float* h_A = (float*)malloc(bytes);
  • 小请求:从 glibc 维护的堆里切(用 brk 系统调用扩展)
  • 大请求(≥128KB):走 mmap 直接向内核要匿名映射

懒分配malloc 只是登记 VMA,不分配物理内存。第一次写入触发缺页中断才分配物理页。

1
2
3
4
5
auto t1 = clock();
float* p = (float*)malloc(100*MB); // 几微秒,不分配物理内存
auto t2 = clock();
for (int i=0; i<N; ++i) p[i] = 0; // 真正的"分配"开销在这里(缺页)
auto t3 = clock();

10.2 cudaMalloc vs malloc

malloc cudaMalloc
内存位置 CPU DRAM GPU 显存
谁能访问 CPU GPU
大请求路径 mmap 驱动管理的显存分配器
分配时机 懒分配(首次访问) 立即分配
速度 几微秒 几十~几百微秒
是否清零

10.3 Pinned Memory(性能优化)

cudaMallocHost / cudaHostAlloc 分配页锁定内存

  • 不会被 swap
  • cudaMemcpy 更快(PCIe 4.0 上 ~25 GB/s vs ~12 GB/s)
  • 支持 cudaMemcpyAsync 异步传输

代价:占真实物理内存,分配慢,分配太多影响系统。


第十一部分:指针与索引

11.1 A[idx] ≡ *(A + idx)

步长由指针类型决定,不是索引类型

1
2
3
4
5
float* A;  int i = 3;
A + i → 地址 = A + 3 × sizeof(float) = A + 12

double* B; int i = 3;
B + i → 地址 = B + 3 × sizeof(double) = B + 24

x86-64 汇编层面:

1
movss xmm0, [rax + rcx*4]   ; *4 是编译器根据 float* 编进指令的

11.2 索引类型只关心溢出

1
2
3
4
5
6
// N 较小,int 够用
int idx = blockIdx.x * blockDim.x + threadIdx.x;

// N 超过 21 亿,要警惕
size_t idx = (size_t)blockIdx.x * blockDim.x + threadIdx.x;
// ↑ 必须在乘法前转换,否则按 32 位算先溢出再转

CUDA 的 threadIdx.x 等都是 unsigned int(32 位)。


第十二部分:浮点数 IEEE 754

12.1 float 结构

32 位 = 1 符号 + 8 指数(偏移 127)+ 23 尾数(隐含前导 1,有效 24 位 ≈ 7 位十进制)

12.2 浮点加法不结合

(a + b) + c ≠ a + (b + c) 一般成立。

例:(1e20 + (-1e20)) + 1 = 1,但 1e20 + (-1e20 + 1) = 0

对向量加法无影响(每位置独立)。对并行归约有影响——并行求和顺序与串行不同,结果会有微小差异,不是 bug。

12.3 GPU 浮点

NVIDIA GPU 默认 IEEE 754 兼容。-use_fast_math 选项会牺牲精度换速度(用近似 sin/cos/exp、融合 FMA 等)。


第十三部分:常见陷阱清单

  1. 混用 host/device 指针*d_A 在 CPU 代码里、*h_A 在 kernel 里都会崩
  2. 忘了边界检查 if (idx < N):最后一个 block 的多余线程越界
  3. 启动配置算错:用 (N + tpb - 1) / tpb 上取整,不是 N / tpb
  4. threadsPerBlock 不是 32 的倍数:有死线程,浪费 occupancy
  5. 忘了同步:异步 kernel 还没跑完就读结果
  6. 不检查 CUDA 错误:错误悄无声息,用 CUDA_CHECK
  7. 整数索引溢出:N 大于 21 亿要用 size_t
  8. 二维访存方向错threadIdx.x 应对应连续维度
  9. 依赖多次 cudaMalloc 返回的地址相邻:实际不保证

第十四部分:心智模型总结

每写一个 kernel 前,先回答四个问题:

  1. 数据怎么分? 多少 block、多少线程?每个线程负责哪些数据?
  2. 数据放哪? Global memory 够,还是要搬到 shared memory?
  3. 线程要不要协作? 独立干活(向量加法)还是协作(归约、矩阵乘法)?
  4. 数据怎么进出 GPU? 一次上传多次算,还是反复来回?

这四个问题贯穿所有 CUDA 题目。每道题的不同就在于答案不同。


Author: YANG
Reprint policy: All articles in this blog are used except for special statements CC BY 4.0 reprint polocy. If reproduced, please indicate source YANG !
  TOC