快速沃尔什变换（FWT）

参考资料

• 快速沃尔什变换 - OI Wiki
• 沃尔什转换 - 维基百科

思路

详见 题解：P4717 【模板】快速莫比乌斯/沃尔什变换 (FMT/FWT)

实现

分治递归

按位或（or）

void fmt_or(ll *a,int n,int type)
{
	if(n==1)return;
	int mid=n>>1;
	ll *a0=a,*a1=a+mid;
	fmt_or(a0,mid,type);
	fmt_or(a1,mid,type);
	for(int i=0;i<mid;i++)a1[i]=(a1[i]+a0[i]*type+mod)%mod;
}

按位与（and）

void fmt_and(ll *a,int n,int type)
{
	if(n==1)return;
	int mid=n>>1;
	ll *a0=a,*a1=a+mid;
	fmt_and(a0,mid,type);
	fmt_and(a1,mid,type);
	for(int i=0;i<mid;i++)a0[i]=(a0[i]+a1[i]*type+mod)%mod;
}

按位异或（xor）

void fwt_xor(ll *a,int n,int type)
{
	if(n==1)return;
	int mid=n>>1;
	ll *a0=a,*a1=a+mid;
	fwt_xor(a0,mid,type);
	fwt_xor(a1,mid,type);
	for(int i=0;i<mid;i++)
	{
		ll x=a0[i],y=a1[i],t=type==1?1:mod+1>>1;
		a0[i]=(x+y)*t%mod;
		a1[i]=(x-y+mod)*t%mod;
	}
}

按位同或（xnor）

void fwt_xnor(ll *a,int n,int type)
{
	if(n==1)return;
	int mid=n>>1;
	ll *a0=a,*a1=a+mid;
	fwt_xnor(a0,mid,type);
	fwt_xnor(a1,mid,type);
	for(int i=0;i<mid;i++)
	{
		ll x=a0[i],y=a1[i],t=type==1?1:mod+1>>1;
		a0[i]=(y-x+mod)*t%mod;
		a1[i]=(y+x)*t%mod;
	}
}

倍增迭代

按位或（or）

void fmt_or(ll *a,int n,int type)
{
	for(int k=1;k<n;k<<=1)
	{
		for(int i=0;i<n;i+=k<<1)
		{
			for(int j=0;j<k;j++)
			{
				a[i+j+k]=(a[i+j+k]+a[i+j]*type+mod)%mod;
			}
		}
	}
}

按位与（and）

void fmt_and(ll *a,int n,int type)
{
	for(int k=1;k<n;k<<=1)
	{
		for(int i=0;i<n;i+=k<<1)
		{
			for(int j=0;j<k;j++)
			{
				a[i+j]=(a[i+j]+a[i+j+k]*type+mod)%mod;
			}
		}
	}
}

按位异或（xor）

void fwt_xor(ll *a,int n,int type)
{
	for(int k=1;k<n;k<<=1)
	{
		for(int i=0;i<n;i+=k<<1)
		{
			for(int j=0;j<k;j++)
			{
				ll x=a[i+j],y=a[i+j+k],t=type==1?1:mod+1>>1;
				a[i+j]=(x+y)*t%mod;
				a[i+j+k]=(x-y+mod)*t%mod;
			}
		}
	}
}

按位同或（xnor）

void fwt_xnor(ll *a,int n,int type)
{
	for(int k=1;k<n;k<<=1)
	{
		for(int i=0;i<n;i+=k<<1)
		{
			for(int j=0;j<k;j++)
			{
				ll x=a[i+j],y=a[i+j+k],t=type==1?1:mod+1>>1;
				a[i+j]=(y-x+mod)*t%mod;
				a[i+j+k]=(y+x)*t%mod;
			}
		}
	}
}

高维前缀和

按位或（or）

void fmt_or(ll *a,int n,int type)
{
	for(int i=1;i<n;i<<=1)
	{
		for(int j=0;j<n;j++)
		{
			if(i&j)a[j]=(a[j]+a[i^j]*type+mod)%mod;
		}
	}
}

按位与（and）

void fmt_and(ll *a,int n,int type)
{
	for(int i=1;i<n;i<<=1)
	{
		for(int j=0;j<n;j++)
		{
			if(!(i&j))a[j]=(a[j]+a[i^j]*type+mod)%mod;
		}
	}
}

例题

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

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

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

$$C_k=\sum_{i\oplus j=k}A_i\times B_j$$

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

题解

原题链接

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

参考资料

• 快速沃尔什变换 - OI Wiki
• 前缀和 & 差分 - OI Wiki
• 容斥原理 - OI Wiki
• Fast Walsh–Hadamard transform - Wikipedia
• Prefix sum - Wikipedia
• 容斥原理与子集卷积（四） - 好地方bug - 知乎
• FMT（快速莫比乌斯变换） - 秋钧 - 知乎
• FWT(快速沃尔什变换)零基础详解qaq（ACM/OI） - 月下桃子树 - 知乎

题意简述

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

$$C_k=\sum_{i\oplus j=k}A_i\times B_j$$

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

解题思路

约定

1. 为统一记号，设数据规模为：

$$m=2^n$$ $$O(m\log m)=O(n2^n)$$

2. 为避免歧义，位运算不再使用 $\mid,\&,\oplus$，统一记为：

• 按位或：$\operatorname{or}$（或 $\cup$）
• 按位与：$\operatorname{and}$（或 $\cap$）
• 按位异或：$\operatorname{xor}$
• 按位同或：$\operatorname{xnor}$

引入

快速傅立叶变换（FFT，Fast Fourier Transform）可以 $O(m\log m)$ 计算 加法卷积（多项式乘法）：

$$C_k=\sum_{i+j=k}A_i\times B_j$$

若将求和符号中的 加法（$+$）替换为 位运算（$\operatorname{or},\operatorname{and},\operatorname{xor},\operatorname{xnor}$），变成 位运算卷积 该如何处理呢？

直接枚举的时间复杂度为 $O(m^2)$，效率较低，因此需要使用 FMT 和 FWT。

简介

快速莫比乌斯变换（FMT，Fast Möbius Transform）可以 $O(m\log m)$ 计算 或卷积 和 与卷积：

$$C_k=\sum_{i\operatorname{or} j=k}A_i\times B_j$$ $$C_k=\sum_{i\operatorname{and} j=k}A_i\times B_j$$

快速沃尔什变换（FWT，Fast Walsh Transform）可以 $O(m\log m)$ 计算 异或卷积 和 同或卷积：

$$C_k=\sum_{i\operatorname{xor} j=k}A_i\times B_j$$ $$C_k=\sum_{i\operatorname{xnor} j=k}A_i\times B_j$$

有些地方会将 FMT 和 FWT 统称为 FWT。

思想

FMT 和 FWT 的思想与 FFT 类似。

首先 $O(m\log m)$ 进行 正变换：

$$A\xrightarrow{\text{FMT/FWT}}A'$$ $$B\xrightarrow{\text{FMT/FWT}}B'$$

然后 $O(m)$ 逐位相乘（点积）：

$$C'=A'\cdot B'$$

最后 $O(m\log m)$ 通过 逆变换 得到结果：

$$C'\xrightarrow{\text{IFMT/IFWT}}C$$

这样整体时间复杂度为 $O(m\log m)$。

基础知识

位运算 & 集合

每个二进制数可以视为一个 集合，二进制位为 $1$ 的位置表示对应元素属于该集合。

$$(10010110)_2\iff \set{a_1,a_2,a_4,a_7}$$

两个数 按位或 对应集合的 并集，而 按位与 对应 交集。

$$(0101)_2\operatorname{or}(1100)_2=(1101)_2\iff \set{a_0,a_2}\cup\set{a_2,a_3}=\set{a_0,a_2,a_3}$$ $$(0101)_2\operatorname{and}(1100)_2=(0100)_2\iff \set{a_0,a_2}\cap\set{a_2,a_3}=\set{a_2}$$

超集是补集的子集的补集，超集与子集在补集意义上互为对偶。

$$T\supseteq S\iff \overline{T}\subseteq \overline{S}$$

简单来说，就是将 $0/1$ 取反，因此子集和超集是 对称 的。

变换 & 反演

变换（Transform）是指根据特定规则，将一个对象（如数、向量、函数、图形等）映射为另一个对象的过程。

$$f\to g$$

变换的逆过程称为 逆变换，亦称 反演（Inversion）。

$$g\to f$$

在数学和算法中，我们会遇到许多变换，如傅里叶变换、数论变换、莫比乌斯变换和沃尔什变换等。

这里解释一下常见变换的命名方式：

• 不带“快速”前缀的名称，指的是变换本身，并不涉及具体的实现方法。
• 而带“快速”前缀的名称，表示使用一种高效实现该变换的算法。

Zeta 变换

考虑一个集合函数 $f:2^V\to \mathbb{Z}$，其中 $2^V$ 表示 $V$ 的幂集，即 $V$ 的所有子集构成的集合。

定义 子集和（Sum Over Subsets）和 超集和（Sum Over Supersets）：

$$g_0(S)=\sum_{T\subseteq S}f_0(T)$$ $$g_1(S)=\sum_{T\supseteq S}f_1(T)$$

计算 子集和 和 超集和 的过程（$f\to g$）称为 Zeta 变换（Zeta Transform）。

莫比乌斯反演

Zeta 的逆变换（$g\to f$）称为 莫比乌斯反演（Möbius Inversion）。

$$f_0(S)=\sum_{T\subseteq S}(-1)^{|S|-|T|}g_0(T)$$ $$f_1(S)=\sum_{T\supseteq S}(-1)^{|T|-|S|}g_1(T)$$

可以通过 容斥原理 证明：

$$\begin{aligned} \sum_{T\subseteq S}(-1)^{|S|-|T|}g_0(T) &= \sum_{T\subseteq S}(-1)^{|S|-|T|}\sum_{V\subseteq T}f_0(V) \\ &= \sum_{V\subseteq S}f_0(V)\sum_{V\subseteq T\subseteq S}(-1)^{|S|-|T|} \\ &= \sum_{V\subseteq S}f_0(V)\sum_{k=0}^{|S|-|V|}(-1)^k\binom{|S|-|V|}{k} \\ &= \sum_{V\subseteq S}f_0(V)[V=S] \\ &= f_0(S) \end{aligned}$$

以上为子集的情形，超集的证明同理。

莫比乌斯反演也广泛出现于其他数学和算法领域，如 数论中的莫比乌斯反演。

它们的共同点是在某种偏序上逆转累加关系以恢复原函数。

快速莫比乌斯变换（FMT）

原理

考虑如何构造变换规则。

若 $i\operatorname{or} j=k$（并集为 $k$），则 $i$ 和 $j$ 必为 $k$ 的子集：

$$i\operatorname{or} j=k\Rightarrow i\subseteq k\land j\subseteq k$$

若 $i\operatorname{and} j=k$（交集为 $k$），则 $i$ 和 $j$ 必为 $k$ 的超集：

$$i\operatorname{and} j=k\Rightarrow i\supseteq k\land j\supseteq k$$

我们可以构造：

$$\operatorname{FMT_{or}}[A]_k=\sum_{i\subseteq k}A_{i}$$ $$\operatorname{FMT_{and}}[A]_k=\sum_{i\supseteq k}A_{i}$$

则有：

$$\begin{aligned} \operatorname{FMT_{or}}[A]_k\cdot \operatorname{FMT_{or}}[B]_k &= \left(\sum_{i\subseteq k} A_i\right)\times \left(\sum_{j\subseteq k} B_j\right) \\ &= \sum_{i\subseteq k}\sum_{j\subseteq k}A_i\times B_j \\ &= \sum_{i\operatorname{or} j\subseteq k}A_i\times B_j \\ &= \operatorname{FMT_{or}}[C]_k \end{aligned}$$ $$\begin{aligned} \operatorname{FMT_{and}}[A]_k\cdot \operatorname{FMT_{and}}[B]_k &= \left(\sum_{i\supseteq k} A_i\right)\times \left(\sum_{j\supseteq k} B_j\right) \\ &= \sum_{i\supseteq k}\sum_{j\supseteq k}A_i\times B_j \\ &= \sum_{i\operatorname{and} j\supseteq k}A_i\times B_j \\ &= \operatorname{FMT_{and}}[C]_k \end{aligned}$$

以上两式分别对应 子集和 和 超集和，即 Zeta 变换 和 莫比乌斯反演。

分治递归

正变换

考虑如何计算 Zeta 变换。

直接枚举的时间复杂度为 $O(m^2)$，效率较低。借鉴 FFT，我们可以考虑分治。

设序列 $A$ 的前一半为 $A_0$，后一半为 $A_1$。不难发现，$A_0$ 最高位都为 $0$，$A_1$ 的最高位都为 $1$。

那么，$A_0$ 的子集就是它本身的子集，而 $A_1$ 满足条件的子集是它本身和 $A_0$ 的子集：

$$\operatorname{FMT_{or}}[A]=\operatorname{merge}(\operatorname{FMT_{or}}[A_0],\operatorname{FMT_{or}}[A_1]+\operatorname{FMT_{or}}[A_0])$$

其中，$\operatorname{merge}$ 表示将两个序列首尾相连（类似字符串拼接），加法（$+$）表示逐项相加。

这样我们就能在 $O(m\log m)$ 的时间复杂度得到 $\operatorname{FMT_{or}}[A]$ 和 $\operatorname{FMT_{and}}[A]$。

逆变换

既然 $A_0$ 的子集是它本身，$A_1$ 的子集是 $\operatorname{FMT}[A_1]+\operatorname{FMT}[A_0]$，即可得到反演的递推式：

$$\operatorname{IFMT_{or}}[A]=\operatorname{merge}(\operatorname{IFMT_{or}}[A_0],\operatorname{IFMT_{or}}[A_1]-\operatorname{IFMT_{or}}[A_0])$$

按位与

同理可得：

$$\operatorname{FMT_{and}}[A]=\operatorname{merge}(\operatorname{FMT_{and}}[A_0]+\operatorname{FMT_{and}}[A_1],\operatorname{FMT_{and}}[A_1])$$ $$\operatorname{IFMT_{and}}[A]=\operatorname{merge}(\operatorname{IFMT_{and}}[A_0]-\operatorname{IFMT_{and}}[A_1],\operatorname{IFMT_{and}}[A_1])$$

代码实现

正变换和逆变换可以合并为一个函数：正变换取 $\text{type}=1$，逆变换取 $\text{type}=-1$。

$$\operatorname{T_{or}}[A]=\operatorname{merge}(\operatorname{T_{or}}[A_0],\operatorname{T_{or}}[A_1]+\operatorname{T_{or}}[A_0]\times\text{type})$$ $$\operatorname{T_{and}}[A]=\operatorname{merge}(\operatorname{T_{and}}[A_0]+\operatorname{T_{and}}[A_1]\times\text{type},\operatorname{T_{and}}[A_1])$$

综上，可用分治写出 $\operatorname{or}$ 和 $\operatorname{and}$ 的递归版本：

void fmt_or(ll *a,int n,int type)
{
	if(n==1)return;
	int mid=n>>1;
	ll *a0=a,*a1=a+mid;
	fmt_or(a0,mid,type);
	fmt_or(a1,mid,type);
	for(int i=0;i<mid;i++)a1[i]=(a1[i]+a0[i]*type+mod)%mod;
}

void fmt_and(ll *a,int n,int type)
{
	if(n==1)return;
	int mid=n>>1;
	ll *a0=a,*a1=a+mid;
	fmt_and(a0,mid,type);
	fmt_and(a1,mid,type);
	for(int i=0;i<mid;i++)a0[i]=(a0[i]+a1[i]*type+mod)%mod;
}

倍增迭代

类似 FFT，可用倍增写出 $\operatorname{or}$ 和 $\operatorname{and}$ 的更高效（常数更小）的迭代版本：

void fmt_or(ll *a,int n,int type)
{
	for(int k=1;k<n;k<<=1)
	{
		for(int i=0;i<n;i+=k<<1)
		{
			for(int j=0;j<k;j++)
			{
				a[i+j+k]=(a[i+j+k]+a[i+j]*type+mod)%mod;
			}
		}
	}
}

void fmt_and(ll *a,int n,int type)
{
	for(int k=1;k<n;k<<=1)
	{
		for(int i=0;i<n;i+=k<<1)
		{
			for(int j=0;j<k;j++)
			{
				a[i+j]=(a[i+j]+a[i+j+k]*type+mod)%mod;
			}
		}
	}
}

高维前缀和

前缀和

一维序列 $A$ 的 前缀和（Prefix Sum）序列 $S$，表示索引 不超过 它的元素之和。

$$S_i=\sum_{j=1}^i A_j$$

高维前缀和

假设有 $n$ 个维度，每个维度的索引范围为 $[0,p)$。

扩展到 高维前缀和，就是各维度索引都 不超过 它的元素之和。

$$S_{i_1,\dots,i_n}=\sum_{j_1\le i_1,\dots,j_n\le i_n}A_{j_1,\dots,j_n}$$

计算高维前缀和的常用方法有两种：容斥原理 与 逐维前缀和。

容斥原理

根据 容斥原理，可得：

$$S_{i_1,\dots,i_n}=\sum_{T\subseteq \set{1,2,\dots,n}}(-1)^{|T|}A_T$$

但该方法的时间复杂度为 $O(2^np^n)$，效率较低。

逐维前缀和

注意到，$n$ 维前缀和可视为 $n$ 次逐维求和：

$$S_{i_1,\dots,i_n}=\sum_{j_1\le i_1,\dots,j_n\le i_n}A_{j_1,\dots,j_n}=\sum_{j_1\le i_1}\dots\sum_{j_n\le i_n} A_{j_1,\dots,j_n}$$

因此，每次只考虑一个维度，固定其余维度，求若干个一维前缀和。

对所有 $n$ 个维度分别求和后，即得到 $n$ 维前缀和。

该方法的时间复杂度为 $O(np^n)$。

正变换

当 $p=2$，每维只有 $0/1$ 两种索引时，高维前缀和 就是 子集和。

同理，也可以通过 高维后缀和（K-Dimensional Suffix Sum）得到 超集和。

所以，我们可以通过 高维前缀和 求出 Zeta 变换。

时间复杂度为 $O(n2^n)=O(m\log m)$。

逆变换

类似地，对所有 $n$ 个维度分别做差分后，可通过 高维差分 得到莫比乌斯反演。

代码实现

我们可以写出 $\operatorname{or}$ 和 $\operatorname{and}$ 的高维前缀和版本：

void fmt_or(ll *a,int n,int type)
{
	for(int i=1;i<n;i<<=1)
	{
		for(int j=0;j<n;j++)
		{
			if(i&j)a[j]=(a[j]+a[i^j]*type+mod)%mod;
		}
	}
}

void fmt_and(ll *a,int n,int type)
{
	for(int i=1;i<n;i<<=1)
	{
		for(int j=0;j<n;j++)
		{
			if(!(i&j))a[j]=(a[j]+a[i^j]*type+mod)%mod;
		}
	}
}

子集卷积

详见 洛谷 P6097 【模板】子集卷积

FMT 还可以计算 子集卷积：

$$C_k=\sum_{\substack{i\operatorname{or} j=k\\i\operatorname{and} j=0}}A_i\times B_j$$

第一个限制 $i\operatorname{or} j=k$，可以通过 或卷积 轻松处理。

$$(A*_{\operatorname{or}}B)_k=\sum_{i\operatorname{or} j=k}A_i\times B_j$$

第二个限制 $i\operatorname{and} j=0$，即没有公共元素，等价于 $|i|+|j|=|i\operatorname{or} j|$。

$$i\operatorname{and} j=0\iff |i|+|j|=|i\operatorname{or} j|$$

我们可以多开一个维度，记录 集合大小。具体令：

$$F_{i,j}=\begin{cases} A_j & |j|=i \\ 0 & |j|\ne i \end{cases}, G_{i,j}=\begin{cases} B_j & |j|=i \\ 0 & |j|\ne i \end{cases}$$

则有：

$$H_{i,j}=\sum_{k=0}^i(F_{k}*_{\operatorname{or}}G_{i-k})_j$$

答案为：

$$C_i=H_{|i|,i}$$

时间复杂度为 $O(n^22^n)$。

快速沃尔什变换（FWT）

原理

考虑构造异或运算的变换。

若令 $i\circ j$ 表示 $i\cap j$ 中 $1$ 数量的奇偶性，即：

$$i\circ j=\operatorname{popcount}(i\cap j)\bmod 2$$

其中 $\operatorname{popcount}(x)$ 表示 $x$ 二进制中 $1$ 的个数。

则有：

$$(i\circ k)\operatorname{xor} (j\circ k)=(i\operatorname{xor} j)\circ k$$

我们可以构造：

$$\operatorname{FWT_{xor}}[A]_k=\sum_{i\circ k=0}A_i-\sum_{i\circ k=1}A_i$$

则有：

$$\begin{aligned} \operatorname{FWT_{xor}}[A]_k\cdot \operatorname{FWT_{xor}}[B]_k &= \left(\sum_{i\circ k=0}A_i-\sum_{i\circ k=1}A_i\right)\times\left(\sum_{j\circ k=0}B_j-\sum_{j\circ k=1}B_j\right) \\ &= \left(\sum_{i\circ k=0}A_i\sum_{j\circ k=0}B_j+\sum_{i\circ k=1}A_i\sum_{j\circ k=1}B_j\right) \\ &- \left(\sum_{i\circ k=0}A_i\sum_{j\circ k=1}B_j+\sum_{i\circ k=1}A_i\sum_{j\circ k=0}B_j\right) \\ &= \sum_{(i\operatorname{xor} j)\circ k=0}A_iB_j-\sum_{(i\operatorname{xor} j)\circ k=1}A_iB_j \\ &= \operatorname{FWT_{xor}}[C]_k \end{aligned}$$

而以上式子就是 沃尔什变换。

分治递归

正变换

我们依然考虑分治。

设序列 $A$ 的前一半为 $A_0$，后一半为 $A_1$。不难发现，$A_0$ 最高位都为 $0$，$A_1$ 的最高位都为 $1$。

当索引 $i$ 的最高位为 $0$ 时，无论与最高位为 $0$ 还是 $1$ 的 $j$ 进行 $\circ$ 运算，因 $0\cap 0=0,\,0\cap 1=0$，奇偶性不变，因此这一部分为两者之和：

$$\mathrm{FWT_{xor}}[A]_i=\mathrm{FWT_{xor}}[A_0]_i+\mathrm{FWT_{xor}}[A_1]_i$$

当索引 $i$ 的最高位为 $1$ 时，与最高位为 $0$ 的 $j$ 运算结果为 $1\cap 0=0$（奇偶性不变），与最高位为 $1$ 的 $j$ 运算结果为 $1\cap 1=1$（奇偶性翻转），因此这一部分为两者之差：

$$\mathrm{FWT_{xor}}[A]_i=\mathrm{FWT_{xor}}[A_0]_i-\mathrm{FWT_{xor}}[A_1]_i$$

综上可得：

$$\operatorname{FWT_{xor}}[A]=\operatorname{merge}(\operatorname{FWT_{xor}}[A_0]+\operatorname{FWT_{xor}}[A_1],\operatorname{FWT_{xor}}[A_0]-\operatorname{FWT_{xor}}[A_1])$$

逆变换

逆变换易得：

$$\operatorname{IFWT_{xor}}[A]=\operatorname{merge}\left(\frac{\operatorname{IFWT_{xor}}[A_0]+\operatorname{IFWT_{xor}}[A_1]}{2},\frac{\operatorname{IFWT_{xor}}[A_0]-\operatorname{IFWT_{xor}}[A_1]}{2}\right)$$

按位同或

同理，将 $\circ$ 定义中的 $\cap$ 替换为 $\cup$，即可得到同或：

$$\operatorname{FWT_{xnor}}[A]=\operatorname{merge}(\operatorname{FWT_{xnor}}[A_1]-\operatorname{FWT_{xnor}}[A_0],\operatorname{FWT_{xnor}}[A_1]+\operatorname{FWT_{xnor}}[A_0])$$ $$\operatorname{IFWT_{xnor}}[A]=\operatorname{merge}\left(\frac{\operatorname{IFWT_{xnor}}[A_1]-\operatorname{IFWT_{xnor}}[A_0]}{2},\frac{\operatorname{IFWT_{xnor}}[A_1]+\operatorname{IFWT_{xnor}}[A_0]}{2}\right)$$

代码实现

正变换和逆变换可以合并为一个函数：正变换取 $\text{type}=1,t=1$，逆变换取 $\text{type}=-1,t=\frac{1}{2}$。

$$\operatorname{T_{xor}}[A]=\operatorname{merge}((\operatorname{T_{xor}}[A_0]+\operatorname{T_{xor}}[A_1])\times t,(\operatorname{T_{xor}}[A_0]-\operatorname{T_{xor}}[A_1])\times t)$$ $$\operatorname{T_{xnor}}[A]=\operatorname{merge}((\operatorname{T_{xnor}}[A_1]-\operatorname{T_{xnor}}[A_0])\times t,(\operatorname{T_{xnor}}[A_1]+\operatorname{T_{xnor}}[A_0])\times t)$$

注意：模意义下的 $\frac{1}{2}$ 指 $2$ 的乘法逆元，即 $\frac{mod+1}{2}$。

综上，可用分治写出 $\operatorname{xor}$ 和 $\operatorname{xnor}$ 的递归版本：

void fwt_xor(ll *a,int n,int type)
{
	if(n==1)return;
	int mid=n>>1;
	ll *a0=a,*a1=a+mid;
	fwt_xor(a0,mid,type);
	fwt_xor(a1,mid,type);
	for(int i=0;i<mid;i++)
	{
		ll x=a0[i],y=a1[i],t=type==1?1:mod+1>>1;
		a0[i]=(x+y)*t%mod;
		a1[i]=(x-y+mod)*t%mod;
	}
}

void fwt_xnor(ll *a,int n,int type)
{
	if(n==1)return;
	int mid=n>>1;
	ll *a0=a,*a1=a+mid;
	fwt_xnor(a0,mid,type);
	fwt_xnor(a1,mid,type);
	for(int i=0;i<mid;i++)
	{
		ll x=a0[i],y=a1[i],t=type==1?1:mod+1>>1;
		a0[i]=(y-x+mod)*t%mod;
		a1[i]=(y+x)*t%mod;
	}
}

倍增迭代

同理，可用倍增写出 $\operatorname{xor}$ 和 $\operatorname{xnor}$ 更高效（常数更小）的迭代版本：

void fwt_xor(ll *a,int n,int type)
{
	for(int k=1;k<n;k<<=1)
	{
		for(int i=0;i<n;i+=k<<1)
		{
			for(int j=0;j<k;j++)
			{
				ll x=a[i+j],y=a[i+j+k],t=type==1?1:mod+1>>1;
				a[i+j]=(x+y)*t%mod;
				a[i+j+k]=(x-y+mod)*t%mod;
			}
		}
	}
}

void fwt_xnor(ll *a,int n,int type)
{
	for(int k=1;k<n;k<<=1)
	{
		for(int i=0;i<n;i+=k<<1)
		{
			for(int j=0;j<k;j++)
			{
				ll x=a[i+j],y=a[i+j+k],t=type==1?1:mod+1>>1;
				a[i+j]=(y-x+mod)*t%mod;
				a[i+j+k]=(y+x)*t%mod;
			}
		}
	}
}

性能分析

• 测试机器：CPU Apple M3 Max，RAM 36GB
• 测试环境：macOS Sequoia 15.6，clang 17.0.0
• 编译选项：-std=c++17 -O2
• 计时方法：调用 std::chrono::high_resolution_clock，单位为 $\mathrm{ms}$
• 数据生成：调用 std::mt19937_64，生成 $2^{24}\approx 1.68\times 10^7$ 个随机整数
• 测试方法：每个算法连续运行 $100$ 次，记录耗时并计算平均值

类型/算法	倍增迭代	分治递归	高维前缀和
或卷积（$\operatorname{or}$）	$575.31$	$665.67(+15.7\%)$	$966.18(+67.9\%)$
与卷积（$\operatorname{and}$）	$559.65$	$693.12(+23.8\%)$	$961.38(+71.8\%)$
异或卷积（$\operatorname{xor}$）	$939.02$	$1167.31(+24.3\%)$	$\text{N/A}$
同或卷积（$\operatorname{xnor}$）	$937.16$	$1164.39(+24.2\%)$	$\text{N/A}$

参考代码

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

using ll=long long;
const int mod=998244353;
const int N=20;
void fmt_or(ll *a,int n,int type)
{
	for(int k=1;k<n;k<<=1)
	{
		for(int i=0;i<n;i+=k<<1)
		{
			for(int j=0;j<k;j++)
			{
				a[i+j+k]=(a[i+j+k]+a[i+j]*type+mod)%mod;
			}
		}
	}
}
void fmt_and(ll *a,int n,int type)
{
	for(int k=1;k<n;k<<=1)
	{
		for(int i=0;i<n;i+=k<<1)
		{
			for(int j=0;j<k;j++)
			{
				a[i+j]=(a[i+j]+a[i+j+k]*type+mod)%mod;
			}
		}
	}
}
void fwt_xor(ll *a,int n,int type)
{
	for(int k=1;k<n;k<<=1)
	{
		for(int i=0;i<n;i+=k<<1)
		{
			for(int j=0;j<k;j++)
			{
				ll x=a[i+j],y=a[i+j+k],t=type==1?1:mod+1>>1;
				a[i+j]=(x+y)*t%mod;
				a[i+j+k]=(x-y+mod)*t%mod;
			}
		}
	}
}
ll A[1<<N],B[1<<N],a[1<<N],b[1<<N];
template<typename F>
void calc(F f,int n)
{
	for(int i=0;i<n;i++){a[i]=A[i];b[i]=B[i];}
	f(a,n,1);
	f(b,n,1);
	for(int i=0;i<n;i++)a[i]=a[i]*b[i]%mod;
	f(a,n,-1);
	for(int i=0;i<n;i++)cout<<a[i]<<' ';
	cout<<'\n';
}
int main()
{
	ios::sync_with_stdio(false);
	cin.tie(nullptr);
	int n;
	cin>>n;
	n=1<<n;
	for(int i=0;i<n;i++)cin>>A[i];
	for(int i=0;i<n;i++)cin>>B[i];
	calc(fmt_or,n);
	calc(fmt_and,n);
	calc(fwt_xor,n);
	return 0;
}

洛谷 P6097 【模板】子集卷积

洛谷 P6097 【模板】子集卷积

给定两个长度为 $2^n$ 的序列 $a_0,a_1,\dots,a_{2^n-1}$ 和 $b_0,b_1,\dots,b_{2^n-1}$，你需要求出一个序列 $c_0,c_1,\dots,c_{2^n-1}$，其中 $c_k$ 满足：

$$c_k=\sum_{\substack{{i \& j=0}\\{i~\mid~ j=k}}} a_i b_j$$

其中$~\mid~$表示按位或，$\&$表示按位与。

答案对 $10^9+9$ 取模。

代码（3）

倍增迭代

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

using ll=long long;
const int mod=1e9+9;
const int N=20;
void fmt(ll *a,int n,int type)
{
	for(int k=1;k<n;k<<=1)
	{
		for(int i=0;i<n;i+=k<<1)
		{
			for(int j=0;j<k;j++)
			{
				a[i+j+k]=(a[i+j+k]+a[i+j]*type+mod)%mod;
			}
		}
	}
}
ll a[N+1][1<<N],b[N+1][1<<N],c[N+1][1<<N];
int main()
{
	ios::sync_with_stdio(false);
	cin.tie(nullptr);
	int n;
	cin>>n;
	int m=1<<n;
	for(int i=0;i<m;i++)
	{
		cin>>a[__builtin_popcount(i)][i];
	}
	for(int i=0;i<m;i++)
	{
		cin>>b[__builtin_popcount(i)][i];
	}
	for(int i=0;i<=n;i++)
	{
		fmt(a[i],m,1);
		fmt(b[i],m,1);
	}
	for(int i=0;i<=n;i++)
	{
		for(int j=0;j<=i;j++)
		{
			for(int k=0;k<m;k++)
			{
				c[i][k]=(c[i][k]+a[j][k]*b[i-j][k])%mod;
			}
		}
	}
	for(int i=0;i<=n;i++)
	{
		fmt(c[i],m,-1);
	}
	for(int i=0;i<m;i++)
	{
		cout<<c[__builtin_popcount(i)][i]<<' ';
	}
	return 0;
}

分治递归

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

using ll=long long;
const int mod=1e9+9;
const int N=20;
void fmt(ll *a,int n,int type)
{
	if(n==1)return;
	int mid=n>>1;
	ll *a0=a,*a1=a+mid;
	fmt(a0,mid,type);
	fmt(a1,mid,type);
	for(int i=0;i<mid;i++)a1[i]=(a1[i]+a0[i]*type+mod)%mod;
}
ll a[N+1][1<<N],b[N+1][1<<N],c[N+1][1<<N];
int main()
{
	ios::sync_with_stdio(false);
	cin.tie(nullptr);
	int n;
	cin>>n;
	int m=1<<n;
	for(int i=0;i<m;i++)
	{
		cin>>a[__builtin_popcount(i)][i];
	}
	for(int i=0;i<m;i++)
	{
		cin>>b[__builtin_popcount(i)][i];
	}
	for(int i=0;i<=n;i++)
	{
		fmt(a[i],m,1);
		fmt(b[i],m,1);
	}
	for(int i=0;i<=n;i++)
	{
		for(int j=0;j<=i;j++)
		{
			for(int k=0;k<m;k++)
			{
				c[i][k]=(c[i][k]+a[j][k]*b[i-j][k])%mod;
			}
		}
	}
	for(int i=0;i<=n;i++)
	{
		fmt(c[i],m,-1);
	}
	for(int i=0;i<m;i++)
	{
		cout<<c[__builtin_popcount(i)][i]<<' ';
	}
	return 0;
}

高维前缀和

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

using ll=long long;
const int mod=1e9+9;
const int N=20;
void fmt(ll *a,int n,int type)
{
	for(int i=1;i<n;i<<=1)
	{
		for(int j=0;j<n;j++)
		{
			if(i&j)a[j]=(a[j]+a[i^j]*type+mod)%mod;
		}
	}
}
ll a[N+1][1<<N],b[N+1][1<<N],c[N+1][1<<N];
int main()
{
	ios::sync_with_stdio(false);
	cin.tie(nullptr);
	int n;
	cin>>n;
	int m=1<<n;
	for(int i=0;i<m;i++)
	{
		cin>>a[__builtin_popcount(i)][i];
	}
	for(int i=0;i<m;i++)
	{
		cin>>b[__builtin_popcount(i)][i];
	}
	for(int i=0;i<=n;i++)
	{
		fmt(a[i],m,1);
		fmt(b[i],m,1);
	}
	for(int i=0;i<=n;i++)
	{
		for(int j=0;j<=i;j++)
		{
			for(int k=0;k<m;k++)
			{
				c[i][k]=(c[i][k]+a[j][k]*b[i-j][k])%mod;
			}
		}
	}
	for(int i=0;i<=n;i++)
	{
		fmt(c[i],m,-1);
	}
	for(int i=0;i<m;i++)
	{
		cout<<c[__builtin_popcount(i)][i]<<' ';
	}
	return 0;
}