代码

#include<bits/stdc++.h>
#define rep(i,l,r) for (int i=l; i<=r; ++i)
#define rep2(i,l,r) for (int i=l; i<r; ++i)

const int N=(1<<21),n2=(1<<20); int rev[N];
struct C { double Re,Im; }w,w2,tmp,f[N];
inline C operator + (C a,C b) { return {a.Re+b.Re,a.Im+b.Im}; }
inline C operator - (C a,C b) { return {a.Re-b.Re,a.Im-b.Im}; }
inline C operator * (C a,C b) { return {a.Re*b.Re-a.Im*b.Im,a.Re*b.Im+a.Im*b.Re}; }

void FFT() {
    rep2(i,0,N) if (i<rev[i]) std::swap(f[i],f[rev[i]]);
    double t=acos(-1);
    for (int i=2; i<=N; i<<=1,t/=2) {
        w={cos(t),sin(t)};
        for (int j=0; j<N; j+=i) {
            w2={1,0}; const int mid=j+(i>>1);
            for (int k1=j,k2=mid; k1<mid; ++k1,++k2,w2=w2*w)
                tmp=w2*f[k2],f[k2]=f[k1]-tmp,f[k1]=f[k1]+tmp;
        }
    }
}

void poly_multiply(unsigned *a,int n,unsigned *b,int m,unsigned *c) {
    rep(i,0,1000000) f[i]={a[i],b[i]};
    rep2(i,0,n2) rev[i<<1|1]=(rev[i<<1]=rev[i]>>1)|n2;
    FFT(); rep2(i,0,N) f[i]=f[i]*f[i];
    FFT(); std::reverse(f+1,f+N);
    rep(i,0,2000000) c[i]=f[i].Im/N/2+0.5;
}

评测结果