CUDA C++ 不是“C++ 加几个关键字”

一个 .cu 翻译单元同时面对 host 与 device 两套执行环境。它们可以共享模板和类型，但地址空间、可调用函数、异常机制、标准库实现和编译目标并不完全相同。

同一份 .cu
├─ host pass   → 系统 C++ 编译器 → x86-64 / Arm CPU 代码
└─ device pass → CUDA 前端/PTXAS → PTX + GPU SASS

同一个函数名也可能生成：
├─ host 版本：CPU 调用
└─ device 版本：GPU 线程调用

因此，遇到模板报错、重载选择异常或链接错误时，先问：这段代码由哪一侧编译、从哪一侧调用、符号最终应落在哪个二进制中？

函数与变量限定符

最常见的误区是给复杂 C++ 工具函数随手加上 __host__ __device__。这会让函数体、它调用的所有函数以及涉及的类型都接受两套编译约束。更稳妥的做法是把纯计算核心缩小为双端函数，把 I/O、容器和异常留在 host 侧。

struct Affine {
    float a;
    float b;

    __host__ __device__ constexpr float operator()(float x) const {
        return a * x + b;
    }
};

__global__ void transform(float* out, const float* in, int n, Affine op) {
    for (int i = blockIdx.x * blockDim.x + threadIdx.x;
         i < n;
         i += blockDim.x * gridDim.x) {
        out[i] = op(in[i]);
    }
}

指针契约：const、restrict 与对齐

编译器只有在能证明别名关系和对齐时，才容易做重排、向量化和只读路径优化。

__global__ void saxpy(float* __restrict__ y,
                      const float* __restrict__ x,
                      float a,
                      std::size_t n) {
    for (std::size_t i = blockIdx.x * blockDim.x + threadIdx.x;
         i < n;
         i += static_cast<std::size_t>(blockDim.x) * gridDim.x) {
        y[i] = fmaf(a, x[i], y[i]);
    }
}

• const 表示通过该表达式不修改对象，不等于没有其他写入者。
• __restrict__ 承诺相关指针在该作用域内不以冲突方式别名；承诺错误会产生未定义行为。
• reinterpret_cast<float4*> 不会自动让地址满足 16 B 对齐；分配、偏移和长度都必须满足契约。
• 向量化 load/store 是否生成，应检查 SASS 或 Nsight Compute 指标，不能只看源码类型。

标量、向量与低精度类型

低精度的核心问题不是“能否表示输入”，而是整个计算图中的动态范围：输入量化、乘法、累加、归一化、输出回写分别使用什么精度。生产代码通常采用 mixed precision：低精度存储和乘法，高精度累加及关键归约。

浮点结果为什么会变化

浮点加法不满足结合律：

$$(a+b)+c \ne a+(b+c)$$

并行 reduction 改变了求和树；线程调度或原子更新顺序还可能让每次运行的末位略有不同。a*b+c 也可能被融合为一次舍入的 FMA，因此 GPU 与未启用 FMA 的 CPU reference 不一定逐位一致。

误差评估至少同时使用绝对与相对容差：

$$|x-y| \le atol + rtol\cdot|y|$$

bool nearly_equal(float got, float ref,
                  float atol = 1e-6f, float rtol = 1e-5f) {
    return std::abs(got - ref) <= atol + rtol * std::abs(ref);
}

对长 reduction、方差、softmax 等问题，还应采用以下策略之一：

• pairwise/tree reduction，减少误差随长度线性累积；
• 用 float 或 double 累加低精度输入；
• Kahan/Neumaier 补偿求和，用更多指令换稳定性；
• 先减最大值再计算 softmax；
• 对参考结果使用更高精度，并报告最大绝对误差、最大相对误差和 ULP。

快速数学不是免费午餐

--use_fast_math 会改变多个运算的精度、舍入或特殊值语义，并可能把除法、平方根和超越函数替换为更快实现。不要把它作为通用“优化开关”。

正确流程是：

定义业务误差预算
      ↓
只替换热点运算或启用受控编译选项
      ↓
覆盖正常值、极值、NaN、Inf、subnormal 的测试
      ↓
比较端到端质量与吞吐，而不是只看单个 kernel

若只需要某个快速 intrinsic，优先局部使用并写清数值契约，而不是让整个翻译单元都改变语义。

CUDA C++ 内存模型与现代原子

原子操作同时包含三个维度：

1. 操作是否不可分割：避免 read-modify-write 丢失更新。
2. 内存顺序：relaxed、acquire、release 或更强顺序。
3. 线程作用域：thread、block、device、system。

“用了 atomic”不自动意味着其他普通内存写入已经对消费者可见。发布数据时通常需要 release，消费标志时通常需要 acquire；跨 CPU/GPU 还必须确认 system scope 与内存类型支持。

#include <cuda/atomic>

struct Mailbox {
    int payload;
    cuda::atomic<int, cuda::thread_scope_device> ready;
};

__device__ void publish(Mailbox* box, int value) {
    box->payload = value;
    box->ready.store(1, cuda::memory_order_release);
}

__device__ int consume(Mailbox* box) {
    while (box->ready.load(cuda::memory_order_acquire) == 0) {}
    return box->payload;
}

优先使用 libcu++ 的 cuda::atomic、cuda::std::atomic 与同步原语表达意图；传统 atomicAdd 等 intrinsic 仍适合简单计数和归约。任何自旋协议都要证明生产者能够被调度，否则可能死锁或长期占用 SM。

模板、lambda 与编译期分派

模板很适合把 tile 大小、数据类型、算法策略变成编译期参数：

template<int BLOCK_SIZE, class T>
__global__ void scale(T* data, T factor, int n) {
    static_assert(BLOCK_SIZE % 32 == 0);
    int i = blockIdx.x * BLOCK_SIZE + threadIdx.x;
    if (i < n) data[i] *= factor;
}

scale<256><<<(n + 255) / 256, 256>>>(data, 2.0f, n);

收益是常量传播、循环展开和去除分支；代价是实例化数量、编译时间和二进制体积。通常只为基准证明有价值的有限组合实例化，而不是生成所有排列。

设备 lambda、constexpr 和 C++ 标准支持随 Toolkit、host compiler 与语言模式变化。构建系统应显式指定 -std=c++17/20/23，并以当前版本的语言支持表为准。设备代码中优先使用 cuda::std/libcu++，不要假设 host 的 std 实现可在 device 侧调用。

分离编译、设备链接与 LTO

默认 whole-program 模式中，device 调用通常在同一编译单元内解析。跨 .cu 调用设备函数时，需要 relocatable device code：

a.cu ─┐                 ┌─ host link
      ├─ nvcc -dc ──────┤
b.cu ─┘                 └─ nvlink/device link → 合并 device code

常用概念：

• -dc 或 -rdc=true：生成可重定位 device code；
• device link：解析跨翻译单元的 device 符号；
• device LTO：尽可能恢复跨单元优化；
• 静态库必须携带对应的 device object，否则 host 链接成功也可能缺设备符号。

分离编译改善大型工程边界，却可能阻断内联和常量传播。对极致性能的小核心，可将模板实现留在头文件或开启合适的 LTO；对稳定 ABI 的模块，再接受清晰边界带来的少量成本。

资源提示：__launch_bounds__ 与寄存器上限

__global__ __launch_bounds__(256, 2)
void kernel(float* out) {
    // 编译器据此假设 block 不超过 256 线程，并尝试允许每 SM 至少 2 个 block。
}

__launch_bounds__、--maxrregcount、__maxnreg__ 都是在计算与 occupancy 之间施加约束，不是越低越好。压低寄存器可能产生 local-memory spill，比 occupancy 的收益更昂贵。

验证时同时看：

• ptxas -v 报告的 registers、stack、spill load/store；
• 实际 active warps/blocks；
• long scoreboard、依赖 stall 和 local memory 流量；
• kernel 总时间，而非单一 occupancy 百分比。

语言层代码审查清单

• 索引乘法是否在转换到 64 位前已经溢出？
• __host__ __device__ 函数是否误用了 host-only 标准库或异常？
• __restrict__ 与对齐承诺是否真实？
• 低精度输入是否用了足够宽的 accumulator？
• 并行重排后的误差阈值是否按规模设计？
• 原子操作的 scope/order 是否覆盖真正的生产者和消费者？
• 模板组合是否造成代码膨胀？
• 分离编译是否丢失关键内联？
• 资源约束是否以 profile 结果为依据？

练习

1. 对同一组数据实现串行、原子、树形和补偿求和，比较输入长度增大时的误差与时间。
2. 用 float4 重写 copy kernel，分别测试对齐地址、偏移 4 B 地址和尾部非四倍数长度。
3. 将一个模板 device 函数拆到另一 .cu，完成 separable compilation，再比较启用 device LTO 前后的 SASS 和性能。
4. 用 cuda::atomic 实现单生产者/单消费者 mailbox，并说明 acquire/release 各保护了什么。

官方参考

• CUDA Programming Guide: C++ Language Support
• CUDA Programming Guide: C/C++ Language Extensions
• CUDA Programming Guide: Floating-Point Computation
• CUDA Programming Guide: CUDA C++ Memory Model
• CUDA C++ Core Libraries / libcu++
• NVCC: Separate Compilation