3. CUDA开发环境与第一个科学计算内核

本节把前面建立的 GPU 执行模型落到工程实践上：用 nvcc 把 .cu 源文件编译成可执行程序，用 CUDA 运行时 API 申请显存、拷贝数据、启动内核，并以向量加法为例完成一个可编译、可运行、可验证、可计时的最小科学计算内核。向量加法虽然简单，却完整覆盖了 GPU 编程的四个关键环节：主机与设备的内存分工、数据搬运、并行映射和结果校验，是后续所有复杂内核的骨架。

nvcc 编译流程：一次编译，两套代码

CUDA 源文件（.cu）中同时存在主机代码和设备代码。nvcc 的核心工作是分离二者：主机代码交给系统 C++ 编译器（Linux 下通常是 gcc/g++，Windows 下是 cl.exe），设备代码交给 NVIDIA 的 PTX 汇编器与 cicc 前端，最终由 nvlink 或宿主链接器合并成可执行文件。理解这条链路，才能在报错时判断问题出在哪一层。

典型编译命令如下：

nvcc -O2 -arch=sm_80 vector_add.cu -o vector_add
# 查看编译各阶段产物
nvcc -O2 -arch=sm_80 --ptx vector_add.cu -o vector_add.ptx
nvcc -O2 -arch=sm_80 --cubin vector_add.cu -o vector_add.cubin
# 查看内核资源占用（寄存器、共享内存）
nvcc -O2 -arch=sm_80 --resource-usage vector_add.cu -o vector_add

其中 -arch=sm_80 指定目标计算能力（对应 A100）。若省略，nvcc 会为默认架构生成代码，可能无法在新 GPU 上运行或无法使用新指令。更稳妥的写法是 -gencode arch=compute_80,code=sm_80，可同时指定虚拟架构与真实架构，便于生成 PTX 以支持向前兼容。编译报错时，若错误信息指向 __global__ 函数内部，通常是设备代码问题；若指向普通 C++ 语法，则多半是主机侧问题。

运行时 API：内存分配与拷贝

CUDA 运行时 API 是宿主程序与 GPU 交互的主要接口。最常用的四个函数构成数据搬运闭环：cudaMalloc 在设备端分配显存，cudaMemcpy 在主机与设备之间或设备内部拷贝数据，cudaFree 释放显存，cudaGetLastError 检查调用是否成功。它们都以 cudaError_t 作为返回值，必须检查。

一个常见的错误是只分配设备内存却忘记释放，或在拷贝时方向写错。cudaMemcpy 的第四个参数必须是 cudaMemcpyHostToDevice、cudaMemcpyDeviceToHost 或 cudaMemcpyDeviceToDevice 之一。方向错误不会导致编译失败，但会在运行时返回错误码或产生错误结果，因此错误检查不是可选项。

向量加法的完整实现

下面给出一个可直接编译运行的向量加法内核。约定：向量长度为 N，每个线程处理一个元素，线程索引 i = blockIdx.x * blockDim.x + threadIdx.x，并做边界判断防止越界。

#include <cstdio>
#include <cuda_runtime.h>

__global__ void vectorAdd(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];
}

#define CUDA_CHECK(call) \
    do { \
        cudaError_t err = call; \
        if (err != cudaSuccess) { \
            printf("CUDA error %s at %s:%d\n", \
                   cudaGetErrorString(err), __FILE__, __LINE__); \
            exit(1); \
        } \
    } while (0)

int main() {
    const int N = 1 << 20;          // 约 100 万个元素
    const size_t bytes = N * sizeof(float);
    float *h_a = (float*)malloc(bytes);
    float *h_b = (float*)malloc(bytes);
    float *h_c = (float*)malloc(bytes);
    for (int i = 0; i < N; ++i) { h_a[i] = 1.0f; h_b[i] = 2.0f; }

    float *d_a, *d_b, *d_c;
    CUDA_CHECK(cudaMalloc(&d_a, bytes));
    CUDA_CHECK(cudaMalloc(&d_b, bytes));
    CUDA_CHECK(cudaMalloc(&d_c, bytes));
    CUDA_CHECK(cudaMemcpy(d_a, h_a, bytes, cudaMemcpyHostToDevice));
    CUDA_CHECK(cudaMemcpy(d_b, h_b, bytes, cudaMemcpyHostToDevice));

    int threads = 256;
    int blocks = (N + threads - 1) / threads;
    vectorAdd<<<blocks, threads>>>(d_a, d_b, d_c, N);
    CUDA_CHECK(cudaGetLastError());          // 检查启动是否成功
    CUDA_CHECK(cudaDeviceSynchronize());     // 等待内核完成

    CUDA_CHECK(cudaMemcpy(h_c, d_c, bytes, cudaMemcpyDeviceToHost));
    printf("c[0]=%f c[N-1]=%f\n", h_c[0], h_c[N-1]);

    cudaFree(d_a); cudaFree(d_b); cudaFree(d_c);
    free(h_a); free(h_b); free(h_c);
    return 0;
}

这段代码中有三个值得注意的点。第一，blocks 用向上取整计算，保证覆盖全部 N 个元素；第二，内核中的 if (i < n) 是必需的，因为最后一个 block 可能有多余线程；第三，cudaGetLastError 只能捕获启动配置错误（如 block 数超限），而 cudaDeviceSynchronize 才能捕获内核执行期间的错误。两者配合才能完整定位问题。

错误检查与性能计时

错误检查应遵循“每次 API 调用后立即检查”的原则。上面的 CUDA_CHECK 宏把错误码、错误字符串和出错位置一起打印，比裸写 cudaGetLastError() 更易定位。对于内核启动，还要注意 <<<blocks, threads>>> 本身不返回错误码，必须紧跟一次 cudaGetLastError()。

性能计时使用 CUDA 事件，而不是 CPU 时钟。事件记录在 GPU 时间线上，能准确反映内核执行时间：

cudaEvent_t start, stop;
cudaEventCreate(&start); cudaEventCreate(&stop);
cudaEventRecord(start);
vectorAdd<<<blocks, threads>>>(d_a, d_b, d_c, N);
cudaEventRecord(stop);
cudaEventSynchronize(stop);
float ms = 0;
cudaEventElapsedTime(&ms, start, stop);
printf("kernel time: %.3f ms\n", ms);

若把 N=2^20、threads=256 作为示例参数，则 blocks=4096。在 A100 上，该内核的合理预期是几十微秒量级；若测得毫秒级，通常说明计时包含了首次启动开销或未预热。科学计算中应至少运行一次预热，再取多次平均。此外，可用有效带宽评估性能：读取 a、b 并写入 c，共搬运 3×N×4 字节，除以耗时即得 GB/s，与 GPU 峰值带宽对比即可判断是否接近上限。

把上述代码保存为 vector_add.cu，执行 nvcc -O2 -arch=sm_80 vector_add.cu -o vector_add && ./vector_add，若输出 c[0]=3.000000 c[N-1]=3.000000，说明数据搬运、内核映射和结果回传全部正确。若输出错误值，优先检查 cudaMemcpy 方向、内核边界判断和 cudaDeviceSynchronize 是否遗漏。这个最小闭环建立后，后续的矩阵乘法、归约和卷积都可以沿用同一套“分配—拷贝—启动—校验—计时”骨架，只是把内核逻辑替换为更复杂的计算。