FWT学习笔记

---

FWT

快速沃尔什变换，用来解决位运算相关的卷积。常见的有与、或、异或三种。

P4717 【模板】快速莫比乌斯/沃尔什变换 (FMT/FWT)#

给定长度为 $2^n$ 两个序列 $A,B$，设

$$C_i=\sum _{j\oplus k=i}{A}_j\times B_k$$

分别当 $\oplus$ 是 or,and,xor 时求出 $C$。

$n\le 17$。

$\mathrm{merge}(A,B)$ 表示将数列 $B$ 拼接到数列 $A$ 后面形成的新数列。

考虑构造一种可逆的线性变换，记 $A$ 变换后为 $\hat{A}$，那么满足

$$C=A\ast B\iff \hat{C}=\hat{A}\times \hat{B}$$

其中 $\times$ 号表示对应位相乘，$\ast$ 号是卷积。然后求出 $\hat{C}$ 的逆变换即可。

因为变换是线性的，所以 $C=A+B\iff \hat{C}=\hat{A}+\hat{B}$。

然后考虑如何构造这个变换。

or#

设现在数列长度为 $2^n$，$A_0$ 表示 $A$ 的前 $2^{(n-1)}$ 个数，其下标二进制下第 $n$ 位为 $0$。$A_1$ 表示 $A$ 的后 $2^{(n-1)}$ 个数，其下标二进制下第 $n$ 位为 $1$。$B_0,B_1,C_0,C_1$ 类似。

显然 $C_0$ 只需要 $A_0,B_0$ 的信息，但是 $C_1$ 需要 $A_0,A_1,B_0,B_1$ 的信息。不妨构造 $\hat{A}=\mathrm{merge}({\hat{A}}_0,{\hat{A}}_0+{\hat{A}}_1)$，当序列长度为 $1$ 时 $\widehat{A_i}=A_i$。证明这样变换是正确的：

因为有 $C_0=A_0\ast B_0,C_1=A_0\ast B_1+A_1\ast B_0+A_1\ast B_1$，

所以有 ${\hat{C}}_0={\hat{A}}_0\times {\hat{B}}_0,{\hat{C}}_1={\hat{A}}_0\times {\hat{B}}_1+{\hat{A}}_1\times {\hat{B}}_0+{\hat{A}}_1\times {\hat{B}}_1$。

$$\begin{array}{rl}& \mathrm{merge}({\hat{A}}_0,{\hat{A}}_0+{\hat{A}}_1)\times \mathrm{merge}({\hat{B}}_0,{\hat{B}}_0+{\hat{B}}_1)\\ & =\mathrm{merge}({\hat{A}}_0\times {\hat{B}}_0,{\hat{A}}_0\times {\hat{B}}_0+{\hat{A}}_0\times {\hat{B}}_1+{\hat{A}}_1\times {\hat{B}}_0+{\hat{A}}_1\times {\hat{B}}_1)\\ & =\mathrm{merge}({\hat{C}}_0,{\hat{C}}_0+{\hat{C}}_1)\\ & =C\end{array}$$

既然正变换是 $\hat{A}=\mathrm{merge}({\hat{A}}_0,{\hat{A}}_0+{\hat{A}}_1)$，那么逆变换应该是 $A=\mathrm{merge}(A_0,A_1-A_0)$。不论是正，逆变换，都可以使用分治在 $O(n{2}^n)$ 的时间复杂度内完成。

void FWT_OR(int* a, int v) {
    for (int k = 2; k <= (1 << n); k <<= 1)
        for (int len = k >> 1, i = 0; i < (1 << n); i += k)
            for (int j = i; j < i + len; j++) (a[j + len] += (a[j] * v + mod) % mod) %= mod;
}

当 v = 1 时是正变换，v = -1 时是逆变换。

事实上，观察上面的代码（或者手玩几组数据）可以发现最终的 ${\hat{A}}_i=\sum _{j|i=i}{A}_j$。代入到式子里面也是正确的。

可以看出来 or 卷积是高维前缀和。

and#

和 or 操作类似的，构造 $\hat{A}=\mathrm{merge}({\hat{A}}_0+{\hat{A}}_1,{\hat{A}}_0)$，容易证明这样变换是正确的。逆变换是 $A=\mathrm{merge}(A_0-A_1,A_1)$。

容易发现 ${\hat{A}}_i=\sum _{j\mathrm{\&}i=i}{A}_j$。

可以看出来 and 卷积实际上是高维后缀和。

void FWT_AND(int* a, int v) {
    for (int k = 2; k <= (1 << n); k <<= 1)
        for (int len = k >> 1, i = 0; i < (1 << n); i += k)
            for (int j = i; j < i + len; j++) (a[j] += (a[j + len] * v + mod) % mod) %= mod;
}

xor#

因为有 $\mathrm{popcount}(x)+\mathrm{popcount}(y)\equiv \mathrm{popcount}(x\mathrm{xor}y)(mod2)$ 所以 $(-1{)}^{\mathrm{popcount}(x)}(-1{)}^{\mathrm{popcount}(y)}=(-1{)}^{\mathrm{popcount}(x\mathrm{xor}y)}$。

而 $x\mathrm{\&}(y\mathrm{xor}z)=(x\mathrm{xor}y)\mathrm{\&}(x\mathrm{xor}z)$ 所以可以构造 ${\hat{A}}_i=\sum _{j=0}^{2^n-1}(-1{)}^{\mathrm{popcount}(i\mathrm{\&}j)}{A}_j$。也可以写成 $\sum _{j=0}^{2^n-1}(-1{)}^{|i\cap j|}{A}_j$。

根据式子，可以得到变换 $\hat{A}=\mathrm{merge}({\hat{A}}_0+{\hat{A}}_1,{\hat{A}}_0-{\hat{A}}_1)$。这是因为原来 ${\hat{A}}_0,{\hat{A}}_1$ 都没有考虑下标二进制下的第 $n$ 位，而 $\hat{A}$ 要考虑第 $n$ 位了。当原先不考虑第 $n$ 位，现在两个下标第 $n$ 位都是 $1$，那么 $\mathrm{popcount}(i^{\prime }\mathrm{\&}{j}^{\prime })=\mathrm{popcount}(i\mathrm{\&}j)+1$，$(-1{)}^{\mathrm{popcount}(i^{\prime }\mathrm{\&}{j}^{\prime })}$ 符号就反过来了。如果其中有一个下标二进制下第 $n$ 位为 $0$，那么 $\mathrm{popcount}(i^{\prime }\mathrm{\&}{j}^{\prime })=\mathrm{popcount}(i\mathrm{\&}j)$，符号不变。

显然其逆变换是 $A=\mathrm{merge}(\frac{A_0+A_1}{2},\frac{A_0-A_1}{2})$。

看起来挺像 FFT 的，事实上 xor 卷积是高维循环卷积。

void FWT_XOR(int* a, int v) {
    for (int k = 2; k <= (1 << n); k <<= 1)
        for (int len = k >> 1, i = 0; i < (1 << n); i += k)
            for (int j = i; j < i + len; j++) {
                int x = a[j], y = a[j + len];
                a[j] = 1ll * (x + y) * v % mod;
                a[j + len] = 1ll * (x - y + mod) * v % mod;
            }
}

Code of P4717

#include<cstdio>
const int M = 2e5 + 10;
const int mod = 998244353;

int n, a[M], b[M], c[M];

void FWT_OR(int* a, int v) {
    for (int k = 2; k <= (1 << n); k <<= 1)
        for (int len = k >> 1, i = 0; i < (1 << n); i += k)
            for (int j = i; j < i + len; j++) (a[j + len] += (a[j] * v + mod) % mod) %= mod;
}

void FWT_AND(int* a, int v) {
    for (int k = 2; k <= (1 << n); k <<= 1)
        for (int len = k >> 1, i = 0; i < (1 << n); i += k)
            for (int j = i; j < i + len; j++) (a[j] += (a[j + len] * v + mod) % mod) %= mod;
}

void FWT_XOR(int* a, int v) {
    for (int k = 2; k <= (1 << n); k <<= 1)
        for (int len = k >> 1, i = 0; i < (1 << n); i += k)
            for (int j = i; j < i + len; j++) {
                int x = a[j], y = a[j + len];
                a[j] = 1ll * (x + y) * v % mod;
                a[j + len] = 1ll * (x - y + mod) * v % mod;
            }
}

int main() {
    scanf("%d", &n);
    for (int i = 0; i < (1 << n); i++) scanf("%d", &a[i]);
    for (int i = 0; i < (1 << n); i++) scanf("%d", &b[i]);

    FWT_OR(a, 1), FWT_OR(b, 1);
    for (int i = 0; i < (1 << n); i++) c[i] = 1ll * a[i] * b[i] % mod;
    FWT_OR(c, -1), FWT_OR(a, -1), FWT_OR(b, -1);
    for (int i = 0; i < (1 << n); i++) printf("%d ", c[i]);
    putchar('\n');

    FWT_AND(a, 1), FWT_AND(b, 1);
    for (int i = 0; i < (1 << n); i++) c[i] = 1ll * a[i] * b[i] % mod;
    FWT_AND(c, -1), FWT_AND(a, -1), FWT_AND(b, -1);
    for (int i = 0; i < (1 << n); i++) printf("%d ", c[i]);
    putchar('\n');

    FWT_XOR(a, 1), FWT_XOR(b, 1);
    for (int i = 0; i < (1 << n); i++) c[i] = 1ll * a[i] * b[i] % mod;
    FWT_XOR(c, (mod + 1) >> 1);
    for (int i = 0; i < (1 << n); i++) printf("%d ", c[i]);
    return 0;
}

CF662C Binary Table#

link of Luogu

有一个 $n$ 行 $m$ 列的表格，每个元素都是 $0/1$ 。

每次操作可以选择一行或一列，把 $0/1$ 翻转，即把 $0$ 换为 $1$ ，把 $1$ 换为 $0$ 。

请问经过若干次操作后，表格中最少有多少个 $1$ 。

$1\le n\le 20,1\le m\le {10}^5$。

$n$ 很小，可以将每一列看做一个状态。

每行只有翻转和不翻转之分。考虑翻转的行的集合为 $S$，原先这一列的状态为 $x$，翻转之后变成 $x\mathrm{xor}S$。

$A_x$ 为状态为 $x$ 的列的数量，$B_x=max(\mathrm{popcount}(x),n-\mathrm{popcount}(x))$。那么翻转行的集合为 $S$ 时答案就是 $\sum _{x=0}^{2^n-1}{A}_x{B}_{x\mathrm{xor}S}$，即 $\sum _{x\mathrm{xor}y=S}{A}_x{B}_y$，FWT 优化即可。

Code of CF662C

#include<cstdio>
#include<algorithm>
const long long M = 2e6 + 10;

long long A[M], B[M], C[M];
int n, m, a[21][100010];

void FWT_XOR(long long* a, int op) {
    for (int k = 2; k <= (1 << n); k <<= 1)
        for (int len = k >> 1, i = 0; i < (1 << n); i += k)
            for (int j = i; j < i + len; j++) {
                long long x = a[j], y = a[j + len];
                if (op == 1)
                    a[j] = x + y, a[j + len] = x - y;
                if (op == -1)
                    a[j] = (x + y) / 2, a[j + len] = (x - y) / 2;
            }
}

int main() {
    scanf("%d%d", &n, &m);
    for (int i = 0; i < n; i++)
        for (int j = 0; j < m; j++) scanf("%1d", &a[i][j]);
    for (int j = 0; j < m; j++) {
        int S = 0;
        for (int i = 0; i < n; i++) S += a[i][j] * (1 << i);
        A[S]++;
    }
    for (int i = 1; i < (1 << n); i++) B[i] = B[i - (i & -i)] + 1;
    for (int i = 0; i < (1 << n); i++) B[i] = std::min(B[i], n - B[i]);

    FWT_XOR(A, 1), FWT_XOR(B, 1);
    for (int i = 0; i < (1 << n); i++) C[i] = A[i] * B[i];
    FWT_XOR(C, -1);
    long long ans = 10000000000007;
    for (int i = 0; i < (1 << n); i++) ans = std::min(ans, C[i]);
    printf("%lld\n", ans);
    return 0;
}

HDU5823 Color II#

给定一张无向图，设点集 $S$ 至少要 $g(S)$ 种颜色染色，能够使得其中有边相连的两点不同颜色。求 $\sum _{S}g(S)\times {233}^{S}mod{2}^{32}$。

$n\le 18$。

设 $f_{x,S}$ 为 $x$ 种颜色给点集 $S$ 染色的方案数。那么有 $f_{x,S}=\sum _{S_1\cap S_2=S}{f}_{x-1,S_1}{f}_{1,S_2}$，即先用 $x-1$ 种颜色染 $S_1$，然后再将 $S_2$ 中的点涂上新的颜色。$g(S)$ 即为 $f_{x,S}$ 不为 $0$ 的最小的 $x$。

先预处理出 $f_{1,S}$，然后 $n$ 次 FWT 求出每个 $f_{x,S}$ 即可。时间复杂度 $O(n^2{2}^n)$。

Code of HDU5823

#pragma GCC optimize("Ofast")
#include <cstdio>
#include <cstring>
const int N = 1 << 18, M = 20;

typedef unsigned int uint;

int T, n;
bool G[M][M], vis[N];
uint f[M][N], pw233[N];

void FWT_OR(uint* a, int op) {
    for (int k = 2; k <= (1 << n); k <<= 1)
        for (int len = k >> 1, i = 0; i < (1 << n); i += k)
            for (int j = i; j < i + len; j++) a[j + len] += a[j] * op;
}

inline void clear();
int main() {
    pw233[0] = 1;
    for (int i = 1; i < N; i++) pw233[i] = pw233[i - 1] * 233;
    scanf("%d", &T);
    while (T--) {
        clear();
        scanf("%d", &n);
        for (int i = 0; i < n; i++)
            for (int j = 0; j < n; j++) scanf("%1d", &G[i][j]);

        for (int S = 1; S < (1 << n); S++) {
            bool flag = true;
            for (int i = 0; i < n && flag; i++)
                for (int j = 0; j < n && flag; j++)
                    if ((S & (1 << i)) && (S & (1 << j)) && G[i][j]) flag = false;
            if (flag) f[1][S] = 1;
        }

        uint ans = 0;
        for (int i = 1; i < n; i++) {
            for (int S = 1; S < (1 << n); S++)
                if (f[i][S] && !vis[S]) vis[S] = true, ans += i * pw233[S];
            FWT_OR(f[i], 1);
            for (int S = 0; S < (1 << n); S++) f[i + 1][S] = f[1][S] * f[i][S];
            FWT_OR(f[i + 1], -1);
        }
        for (int S = 1; S < (1 << n); S++)
            if (f[n][S] && !vis[S]) ans += n * pw233[S];
        printf("%u\n", ans);
    }
    return 0;
}

inline void clear() {
    for (int i = 0; i < n; i++)
        for (int j = 0; j < n; j++) G[i][j] = 0;
    for (int i = 0; i < n; i++)
        for (int j = 0; j < (1 << n); j++) f[i][j] = 0;
    for (int i = 0; i < (1 << n); i++) vis[i] = false;
}

UOJ310 黎明前的巧克力#

长为 $n$ 的序列 $a$，选出两个互不相交的集合，满足这两个集合内元素异或和相同且不同时为空，问方案数。

$n\le {10}^6,a_i\le {10}^6$。

两个集合元素异或和相同，说明将选出的所有元素异或起来和为 $0$。

设 $f_{i,x}$ 为前 $i$ 个数选出若干个分成两个集合，两个集合元素异或和异或起来为 $x$ 的方案数。

显然有转移 $f_{i,x}=f_{i-1,x}+2f_{i-1,x\mathrm{xor}{a}_i}$，但是 $O(n^2)$ 的。

FWT 优化转移，构造向量 $g_i$ 满足 $g_{i,0}=1,g_{i,a_i}=2$，其余元素都为 $0$。那么 ${\widehat{f_i}}_x={\widehat{f_{i-1}}}_x{\widehat{g_i}}_x$，即 ${\widehat{f_n}}_x=\prod _{i=1}^n{\widehat{g_i}}_x$。

但是这样要做 $n$ 轮异或卷积，会 T 飞。

考虑异或卷积的变换：

$$A_x=\sum _{y=0}^{2^n-1}(-1{)}^{|x\cap y|}{A}_y$$

但是 $g_i$ 里面只有 $g_{i,0}$ 和 $g_{i,a_i}$ 不为 $0$，可以简化：

$$\begin{array}{rl}{A}_x& =(-1{)}^{|0\cap x|}{A}_0+(-1{)}^{|a_i\cap x|}{A}_{a_i}\\ & =1+2\times (-1{)}^{|a_i\cap x|}\end{array}$$

故 $g_{i,x}$ 的取值只有可能为 $-1$ 或 $3$。设 $g_{\ast ,x}$ 共有 $p$ 个 $3$，$n-p$ 个 $-1$，那么最终 $\prod _{i=1}^n{g}_{i,x}=3^p(-1{)}^{n-p}$。

考虑如何快速求得 $p$。因为只有两种取值，可以求出 $\sum _{i=1}^n{\widehat{g_i}}_x$，然后鸡兔同笼就能求出 $p$。而 FWT 是线性变换，所以可以全部加起来之后再变换。最终得到 ${\hat{h}}_x=\sum _{i=1}^n{\widehat{g_i}}_x$，那么 $3p-(n-p)={\hat{h}}_x$ 即 $p=\frac{{\hat{h}}_x+n}{4}$。

于是 $\prod _{i=1}^n{\widehat{g_i}}_x$ 就容易求得。和初始状态 $f_{0,0}=1$ 卷积即可。

其实直接对 $\prod _{i=1}^n{\widehat{g_i}}_x$ 做 IFWT 即可，因为只有 $f_{0,0}=1$ 时，$\widehat{f_0}$ 的每个元素都是 $1$。

答案是 $\hat{(\prod _{i=1}^n{\widehat{g_i}}_0)}$。

Code of UOJ310

#include <cstdio>
#include <cassert>
#include <cstring>
#include <algorithm>
#define int long long

const int N = 1 << 20, M = N + 10;
const int mod = 998244353;

int n, a[M], f[M], g[M], h[M], pw3[M];

void FWT_XOR(int* a, int v, bool flag = true) {
    for (int k = 2; k <= N; k <<= 1)
        for (int l = k >> 1, i = 0; i < N; i += k)
            for (int j = i; j < i + l; j++) {
                int x = a[j], y = a[j + l];
                a[j] = 1ll * (x + y) * v;
                a[j + l] = 1ll * (x - y) * v;
                if (flag) ((a[j] %= mod) += mod) %= mod, ((a[j + l] %= mod) += mod) %= mod;
            }
}

signed main() {
    pw3[0] = 1;
    for (int i = 1; i < M; i++) pw3[i] = 3ll * pw3[i - 1] % mod;

    scanf("%lld", &n);
    for (int i = 1; i <= n; i++) scanf("%lld", &a[i]);
    for (int i = 1; i <= n; i++) h[0]++, h[a[i]] += 2;
    FWT_XOR(h, 1, false);
    int mn = 0, mx = 0;
    for (int i = 0; i < N; i++) mn = std::min(mn, h[i]), mx = std::max(mx, h[i]);

    for (int i = 0; i < N; i++) {
        int _3 = (h[i] + n) / 4;
        h[i] = ((n - _3) & 1) ? (mod - pw3[_3]) : pw3[_3];
    }
    FWT_XOR(h, (mod + 1) >> 1, true);
    printf("%d\n", (h[0] + mod - 1) % mod);
    return 0;
}

其他的题还有 BZOJ4589 Hard Nim | sol

---

其他位运算#

FWT 不止可以应用于 or, and, xor 卷积。

FWT 可以看做一个向量 $A$ 乘上变换矩阵 $G$ 得到 $\hat{A}$。显然矩阵 $G$ 可逆。

$$\begin{array}{rl}{\hat{C}}_x& ={\hat{A}}_x{\hat{B}}_x\\ \sum _{i=0}^n{G}_{x,i}{C}_i& =\sum _{j=0}^n{G}_{x,j}{A}_j\sum _{k=0}^n{G}_{x,k}{B}_k\\ \sum _{i=0}^n{G}_{x,i}\sum _{j\oplus k=i}{A}_j{B}_k& =\sum _{j=0}^n\sum _{k=0}^n{G}_{x,j}{G}_{x,k}{A}_j{B}_k\\ \sum _{j=0}^n\sum _{k=0}^n{G}_{x,i\oplus j}{A}_j{B}_k& =\sum _{j=0}^n\sum _{k=0}^n{G}_{x,j}{G}_{x,k}{A}_j{B}_k\\ G_{x,j\oplus k}& =G_{x,j}{G}_{x,k}\end{array}$$

所以如果矩阵 $G$ 满足上述条件，那么 $G$ 就是一个合法的变换矩阵。

例如 or 卷积，$G=\left[\begin{array}{cc}1& 1\\ 0& 1\end{array}\right]$，$G\left[\begin{array}{c}{A}_0\\ A_1\end{array}\right]=\left[\begin{array}{c}{A}_0\\ A_0+A_1\end{array}\right]$。而 $G^{-1}=\left[\begin{array}{cc}1& 0\\ -1& 1\end{array}\right]$ 恰好是逆变换的矩阵。

HDU6966 I love sequences#

给定三个长为 $n$ 的序列 $a,b,c$。

设 $d_{p,k}=\sum _{k=i\oplus j}{a}_i{b}_j$，其中 $1\le i,j\le \frac{n}{p}$。$\oplus$ 代表三进制下的按位 $gcd$。

求 $\sum _{p=1}^n\sum _{k=1}^{+\mathrm{\infty }}{d}_{p,k}{c}_p^k$。

$n\le 2\times {10}^5$。

枚举 $p$，$\sum \frac{n}{p}=O(n\ln n)$，所以对于每个 $p$ 暴力求出 $d_p$。

虽然是三进制，但是还是位运算卷积，考虑 FWT。

列出关于变换矩阵 $G$ 的方程组：

$$\{\begin{array}{l}{G}_{x,0}^2=G_{x,0}\\ G_{x,1}{G}_{x,0}=G_{x,1}\\ G_{x,1}^2=G_{x,1}\\ G_{x,2}{G}_{x,0}=G_{x,2}\\ G_{x,2}{G}_{x,1}=Gx,1\\ G_{x,2}^2=G_{x,2}\end{array}$$

容易解得 $(G_{x,0},G_{x,1},G_{x,2})=(0,0,0),(1,0,1),(1,1,1),(1,0,0)$。显然不能取 $(0,0,0)$ 因为会导致不存在 $G^{-1}$。故最终 $G=\left[\begin{array}{ccc}1& 0& 0\\ 1& 0& 1\\ 1& 1& 1\end{array}\right]$。

其实可以交换任意两行，仍然满足转移矩阵的性质。

故 $\hat{A}=\mathrm{merge}({\hat{A}}_0,{\hat{A}}_0+{\hat{A}}_2,{\hat{A}}_0+{\hat{A}}_1+{\hat{A}}_2)$，逆变换是 $A=\mathrm{merge}(A_0,A_2-A_1,A_1-A_0)$。

对于每个 $p$，FWT 求出数组 $d_p$ 即可。时间复杂度 $O(n{\log }^{2}n)$。

Code of HDU6966

#pragma GCC optimize("Ofast")
#include <cstdio>
#define int long long

const int mod = 1e9 + 7;
const int M = 6e5 + 10;

int n, a[M], b[M], c[M], d[M], a0[M], b0[M];
int _lg3[M], pw3[M];

inline int read() {
    char ch = getchar();
    int x = 0;
    while (ch < '0' || ch > '9') ch = getchar();
    while (ch >= '0' && ch <= '9') x = x * 10 + ch - '0', ch = getchar();
    return x;
}

inline void FWT(int* a, int op, int L) {
    for (int k = 3; k <= L; k *= 3)
        for (int len = k / 3, i = 0; i < L; i += k)
            for (int j = i; j < i + len; j++) {
                int x = a[j], y = a[j + len], z = a[j + 2 * len];
                if (op == 1) {
                    a[j] = x, a[j + len] = (x + z) % mod, a[j + 2 * len] = (x + y + z) % mod;
                } else {
                    a[j] = x, a[j + len] = (z - y + mod) % mod, a[j + 2 * len] = (y - x + mod) % mod;
                }
            }
}

signed main() {
    _lg3[0] = -1, pw3[0] = 1;
    for (int i = 1; i < M; i++) _lg3[i] = _lg3[i / 3] + 1;
    for (int i = 1; i <= 20; i++) pw3[i] = pw3[i - 1] * 3;

    n = read();
    for (int i = 1; i <= n; i++) a[i] = read();
    for (int i = 1; i <= n; i++) b[i] = read();
    for (int i = 1; i <= n; i++) c[i] = read();

    int ans = 0;
    for (int p = 1; p <= n; p++) {
        int L = 1;
        while (L <= n / p) L *= 3;
        //int L = pw3[_lg3[n / p] + 1];
        for (int i = 1; i <= n / p; i++) a0[i] = a[i], b0[i] = b[i];
        FWT(a0, 1, L), FWT(b0, 1, L);
        for (int i = 0; i < L; i++) d[i] = 1ll * a0[i] * b0[i] % mod;
        FWT(d, -1, L);
        //for (int i = 0; i < L; i++) (ans += 1ll * qpow(c[p], i) * d[i] % mod) %= mod;
        int pwc = 1;
        for (int i = 0; i < L; i++)
            (ans += 1ll * pwc * d[i] % mod) %= mod, pwc = 1ll * pwc * c[p] % mod;

        for (int i = 0; i < L; i++) d[i] = a0[i] = b0[i] = 0;
    }
    printf("%lld\n", ans);
    return 0;
}

子集卷积

P6097 【模板】子集卷积#

给定两个长度为 $2^n$ 的序列 $a,b$，求出序列 $c$：

$$c_k=\sum _{\begin{array}{c}i\mathrm{\&}j=0\\ i\mid j=k\end{array}}{a}_i{b}_j$$

$n\le 20$。

只有 $i\mid j=k$ 的限制的话是好求的。考虑如何加上 $i\mathrm{\&}j=0$ 这个限制。

注意到这两个限制相当于 $|S_1|+|S_2|=|S|\wedge S_1\cup S_2=S$，考虑 dp。

设 $f_{|i|,i}=a_i,g_{|j|,j}=b_j,h_{|k|,k}=c_k$，其余项都为 $0$，那么有：

$$h_{x,S}=\sum _{y+z=x}\sum _{S_1\mid S_2=S}{f}_{y,S_1}{g}_{z,S_2}$$

暴力枚举 $y,z$，然后后面部分用 FWT 优化即可。时间复杂度 $O(n^2{2}^n)$。

Code of P6097

#include <cstdio>
const int M = 21, N = 1 << 20;
const int mod = 1e9 + 9;

int n, a[M][N], b[M][N], c[M][N];
int pcnt[N];

inline void FWT_OR(int* a, int op) {
    for (int k = 2; k <= (1 << n); k <<= 1)
        for (int l = k >> 1, i = 0; i < (1 << n); i += k)
            for (int j = i; j < i + l; j++) (((a[j + l] += 1ll * op * a[j] % mod) %= mod) += mod) %= mod;
}

int main() {
    for (int i = 1; i < N; i++)
        pcnt[i] = pcnt[i - (i & -i)] + 1;
    scanf("%d", &n);
    for (int i = 0; i < (1 << n); i++) scanf("%d", &a[pcnt[i]][i]);
    for (int i = 0; i < (1 << n); i++) scanf("%d", &b[pcnt[i]][i]);
    for (int i = 0; i <= n; i++) FWT_OR(a[i], 1), FWT_OR(b[i], 1);
    for (int i = 0; i <= n; i++)
        for (int j = 0; i + j <= n; j++)
            for (int S = 0; S < (1 << n); S++)
                (c[i + j][S] += 1ll * a[i][S] * b[j][S] % mod) %= mod;
    for (int i = 0; i <= n; i++) FWT_OR(c[i], -1);
    for (int i = 0; i < (1 << n); i++) printf("%d ", c[pcnt[i]][i]);
    return 0;
}