跳转至

位运算卷积（FWT&FMT）

FWT

FWT 是作用于序列的一种线性变换，即

\[ fwt(A)_i=\sum_{j=0}^{2^p-1}A_jc_{i,j} \]

因此其满足线性性和可逆性：

\[ \begin{align*} &fwt(k\cdot A)=k\cdot fwt(A)\\ &fwt(A+B)=fwt(A)+fwt(B)\\ &ifwt(fwt(A))=A \end{align*} \]

并且，为了实现位运算卷积，对于某种运算 \(\oplus\)，FWT 对于任意的序列 \(A,B\)满足（星号表示 \(i,j\to i\oplus j\) 的卷积，点表示对位相乘）

\[ \boxed{ fwt(A)\cdot fwt(B)=fwt(A*B) } \]

对于常见的位运算（按位与，按位或，按位异或），我们给出构造：

Bitwise Or

构造 \(g=fwt(f)\) 满足（高维前缀和）

\[ g_S=\sum_{T\subseteq S}f(T) \]

Bitwise And

构造 \(g=fwt(f)\) 满足（高维后缀和）

\[ g_S=\sum_{T\supseteq S}f(T) \]

Bitwise Xor

构造 \(g=fwt(f)\) 满足

\[ g_S=\sum_{T}(-1)^{\operatorname{popcnt}(S\cap T)}f(T) \]

我们现在来验证上面三个构造满足线性性，可逆性，以及对对应的卷积满足 \(fwt(a)\cdot fwt(b)=fwt(a*b)\)。

线性性显然。对于可逆性，注意到按位或和按位与本质上就是高维前/后缀和，因此其逆运算为高维差分；而对于按位异或，构造 \(ifwt(A)=\frac{1}{n}fwt(A)\) 然后发现是对的。

现在证明以上三种构造可以正确实现位运算卷积，即

\[ w\big(S_1\to T\big)\times w\big(S_2\to T\big)=w\big(S_1\oplus S_2\to T\big) \]

其中 \(w\) 表示提取贡献系数。对于按位或：

\[ w(S_1\to T)w(S_2\to T)=[S_1\subseteq T][S_2\subseteq T]=[S_1\cup S_2\subseteq T]=w(S_1\cup S_2\to T) \]

按位与同理。对于按位异或：

\[ w(S_1\to T)w(S_2\to T)=(-1)^{|S_1\cap T|+|S_2\cap T|}=(-1)^{|(S_1\operatorname{xor}S_2)\cap T|}=w(S_1\operatorname{xor}S_2\to T) \]

运算 \(\langle S,T\rangle=(-1)^{|S\cap T|}\) 对异或运算具有神秘的分配律，分配出去会变成乘法。

实现

那么如何把 FWT 实现到 \(O(k2^k)\) 呢？注意到上面三种变换的每一位都是独立的，因此一位一位考虑，每次加入该位产生的交叉贡献。

示例代码（Xor） 
 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
inline void fwt () {
    for(int k = 1; k < n; k <<= 1) {
        int l = k << 1;
        for(int i = 0; i < n; i += l) {
            for(int j = 0; j < k; j++) {
                int x = a[i + j], y = a[i + j + k];
                a[i + j] = (x + y) % MOD;
                a[i + j + k] = (x - y + MOD) % MOD;
            }
        }
    }
}

FMT

对于 Bitwise Or&And，注意到其 fwt 其实就是在进行高维前缀/后缀和。逐维进行贡献，得到 FMT：

for(int i = 0; i < n; i++) {
    for(int s = 0; s < (1 << n); s++) {
        if((s >> i) & 1) f[s] += f[s ^ (1 << i)];
    }
}

优点：代码短。缺点：比 fwt 慢将近一倍。

CF1906K Deck-Building Game

给定 \(n\) 个数 \(a_{1\sim n}\)，求：

\[ \sum_{S\subseteq \{1,2,3,\cdots,n\}} \left[\bigoplus a[S_i]=0\right]2^{|S|} \]

\(n\le 10^5\)

定义 \(x^p\cdot x^q=x^{p\oplus q}\)（异或卷积），那么答案即

\[ [x^0]\prod_{i=1}^{n}(2x^{a_i}+1) \]

现在问题变为如何快速求出 \(10^5\) 个二项式的异或卷积。注意到在我们钦定了值域的前 \(k\) 位之后，剩下的数字只能组合出 \(2^{17-k}\) 级别个异或和。因为前 \(k\) 位要么不变，要么全是 \(0\)，只和奇偶性有关。而后 \(17-k\) 位显然只有 \(2^{17-k}\) 种组合。

考虑在值域上分治，钦定一个二进制位的前缀。由于前面保证了一段区间能组合出的值是和区间长度同阶的，这里直接模仿分治 ntt 就行了。由于需要讨论当前第 \(k\) 位是多少，因此需要记一下区间内集合大小的奇偶性。时间复杂度 \(O(n\log^2 n)\)。

实际上，由于 fwt 数组的每一位都是 \(-1\) 或 \(3\)，因此可以先把所有二项式加起来，算出每一个位置有几个 \(-1\) 和 \(3\)，然后再用快速幂/预处理解决，可以做到 \(O(n\log n)\)。

代码 
 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
#include<iostream>
#define int long long
using namespace std;
const int N = 1.5e5 + 10;
const int MOD = 998244353;

int n;
int cnt[N];
int res[2][N];

int inv[N], fact[N], ifact[N], pw2[N];

inline int C(int a, int b) { return fact[a] * ifact[a - b] % MOD * ifact[b] % MOD; }

void fwt(int a[], int n, int t = 1) {
    for(int k = 1; k < n; k <<= 1) {
        for(int i = 0; i < n; i += 2 * k) {
            for(int j = 0; j < k; j++) {
                int x = a[i + j], y = a[i + j + k];
                a[i + j] = (x + y % MOD) % MOD;
                a[i + j + k] = (x + MOD - y % MOD) % MOD;
            }
        }
    }
    if(t == -1)
    for(int i = 0; i < n; i++) a[i] = a[i] * inv[n] % MOD;
}

void solve(int l, int r) {
    if(l + 1 == r) {
        for(int i = 0; i <= cnt[l]; i++) (res[i & 1][l] += C(cnt[l], i) * pw2[i] % MOD) %= MOD;
        return;
    }
    int mid = (l + r) >> 1, len = mid - l;
    solve(l, mid);
    solve(mid, r);
    fwt(res[0] + l, len);
    fwt(res[0] + mid, len);
    fwt(res[1] + l, len);
    fwt(res[1] + mid, len);
    for(int i = l; i < mid; i++) {
        int x = res[0][i], y = res[1][i], z = res[0][i + len], w = res[1][i + len];
        res[0][i] = x * z % MOD;
        res[1][i] = y * z % MOD;
        res[0][i + len] = y * w % MOD;
        res[1][i + len] = x * w % MOD;
    }
    fwt(res[0] + l, len, -1);
    fwt(res[0] + mid, len, -1);
    fwt(res[1] + l, len, -1);
    fwt(res[1] + mid, len, -1);
}

signed main() {

    ios::sync_with_stdio(0);
    cin.tie(0);

    inv[1] = fact[0] = ifact[0] = pw2[0] = 1;
    for(int i = 2; i < N; i++) inv[i] = MOD - inv[MOD % i] * (MOD / i) % MOD;
    for(int i = 1; i < N; i++) fact[i] = fact[i - 1] * i % MOD;
    for(int i = 1; i < N; i++) ifact[i] = ifact[i - 1] * inv[i] % MOD;
    for(int i = 1; i < N; i++) pw2[i] = pw2[i - 1] * 2 % MOD;

    cin >> n;
    for(int i = 1; i <= n; i++) {
        int x;
        cin >> x;
        ++cnt[x];
    }

    solve(0, 131072);

    cout << (res[0][0] + res[1][0]) % MOD << '\n';

    return 0;
}