FFT

const double pi=acos(-1.0);
int rev[N];

void FFT(complex<double> *a,int nr,int flag){
    for(int i=0;i<nr;i++){
        if(i<rev[i]) swap(a[i],a[rev[i]]);
    }
    for(int i=1;i<nr;i<<=1){
        complex<double> wn(cos(pi/i),flag*sin(pi/i));
        for(int j=0;j<nr;j+=(i<<1)){
            complex<double> w0(1,0);
            for(int k=0;k<i;k++,w0*=wn){
                complex<double> x=a[j+k],y=w0*a[i+j+k];
                a[j+k]=x+y;
                a[i+j+k]=x-y;
            }
        }
    }
}

多项式乘法

complex<double> f[N];

void multiply(double *a,double *b,int n,int m){
    for(int i=0;i<n;i++) f[i].real(a[i]);
    for(int i=0;i<m;i++) f[i].imag(b[i]);
    int l=max((int)ceil(log2(n+m)),1),len=1<<l;
    for(int i=0;i<len;i++) rev[i]=(rev[i>>1]>>1)|((i&1)<<(l-1));
    FFT(f,len,1);
    for(int i=0;i<=len;i++) f[i]=f[i]*f[i];
    FFT(f,len,-1);
    for(int i=0;i<=n+m;i++) a[i]=f[i].imag()/2.0/len+0.5;
}