位运算卷积
我的 FWT 因为实现原因可能爆负数,但是加个取模是对的喵。
概念
位运算卷积(FWT)
加法卷积的 “加法” 意为第 \(i\) 项和第 \(j\) 项的结果贡献到第 \(i + j\) 项。
类似地,做位运算卷积时多项式的第 \(i\) 项和第 \(j\) 项贡献到第 \(i \oplus j\) 项,其中 \(\oplus\) 是某一位运算。
- 形式化定义:\(S[k] = \sum\limits_{i \oplus j = k} A[i] B[j]\),记作 \(A * B = S\).
位矩阵
FFT 的作用是把多项式转化成点值表示之后,将卷积直接转化成点积。类似地,FWT 的作用也是将位运算卷积转化成直接点积。
设 \(FWT(A)\) 表示幂级数 \(A\) 经过 FWT 变换之后得到的幂级数。
FWT 和 FFT 都是线性变换,考虑设变换系数 \(c(i, j)\),表示 \(A[j]\) 对 \(FWT(A)[i]\) 的贡献系数,也就是 \(FWT(A)[i] = \sum\limits_{j = 0}^{n - 1} c(i, j) A[j]\).
将上式代入 \(FWT(A) \cdot FWT(B) = FWT(C)\),化简可得 \(\sum\limits_{j = 0}^{n - 1} \sum\limits_{k = 0}^{n - 1} c(i, j) c(i, k) A[j] B[k] = \sum\limits_{p = 0}^{n - 1} c(i, p) C[p]\).
又因为 \(A * B = C\),可以得到 \(C[p] = \sum\limits_{t1 \oplus t2 = p} A[t1] B[t2]\).
将上面二式代入分别代入 \(\sum\limits_{p = 0}^{n - 1} c(i, p) C[p]\),可以得到 \(\sum\limits_{j = 0} \sum\limits_{k = 0}^{n - 1} c(i, j) c(i, k) A[j] B[k] = \sum\limits_{p = 0}^{n - 1} \sum\limits_{t1 \oplus t2 = p} A[t1] B[t2] c(i, t1 \oplus t2) = \sum\limits_{t1 = 0}^{n - 1} \sum\limits_{t2 = 0}^{n - 1} A[t1] B[t2] c(i, t1 \oplus t2)\).
对比可得 \(c(i, j) c(i, k) = c(i, j \oplus k)\),只需要考虑构造出符合条件的 \(c\) 即可。
又因为位运算具有独立性,所以可以拆位考虑,只需要求出满足条件的 \(c[0, 1], [0, 1]\)。
\(c\) 可以写作矩阵的形式,称为位矩阵。
一些常用的位矩阵:
-
或卷积:\(\begin{aligned} \begin{bmatrix} 1 &0 \\ 1 &1 \end{bmatrix} \end{aligned}\)
逆矩阵:\(\begin{aligned} \begin{bmatrix} 1 &0 \\ -1 &1 \end{bmatrix} \end{aligned}\)
-
与卷积:\(\begin{aligned} \begin{bmatrix} 1 &1 \\ 0 &1 \end{bmatrix} \end{aligned}\)
逆矩阵:\(\begin{aligned} \begin{bmatrix} 1 &-1 \\ 0 &1 \end{bmatrix} \end{aligned}\)
-
异或卷积:\(\begin{aligned} \begin{bmatrix} 1 &1 \\ 1 &-1 \end{bmatrix} \end{aligned}\)
逆矩阵:\(\begin{aligned} \begin{bmatrix} \frac{1}{2} &\frac{1}{2} \\ \frac{1}{2} &-\frac{1}{2} \end{bmatrix} \end{aligned}\)
FWT 变换
类似 FFT,考虑分类,按位拆半。
\(FWT(A)[i] = \sum\limits_{j = 0}^{n - 1} c(i, j) A[j] = \sum\limits_{j = 0}^{\lfloor \frac{n}{2} \rfloor - 1} c(i, j) A[j] + \sum\limits_{j = \lfloor \frac{n}{2} \rfloor}^{n - 1} c(i, j) A[j]\).
令 \(i^{\prime}\) 为 \(i\) 去掉二进制下首位剩下的数,\(i_0\) 为 \(i\) 二进制表示中最高位上的数(此处均钦定二进制表示的最高位)
根据位矩阵的定义可得原式等价于:
\(\sum\limits_{j = 0}^{\lfloor \frac{n}{2} \rfloor - 1} c(i_0, j_0) c(i^{\prime}, j^{\prime}) A[j] + \sum\limits_{j = \lfloor \frac{n}{2} \rfloor}^{n - 1} c(i_0, j_0) c(i^{\prime}, j^{\prime}) A[j]\).
根据按位拆半可得:
\(\sum\limits_{j = 0}^{\lfloor \frac{n}{2} \rfloor - 1} c(i_0, 0) c(i^{\prime}, j^{\prime}) A[j] + \sum\limits_{j = \lfloor \frac{n}{2} \rfloor}^{n - 1} c(i_0, 1) c(i^{\prime}, j^{\prime}) A[j]\).
注意到 \(c(i^{\prime}, j^{\prime})\) 较 \(c(i, j)\) 规模减半。
令 \(A_0\) 为幂级数下标的二进制表示最高位为 \(0\) 的部分,\(A_1\) 同理,则:
当 \(i_0 = 0\) 时:
\(FWT(A)[i] = c(0, 0) FWT(A_0)[i] + c(0, 1) FWT(A_1)[i]\).
当 \(i_0 = 1\) 时:
\(FWT(A)[i + \frac{n}{2}] = c(1, 0) FWT(A_0)[i] + c(1, 1) FWT(A_1)[i]\).
于是根据上式,可以以 \(O(n)\) 的代价合并两个长度为 \(\frac{n}{2}\) 的子变换。
当 \(n = 2^m\) 时,复杂度为 \(O(m 2^m)\)。看起来很大,但实际上是 \(O(n \log n)\) 级别的 /xk
逆变换 IFWT 就是对 \(c\) 求逆再进行变换。
做一次或卷积相当于求子集信息和,做一次与卷积相当于求超集信息和。
FWT 在一定程度上是可以代替 FMT 和高维前缀和的。
下面模板图一乐就行,正常不会这么写。
P4717 【模板】快速莫比乌斯/沃尔什变换 (FMT/FMT)
#include <cstdio>
#include <cstring>
using namespace std;
typedef long long ll;
const int maxn = 1.5e5 + 5;
const int mod = 998244353;
const int inv2 = 499122177;
const ll
Cor[2][2] = {{1, 0}, {1, 1}},
Cand[2][2] = {{1, 1}, {0, 1}},
Cxor[2][2] = {{1, 1}, {1, mod - 1}},
ICor[2][2] = {{1, 0}, {mod - 1, 1}},
ICand[2][2] = {{1, mod - 1}, {0, 1}},
ICxor[2][2] = {{inv2, inv2}, {inv2, mod - inv2}};
int n;
ll f[maxn], g[maxn], F[maxn], G[maxn];
void FWT(ll *F, const ll C[2][2], int n)
{
for (int len = 2, m = 1; len <= n; m = len, len <<= 1)
{
for (int l = 0, r = len - 1; r <= n; l += len, r += len)
{
for (int p = l; p < l + m; p++)
{
ll val = F[p];
F[p] = (C[0][0] * F[p] + C[0][1] * F[p + m]) % mod;
F[p + m] = (C[1][0] * val + C[1][1] * F[p + m]) % mod;
}
}
}
}
void mul(ll *F, ll *G, const ll C[2][2], const ll IC[2][2], int n)
{
FWT(F, C, n), FWT(G, C, n);
for (int i = 0; i < n; i++) F[i] = F[i] * G[i] % mod;
FWT(F, IC, n);
}
int main()
{
scanf("%d", &n), n = (1 << n);
for (int i = 0; i < n; i++) scanf("%lld", &F[i]);
for (int i = 0; i < n; i++) scanf("%lld", &G[i]);
memcpy(f, F, n * sizeof(ll)), memcpy(g, G, n * sizeof(ll));
mul(f, g, Cor, ICor, n);
for (int i = 0; i < n; i++) printf("%lld ", f[i]); putchar('\n');
memcpy(f, F, n * sizeof(ll)), memcpy(g, G, n * sizeof(ll));
mul(f, g, Cand, ICand, n);
for (int i = 0; i < n; i++) printf("%lld ", f[i]); putchar('\n');
memcpy(f, F, n * sizeof(ll)), memcpy(g, G, n * sizeof(ll));
mul(f, g, Cxor, ICxor, n);
for (int i = 0; i < n; i++) printf("%lld ", f[i]); putchar('\n');
return 0;
}
子集卷积
求 \(c_k = \sum\limits_{i \& j = 0, i \mid j = k} a_i b_j\).
应该可以用来优化一类状压 dp.
令 \(|i|\) 表示 \(i\) 的二进制表示中 \(1\) 的个数。\(i \& j = 0, i \mid j = k\) 等价于 \(i \mid j = s, |i| + |j| = |s|\).
考虑再开一维记录 \(|i|\),对每一维分别卷积就行。
时间复杂度 \(O(m^2 2^m)\)
#include <cstdio>
#include <cstring>
using namespace std;
#define pcnt __builtin_popcount
const int maxn = 1.1e6 + 5;
const int lg_sz = 21;
const int mod = 1e9 + 9;
const int Cor[2][2] = {{1, 0}, {1, 1}}, ICor[2][2] = {{1, 0}, {mod - 1, 1}};
int n;
int f[lg_sz][maxn], g[lg_sz][maxn], h[lg_sz][maxn];
void FWT(int *F, const int C[2][2], int n)
{
for (int len = 2, m = 1; len <= n; m = len, len <<= 1)
{
for (int l = 0, r = len - 1; r <= n; l += len, r += len)
{
for (int p = l; p < l + m; p++)
{
int val = F[p];
F[p] = (C[0][0] * F[p] + C[0][1] * F[p + m]) % mod;
F[p + m] = (1ll * C[1][0] * val + C[1][1] * F[p + m]) % mod;
}
}
}
}
int main()
{
int lim;
scanf("%d", &lim), n = (1 << lim);
for (int i = 0; i < n; i++) scanf("%d", &f[pcnt(i)][i]);
for (int i = 0; i < n; i++) scanf("%d", &g[pcnt(i)][i]);
for (int i = 0; i <= lim; i++) FWT(f[i], Cor, n), FWT(g[i], Cor, n);
for (int i = 0; i <= lim; i++)
for (int j = 0; j <= i; j++)
for (int k = 0; k < n; k++)
h[i][k] = (h[i][k] + 1ll * f[j][k] * g[i - j][k] % mod) % mod;
for (int i = 0; i <= lim; i++) FWT(h[i], ICor, n);
for (int i = 0; i < n; i++) printf("%d ", h[pcnt(i)][i]);
putchar('\n');
return 0;
}
FWT 快速幂
在 FWT 意义下求 \(F^t(x)\).
考虑先进行一次 FWT,然后对于每个点值直接取其 \(t\) 次方,最后 IWFT 回来就行。
例题:[BZOJ4589] Hard Nim