FFT介绍

写在前面

详细内容请参见《算法导论》3rd edition,或者参见博客
FFT算法学习笔记
这里只展示它的迭代形式的代码

IterativeFFT

#include <iostream>
#include<cmath>
#include<cstdio>
#include<complex>
#include<algorithm>


using namespace std;

typedef complex<double> CLDD;
const double PI = acos(-1);

int rev(int id,int len)//向量长度，为2^k
{
    int res = 0;
    for(int i=0 ; (1<<i)<len ; ++i)
    {
        res<<=1;
        if(id& (1<<i))res|=1;
    }
    return res;
}
//DFT为-1，表示逆变换
void FFT(CLDD *a,int len,CLDD *ans,int DFT)
{
    for(int i=0 ; i<len ; ++i)
        ans[rev(i,len)] = a[i];
    for(int s = 1 ; (1<<s)<=len ; ++s)
    {
        int m = (1<<s);
        CLDD wm = CLDD(cos(DFT*2*PI/m),sin(DFT*2*PI/m));
        for(int k =0 ; k<len ; k+=m)
        {
            CLDD w = CLDD(1,0);
            for(int j=0 ; j<(m>>1) ; ++j)
            {
                CLDD u = ans[k+j];CLDD t = w*ans[k+j+(m>>1)];
                ans[k+j] = u+t;
                ans[k+j+(m>>1)] = u-t;
                w = w*wm;
            }
        }
    }
    if(DFT==-1)for(int i=0 ; i<len ; ++i)ans[i]/=len;
}


int main()
{

    CLDD a[4];
    CLDD A[4];
    int len = 4;
    for(int i=0  ; i<len ; ++i)
        a[i]=CLDD(i,0);
    printf("\n----------rev--------------\n");
    for(int i=0 ; i<len ; ++i)
        cout<<rev(i,len)<<endl;;

    FFT(a,len,A,1);
    for(int i=0 ; i<len ; ++i)
        cout<<A[i]<<endl;
    FFT(A,len,a,-1);
    printf("\n-------DFT-1---------\n");
    for(int i=0 ; i<len ;++i)
        cout<<a[i]<<endl;
    return 0;
}

结果。