【多项式】快速沃尔什变换 (FWT)
快速沃尔什变换(Fast Walsh Transform, FWT)是一种用于高效计算序列卷积的算法,特别适用于按位运算(如异或、与、或)定义的卷积。它在组合数学、概率论、信号处理和算法竞赛中有着重要应用。
1. 问题引入:什么是按位卷积?
1.1 普通卷积回顾
在学习FFT(快速傅里叶变换)时,我们处理的是加法卷积:
ck=∑i+j=kai⋅bj
c_k = \sum_{i+j=k} a_i \cdot b_j
ck=i+j=k∑ai⋅bj
FFT可以将时间复杂度从 O(n2)O(n^2)O(n2) 优化到 O(nlogn)O(n \log n)O(nlogn)。
1.2 按位卷积
现在考虑以下三种按位运算卷积:
(1) 按位或卷积(OR Convolution)
ck=∑i∣j=kai⋅bj c_k = \sum_{i \mid j = k} a_i \cdot b_j ck=i∣j=k∑ai⋅bj
(2) 按位与卷积(AND Convolution)
ck=∑i&j=kai⋅bj c_k = \sum_{i \& j = k} a_i \cdot b_j ck=i&j=k∑ai⋅bj
(3) 按位异或卷积(XOR Convolution)
ck=∑i⊕j=kai⋅bj c_k = \sum_{i \oplus j = k} a_i \cdot b_j ck=i⊕j=k∑ai⋅bj
这些卷积无法直接用FFT解决,因为按位运算不满足加法的线性性质。FWT就是为了解决这类问题而设计的。
2. FWT 的核心思想
FWT 的核心思想是:
通过线性变换,将按位卷积转换为逐点乘积,类似于FFT将加法卷积转换为逐点乘积。
即:
FWT(C)=FWT(A)∘FWT(B)
\text{FWT}(C) = \text{FWT}(A) \circ \text{FWT}(B)
FWT(C)=FWT(A)∘FWT(B)
其中 ∘\circ∘ 表示逐点乘积。
然后通过逆变换得到结果:
C=IFWT(FWT(A)∘FWT(B))
C = \text{IFWT}(\text{FWT}(A) \circ \text{FWT}(B))
C=IFWT(FWT(A)∘FWT(B))
3. FWT 的数学基础
3.1 变换的本质
FWT 是一种线性变换,可以表示为矩阵乘法:
A^=M⋅A
\hat{A} = M \cdot A
A^=M⋅A
其中 MMM 是一个 n×nn \times nn×n 的变换矩阵(n=2mn = 2^mn=2m,mmm 是位数)。
对于三种运算,MMM 的结构不同。
3.2 分治结构
FWT 利用分治策略,类似于FFT。假设序列长度为 2n2^n2n,我们可以将其分为前半部分和后半部分(对应最高位为0和1)。
4. 三种 FWT 的详细推导
我们以长度 n=2mn = 2^mn=2m 的序列为例。
4.1 按位或卷积 (OR Convolution)
目标:
ck=∑i∣j=kaibj c_k = \sum_{i \mid j = k} a_i b_j ck=i∣j=k∑aibj
变换定义(FWT):
设 AAA 是原序列,FWT(A)\text{FWT}(A)FWT(A) 定义为:
FWT(A)[i]=∑j⊆iA[j]
\text{FWT}(A)[i] = \sum_{j \subseteq i} A[j]
FWT(A)[i]=j⊆i∑A[j]
即子集和。
分治实现:
将序列分为两部分:
- A0A_0A0:最高位为0的部分
- A1A_1A1:最高位为1的部分
则变换后:
- FWT(A)0=FWT(A0)\text{FWT}(A)_0 = \text{FWT}(A_0)FWT(A)0=FWT(A0)
- FWT(A)1=FWT(A0)+FWT(A1)\text{FWT}(A)_1 = \text{FWT}(A_0) + \text{FWT}(A_1)FWT(A)1=FWT(A0)+FWT(A1)
逆变换(IFWT):
- A0=IFWT(FWT(A)0)A_0 = \text{IFWT}(\text{FWT}(A)_0)A0=IFWT(FWT(A)0)
- A1=IFWT(FWT(A)1−FWT(A)0)A_1 = \text{IFWT}(\text{FWT}(A)_1 - \text{FWT}(A)_0)A1=IFWT(FWT(A)1−FWT(A)0)
变换矩阵:
对于 n=2n=2n=2,矩阵为:
M=[1011]
M = \begin{bmatrix}
1 & 0 \\
1 & 1
\end{bmatrix}
M=[1101]
4.2 按位与卷积 (AND Convolution)
目标:
ck=∑i&j=kaibj c_k = \sum_{i \& j = k} a_i b_j ck=i&j=k∑aibj
变换定义:
FWT(A)[i]=∑i⊆jA[j]
\text{FWT}(A)[i] = \sum_{i \subseteq j} A[j]
FWT(A)[i]=i⊆j∑A[j]
即超集和。
分治实现:
- FWT(A)0=FWT(A0)+FWT(A1)\text{FWT}(A)_0 = \text{FWT}(A_0) + \text{FWT}(A_1)FWT(A)0=FWT(A0)+FWT(A1)
- FWT(A)1=FWT(A1)\text{FWT}(A)_1 = \text{FWT}(A_1)FWT(A)1=FWT(A1)
逆变换:
- A0=IFWT(FWT(A)0−FWT(A)1)A_0 = \text{IFWT}(\text{FWT}(A)_0 - \text{FWT}(A)_1)A0=IFWT(FWT(A)0−FWT(A)1)
- A1=IFWT(FWT(A)1)A_1 = \text{IFWT}(\text{FWT}(A)_1)A1=IFWT(FWT(A)1)
变换矩阵:
M=[1101] M = \begin{bmatrix} 1 & 1 \\ 0 & 1 \end{bmatrix} M=[1011]
4.3 按位异或卷积 (XOR Convolution)
这是最复杂也最常用的一种。
目标:
ck=∑i⊕j=kaibj c_k = \sum_{i \oplus j = k} a_i b_j ck=i⊕j=k∑aibj
变换定义(Walsh-Hadamard Transform):
FWT(A)[i]=∑j=0n−1(−1)popcount(i&j)A[j] \text{FWT}(A)[i] = \sum_{j=0}^{n-1} (-1)^{\text{popcount}(i \& j)} A[j] FWT(A)[i]=j=0∑n−1(−1)popcount(i&j)A[j]
分治实现:
- FWT(A)0=FWT(A0)+FWT(A1)\text{FWT}(A)_0 = \text{FWT}(A_0) + \text{FWT}(A_1)FWT(A)0=FWT(A0)+FWT(A1)
- FWT(A)1=FWT(A0)−FWT(A1)\text{FWT}(A)_1 = \text{FWT}(A_0) - \text{FWT}(A_1)FWT(A)1=FWT(A0)−FWT(A1)
逆变换:
- A0=IFWT(FWT(A)0+FWT(A)12)A_0 = \text{IFWT}\left(\frac{\text{FWT}(A)_0 + \text{FWT}(A)_1}{2}\right)A0=IFWT(2FWT(A)0+FWT(A)1)
- A1=IFWT(FWT(A)0−FWT(A)12)A_1 = \text{IFWT}\left(\frac{\text{FWT}(A)_0 - \text{FWT}(A)_1}{2}\right)A1=IFWT(2FWT(A)0−FWT(A)1)
注意:逆变换需要除以2。
变换矩阵:
M=[111−1] M = \begin{bmatrix} 1 & 1 \\ 1 & -1 \end{bmatrix} M=[111−1]
5. C++ 实现
以下是三种 FWT 的完整 C++ 实现:
#include <iostream>
#include <vector>
using namespace std;typedef long long ll;
const int MOD = 998244353; // 常用模数
const int INV2 = (MOD + 1) / 2; // 2的逆元,因为 (MOD+1) % 2 == 0class FWT {
public:// 快速幂(用于求逆元)static ll pow_mod(ll a, ll b) {ll res = 1;while (b) {if (b & 1) res = res * a % MOD;a = a * a % MOD;b >>= 1;}return res;}// 按位或卷积 FWTstatic void fwt_or(vector<ll>& a, bool inv) {int n = a.size();for (int len = 1; 2 * len <= n; len <<= 1) {for (int i = 0; i < n; i += 2 * len) {for (int j = 0; j < len; j++) {ll u = a[i + j];ll v = a[i + j + len];if (!inv) {// 正变换: a[i+j+len] += a[i+j]a[i + j + len] = (u + v) % MOD;} else {// 逆变换: a[i+j+len] -= a[i+j]a[i + j + len] = (v - u + MOD) % MOD;}}}}}// 按位与卷积 FWTstatic void fwt_and(vector<ll>& a, bool inv) {int n = a.size();for (int len = 1; 2 * len <= n; len <<= 1) {for (int i = 0; i < n; i += 2 * len) {for (int j = 0; j < len; j++) {ll u = a[i + j];ll v = a[i + j + len];if (!inv) {// 正变换: a[i+j] += a[i+j+len]a[i + j] = (u + v) % MOD;} else {// 逆变换: a[i+j] -= a[i+j+len]a[i + j] = (u - v + MOD) % MOD;}}}}}// 按位异或卷积 FWTstatic void fwt_xor(vector<ll>& a, bool inv) {int n = a.size();for (int len = 1; 2 * len <= n; len <<= 1) {for (int i = 0; i < n; i += 2 * len) {for (int j = 0; j < len; j++) {ll u = a[i + j];ll v = a[i + j + len];a[i + j] = (u + v) % MOD;a[i + j + len] = (u - v + MOD) % MOD;// 如果是逆变换,则除以2if (inv) {a[i + j] = a[i + j] * INV2 % MOD;a[i + j + len] = a[i + j + len] * INV2 % MOD;}}}}}// 三种卷积的完整接口// A[i] = sum_{j|k=i} a[j] * b[k]static vector<ll> or_convolution(vector<ll> a, vector<ll> b) {fwt_or(a, false);fwt_or(b, false);int n = a.size();for (int i = 0; i < n; i++) {a[i] = a[i] * b[i] % MOD;}fwt_or(a, true);return a;}// A[i] = sum_{j&k=i} a[j] * b[k]static vector<ll> and_convolution(vector<ll> a, vector<ll> b) {fwt_and(a, false);fwt_and(b, false);int n = a.size();for (int i = 0; i < n; i++) {a[i] = a[i] * b[i] % MOD;}fwt_and(a, true);return a;}// A[i] = sum_{j^k=i} a[j] * b[k]static vector<ll> xor_convolution(vector<ll> a, vector<ll> b) {fwt_xor(a, false);fwt_xor(b, false);int n = a.size();for (int i = 0; i < n; i++) {a[i] = a[i] * b[i] % MOD;}fwt_xor(a, true);return a;}
};// 辅助函数:扩展序列到2的幂
vector<ll> extend_to_power_of_two(const vector<ll>& a) {int n = 1;while (n < (int)a.size()) n <<= 1;vector<ll> res(n, 0);for (int i = 0; i < (int)a.size(); i++) {res[i] = a[i];}return res;
}// 测试函数
int main() {// 示例:计算异或卷积vector<ll> a = {1, 2, 3, 4};vector<ll> b = {1, 1, 1, 1};// 扩展到2的幂a = extend_to_power_of_two(a);b = extend_to_power_of_two(b);cout << "Original a: ";for (ll x : a) cout << x << " ";cout << endl;cout << "Original b: ";for (ll x : b) cout << x << " ";cout << endl;// 计算异或卷积vector<ll> c = FWT::xor_convolution(a, b);cout << "XOR Convolution result: ";for (ll x : c) cout << x << " ";cout << endl;// 验证手动计算// c[0] = a[0]*b[0] + a[1]*b[1] + a[2]*b[2] + a[3]*b[3] = 1*1 + 2*1 + 3*1 + 4*1 = 10// c[1] = a[0]*b[1] + a[1]*b[0] + a[2]*b[3] + a[3]*b[2] = 1*1 + 2*1 + 3*1 + 4*1 = 10// c[2] = a[0]*b[2] + a[1]*b[3] + a[2]*b[0] + a[3]*b[1] = 1*1 + 2*1 + 3*1 + 4*1 = 10// c[3] = a[0]*b[3] + a[1]*b[2] + a[2]*b[1] + a[3]*b[0] = 1*1 + 2*1 + 3*1 + 4*1 = 10// 实际上所有位置都是10,因为b全为1return 0;
}
6. 算法复杂度分析
- 时间复杂度:O(nlogn)O(n \log n)O(nlogn),其中 nnn 是序列长度(必须是 2m2^m2m)。
- 空间复杂度:O(n)O(n)O(n)。
与暴力算法 O(n2)O(n^2)O(n2) 相比,FWT 在 nnn 较大时优势明显。
7. 应用场景
-
算法竞赛:
- 计算子集卷积
- 概率DP中状态转移
- 生成函数优化
-
组合数学:
- 莫比乌斯反演的快速实现
- 子集和问题
-
信号处理:
- 沃尔什-阿达玛变换用于图像压缩
-
密码学:
- 分析布尔函数的非线性性质
8. 注意事项
- 序列长度:必须是 2m2^m2m,如果不是,需要补0。
- 模运算:通常在模意义下进行,注意逆元的计算。
- 精度:整数运算不会出现精度问题。
- 变换顺序:先变换,再逐点乘,最后逆变换。
9. 总结
FWT 是处理按位运算卷积的强大工具,其核心思想是线性变换 + 分治。通过将卷积转换为逐点乘积,实现了 O(nlogn)O(n \log n)O(nlogn) 的高效计算。
掌握 FWT 需要理解:
- 三种变换的数学定义
- 分治实现的递推关系
- 正逆变换的对称性
- 模运算下的实现细节