实现：从本原元到查表法)
我第一次接触“伽罗华域”这个概念是因为调一个 Reed-Solomon 纠错模块。当时数据算出来总是差那么几位翻来覆去查不出问题最后才发现是我把有限域里的乘法当成了普通整数乘法去算。这件事给我的教训很直接想在密码学、纠错码或者 AI 加速里不被人绕晕有限域也就是伽罗华域Galois Field的“元素生成”和“运算原理”必须从头亲手推一遍。这篇博客我就用工程视角把这件事讲透它到底是什么、为什么元素个数必须是素数的幂、本原元生成元怎么找、GF(p) 和 GF(2^m) 的运算怎么落地以及怎么用查表法把乘法、除法、逆元全部加速到 O(1)。读完之后你可以直接照着代码搭一套自己的有限域工具而不是停留在只会调库的阶段。1. 有限域到底是个什么东西1.1 一切从钟表算术说起先别被“群环域”这些抽象代数的名词吓住。有限域最贴近生活的模型就是钟表12 小时制里11 点再过 3 小时是 2 点而不是 14 点。这就是一种“模运算”。你把钟面上的数字看作一个集合 {0, 1, 2, ..., 11}定义加法和乘法都在“超过 12 就回绕”的规则下进行。比如模 12 下7 × 9 6363 mod 12 3。集合有限加法和乘法结果也永远落在集合内部这叫“封闭性”。但是钟表算术有个麻烦你没法保证每个非零元素都有乘法逆元。比如模 12 下2 × ? 1 (mod 12) 是无解的因为 2 和 12 不互质。所以 {0...11} 在模 12 运算下只是一个环不是一个域。有限域要求更严格集合里的元素数量必须有限而且每个非零元素都得有乘法逆元。换句话说在这个世界里“除以 a”永远可以等价成“乘以 a 的逆元”不会出现除不下去的情况。1.2 严格定义与封闭性约束严谨地说一个域是一个集合 F上面定义了加法和乘法两种运算满足加法构成交换群乘法在去掉加法单位元之后也构成交换群并且乘法对加法满足分配律。所谓“伽罗华域 GF(q)”就是指元素个数为 q 的有限域。这里有一个非常重要的结论有限域的元素个数只能是某个素数 p 的幂即 q p^m。p 叫域的特征m 叫扩张次数。这个结论不是拍脑袋规定的它来自域的特征理论。任何一个有限域如果往里不断地加 1必定会在某个时刻回到 0。最早回来时经过的次数就是一个素数 p。也就是说p × 1 0 在这个域里成立。再结合向量空间维度理论整个域可以看成是某个基本功为 p 的元素集合上的 m 维向量空间于是元素总数自然是 p^m。这个约束带来一个实际判断标准GF(10)、GF(15) 这种“想当然”的域根本不存在。常见的 GF(2)、GF(256)、GF(2^128) 都有严格的数学依据而 GF(12) 这种域设计从第一行代码开始就是错的。1.3 为什么工程上爱用二进制域工程领域尤其是密码学和纠错编码里最常出现的是 GF(2^m)。原因非常朴素计算机内部数据就是二进制每个 0/1 刚好对应 GF(2) 里的元素。GF(2) 自身只有两个元素 0 和 1规则很简单加法0 0 00 1 11 1 0注意这里 1 1 0不是 2。乘法0 × 0 00 × 1 01 × 1 1。进一步扩展到 GF(2^m)一个字节8 bit就能天然表示 GF(2^8) 里的一个元素。AES 密码算法的 S 盒、Reed-Solomon 纠错码、BCH 码、CRC 校验底层全在跑这套运算。所以你把 GF(2^8) 研究透等于把一大票密码和通信算法的基础课补上了。2. 素域 GF(p) 的构建与元素生成2.1 模运算就是全部构建规则最基础的有限域是 GF(p)其中 p 是素数。它的元素就是 {0, 1, 2, ..., p-1}加法和乘法都定义为模 p 的结果a b (a b) mod p a * b (a * b) mod p为什么 p 必须是素数因为只有模素数的剩余类集合才保证每个非零元素都有逆元。这是数论里“互质才有逆元”的直接推论p 是素数时1 到 p-1 每个数都和 p 互质所以每个非零元素都能找到自己的乘法逆元。举个最简单的例子GF(5) 里2 4 6 mod 5 13 × 4 12 mod 5 22 的逆元是 3因为 2 × 3 6 mod 5 1求逆元的标准算法是扩展欧几里得算法要解 a × x ≡ 1 (mod p)本质就是找整数 x 和 y 满足 a × x p × y 1。在代码里可以用一行递归或者循环搞定。2.2 本原元元素生成的发动机“元素生成”这个词在有限域里有一个专门含义找一组幂运算让一个数的所有幂次恰好把所有非零元素全部遍历一遍。这样的数叫本原元primitive element也叫生成元。如果 g 是本原元那么{g^0, g^1, g^2, ..., g^(p-2)} {1, 2, ..., p-1}这里 g^0 1。因为非零元素一共有 p-1 个所以 g 的幂次周期刚好是 p-1。从费马小定理我们知道任何非零元素 a 都满足 a^(p-1) ≡ 1 (mod p)所以周期一定是 p-1 的因数。而“本原元”的特殊之处就是它的周期恰好达到最大值 p-1一天不多一天不少。找本原元的方法也固定把 p-1 做素因子分解得到一组素因子 q1, q2, ..., qk。对于候选 g只需要检查对所有 qi是否满足g^((p-1)/qi) ≠ 1 (mod p)如果全部成立g 就是本原元。如果任何一个该幂次等于 1说明 g 的阶被压缩在了某个真子群里它只能生成部分元素。我来手算一遍 GF(7) 验证这个过程。p 7p-1 6素因子是 {2, 3}。候选 g 2检查 2^(6/2) 2^3 8 mod 7 1不合格因为 2 的阶其实只有 3。候选 g 33^(6/2) 3^3 27 mod 7 6 ≠ 13^(6/3) 3^2 9 mod 7 2 ≠ 1。全部通过3 是 GF(7) 的本原元。于是生成序列为幂次 i0123453^i mod 7132645这张表就是 GF(7) 的“离散对数表”。任何非零元素 a 都能写成一个唯一的指数 log(a)反过来指数 i 也能还原成元素 3^i。有了这张表乘除法变成了指数加减法这是整个有限域工程优化的起点。2.3 用幂表统一乘除有了本原元GF(7) 里的乘除法可以这样做查表把元素转成指数指数相加或相减再查表还原。比如算 2 × 6log(2) 2log(6) 3指数和 5查表 3^5 5所以 2 × 6 5验证一下2 × 6 12 mod 7 5正确。再算 4 ÷ 5log(4) 4log(5) 5指数差 4 - 5 -1在周期 6 下等价于 5查表 3^5 5所以 4 ÷ 5 5验证4 × 5 的逆元5 的逆元是 35 × 3 15 mod 7 14 × 3 12 mod 7 5正确。这就是“元素生成”和“运算原理”结合的经典套路先通过本原元把所有非零元素按指数顺序生成一遍然后用指数运算统一乘除。3. 扩展域 GF(p^m)多项式才是主角3.1 为什么要扩展域GF(p) 只有 p 个元素但很多场景需要更大的域。比如一个字节 8 bit有 256 种可能你想让每个字节都能参与运算并且保持“域”的性质就得构造 GF(2^8)而不是生造一个 GF(256)256 不是素数。构造扩展域 GF(p^m) 的思路是把元素看成最高次数不超过 m-1 的多项式系数都在 GF(p) 里取。两个元素相加就是按对应系数在 GF(p) 下相加两个元素相乘要把乘积结果对某个“不可约多项式”取模。这个“不可约多项式”有点类似整数里的素数。它在 GF(p) 上不能被分解成两个更低次数的多项式之积。如果选了一个可约多项式当模两个非零多项式乘起来可能变成零多项式域的结构就被破坏了等于在整数环里混进了零因子。3.2 为什么二进制域里多项式可以当数字用工程上最常用的是 GF(2^m)因为系数只有 0 和 1。一个多项式x^7 x^5 x^3 x 1每个“有无 x^i 项”刚好对应二进制位。这个例子写成二进制是 10101011也就是十六进制 0xAB。位 7、5、3、1、0 为 1其余为 0。因此 GF(2^8) 里的一个元素本质是一个最高 8 位的数值。只是它的加法和乘法规则不再等于整数加和整数乘。多项式加法因为系数只在 GF(2) 里1 1 0恰好就是按位异或XOR。多项式乘法按位相乘再相加但相加是 XOR所以相当于“先做普通乘法再把超过 8 位的部分通过模多项式约回来”。3.3 不可约多项式怎么选以 AES 的 0x11B 为例GF(2^8) 最常见的不可约多项式是x^8 x^4 x^3 x 1对应的二进制是 1 0001 1011十六进制就是 0x11B。注意这里多出来的最高次项 x^8 占 9 位所以用 0x11B 表示。AES 标准选中它就是看中了它的性质不可约、本原且系数尽量少硬件实现和软件查表都方便。虽然 GF(2^8) 也可以选别的不可约多项式比如 x^8 x^4 x^3 x^2 1但“域”的结构只取决于你选哪个模多项式。同一个十六进制数 0x57在不同不可约多项式下乘出来的结果可能完全不同。所以做 AES、Reed-Solomon 时必须完整对齐协议规定的那个不可约多项式一个字节都不能错。4. 元素生成与核心运算的落地实现4.1 加法XOR 一锤子买卖GF(2^m) 里加法最简单没有任何条件分支def gf_add(a, b): return a ^ b多项式加法按系数逐位 XOR。因为每一项系数本身就是 1 或 0相同次数的项相遇1 1 0。所以一个字节的加法和乘法逆元在这个意义上很“廉价”你甚至不需要担心进位因为这里根本没有进位。这给硬件设计带来巨大优势。你去看 AES 的 S 盒高速实现加法电路就是一个按位异或门阵列几纳秒就出结果。4.2 乘法多项式相乘再对不可约多项式取模乘法是有限域里最核心、也最容易犯错的地方。以 GF(2^8) 和 AES 不可约多项式 0x11B 为例先看一个标准参考结果0x57 · 0x13 0xFE这个等式在 AES 官方文档里经常出现用来验证实现是否写对。它的原理是把 0x57 和 0x13 都展开成多项式做普通多项式乘法得到的结果次数可能高达 14。之后不断用 x^8 x^4 x^3 x 1 把高于 7 次的项降下来直到最高次数小于 8。软件里最直观的做法是俄罗斯农民乘法也叫 shift-and-adddef gf_mul_bit(a, b, mod0x11B): r 0 while a: if a 1: r ^ b a 1 b 1 if b 0x100: # 超过了第 8 位 b ^ mod # 用不可约多项式约掉 x^8 return r这个算法我很推荐初学者先跑一遍。它逻辑清晰a 的每一位对应要不要把当前的 b 加进结果每次 b 左移等于多项式乘 x一旦 b 出现第 8 位就用 0x11B 做一次模约减。不过这个按位循环在性能敏感场景下不太够看。一次乘法平均要循环 4 次左右虽然常数不大但 AES 一个分组就要做十几轮Reed-Solomon 编码更是海量调用。所以工程上更常上查表法。4.3 查表法用离散对数统一乘除查表法的思路来自第 2 节的幂表。GF(2^8) 的非零元素一共 255 个它们在本原元 g 的幂次下构成一个循环群{g^0, g^1, g^2, ..., g^254} {1, 2, ..., 255}这里注意指数范围是 0 到 254因为非零元素有 255 个。只要选对本原元 g就能把每个非零元素映射到一个指数也就是离散对数 log(a)。于是乘法变成a · b g^(log(a) log(b))除法变成a / b g^(log(a) - log(b))指数相加或相减之后必须对 255 取模因为非零元素群的周期是 255。查表实现可以把指数表 exp 扩展成 512 个元素让 log(a) log(b) 的最大值 508 直接落到扩展索引里省掉一次取模。4.4 逆元与除法指数表带来的便利AES 的 S 盒计算需要求一个字节在 GF(2^8) 里的乘法逆元。如果不用查表法可以跑扩展欧几里得算法也可以直接用费马小定理把逆元算成 a 的 254 次方。但这些都比较慢。有了对数表逆元就是一行代码def gf_inv(a): if a 0: raise ZeroDivisionError(0 没有乘法逆元) return exp[255 - log[a]]为什么是 255 - log[a]因为 g^255 g^0 1所以a^{-1} g^(-log(a)) g^(255 - log(a))这个公式对 log(a) 0 也成立255 - 0 255exp[255] exp[0] 1恰好 a 1 的逆元是 1。除法同样只需要指数差def gf_div(a, b): if b 0: raise ZeroDivisionError(除数不能为 0) if a 0: return 0 return exp[(log[a] - log[b]) % 255]这里% 255对负数要小心。Python 的取模对负数会返回一个非负结果自然没问题但如果你用 C/C(log[a] - log[b]) % 255在指数差为负时可能得到负数索引必须写成(log[a] - log[b] 255) % 255这是跨语言移植最容易踩的坑之一。5. 实操从零搭一套 GF(2^8) 运算表5.1 生成元怎么选0x03 带来的意外在 GF(2^8) 模 0x11B 下很多初学者会下意识认为“生成元是 x也就是 0x02”因为普通多项式里 x 看起来最自然。但实际上 0x02 的周期不是 255它只能生成一部分非零元素不是一个合格的本原元。在 AES 常见实现里本原元通常取 0x03也就是多项式 x 1。为什么 0x03 行而 0x02 不行这只能通过逐个计算乘法阶来验证。0x02^51 会回到 1而 0x03 的幂次能跑遍全部 255 个非零元素。所以在搭建运算表前必须确认两件事不可约多项式确定了吗本文默认 0x11B。本原元确定了吗本文默认 0x03。这两个参数一旦选错后面所有表都是错的。而且这种错误很隐蔽因为表本身看起来完全正常2 的幂次序列也很“有规律”只是别人对表时必然对不上。5.2 Python 代码指数表和对数表有了本原元可以先写一个按位乘法函数用它循环生成本原元的 255 个幂次同时记录每个元素对应的指数MOD 0x11B GEN 0x03 def gf_mul_bit(a, b, modMOD): r 0 while a: if a 1: r ^ b a 1 b 1 if b 0x100: b ^ mod return r exp [0] * 512 log [0] * 256 x 1 for i in range(255): exp[i] x log[x] i x gf_mul_bit(x, GEN) # 验证生成元的阶是否为 255 assert x 1, 生成元不是本原元循环没有在 255 次后闭合 # 扩展指数表让乘法里避免取模 for i in range(255, 512): exp[i] exp[i - 255]这段代码的关键点在于循环里的x gf_mul_bit(x, GEN)每一步都在让 x 乘一个 0x03而不是普通左移。如果你拿 0x02 当生成元用移位加约减确实也能得到一张表但那张表只覆盖 51 阶子群不是完整域。构建完成后def gf_mul(a, b): if a 0 or b 0: return 0 return exp[log[a] log[b]] def gf_inv(a): if a 0: raise ZeroDivisionError(0 没有乘法逆元) return exp[255 - log[a]] def gf_div(a, b): if b 0: raise ZeroDivisionError(除数不能为 0) if a 0: return 0 return exp[(log[a] - log[b]) % 255]5.3 验证正确性的关键检查点表建好之后不要急着跑业务逻辑。先做几组 sanity check健全性检查能省下后面大量排错时间。第一个检查点乘法的已知参考值。拿 0x57 和 0x13 试assert gf_mul(0x57, 0x13) 0xFE这个等式来自 AES 标准文档是验证 GF(2^8) 乘法实现的金标准。如果这一步不对说明不可约多项式、生成元或者表构建有一个地方错了。第二个检查点逆元。0x53 的乘法逆元是 0xCA所以assert gf_mul(0x53, 0xCA) 0x01 assert gf_inv(0x53) 0xCA第三个检查点指数周期闭合。构建表之后验证assert exp[255] 1 assert all(log[exp[i]] i for i in range(255))第四个检查点代数定律抽查。随便抓几组随机数验证交换律、结合律和分配律import random for _ in range(1000): a, b, c random.randrange(256), random.randrange(256), random.randrange(256) assert gf_mul(a, b) gf_mul(b, a) assert gf_mul(gf_mul(a, b), c) gf_mul(a, gf_mul(b, c)) assert gf_mul(a, b ^ c) (gf_mul(a, b) ^ gf_mul(a, c))这些检查全过基本可以放心你的表是对的。6. 避坑指南与调试实录6.1 我踩过的几个典型坑第一个坑把有限域乘法直接当整数乘法。这是所有新手都会犯的错误。0x53 × 0xCA 用普通整数乘法是 0x43AE完全不是 1。有限域乘法不是整数进位乘法而是多项式乘法后再做不可约多项式约减。只要记住“一次进位都不允许出现”很多误解都能消除。第二个坑忽略了不可约多项式约减的时机。按位移乘法里如果 b 左移后超过了 8 位必须立刻对 0x11B 做 XOR 约减。如果等到最后才约减会导致中间结果溢出的位数越堆越多最终结果看起来“有点对”但实际全错。第三个坑对数表里把 0 当成合法输入。离散对数的定义域是非零元素0 没有对数。很多 bug 不是出现在表构建阶段而是出现在边界输入某个数据块里恰好出现一字节的 0直接丢进 log 数组取下标取到 0 之后算出一个看似合理的错值。正确处理方式是在乘法、除法、逆元入口处显式检查 0。第四个坑不同不可约多项式之间混用。同一字节 0x57在 AES 的 0x11B 模多项式下乘 0x13 得到 0xFE但换一个不可约多项式结果就不是 0xFE。如果你做的是标准协议比如 AES、Reed-Solomon协议里规定哪个多项式就用哪个不要为了“省事”自定义。6.2 快速排错 checklist我把排查流程整理成一张表遇到问题按顺序过一遍就能定位现象优先排查点乘法结果跟标准参考值对不上不可约多项式是否是 0x11B按位乘法是否有按位右移左移只用 0x02 当生成元结果覆盖不全换 0x03 或其他本原元验证阶是否为 255exp[255] ≠ 1生成元不是本原元或者迭代次数写错除法结果为负数索引C/C 里% 255对负数的行为改为(log[a]-log[b]255)%2550 参与运算时结果诡异检查乘法/逆元/除法入口是否单独处理 0查表乘法在指数相加超过 255 时越界exp 表是否扩展到 512或者显式% 2556.3 优化方向从查表到硬件友好查表法在 8 位微控制器和普通 CPU 上已经很快但放到更大域比如 GF(2^128)、GF(2^256)时指数表和对数表会爆炸式增长——表的大小是 2^m 级别完全不可行。这时候要用更高级的算法比如分块乘法、多项式基下的 Karatsuba 乘法或者用位切片bitslice技术在通用寄存器里并行计算多个元素。另外如果你在写密码学库还有个细节不能忽略查表法依赖数据索引访问而索引本身和运算值强相关容易成为侧信道攻击的突破口。AES 的查表实现就需要配合随机化掩码或者恒定时间逻辑来防 Cache 时序攻击。这个属于更深一层的话题但至少要知道“快的实现”和“安全且快的实现”之间还有一道坎。7. 一点个人体会把伽罗华域从定义推到一张可以用的运算表这个过程最大的价值不在于你会写一个乘法函数而在于你建立了“域”的直觉。之后再去看 AES S 盒、Reed-Solomon 生成多项式、椭圆曲线点运算会发现它们都在同一套逻辑下运转先选域再找本原元然后用适当的表示方式把元素和运算映射到硬件友好的操作上。我个人经验里最值得刻意练习的是手算一遍小域上的幂表比如 GF(7) 或者 GF(11)。虽然规模小但你能直观看到本原元的幂次如何遍历所有非零元素也能理解为什么指数加减就能等价乘除。有了这个直觉再上 GF(2^8)、GF(2^16) 就是顺水推舟的事了。最后再分享一个小技巧如果你在实现 AES 或者 Reed-Solomon 时检查不出 bug先别急着怀疑协议逻辑回到有限域层做 0x57 · 0x13 0xFE 这种标准验算。有限的域最大的好处就是可以穷举验证花十分钟跑完生成元周期、乘法结合性、分配律基本能把问题范围从“整段代码”缩小到“某个具体表项”。这套排查思路我在不同领域反反复复用从来没有失效过。
有限域运算与GF(2^8)实现：从本原元到查表法
我第一次接触“伽罗华域”这个概念是因为调一个 Reed-Solomon 纠错模块。当时数据算出来总是差那么几位翻来覆去查不出问题最后才发现是我把有限域里的乘法当成了普通整数乘法去算。这件事给我的教训很直接想在密码学、纠错码或者 AI 加速里不被人绕晕有限域也就是伽罗华域Galois Field的“元素生成”和“运算原理”必须从头亲手推一遍。这篇博客我就用工程视角把这件事讲透它到底是什么、为什么元素个数必须是素数的幂、本原元生成元怎么找、GF(p) 和 GF(2^m) 的运算怎么落地以及怎么用查表法把乘法、除法、逆元全部加速到 O(1)。读完之后你可以直接照着代码搭一套自己的有限域工具而不是停留在只会调库的阶段。1. 有限域到底是个什么东西1.1 一切从钟表算术说起先别被“群环域”这些抽象代数的名词吓住。有限域最贴近生活的模型就是钟表12 小时制里11 点再过 3 小时是 2 点而不是 14 点。这就是一种“模运算”。你把钟面上的数字看作一个集合 {0, 1, 2, ..., 11}定义加法和乘法都在“超过 12 就回绕”的规则下进行。比如模 12 下7 × 9 6363 mod 12 3。集合有限加法和乘法结果也永远落在集合内部这叫“封闭性”。但是钟表算术有个麻烦你没法保证每个非零元素都有乘法逆元。比如模 12 下2 × ? 1 (mod 12) 是无解的因为 2 和 12 不互质。所以 {0...11} 在模 12 运算下只是一个环不是一个域。有限域要求更严格集合里的元素数量必须有限而且每个非零元素都得有乘法逆元。换句话说在这个世界里“除以 a”永远可以等价成“乘以 a 的逆元”不会出现除不下去的情况。1.2 严格定义与封闭性约束严谨地说一个域是一个集合 F上面定义了加法和乘法两种运算满足加法构成交换群乘法在去掉加法单位元之后也构成交换群并且乘法对加法满足分配律。所谓“伽罗华域 GF(q)”就是指元素个数为 q 的有限域。这里有一个非常重要的结论有限域的元素个数只能是某个素数 p 的幂即 q p^m。p 叫域的特征m 叫扩张次数。这个结论不是拍脑袋规定的它来自域的特征理论。任何一个有限域如果往里不断地加 1必定会在某个时刻回到 0。最早回来时经过的次数就是一个素数 p。也就是说p × 1 0 在这个域里成立。再结合向量空间维度理论整个域可以看成是某个基本功为 p 的元素集合上的 m 维向量空间于是元素总数自然是 p^m。这个约束带来一个实际判断标准GF(10)、GF(15) 这种“想当然”的域根本不存在。常见的 GF(2)、GF(256)、GF(2^128) 都有严格的数学依据而 GF(12) 这种域设计从第一行代码开始就是错的。1.3 为什么工程上爱用二进制域工程领域尤其是密码学和纠错编码里最常出现的是 GF(2^m)。原因非常朴素计算机内部数据就是二进制每个 0/1 刚好对应 GF(2) 里的元素。GF(2) 自身只有两个元素 0 和 1规则很简单加法0 0 00 1 11 1 0注意这里 1 1 0不是 2。乘法0 × 0 00 × 1 01 × 1 1。进一步扩展到 GF(2^m)一个字节8 bit就能天然表示 GF(2^8) 里的一个元素。AES 密码算法的 S 盒、Reed-Solomon 纠错码、BCH 码、CRC 校验底层全在跑这套运算。所以你把 GF(2^8) 研究透等于把一大票密码和通信算法的基础课补上了。2. 素域 GF(p) 的构建与元素生成2.1 模运算就是全部构建规则最基础的有限域是 GF(p)其中 p 是素数。它的元素就是 {0, 1, 2, ..., p-1}加法和乘法都定义为模 p 的结果a b (a b) mod p a * b (a * b) mod p为什么 p 必须是素数因为只有模素数的剩余类集合才保证每个非零元素都有逆元。这是数论里“互质才有逆元”的直接推论p 是素数时1 到 p-1 每个数都和 p 互质所以每个非零元素都能找到自己的乘法逆元。举个最简单的例子GF(5) 里2 4 6 mod 5 13 × 4 12 mod 5 22 的逆元是 3因为 2 × 3 6 mod 5 1求逆元的标准算法是扩展欧几里得算法要解 a × x ≡ 1 (mod p)本质就是找整数 x 和 y 满足 a × x p × y 1。在代码里可以用一行递归或者循环搞定。2.2 本原元元素生成的发动机“元素生成”这个词在有限域里有一个专门含义找一组幂运算让一个数的所有幂次恰好把所有非零元素全部遍历一遍。这样的数叫本原元primitive element也叫生成元。如果 g 是本原元那么{g^0, g^1, g^2, ..., g^(p-2)} {1, 2, ..., p-1}这里 g^0 1。因为非零元素一共有 p-1 个所以 g 的幂次周期刚好是 p-1。从费马小定理我们知道任何非零元素 a 都满足 a^(p-1) ≡ 1 (mod p)所以周期一定是 p-1 的因数。而“本原元”的特殊之处就是它的周期恰好达到最大值 p-1一天不多一天不少。找本原元的方法也固定把 p-1 做素因子分解得到一组素因子 q1, q2, ..., qk。对于候选 g只需要检查对所有 qi是否满足g^((p-1)/qi) ≠ 1 (mod p)如果全部成立g 就是本原元。如果任何一个该幂次等于 1说明 g 的阶被压缩在了某个真子群里它只能生成部分元素。我来手算一遍 GF(7) 验证这个过程。p 7p-1 6素因子是 {2, 3}。候选 g 2检查 2^(6/2) 2^3 8 mod 7 1不合格因为 2 的阶其实只有 3。候选 g 33^(6/2) 3^3 27 mod 7 6 ≠ 13^(6/3) 3^2 9 mod 7 2 ≠ 1。全部通过3 是 GF(7) 的本原元。于是生成序列为幂次 i0123453^i mod 7132645这张表就是 GF(7) 的“离散对数表”。任何非零元素 a 都能写成一个唯一的指数 log(a)反过来指数 i 也能还原成元素 3^i。有了这张表乘除法变成了指数加减法这是整个有限域工程优化的起点。2.3 用幂表统一乘除有了本原元GF(7) 里的乘除法可以这样做查表把元素转成指数指数相加或相减再查表还原。比如算 2 × 6log(2) 2log(6) 3指数和 5查表 3^5 5所以 2 × 6 5验证一下2 × 6 12 mod 7 5正确。再算 4 ÷ 5log(4) 4log(5) 5指数差 4 - 5 -1在周期 6 下等价于 5查表 3^5 5所以 4 ÷ 5 5验证4 × 5 的逆元5 的逆元是 35 × 3 15 mod 7 14 × 3 12 mod 7 5正确。这就是“元素生成”和“运算原理”结合的经典套路先通过本原元把所有非零元素按指数顺序生成一遍然后用指数运算统一乘除。3. 扩展域 GF(p^m)多项式才是主角3.1 为什么要扩展域GF(p) 只有 p 个元素但很多场景需要更大的域。比如一个字节 8 bit有 256 种可能你想让每个字节都能参与运算并且保持“域”的性质就得构造 GF(2^8)而不是生造一个 GF(256)256 不是素数。构造扩展域 GF(p^m) 的思路是把元素看成最高次数不超过 m-1 的多项式系数都在 GF(p) 里取。两个元素相加就是按对应系数在 GF(p) 下相加两个元素相乘要把乘积结果对某个“不可约多项式”取模。这个“不可约多项式”有点类似整数里的素数。它在 GF(p) 上不能被分解成两个更低次数的多项式之积。如果选了一个可约多项式当模两个非零多项式乘起来可能变成零多项式域的结构就被破坏了等于在整数环里混进了零因子。3.2 为什么二进制域里多项式可以当数字用工程上最常用的是 GF(2^m)因为系数只有 0 和 1。一个多项式x^7 x^5 x^3 x 1每个“有无 x^i 项”刚好对应二进制位。这个例子写成二进制是 10101011也就是十六进制 0xAB。位 7、5、3、1、0 为 1其余为 0。因此 GF(2^8) 里的一个元素本质是一个最高 8 位的数值。只是它的加法和乘法规则不再等于整数加和整数乘。多项式加法因为系数只在 GF(2) 里1 1 0恰好就是按位异或XOR。多项式乘法按位相乘再相加但相加是 XOR所以相当于“先做普通乘法再把超过 8 位的部分通过模多项式约回来”。3.3 不可约多项式怎么选以 AES 的 0x11B 为例GF(2^8) 最常见的不可约多项式是x^8 x^4 x^3 x 1对应的二进制是 1 0001 1011十六进制就是 0x11B。注意这里多出来的最高次项 x^8 占 9 位所以用 0x11B 表示。AES 标准选中它就是看中了它的性质不可约、本原且系数尽量少硬件实现和软件查表都方便。虽然 GF(2^8) 也可以选别的不可约多项式比如 x^8 x^4 x^3 x^2 1但“域”的结构只取决于你选哪个模多项式。同一个十六进制数 0x57在不同不可约多项式下乘出来的结果可能完全不同。所以做 AES、Reed-Solomon 时必须完整对齐协议规定的那个不可约多项式一个字节都不能错。4. 元素生成与核心运算的落地实现4.1 加法XOR 一锤子买卖GF(2^m) 里加法最简单没有任何条件分支def gf_add(a, b): return a ^ b多项式加法按系数逐位 XOR。因为每一项系数本身就是 1 或 0相同次数的项相遇1 1 0。所以一个字节的加法和乘法逆元在这个意义上很“廉价”你甚至不需要担心进位因为这里根本没有进位。这给硬件设计带来巨大优势。你去看 AES 的 S 盒高速实现加法电路就是一个按位异或门阵列几纳秒就出结果。4.2 乘法多项式相乘再对不可约多项式取模乘法是有限域里最核心、也最容易犯错的地方。以 GF(2^8) 和 AES 不可约多项式 0x11B 为例先看一个标准参考结果0x57 · 0x13 0xFE这个等式在 AES 官方文档里经常出现用来验证实现是否写对。它的原理是把 0x57 和 0x13 都展开成多项式做普通多项式乘法得到的结果次数可能高达 14。之后不断用 x^8 x^4 x^3 x 1 把高于 7 次的项降下来直到最高次数小于 8。软件里最直观的做法是俄罗斯农民乘法也叫 shift-and-adddef gf_mul_bit(a, b, mod0x11B): r 0 while a: if a 1: r ^ b a 1 b 1 if b 0x100: # 超过了第 8 位 b ^ mod # 用不可约多项式约掉 x^8 return r这个算法我很推荐初学者先跑一遍。它逻辑清晰a 的每一位对应要不要把当前的 b 加进结果每次 b 左移等于多项式乘 x一旦 b 出现第 8 位就用 0x11B 做一次模约减。不过这个按位循环在性能敏感场景下不太够看。一次乘法平均要循环 4 次左右虽然常数不大但 AES 一个分组就要做十几轮Reed-Solomon 编码更是海量调用。所以工程上更常上查表法。4.3 查表法用离散对数统一乘除查表法的思路来自第 2 节的幂表。GF(2^8) 的非零元素一共 255 个它们在本原元 g 的幂次下构成一个循环群{g^0, g^1, g^2, ..., g^254} {1, 2, ..., 255}这里注意指数范围是 0 到 254因为非零元素有 255 个。只要选对本原元 g就能把每个非零元素映射到一个指数也就是离散对数 log(a)。于是乘法变成a · b g^(log(a) log(b))除法变成a / b g^(log(a) - log(b))指数相加或相减之后必须对 255 取模因为非零元素群的周期是 255。查表实现可以把指数表 exp 扩展成 512 个元素让 log(a) log(b) 的最大值 508 直接落到扩展索引里省掉一次取模。4.4 逆元与除法指数表带来的便利AES 的 S 盒计算需要求一个字节在 GF(2^8) 里的乘法逆元。如果不用查表法可以跑扩展欧几里得算法也可以直接用费马小定理把逆元算成 a 的 254 次方。但这些都比较慢。有了对数表逆元就是一行代码def gf_inv(a): if a 0: raise ZeroDivisionError(0 没有乘法逆元) return exp[255 - log[a]]为什么是 255 - log[a]因为 g^255 g^0 1所以a^{-1} g^(-log(a)) g^(255 - log(a))这个公式对 log(a) 0 也成立255 - 0 255exp[255] exp[0] 1恰好 a 1 的逆元是 1。除法同样只需要指数差def gf_div(a, b): if b 0: raise ZeroDivisionError(除数不能为 0) if a 0: return 0 return exp[(log[a] - log[b]) % 255]这里% 255对负数要小心。Python 的取模对负数会返回一个非负结果自然没问题但如果你用 C/C(log[a] - log[b]) % 255在指数差为负时可能得到负数索引必须写成(log[a] - log[b] 255) % 255这是跨语言移植最容易踩的坑之一。5. 实操从零搭一套 GF(2^8) 运算表5.1 生成元怎么选0x03 带来的意外在 GF(2^8) 模 0x11B 下很多初学者会下意识认为“生成元是 x也就是 0x02”因为普通多项式里 x 看起来最自然。但实际上 0x02 的周期不是 255它只能生成一部分非零元素不是一个合格的本原元。在 AES 常见实现里本原元通常取 0x03也就是多项式 x 1。为什么 0x03 行而 0x02 不行这只能通过逐个计算乘法阶来验证。0x02^51 会回到 1而 0x03 的幂次能跑遍全部 255 个非零元素。所以在搭建运算表前必须确认两件事不可约多项式确定了吗本文默认 0x11B。本原元确定了吗本文默认 0x03。这两个参数一旦选错后面所有表都是错的。而且这种错误很隐蔽因为表本身看起来完全正常2 的幂次序列也很“有规律”只是别人对表时必然对不上。5.2 Python 代码指数表和对数表有了本原元可以先写一个按位乘法函数用它循环生成本原元的 255 个幂次同时记录每个元素对应的指数MOD 0x11B GEN 0x03 def gf_mul_bit(a, b, modMOD): r 0 while a: if a 1: r ^ b a 1 b 1 if b 0x100: b ^ mod return r exp [0] * 512 log [0] * 256 x 1 for i in range(255): exp[i] x log[x] i x gf_mul_bit(x, GEN) # 验证生成元的阶是否为 255 assert x 1, 生成元不是本原元循环没有在 255 次后闭合 # 扩展指数表让乘法里避免取模 for i in range(255, 512): exp[i] exp[i - 255]这段代码的关键点在于循环里的x gf_mul_bit(x, GEN)每一步都在让 x 乘一个 0x03而不是普通左移。如果你拿 0x02 当生成元用移位加约减确实也能得到一张表但那张表只覆盖 51 阶子群不是完整域。构建完成后def gf_mul(a, b): if a 0 or b 0: return 0 return exp[log[a] log[b]] def gf_inv(a): if a 0: raise ZeroDivisionError(0 没有乘法逆元) return exp[255 - log[a]] def gf_div(a, b): if b 0: raise ZeroDivisionError(除数不能为 0) if a 0: return 0 return exp[(log[a] - log[b]) % 255]5.3 验证正确性的关键检查点表建好之后不要急着跑业务逻辑。先做几组 sanity check健全性检查能省下后面大量排错时间。第一个检查点乘法的已知参考值。拿 0x57 和 0x13 试assert gf_mul(0x57, 0x13) 0xFE这个等式来自 AES 标准文档是验证 GF(2^8) 乘法实现的金标准。如果这一步不对说明不可约多项式、生成元或者表构建有一个地方错了。第二个检查点逆元。0x53 的乘法逆元是 0xCA所以assert gf_mul(0x53, 0xCA) 0x01 assert gf_inv(0x53) 0xCA第三个检查点指数周期闭合。构建表之后验证assert exp[255] 1 assert all(log[exp[i]] i for i in range(255))第四个检查点代数定律抽查。随便抓几组随机数验证交换律、结合律和分配律import random for _ in range(1000): a, b, c random.randrange(256), random.randrange(256), random.randrange(256) assert gf_mul(a, b) gf_mul(b, a) assert gf_mul(gf_mul(a, b), c) gf_mul(a, gf_mul(b, c)) assert gf_mul(a, b ^ c) (gf_mul(a, b) ^ gf_mul(a, c))这些检查全过基本可以放心你的表是对的。6. 避坑指南与调试实录6.1 我踩过的几个典型坑第一个坑把有限域乘法直接当整数乘法。这是所有新手都会犯的错误。0x53 × 0xCA 用普通整数乘法是 0x43AE完全不是 1。有限域乘法不是整数进位乘法而是多项式乘法后再做不可约多项式约减。只要记住“一次进位都不允许出现”很多误解都能消除。第二个坑忽略了不可约多项式约减的时机。按位移乘法里如果 b 左移后超过了 8 位必须立刻对 0x11B 做 XOR 约减。如果等到最后才约减会导致中间结果溢出的位数越堆越多最终结果看起来“有点对”但实际全错。第三个坑对数表里把 0 当成合法输入。离散对数的定义域是非零元素0 没有对数。很多 bug 不是出现在表构建阶段而是出现在边界输入某个数据块里恰好出现一字节的 0直接丢进 log 数组取下标取到 0 之后算出一个看似合理的错值。正确处理方式是在乘法、除法、逆元入口处显式检查 0。第四个坑不同不可约多项式之间混用。同一字节 0x57在 AES 的 0x11B 模多项式下乘 0x13 得到 0xFE但换一个不可约多项式结果就不是 0xFE。如果你做的是标准协议比如 AES、Reed-Solomon协议里规定哪个多项式就用哪个不要为了“省事”自定义。6.2 快速排错 checklist我把排查流程整理成一张表遇到问题按顺序过一遍就能定位现象优先排查点乘法结果跟标准参考值对不上不可约多项式是否是 0x11B按位乘法是否有按位右移左移只用 0x02 当生成元结果覆盖不全换 0x03 或其他本原元验证阶是否为 255exp[255] ≠ 1生成元不是本原元或者迭代次数写错除法结果为负数索引C/C 里% 255对负数的行为改为(log[a]-log[b]255)%2550 参与运算时结果诡异检查乘法/逆元/除法入口是否单独处理 0查表乘法在指数相加超过 255 时越界exp 表是否扩展到 512或者显式% 2556.3 优化方向从查表到硬件友好查表法在 8 位微控制器和普通 CPU 上已经很快但放到更大域比如 GF(2^128)、GF(2^256)时指数表和对数表会爆炸式增长——表的大小是 2^m 级别完全不可行。这时候要用更高级的算法比如分块乘法、多项式基下的 Karatsuba 乘法或者用位切片bitslice技术在通用寄存器里并行计算多个元素。另外如果你在写密码学库还有个细节不能忽略查表法依赖数据索引访问而索引本身和运算值强相关容易成为侧信道攻击的突破口。AES 的查表实现就需要配合随机化掩码或者恒定时间逻辑来防 Cache 时序攻击。这个属于更深一层的话题但至少要知道“快的实现”和“安全且快的实现”之间还有一道坎。7. 一点个人体会把伽罗华域从定义推到一张可以用的运算表这个过程最大的价值不在于你会写一个乘法函数而在于你建立了“域”的直觉。之后再去看 AES S 盒、Reed-Solomon 生成多项式、椭圆曲线点运算会发现它们都在同一套逻辑下运转先选域再找本原元然后用适当的表示方式把元素和运算映射到硬件友好的操作上。我个人经验里最值得刻意练习的是手算一遍小域上的幂表比如 GF(7) 或者 GF(11)。虽然规模小但你能直观看到本原元的幂次如何遍历所有非零元素也能理解为什么指数加减就能等价乘除。有了这个直觉再上 GF(2^8)、GF(2^16) 就是顺水推舟的事了。最后再分享一个小技巧如果你在实现 AES 或者 Reed-Solomon 时检查不出 bug先别急着怀疑协议逻辑回到有限域层做 0x57 · 0x13 0xFE 这种标准验算。有限的域最大的好处就是可以穷举验证花十分钟跑完生成元周期、乘法结合性、分配律基本能把问题范围从“整段代码”缩小到“某个具体表项”。这套排查思路我在不同领域反反复复用从来没有失效过。
