快速沃尔什变换浅析
问题描述
给定长度为 \(2^n\) 的序列 \(A,B\)。设序列 \(C\) 满足 \(C_i=\sum_{j⊕k=i}A_j \times B_j\)。当二元运算 \(⊕\) 分别为 \(\text{or},\text{and},\text{xor}\) 时求序列 \(C\)。
关于生成函数与卷积
生成函数
简而言之，生成函数就是将长度为的序列 \(A\) 中的每一个元素按下标赋值到一个形式幂级数 \(F(x)=\sum_{i=0}^\infty A_ix^i\) 中。
卷积
卷积通常是对于两个多项式 \(F(x),G(x)\)（此处只讨论序列中的有限部分，故将形式幂级数看作多项式）生成的第三个数学算子。
对于常见的和卷积（即 \(\sum_{i+j=k}A_ix^i\)），常用快速傅里叶变换解决。而对于位运算的卷积就需要引入此处讨论的算法——快速沃尔什变换。
从FFT中得到的算法思想
变换与逆变换
对FFT有所了解的读者不难知道，FFT的算法流程是将系数表示法下的多项式变换为点值表示法，在 \(O(n)\) 的运算求出点值乘法结果后再通过逆变换将点值变换为系数。同样地，FWT也会有相似的算法流程。
为后文方便表述，此处定义 \(\text{FWT}(A)\) 表示对序列 \(A\) 进行快速沃尔什变换的结果。
系数分治
FFT时间复杂度能降低至 \(O(n\log n)\) 的主要原因，是利用折半引理与消去引理，将系数按奇偶性分类后用分治法将求点值的时间复杂度优化为对数时间复杂度。
因此，在FWT中，我们也希望出现类似的性质。即同样将一个序列分为可以在 \(O(1)\) 时间复杂度内相互转换的两个子序列。
卷积的部分性质
在卷积运算中，不同算子之间的交换律、结合律、分配率仍然成立。简略证明如下文（此处以定义 \(⊗\) 为序列间的卷积运算， \(⊕\) 为要求的系数二元运算）。
交换律：
结合律：
分配率：
快速沃尔什变换
与或运算卷积
此处先讨论或运算卷积的处理，与运算卷积同理。
对于三个整数 \(i,j,k\)，不难推导出有 \(j|i=i,k|i=i⇒(j|k)|i=i\)。所以得到如下等式：
对于三个序列的变换与逆变换只需要用类似FFT的方法系数分治即可，即如下过程（\(merge(A,B)\) 表示将序列合并，\(len\) 为序列长度）：
将序列分为下标最高位为 \(0\) 和 \(1\) 的两个子序列，并不断往下递归，就完成了变换过程；再将该过程中的符号取反，就得到逆变换过程。
异或运算卷积
首先引入一个结论：
其中 \(\text{popcount}(x)\) 表示整数 \(x\) 所对应二进制数中 \(1\) 的数量。
通俗一点说，就是异或运算不改变两个整数对应二进制数中 \(1\) 数量的奇偶性。
由于在异或运算中需要计算 \(C_i=\sum_{j⊕k=i}A_jB_k\)，所以需要构造 \(\text{FWT}(C)=\text{FWT}(A)⊗\text{FWT}(B)\)。所以引入一种构造方式 \(\text{FWT}(A)_{i=0}^{n}=g(x,i)A_i\)。而该构造方式需要满足如下：
易得 \(g(x,j⊕k)=g(x,j)g(x,k)\)。
这条性质刚好与引入的结论相似，所以不难推知 \(g(x,i)=(-1)^{\text{popcount}(i\&x)}\)。
因此我们就得到了异或运算下的变换方式：
异或运算下的快速沃尔什变换表达式为：
算法总结与补充
一些注意事项
由于所有变换过程中 \(\text{FWT}(A)\) 与 \(\text{IFWT}(A)\) 只有变换符号的区别，所以实现过程中可以添加一个操作命令符将变换与逆变换整合为一个函数。
洛谷的模板题数据有点大，需要开long long。
后记
由于笔者数学功底有限（加之期末时间较为紧张），算法推导过程并不尽详细。笔者拟在假期后对推导和证明部分作部分修改。
实现代码
#include<iostream>
#include<cstdio>
#include<cstdlib>
#include<algorithm>
#include<cmath>
#include<cstring>
#include<queue>
#include<set>
#include<ctime>
#define ll long long
#define RL register long long
using namespace std;
const ll maxn=2e5+10,mod=998244353ll;
ll lev,n;
ll A[maxn],B[maxn],aph[maxn],bet[maxn];
inline ll gmd(ll x);
inline ll qpow(ll x,ll y);
inline ll inv(ll x);
inline void copy(ll len,ll *ori,ll *ret);/*数组赋值*/
inline void output(ll len,ll *f);/*输出模块*/
inline void calc(ll len,ll *a,ll *b);/*点值运算*/
inline void OR(ll len,ll *f,ll opt);/*对或运算的FWT*/
inline void AND(ll len,ll *f,ll opt);/*对与运算的FWT*/
inline void XOR(ll len,ll *f,ll opt);/*对异或运算的FWT*/
int main()
{
ios::sync_with_stdio(false);cin.tie(0);cout.tie(0);
cin >> lev;n=1<<lev;
for(RL i=1;i<=n;++i) cin >> A[i];
for(RL i=1;i<=n;++i) cin >> B[i];
copy(n,A,aph);copy(n,B,bet);OR(n,aph,1);OR(n,bet,1);calc(n,aph,bet);OR(n,aph,-1);output(n,aph);
copy(n,A,aph);copy(n,B,bet);AND(n,aph,1);AND(n,bet,1);calc(n,aph,bet);AND(n,aph,-1);output(n,aph);
copy(n,A,aph);copy(n,B,bet);XOR(n,aph,1);XOR(n,bet,1);calc(n,aph,bet);XOR(n,aph,inv(2));output(n,aph);
return 0;
}
inline ll gmd(ll x)
{
return (x%mod+mod)%mod;
}
inline ll qpow(ll x,ll y)
{
RL ret=1;
for(;y;x=gmd(x*x),y>>=1)
if(y&1) ret=gmd(ret*x);
return ret;
}
inline ll inv(ll x)
{
return qpow(x,mod-2);
}
inline void copy(ll len,ll *ori,ll *ret)
{
for(RL i=1;i<=len;++i) ret[i]=ori[i];
return;
}
inline void output(ll len,ll *f)
{
for(RL i=1;i<=len;++i) cout << f[i] << " ";
cout << endl;
return;
}
inline void calc(ll len,ll *a,ll *b)
{
for(RL i=1;i<=len;++i) a[i]=gmd(a[i]*b[i]);
return;
}
inline void OR(ll len,ll *f,ll opt)
{
for(RL k=1;(k<<1)<=len;k<<=1)
for(RL i=1;i<=len;i+=(k<<1))
for(RL j=0;j<k;++j)
f[i+j+k]=gmd(f[i+j+k]+gmd(f[i+j]*opt));
return;
}
inline void AND(ll len,ll *f,ll opt)
{
for(RL k=1;(k<<1)<=len;k<<=1)
for(RL i=1;i<=len;i+=(k<<1))
for(RL j=0;j<k;++j)
f[i+j]=gmd(f[i+j]+gmd(f[i+j+k]*opt));
return;
}
inline void XOR(ll len,ll *f,ll opt)
{
for(RL k=1;(k<<1)<=len;k<<=1)
for(RL i=1;i<=len;i+=(k<<1))
for(RL j=0;j<k;++j)
{
f[i+j]=f[i+j]+f[i+j+k];
f[i+j+k]=f[i+j]-(f[i+j+k]<<1);
f[i+j]=gmd(f[i+j]*opt);f[i+j+k]=gmd(f[i+j+k]*opt);
}
return;
}
