[学习笔记]FFT——快速傅里叶变换

最新推荐文章于 2024-09-06 16:54:11 发布
最新推荐文章于 2024-09-06 16:54:11 发布
阅读量146 收藏
点赞数
原文链接：http://www.cnblogs.com/Miracevin/p/9994262.html
版权

大力推荐博客：

傅里叶变换(FFT)学习笔记

 

一、多项式乘法：


我们要明白的是：

FFT利用分治，处理多项式乘法，达到O(nlogn)的复杂度。（虽然常数大）

FFT=DFT+IDFT

DFT:

本质是把多项式的系数表达转化为点值表达。因为点值表达，y可以直接相乘。点值表达下相乘的复杂度是O(n)的。

我们分别对两个多项式求x为$\omega_n^i$时的y值。

然后可以O(n)求出乘积多项式x为$\omega_n^i$时的y值。

求法：

把F(x)奇偶分类。

$FL(x)=a_0+a_2x+...+a_{n-2}x^{n/2-1}$

$FR(x)=a_1+a_3x+...+a_{n-1}x^{n/2-1}$

$F(x)=FL(x^2)+xFR(x^2)$


带入那些神奇的单位根之后，
发现有：

$0<=k<n/2$

$F(\omega_n^k)=Fl(\omega_{n/2}^k)+\omega_{n}^kFR(\omega_{n/2}^k)$


$F(\omega_n^{k+n/2})=Fl(\omega_{n/2}^k)-\omega_{n}^kFR(\omega_{n/2}^k)$

我们只要知道Fl、FR多项式在那n/2个位置的点值，那么就可以知道F那n个位置的点值了。

分治就可以处理出来。



IDFT：

经过一系列矩阵的运算之后，，，，

可以得到：

$b_k=[(ω_n^{-k})^0y_0+(ω_n^{-k})^1y_1+(ω_n^{-k})^2y_2+...+(ω_n^k)^{n-1}y_{n-1}]/n$

可以把y当做系数，

只要知道，当x是一系列w的时候，值是多少。

那么就求出来了$b_k$

FFT再写一遍。

注意这里带入的是$ω_n^{-k}$

开始的$tmp$有所不同

 

// luogu-judger-enable-o2
#include<bits/stdc++.h>
#define reg register int
#define il inline
#define numb (ch^'0')
using namespace std;
typedef long long ll;
il void rd(int &x){
    char ch;x=0;bool fl=false;
    while(!isdigit(ch=getchar()))(ch=='-')&&(fl=true);
    for(x=numb;isdigit(ch=getchar());x=x*10+numb);
    (fl==true)&&(x=-x);
}
namespace Miracle{
const int N=1e6+5;
const double Pi=acos(-1);
struct node{
    double x,y;
    node(){}
    node(double xx,double yy){x=xx,y=yy;}
    node operator +(const node &b){
        return node(x+b.x,y+b.y);
    }
    node operator -(const node &b){
        return node(x-b.x,y-b.y);
    }
    node operator *(const node &b){
        return node(x*b.x-y*b.y,x*b.y+y*b.x);
    }
}a[4*N],b[4*N];
int n,m;
int r[4*N];
void FFT(node *f,short op){
    for(reg i=0;i<n;++i){
        if(i<r[i]){
            node tmp=f[i];
            f[i]=f[r[i]];
            f[r[i]]=tmp;
        }
    }
    for(reg p=2;p<=n;p<<=1){
        int len=p/2;
        node tmp(cos(Pi/len),op*sin(Pi/len));
        for(reg k=0;k<n;k+=p){
            node buf(1,0);
            for(reg l=k;l<k+len;++l){
                node tt=buf*f[l+len];
                f[l+len]=f[l]-tt;
                f[l]=f[l]+tt;
                buf=buf*tmp;
            }
        }
    }
}
int main(){
    scanf("%d%d",&n,&m);
    for(reg i=0;i<=n;++i) scanf("%lf",&a[i].x);
    for(reg i=0;i<=m;++i) scanf("%lf",&b[i].x);
    for(m=n+m,n=1;n<=m;n<<=1);
    for(reg i=0;i<n;++i){
        r[i]=r[i>>1]>>1|((i&1)?n>>1:0);
    }
    FFT(a,1);FFT(b,1);
    for(reg i=0;i<n;++i) b[i]=a[i]*b[i];
    FFT(b,-1);
    for(reg i=0;i<=m;++i) printf("%.0lf ",fabs(b[i].x)/n);
    return 0;
}

}
int main(){
    Miracle::main();
    return 0;
}

/*
   Author: *Miracle*
   Date: 2018/11/21 8:05:13
*/
多项式乘法

 

 

关键点就是在于，用了单位根这个东西，可以避免平方、避免爆long long 以及精度损失的情况下，再利用乘法分配律，可以O(nlogn)得到多项式的点值表达。

 

例题：

P3338 [ZJOI2014]力

 

思路：要用FFT，必然要化成多项式卷积的形式

即形如：$h[j]=\sum_{i=0}^j f[i]*g[j-i]$

这样的话，我们把f，g分别作为两个多项式的系数，那么，发现，h[j]的值，恰好是f,g两个多项式乘积得到的多项式的第j+1项的系数。（考虑次数j是怎么来的）

就可以FFT优化这个n^2的算式了。

 

对于这个题：

令$f[i]=q[i]$,$g[i]=\frac{1}{i*i}$

特别的；有$g[0]=0$

则有$E[j]=\sum_{i=0}^jf[i]*g[j-i]-\sum_{i=j}^nf[i]*g[i-j]$

我们可以分开算，

后面的减法部分类似一个后缀，把$f$数组$reverse$一下，就变成了前缀了。$g$数组不用，因为距离要保持这样。

于是；

$E[j]=\sum_{i=0}^jf[i]*g[j-i]-\sum_{i=0}^{n-j}f'[i]*g[n-j-i]$

两次$FFT$即可

 

值得注意的是：

1.g数组赋值的时候，i*i可能会爆int，导致精度误差。所以，写成1/i/i比1/(i*i)要好得多。（30pts->100pts）

2.乘积多项式一定要n+n项都算出来，因为最后的插值和每一项的点值都有关系。即使我们只关心前n项。

 

代码：

#include<bits/stdc++.h>
#define reg register int
#define il inline
#define numb (ch^'0')
using namespace std;
typedef long long ll;
il void rd(int &x){
    char ch;x=0;bool fl=false;
    while(!isdigit(ch=getchar()))(ch=='-')&&(fl=true);
    for(x=numb;isdigit(ch=getchar());x=x*10+numb);
    (fl==true)&&(x=-x);
}
namespace Miracle{
const int N=200000+5;
const double Pi=acos(-1);
struct node{
    double x,y;
    node(){}
    node(double xx,double yy){
        x=xx;y=yy;
    }
}f[2*N],g[2*N],h[2*N];
double q[2*N];
int r[2*N];
int n,m;
node operator+(const node &a,const node &b){
    return node(a.x+b.x,a.y+b.y);
}
node operator-(const node &a,const node &b){
    return node(a.x-b.x,a.y-b.y);
}
node operator*(const node &a,const node &b){
    return node(a.x*b.x-a.y*b.y,a.x*b.y+a.y*b.x);
}
void FFT(node *f,short op){
    for(reg i=0;i<n;++i){
        if(i<r[i]){
            node tmp=f[i];
            f[i]=f[r[i]];
            f[r[i]]=tmp;
        }
    }
    for(reg p=2;p<=n;p<<=1){
        int len=p/2;
        node tmp=node(cos(Pi/len),op*sin(Pi/len));
        for(reg k=0;k<n;k+=p){
            node buf=node(1,0);
            for(reg l=k;l<k+len;++l){
                node tt=buf*f[l+len];
                f[l+len]=f[l]-tt;
                f[l]=f[l]+tt;
                buf=buf*tmp;
            }
        }
    }
}
int main(){
    scanf("%d",&m);
    for(reg i=1;i<=m;++i){
        scanf("%lf",&q[i]);
        if(i)g[i]=node((double)1/(double)i/(double)i,0);
    }
    for(n=1;n<=2*m;n<<=1);
    //cout<<" nn "<<n<<endl;
    for(reg i=0;i<n;++i){
        f[i]=node(q[i],0);
        //cout<<f[i].x<<" ";
    }
    
    //g[0]=node(0,0);
    for(reg i=0;i<n;++i){
        r[i]=(r[i>>1]>>1)|((i&1)?(n>>1):0);
    }
    
    FFT(f,1);
    FFT(g,1);
    for(reg i=0;i<n;++i) f[i]=g[i]*f[i];
    FFT(f,-1);
    
    
    
    reverse(q+1,q+n);
    for(reg i=0;i<n;++i){
        h[i]=node(q[i],0);
    }
    FFT(h,1);
    for(reg i=0;i<n;++i) h[i]=h[i]*g[i];
    FFT(h,-1);
    
    for(reg i=1;i<=m;++i){
        node ans=f[i]-h[n-i];
        printf("%lf\n",ans.x/n);
    }
    return 0;
}

}
int main(){
    Miracle::main();
    return 0;
}

/*
   Author: *Miracle*
   Date: 2018/11/21 10:17:15
*/
力

 

FFT优化高精乘法：

把数字看成系数，把10^k看做是x^k，那么就可以得到多项式。

这两个多项式相乘，得到的多项式，各个系数通过进位变成个位数之后，直接输出系数即可。

值得注意的是：

浮点数四舍五入赋值：

$a=floor(b+0.5);$

#include<bits/stdc++.h>
#define reg register int
#define il inline
#define numb (ch^'0')
using namespace std;
typedef long long ll;
il void rd(int &x){
    char ch;x=0;bool fl=false;
    while(!isdigit(ch=getchar()))(ch=='-')&&(fl=true);
    for(x=numb;isdigit(ch=getchar());x=x*10+numb);
    (fl==true)&&(x=-x);
}
namespace Miracle{
const int N=60000+5;
const double Pi=acos(-1);
struct node{
    double x,y;
    node(){}
    node(double xx,double yy){
        x=xx;y=yy;
    }
}a[4*N],b[4*N];
char p[N],q[N];
int c[4*N];
int n,m;
int r[4*N];
node operator+(node a,node b){
    return node(a.x+b.x,a.y+b.y);
}
node operator-(node a,node b){
    return node(a.x-b.x,a.y-b.y);
}
node operator*(node a,node b){
    return node(a.x*b.x-a.y*b.y,a.x*b.y+a.y*b.x);
}
void FFT(node *f,short op){
    for(reg i=0;i<n;++i){
        if(i<r[i]){
            node tmp=f[i];
            f[i]=f[r[i]];
            f[r[i]]=tmp;
        }
    }
    for(reg p=2;p<=n;p<<=1){
        int len=p/2;
        node tmp=node(cos(Pi/len),op*sin(Pi/len));
        for(reg k=0;k<n;k+=p){
            node buf=node(1,0);
            for(reg l=k;l<k+len;++l){
                node tt=buf*f[l+len];
                f[l+len]=f[l]-tt;
                f[l]=f[l]+tt;
                buf=buf*tmp;
            }
        }
    }
}
int main(){
    scanf("%d",&m);
    scanf("%s",p);scanf("%s",q);
    for(reg i=0;i<m;++i){
        a[m-i-1].x=p[i]-'0';
        b[m-i-1].x=q[i]-'0';
    }
    for(m=m+m,n=1;n<m;n<<=1);
    for(reg i=0;i<n;++i){
        r[i]=r[i>>1]>>1|((i&1)?n>>1:0);
    }
    FFT(a,1);FFT(b,1);
    
    for(reg i=0;i<n;++i) b[i]=a[i]*b[i];
    FFT(b,-1);
    for(reg i=0;i<n;++i){
        c[i]=floor(b[i].x/n+0.5);
    }
    
    int x=0;
    for(reg i=0;i<n;++i){
        c[i]+=x;
        x=(int)c[i]/10;
        c[i]%=10;
    }
    while(c[n-1]==0&&n>=1) --n;
    for(reg i=n-1;i>=0;--i){
        printf("%d",c[i]);
    }
    return 0;
}

}
int main(){
    Miracle::main();
    return 0;
}

/*
   Author: *Miracle*
   Date: 2018/11/21 16:30:14
*/
FFT高精

 

转载于:https://www.cnblogs.com/Miracevin/p/9994262.html

weixin_33817333
关注
fft变换理解
youzjuer的博客
01-02 648
1、N点的离散傅里叶变换X（k）（其取值范围是k=,0，1，...N-1）当中每一个均可由N个原始序列 x(n)与旋转因子进行一系列的乘累加操作来获得。具体过程如下： 其中：
快速傅里叶变换（FFT）学习笔记
Jun
07-27 605
Step1 模板理解: 1.背景：完成两个多项式的卷积 即：两多项式形如∑i=0naixi\sum_{i=0}^n a_{i}x^{i}∑i=0n​ai​xi ，∑i=0mbixi\sum_{i=0}^m b_{i}x^{i}∑i=0m​bi​xi 卷积后为 ∑i=0n+m(ci=∑j=0iajbi−j)xi\sum_{i=0}^{n+m}( c_{i}=\sum_{j=0}^ia_{j}b_{i...
【学习笔记】快速傅里叶变换(FFT)
weixin_30432179的博客
03-17 303
【学习笔记】快速傅里叶变换 学习之前先看懂这个 浅谈范德蒙德(Vandermonde)方阵的逆矩阵的求法以及快速傅里叶变换(FFT)中IDFT的原理——gzy hhh开个玩笑。 讲一下\(FFT\) 的流程，我也不准备长篇大论地分析\(FFT...\) 将系数表示法转换为点值表示法 \(O(n \log n)​\) 对于点值表示法直接进行操作 \(O(n)\) 将点值表示法转换为系数表示法 \(...
FFT(快速傅里叶变换)学习小记
Cold_Chair的博客
03-28 615
磕了三天，也算是磕了出来，还需要强化一下。谢谢yl的ppt。问题引入：51nod 大数乘法 V2 题目就是大数乘法，n&lt;=10^5。一个科学的做法：把a,b视作多项式。先分别做点值运算，再做插值运算搞回来就行了。但是我们知道正常的点值和插值都需要O(n^2)的时间去完成。FFT就是利用特殊的x去优化这个过程。单位复数根：n次单位复数根是满足wn=1w^n=1的复数。n次单位复数根恰好有n个，它们均
FFT变换
u012632105的博客
12-30 150
FFT变换
【学习笔记】极其美妙的算法——FFT（快速傅里叶变换）
LeeCongWei的博客
07-20 4477
文章目录前言FFT是什么？多项式表示法系数表示法点值表示法初步想法问题解决问题2复数定义运算复平面单位圆 前言 最近学会了FFT，一个极其美妙的算法！ 众所周知FFT英文全称是 Fast Fast TLE Fast Fourier Transform FFT是什么？ FFT（快速傅里叶变换）是DFT（离散傅里叶变换）的优化，用于加速多项式乘法。 即： 设 F(x)=∑i=0nfixiF(x)=\sum_{i=0}^nf_ix^iF(x)=∑i=0n​fi​xi , G(x)=∑i=0ngixiG(x)
【学习笔记】超简单的快速傅里叶变换（FFT）（含全套证明）
繁凡さん的博客
12-25 1万+
整理的算法模板合集： ACM模板 目录多项式复数复数的单位根 / 单位向量离散傅里叶变换（DTF）离散傅里叶逆变换（IDTF）FFT算法整体流程代码实现递归版迭代版“三步变两步”优化P1919 【模板】FFT快速傅里叶变换 学习笔记（用ipad记录的笔记）： 超级简单的快速傅里叶变换！只要基础够扎实，顺着推一遍没有什么难以理解的，我学的整个过程没有一点卡壳，真的很爽，整整写了三个小时，写了满满6页的笔记。主要是内容太多了。 首先是一些基础概念： 多项式 复数 复数的单位根 / 单位向量 离散傅里叶变
【数字信号处理 | 学习笔记】二、离散傅里叶变换及其快速算法
qq_56734650的博客
07-06 3321
1. 傅里叶变换的几种可能形式：2. 离散傅里叶级数3. 周期序列的离散傅里叶级数变换对其中是周期为N的离散周期信号。4. DFS与Z变换的关系周期序列DFS可以看作是的一个周期x(n)做Z变换，再将Z变换在Z平面单位圆上按等间隔角抽样而得到。 周期卷积与线性卷积的不同之处在于周期卷积仅在一个周期内求和。离散傅里叶变换DFT实际上是DFS的主值，DFT隐含了周期性。x(n)的N点DFT对应x(n)的Z变换在单位圆上N个等间隔点的采样。 计算过程：补零，周期延拓，翻褶，移位取主值序列，相乘相加。循环卷积是周期
学习笔记 - 位逆序置换 / 非递归FFT / 快速数论变换(NTT)
qq_26087481的博客
03-28 2308
文章讨论了快速傅立叶变换的非递归实现以及快速数论变换的实现.
快速傅里叶变换 FFT 学习笔记
weach的博客
08-24 393
概述 离散傅里叶变换（Discrete Fourier Transform，缩写为 DFT），是傅里叶变换在时域和频域上都呈离散的形式，将信号的时域采样变换为其 DTFT 的频域采样。 FFT 是一种高效实现 DFT 的算法，称为快速傅立叶变换（Fast Fourier Transform，FFT）。它对傅里叶变换的理论并没有新的发现，但是对于在计算机系统或者说数字系统中应用离散傅立叶变换，可以说是进了一大步。快速数论变换 (NTT) 是快速傅里叶变换（FFT）在数论基础上的实现。 在 1965 年，Coo
快速傅里叶变换算法与应用
08-21
快速傅里叶变换算法与应用，电子书，好不容易找到，希望对你的学习工作有帮助，如果觉得不错，请支持正版。
第一次邂逅快速傅立叶变换(FFT)
海纳百川
04-19 3万+
 为了毕业设计，我要学习JPEG，还有视频压缩技术，在JPEG的时候，我就被前面的DCT给挡住了，现如今我终于写了一个FFT程序，发了我好长的时间。如果说是因为我的无知，还是什么，我对学习这类有关数学的东西，总是显得那么的迟钝，也许是因为人老了吧。其它我还像个小孩子一样，唉，这年头，还真是搞不懂自己了。进入正题吧，我对FFT的完全不了解，到最后，实现FFT正变换与反变换，其中有太多的细节，如果我现
信号处理之快速傅里叶变换(FFT)-通俗易懂
热门推荐
红叶幽香的博客
03-19 4万+
历史溯源 欧拉公式 傅里叶级数(FS) 傅里叶变换(FT) 离散傅里叶级数(DFS) 离散时间傅里叶变换(DTFT) 离散傅里叶变换(DFT) 快速傅里叶变换(FFT) MATLAB中常用的FFT工具 FFT中常见的问题(频率混叠，频谱泄露，栅栏效应，旁瓣效应)
快速傅里叶变换（FFT）及其在多项式乘法中的应用 —— 深入分析与 Python 实现
最新发布
qq_52010229的博客
09-06 1577
单位根是指在复平面上均匀分布在单位圆上的点，满足 x^n = 1 的复数点。它们可以表示为：这些点的对称性和周期性使得我们可以通过分治法高效计算傅里叶变换。通过本文，我们深入探讨了如何利用快速傅里叶变换FFT进行高效的多项式乘法。我们从最基本的多项式操作出发，介绍了如何通过单位根和分治法将乘法复杂度从 O(n^2) 降低到 O(n log n)。我们还通过 Python 代码展示了 FFT 和 IFFT 的具体实现，并进行了多项式乘法的示例。FFT 在多项式乘法中的应用只是其众多用途之一。
简单的理解FFT变换
wangyu088的专栏
07-05 3644
任何一个复杂的信号都可以分解成若干个正弦信号，分解的过程采用的是FFT变换。 一个正弦信号的表达式如下： ），其中A表示幅值，f表示频率，表示相偏移。 FFT变换通过对模拟信号进行采样，得到该模拟信号组成的每个正弦信号的幅值，频率，相位，从而复原该模拟信号。 FFT变换对每个采样点的值采用的是复数表示，也即，对应局部较大值的点表示的是一个正弦信号的频率，例如某个信号采样的第N个点是局部最...
FFT （快速傅里叶变换）
海边拾贝，沧海一粟
07-28 7225
FFT 标签（空格分隔）： FFT 模板 FFT FFT的全称是 FAST Fourier Transform 即快速傅里叶变换 傅里叶变换是复变函数的重要内容，傅里叶变换分为离散和连续傅里叶变换 傅里叶变换实现从时域到频域的转换，是信号与系统重要的分析工具 连续傅里叶变换 离散傅里叶变换(DFT discrete Fourier transform) 将离散...
FFT变换-学习笔记
qq_42640275的博客
04-09 3690
1.多项式相乘 2.用d+1个点来表示d阶多项式，怎么取这些点 3.利用奇偶性取点，并把这些点取成复数，然后递归 4.将需求的点比作1开n次方的根，也就是复平面单位圆上的n等分的点，引入旋转因子(欧拉公式可推)。 5.将多项式中的参数换成单位圆上的旋转因子，然后递归。 会发现单位圆上成相反数的点平方一次后还是一对相反数的点，只是点数减少一半，所以可以一直递归，直到最后剩一个点，那么反过来由一个点往n个点推出来的就是1开n次方的根，也就是FFT的结果。 为了计算快速，对于FFT尽量取2的幂次个点，
FFT（快速傅里叶变换）
fzyyyyyyyy的博客
05-07 8308
FFT
评论
添加红包

请填写红包祝福语或标题

个

红包个数最小为10个

元

红包金额最低5元

当前余额3.43元 前往充值 >
需支付：10.00元
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付元
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明：

1.余额是钱包充值的虚拟货币，按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载，可以购买VIP、付费专栏及课程。

余额充值