📚 课程目录第4章向量加法核函数实战

4. 向量加法核函数实战

向量加法是理解CUDA线程模型最合适的起点：计算本身足够简单，不会掩盖并行映射的细节；而它又完整包含了主机端内存分配、设备端内存拷贝、核函数启动配置、线程索引映射、边界判断和结果校验这一整套流程。本节以长度 N = 1,000,000 的单精度浮点向量 C = A + B 为例，把CPU循环版本逐步改造成GPU并行版本，并说明如何用 nvprof / Nsight Systems 观察核函数的实际耗时。

一、CPU版本：一个串行循环

CPU版本的逻辑非常直接：对每个下标 i，读取 A[i] 与 B[i]，相加后写回 C[i]。在单线程下，这 100 万次加法按顺序执行，循环变量 i 同时承担“访问哪个元素”和“执行到第几步”两个角色。GPU版本的关键转变在于：把 i 的计算交给线程本身完成，让每个线程只负责一个元素，从而把“循环”变成“并行展开”。

// CPU 版本
void vecAddCPU(const float* A, const float* B, float* C, int N) {
    for (int i = 0; i < N; ++i) {
        C[i] = A[i] + B[i];
    }
}

二、核函数与一维线程索引映射

GPU版本中，每个线程需要自己算出全局下标。CUDA 提供两级内置变量：blockIdx.x 是当前线程块在网格中的编号，threadIdx.x 是线程在块内的编号，blockDim.x 是每块线程数。一维映射公式为：

int i = blockIdx.x * blockDim.x + threadIdx.x;

这个公式的含义是“先数块、再数块内线程”：第 0 块覆盖下标 0～blockDim.x-1，第 1 块紧接着覆盖下一段，依此类推。若取 blockDim.x = 256，则 100 万元素需要 ceil(1000000/256) = 3907 个线程块，最后一个块只有 1000000 - 3906×256 = 64 个有效线程，其余 192 个线程的下标会越界。这正是边界条件必须处理的原因。

三、核函数实现与边界条件

把映射公式和边界判断写进核函数，就得到最基础的一维向量加法实现。注意核函数没有返回值，所有结果通过指针写回全局内存；参数中的 N 必须传入，因为设备端无法自动知道数组长度。

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

// 主机端启动
int threadsPerBlock = 256;
int blocksPerGrid = (N + threadsPerBlock - 1) / threadsPerBlock;
vecAddKernel<<<blocksPerGrid, threadsPerBlock>>>(d_A, d_B, d_C, N);

这里有两个容易出错的细节。第一，blocksPerGrid 必须向上取整，写成 (N + 255) / 256 而不是 N / 256，否则末尾不足一块的元素会被漏算。第二，if (i < N) 不能省略：末块中越界线程若继续访问 A[i]、B[i]，会读到非法地址，轻则结果错误，重则触发 “an illegal memory access was encountered” 并让后续所有 CUDA 调用失败。边界判断是正确性的一部分，不是可选的优化。

四、主机端完整流程与结果校验

一次完整的GPU向量加法包含六个动作：分配主机内存并初始化、分配设备内存、把 A、B 拷贝到设备、启动核函数、把 C 拷回主机、释放设备内存。核函数启动是异步的，若紧接着用 CPU 读取拷贝结果，必须用 cudaMemcpy 的同步语义或 cudaDeviceSynchronize() 等待完成。校验时可将GPU结果与CPU结果逐元素比较，统计最大绝对误差；对单精度加法，误差应为 0（同一运算顺序下结果完全一致），若出现非零偏差，优先检查索引映射和边界条件。

cudaMalloc(&d_A, bytes); cudaMalloc(&d_B, bytes); cudaMalloc(&d_C, bytes);
cudaMemcpy(d_A, h_A, bytes, cudaMemcpyHostToDevice);
cudaMemcpy(d_B, h_B, bytes, cudaMemcpyHostToDevice);
vecAddKernel<<<blocksPerGrid, threadsPerBlock>>>(d_A, d_B, d_C, N);
cudaMemcpy(h_C, d_C, bytes, cudaMemcpyDeviceToHost);
cudaFree(d_A); cudaFree(d_B); cudaFree(d_C);

五、用 nvprof / Nsight Systems 做时间测量

性能对比不能只看“感觉更快”，而要拿到核函数自身的耗时。nvprof 的常用命令是 nvprof ./vecAdd，输出中会列出每个核函数的调用次数、耗时和吞吐；也可用 nvprof --print-gpu-trace ./vecAdd 查看单次启动的时间线。较新版本的 CUDA 已用 Nsight Systems 替代 nvprof，对应命令为 nsys profile -o report ./vecAdd，再用 nsys stats report.nsys-rep 汇总核函数耗时。

测量时要注意三点：其一，首次启动包含上下文初始化开销，应预热若干次再计时；其二，主机端计时需用 cudaEventRecord 记录设备端事件，避免把异步启动的返回时间误当成计算时间；其三，向量加法是访存密集型任务，100 万元素约 12 MB 数据量，核函数耗时通常在几十微秒量级，此时 cudaMemcpy 的传输时间可能占比很高，因此对比时应分别报告“核函数时间”和“含拷贝的端到端时间”，否则容易把 PCIe 传输开销误判为核函数性能。

判断优化是否有效的基准：先确认结果与CPU版本一致，再比较核函数耗时随 blockDim.x 的变化（如 128、256、512）。若耗时基本不变，说明已接近带宽上限；若明显变差，检查是否因块过大导致占用率下降或末块浪费过多。

把上述流程跑通后，读者就拥有了一个可复用的最小GPU程序骨架：索引映射决定“谁算哪个元素”，边界判断保证“不算不该算的”，事件计时与 profiler 回答“到底快了多少”。后续更复杂的核函数，都是在这个骨架上增加计算密度和访存策略。