FFT模板

最近很忙,有空了补描述.


FFT递归形式

const int maxn=1e6;
void FFT(complex<double> *a,int d,int inv)
{
    if(d==10) return;
    int mid=d/2;
    complex<double> p1[maxn],p2[maxn];
    for(int i=0;i<mid;i++)
    {
        p1[i]=a[i*2];
        p2[i]=a[i*2+1];
    }
    FFT(p1,mid,inv);
    FFT(p2,mid,inv);
    double sitar=(2*acos(-1.0))/d;
    complex<double> w0(1,0),w(cos(sitar),inv*sin(sitar));
    for(int k=0;k<mid;k++,w0*=w)
    {
        a[k]=p1[k]+w0*p2[k];
        a[k+d/2]=p1[k]-w0*p2[k];
    }
}

FFT的非递归形式

const int maxn=1e6;
int rev[maxn];
void FFT(complex<double> *a,int d,int inv)
{
    for(int i=0;i<d;i++)
    {
        rev[i]=(rev[i>>1]>>1)|((i&1)<<(int(log2(d)-1)));
        if(i<rev[i])
        {
            complex<double> temp=a[i];
            a[i]=a[rev[i]];
            a[rev[i]]=temp;
        }
    }
    complex<double> p1,p2;
    for(int mid=1;mid<d;mid*=2)
    {
        double sitar=acos(-1.0)/mid;
        complex<double> w(cos(sitar),inv*sin(sitar));
        for(int j=0;j<d;j+=(mid*2))
        {
            complex<double> w0(1,0);
            for(int k=0;k<mid;k++,w0*=w)
            {
                p1=a[j+k],p2=a[j+mid+k];
                a[j+k]=p1+w0*p2;
                a[j+mid+k]=p1-w0*p2;
            }
        }
    }
}

img_show