RainAir
My OI Blog
RainAir
FFT 学习笔记

快速傅里叶变换 (Fast Fourier Transform, FFT) 是一种 $ O(n\ log_2\ n) $ 的时间复杂度内完成离散傅里叶变换 (DFT) 的算法,OI 中通常用来优化多项式乘法。

常见题目类型链接
给定一个 $n$ 次多项式 $F(x)$ 和一个 $m$ 次多项式 $G(x) $。求出它们的卷积。

FFT 就可以在 $ n\ log_2\ n $ 的时间复杂度内解决该类问题。

先来学习一些 FFT 前备知识。

Q:我是初中的啊…看不懂怎么办?

A:本人尽量写到初中生理解水平…因为我就是一个没学过数学奥赛的初中生…大可放心食用

多项式

系数表示法

设 $F(x)$ 表示一个 $n-1$ 次的多项式,则 $$ F(x) = \Sigma_{i=0}^{n-1} a_i * x_i $$

利用系数直接运算多项式乘法的时间复杂度是 $ O(n^2) $ 。(即每个系数都要枚举相乘)

点值表示法

将这个多项式 $F(x) $ 看成一个函数,就是说代入 $ n $ 个互不相同的 $x$ ,会得到 $n$ 个不同的取值 $y$。

那么这个多项式就被这 $n$ 个点 $(x_1,y_1),(x_2,y_2) \cdots (x_n,y_n) $ 唯一确定。

当然这样计算也是 $O(n^2)$ 的。

当然这些都需要优化。对于第一种好像不会优化,于是我们考虑第二种表示法的优化。


我们为了照顾初中选手,先介绍一些简单的东西。

向量(矢量)

有方向也有大小的量,与标量(无方向)相对。

在集合中我们用有箭头的直线来表示。

圆的弧度制

我们之前学的角度制的定义是将圆 $ 360 $ 等分。

弧度制的定义是等于半径长的圆弧对应的圆心角叫做 1 弧度的角。用符号 rad 表示,读作弧度,用 rad 作为单位来度量角的制度叫做弧度制(字面意思)

显然有:$ 1 ^{ ∘ } = \frac{\pi}{180}rad $

因为设这个圆的半径为 r ,我们发现周长是 $2\pi r$ ,那么在这个周长上取长度为 r 的弧,圆心角显然是 $ ( \frac{r}{2\pi r} * 360^∘ )rad= \frac{\pi}{180}rad $

平行四边形定则

有一个平行四边形 $ABCD $ ,那么有 $ AB + AD = AC$。

可以表示成向量加法。

等比数列求和

设一个等比数列中:首项为 $ a_1 $,公比是 $q$。

那么前 $n$ 项的和 $S = \frac{a_1(1-q^n)}{1-q}$

证明可以通过错位相减来实现。

复数

定义

设实数 $a,b$,定义虚数单位 $ i $,有 $i^2 = -1 $ ,将形如 $a + bi$ 的数叫做复数。复数域是目前已知最大的域。

在复平面里(x 轴表示实数,y 轴表示虚数)。从 $ O(0,0) $ 到 $A(a,b)$ 的向量能表示复数 $a + bi$。

模长:$OA$ 的长度,即 $\sqrt{a^2 + b^2}$

幅角:假设以逆时针为正方向,从 $x$ 正半轴旋转至向量的转角的有向角叫做幅角。

运算法则

加法

向量相加,满足平行四边形定则。即:$ (a + bi) + (c + di) = (a + c) + (b+d)i $

复数的加法是封闭的。

乘法

向量的乘法几何定义:复数相乘,模长相乘,幅角相加。

代数定义:

$$(a + bi) * (c + di) =ac + adi + bci + bdi^2 = ac + adi + bci – bd = (ac-bd) + (bc + ad)i$$

单位根

下文若不做特别的声明,均默认 $n$ 为 $2$ 的整数次幂。

在复平面上以原点为圆心,1 为半径作出的圆叫做单位圆。以圆点为起点,圆的 n 等分点为终点作出 n 个向量。设幅角为正且最小的向量对应的复数为 $\omega_n$ ,称为 $n$ 次单位根。

其余的 $n-1$ 个复数显然是 $ \omega_n^2,\omega_n^3 \cdots \omega_n^n$

注意到有 $\omega_n^0 = \omega_n^n = 1$ (对应方向为 x 轴正半轴的向量)

单位根的幅角是周角的 $\frac{1}{n}$ (显然 $n$ 等分)。

一些性质:

性质1

$$ \omega_n^k = cos\ k * \frac{2\pi}{n} + i\ sin\ k * \frac{2\pi}{n} \Rightarrow \omega_{2n}^{2k} = \omega_n^k$$

证明:

$$\omega_{2n}^{2k} = \cos 2k * \frac{2 \pi }{2n} + i \sin 2k* \frac{2 \pi }{2n} = \cos k* \frac{2 \pi }{n} + i \sin k * \frac{2\pi}{n} =\omega_n^k$$

性质2

$ \omega_n^{k+\frac{n}{2}} = -\omega_n^k $

证明:

$$ w_n^{ \frac{n}{2} } = \cos \frac{n}{2} * \frac{2 \pi }{n} + i \sin \frac{n}{2} * \frac{2 \pi }{n} = \cos \pi + i \sin \pi = -1$$

我们终于把预备知识讲完了,那么接下来开始正题。

快速傅里叶变换

朴素的求一个多项式的点值表达法的复杂度是 $ O(n^2) $,我们先来探究一些性质。

设多项式 $F(x)$ 的系数为 $$ {a_0,a_{1},a_{2},\cdots,a_{n-1}}$$

那么

$ F(x) = \Sigma_{i=0}^{n-1} a_i * x^i $

将其按照下标分成两个多项式:

$$ F_1 (x) = \Sigma_{i=0}^{ \frac{n}{2} – 1 } a_{i * 2} x^{i * 2} $$
$$F_2(x) = \Sigma_{i=0}^{ \frac{n}{2} – 1 }a_{i * 2 + 1 } x^{i * 2 + 1 }$$

那么可以得到:

$$ F(x) = F_1(x^2) + xF_2(x^2) $$

我们代入 $\omega_n^k(k < \frac{n}{2})$ 得到

$$ F( \omega^k_n ) = F_1( \omega^{2k}_ n ) + \omega^k_n F_2( \omega_n^{2k} ) $$

$$ = F_1(\omega^k_{ \frac{n}{2} }) + \omega^k_nF_2( \omega_{ \frac{n}{2}}^{k})$$

同理,代入 $w_n^{k+\frac{n}{2}}$得到

$$ F( \omega_n^{k + \frac{n}{2} } ) = F_1({ \omega^{2k+n}_ n }) + \omega_n^{k+ \frac{n}{2} }F_2(\omega_n^{2k+n}) $$
$$ = F_1(\omega_n^{2k}) – \omega_n^kF_2(\omega_n^{2k} * \omega_n^n) $$
$$= F_1(\omega_n^{2k}) – \omega_n^kF_2(\omega_n^{2k}) $$
$$ = F_1(\omega^k_{\frac{n}{2}}) – \omega^k_nF_2(\omega_{\frac{n}{2}}^{k})$$

发现这两个式子仅由一个常数项不同(-号)

所以我们只需要算出第一个式子,第二个式子就可以 $O(1)$ 求。

所以说我们可以递归二分去计算这个东西,遇到常数项就返回。

时间复杂度:$ T(n) = 2T(n/2) + O(n) = O(n\ log_2\ n)$

快速傅里叶逆变换

还没有结束…

我们要考虑怎么输出答案啊。

一般来说没有人会让你输出点值表示法的,所以我们需要用逆变换来换回来。

首先我们可以将 FFT 后的结果看做一个向量 $$ (y_0,y_1,y_2,\cdots,y_{n-1}) $$。

我们设另一个向量 $$ (c_0,c_1,c_2,\cdots,c_{n-1}) $$ 令其满足

$$ c_k=\Sigma_{i=0}^{n-1}y_i(\omega^{-k}_n)^i $$

浅显易懂的讲就是将结果看做一个多项式,求出这个多项式的 FFT 的结果。

推一波式子:

$$ c_k=\Sigma_{i=0}^{n-1}y_i(\omega^{-k}_ n)^i $$
$$=\Sigma_{i=0}^{n-1}(\Sigma_{j=0}^{n-1}a_j(\omega^i_n)^j)(\omega_n^{-k})^i $$
$$=\Sigma_{i=0}^{n-1}(\Sigma_{j=0}^{n-1}a_j(\omega_n^j)^i)(\omega_n^{-k})^i $$
$$= \Sigma_{i=0}^{n-1}(\Sigma_{j=0}^{n-1}a_j(\omega_n^j)^i(\omega_n^{-k})^i $$
$$= \Sigma_{i=0}^{n-1}\Sigma_{j=0}^{n-1}a_j(\omega_n^j)^i(\omega_n^{-k})^i $$
$$=\Sigma_{i=0}^{n-1}\Sigma_{j=0}^{n-1}a_j(\omega_n^{j-k})^i $$
$$ = \Sigma_{j=0}^{n-1}a_j \Sigma_{i=0}^{n-1}(\omega_n^{j-k})^i $$

设 $ S(x) = \Sigma_{i=0}^{n-1} x^i$

代入 $ \omega_n^k $ 得:

$$ S(\omega_n^k) = \Sigma_{i=0}^{n-1}w_n^i $$

直接套用等比数列求和(错位相减法),得:

$ S(\omega_n^k) = \frac{1-1}{\omega_n^k-1} $

不难看出分子为 0 ,分母不为 0 。

分类讨论一下:

$ K \neq 0$ :$ S(\omega_n^k) = 0 $

$ K = 0$:$S(\omega_n^0) = n$

我们继续考虑式子 $$c_k = \Sigma_{j=0}^{n-1}a_j\Sigma_{i=0}^{n-1}(\omega_n^{j-k})^i $$

同理:

$ j \neq k \to (j-k) \neq 0 $:$c_k = 0$

$ j = k \to (j – k) = 0$:$c_k = n $

这样我们就得到点值与系数之间的表示方法:

$$ c_k = na_k \Longleftrightarrow a_k = \frac{c_k}{n} $$

递归实现的代码

蝴蝶操作

稍微优化一下小常数。

发现 $w*a2[i]$ 被计算了两次,开一个 $t$ 来存储仅需计算一次。

迭代实现

递归的效率过于低下,我们需要更快的解决这个问题。

我们观察原序列的下标和分类后的序列下标。

发现二进制翻转可以就直接解决这个问题。

这样我们就可以 $O(n)$ 的利用某种操作来求得序列,然后依照规律合并即可。

赞赏
知识共享许可协议
本文链接: https://blog.aor.sd.cn/archives/74
如文中无特殊声明,本文采用 CC BY-NC-SA 4.0 进行许可,转载请说明出处!
希望 CSP 不要翻车,希望省选不要翻车
https://secure.gravatar.com/avatar/97c17c68a1e55e11bb5558bc0f10cc0d?s=256&d=mm&r=g

RainAir

文章作者

一个OIer。

推荐文章

发表评论

textsms
account_circle
email

RainAir

FFT 学习笔记
快速傅里叶变换 (Fast Fourier Transform, FFT) 是一种 $ O(n\ log_2\ n) $ 的时间复杂度内完成离散傅里叶变换 (DFT) 的算法,OI 中通常用来优化多项式乘法。 常见题目类型链接: 给定…
扫描二维码继续阅读
2018-09-05
标签
近期评论