FFT 快速傅里叶变换

来源:互联网 发布:pullandbear怎么样知乎 编辑:程序博客网 时间:2024/06/06 03:44

http://blog.miskcoo.com/2015/04/polynomial-multiplication-and-fast-fourier-transform#i-10
唔看了算导和一个讲的很清楚的博客
我只是做了下笔记,尽量简洁一点写吧。

这样复杂度就是

T(n)=2T(n/2)+O(n)=O(n log n)

好像发现了惊天大秘密QwQ 它们的位置是颠倒的诶!
怎么计算rev数组呢?

for (int i=1;i<n;i++)    rev[i]=rev[i>>1]>>1 | ((i&1) << l-1);

再看一下迭代的过程

void fft(complex<double> *a){    for (int i=1;i<n;i++)        if (i<rev[i])            swap(a[i],a[rev[i]]);    for (int len=2;len<=n;len<<=1){        complex<double> wn(cos(2*PI/len),sin(2*PI/len));        for (int j=0;j<n;j+=len){            complex<double> w(1,0),x;            for (int k=0;k<(len>>1);k++,w*=wn)                x=w*a[j+k+(len>>1)],                a[j+k+(len>>1)]=a[j+k]-x,                a[j+k]+=x;        }    }}

最后贴一下代码(UOJ 34 模板题)

#include <bits/stdc++.h>#define N 300010#define PI acos(-1)using namespace std;int n,m,l,rev[N];complex<double> a[N],b[N],c[N];template <class Aqua>inline void read(Aqua &s){    s=0; char c=getchar();    while (!isdigit(c)) c=getchar();    while (isdigit(c)) s=s*10+c-'0',c=getchar();}void pre(){    for (m+=n,n=1,l=0;n<=m;n<<=1,l++);    for (int i=1;i<n;i++)        rev[i]=rev[i>>1]>>1 | ((i&1) << l-1);}void fft(complex<double> *a,int f){    for (int i=1;i<n;i++)        if (i<rev[i])            swap(a[i],a[rev[i]]);    for (int len=2;len<=n;len<<=1){        complex<double> wn(cos(2*PI/len),f*sin(2*PI/len));        for (int j=0;j<n;j+=len){            complex<double> w(1,0),x;            for (int k=0;k<(len>>1);k++,w*=wn)                x=w*a[j+k+(len>>1)],                a[j+k+(len>>1)]=a[j+k]-x,                a[j+k]+=x;        }    }}int main(){    read(n),read(m);    for (int i=0;i<=n;i++)        read(a[i].real());    for (int i=0;i<=m;i++)        read(b[i].real());    pre();    fft(a,1); fft(b,1);    for (int i=0;i<=n;i++)        c[i]=a[i]*b[i];    fft(c,-1);    for (int i=0;i<m;i++)        printf("%d ",(int)(c[i].real()/n+0.5));    printf("%d\n",(int)(c[m].real()/n+0.5));    return 0;}

讲完啦
运用下次再写吧QwQ

原创粉丝点击