快速傅里叶变换FFT

快速傅里叶变换FFT

作为我5.1出来培训这几天的原因,学会了之后当然要记录一下awa

虽然这是我唯一学会的内容qwq

首先,FFT是干什么的:

如果有这么两个东西:

A(x)=i=0n1aixi

B(x)=i=0n1bixi

让你求出 A(x)B(x)

啊,我会我会!直接 n2 枚举,进行多项式乘法就好了awa

是的,但是如果 n<=1e6 呢?

你会发现 n2 就解决不了了,而且,枚举的方法显然已经无法继续优化了,所以我们需要另辟蹊径。

首先我们设 C(x)=A(x)B(x)

然后如果我们能快速的求出 A(x)B(x) 的点值,那么我们就可以把他们乘在一起,之后就可以得到 C(x) 的点值,然后就可以逆向唯一确定 C(x) 的表达式了

(原理:n 次多项式唯一的由 n + 1 个不同点的值确定。)

(求点值就是随便带入一个 x ,然后求出它在函数图像上的 y

所以FFT就是干这个的:

(分为以下三步)

1.快速求出 A(x)B(x) 的点值。(DFT)

2.将求出的 A(x)B(x) 相乘

3.通过得到的 C(x)n 个点值求出 C(x)n 个系数 (IDFT)

(可以发现第三步求出的数就是答案)

开始求解:

我们可以记单位根 ωn=e2πin

i 就是虚数,一般为 1

(选择 e2πin 的原因是 e2πi=1 ,一会方便计算awa)

那么我们就可以将选的 n 个点值设定为:

x=ωn0,ωn1,ωn2ωnn1

之后我们将 A(x)=i=0n1aixi 拆成 A1(x2)+xA2(x2)

其中:

A1(x)=a0+a2x+a4x2++an2xn21

A2(x)=a1+a3x+a5x2++an1xn21

那我们再把我们设的 x=ωnk 带入:

A(ωnk)=A1((ωnk)2)+ωnkA2((ωnk)2)

=A1((e2πi2n)k)+ωnkA2((e2πi2n)k)

=A1((e2πin2)k)+ωnkA2((e2πin2)k)

=A1(ωn2k)+ωnkA2(ωn2k)

因为我们设的是x=ωn0,ωn1,ωn2ωnn1,所以当 ωin2 阶的时候, k 一定是小于 n2 的,所以这只是前一半。

然后对于另一半,我们可以带入带入 x=ωnk+n2 ,可以得出:

A(ωnk+n2)=A1((ωnk)2ωnn)+ωnk+n2A2((ωnk)2ωnn)

=A1((e2πi2n)k1)+ωnn2ωnkA2((e2πi2n)k1)

=A1((e2πin2)k)ωnkA2((e2πin2)k)

=A1(ωn2k)ωnkA2(ωn2k)

也就是说,如果我们已经求出了 A1,A2 在阶为 n2 时,某个单位根的点值,那么我们可以通过 O(n) 计算出 A(x) 在阶为 n 时,相应单位根的点值。

然后一直递归下去,当 n=1 时,ω1=e2πi1=e2πi=1 。答案就是 A(x) 的系数,然后我们再一层一层的返回,就可以求出最终答案。

这样就会发现,复杂度就是 T(n)=2T(n/2)+O(n)=O(nlogn)awa

然后我们发现当上式可做时,最好需要 n=2k , 那我们就不妨让 n=2k ,然后多出来的位数就让它的系数为 0 就好了。

所以我们就通过分治的方法,在 O(logn) 的时间内,解决了第一步(DFT)awa

之后我们就分别对 A(x)B(x) 求一遍,再将它们对应的点值相乘,就求出了 C(x) 的点值( 复杂度: O(n)

第二步解决awa

然后是第三步。

乍一看第三步十分的困难 qwq ,但是我们从另一个角度去看第一步(DFT)

可以发现:我们实际上是把一个单位根矩阵 M=[ωij] 左乘在了系数向量上, 那么实现 IDFT 实际上我们需要找到这个单位根矩阵的逆。

(上面一句一听就不是人话,让我们来形象的翻译一下)

首先,我们可以发现, DFT 其实干了下面这个事情:

[(ωn0)0(ωn0)1(ωn0)n2(ωn0)n1(ωn1)0(ωn1)1(ωn1)n2(ωn1)n1(ωnn2)0(ωnn2)1(ωnn2)n2(ωnn2)n1(ωnn1)0(ωnn1)1(ωnn1)n2(ωnn1)n1]×[a0a1an2an1]=[A(ωn0)A(ωn1)A(ωnn2)A(ωnn1)]

对应上面的:“我们实际上是把一个单位根矩阵 M=[ωij] 左乘在了系数向量上”

只不过我们用了一些比较魔幻的方法将它加速到了 O(nlogn) awa

为了方便表示,我们设:

Mω=[(ωn0)0(ωn0)1(ωn0)n2(ωn0)n1(ωn1)0(ωn1)1(ωn1)n2(ωn1)n1(ωnn2)0(ωnn2)1(ωnn2)n2(ωnn2)n1(ωnn1)0(ωnn1)1(ωnn1)n2(ωnn1)n1]

Ma=[a0a1an2an1]

MA=[A(ωn0)A(ωn1)A(ωnn2)A(ωnn1)]

相应的,McMC 同理

然后现在的情况是我们通过某些奇妙的方法求出了 MC ,并且我们已经知道了 Mω 然后现在我们就是要求出 Mc

怎么求呢 qwq?

可以发现其实就是在 MC 上乘上 Mω 的逆矩阵(不知道为什么的可以先去学一下矩阵相关的芝士)

我们先设 MωMω 的逆矩阵。

那么 Mω 是什么呢?

我们先来浅猜一手:

Mω=[ωij]

(就是把 Mω 的每一项的次数都 ×(1)

然后我们发现 Mω×Mω=n×I0

I0 就是单位矩阵,即除了对角线是 1 ,其他都是 0 该不会有人不知道吧awa

所以我们就找到了 Mω 的逆矩阵 Mω=[ωijn]

第三步解……

Q:哎!不对!!还没解决呢!!!

Q:你这样进行矩阵乘法,不就是 O(n2) 的了吗?

Q:你这复杂的有问题啊。

不不不,我们可以发现,其实 IDFT 就是把 DFT 里面的 Mω 换成了 Mω , Ma 换成了 MC

其他的没有任何的改变。

所以可以如法炮制,直接向上面一样 O(nlogn) 解决 awa

所以第三步(IDFT)解决.

所以,我宣布,FFT,搞定 awa

posted @   YT0104  阅读(10)  评论(0编辑  收藏  举报
相关博文:
阅读排行:
· 一个适用于 .NET 的开源整洁架构项目模板
· API 风格选对了,文档写好了,项目就成功了一半!
· 【开源】C#上位机必备高效数据转换助手
· .NET 9.0 使用 Vulkan API 编写跨平台图形应用
· MyBatis中的 10 个宝藏技巧!
点击右上角即可分享
微信分享提示