📚 课程目录第3章 CUDA编程入门

CUDA 编程入门：核函数编写、线程层次结构、内存管理基础

好，咱们正式开始接触 CUDA。说实话，我刚学 CUDA 那会儿最头疼的不是语法，而是思维转换——从串行 CPU 代码跳到并行 GPU 代码，脑子得转个弯。不过别怕，这一节我把最核心的三个东西给你讲透：核函数怎么写、线程怎么排、内存怎么管。这三样拿下了，你就能看懂大部分单细胞分析的 CUDA 加速代码了。

核函数编写——GPU 上跑的第一个函数

核函数就是跑在 GPU 上的函数。它和普通 C++ 函数最大的区别是前面加了个 __global__，而且返回值必须是 void。我个人习惯把核函数命名加个 _kernel 后缀，方便阅读。

看个最简例子：

// 每个线程把一个数组元素加 1
__global__ void add_one_kernel(float* d_data, int n) {
    int idx = blockIdx.x * blockDim.x + threadIdx.x;
    if (idx < n) {
        d_data[idx] += 1.0f;
    }
}

注意那个 if (idx < n) 判断。我在项目里遇到过数组长度不是线程块倍数的情况，漏掉这个判断就会越界写坏内存，debug 能让你抓狂。记住：核函数里每个线程都要做边界检查。

调用核函数用三尖括号语法：

int n = 1000000;
float* d_data;
cudaMalloc(&d_data, n * sizeof(float));
// 启动 256 个线程每块，共 ceil(n/256) 块
int threads_per_block = 256;
int blocks = (n + threads_per_block - 1) / threads_per_block;
add_one_kernel<<<blocks, threads_per_block>>>(d_data, n);
cudaDeviceSynchronize();

核函数执行是异步的，cudaDeviceSynchronize() 可以等 GPU 算完。我在单细胞分析 pipeline 里经常在关键 kernel 后加一个同步，确保数据就绪。

线程层次结构——Grid、Block、Thread 怎么搭

CUDA 的线程组织分三层：Grid（网格） → Block（线程块） → Thread（线程）。你想想看，一个核函数启动时可以指定 grid 和 block 的维度，dim3 类型支持 x/y/z 三个方向。

为什么要有层次？说白了，是为了适配硬件和方便索引。同一个 block 里的线程可以共享数据（用 shared memory），不同 block 之间线程一般只能通过全局内存通信。单细胞分析中，我们通常一个 block 处理一个细胞或一个基因，这样 block 内线程协作读/写共享内存会非常高效。

dim3 grid(32, 16);       // 32*16 = 512 个 block
dim3 block(1024, 1, 1);  // 每个 block 1024 个线程
kernel_name<<<grid, block>>>(args);

计算线程全局索引的公式要烂熟于心：

int tid = blockIdx.x * blockDim.x + threadIdx.x;     // 1D 情况
// 2D 情况（处理图像或矩阵时常用）
int ix = blockIdx.x * blockDim.x + threadIdx.x;
int iy = blockIdx.y * blockDim.y + threadIdx.y;
int linear_idx = iy * width + ix;

我曾经在写空间转录组邻域计算时用了 2D grid，每个线程对应一个 spot 坐标，然后对每个 spot 的邻居做归约。那种场景下线程索引和空间坐标一一对应，代码可读性特别好。

  硬件约束提醒：每个 block 最多 1024 个线程，但受共享内存和寄存器限制实际值可能更小。调试时先用 256 线程/块起步，性能优化时再压榨。

这里我画一张线程层次结构图，把 grid、block、thread 的关系说清楚：

内存管理基础——数据怎么搬到 GPU 上

GPU 不能直接用 CPU 的内存指针。你得显式地在设备上分配内存，然后把数据拷过去。这套流程我刚开始写了不下上百次，闭着眼睛都能敲出来。

基础三件套：

函数	作用	参数形式
`cudaMalloc`	在 GPU 上分配全局内存	`cudaMalloc((void**)&d_ptr, size)`
`cudaMemcpy`	CPU ↔ GPU 或 GPU ↔ GPU 数据拷贝	`cudaMemcpy(dst, src, size, direction)`
`cudaFree`	释放 GPU 内存	`cudaFree(d_ptr)`

一定要检查 cudaMalloc 和 cudaMemcpy 的返回值。我曾经在单细胞分析中分配一个很大的表达矩阵（基因数×细胞数），结果显存不够直接崩了，没有任何报错提示——后来加了错误检查才抓出来。

一个典型的数据搬运模式：

float* h_data = (float*)malloc(n * sizeof(float));  // 主机端数据
// ... 填充 h_data ...

float* d_data;
cudaMalloc(&d_data, n * sizeof(float));             // 分配设备内存
cudaMemcpy(d_data, h_data, n * sizeof(float), 
           cudaMemcpyHostToDevice);                  // 拷贝到 GPU

// 调用核函数运算...

cudaMemcpy(h_data, d_data, n * sizeof(float),
           cudaMemcpyDeviceToHost);                  // 结果拷回 CPU

cudaFree(d_data);
free(h_data);

CUDA 还提供了统一定址内存（Unified Memory，cudaMallocManaged），它能自动在 CPU/GPU 间迁移数据。我在快速原型验证时喜欢用，但生产环境一般还是手动管理，因为性能更可控。

内存类型速查表

除了全局内存，CUDA 还有好几种内存类型，单细胞分析里用得比较多的是共享内存（shared memory）和常量内存（constant memory）。下面这个表帮你快速区分：

内存类型	访问位置	读写权限	生命周期	典型用途
全局内存	所有线程	读/写	程序结束	大规模数据
共享内存	同一 block	读/写	block 结束	block 内协作、归约
常量内存	所有线程	只读	程序结束	核函数参数、查找表
寄存器	单个线程	读/写	线程结束	局部变量

共享内存的延迟比全局内存低一个数量级。我在做单细胞基因表达归约时，每个 block 先把数据拉到共享内存里，然后原子操作累加，速度直接翻了 3 倍。

避坑指南——我当时犯过的错

忘记同步：多 kernel 串联时中间不 cudaDeviceSynchronize()，后一个 kernel 读到的是脏数据。
索引越界：核函数里索引计算超出数组边界，写坏了其他数据。调试时用 cuda-memcheck 跑一下。
块/线程配比不当：block 内线程数不是 32 的倍数（warp 大小）会浪费计算资源。习惯用 128、256 这种。
内存泄漏：每 cudaMalloc 必配 cudaFree，否则显存越用越少。我在分析 pipeline 里统一用 RAII 管理。