­

集合冪級數相關

CHANGE LOG

NOI 大綱里沒有把位運算卷積如 FMT,FWT,子集卷積等知識點單獨列出,但高維前綴和(SOSDP)是應用比較廣泛的重要演算法。

學習上述演算法,首先要理解什麼是集合冪級數。

1. 集合冪級數

1.1 定義

集合冪級數最初由呂凱風在他的 2015 年集訓隊論文《集合冪級數的性質與應用及其快速演算法》中提出。這個名字聽起來非常高級,但實際上不難理解。如果讀者能夠有一點點形式冪級數的基礎知識會更容易理解。

為方便說明,我們引入一些記號:

  • 定義 \(2 ^ X\) 表示 \(X\) 的所有冪集,即 \(X\) 的所有子集組成的集合。
  • 定義 \(U\) 為全集 \(\{1, 2, \cdots, n\}\),其中 \(n = |U|\),即 \(U\) 的大小。用數字來表示元素是因為我們並不關心集合中的元素具體是什麼,數字只是一個易於理解的符號。
  • 若無特殊說明,\(S\) 均為 \(U\) 的子集。

形式化地,域 \(F\) 上的集合冪級數 \(f\)\(2 ^ U \to F\) 的函數,其中對於每個 \(S \subseteq U\),函數 \(f\) 均有對應的取值 \(f_S \in F\)

換句話說,集合冪級數就是從全集的所有子集(全集的冪集)到某個域 \(F\)(OI 中一般研究模某個大質數 \(p\) 意義下的整數域 \(\mathbb Z_p\) 或實數域 \(\mathbb R\))上的 映射。把 \(U\) 的每個子集 \(S\) 帶入函數,均有對應的取值 \(f_S\)

通俗地(不嚴謹),集合冪級數就是把函數的自變數值類型從一個數 \(x\) 換成了一個集合 \(S\),且函數的定義域為全集的所有子集。

  • 對比:形式冪級數由序列導出,而集合冪級數由集合導出。

1.2 表示

有了定義,我們考慮如何表示一個集合冪級數。由於它的定義實在有些抽象,我們需要一個方法清晰地描述集合冪級數 \(f\) 在每個集合 \(S\) 處的取值(總不能說 \(f\)\(\varnothing\) 映射到 \(a\),把 \(\{1\}\) 映射到 \(b\),把 \(\{\cdots\}\) 映射到 \(\cdots\) 吧?這樣太麻煩了,而且不方便形式化表述)。

類比形式冪級數,我們可以給每個集合 \(S\) 對應的取值 \(f_S\) 後添加 佔位符 \(x ^ S\)。即 \(cx ^ S\) 表示一個將 \(S\) 映射到 \(c\ (f_S = c)\) 的集合冪級數,且對於不等於 \(S\)\(U\) 的子集 \(T\)\(f_T\) 的取值為 \(0\)(因為 \(0x ^ T\) 可以忽略)。因此,\(f = \sum\limits_{S \in 2 ^ U} f_S x ^ S\)

例如,\(f = 3x ^ {\varnothing} + 4x ^ {\{1\}} + 9 ^ {\{1, 2\}}\) 對應一個 \(f_{\varnothing} = 3\)\(f_{\{1\}} = 4\)\(f_{\{2\}} = 0\)\(f_{\{1, 2\}} = 9\) 的集合冪級數 \(f\)

  • 注意,\(cx ^ S\) 是一個整體而非 \(c\times x ^ S\),因為 \(x ^ S\) 本身沒有意義,它是佔位符。只有它和它前面的係數作為整體時有意義。\(x ^ S\) 表示集合為 \(S\),而 \(c\) 表示 \(f_S\) 對應的值為 \(c\)
  • 下文有些地方會將 \(S\) 看成一個 \(0 \sim 2 ^ n – 1\) 的二進位數 \(i\)。這表示 \(i = \sum\limits_{j \in S} 2 ^ {j – 1}\),或者說 \(S\)\(i\) 在二進位下值為 \(1\) 的位。讀者需要在集合 \(S\) 和二進位數 \(i\) 之間建立起快速的對應關係。

1.3 加法與乘法

形式冪級數有加法(對應係數相加)和乘法(加法卷積,相關演算法有著名的 FFT 和 NTT),集合冪級數自然也需要定義加法和乘法。

加法比較容易,只要將對應係數相加即可。若集合冪級數 \(h = f + g\),則 \(h_S = f_S + g_S\)。注意,\(f, g\) 對應的全集需相同。

乘法的定義就稍微複雜了一點。考慮集合冪級數 \(h = f \cdot g = \left(\sum\limits_{L \in 2 ^ U} f_L x ^ L \right) \cdot \left(\sum\limits_{R \in 2 ^ U} g_R x ^ R \right)\)。我們希望集合冪級數的乘法對加法有分配律,這樣才能推導出更多性質,擴展集合冪級數的應用。因此,\(h = \sum\limits_{L \in 2 ^ U}\sum\limits_{R \in 2 ^ U} f_L x ^ L \cdot g_R x ^ R\)

考慮如何定義 \(f_L x ^ L \cdot g_R x ^ R\)。設其結果為 \((f_L \times g_R) x ^ {L * R}\),其中 \(*\)\(2 ^ U\) 上的二元運算。

  • 為滿足集合冪級數的乘法對加法的分配律,\(f_L\)\(g_R\) 應當以域 \(F\) 中的乘法相結合(存在域 \(F\) 不滿足乘法對加法的分配律嗎?至少筆者沒見過,而且這不重要 = . =),即 \(\times\) 是域 \(F\) 中的乘法。一般即模 \(p\) 意義下的整數域 \(\mathbb Z_p\)\(2 \times 3 \equiv 1 \pmod 5\))或有理數域 \(\mathbb R\) 上的乘法(\(2.5 \times 2.5 = 6.25\))。
  • 為滿足集合冪級數的乘法交換律,\(*\) 應滿足交換律。
  • 為滿足集合冪級數的乘法結合律,\(*\) 應滿足結合律。
  • 當然,\(*\) 還應該有單位元 \(\varnothing\)。這符合集合之間二元運算關係的基本直覺。

常見的滿足上述性質的集合二元運算並不唯一。因此,集合冪級數的乘法類型也並不唯一,常見的有 並卷積對稱差卷積 以及 子集卷積。如果將集合的所有子集用二進位數表示,則它們分別對應了位運算卷積中的 或卷積異或卷積子集卷積。相關演算法將在接下來介紹。

2. 集合併卷積

2.1 定義

令集合之間的二元運算 \(*\) 為集合併運算即可得到集合併卷積的形式

\[h_S = \sum_{L \in 2 ^ U} \sum_{R \in 2 ^ U} [L \cup R = S] f_L g_R
\]

即下標分別為 \(L\)\(R\) 的兩個數 \(f_L, g_R\) 相乘後會貢獻到下標為 \(L \cup R\) 的值。如果從位運算的角度理解,它就是或卷積,即 \(h_i = \sum\limits_{j | k = i} f_j g_k\)

最暴力地執行卷積,時間複雜度 \(\mathcal{O}(4 ^ n)\),非常劣。

2.2 分治解法

根據 位運算在每一位的獨立性,我們考慮分治求解卷積。

\(f = f ^ – + x ^ {\{n\}} \cdot f ^ +\)\(g = g ^ – + x ^ {\{n\}} \cdot g ^ +\)(注意 \(\cdot\) 表示集合併卷積,即兩個集合貢獻到它們的並)。即我們將第 \(n\) 位單獨提出來,\(f ^ -\) 表示 \(f\) 不包含元素 \(n\) 的所有集合對應的集合冪級數,或者說 \(f_0x ^ 0 \sim f_{2 ^ {n – 1} – 1}x ^ {2 ^ {n – 1} – 1}\) 單獨提取出來,即 \(f ^ – = \sum\limits_{n \notin S} f_Sx ^ S = \sum\limits_{0\leq i < 2 ^ {n – 1}} f_ix ^ i\)。而 \(f ^ +\) 表示 \(f\) 包含元素 \(n\) 的所有集合對應的集合冪級數,同時將這些集合中的元素 \(n\) 去掉,但係數保留,即 \(f ^ + = \sum\limits_{n\in S} f_Sx ^ {S\backslash\{n\}} = \sum\limits_{2 ^ {n – 1} \leq i < 2 ^ n} f_ix ^ {i – 2 ^ {n – 1}}\)

顯然,\(f \cdot g = (f ^ – + x ^ {\{n\}} \cdot f ^ +) \cdot (g ^ – + x ^ {\{n\}} \cdot g ^ +)\)。由於 \(x ^ {\{n\}} \cdot x ^ {\{n\}} = x ^ {\{n\}}\),因此

\[f \cdot g = f ^ – \cdot g ^ – + x ^ {\{n\}} \cdot ((f ^ – + f ^ +) \cdot (g ^ – + g ^ +) – f ^ – \cdot g ^ -)
\]

注意到我們將問題分成了兩次不包含 \(n\) 的兩個集合冪級數的集合併卷積。求解 \(\mathcal{T}(n) = 2\mathcal{T}(\frac n 2) + \mathcal{O}(n)\),得到 \(\mathcal{T}(n) = \mathcal{O}(n \log n)\)。因此,上述分治做法的時間複雜度為 \(\mathcal{O}(n 2 ^ n)\)

vector <int> or_conv(int n, vector <int> &f, vector <int> &g) {
	vector <int> h(1 << n);
	if(n == 0) return h[0] = 1ll * f[0] * g[0] % mod, h;
	vector <int> a(1 << n - 1), b(1 << n - 1), c(1 << n - 1), d(1 << n - 1);
	for(int i = 0; i < 1 << n - 1; i++) {
		a[i] = f[i], b[i] = g[i];
		c[i] = (f[i] + f[i + (1 << n - 1)]) % mod;
		d[i] = (g[i] + g[i + (1 << n - 1)]) % mod; 
	}
	vector <int> l = or_conv(n - 1, a, b), r = or_conv(n - 1, c, d);
	for(int i = 0; i < 1 << n - 1; i++) h[i] = l[i];
	for(int i = 0; i < 1 << n - 1; i++) h[i + (1 << n - 1)] = (r[i] + mod - l[i]) % mod;
	return h;
}

2.3 莫比烏斯變換

上述分治做法已經達到了相當優秀的複雜度,但是由於其涉及到遞歸和數組複製,導致常數並不理想。對分治做法加以改進,我們得到了不需要遞歸的快速莫比烏斯變換(FMT)和快速沃爾什變換(FWT,這將在第 4 章提到)。

觀察分治做法,可以發現在後半部分,我們相當於先將所有 \(f_i (2 ^ {n – 1} \leq i < 2 ^ n)\) 加上 \(f_{i – 2 ^ {n – 1}}\),相乘後再減去 \(f_{i-2 ^ {n – 1}}\) 相乘之後的結果。前面一部分給我們的感覺是 將所有第 \(n\) 位為 \(0\) 的項前的係數加到其對應的第 \(n\) 位為 \(1\) 的項前的係數上,後面一部分則是 進行了這種變換的兩個集合冪級數的對應係數相乘(因為遞歸到最底層是對應係數相乘)後,將所有第 \(n\) 位為 \(1\) 的項前的係數減去其對應的第 \(n\) 位為 \(0\) 的項前的係數

因此,我們考慮 模擬遞歸的過程從大到小 枚舉每一位 \(j(1\leq j \leq n)\),對所有第 \(j\) 位為 \(1\) 的數 \(i\),執行 \(f_i\gets f_i + f_{i – 2 ^ {j – 1}}\)\(g_i\gets g_i + g_{i – 2 ^ {j – 1}}\)。這樣,我們得到了兩個新的集合冪級數 \(\hat f\)\(\hat g\)。令 \(\hat h\)\(\hat f, \hat g\) 對應位置係數相乘後得到的集合冪級數,從小到大 枚舉每一位 \(j(1 \leq j \leq n)\),對所有第 \(j\) 位為 \(1\) 的數 \(i\),執行 \(\hat h_i \gets \hat h_i – \hat h_{i – 2 ^ {j – 1}}\),得到 \(h\)\(h\) 即為 \(f\cdot g\)

更進一步地,注意到枚舉 \(j\) 的順序實際上並沒有影響,因為最終得到的形式均為

\[\begin{aligned}
\hat f_S = \sum_{T\subseteq S} f_T \\
\hat f_i = \sum_{j \mid i = i} f_j
\end{aligned}
\]
  • 註:關於上述操作得到該形式的原因見 2.4 小節。關於該形式的正確性,見本小節最後。

對於 集合冪級數 \(f\),我們稱形如 \(\hat f_i = \sum\limits_{j \mid i = i} f_j\) 的變換為 莫比烏斯變換\(\hat f\) 又記作 \(\mathrm{FMT}(f)\),因為我們通常用快速莫比烏斯變換(Fast Mobius Transform)求解莫比烏斯變換。

考慮 \(\hat h_i \gets \hat h_i – \hat h_{i – 2 ^ {j – 1}}\) 的過程,不難發現對於 \(T\subseteq S\)\(\hat h_T\)\(h_S\) 的貢獻的係數為 \((-1) ^ {|S| – |T|}\),這是因為 \(T\) 每減少一個元素,貢獻就會乘上 \(-1\)。根據 容斥原理 也可以推導得到該結果。對於集合冪級數 \(f\),我們稱形如 \(\hat f_i = \sum\limits_{j\mid i = i} (-1) ^ {\mathrm{popcount}(i) – \mathrm{popcount}(j)}f_j\) 的變換為 莫比烏斯反演(注意與數論的莫比烏斯反演進行區分)。\(\hat f\) 又記作 \(\mathrm{FMI}(f)\)

  • FMT 和 FMI 互為逆變換。

快速莫比烏斯變換和快速莫比烏斯反演的具體演算法在上文已經進行了講解。觀察式子,可以總結出一般形式:枚舉每一位 \(j\),令第 \(j\) 位為 \(1\)\(i\) 執行 \(f_j\gets f_j + c \times f_{j – 2 ^ {j – 1}}\)。容易發現 FMT 即 \(c = 1\),FMI 即 \(c = -1\)。因此,我們有如下程式碼(P4717 【模板】快速莫比烏斯/沃爾什變換 (FMT/FWT) 部分程式碼)。

#include <bits/stdc++.h>
using namespace std;

const int N = 1 << 17;
const int mod = 998244353;

int n, a[N], b[N], A[N], B[N], C[N];
void FMT(int *f, int c) {
	for(int i = 1; i < 1 << n; i <<= 1)
		for(int j = 0; j < 1 << n; j++)
			if(i & j)
				f[j] = (f[j] + 1ll * c * f[j ^ i] + mod) % mod;
}
int main() {
	cin >> 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]);
	memcpy(A, a, N << 2), memcpy(B, b, N << 2);
	FMT(A, 1), FMT(B, 1);
	for(int i = 0; i < 1 << n; i++) C[i] = 1ll * A[i] * B[i] % mod;
	FMT(C, mod - 1);
	for(int i = 0; i < 1 << n; i++) printf("%d ", C[i]);
	return 0;
}

接下來我們證明 FMT 進行集合併卷積的正確性。

據定義,\(h_S = \sum\limits_L \sum\limits_R [L \cup R = S] f_L g_R\),因此 \(\hat h_S = \sum\limits_L \sum\limits_R [L \cup R \subseteq S] f_L g_R\)。由於 \([L \subseteq S][R \subseteq S]\) 等價於 \([L \cup R \subseteq S]\),故 \(\hat h_S = \sum\limits_{L \subseteq S} f_L \sum\limits_{R \subseteq S} g_R = \hat f_L \hat g_R\),證畢。

2.4 高維前綴和

我們將說明 FMT 操作為什麼能得到一個集合冪級數 \(f\) 的莫比烏斯變換形式 \(\hat f\)。它和 高維前綴和(SOSDP,Sum over Subsets DP) 有著密不可分的關係。

回憶二維前綴和的做法,我們利用容斥原理得到遞推式 \(s_{i, j} = s_{i – 1, j} + s_{i, j – 1} – s_{i – 1, j – 1}\)。但這對於更高維數的前綴和並沒有擴展價值,因為對於 \(n\) 維前綴和,利用容斥原理得到的 \(s_i\) 的遞推式共有 \(2 ^ n\) 項(包括 \(s_i\) 本身)。此時時間複雜度驟然上升至 \(2 ^ n \times \prod V_i\),其中 \(V_i\) 是第 \(i\) 個維度的大小。

事實上,求解高維前綴和有更一般的做法。即首先枚舉每一維,對所有數 僅關於這一維 做前綴和。例如,對於三維前綴和 \(s_{n, m, t}\)。首先對於 \(n\) 這一維枚舉所有 \(s_{i, j, k}\),然後令 \(s_{i + 1, j, k} \gets s_{i + 1, j, k} + s_{i, j, k}\),注意 \(i\) 必須升序枚舉(顯然,做一維前綴和的時候沒有人會按下標從大到小枚舉),其次再對 \(m\) 這一維和 \(t\) 這一維做類似的操作。容易證明 \(s_{n, m, t} = \sum\limits_{i \in [1, n]}\sum\limits_{j \in [1, m]}\sum\limits_{k \in [1, t]} s_{i, j, k}\)

對於更高的維度同理。這樣做高維前綴和的時間複雜度是 \(\sum V_i\times \prod V_i\),非常優秀。

關於 FMT,令第 \(j\) 位為 \(1\)\(i\) 執行 \(f_i\gets f_i + f_{i – 2 ^ {j – 1}}\) 可以看做 \(j\) 這一位執行前綴和,因此變換結束後得到的就是 \(f\) 的高維前綴和。故有 \(\hat f_i = \sum\limits_{j\mid i = i} f_j\) 這一形式。

2.5 FMT 的性質

注意到 FMT 本質上是一個形如 \(\hat f_i = \sum\limits_{0 \leq j < 2 ^ n} c(i, j) f_j\)線性變換,其標準矩陣 \(C\)\(C_{i, j} = [j \mid i = i]\)。故 FMT 具有 線性性。對於 FMI 同理。根據 FMT 和 FMI 互逆,它們的標準矩陣也互逆。

因此 \(\mathrm{FMT}(f + g) = \mathrm{FMT}(f) + \mathrm{FMT}(g)\),且對於任意正整數 \(c\),均有 \(\mathrm{FMT}(cf) = c\mathrm{FMT}(f)\)。更一般地,由 疊加原理,對於若干集合冪級數 \(f_i\)\(\mathrm{FMT}(\sum c_if_i) = \sum c_i\mathrm{FMT}(f_i)\)。這裡 \(i\) 表示第 \(i\) 個集合冪級數而非集合冪級數 \(f\)\(x ^ i\) 項係數。

FMT 還有一個重要的性質,那就是 \(\prod f_i = \mathrm{FMI}(\prod \mathrm{FMT}(f_i))\)。等式左邊的求積符號表示 集合併卷積,右邊的求積符號表示 對應項係數相乘。同樣的,\(i\) 表示第 \(i\) 個集合冪級數而非集合冪級數 \(f\)\(x ^ i\) 項係數。

這使得我們在計算若干集合冪級數的集合併卷積時,可以先全部求莫比烏斯變換,對應位置相乘後再莫比烏斯反演回來,而不用每計算兩個集合冪級數的卷積就要莫比烏斯反演一遍。顯然,一個集合冪級數 FMI 之後再 FMT,所得結果不變(互逆性)。

注意,這和若干形式冪級數求加法卷積通常使用的分治 + FFT(而非分治 FFT!)不同。這是因為兩個 \(n\) 次多項式相乘會得到 \(2n\) 次多項式(除非對 \(x ^ {n + 1} – 1\) 取模,這相當於循環卷積),次數快速增長,所以求 DFT 相乘之後要趕緊 IDFT 回去。但兩個 \(n\) 維集合冪級數相乘,結果仍然是 \(n\) 維。

3. 集合對稱差卷積

3.1 定義

定義二元集合運算 對稱差 表示 \(A\oplus B = \{x \mid (x \in A) \oplus (x \in B)\}\),即只屬於其中一個集合的所有元素構成的集合。它相當於邏輯運算中的 異或 運算符。

令集合之間的二元運算 \(*\) 為集合對稱差運算即可得到集合對稱差卷積的形式

\[h_S = \sum_{L \in 2 ^ U} \sum_{R \in 2 ^ U} [L \oplus R = S] f_L g_R
\]

即下標分別為 \(L\)\(R\) 的兩個數 \(f_L, g_R\) 相乘後會貢獻到下標為 \(L \oplus R\) 的值。如果從位運算的角度理解,它就是異或卷積,即 \(h_i = \sum\limits_{j \oplus k = i} f_j g_k\)

最暴力地執行卷積,時間複雜度 \(\mathcal{O}(4 ^ n)\),非常劣。

3.2 分治解法

類似地,利用位運算在每一維獨立的性質,考慮分治。

\(f = f ^ – + x ^ {\{n\}} \cdot f ^ +\)\(g = g ^ – + x ^ {\{n\}} \cdot g ^ +\),類似集合併卷積的推導,我們有

\[\begin{aligned}
fg & = (f ^ – + x ^ {\{n\}}f ^ +)(g ^ – + x ^ {\{n\}}g ^ +) \\
& = (f ^ – g ^ – + f ^ + g ^ +) + x ^ {\{n\}}(f ^ – g ^ + + f ^ + g ^ -)
\end{aligned}
\]

如果直接拆成四個子卷積,時間複雜度仍然是 \(\mathcal{O}(4 ^ n)\)。考慮計算 \(p = (f ^ – + f ^ +)(g ^ – + g ^ +)\) 以及 \(q = (f ^ – – f ^ +)(g ^ – – g ^ +)\)(當然,\((f ^ + – f ^ -)(g ^ + – g ^ -)\) 也可以,可以說明兩者是等價的)。那麼 \(fg = \dfrac {p + q} 2 + x ^ {\{n\}} \dfrac {p – q} 2\)。時間複雜度 \(\mathcal{O}(n 2 ^ n)\)

3.3 沃爾什變換

以下內容均來自於呂凱風的論文。

我們進行一步非常神仙的轉化:\([S = \varnothing] = \dfrac 1 {2 ^ n}\sum\limits_{T \subseteq 2 ^ U} (-1) ^ {|S \cap T|}\)。若 \(S \neq \varnothing\),則選出 \(S\) 中任意一個元素 \(x\),包含 \(x\)\(T\) 和不包含 \(x\)\(T\) 一一對應且貢獻相反。若 \(S = \varnothing\) 顯然成立。

注意到 \(L \oplus R = S\) 等價於 \(L \oplus R \oplus S = \varnothing\),且 \((L \oplus R) \cap S = (L \cap S) \oplus (R \cap S)\),所以

\[\begin{aligned}
h_S & = \sum_{L \in 2 ^ U} \sum_{R \in 2 ^ U} [L \oplus R = S] f_L g_R \\
& = \dfrac 1 {2 ^ n} \sum_{L \in 2 ^ U} \sum_{R \in 2 ^ U} \sum_{T \in 2 ^ U} (-1) ^ {|(S \oplus L \oplus R) \cap T|} f_L g_R \\
& = \dfrac 1 {2 ^ n} \sum_{T \in 2 ^ U} (-1) ^ {|S \cap T|} \left(\sum_{L \in 2 ^ U} (-1) ^ {|L \cap T|} f_L \right) \left(\sum_{R \in 2 ^ U} (-1) ^ {|R \cap T|} g_R \right) \\
\end{aligned}
\]

因此,考慮令 \(\hat f_S = \sum\limits_{T \in 2 ^ U} f_T(-1) ^ {|S \cap T|}\),它被稱為 \(f\)沃爾什變換\(\hat f\) 也記作 \(\mathrm{FWT}(f)\),因為我們用快速沃爾什變換(Fast Walsh-Hadamard Transform)求解沃爾什變換。

同理,根據 \([S = \varnothing] = \dfrac 1 {2 ^ n}\sum\limits_{T \subseteq 2 ^ U} (-1) ^ {|S \cap T|}\),我們定義 \(\hat f\)沃爾什逆變換\(f_S = \dfrac 1 {2 ^ n} \sum\limits_{T \subseteq 2 ^ U} \hat f_T(-1) ^ {|S\cap T|}\)\(f\) 也記作 \(\mathrm{IFWT}(\hat f)\),因為我們用快速沃爾什逆變換(Inverse Fast Walsh-Hadamard Transform)求解沃爾什你變換。

\[\begin{aligned}
f_S & = \dfrac 1 {2 ^ n} \sum\limits_{T \subseteq 2 ^ U} \hat f_T(-1) ^ {|S\cap T|} \\
& = \dfrac 1 {2 ^ n} \sum\limits_{T \subseteq 2 ^ U} \sum\limits_{X \subseteq 2 ^ U} (-1) ^ {|S\cap T|} f_X (-1) ^ {|T \cap X|} \\
& = \sum\limits_{X \subseteq 2 ^ U} f_X \dfrac 1 {2 ^ n} \sum\limits_{T \subseteq 2 ^ U} (-1) ^ {|(S\oplus X)\cap T|} \\
& = \sum\limits_{X \subseteq 2 ^ U} f_X [S\oplus X = \varnothing] \\
& = f_S
\end{aligned}
\]

因此沃爾什變換和沃爾什逆變換互逆。


考慮如何計算沃爾什(逆)變換。直接枚舉子集,時間複雜度 \(\mathcal{O}(3 ^ n)\),比一開始的 \(4 ^ n\) 更優,但仍然不夠快。我們嘗試一位一位地構造。

設當前考慮第 \(j\) 位,則對於兩個除了第 \(j\) 位以外均相同的位置 \(p, q(p + 2 ^ {j – 1} = q)\),由於 \(p\) 這一位是 \(0\),由於 \(0\)\(0 / 1\) 的與均為 \(0\),因此貢獻了 \(0\)\(-1\),故 \(f_p \gets f_p + f_q\)。同理,因為 \(1 \land 0 = 0\)\(1\land 1 = 1\),可以推出 \(f_q \gets f_p – f_q\)。因此對每一位的每兩對這樣的 \(p, q\) 做如上變換,即可得到 \(f\) 的沃爾什變換。這相當於依次把每一位 \(j\)(每一個元素)放到 \(f_i = \sum\limits_{0\leq k < n} (-1) ^ {\mathrm{popcount}(i\ \mathrm{and}\ k)} f_k\) 的考慮範圍內。

  • 方便起見,我們將 \(\mathrm{popcount}(i\ \mathrm{and}\ k)\)\(i\)\(k\) 的按位與的 \(1\) 的個數的奇偶性記作 \(i\odot k\)

沒聽懂?沒關係,我們換一種理解方式。在計算 \(f\) 的沃爾什變換時,設 $f = f ^ – + x ^ {{n}} \cdot f ^ + $,先計算 \(f ^ -\)\(f ^ +\) 的沃爾什變換,再考慮將它們合併。對於 \(f_i\),若 \(i\) 的第 \(n\) 位為 \(0\),由於 \(0 \land 0 = 0\)\(0\land 1 = 0\),所以合併

\[\begin{aligned}
\hat f_i & = \left(\sum\limits_{0\leq j < 2 ^ {n – 1}} f_j \times (-1) ^ {i\odot j} \right) + \left(\sum\limits_{2 ^ {n – 1}\leq j < 2 ^ n} f_j \times (-1) ^ {i \odot j} \right) \\
& = \left(\sum\limits_{0\leq j < 2 ^ {n – 1}} f_j \times (-1) ^ {i\odot j} \right) + \left(\sum\limits_{2 ^ {n – 1}\leq j < 2 ^ n} f_j \times (-1) ^ {i\odot (j – 2 ^ {n – 1})} \right) \\
& = \hat f ^ -_j + \hat f ^ +_j
\end{aligned}
\]

同理,當 \(i\) 的第 \(n\) 位為 \(1\) 時,由於 \(1 \land 0 = 0\)\(1\land 1 = 1\),故

\[\begin{aligned}
\hat f_i & = \left(\sum\limits_{0\leq j < 2 ^ {n – 1}} f_j \times (-1) ^ {i\odot j} \right) + \left(\sum\limits_{2 ^ {n – 1}\leq j < 2 ^ n} f_j \times (-1) ^ {i \odot j} \right) \\
& = \left(\sum\limits_{0\leq j < 2 ^ {n – 1}} f_j \times (-1) ^ {(i – 2 ^ {j – 1})\odot j} \right) {\color{red}-} \left(\sum\limits_{2 ^ {n – 1}\leq j < 2 ^ n} f_j \times (-1) ^ {(i – 2 ^ {n – 1})\odot (j – 2 ^ {n – 1})} \right) \\
& = \hat f ^ -_j – \hat f ^ +_j
\end{aligned}
\]

因此,我們可以考慮化遞歸為遞推,直接從小到大枚舉每一位執行 \(p \gets p + q\)\(q\gets p – q\) 即可。實際上按任意順序枚舉位均可,因為各位獨立且等價。程式碼實現類似 DFT 的寫法。

對於沃爾什逆變換,只需在最後對每一項乘以 \(2 ^ n\) 的逆元即可。需要保證 \(2\) 的逆元存在。

  • 上述快速沃爾什變換的演算法和遞歸演算法本質相同,注意我們令 \(p = f ^ – + f ^ +\)\(q = f ^ – – f ^ +\) 相當於沃爾什變換,而 \(f \cdot g = \dfrac{p + q} 2 + x ^ {\{n\}} \dfrac {p – q} 2\) 相當於逆沃爾什變換,並將 \(2 ^ n\) 的逆元分攤在了所有 \(n\) 層當中。

P4717 【模板】快速莫比烏斯/沃爾什變換 (FMT/FWT) 部分程式碼。

#include <bits/stdc++.h>
using namespace std;

const int N = 1 << 17;
const int mod = 998244353;

int n, a[N], b[N], A[N], B[N], C[N];
void FWT(int *f, int n, int coef) {
	for(int k = 1; k < 1 << n; k <<= 1)
		for(int i = 0; i < 1 << n; i += k << 1)
			for(int j = 0, x, y; j < k; j++)
				x = f[i | j], y = f[i | j | k],
				f[i | j] = (x + y) % mod, f[i | j | k] = (x - y + mod) % mod;
	for(int i = 0; i < 1 << n; i++) f[i] = 1ll * f[i] * coef % mod;
}
int main() {
	cin >> 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]);
	memcpy(A, a, N << 2), memcpy(B, b, N << 2);
	FWT(A, n, 1), FWT(B, n, 1);
	for(int i = 0; i < 1 << n; i++) C[i] = 1ll * A[i] * B[i] % mod;
	int inv = mod + 1 >> 1;
	for(int i = 1; i < n; i++) inv = 1ll * inv * (mod + 1 >> 1) % mod;
	FWT(C, n, inv);
	for(int i = 0; i < 1 << n; i++) printf("%d ", C[i]);
	cout << endl;
	return 0;
}

4. 快速沃爾什變換

本章將從另一個角度詳細講解快速沃爾什變換的原理和推導過程,需要一點點線性代數的知識(相關資訊見 線性代數學習筆記 第一章)。內容部分來自於 command_block 的部落格 位運算卷積(FWT) & 集合冪級數

4.1 核心思想

依據 DFT 的思想,我們 化卷積為乘積,考慮對集合冪級數進行 線性變換 \(T\),設其標準矩陣為 \(C\),即令 \(\hat f_i = \sum\limits_{j = 0} ^ {2 ^ n – 1} C_{i, j}f_j\)。若要求 \(\hat f\)\(\hat g\) 各位相乘後等於 \(\hat f\),不難發現需要有 \(C(i, j) C(i, k) = C(i, j\oplus k)\)。這是因為

\[\begin{aligned}
\hat h_i & = \sum\limits_{j = 0} ^ {2 ^ n – 1} C_{i, j} h_j \\
& = \sum\limits_{j = 0} ^ {2 ^ n – 1} \sum\limits_{k = 0} ^ {2 ^ n – 1}\sum\limits_{l = 0} ^ {2 ^ n – 1} [k \oplus l = j]C_{i, j} f_kg_l \\
& = \sum\limits_{k = 0} ^ {2 ^ n – 1}\sum\limits_{l = 0} ^ {2 ^ n – 1} C_{i, k\oplus l} f_k g_l
\end{aligned}
\]

\[\begin{aligned}
\hat h_i & = \hat f_i \hat g_i \\
& = \left(\sum\limits_{k = 0} ^ {2 ^ n – 1} C_{i, k} f_k\right) \left(\sum\limits_{l = 0} ^ {2 ^ n – 1} C_{i, l} g_l \right)\\
& = \sum\limits_{k = 0} ^ {2 ^ n – 1}\sum\limits_{l = 0} ^ {2 ^ n – 1} C_{i, k} C_{i, l} f_k g_l\\
\end{aligned}
\]

我們利用位運算的一個非常重要的性質,即 各位獨立性,將 \(i, j\) 二進位分解,從而將 \(C(i, j)\) 表示成若干個 \(c(i_p, j_p)\) 相乘,其中 \(i_p, j_p\in [0, 1]\)。即設 \(i\) 的各位分別為 \(i_1, i_2, \cdots, i_p,\cdots,i_n\)\(j\) 同理,則 \(C(i, j) = \prod\limits_{p = 1} ^ n c(i_p, j_p)\)。根據各位獨立性,\(C(i, j) C(i, k) = C(i, j\oplus k)\) 等價於 \(c(i_p, j_p) c(i_p, k_p) = c(i_p, j_p \oplus k_p)\),我們可以從這一結論出發,反過來推導 \(c\) 的具體值。

  • 注意,\(\oplus\) 可以是 任何位運算符,即按位與,按位或和按位異或。
  • 注意,矩陣 \(c\)存在逆元,否則沒有逆變換。

考慮得到矩陣 \(c\) 之後如何快速實現變換。繼續使用按位討論的思路,設 $f = f ^ – + x ^ {{n}} \cdot f ^ + $,當 \(i\) 的第 \(n\) 位為 \(0\) 時,我們有

\[\begin{aligned}
\hat f_i & = \left(c_{0, 0} \times \sum\limits_{0\leq j < 2 ^ {n – 1}}C_{i, j} f_j\right) + \left(c_{0, 1} \times \sum\limits_{2 ^ {n – 1} \leq j < 2 ^ n} C_{i, j – 2 ^ {n – 1}}f_j \right) \\
& = c_{0, 0} \hat f ^ -_j + c_{0, 1} \hat f ^ +_j
\end{aligned}
\]

同理,當 \(i\) 的第 \(n\) 位為 \(1\) 時,我們有 \(\hat f_i = c_{1, 0} \hat f ^ -_j + c_{1, 1} \hat f ^ +_j\)。因此,我們只需遞歸計算 \(\hat f ^ -\)\(\hat f ^ +\) 即可 \(\mathcal{O}(2 ^ n)\) 合併,時間複雜度 \(\mathcal{O}(n 2 ^ n)\)。可以通過從低到高按位枚舉化為遞推。

4.2 或卷積

\(\oplus\) 是按位或時,我們有如下構造:

  • \(0 \lor 0 = 0\),因此 \(c_{0, 0}c_{0, 0} = c_{0, 0}\),以及 \(c_{1, 0}c_{1, 0} = c_{1, 0}\)。這說明 \(c_{0, 0}\)\(c_{1, 0}\) 均等於 \(0\)\(1\),但不同時等於 \(0\),否則 \(c\) 無逆。
  • \(0 \lor 1 = 1\),因此 \(c_{0, 0}c_{0, 1} = c_{0, 1}\),以及 \(c_{1, 0}c_{1, 1} = c_{1, 1}\)。這說明 \(c_{0, 0} = 1\)\(c_{0, 1} = 0\),並且 \(c_{1, 0} = 1\)\(c_{1, 1} = 0\)
  • \(1\lor 0 = 1\),這和上面一條限制等價。
  • \(1 \lor 1 = 1\),因此 \(c_{0, 1}c_{0, 1} = c_{0, 1}\),以及 \(c_{1, 1}c_{1, 1} = c_{1, 1}\)。這說明 \(c_{0, 1}\)\(c_{1, 1}\) 均等於 \(0\)\(1\),但不同時等於 \(0\),否則 \(c\) 無逆。

根據上述限制,我們容易枚舉得到所有合法的矩陣 \(c\)\(\begin{bmatrix} 1 & 1 \\ 1 & 0\end{bmatrix}\) 以及 \(\begin{bmatrix} 1 & 0 \\ 1 & 1\end{bmatrix}\)。其中第二個矩陣相當於子集求和,因為當且僅當 \(i_p = 0\)\(j_p = 1\)\(i_p \lor j_p \neq i_p\),此時 \(c(i_p, j_p) = c_{0, 1} = 0\),即不屬於 \(i\) 的子集的集合對位置 \(i\) 的貢獻為 \(0\),其它情況均為 \(1\)

上述兩個矩陣對應的逆矩陣分別為 \(\begin{bmatrix} 0 & 1 \\ -1 & 1\end{bmatrix}\)\(\begin{bmatrix} 1 & 0 \\ -1 & 1\end{bmatrix}\),後者相當於求子集和的逆變換(互逆的標準矩陣顯然對應著互逆的線性變換)。

因此,通過 FWT 實現的或卷積和通過 FMT 實現的或卷積 本質相同

對於與卷積同理,讀者可自行推導,這裡給出結論:常用的矩陣是 \(c = \begin{bmatrix} 1 & 1 \\ 0 & 1\end{bmatrix}\),這相當於對 \(i\) 的超集求和。

4.3 異或卷積

\(\oplus\) 是按位異或時,我們有如下構造:

  • 對於任意 \(i, j\),均有 \(c_{i, 0}c_{i, j} = c_{i, j}\)。所以 \(c_{0, 0} = c_{1, 0} = 1\)
  • \(0 \oplus 1 = 1\),因此 \(c_{0, 0}c_{0, 1} = c_{0, 1}\),以及 \(c_{1, 0}c_{1, 1} = c_{1, 1}\)。由於上面一條限制,該限制自然滿足。
  • \(1 \oplus 1 = 0\),因此 \(c_{0, 1}c_{0, 1} = c_{0, 0}\),以及 \(c_{1, 1}c_{1, 1} = c_{1, 0}\)。這說明 \(c_{0, 1}\)\(c_{1, 1}\) 等於 \(1\)\(-1\),但不相等,否則 \(c\) 無逆。

綜上,矩陣 \(\begin{bmatrix} 1 & -1 \\ 1 & 1 \end{bmatrix}\)\(\begin{bmatrix} 1 & 1 \\ 1 & -1 \end{bmatrix}\) 合法。其中第二個矩陣所對應的的變換和 3.3 小節所講解的快速沃爾什變換相等,因為 \(\begin{bmatrix} 1 & 1 \\ 1 & -1 \end{bmatrix} \begin{bmatrix} p\\ q \end{bmatrix} = \begin{bmatrix} p + q \\ p – q \end{bmatrix}\)

對第二個矩陣求逆(可以直接手算),得到 \(\begin{bmatrix} 0.5 & 0.5 \\ 0.5 & -0.5 \end{bmatrix}\),恰對應著快速沃爾什逆變換。根據上文的分析,我們可以將 \(\dfrac 1 2\) 提出,放到最後乘以 \(\dfrac 1 {2 ^ n}\),以減少取模次數,減小常數。

5. 子集卷積

5.1 定義

將集合之間的二元運算定義為 不相交集合併,即 \(x ^ L \cdot x ^ R = \begin{cases} 0 & (L \cap R \neq \varnothing) \\ x ^ {L \cup R} & (L \cap R = \varnothing)\end{cases}\),我們得到了 子集卷積 的形式

\[h_S = \sum_{L \in 2 ^ U} \sum_{R \in 2 ^ U} [L \cup R = S][L \cap R = \varnothing] f_L g_R
\]

它和集合併卷積的唯一區別在於兩個集合不能有交。暴力枚舉子集,時間複雜度 \(\mathcal{O}(3 ^ n)\)

5.2 快速演算法

注意到當 \(L \cup R = S\) 時,\(L \cap R = \varnothing\)充要條件\(|L| + |R| = |S|\)。因此,我們設 \(f_{i, S} = [|S| = i]f_S\)\(g\) 同理,令 \(H_i = \sum\limits_{j + k = i} f_j \cdot g_k\)。這裡 \(H_i\) 表示一個集合冪級數而非集合冪級數 \(H\) 的第 \(i\) 項,\(\cdot\) 表示集合併卷積。

根據性質我們有 \(h_S = H_{|S|, S}\),因為只有所有大小之和為 \(|S|\)\(L, R\) 才會貢獻到 \(H_{|S|}\) 上,而取 \(h_{|S|}\) 位置 \(S\) 上的值保證了貢獻到該位置上的集合 \(L, R\) 的並為 \(S\),因此符合限制。

將等式兩邊同時取 FMT,得到 \(\hat H_i = \sum\limits_{j + k = i} \hat f_j \cdot \hat g_k\),這裡 \(\cdot\) 表示對應位置相乘。因此,我們首先 \(\mathcal{O}(n ^ 2 2 ^ n)\) 計算出每個 \(f_j\)\(g_k\) 的莫比烏斯變換(共有 \(n\) 個集合冪級數,每做一次 FMT 需要 \(\mathcal{O}(n 2 ^ n)\) 的時間),然後 \(\mathcal{O}(n ^ 2 2 ^ n)\) 按式子計算出 \(\hat H_i\)(共要做 \(n ^ 2\) 次對應位置相乘,每次需要 \(\mathcal{O}(2 ^ n)\) 的時間),最後再將每個 \(\hat H_i\) 逆變換回來即可。三個部分時間複雜度均為 \(\mathcal{O}(n ^ 2 2 ^ n)\),因此這也是總時間複雜度。

P6097【模板】子集卷積 程式碼如下。

#include <bits/stdc++.h>
using namespace std;

inline int read() {
	int x = 0; char s = getchar();
	while(!isdigit(s)) s = getchar();
	while(isdigit(s)) x = x * 10 + s - '0', s = getchar();
	return x;
}
inline void print(int x) {
	if(x >= 10) print(x / 10);
	putchar(x % 10 + '0');
}

const int N = 20, mod = 1e9 + 9;
long long f[N + 1][1 << N], g[N + 1][1 << N], h[1 << N], ans[1 << N];
int n;
template <class T> void FMT(T *f, int op) {
	for(int j = 1; j < 1 << n; j <<= 1)
		for(int i = 0; i < 1 << n; i++)
			if(i & j)
				f[i] += op ? f[i ^ j] : -f[i ^ j];
	for(int i = 0; i < 1 << n; i++) f[i] = (f[i] % mod + mod) % mod;
}
int main() {
	cin >> n;
	for(int i = 0; i < 1 << n; i++) f[__builtin_popcount(i)][i] = read();
	for(int i = 0; i < 1 << n; i++) g[__builtin_popcount(i)][i] = read();
	for(int i = 0; i <= n; i++) FMT(f[i], 1), FMT(g[i], 1);
	ans[0] = f[0][0] * g[0][0] % mod;
	for(int i = 1; i <= n; i++) {
		memset(h, 0, sizeof(h));
		for(int j = i; ~j; j--) {
			for(int p = 0; p < 1 << n; p++) h[p] += f[j][p] * g[i - j][p];
			if(j % 9 == 0) for(int p = 0; p < 1 << n; p++) h[p] %= mod;
		}
		FMT(h, 0);
		for(int p = 0; p < 1 << n; p++) if(__builtin_popcount(p) == i) ans[p] = h[p];
	}
	for(int p = 0; p < 1 << n; p++) print(ans[p]), putchar(' ');
	return 0;
}

6. ln 和 exp

該部分內容過於高深,筆者能力有限(筆者對形式冪級數的 ln 和 exp 都不甚了解),因此挖坑待填。

7. 例題

I. CF1530F Bingo

非正解。考慮補集轉化,求沒有任何一行 / 列(將對角線看成兩個特殊的列)的事件全部發生的概率。設 \(f_{i, S}\) 表示前 \(i\) 行已經有集合 \(S\) 中的列和對角線存在一個事件沒有發生。轉移則考慮枚舉當前行的 \(2 ^ n\) 個狀態,設沒有發生事件的列和對角線為 \(S\),需滿足 \(S \neq \varnothing\) 否則這一行的事件就全部發生了,不符合限制。求出 \(g_S\) 表示沒有發生事件的列和對角線恰好為 \(S\) 的可能性(儘管 \(S\)\(2 ^ {n + 2}\) 種可能,但 \(g\) 只有 \(2 ^ n – 1\) 個位置上有值,因為它依賴於第 \(i\) 行的 \(2 ^ n – 1\) 個合法狀態),轉移 \(f_i \cdot g \to f_{i + 1}\) 是 OR 卷積的形式。最終答案即 \(1 – f_{n, 2 ^ {n + 2} – 1}\),因為所有列和對角線都需要存在一個事件沒有發生。

時間複雜度 \(\mathcal{O}(2^nn^2)\)\(4\) 倍常數。把 long long 換成 int 加上 7s 時限,勉強能過。程式碼

II. ABC212H Nim Counting

\(f_i\) 表示 \(i\)\(A\) 中出現的次數。根據 NIM 遊戲的結論,先手獲勝當且僅當初始石子異或和個數異或和不為 \(0\)。故答案為 \(\sum\limits_{i = 1} ^ n \mathrm{IFWT} (\mathrm{FWT} ^ i(f))\)\(x ^ 1 \sim x ^ {2 ^ {16} – 1}\) 處的取值。根據 FWT 的性質,變成 \(\mathrm{IFWT} \left(\sum\limits_{i = 1} ^ n\mathrm{FWT} ^ i (f)\right)\) 即可用等比數列求和計算,注意特判 \(f_i = 1\) 的情況等比數列求和不管用。時間複雜度線性對數。程式碼

*III. 某聯考題 /kk

\(k\)\(A_{i,j}\gets A_{i,j}\oplus A_{i,j+1}\oplus A_{i+1,j}\oplus A_{i+1,j+1}\ (0\leq i<n,0\leq j<m)\) 的變換,求 \(A_{0,0}\)\(q\) 組詢問。

\(nm \leq 5\times 10 ^ 6\)\(a_{i, j}, k < 2 ^ {31}\)\(q \leq 5 \times 10 ^ 7\)。時間限制 2.5s,空間限制 1GB。

\((i, j)\) 對於 \((0, 0)\) 的貢獻次數是 \(\dbinom k i \dbinom k j\),因為 \(k\) 次操作里有 \(i\) 次選行 \(+1\)(剩下選 \(0\)),\(j\) 次選列 \(+1\)。根據 Lucas 定理,當且僅當 \(k \& i = i\)\(k \& j = j\)\(k \& (i | j) = i | j\) 時上式為奇數。高維前綴異或和即可。時間複雜度 \(\mathcal{O}(\max(n, m) \log\max(n, m) + q)\)

IV. AT4168 [ARC100C] Or Plus Max

ARC 也會出裸題?

記錄 \(f_S, g_S\) 表示 \(a_i\ (i\in S)\) 的最大值和次大值,用高維前綴和做。\(K = i\) 時的答案 \(ans_i\)\(\max(ans_{i – 1}, f_i + g_i)\)程式碼

*V. P4221 [WC2018]州區劃分

首先判斷一個點集 \(S\) 是否合法:若不連通顯然合法,否則根據歐拉定理可知,若原圖以 \(S\) 為點集的點導出子圖中存在度數為奇數的點則合法,反之不合法。

\(q_S\) 表示 \(\sum_{i \in S} w_i\)\(f_S\) 表示點集為 \(S\) 時的滿意度之和,則轉移方程可寫為

\[f_S=\sum_{T\subset S}f_T\times\left(\dfrac{q_{S\setminus T}}{q_S}\right)^p=\dfrac{1}{q_S^p}\sum_{T\subset S}f_Tq_{S\setminus T}^p
\]

一個顯然的子集卷積形式。

這個 \(\dfrac{1}{q_S^p}\) 不太好處理,我們可以這樣:子集卷積需要額外記錄一維 \(i\) 表示集合大小,而 \(i\) 只會從比它小的 \(j\) 轉移過來,因此我們可以按 \(i\) 從小到大 DP,每 DP 完一層就將該層的 \(f_S\ (|S| = i)\) 做一遍 FMT。這樣就不會出現係數的問題了。

時間複雜度 \(\mathcal{O}(2^nn^2)\)程式碼

VI. CF914G Sum the Fibonacci

三合一題。

注意到 \(s_a | s_b = i\)\(s_a\&s_b = 0\) 是子集卷積,\(s_d \oplus s_e = i\) 是異或卷積,外層是與卷積。於是暴力捲起來即可。時間複雜度 \(\mathcal{O}(n + s\log ^ 2 s)\)程式碼

VII. P5387 [Cnoi2019]人形演舞

打表求出每個數的 SG 值,發現是 \(i – 2 ^ {\lfloor \log_2 i\rfloor} + 1\)。根據 nim 遊戲 SG 值非零先手必勝的結論,直接 FWT 即可,時間複雜度 \(m(\log p + \log m)\)

VIII. CC21JAN ORAND

將存在性問題轉化為計數問題,直接 FWT 即可。一直 FWT 到答案不再更新為止,複雜度為 \(\mathcal{O}(TV \log ^ 2 V)\),因為最多進行 \(\log V\) 次位運算卷積(考慮每次卷積必然更新至少一位)。實際遠遠達不到上界,跑得很快。程式碼

8. 參考文章