FFT-快速傅里叶变换
一,时间复杂度
时间复杂度: (mathcal{O(nlogn)})
二，限制
- 计算过程(即 多项式相乘过程)不能取模；
- 常数较大；
- 会存在精度差。
三,模板
计算多项式 (C=A imes B), 其中 多项式 (A) 的长度是 (n+1), 多项式 (B) 的长度是 (m+1) ,得到多项式 (C) 的长度是 (n+m+1)。
#include <bits/stdc++.h>
using namespace std;
const int maxn=2e6+5;
const double PI=acos(-1);
struct comp{
double x,y;
comp operator +(comp b){return comp{x+b.x,y+b.y};}
comp operator -(comp b){return comp{x-b.x,y-b.y};}
comp operator *(comp b){return comp{x*b.x-y*b.y,x*b.y+y*b.x};}
};
int up,L,R[maxn<<1];
void fft_init(int n,int m){
up=1;L=0;
while(up<n+m+2)up<<=1,L++;
for(int i=0;i<up;i++)
R[i]=(R[i>>1]>>1)|((i&1)<<(L-1));
}
comp a[maxn<<1],b[maxn<<1],c[maxn<<1];
void fft(comp a[],int type)
{
for(int i=0;i<up;i++)
if(i<R[i])swap(a[i],a[R[i]]);
for(int mid=1;mid<up;mid<<=1)
{
comp wn=comp{cos(PI/mid),type*sin(PI/mid)};
for(int r=mid<<1,j=0;j<up;j+=r)
{
comp w{1,0};
for(int k=0;k<mid;k++,w=w*wn)
{
comp x=a[j+k],y=w*a[j+mid+k];
a[j+k]=x+y;
a[j+mid+k]=x-y;
}
}
}
}
int main()
{
int n,m;
scanf("%d%d",&n,&m);
for(int i=0;i<=n;i++)scanf("%lf",&a[i].x);
for(int i=0;i<=m;i++)scanf("%lf",&b[i].x);
fft_init(n,m);
fft(a,1);fft(b,1);
for(int i=0;i<=up;i++)c[i]=a[i]*b[i];
fft(c,-1);
for(int i=0;i<=n+m;i++)printf("%d ",(int)(c[i].x/up+0.5));
}
NTT-快速数论变换
一, 时间复杂度
时间复杂度: (mathcal{O(nlogn)})
二,限制
-
系数必须是整数；
-
系数取模时, 模数有限制, 需要知道模数的原根。
-
常见模数对应原根：
998244353 ( o) 3
100000007 ( o) 5
-
三,模板
计算多项式 (C=A imes B), 其中 多项式 (A) 的长度是 (n+1), 多项式 (B) 的长度是 (m+1) ,得到多项式 (C) 的长度是 (n+m+1)。
#include <bits/stdc++.h>
using namespace std;
typedef long long ll;
const int maxn=2e6+5;
const int mod=998244353;
const int g=3;
const int gi=332748118;//invg
ll kpow(ll a,ll b){
ll ans=1;
while(b){
if(b&1)ans=ans*a%mod;
a=a*a%mod;
b>>=1;
}
return ans;
}
int up,L,R[maxn<<1];
void ntt_init(int n,int m){
up=1;L=0;
while(up<n+m+2)up<<=1,L++;
for(int i=0;i<up;i++)
R[i]=(R[i>>1]>>1)|((i&1)<<(L-1));
}
int a[maxn<<1],b[maxn<<1],c[maxn<<1];
void ntt(int a[],int type)
{
for(int i=0;i<up;i++)
if(i<R[i])swap(a[i],a[R[i]]);
for(int mid=1;mid<up;mid<<=1){
ll wn=kpow(type==1?g:gi,(mod-1)/(mid<<1));
for(int r=mid<<1,j=0;j<up;j+=r)
{
ll w=1;
for(int k=0;k<mid;k++,w=w*wn%mod)
{
int x=a[j+k],y=w*a[j+k+mid]%mod;
a[j+k]=(x+y)%mod;
a[j+k+mid]=(x-y+mod)%mod;
}
}
}
}
int main()
{
int n,m;
scanf("%d%d",&n,&m);
for(int i=0;i<=n;i++)scanf("%d",&a[i]);
for(int i=0;i<=m;i++)scanf("%d",&b[i]);
ntt_init(n,m);
ntt(a,1);ntt(b,1);
for(int i=0;i<up;i++)c[i]=1ll*a[i]*b[i]%mod;
ntt(c,-1);
int inv=kpow(up,mod-2);
for(int i=0;i<=n+m;i++)printf("%lld ",1ll*c[i]*inv%mod);
}
FWT-快速沃尔什变换
一,作用
FFT可以解决多项式卷积, 即:
[C_k=sum_{k=i+j}\,A_i*B_j
]
FWT可以解决 或/与/异或 卷积, 即:
[C_k=sum_{k=i|j}\,A_i*B_j\
C_k=sum_{k=i&j}\,A_i*B_j\
C_k=sum_{k=iigotimes j}\,A_i*B_j
]
二,时间复杂度
时间复杂度: (mathcal{O}(nlogn))
三,做法
1, 将 (A) 转换成 (FWT[A]) , 将 (B) 转换成 (FWT[B])
2, 计算
[FWT[C][i]=FWT[A][i]*FWT[B][i]
]
3, 将 (FWT[C]) 转换成 (C) .
四,模板
(code:)
void FWT_or(long long a[],int len)//A -> FAT[A]
{
for(int mid=2;mid<=len;mid<<=1)
for(int i=0;i<len;i+=mid)
for(int j=i;j<i+(mid>>1);j++)a[j+(mid>>1)]+=a[j];
}
void IFWT_or(long long a[],int len)// FAT[A] -> A
{
for(int mid=2;mid<=len;mid<<=1)
for(int i=0;i<len;i+=mid)
for(int j=i;j<i+(mid>>1);j++)a[j+(mid>>1)]-=a[j];
}
void FWT_and(long long a[],int len)
{
for(int mid=2;mid<=len;mid<<=1)
for(int i=0;i<len;i+=mid)
for(int j=i;j<i+(mid>>1);j++)a[j]+=a[j+(mid>>1)];
}
void IFWT_and(long long a[],int len)
{
for(int mid=2;mid<=len;mid<<=1)
for(int i=0;i<len;i+=mid)
for(int j=i;j<i+(mid>>1);j++)a[j]-=a[j+(mid>>1)];
}
void FWT_xor(long long a[],int len)
{
for(int mid=2;mid<=len;mid<<=1)
for(int i=0;i<len;i+=mid)
for(int j=i;j<i+(mid>>1);j++)
{
long long x=a[j],y=a[j+(mid>>1)];
a[j]=x+y,a[j+(mid>>1)]=x-y;
}
}
inline void IFWT_xor(long long a[],int len)
{
for(int mid=2;mid<=len;mid<<=1)
for(int i=0;i<len;i+=mid)
for(int j=i;j<i+(mid>>1);j++)
{
long long x=a[j],y=a[j+(mid>>1)];
a[j]=(x+y)>>1,a[j+(mid>>1)]=(x-y)>>1;
}
}
五,例题
题意:
给定长度为 (n) 的数组 (A) ，其中 (1le nle 2 imes 10^5,\,0le A_i<2^{18}) 。
对于所有的 (1le ile n) ，要求：
从数组 (A) 中取出 (i) 个数(可重复取同一个数) , 计算出 (a_1igotimes a_2igotimes...igotimes a_i) 的值。
题解:
当 (i> 19) 之后 , 有 (ans_i=ans_{i-2}) 。
只需用异或卷积求出前 (19) 个答案。
代码:
#include<bits/stdc++.h>
#define mem(a,b) memset(a,b,sizeof(a))
using namespace std;
typedef long long ll;
const int mod=1e9+7;
const int inf=0x3f3f3f3f;
const int maxn=1e6+5;
void read(int&x)
{
char c;
while(!isdigit(c=getchar()));x=c-'0';
while(isdigit(c=getchar()))x=(x<<3)+(x<<1)+c-'0';
}
void FWT_xor(int a[],int n){
for(int mid=2;mid<=n;mid<<=1)
for(int i=0;i<n;i+=mid)
for(int j=i,x,y;j<i+(mid>>1);j++)
{
x=a[j];y=a[j+(mid>>1)];
a[j]=x+y;a[j+(mid>>1)]=x-y;
}
}
void IFWT_xor(int a[],int n){
for(int mid=2;mid<=n;mid<<=1)
for(int i=0;i<n;i+=mid)
for(int j=i,x,y;j<i+(mid>>1);j++)
{
x=a[j];y=a[j+(mid>>1)];
a[j]=(x+y)>>1;a[j+(mid>>1)]=(x-y)>>1;
}
}
int a[maxn],b[maxn],c[maxn];
int ans[maxn];
int main()
{
int n,up,ed=1<<18;
read(n);
for(int i=1,x;i<=n;i++){
read(x);
a[x]=1;ans[1]=max(ans[1],x);
}
up=min(n,19);
FWT_xor(a,ed);
for(int i=0;i<ed;i++)b[i]=a[i];
for(int i=2;i<=up;i++){
for(int j=0;j<ed;j++)
c[j]=a[j]*b[j];
IFWT_xor(c,ed);
for(int j=0;j<ed;j++)
if(c[j])ans[i]=j,c[j]=1;
FWT_xor(c,ed);
for(int j=0;j<ed;j++)b[j]=c[j];
}
for(int i=up+1;i<=n;i++)ans[i]=ans[i-2];
for(int i=1;i<=n;i++)printf("%d ",ans[i]);putchar(10);
}