做算子库这些年,“泛化能力”这四个字特别容易变成口号。真正落到代码层面,无非是三件事:shape 对不上的时候能不能算、标量和张量混在一起的时候怎么算、dtype 五花八门的时候按什么规则算。我的 ops-math 库在这三块上花的时间远超内核本身,因为它们看着是独立议题,实际纠缠得很深。广播要考虑标量的退化形态,标量要参与类型提升,类型提升又会影响广播后的中间精度。这篇就当踩坑实录,把我从设计到实现、从 CPU 到 GPU kernel 的完整过程捋一遍。想自己写算子库的朋友,或者想知道框架底层到底帮你处理了什么的朋友,这篇应该能给你省点时间。
1. 泛化能力到底泛化在哪里
1.1 一个加法算子的三种用法
我们最早实现add算子时,接口写得很“老实”:只接受两个 shape 完全一致、dtype 完全一致的张量。这当然好实现,循环嵌套一摆,直接开算,两小时就能跑通。但真实调用场景很快就给了一巴掌。
第一种用法是张量加标量,比如给整个图片矩阵每个像素加一个亮度偏移。第二种用法是张量加低维张量,比如给一批向量组成的矩阵加上一个公共向量,shape 从[4, 3]和[1, 3]变成输出[4, 3]。第三种场景更隐蔽,同一个深度学习模型里,前向计算是 float32,反向传播梯度是 float16,两个张量做加法到底按什么精度算?
这三种用法在数学上都能用“逐元素操作”一句话概括,但落到实现层是完全不同的路径。我们最初的代码只能覆盖第一种中的一小部分——同 shape 同 dtype。想补上缺口,有三种补法:给每个算子写三个特化版本;在算子内部做一堆 if 分支处理不同情况;或者抽象出一套通用机制。ops-math 选的是最后一种,也就是标题里“泛化能力”的实质:广播机制管 shape,标量操作管低维输入,类型系统管 dtype 组合。
1.2 三件事其实是一件事
很多开发者会误以为广播、标量、类型系统是三个可以独立开工的模块。我一开始也这么规划,真正动手才发现,它们的耦合点无处不在。
先说广播和标量的关系。把标量看成一个 shape 为空的 0 维张量,它在内存里只有一个值。如果把这个值“复制”到输出张量的每一个位置上,那它的每个维度步长自然都是 0——步长为 0 意味着无论坐标怎么变,都指向同一个内存地址。这跟广播张量里“某维度大小为 1”的处理方式完全一致。换句话说,标量就是广播机制的最极端退化情形。
类型系统和这两个还纠缠得更深。标量携带自己的 dtype,张量也携带 dtype,两者相加前必须先确定结果类型。而类型提升规则一旦定下来,又会影响广播后的 kernel 里到底用什么寄存器和指令。我们最终的开发顺序是:先定类型系统,再做广播机制,最后做标量封装。这个顺序是有道理的,后面详述。先抛开这些纠缠,把广播机制这一层单独讲清楚。
2. 广播机制:不是复制数据,而是重新解释内存
2.1 五条规则就能覆盖所有情况
广播的数学规则很简单,我把它们归纳成五条,开发时反复对照:
- 从右往左对齐两个张量的 shape。
- 对应维度大小相等,或者其中一个是 1,或者其中一个维度不存在。
- 如果某个维度大小为 1,就沿该方向“扩展”到另一方的大小。
- 如果某个维度不存在,按大小为 1 处理。
- 两个维度大小都大于 1 且不相等,直接报错。
举个例子:[4, 1]加[1, 3],右对齐后第一维 4 对 1,扩展到 4;第二维 1 对 3,扩展到 3,最终输出[4, 3]。再比如[5]加[3, 5],缺失维度按 1 处理,输出[3, 5]。这套规则基本是行业通用标准,用起来不需要发明新东西。
但规则的实现细节才是关键:一个[4, 1]的张量加[1, 3],如果真把[4, 1]展开成[4, 3]然后再计算,内存开销是 4 倍,而且没必要。业界通用的做法是“视图广播”,也就是不实际复制数据,只改步长。
2.2 零拷贝广播背后的 stride 算术
理解广播底层机制,必须懂 shape 之外的另一个概念:stride。stride 表示沿某个维度前进一个索引时,在内存中要跳过的元素个数。比如一个连续存储的[4, 3]矩阵,行主序下第一维 stride 是 3,第二维 stride 是 1。
广播视图要做的事情非常暴力:把某个大小为 1 的维度的 stride 直接改成 0。stride 为 0 意味这个维度无论下标是多少,指向的都是同一块内存。这样[4, 1]和[1, 3]不需要任何数据搬运,只需要构造一组新的 strides,就能在 kernel 里像访问[4, 3]一样去读取它。
这就是广播机制的第一个关键认知:广播不是数据的扩展,而是对数据寻址方式的重新解释。我在给团队新人讲的时候常用一个比喻:同样一箱矿泉水,你可以把它说成“4 排每排 1 瓶”,也可以说成“4 排每排 3 瓶”并且允许多次拿到同一瓶。后者没有多生产水,只是规则变了。
2.3 CPU 路径:快速路径和通用路径
CPU 上的实现会比 GPU 稍微简单,因为没有线程组织的问题。但性能上要注意两种路径的区分。
理想情况是输出张量在内存里连续,且输入本身不需要广播,此时直接用扁平索引从 0 遍历到total_elements,每次做两次连续内存读取即可。这条路径编译器能很好地在底层自动向量化。
有广播时就不能这么干了。必须对输出逻辑索引坐标做解算,再通过各自 strides 映射到物理偏移。一个经典的 kernel 长这样:
for (int64_t flat_idx = 0; flat_idx < out_total; ++flat_idx) { int64_t offset_a = 0; int64_t offset_b = 0; int64_t tmp = flat_idx; // 从高维到低维解坐标,逐维累加步长贡献 for (int d = ndim - 1; d >= 0; --d) { int64_t coord = tmp % shape_a[d]; tmp /= shape_a[d]; offset_a += coord * strides_a[d]; offset_b += coord * strides_b[d]; // 注意如果用同一个coord去索引b,隐含了shape对齐逻辑 } out_ptr[flat_idx] = a_ptr[offset_a] + b_ptr[offset_b]; }这里面有两点容易踩坑。第一,坐标解算用的是右边较小的那个 shape 维度,否则越界。第二,两个输入张量的 ndim 可能不同,实现时要做 pad 对齐。我们在 ops-math 里会把缺维度的输入补成一个逻辑 shape 和全 0 stride 数组,统一 ndim,kernel 里就不用再分情况。
2.4 GPU 路径:让每个线程都能正确取数
GPU kernel 和 CPU 在逻辑上很像,只是多了线程组织和寄存器压力。我们的实现思路通常是每个线程负责一个或几个输出元素,通过全局线程编号解算逻辑坐标,再用 strides 映射物理偏移。调优的关键在减少坐标解算的开销。
坐标解算需要做整数除法和取模,在 GPU 上这些操作比乘法慢一个数量级。因此我强烈建议做“广播轴合并”。思路很简单:把相邻的广播维(stride 为 0 或 size 为 1)和相邻的非广播维合并成一组,用乘法次数更少的组合坐标去解算。举个例子,shape 为[1, 2, 1, 5]和[4, 1, 3, 1],合并后实际只剩两组有效轴:第一组由第 0、1、2 维合并,第二组由第 3 维单独构成。这样 kernel 里的取模运算从 4 次降到 2 次,实测整体性能能提升 20% 到 40%,在 kernel 本来就简单、带宽不高的逐元素算子场景里非常可观。
另一件重要的事情是索引宽度。广播后输出元素总量可能非常大,达到几十亿级别。如果你的索引用的是 32 位整数,就会溢出绕回,数据错乱且毫无提示。ops-math 内部从头到尾统一用int64_t计算扁平索引,宁愿牺牲一点寄存器效率也不愿在某个深夜排查一个莫名奇妙的越界问题。
2.5 广播里的三个边界案例
实测中我们被三个边界情况卡过,都值得单独记录。
第一个是 0 长度维度。shape 为[0]和[1]广播,按规则输出是[0],理论上 kernel 一个元素都不该算。但如果循环条件写的不是<而是<=,就会故意访问 out_ptr[0],导致非法写入。这个看似弱智的 bug 在向量化时尤其隐蔽。我在代码里加了一行防御:total_elements == 0时直接 return。
第二个是大 shape 和小 shape 的异步广播。比如[1, 1, 1, 4096]广播到[128, 128, 128, 4096],内存占用巨大但每个输入元素会被复用海量次。这种情况下输出写入会形成严重的 cache 争用,因为不同线程同时往同一块地址附近写。此时“广播视图偷懒”反而变成了性能陷阱,我们后来会对广播输出做写合并,让每个 inactive 线程写连续块。
第三个是索引溢出。前面说的 int64 方案解决了这问题,但还连带着一个隐患:坐标解算过程中临时变量可能溢出。我们的策略是做完整测试覆盖[0xFFFFFFFE]这种极端 shape。
3. 标量操作:当“零维输入”碰到硬件指令
3.1 标量在库里的正确定位
先明确一个概念:标量在算子库里应该是一种独立类型,而不是一个 shape 为空的张量。张量需要管内存生命周期、引用计数、设备归属,而标量只是一个值。把标量硬塞进张量体系,会让底层代码到处检查“这个张量是不是标量”,烦不胜烦。
我们在 ops-math 里定义了一个极简的Scalar:
struct Scalar { DataType dtype; union { bool b; int8_t i8; int16_t i16; int32_t i32; int64_t i64; uint8_t u8; uint16_t u16; uint32_t u32; uint64_t u64; float f32; double f64; } v; };每个算子会提供一个重载,接受Scalar,最终和广播张量走同一个 kernel,只不过程序内部预先知道“这是一个标量输入”,从内存布局上直接规避了张量开销。
3.2 CPU 上别给标量分配缓冲区
我见过不少实现,图省事把标量包装成 1 元素张量,然后走广播路径。这在功能上没错,但性能上是灾难。1 元素张量意味着标量数据在内存里,kernel 每次读取都要做一次 load;如果输出很大,这份标量还会因为 cache 逐出而被反复从主存取,拖慢整个循环。
正确做法是让标量直接参与计算表达式。CPU 端的编译器对常量表达式优化很激进,像out[i] = in[i] * Scalar(2.0f)这种操作,标量会被搬进 SIMD 寄存器,通过向量立即数指令完成。数据根本没机会进 cache,因为压根没落到内存。
3.3 GPU 上标量应该走 uniform 寄存器
GPU 侧类似,但多了一层硬件特殊性。标量作为 kernel 参数传入后,在所有线程里都是同一个值,硬件会把它放在普通寄存器或者 uniform/constant bank。所有线程同时读同一个 uniform 地址时带宽极高,因为硬件做的是 broadcast。
反例是把标量塞进显存里的某个全局缓冲区,然后每个线程__ldg(&scalar_ptr[0])。体验过就知道,这不仅浪费带宽,还可能因为不同线程读同一地址造成不必要的 L2 压力。所以我们内部约定:标量一律作为 kernel 形参,绝不分配张量缓冲区。这条约定甚至写进了开发文档,我见过太多次图方便而犯这个错的案例。
GPU kernel 里标量和非标量输入统一处理有一个技巧:我们可以构造一个“全 stride 为 0”的输入描述,让标量和张量走同一套索引寻址逻辑。这样从代码结构上看,标量只是广播机制的一个特例,不需要单独写一套标量 kernel。
3.4 标量与类型系统交汇处的坑
标量带来的最大麻烦不是性能,而是语义。你在用户层写tensor / 2,这个2在脚本语言里通常是一个整数。如果 tensor 是 float32,我们很容易直接提升成 float32 后运算。如果 tensor 是 int32,int32 / 2按 C++ 语义是整除,结果可能丢掉小数。如果用 Python 的思维,期望的是精确除法。两种行为没法同时满足。
ops-math 的决策是:默认严格跟随 C/C++ 的算术语义。用户想要精确除法,应该显式地把标量转成浮点。这个决策在文档里显著标注,还在 API 层面区分了divide_floor和divide_trunc。宁可对外多暴露几个方法,也不要让同一个算子在不同 dtype 下表现飘忽。
另一个隐蔽问题是默认标量类型。常见调用场景里,用户不会用一个Scalar结构体去传值,而是直接传一个编程语言原生数字。脚本语言里的0.1通常是 double 类型的。double 和 float32 张量做运算时,如果提升规则把结果提升成 double,那么中间计算精度会高,但后半段可能还得转回 float32,造成“提升-回降”的精度损耗。更麻烦的是,一旦标量是 double 而张量是 float16,常见的提升规则会先统一到 double,再计算,再转回 float16,中间过程浪费了寄存器且没有任何收益。我们后续在类型系统章节会专门讨论如何用“公共类型但保留中间精度”的规则来解决。
4. 类型系统:把“能不能加”变成编译期事实
4.1 dtype 基础设施:枚举只是起点
类型系统的第一步是定义完整的 dtype 集合。ops-math 支持的包括 bool、uint8/16/32/64、int8/16/32/64、float16、bfloat16、float32、float64。每个 dtype 至少需要三个信息:枚举值、字节数、类别标志(整型/浮点/无符号/有符号)。
千万别小看“字节数”和“类别标志”这两个看似平凡的东西。kernel 类型分派时需要根据这些信息判断能不能简化寻址、需不需要特殊处理符号扩展。我们早期把这两个属性散落在各处判断,后来集中放进了traits表,每个 dtype 一行:
| dtype | 字节数 | 类别 | 说明 |
|---|---|---|---|
| f32 | 4 | 浮点 | 默认计算类型 |
| f16 | 2 | 浮点 | 存储优化类型 |
| i32 | 4 | 有符号整型 | 常见索引类型 |
| u8 | 1 | 无符号整型 | 图像/掩码场景 |
这看起来简单,但定义 bool 的类别时就出过分歧:bool 是逻辑类型,不能参与常规算术,和 uint8 看似一样但不能直接互操作。我们在traits里单独给 bool 一个is_logical标记,所有算术算子入口先把它挡在外面。
4.2 类型提升规则:一张表还是编译期分支
类型提升规则本身不复杂。基本原则是:整型按位宽升级,无符号和有符号同宽度时结果按无符号处理,整型和浮点混合时按浮点处理,浮点按位宽升级。真正难的是实现方式。
最直观的方案是类型提升查表法。两张输入 dtype 共上百种组合,定义一张二维静态表,每个单元格写结果 dtype。查表法对代码生成友好,改规则也简单,但运行时有查表开销。这点开销在逐元素算子动辄上千万次循环时会被放大,同时也没法让编译器在编译期就确定最终计算类型,进而影响算子内部是否能用常量折叠。
ops-math 的最终选择是编译期类型分派。我们把每个 dtype 映射为一个编译期标签,比如DTypeTag<int32_t>。当两个操作数类型确定后,调用一个promote_t<A, B>模板函数,由模板特化在编译期完成提升。核心代码大约长这样:
template <typename A, typename B> struct PromoteType { using type = std::conditional_t< std::is_floating_point_v<A> || std::is_floating_point_v<B>, PromoteFloat<A, B>, PromoteInt<A, B> >; };这样所有提升判定都在编译期完成,运行时只对应一组已经编译好的指令。代价是模板组合数量膨胀,这点放到下一节讲。
4.3 类型分派:模板实例化和代码膨胀的平衡术
编译期方案最大的敌人是实例化数量。ops-math 有 14 种 dtype,两两组合就有接近 100 个类型对,每个算子都要为每一对生成一套指令。如果 100 个算子都这么搞,编译时间能让人怀疑人生,二进制体积也会变成灾难。
我们的平衡方案分三层。第一层,对于热点算子(add、multiply、scale 等),完整特化所有类型组合,保证最优性能。第二层,对于冷门算子,先做一个轻量类型转换,把两个输入统一提升到某个公共类型,再调用该公共类型对应的一个通用内核。第三层,用一种DTypeVisitor辅助生成分派代码,把手工写的 switch 次数降到最低。
有一类组合我特意不特化:f16 + f16或者bf16 + bf16。这类低精度加低精度在硬件上通常没有直接指令,硬算反而慢。我们把它提升到 f32 计算,再降回原类型,实测比从头用 f16 指令攒一套快得多。类型提升规则虽然在语义上保持不变,但底层对中等精度组合做了优化,这个细节对外部是不可见的。
4.4 提升和广播,先算哪个
类型提升和广播的先后顺序不只是代码风格问题,它会直接影响结果精度。ops-math 的选择是:先做类型提升,再做广播,所有 kernel 使用统一提升后的类型计算。
这个顺序在大多数情况下没有问题,但有一个场景要特别注意:低精度张量做广播时,如果半路提升到高精度,那么输出是有高精度中间结果的。这个“多余”的高精度空间实际上保护了精度不损失。反过来,如果先广播再提升,输入张量里的低精度值已经被复制成海量元素,整数算术里这一步可能产生大量截断错误。
遇到 f16 和标量 0.1 的经典场景时,我们的处理是:先看标量是 double,张量是 f16,查到提升类型是 double,但 kernel 里把 double 中间结果直接转成 float16 参与最终运算,而不是把整个循环都跑成 double。这相当于做了一个“优先级微调”:保精度但不浪费带宽。实现方法是在类型系统层额外提供一个promote_result,在线性 kernel 初始化时查出来,然后传入模板参数。
5. 三块能力耦合后的实战排坑
5.1 坑一:广播轴合并后索引错乱
我们在最初实现广播轴合并优化时,犯过和很多团队一样的错误:合并维度时拿原始 shape 去算新坐标,而不是用合并后的 shape。结果当一个维度大小为 1 的轴被合并进另一个维度后,坐标解算的模运算维度对不上,输出数据随机错乱。
排查这类问题不能靠肉眼观察。我的经验是先写一个只做索引映射不依赖算子的调试内核,输入输出都是假数据,把每个输出位置的 offset 打出来和 CPU 逐元素对照。把对比结果用断言写在测试里,后面回归时能立刻暴露问题。现在 ops-math 的 CI 里有专门的广播索引一致性测试,覆盖各种形状组合,包括那个典型的[1, 2, 1, 5]对[4, 1, 3, 1]。
5.2 坑二:标量默认 double 带来的精度幻觉
另一个高频问题出在用户层默认标量类型。脚本语言里的原生浮点数字通常是 double,而模型里的张量经常是 float32 甚至 float16。用户写一句tensor * 0.1,我们如果把 0.1 按 double 传入 Scalar,那么类型提升结果就是 double。kernel 算完输出后,再被隐式转回 float32。整个过程浪费了 double 的存储和寄存器,还让用户误以为“底层用高精度算了会更准”。
我们的解决手段是提供“严格模式”和“兼容模式”。兼容模式下,标量的 dtype 会同张量统一到同一个“默认计算类型”,比如 float32。严格模式下,遵循标准提升规则,float32 张量乘 double 0.1 结果就是 double。大部分上层框架调用会选择兼容模式,因为能省内存且性能更好,而严格模式留给对数值精度敏感的底层用户。
5.3 坑三:0 长度维度和广播的“巧合”
前面提过[0]与[1]广播输出[0]。更刁钻的组合是[0, 5]和[1, 5]广播输出[0, 5],这时 kernel 会正常跑,但坐标解算时某维度被广播到 0 后,循环体一次都进不去。如果代码里有一条初始化的路径,可能出现“明明没跑任何元素但某些变量还是被更新了”的状态。
我的建议是在广播预处理阶段直接特殊处理:只要输出某个维度大小是 0,整个 kernel 直接不调度。这样能同时规避大量晦涩的边界 bug。为了保险,还在 debug 模式里对输出张量做了预填充,一旦某个输出位没被写入就能在断言中暴露。
5.4 常用的调试三板斧
排查这三块耦合问题时,我有一套固定流程,效率非常高。
第一板斧是统一日志。在算子入口打一条简短的日志:输入 shape、strides、dtype,输出 shape、dtype。配合一个能解析这些日志的小工具,一眼就能看出广播是否如期进行、提升类型是否和预期一致。
第二板斧是 CPU 参考实现对照。GPU 和 CPU kernel 各跑一遍相同的数据,把输出张量的字节序列拉出来做 diff。这一步能覆盖大量索引错位和类型不匹配问题。
第三板斧是精心设计的 golden test 用例。我把所有 dtype 组合和广播形态的交叉测试,用一个参数化测试框架统一生成,每增加一种 dtype 就自动补全测试矩阵。这看起来笨重,但实际是唯一能长期防止回归的办法。
5.5 常见问题速查表
| 问题 | 根因 | 排查与解法 |
|---|---|---|
| 输出数据偶尔错乱 | 广播维度合并后坐标解算模数用错 | 用索引一致性测试对比 offset |
| 标量参与运算后精度异常 | 标量默认 double 与张量类型不一致 | 设置默认计算类型或显式 cast 标量 |
| 空 shape 张量崩溃 | 0 长度维度没有提前返回 | 预处理阶段检测并跳过 kernel |
| 输出元素总量巨大时越界 | 使用了 32 位扁平索引 | 全链路改用 int64_t 索引 |
| GPU kernel 变慢 | 坐标解算取模太多 | 合并广播轴,减少取模次数 |
| 类型组合太多编译时间爆炸 | 所有算子全组合特化 | 冷门算子转公共类型后走通用内核 |
6. 泛化能力给我们的实际回报
6.1 上层接口从三个函数变成一个
改造前我们的 API 里有add_tensor、add_scalar、broadcast_add三个入口,使用者需要仔细辨别自己属于哪种情况。改造后只剩一个add,内部根据输入是否为 Scalar、shape 是否需广播、dtype 是否需提升,自动调度。
这不是代码行数的简单削减。上层框架调用我们库时,不再需要自己拼凑一堆形状检查逻辑。模型导出和部署流程也跟着简化,因为所有组合都能在一个普通调用里完成。某位对接的同事原话是“以前猜 API 要猜半天,现在直接传参就能跑”。
6.2 同一套 kernel 在多种场景复用
从性能上看,统一机制并没有带来额外开销。相反,因为所有路径都最终落到了少数几个模板化 kernel,指令缓存友好度反而变高。GPU 上的实测里,广播合并后性能比之前手动写的特化版本快不少。
6.3 最后的经验总结
做算子库这几年,最大的体感是:泛化能力不是靠堆特例代码实现的,而是靠恰当的抽象,把“不同情况”统一成一个模型。广播机制用 stride 为 0 统一了内存寻址;Scalar 结构统一了标量输入;类型系统用模板元编程统一了 dtype 组合。三者结合,才真正把“逐元素序列”抽象成了一个通用算子族。
如果让我重新来一遍,我会把时间更多花在测试上面。边界情况的覆盖面远比多写几个内核重要。你已经不会在一开始就处理[0]维度加标量这种组合,但迟早会遇到,那时你就会明白,当初在类型分派上多花的那点心思,全都会还回来。