博客
关于我
强烈建议你试试无所不能的chatGPT,快点击我
FFT
阅读量:5118 次
发布时间:2019-06-13

本文共 2222 字,大约阅读时间需要 7 分钟。

 这里是板子,虽然看懂了原理,但是代码还是好难理解哦
void fft(int n,complex<double>*buffer,int offset,int step,complex<double>* epsilon)     
{                                                                                       
if(n==1) return;                                                                       
int m=n>>1;                                                                            
fft(m,buffer,offset,step<<1,epsilon);                                                  
fft(m,buffer,offset+step,step<<1,epsilon);                                             
for(int k=0;k!=m;++k)                                                                  
{                                                                                      
int pos=2*step*k;                                                                     
temp[k]=buffer[pos+offset]+epsilon[k*step]*buffer[pos+offset+step];                   
temp[k+m]=buffer[pos+offset]-epsilon[k*step]*buffer[pos+offset+step];                 
}                                                                                      
for(int i=0;i!=n;++i)                                                                  
buffer[i*step+offset]=temp[i];                                                        
}                                                                                       
void init_epsilon(int n)                                                                
{                                                                                       
double pi=acos(-1);                                                                    
for(int i=0;i!=n;++i)                                                                  
{                                                                                      
epsilon[i]=complex<double>(cos(2.0*pi*i/n),sin(2.0*pi*i/n));                          
arti_epsilon[i]=conj(epsilon[i]);                                                     
}                                                                                      
}                                                                                       
int reverse_add(int x)                                                                  
{                                                                                       
for(int l=1<<bit_length;(x^=l)<l;l>>=1);                                               
return x;                                                                              
}                                                                                       
/* 这时候 n 已经补齐到 2 的幂次 */                                                      
void bit_reverse(int n, complex_t *x)                                                   
{                                                                                       
for(int i=0,j=0;i!=n;++i)                                                              
{                                                                                      
if(i>j) swap(x[i],x[j]);                                                              
for(int l=n>>1;(j^=l)<l;l>>=1);                                                       
}                                                                                      
}                                                                                       
void transform(int n,complex_t *x,complex_t *w)                                         
{                                                                                       
bit_reverse(n, x);                                                                     
for(int i=2;i<=n;i<<=1)                                                                
{                                                                                      
int m=i>>1;                                                                           
for(int j=0;j<n;j+=i)                                                                 
{                                                                                     
for(int k=0;k!=m;++k)                                                                
{                                                                                    
complex_t z=x[j+m+k]*w[n/i*k];                                                      
x[j+m+k]=x[j+k]-z;                                                                  
x[j+k]+=z;                                                                          
}                                                                                    
}                                                                                     
}                
}
 这个是带注释的版本
#include<cstdio>
#include<cstdlib>
#include<cstring>
#include<cmath>
#include<iostream>
#include<algorithm>
using namespace std;
#define N 301000
const double pi=acos(-1);
struct node
{
    double x,y;
    node(){x=y=0;}
    node(double x,double y):x(x),y(y){}
}a[N],b[N];
node operator + (node x,node y) {return node(x.x+y.x,x.y+y.y);}
node operator - (node x,node y) {return node(x.x-y.x,x.y-y.y);}
node operator * (node x,node y) {return node(x.x*y.x-x.y*y.y,x.x*y.y+x.y*y.x);}
void fft(node *s,int n,int t)
{
    if (n==1) return;
    node a0[n>>1],a1[n>>1];
    for (int i=0;i<=n;i+=2) 
     a0[i>>1]=s[i],a1[i>>1]=s[i+1];
    fft(a0,n>>1,t);fft(a1,n>>1,t);
    node wn(cos(2*pi/n),t*sin(2*pi/n)),w(1,0);
    for (int i=0;i<(n>>1);i++,w=w*wn) 
     s[i]=a0[i]+w*a1[i],s[i+(n>>1)]=a0[i]-w*a1[i];
}
int main()
{
    int n,m,fn,i;
    scanf("%d%d",&n,&m);
    for (i=0;i<=n;i++) scanf("%lf",&a[i].x);
    for (i=0;i<=m;i++) scanf("%lf",&b[i].x);
    fn=1;while (fn<=n+m) fn<<=1;
    fft(a,fn,1);fft(b,fn,1);
    for (i=0;i<=fn;i++) a[i]=a[i]*b[i];
    fft(a,fn,-1);
    for (i=0;i<=n+m;i++) printf("%d ",(int)(a[i].x/fn+0.5));
    printf("\n");
    return 0;
}

转载于:https://www.cnblogs.com/OcahIBye/p/6729320.html

你可能感兴趣的文章
Java泛型的基本使用
查看>>
1076 Wifi密码 (15 分)
查看>>
rsync
查看>>
noip模拟赛 党
查看>>
bzoj2038 [2009国家集训队]小Z的袜子(hose)
查看>>
Java反射机制及其Class类浅析
查看>>
Postman-----如何导入和导出
查看>>
移动设备显示尺寸大全 CSS3媒体查询
查看>>
图片等比例缩放及图片上下剧中
查看>>
【转载】Linux screen 命令详解
查看>>
background-clip,background-origin
查看>>
Android 高级UI设计笔记12:ImageSwitcher图片切换器
查看>>
Blog文章待看
查看>>
【Linux】ping命令详解
查看>>
对团队成员公开感谢博客
查看>>
java学习第三天
查看>>
python目录
查看>>
django+uwsgi+nginx+sqlite3部署+screen
查看>>
Andriod小型管理系统(Activity,SQLite库操作,ListView操作)(源代码下载)
查看>>
在Server上得到数据组装成HTML后导出到Excel。两种方法。
查看>>