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)
上一篇 2024年 6月 19日
下一篇 2024年 6月 19日

相关推荐

关注微信