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

lailai0916 · 2025-08-06 19:00:55 · 题解

参考资料

快速沃尔什变换 - 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。

解题思路

约定

数据规模统一记为：

m=2^n

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

位运算符号统一记为：

按位或：\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

变换的逆过程称为 逆变换（Inverse Transform），亦称反演（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)=O(m\log^2 m)。

快速沃尔什变换（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;
}