在多项式卷积的处理中,我们实际上实现的是下面的一个式子
然而事实上有些和(sang)蔼(xin)可(bing)亲(kuang)的出题人,并不会让你直接求这样的式子,比如把+换成-什么的
但是更多时候,你拿到手上的却是这样一个式子
其中(igoplus)可能是(or,and,xor)中的一种或多种
那么现在我们就会使用FWT(快速沃尔什变化)来优化时间复杂度
前置芝士
FFT
理解思想即可,不在赘述
正文
序
(p.s. 注意下文的乘号应结合具体语境理解,是表示按位相乘还是位运算卷积)
考虑我们在实现FFT的时候干了什么?
我们搞出了一个点值表示法,然后直接对点值进行求乘积,再把它转回去
那么这对位运算是否成立呢?
即我们希望有(FWT(C)=FWT(A) imes FWT(B)),此时可以直接按位相乘
直觉告诉我们是有的,但是是什么样的呢?
or
你要问我是怎么找的我也无从说起,十分诡异
想要了解更多的可以参考:https://www.cnblogs.com/ACMLCZH/p/8022502.html
这里直接给出结论及其正确性的证明
对于一个有(2^n)项的(A_i)
什么是(A_0)和(A_1)?
我们把(A)想象成一个(2^n)维的向量,那么前(2^{n-1})项是(A_0),后(2^{n-1})项是(A_1)
这之后的括号表示的意义同样(在之后会阐述(FWT(A))的本质)
如果只是感性理解这玩意还是比较容易的:对于(A_0)中的每一项,它们的二进制表示最高位均是(0),当它们与(A_1)的某一项进行了(or)卷积的话,那么最高位就变成了1,对后面的(2^{n-1})项产生了贡献
那么问题就变成了如何证明这时的(FWT(C)=FWT(A) imes FWT(B))
我们考虑使用数学归纳法证明,设此时(A,B,C)均为(2^n)项
若(n=0),结论显然成立
若(n>0),假设结论在(n=k)时成立
那么此时就应该有(注意这里用(*)表示(or)卷积)
(理解一下,在括号里面的表示位运算卷积,外面的则为按位相乘)
同时显然(FWT(A+B)=FWT(A)+FWT(B))(+表示按位相加)
首先是等式左边
这也很好理解,和上面(FWT)的公式一样理解
再来看等式左边
证毕
好像结束了?让我们来看一下(FWT(A))的本质
回到最开始的式子(c_i=sum_{j|k=i}a_jb_k),假设我们记(C_i=sum _{jsubseteq i}c_j),则
我们发现(A,B,C)就起到了上文中(FWT)的作用,于是当我们要求某个集合的子集中的元素和时,除了(O(3^n))的枚举子集以外,还有(O(n2^n))的求出数组(A)的(FWT(A))的方法,(A_i)的(i)的子集和就是(FWT(A)_i)
and
证明和上面是类似的,在这里直接给出结论
以及我们注意到,当求出(FWT(A))时,我们有(FWT(A)_i=sum_{isubseteq j}A_j),它和上面的(or)卷积在其他除了FWT的题目中也有很好的应用
xor
依然是只给出结论,证明留给读者
IDFT相关
我们实现了FFT中类似于DFT的一步,那么IDFT呢?
直接按照上面的逆运算做回去就可以了,直接上结论:
对于(or)卷积
对于(and)卷积
对于(xor)卷积
code
模板题:luogu 4717
#include<iostream>
#include<string.h>
#include<string>
#include<stdio.h>
#include<algorithm>
#include<math.h>
#include<vector>
#include<queue>
#include<map>
using namespace std;
typedef long long ll;
const ll maxd=998244353,inv2=499122177;
const double pi=acos(-1.0);
int n,N;
ll a[1001000],b[1001000],c[1001000];
int read()
{
int x=0,f=1;char ch=getchar();
while ((ch<'0') || (ch>'9')) {if (ch=='-') f=-1;ch=getchar();}
while ((ch>='0') && (ch<='9')) {x=x*10+(ch-'0');ch=getchar();}
return x*f;
}
void fwt_or(int lim,ll *a,int typ)
{
int mid;
for (mid=1;mid<lim;mid<<=1)
{
int len=(mid<<1),sta,j;
for (sta=0;sta<lim;sta+=len)
{
for (j=0;j<mid;j++)
{
if (typ==1) a[sta+mid+j]=(a[sta+j+mid]+a[sta+j])%maxd;
else a[sta+mid+j]=(a[sta+j+mid]+maxd-a[sta+j])%maxd;
}
}
}
}
void fwt_and(int lim,ll *a,int typ)
{
int mid;
for (mid=1;mid<lim;mid<<=1)
{
int len=(mid<<1),sta,j;
for (sta=0;sta<lim;sta+=len)
{
for (j=0;j<mid;j++)
{
if (typ==1) a[sta+j]=(a[sta+j]+a[sta+mid+j])%maxd;
else a[sta+j]=(a[sta+j]+maxd-a[sta+j+mid])%maxd;
}
}
}
}
void fwt_xor(int lim,ll *a,int typ)
{
int mid;
for (mid=1;mid<lim;mid<<=1)
{
int len=(mid<<1),sta,j;
for (sta=0;sta<lim;sta+=len)
{
for (j=0;j<mid;j++)
{
int x=a[sta+j],y=a[sta+j+mid];
a[sta+j]=(x+y)%maxd;
a[sta+j+mid]=(x+maxd-y)%maxd;
if (typ==-1) {a[sta+j]=(a[sta+j]*inv2)%maxd;a[sta+j+mid]=(a[sta+j+mid]*inv2)%maxd;}
}
}
}
}
int main()
{
n=read();N=1<<n;int i;
for (i=0;i<N;i++) a[i]=read();
for (i=0;i<N;i++) b[i]=read();
fwt_or(N,a,1);fwt_or(N,b,1);
for (i=0;i<N;i++) c[i]=(a[i]*b[i])%maxd;
fwt_or(N,a,-1);fwt_or(N,b,-1);fwt_or(N,c,-1);
for (i=0;i<N;i++) printf("%lld ",c[i]);printf("
");
fwt_and(N,a,1);fwt_and(N,b,1);
for (i=0;i<N;i++) c[i]=(a[i]*b[i])%maxd;
fwt_and(N,a,-1);fwt_and(N,b,-1);fwt_and(N,c,-1);
for (i=0;i<N;i++) printf("%lld ",c[i]);printf("
");
fwt_xor(N,a,1);fwt_xor(N,b,1);
for (i=0;i<N;i++) c[i]=(a[i]*b[i])%maxd;
fwt_xor(N,a,-1);fwt_xor(N,b,-1);fwt_xor(N,c,-1);
for (i=0;i<N;i++) printf("%lld ",c[i]);printf("
");
return 0;
}