1.基本介绍
代码语言:javascript
复制//495ms
#include <bits/stdc .h>
#define rep(i,l,r) for(int i=l,ed=r;i<ed; i)
typedef long long ll;
const double PI = acos(-1);
const int N = 1<<20;
const int BUF_SIZE=33554431;
using namespace std;
struct buf{
char a[BUF_SIZE],b[BUF_SIZE],*s,*t;
buf():s(a),t(b){a[fread(a,1,sizeof a,stdin)]=0;}
~buf(){fwrite(b,1,t-b,stdout);}
operator int(){
int x=0;
while(*s<48) s;
while(*s>32)
x=x*10 *s -48;
return x;
}
void out(int x){
static char c[12];
char*i=c;
if(!x)*t =48;
else{
while(x){
int y=x/10;
*i =x-y*10 48,x=y;
}
while(i!=c)*t =*--i;
}
*t =10;
}
}it;
struct cp{
double x,y;
cp(double _x=0,double _y=0):x(_x),y(_y){}
cp operator (const cp& b)const{return cp(x b.x,y b.y);}
cp operator -(const cp& b)const{return cp(x-b.x,y-b.y);}
cp operator *(const cp& b)const{return cp(x*b.x-y*b.y,x*b.y y*b.x);}
cp operator !()const{return cp(x,-y);}
}w[N];
void fft(cp p[],int n){
for(int i=0,j=0;i<n; i){
if(i>j)swap(p[i],p[j]);
for(int l=n>>1;(j^=l)<l;l>>=1);
}
for(int i=2;i<=n;i<<=1)
for(int j=0,m=i>>1;j<n;j =i)
rep(k,0,m){
cp b=w[n/i*k]*p[j m k];
p[j m k]=p[j k]-b;
p[j k]=p[j k] b;
}
}
void conv(int n,ll *x,ll *y,ll *z){
static cp p[N],q[N],h(0,-0.25);
rep(i,0,n){
w[i]=cp(cos(2*PI*i/n),sin(2*PI*i/n));
p[i]=cp(x[i],y[i]);
}
fft(p,n);
rep(i,0,n){
int j=i?(n-i):0;
q[j]=(p[i]*p[i]-!p[j]*!p[j])*h;
}
fft(q,n);
rep(i,0,n)z[i]=q[i].x/n 0.5;
}
int n,m,p;
ll a[N],b[N],c[N];
int main(){
n=it 1;m=it 1;
rep(i,0,n) a[i]=it;
rep(i,0,m) b[i]=it;
for(n =m-1,p=1;p<n;p<<=1);
conv(p,a,b,c);
rep(i,0,n)it.out(c[i]);
return 0;
}
2.更快的卷积
代码语言:javascript
复制//325ms
#include <bits/stdc .h>
#define rep(i,l,r) for(int i=l,ed=r;i<ed; i)
typedef long long ll;
const double PI = acos(-1);
const int N = 1<<20;
const int BUF_SIZE=33554431;
using namespace std;
struct buf{
char a[BUF_SIZE],b[BUF_SIZE],*s,*t;
buf():s(a),t(b){a[fread(a,1,sizeof a,stdin)]=0;}
~buf(){fwrite(b,1,t-b,stdout);}
operator int(){
int x=0;
while(*s<48) s;
while(*s>32)
x=x*10 *s -48;
return x;
}
void out(int x){
static char c[12];
char*i=c;
if(!x)*t =48;
else{
while(x){
int y=x/10;
*i =x-y*10 48,x=y;
}
while(i!=c)*t =*--i;
}
*t =10;
}
}it;
struct cp{
double x,y;
cp(double _x=0,double _y=0):x(_x),y(_y){}
cp operator (const cp& b)const{return cp(x b.x,y b.y);}
cp operator -(const cp& b)const{return cp(x-b.x,y-b.y);}
cp operator *(const cp& b)const{return cp(x*b.x-y*b.y,x*b.y y*b.x);}
cp operator *(double b)const{return cp(b*x,b*y);}
cp operator !()const{return cp(x,-y);}
}w[N];
void fft(cp *p,int n){
for(int i=0,j=0;i<n; i){
if(i>j)swap(p[i],p[j]);
for(int l=n>>1;(j^=l)<l;l>>=1);
}
for(int i=2;i<=n;i<<=1)
for(int j=0,m=i>>1;j<n;j =i)
rep(k,0,m){
cp b=w[n/i*k]*p[j m k];
p[j m k]=p[j k]-b;
p[j k]=p[j k] b;
}
}
void conv(int n,ll *x,ll *y,ll *z){
static cp p[N],q[N],a[N];
rep(i,0,n){
(i&1?p[i>>1].y:p[i>>1].x)=x[i];
(i&1?q[i>>1].y:q[i>>1].x)=y[i];
}
rep(i,0,n>>=1)w[i]=cp(cos(2*PI*i/n),sin(2*PI*i/n));
fft(p,n);fft(q,n);
rep(i,0,n){
int j=i?n-i:0;
a[j]=p[i]*q[i]-((cp(1,0) w[i])*(p[i]-!p[j])*(q[i]-!q[j]))*0.25;
}
fft(a,n);
rep(i,0,n)z[i<<1]=a[i].x/n 0.5,z[i<<1|1]=a[i].y/n 0.5;
}
int n,m,p;
ll a[N],b[N],c[N];
int main(){
n=it 1;m=it 1;
rep(i,0,n) a[i]=it;
rep(i,0,m) b[i]=it;
for(n =m-1,p=2;p<n;p<<=1);
conv(p,a,b,c);
rep(i,0,n)it.out(c[i]);
return 0;
}
3.拆系数FFT
代码语言:javascript
复制//933ms
#include <bits/stdc .h>
#define rep(i,l,r) for(int i=l,ed=r;i<ed; i)
typedef long long ll;
const double PI = acos(-1);
const int N = 1<<20;
const ll mod = 1e9 7;
const int BUF_SIZE=33554431;
using namespace std;
struct buf{
char a[BUF_SIZE],b[BUF_SIZE],*s,*t;
buf():s(a),t(b){a[fread(a,1,sizeof a,stdin)]=0;}
~buf(){fwrite(b,1,t-b,stdout);}
operator int(){
int x=0;
while(*s<48) s;
while(*s>32)
x=x*10 *s -48;
return x;
}
void out(int x){
static char c[12];
char*i=c;
if(!x)*t =48;
else{
while(x){
int y=x/10;
*i =x-y*10 48,x=y;
}
while(i!=c)*t =*--i;
}
*t =10;
}
}it;
struct cp{
double x,y;
cp(double _x=0,double _y=0):x(_x),y(_y){}
cp operator (const cp& b)const{return cp(x b.x,y b.y);}
cp operator -(const cp& b)const{return cp(x-b.x,y-b.y);}
cp operator *(const cp& b)const{return cp(x*b.x-y*b.y,x*b.y y*b.x);}
cp operator !()const{return cp(x,-y);}
}w[N];
void fft(cp p[],int n){
for(int i=0,j=0;i<n; i){
if(i>j)swap(p[i],p[j]);
for(int l=n>>1;(j^=l)<l;l>>=1);
}
for(int i=2;i<=n;i<<=1)
for(int j=0,m=i>>1;j<n;j =i)
rep(k,0,m){
cp b=w[n/i*k]*p[j m k];
p[j m k]=p[j k]-b;
p[j k]=p[j k] b;
}
}
void conv(int n,ll *x,ll *y,ll *z){
static cp p[N],q[N],a[N],b[N],c[N],d[N];
static cp r(0.5,0),h(0,-0.5),o(0,1);
rep(i,0,n){
w[i]=cp(cos(2*PI*i/n),sin(2*PI*i/n));
x[i]=(x[i] mod)%mod,y[i]=(y[i] mod)%mod;
p[i]=cp(x[i]>>15,x[i]&32767),q[i]=cp(y[i]>>15,y[i]&32767);
}
fft(p,n);fft(q,n);
rep(i,0,n){
int j=i?(n-i):0;
static cp ka,ba,kb,bb;
ka=(p[i] !p[j])*r;
ba=(p[i]-!p[j])*h;
kb=(q[i] !q[j])*r;
bb=(q[i]-!q[j])*h;
a[j]=ka*kb;b[j]=ka*bb;
c[j]=kb*ba;d[j]=ba*bb;
}
rep(i,0,n){
p[i]=a[i] b[i]*o;
q[i]=c[i] d[i]*o;
}
fft(p,n);fft(q,n);
rep(i,0,n){
ll a,b,c,d;
a=(ll)(p[i].x/n 0.5)%mod;
b=(ll)(p[i].y/n 0.5)%mod;
c=(ll)(q[i].x/n 0.5)%mod;
d=(ll)(q[i].y/n 0.5)%mod;
z[i]=((a<<30) ((b c)<<15) d)%mod;
}
}
int n,m,p;
ll a[N],b[N],c[N];
int main(){
n=it 1;m=it 1;
rep(i,0,n) a[i]=it;
rep(i,0,m) b[i]=it;
for(n =m-1,p=1;p<n;p<<=1);
conv(p,a,b,c);
rep(i,0,n)it.out((c[i] mod)%mod);
return 0;
}