快速 Fourier 变换 (FFT)

`ZZ_N` 上的函数空间 对固定的正整数 `N`, 设 `f(k)` 是定义在 `ZZ_N = 0, 1, 2, cdots, N-1` 上的复值函数. 这个函数也可以写成向量形式: `f = (f_0, f_1, cdots, f_(N-1))`. 它的生成函数定义为一个多项式 `F(x) = sum_(0 le k lt N) f_k x^k`. `ZZ_n` 上两个函数 (即向量) 的内积定义为 `(f, g) = sum_(0 le r lt N) f(r) bar(g(r))`.

正交基底 把单位圆周平均分为 `N` 份, 第 `j` 个分点记为 `omega^j = "e"^(2pi"i"j//N)`, 取 `ZZ_N` 上的一组函数 `epsi_j = (1, omega^j, cdots, omega^((N-1)j))`, `quad j in ZZ_N`, 则它们两两正交: `(epsi_j, epsi_k) = { N, if j = k; 0, otherwise :}` 且 `epsi_0, epsi_1, cdots, epsi_(N-1)` 构成 `ZZ_N` 上的函数空间的一组正交基底.

当 `j = k` 时 `(epsi_j, epsi_k) = sum_(0 le r lt N) omega^(j r) omega^(-j r) = N`.
`j != k` 时, 记 `omega^(j-k) = q`, 有 `q != 1` 和 `q^N = 1`. 于是 `(epsi_j, epsi_k) = sum_(0 le r lt N) omega^(j r) omega^(-k r)` `= sum_(0 le r lt N) q^r` `= (1-q^N)/(1-q) = 0`, 最后, 任意 `f = (f_0, f_1, cdots, f_(N-1))` 都能被 `epsi_0, cdots, epsi_(N-1)` 线性表示, 这是因为假设 `f = sum_(r=0)^(N-1) k_r epsi_r`, 则方程的系数行列式为 Vandermonde 行列式, 必有唯一解.

离散 Fourier 变换 `f` 在 `ZZ_N` 上的 Fourier 系数定义为 `a_N^j(f) ==^"简记" a_j` `:= 1/N (f, epsi_j) = 1/N sum_(0 le k lt N) f_k omega^(-j k)`. 向量 `(a_0, a_1, cdots, a_(N-1))` 称为 `f` 的 Fourier 变换, 记为 `a = hat f`. 如果已知 `a`, 要求 `f`, 只需对 `a_j` 施行 Fourier 逆变换: `f_k = sum_(0 le j lt N) a_j omega^(j k)`.

由于 `{epsi_j//sqrt N}_(j=0)^(N-1)` 构成 `ZZ_N` 上函数空间的标准正交基, 有 `f = sum_(0 le j lt N) (f, epsi_j/sqrt N) epsi_j/sqrt N`, 取第 `k` 个分量得 `f_k = sum_(0 le j lt N) 1/N (f, epsi_j) epsi_j(k)` `= sum_(0 le j lt N) a_j omega^(j k)`.

从生成函数的角度看, `f_k` 等于多项式 `A(x) = sum_(0 le j lt N) a_j x^j` 在 `x = omega^k` 处的值, 而 `a_j` 等于多项式 `F(x) = sum_(0 le k lt N) f_k x^k` 在 `x = omega^-j` 处的值除以 `N`.

已知 `hat f = (4, 3, 2, 1)`, 求 `f`.

由 Fourier 逆变换公式, `f_k = 4 * 1^k + 3 * "i"^k + 2 * (-1)^k + 1 * (-"i")^k`, 代入 `k = 0, 1, 2, 3` 得 `f = (10, 2+2"i", 2, 2-2"i")`.

递推公式 设 `N = 2^n`, 若 `f` 是 `ZZ_(2N)` 上的函数, 则 `g(r) := f(2r)`, `h(r) := f(2r+1)` 是 `ZZ_N` 上的函数. 有如下递推公式 `a_(2N)^j(f) = (a_N^j(g) + a_N^j(h) * omega_(2N)^-j) // 2`,
`a_(2N)^(N+j)(f) = (a_N^j(g) - a_N^j(h) * omega_(2N)^-j) // 2`. 这里 `omega_N^j = "e"^(2pi"i"j//N)`. 对于逆变换, 也有完全类似的公式, 只需把上式 `f, a` 的地位对调, 再用 `omega` 的共轭代替 `omega`.

由 `a_(2N)^j(f)`
`= 1/(2N) sum_(0 le r lt N) (f(2r) omega_(2N)^(-2r j) + f(2r+1) omega_(2N)^(-(2r+1)j))`
`= 1/(2N) sum_(0 le r lt N) (g(r) omega_N^(-r j) + h(r) omega_N^(-r j) * omega_(2N)^-j)` 即得结论.

从生成函数的角度看, `F(x)` 是 `2N-1` 次多项式, `G(x) := f_0 + f_2 x + cdots + f_(2N-2) x^(N-1)`,
`H(x) := f_1 + f_3 x + cdots + f_(2N-1) x^(N-1)` 是两个 `N-1` 次多项式, 满足关系 `F(x) = G(x^2) + x H(x^2)`. 代入 `x = omega_(2N)^-j` 有 `F(omega_(2N)^-j) = G(omega_(2N)^(-2j)) + omega_(2N)^-j H(omega_(2N)^(-2j))` `= G(omega_N^-j) + omega_(2N)^-j H(omega_N^-j)`, 上式两边同除以 `2N` 就得到结论. 此外还有 `F(omega_(2N)^(-N-j)) = F(-omega_(2N)^-j)` `= G(omega_N^-j) - omega_(2N)^-j H(omega_N^-j)`.

快速 Fourier 变换 (FFT) 由递推公式知道, 欲求 `ZZ_(2N)` 上的函数 `f` 的 Fourier 变换, 可以将它拆成两个 `ZZ_N` 上的函数 `g, h`, 分别求得 Fourier 变换后, 再合并得到 `hat f`. 由于 `hat f` 共有 `2N` 个分量, 所以合并操作的时间复杂度为 `O(N)`. 设 FFT 的时间复杂度为 `T(N)`, 则 `T(2N) = 2 T(N) + O(N)`, 解得 `T(N) = O(N log N)`.

校验码

CRC 校验

CRC 校验 (Cyclic Redundancy Check, 循环冗余校验) 是用于通信的校验算法. 该算法的输入是一个二进制串 `b` 和 一个 `bbb F_2` 上的 `n` 次既约多项式 `f`, 输出是一个次数小于 `n` 的多项式.

`bbb F_(2,n)[x]` 上的运算

所谓 `bbb F_2` 上的 `n` 次多项式是指, 系数为 0 或 1, 系数上的加法和乘法是模 2 的 `n` 次多项式. `bbb F_2` 上所有次数不超过 `n` 的多项式共有 `2^n` 个, 它们构成一个环 `bbb F_(2, n)[x]`, 其上的乘法为 `(x^3 + x^2 + 1)(x^3 + x + 1)` `= x^6 + x^5 + x^4 + 3x^3 + x^2 + x + 1` `= x^6 + x^5 + x^4 + x^3 + x^2 + x + 1`. 以上过程简记为 1101 × 1011 = 111 1111. 又, `bbb F_(2, n)[x]` 上的加减法恰好是异或运算: 1101 + 1011 = 1101 - 1011 = 0110. 所谓 `f` 是既约的, 是指 `f` 在 `bbb F_2` 上没有根, 这里就是说 0 和 1 都不是 `f` 的根. 这容易办到, 只要保证 `f` 有奇数个系数为 1 的项, 且常数项等于 1 即可.
当我们计算乘法所得的结果超过 `n` 次, 这个结果需要模 `f`, 使它落在 `bbb F_(2, n)[x]` 中. 比如, 当 `n = 4`, 取既约多项式 `f = x^4 + x + 1`, 则上面的结果化为 `x^6 + x^5 + x^4 + x^3 + x^2 + x + 1` `= (x^4+x+1)(x^2+x+1) + x^2 + x` 即, 111 1111 模 `f =` 1 0011 的最终结果是 `x^2 + x =` 110. 这个求余数的过程很重要, 我们用竖式再算一遍 (注意减法就是异或运算, 不需要借位): 
              111 --- 商 x^2 + x + 1 不重要, 可以不写
        ---------
1 0011 / 111 1111
         100 11
         -------
          11 001
          10 011
          -------
           1 0101
           1 0011
           ------
              110 --- 余数 x^2 + x

说明: mul2 表示 `bbb F_2` 上的乘法, mod2 表示 `bbb F_2` 上的取余.
例3 表明, 若 mod2(x, y) = r, 则 mod2(x + r, y) = 0. 这里的 + 代表异或.
例4 的结果与下文的 CRC-16 的例2 结果相同.
例5 和例6 展示了分步求余的技巧: 

 310000          31320000
% 11021          %  11021
-------          --------
   2672              20b5
+  3200
---------
   147200
%   11021
---------
     20b5

310000 = 11021 * q + 2672
31000000 = 11021 * q * 100 + 267200
31320000 = 11021 * q * 100 + 267200 + 320000

故 mod2(31320000, 11021) 等于 mod2(267200 + 320000, 11021)

CRC 的计算

现在可以介绍 CRC 的计算过程了: `CRC(b, f) = x^n b(x) mod f(x) + x^n b(x)`. 将字节串 `b` 视作多项式, 首先将它乘以 `x^n`, 这相当于在串尾补 `n` 个 0. 然后计算余数 `r(x) = x^n b(x) mod f(x)`, 这是一个次数不超过 `n` 的多项式. 最后将 `x^n b(x)` 与 `r(x)` 相加, 即用 `r(x)` 替换串尾的 `n` 个 0 即可.
将数据 `b(x)` 和校验帧 `r(x)` 相连, 一起发送到接收方, 接收方只需验证 `CRC(b, f) mod f(x) -= 0`, 就通过了 CRC 校验. 若余数不为零, 则说明通信传输出现了问题.

数字 `n` 称为 CRC 算法的位宽 (width), 比如 CRC-8 的位宽等于 8, CRC-32 的位宽等于 32 等. 多项式 `f` 称为 CRC 算法的生成多项式. 如 CRC-8 的生成多项式是 `x^8 + x^2 + x + 1`, CRC-32 的生成多项式是 `x^32 + x^26 + x^23 + x^22 + x^16 + x^12 + x^11 + x^10 + x^8 + x^7 + x^5 + x^4 + x^2 + x + 1`. `f` 用二进制表示时, 由于首项系数必为 1, 只记录 `n` 个低次项即可. 如 CRC-8 和 CRC-32 的生成多项式分别记为 0x07 和 0x04C1 1DB7

下面这个计算 CRC-16 的例子来自 CSDN, 其中运用了分步求余的技巧.
例2、例3 表明 bytes 可以分成多段计算, 把上一次的结果作为 init 传入下一次的计算中.
例4 表明把 CRC 附在数据后面一起计算时, 结果为 0.

这里是 CRC 在线计算.

Hamming 码

Hamming 距离
是指一个数据要变为另一个数据所需要翻转的 bit 数目.
冗余位
为避免数据传输过程出错, 每 `n` bit 的数据需要 `r` bit 的冗余位. 它们满足不等式: `2^r ge n + r + 1`. 比如 `n = 7` 时, 有 `r ge 4`.

Hamming 码 规定数据的第 `2^i` 位 (`i = 1, 2, 4, 8, cdots`) 作为校验位, 其余位置作为数据位. 当数据传输时出现 1bit 的错误时, Hamming 码可以发现并修复错误. 当校验位本身出错时, Hamming 码不能正确修复错误. 以下是一个例子: 

data: 1001 1010 // 待传输的数据, 低位在前: 0x59

encode: __1_ 001_ 1010 // 1 所在的位置是 3 7 9 11

parity: 3^7^9^11 = 0110 // 校验位

hamming: 0111 0010 1010 // 从右向左填充校验位

received: 0111 0010 1110 // 假设第 10 位在传输中出错

parity': 3^7^9^10^11 = 1100 // 接收方计算的校验位

parity xor: 1100 ^ 0110 = 1010 // 即十进制的 10, 表明第 10 位在传输中出错