1.4.3 CUDA C++ 入门

改 kernel 比写 kernel 更常见。

1. CUDA 编程模型:Grid / Block / Thread 与 nvcc

CUDA 把 GPU 当成一台”高度并行的小型 NUMA 机器”用,执行单元按三层调度。nvcc 是异构编译器,前端是标准 C++,device 代码被分拣出来交给 ptxas 编译成 SASS 机器码。

概念	层级	典型规模	同步域	索引 API
Thread	一个 CUDA core	数万	warp(32)内 SIMT	threadIdx.{x,y,z}
Block	一个 SM 上的协程组	≤1024 / block	__syncthreads()	blockIdx.{x,y,z} + blockDim
Grid	一次 launch 的全部 block	≤2^31 / dim	不同 block 无序	全局 ID:blockIdx*blockDim + threadIdx
Kernel	__global__ 函数	由 launch 决定	launch 间通过 stream	<<<grid, block, smem, stream>>>

// hello.cu —— 编译:nvcc hello.cu -o hello && ./hello
#include <cstdio>
__global__ void hello() {
    printf("block %d thread %d\n", blockIdx.x, threadIdx.x);
}
int main() {
    hello<<<2, 4>>>();
    cudaDeviceSynchronize();        // 否则主机退出,设备输出丢失
}

调研依据:SM(Streaming Multiprocessor) 是硬件调度单位,A100 有 108 个 SM;CUDA 12.x 起 nvcc 支持 --threads N --gpu-architecture=sm_80 多线程编译,CUDA C++ 不再是”另一门语言”而是带 __host__ / __device__ / __global__ 限定符的 C++。

2. 内存层级:Global / Shared / Constant / Local 与合并访问

GPU 内存按”延迟带宽漏斗”排开,优化本质就是让数据待在它该待的地方。合并访问(Coalesced Access)是指同一 warp 的 32 个线程访问同一 cache line 的连续 32×4=128 字节,硬件一次性满足。

内存	位置	作用域	生命周期	典型延迟(cycles)	申请方式
Register	SM 寄存器文件	单 thread	kernel 内	0	局部变量
Shared	SM 片上 SRAM	block 内线程	kernel/block	~20	__shared__ float tile[16][16]
Constant	设备只读缓存	所有 thread	应用周期	~5 (命中)	__constant__ + cudaMemcpyToSymbol
Local	DRAM(逻辑私有)	单 thread	kernel 内	~400	寄存器溢出
Global	DRAM	所有 thread	应用周期	~400	cudaMalloc / 统一内存

__global__ void saxpy(float* y, const float* x, float a, int n) {
    int i = blockIdx.x * blockDim.x + threadIdx.x;  // 合并:y[i], x[i] 连续
    if (i < n) y[i] = a * x[i] + y[i];
}
// 启动:saxpy<<<(n+255)/256, 256>>>(d_y, d_x, 2.0f, n);

调研依据:Ampere 架构(A100)全局内存 HBM2e 带宽 1.5~2 TB/s,共享内存 164 KB/SM;非合并访问会把 32 次请求拆成 32 次事务,带宽直降 32 倍——这是新手 kernel 慢的最常见原因。

3. Kernel 优化的 5 大技巧:occupancy / coalescing / branch / shared / async

改 kernel 比写 kernel 更常见,90% 的性能提升来自这五条。互相制约,需要 nvprof / Nsight Compute 量化后再取舍。

技巧	解决的问题	关键 API / 配置	代价 / 副作用
Occupancy 占用率	SM warp 数太少隐藏不了延迟	--ptxas-options=-v、maxThreadsPerBlock、shared 用量	注册器/spill
Coalescing 合并访问	32 lane 跨 cache line	改访问 stride、用 SoA 替 AoS	重构数据结构
Branch 分支收敛	warp 内 if/else 串行	__ballot_sync、__shfl_xor_sync、predication	代码复杂度
Shared Memory 复用	重复读 global	__shared__ tile[16][16] + __syncthreads()	bank conflict
Async 异步流水线	全局加载与计算重叠	cuda::pipeline、cp.async (sm_80+)	资源占用翻倍

// 分支收敛:32 个 lane 只算偶数索引,避免串行 if
__global__ void evens_only(float* x, int n) {
    int i = blockIdx.x * blockDim.x + threadIdx.x;
    if ((i & 31) < 16) {                 // half-warp 不串行
        x[i] *= 2.0f;
    }
}
// bank conflict: shared[16][32] 列访问 → 改 padding [16][32+1]

调研依据:NVIDIA 在 Hopper 引入 Thread Block Cluster 与 Distributed Shared Memory;Volta 起独立线程调度(Independent Thread Scheduling)让 warp 内分支不再像 Maxwell 那样锁步,这是理解 __syncwarp() 必要性的前提。

4. 矩阵乘法从 naive 到 50 倍:一个完整优化案例

N×N 矩阵乘 C = A × B 是 CUDA 教学标杆,naive 版本在 A100 上只能跑到 ~1 TFLOPS,优化版可达 19.5 TFLOPS(A100 FP16 峰值 312 TFLOPS 的 ~6%)。下面以 FP32 为例,跑过 1024×1024 案例。

版本	核心改动	关键参数	1024² 实测 GFLOPS	加速比
v1 naive	一个 thread 算一元素	<<<N, 1>>>	~2	1×
v2 tiled	16×16 tile + shared mem	blockDim=(16,16)	~30	15×
v3 + register	每 thread 算 4×4 子块	tile[4][4] 寄存器	~80	40×
v4 + vectorized	float4 加载 + 边界处理	reinterpret_cast<float4*>	~95	47×
v5 cuBLAS	cublasSgemm	内核工厂实现	~14000	~50×(理论上限)

// v3 register tile 骨架:每线程算 4×4,block 16×16 → 64×64 输出/launch
constexpr int BM = 16, BN = 16, BK = 16, TM = 4, TN = 4;
__global__ void sgemm_v3(const float* A, const float* B, float* C, int N) {
    __shared__ float sA[BK][BM], sB[BK][BN];
    int bx = blockIdx.x * BN, by = blockIdx.y * BM;
    int tx = threadIdx.x, ty = threadIdx.y;
    float reg[TM][TN] = {{0}};
    for (int k = 0; k < N; k += BK) {
        sA[tx][ty] = A[(by+ty)*N + k+tx];
        sB[tx][ty] = B[(k+tx)*N + bx+ty];
        __syncthreads();
        #pragma unroll
        for (int i = 0; i < BK; ++i)
            #pragma unroll
            for (int m = 0; m < TM; ++m)
                #pragma unroll
                for (int n = 0; n < TN; ++n)
                    reg[m][n] += sA[i][ty*TM+m] * sB[i][tx*TN+n];
        __syncthreads();
    }
    // 写回 C[by..by+TM][bx..bx+TN]
}

调研依据:Scott Gray 等 2012 年在博客《CUDA Pro Tip: Write a SGEMM kernel》给出上述 5 步路径,被 NVIDIA 官方培训引用;v5 cuBLAS 用 Tensor Core + 双缓冲 + warp-specialization,实测是 v1 的 6000 倍而非 50 倍——50× 是”手写 FP32 不用 Tensor Core”的天花板。

5. 现代 AI 推理的 CUDA 集成:Triton / vLLM / Flash Attention

直接写 CUDA 的场景在缩小——推理框架把 kernel 抽象掉,改用更高层的 DSL。但性能护城河依然是 CUDA:这几层工具最终还是落到 nvcc + cuBLAS/cuDNN/cutlass。

工具	层级	用途	与 CUDA 的关系	学习曲线
Triton	Python DSL → PTX	自定义 elementwise/reduction/matmul kernel	OpenAI 出品,直接生成 PTX	低,NumPy 风
vLLM	推理引擎	LLM 服务的 PagedAttention + continuous batching	内部用 CUDA + flash-attn 库	中,Python 部署
Flash Attention-2	CUDA kernel 库	O(N²) attention 但 IO 极低(HBM 读写)	Tri Dao 开源,flash_attn Python 包	高,需懂 attention
CUTLASS	C++ 模板	可组合的 GEMM/卷积 template	NVIDIA 官方,生成 SASS	高,C++ template
cuDNN	闭源库	训练/推理标准算子(conv/attention/norm)	厂商优化	无,API 调用

# Triton 版向量化 add —— 编译时自动生成 PTX
import triton, triton.language as tl
@triton.jit
def add_kernel(x_ptr, y_ptr, out_ptr, n, BLOCK: tl.constexpr):
    pid = tl.program_id(0)
    offs = pid * BLOCK + tl.arange(0, BLOCK)
    mask = offs < n
    x = tl.load(x_ptr + offs, mask=mask)
    y = tl.load(y_ptr + offs, mask=mask)
    tl.store(out_ptr + offs, x + y, mask=mask)
# 比手写 cudaMemcpy + kernel 慢约 1.2×,但代码量 1/10

调研依据:Flash Attention-2 论文(Tri Dao, 2023)在 A100 上把 GPT-2 训练加速 3×;vLLM 论文(Kwon et al., SOSP’23)用 PagedAttention 把 throughput 提到 HuggingFace 的 14–24×;Triton 2.0 引入 tl.dot 直接调用 Tensor Core,使自定义 kernel 也能达到 cuBLAS 80% 性能,门槛进一步降低。

说明 · 本站内容均为学习笔记与经验总结,所有菜谱与技法请结合实际食材、季节与个人口味灵活调整。涉及生食、营养与健康的内容仅供参考,特殊体质或疾病请咨询专业营养师/医生。