1. 激活谷首页
  2. 激活谷笔记

fft 应用领域_fft应用领域

• 激活谷笔记

fft 应用领域_fft应用领域FFT(快速傅里叶变化) 介绍&相关应用1.前言本人有幸在算法学习过程中接触该算法 , 在了解算法背后的原理后 , 不禁为先人的巧妙设计感到倾佩.不管以后是算法竞赛还是项目研究 , FFT是一种多领域的、十分常用的多项式乘积优化算法.就本人

FFT(快速傅里叶变化) 介绍&相关应用 　　1.前言 　　本人有幸在算法学习过程中接触该算法 , 在了解算法背后的原理后 , 不禁为先人的巧妙设计感到倾佩. 　　不管以后是算法竞赛还是项目研究 , FFT是一种多领域的、十分常用的多项式乘积优化算法. 　　就本人的理解 , 在这里写下一篇文章好好介绍这种算法. 　　2.问题的产生 　　我们思考这么一个问题： 　　对于给定的两个多项式
f(x)=\sum_{i=0}^{n}{a_i x^i} 与
g(x)=\sum_{i=0}^{m}{b_ix^i} , 试计算
h(x)=f(x)· g(x) 　　显然 , 我们可以计算出
h(x) 中
x^i 项的系数
c_i , 就可以确定
h(x) 的表达式 　　根据幂运算 ,
c_i=\sum_{j=0}^{i}{a_jb_{i-j}} 　　我们可以用一层循环枚举
i , 然后再用一层循环枚举
j , 即可确定所有
c_i 　　显然这种算法的时间复杂度为
O(n·m) (因为每个
a_i 要与所有
b_j 乘一次) 　　当
n,m 的数量级很大时,这种算法的计算效率显然会比较低. 　　有没有任何方法降低计算的时间复杂度呢? 　　3.重新整理思考方向 　　我们不妨从确定函数的方式上做文章. 　　在数值分析这门课程上 , 我们知道 , 任何一个
n 次多项式函数都可以由平面中
n+1 个横坐标互异的点个唯一确定. 　　因此 , 如果我们可以找到函数
h(x) 上任意
n+1 个互不重合的点,那么我们就可以确定唯一的
h(x) 　　但是 , 如果我们仅仅就是随便选取互不相等的
x_0,x_1...x_n , 然后计算
h(x_i)=f(x_i)·g(x_i) ,并用
(x_0,h(x_0)),(x_1,h(x_1))···(x_n,h(x_n)) 来间接表示
h(x) ,我们的时间复杂度还是保持在
O((n+1)·max(n,m)) \approx O(n·m) 的水平,仿佛问题又回到了原点. 　　真的回到了原点嘛?我们不妨思考这样一个问题: 　　现在有一个程序内载
n 次多项式偶函数
h(x) ,它能够输出每次输入的值对应的函数值,那么,你至少需要输入几次数据,根据它的输出,才能保证得出具体的表达式? 　　答案是
\frac{n}{2} 次,这是因为如果我们算出了
h(\alpha) 的值,那么根据偶函数的性质,可以得到图像上必有
(\alpha,h(\alpha)),(-\alpha,h(\alpha)) 两点 　　显然这点对于奇函数也是成立的,对于奇函数
h(x) ,如果我们算出了
h(\alpha) ,那么图像上必有
(\alpha,h(\alpha)),(-\alpha,-h(\alpha)) 两点 　　那我们可不可以利用这种奇偶性质来简化原问题中
h(x_i) 计算过程呢？ 　　4.深入思考与推导 　　对于待确定函数
h(x)=\sum_{i=0}^{N}{c_ix^i} (
N 为奇数) 　　可以改写为
h(x) = \sum_{i=0}^{(N-1)/2}{c_{2i}x^{2i}} + \sum_{i=0}^{(N-1)/2}{c_{2i+1}x^{2i+1}} 　　即
h(x)=\sum_{i=0}^{(N-1)/2}{c_{2i}x^{2i}} +x \sum_{i=0}^{(N-1)/2}{c_{2i+1}x^{2i}} 　　上式中,
\sum_{i=0}^{(N-1)/2}{c_{2i}x^{2i}},\sum_{i=0}^{(N-1)/2}{c_{2i+1}x^{2i}} 都可以看作自变量
z=x^2 的多项式函数,即 　　
\sum_{i=0}^{(N-1)/2}{c_{2i}x^{2i}}=\sum_{i=0}^{(N-1)/2}{c_{2i}z^{i}} 　　
\sum_{i=0}^{(N-1)/2}{c_{2i}x^{2i}}=\sum_{i=0}^{(N-1)/2}{c_{2i+1}z^{i}} 　　这意味着,如果我们能算出两个子多项式
h_1(z)=\sum_{i=0}^{(N-1)/2}{c_{2i}z^{i}} ,
h_2(z)=\sum_{i=0}^{(N-1)/2}{c_{2i+1}z^{i}} 的在
z=\alpha^2 处的函数值
h_1(\alpha^2),h_2(\alpha^2),那么原多项式
h(x) 在
x=\alpha,x=-\alpha 可以这样导出: 　　
h(\alpha)=h_1(\alpha^2)+\alpha h_2(\alpha^2) 　　
h(-\alpha)=h'_1(\alpha^2)-\alpha h'_2(\alpha^2) 　　注意到
h_1(\alpha^2),h_2(\alpha^2)都是关于
\alpha 的偶函数,因此,我们通过计算出两个子多项式上的两点
(\alpha,h_1(\alpha^2)),(\alpha,h_2(\alpha^2)) ,实际上可以通过上式子导出
h(x) 上的两个点: 　　
(\alpha,h(\alpha))=(\alpha,h_1(\alpha^2)+\alpha h_2(\alpha^2) ) 　　
(-\alpha,h(-\alpha))=(-\alpha,h_1(\alpha^2)-\alpha h_2(\alpha^2) ) 　　将子多项式看作
z=\alpha^2 的函数,两个子多项式的次数都是
(N-1)/2 ,这意味着我们只需要解决两个规模更小的子问题即可以解决原问题 　　仔细观察,两个子多项式仍然是多项式,除了次数规模与原问题不一样外基本一致,这意味着我们可以按照上述方法继续递归求解两个子多项式. 　　形式化地说,整个算法的解决思路如下: 　　
\bullet 要解决
n 次多项式问题
Q(n) ,相当于先解决两个
Q(n/2) 的子问题,然后利用子问题的计算结果计算原问题. 　　这样做的时间复杂度是多少呢? 　　设求解
n-1 次多项式的时间规模为
T(n) 　　那么由上述思路,必有
T(n)=2T(\frac{n}{2})+n (这里加
n 是因为有
n 个点由两个子问题计算结果推导得出) 　　显然
T(1)=1 (
fft 应用领域_fft应用领域0 次多项式是一条水平直线,其值恰好等于唯一点的纵坐标) 　　所以
T(n)=n+2·(n/2)+4·(n/4)+8·(n/8)+···+n·(n/n) 　　由于
n=2^{log_2n} ,即
T(n)=n·(1+2\times(1/2)+4\times(1/4)+···) 后面的括号中有
log_2n 项. 　　所以
T(n)\approx nlog_2n 　　也就是说,这种做法的时间复杂度应该是
\bold {O(nlog_2n)} ，明显优于我们之前的朴素暴力解法. 　　5.新问题的产生&解决 　　貌似我们的问题得以解决了,但是新的问题产生了. 　　对于多项式
h(x) 拆分得到的
h'_1(x^2),h'_2(x^2) 不能利用奇偶性质选取正负点对,这是因为
\forall x\in R,x^2\ge0 　　这意味着,如果我们要将上述算法继续递归进行的话,我们有必要将自变量从实数范围扩展到复数范围内 　　拿待确定多项式
f(x)=c_3x^3+c_2x^2+c_1x+c_0 举例,将其转化为
f(x)=(c_2x^2+c_0)+x(c_3x^2+c_1) 　　令
y=x^2 ,也就是要确定子函数
g(y)=c_0+c_2y,h(y)=c_1+c_3y 上的点 　　继续进行拆分: 　　
g(y)=g_1(y)+yg_2(y),其中 g_1(y)=c_0,g_2(y)=c_2 　　
h(y)=h_1(y)+yh_2(y),其中 h_1(y)=c_1,h_2(y)=c_3 　　根据算法,对于
g_1(y),g_2(y),h_1(y),h_2(y) 这四个
fft 应用领域_fft应用领域0 次多项式,此时可以直接任意取横坐标求得点. 　　但如果我们取横坐标形如
y=-1 这样的点,那么就会出现一个问题: 　　
\bullet 我们当初设
y=x^2 ,当
y=-1 时,
x 在实数范围内无解！ 　　因此我们必须将数域扩展到复数并选取复数正负点对,以保证算法的正确进行,在这个例子中我们在四个
g_1(y),g_2(y),h_1(y),h_2(y) 函数上分别取横坐标
y=1,-1,i,-i 　　这样可以算得 　　
·点(\pm1,c_0),(\pm i,c_0)在g_1(y)上 　　
·点(\pm1,c_2),(\pm i,c_2)在g_2(y)上 　　
·点(\pm1,c_1),(\pm i,c_1)在h_1(y)上 　　
·点(\pm1,c_3),(\pm i,c_3)在h_2(y)上 　　根据上述方程，推导出
g(y),h(y) 上的点: 　　
·点(1,c_0+c_2)、(-1,c_0-c_2)、(i,c_0+c_2i)、(-i,c_0-c_2i)在g(y)上 　　
·点(1,c_1+c_3)、(-1,c_1-c_3)、(i,c_1+c_3i)、(-i,c_1-c_3i)在h(y)上 　　最后导出
f(x) 上的点: 　　
·当x=1或-1时,y=1,可导出点(1,c_0+c_1+c_2+c_3)、(-1,c_0-c_1+c_2-c_3)在f(x)上 　　
·当x=i或-i时,y=-1,可导出点(i,(c_0-c_2)+(c_1-c_3)i)、(-i,(c_0-c_2)+(c_3-c_1)i)在f(x)上 　　这个算法对于求解
k 次多项式
f(x) ,更为一般性的方法应该是这样的: 　　·找到最小的整数
N ,满足
2^N\gt k 　　·对于每次取点计算,自变量选取
x=w_0,w_1,w_2···w_{2^N-1} ,其中
w_j=e^{\frac{2j\pi}{2^N}i} ,即方程
x^{2^N}=1 在复平面内按辐角排列的第
j+1 个解 　　·按照递归思路,利用上述方法可以求出原函数
f(x) 上的
2^N 个不同的点 　　6.关于单位根的疑问与解答 　　看了,上面的算法描述,好像确实是这么回事,但为什么要这样设计？ 　　我们从横坐标的选取入手,不难看出,样例中换
y=x^2之后,如果不对数域进行扩展,那么解决两个子问题产生的点会有一半用不上.更为形式化的说,我们必须找到一个极小横坐标集合
S ,使得对于
\forall x \in S,x^2\in S 成立,才能保证所有子多项式的计算结果可以被重复利用,达到优化效果. 　　
x^{(2^N)}=1 在复数域内的解集合是恰好满足这个性质的. 　　这是因为
\forall j \in[0,2^N-1],w_j^2=(e^{\frac{2j\pi}{2^N}i})^2=e^{\frac{4j\pi}{2^N}i} 　　而
e^{\pi i}=-1 ,所以
e^{2 \pi i}=1 　　利用上面这一点性质可推出
w_j^2=w_{2j \% 2^N} (
\% 的意思表示求模),而
2j\%2^N \in [0,2^N-1] 　　因此满足我们需求的性质. 　　7.问题剩下的一半？ 　　那我们现在求出了
N 个不同的点,那我们该如何求出待求
N-1 次多项式
f(x) 呢? 　　经过上述算法操作,对于多项式
f(x)=\sum_{i=0}^{N-1}c_ix^i,我们找到了它上面的
N 个点
(w_0,f(w_0)),(w_1,f(w_1))···(w_{N-1},f(w_{N-1})) 　　相当于确定了一个
N 一次方程组. 　　
\begin{cases} f(w_0) =\sum_{i=0}^{N-1}c_iw_0^i \\ f(w_1) =\sum_{i=0}^{N-1}c_iw_1^i \\ ··· \\ f(w_{N-1}) =\sum_{i=0}^{N-1}c_iw_{N-1}^i \end{cases} 　　根据线性代数的知识,可以写成如下形式: 　　
\bold c \bold X=\bold f 　　其中: 　　
\bold c=[c_0,c_1,···,c_{N-1}] 　　
\bold X^T= \begin {bmatrix} 1 & w_0 &w_0^2 & ··· &w_0^{N-1} \\ 1& w_1 &w_1^2 & ··· &w_1^{N-1} \\ 1 & w_2 &w_2^2 & ··· &w_2^{N-1} \\ ··· & ··· &···&···&···\\ 1&w_{N-1} &w_{N-1}^2 & ··· &w_{N-1}^{N-1} \\ \end{bmatrix} ,显然这也是范德蒙矩阵 　　
\bold f=[f(w_0),f(w_1)···,f(w_{N-1})] 　　我们相当于要求出
\bold c ,即求出
\bold c= \bold f \bold X^{-1} 　　到这里,FFT最为巧妙的一点来了,我们只需要把所求
N 个点的纵坐标构造函数
f'(x)=\sum_{i=0}^{N-1}f(w_i)x^i ,然后以横坐标点集
S'=\{ \frac{w^{-1}_0}{N} ,\frac{w^{-1}_1}{N},···,\frac{w^{-1}_{N-1}}{N}\},求出
f'(x) 上的
N 个点即可,
c_i=f'( \frac{w^{-1}_i}{N}) 　　(这样做的原因涉及到离散数学中的对偶、线性代数中范德蒙矩阵的逆等问题,这里不再细究.) 　　这意味着我们算法思路与大体框架不变,可以仅仅通过修改参数的方式实现多项式系数的求解. 　　8.算法总结&代码实现 　　综上分析,我们来总结一下解决问题的算法: 　　
\bullet 系数
\rightarrow 点:(也叫FFT正变换) 　　·对于两个多项式
A(x)=\sum_{i=0}^{N-1}a_ix^i,B(x)=\sum_{i=0}^{M-1}b_ix^i ,找到最小的整数
K ,使得
fft 应用领域_fft应用领域M+N-2″ eeimg=”1″> 　　·利用第5点中的取点算法,分别求出
A(x),B(x) 上的
2^K 个点坐标（横坐标取值依次为
w_0,w_1,w_2···w_{2^K-1} ）,时间复杂度为
\bold {O(2^Klog_2(2^K))=O(K·2^K)} 　　
\bullet 中间过渡 　　·由
A(x),B(x) 上的
2^K 个点坐标,算出目标函数
C(x)=A(x)·B(x) 上对应的
2^K 个点坐标,即
(w_i,A(w_i)·B(w_i)) ,时间复杂度为
\bold {O(2^K)} 　　
\bullet 点
\rightarrow 系数:(也叫FFT逆变换) 　　·构造
f'(x)=\sum_{i=0}^{N-1}C(w_i)·x^i ,按照第7点中的方法逆向求出
C(x) 中的系数,时间复杂度为
\bold {O(K·2^K)} 　　综上分析,时间复杂度为
\bold {O(K·2^K)=O((M+N)log_2(M+N))} 　　代码(C++) 　　你可以自己测试一下这段代码,输入标准如下: 　　第一行是两个整数
N,M ,分别表示
A(x),B(x) 的最高次项. 　　第二行有
N+1 个实数,第
i 个数表示
A(x) 的
i-1 次项系数
a_{i-1} 　　第三行有
M+1 个实数,第
i 个数表示
B(x) 的
i-1 次项系数
b_{i-1} 　　如果你的输入无误,那么应该会得到
N+M+1 个实数输出,从左到右第
i 个数表示目标函数
C(x) 的
i-1 次项系数
c_{i-1} 　　例如: 　　输入: 输出:
fft 应用领域_fft应用领域
fft 应用领域_fft应用领域 　　表示函数
A(x)=1+2x 与
B(x)=1+2x+x^2 的乘积为
C(x)=1+4x+5x^2+2x^3 　　9. 应用领域 　　FFT算法的应用领域十分广阔. 　　·大数乘法 　　对于任何一个数字,我们都可以将其写成一个多项次形式，例如
426 =6\times10^0+2\times10^1+4\times10^2 　　因此对于两个大数乘法,我们可以先将它们对应的多项式形式写出来,然后利用FFT算法进行函数乘积,最后令目标函数中的自变量为
10 即可得到答案,这种做法比朴素竖式乘法快很多. 　　例如,要计算
125 \times903 　　可以将
125\rightarrow 5\times10^0+2\times10^1+1\times10^2\rightarrow A(x)=5+2x+x^2 　　同理
903\rightarrow B(x)=3+9x^2 　　计算出
C(x)=15+6x+ 48x^2+ 18x^3+ 9x^4 ,
C(10) 即是乘积 　　目前这种做法被运用于大数相乘,如Java BigInt类 与Python 自带的高精度. 　　·波形分解 　　无论是地震波还是声波,它们在一定时间内可以看作是周期函数,利用FFT算法,我们可以从总分解得到这些波形的构成,或者更为专业地说,一组正弦函数正交基.以便我们用于对波形进行分析、预测、处理. 　　详情可见：【官方双语】形象展示傅里叶变换_哔哩哔哩_bilibili 　　·图像处理 　　利用FFT算法处理图像,可以获知图中每个点在灰度空间地变化频率，图像上某一点与邻域点灰度值差异的强弱，即梯度的大小，也即该点的频率的大小（差异/梯度越大，频率越高，能量越低，在频谱图上就越 暗。差异/梯度越小，频率越低，能量越高，在频谱图上就越 亮。换句话说，频率谱上越亮能量越高，频率越低，图像差异越小/平缓）。一般来讲，梯度大则该点的亮度强，否则该点亮度弱。 　　多用于训练集的提前处理.

2024最新激活全家桶教程，稳定运行到2099年，请移步至置顶文章：https://sigusoft.com/99576.html

版权声明：本文内容由互联网用户自发贡献，该文观点仅代表作者本人。本站仅提供信息存储空间服务，不拥有所有权，不承担相关法律责任。如发现本站有涉嫌侵权/违法违规的内容， 请联系我们举报，一经查实，本站将立刻删除。 文章由激活谷谷主-小谷整理,转载请注明出处：https://sigusoft.com/89942.html

赞 (0)
0
上一篇 2024年 6月 19日
下一篇 2024年 6月 19日

相关推荐

关注微信