#include <math.h>
#define N 1048576
#define pi 3.1415926535897932384626
#define op cplx operator
typedef double fp;
struct cplx{
fp c,s;
cplx(fp _c=0,fp _s=0):c(_c),s(_s){}
op+(const cplx&w)const{return cplx(c+w.c,s+w.s);}
op-(const cplx&w)const{return cplx(c-w.c,s-w.s);}
op*(const cplx&w)const{return cplx(c*w.c-s*w.s,s*w.c+c*w.s);}
op-()const{return cplx(-c,-s);}
cplx mi()const{return cplx(-s,c);}
}A[N],B[N],C[N],E[N];
typedef cplx*const cc;
int rev[N];
void dft(cc X){
for(int j=0;j<N;j+=2){
cplx x=X[j],y=X[j+1];
X[j]=x+y;
X[j+1]=x-y;
}
for(int j=0;j<N;j+=4){
cplx x=X[j],y=X[j+2];
X[j]=x+y;
X[j+2]=x-y;
x=X[j+1],y=X[j+3].mi();
X[j+1]=x+y;
X[j+3]=x-y;
}
for(int j=0;j<N;j+=8){
cplx x=X[j],y=X[j+4];
X[j]=x+y;
X[j+4]=x-y;
x=X[j+1],y=X[j+5]*E[5];
X[j+1]=x+y;
X[j+5]=x-y;
x=X[j+2],y=X[j+6].mi();
X[j+2]=x+y;
X[j+6]=x-y;
x=X[j+3],y=X[j+7]*E[7];
X[j+3]=x+y;
X[j+7]=x-y;
}
for(int i=8;i<N;i<<=1){
for(int j=0;j<N;j+=i<<1){
cc f=X+j,g=f+i,e=E+i;
for(int k=0;k<i;++k){
cplx x=f[k],y=g[k]*e[k];
f[k]=x+y;
g[k]=x-y;
++k;
x=f[k],y=g[k]*e[k];
f[k]=x+y;
g[k]=x-y;
}
}
}
}
void cal(unsigned*a,int n,cc A){
for(int i=0;i<=n;i+=2)A[rev[i>>1]]=cplx(a[i],i==n?0:a[i+1]);
dft(A);
A[N]=A[0];
}
void trans(cplx&a0,cplx&a1,const cplx&v0,const cplx&v1){
a0=cplx(v0.c+v1.c,v0.s-v1.s);
a1=cplx(-v0.s-v1.s,v0.c-v1.c);
}
void poly_multiply(unsigned*a,int n,unsigned*b,int m,unsigned*c){
for(int i=1;i<N;++i)rev[i]=rev[i>>1]>>1|(i&1)<<19;
E[1]=cplx(1,0);
for(int i=2;i<N;i<<=1){
cc e0=E+i/2,e1=E+i;
cplx w(cos(pi/i),sin(pi/i));
for(int j=0;j<i;j+=2)e1[j]=e0[j>>1],e1[j+1]=e1[j]*w;
}
cal(a,n,A);
cal(b,m,B);
for(int i=0;i<N;++i){
cplx A0,A1,B0,B1;
trans(A0,A1,A[i],A[N-i]);
trans(B0,B1,B[i],B[N-i]);
cplx e=i<N/2?E[N/2+i]:-E[i];
C[rev[i]]=A0*B0+A1*B1*e+(A0*B1+A1*B0).mi();
}
for(int i=0;i<N;++i)C[i].s*=-1;
dft(C);
for(int i=0;i<=n+m;i+=2){
cplx x=C[i>>1]-C[N-1];
c[i]=int(x.c/(N*4)+.5);
if(i==n+m)break;
c[i+1]=int(x.s/(N*4)+.5);
}
}
Compilation | N/A | N/A | Compile OK | Score: N/A | 显示更多 |
Testcase #1 | 163.891 ms | 75 MB + 656 KB | Accepted | Score: 100 | 显示更多 |