快速傅里叶变换的相关概念

快速傅里叶变换： 快速傅里叶变换算法是用来快速求两个多项式乘积的算法，又称为 $FFT$。设 $A(x)=a_0+a_1x+a_2x^2+…+a_nx^n,B(x)=b_0+b_1x+b_2x^2+…+b_mx^m$，则快速傅里叶变换算法就能快速求出 $A(x)B(x)$ 的结果。决定乘积的关键在于乘积的所有系数是多少，正常需要用 $O(n^2)$ 的时间才能求出来，而快速傅里叶变换算法则可以在 $O(nlogn)$ 的时间内求出乘积的所有系数。

快速傅里叶变换的原理：

在进行快速傅里叶变换之前，我们需要知道一些多项式的前置知识。

设 $A(x) = a_0+a_1x+…+a_nx^n$。

多项式的一个性质就是，我们可以用任意 $n+1$ 个不同点来唯一确定一个 $n$ 次多项式。设这 $n+1$ 个点分别是 $(x_1,y_1) \sim (x_{n+1},y_{n+1})$，则代入就能得到一个方程组：

$$ \begin{cases} a_0+a_1x_1+…+a_nx_1^n=y_1 \newline a_0+a_1x_2+…+a_nx_2^n=y_2 \newline \quad \quad \quad \quad \vdots \newline a_0+a_1x_{n+1}+…+a_nx_{n+1}^n=y_{n+1} \end{cases} $$

这个方程组要想有唯一解，就意味着它的系数矩阵的行列式不为 $0$。

$$ \begin{cases} 1+x_1+…+x_1^n \newline 1+x_2+…+x_2^n \newline \quad \quad \quad \vdots \newline 1+x_{n+1}+…+x_{n+1}^n \end{cases} $$

而这个系数矩阵是和范德蒙行列式类似的，两个行列式是转置的关系，而行列式转置后的值是相同的，因此这个系数矩阵的值就等于范德蒙行列式的值，范德蒙行列式的值就是 $\prod_{1 \leq i<j \leq n+1} (x_i-x_j)$，由于 $x_1 \sim x_{n+1}$ 是不相同的，因此这个公式中的每一项 $x_i-x_j$ 都是不为 $0$ 的，所以系数矩阵的行列式就一定不是 $0$，这就意味着这个方程组是由唯一解的，因此任意 $n+1$ 个不同的点都可以唯一确定一个 $n$ 次多项式。

这里我们在求乘积的过程中用的不是系数表示法来求乘积，而是用点表示法求乘积，由于 $A(x)B(x)$ 一定是 $n+m$ 次，因此我们需要 $n+m+1$ 个点，如果这些点在 $A(x)$ 中分别是 $(x_1,A(x_1)),(x_2,A(x_2)),…,(x_{n+m+1},A(x_{n+m+1}))$，在 $B(x)$ 中分别是 $(x_1,B(x_1)),(x_2,B(x_2)),…,(x_{n+m+1},B(x_{n+m+1}))$，那么对应的在 $A(x)B(x)$ 中的所有点就是 $(x_1,A(x_1)B(x_1)),(x_2,A(x_2)B(x_2)),…,(x_{n+m+1},A(x_{n+m+1})B(x_{n+m+1}))$，可以发现这样我们只需要用 $O(n+m)$ 的级别就能求出多项式的乘积，但是如果用系数表示法，我们需要到对应的系数矩阵中去求，那将是 $O(n^2)$ 级别的，可以发现这两种方式的在时间上的差距非常明显。

而我们经常会得到一个系数表示法的多项式，那么我们就需要快速的将一个系数表示法的多项式转化成一个点表示法的多项式，同时还能快速转化回去。这样我们在得到两个系数表示法的多项式 $A(x),B(x)$ 之后，我们可以将他们分别转化成点表示法的形式，然后就能用 $O(n+m)$ 的时间复杂度求出乘积的点表示法，然后再将乘积的点表示法转化回系数表示法，这样我们就能得到两个多项式的乘积的系数表示法。这里将系数表示法和点表示法进行快速转换的方法就是快速傅里叶变换，整个算法的时间复杂度瓶颈就是快速傅里叶变换，是 $O(nlogn)$ 的时间，那么如何进行快速傅里叶变换呢，这里为了方便我们特定项数从 $0$ 开始，即 $A(x)=a_0+a_1x+…+a_{n-1}x^{n-1}$，$B(x)$ 同理，这样我们就是要任取 $n$ 个不同点，然后唯一确定一个 $n-1$ 次多项式，这里要想快速转换，就需要在这 $n$ 个不同点上下功夫。

这里我们的 $n$ 个点取的是复数意义上的单位根，对于复数，任意的复数都可以表示成 $a+bi$ 的形式，其中 $i=\sqrt{-1}$，任意一个复数都可以转化成坐标系上的一个点 $(a,b)$，这个复数也可以表示成向量，复数的模长就是向量 $(a,b)$ 的长度 $\sqrt{a^2+b^2}$，复数和 $x$ 轴在正方向上的夹角称为复角。复数还存在一些运算性质，如两个复数 $a+bi$ 和 $c+di$ 相加，就是 $(a+bi)+(c+di)=(a+c)+(b+d)i$，这个运算和向量类似，并且复数相加的几何意义就是向量相加的几何意义。而两个复数相乘就是 $(a+bi)+(c+di)=(ac-bd)+(ad+bc)i$，而复数相乘有一个很好的几何意义，相乘后的复数的长度是原两个复数长度的乘积，相乘后的复数的复角是原两个复数的复角的和。

而对于复数意义上的单位根，其实就是一个单位圆，然后将这个单位圆等分成 $n$ 份，此时 $\omega_n^k$ 表示取其中的 $k$ 份，$\omega_n^k$ 的长度就是单位圆的长度，也就是 $1$，夹角就是一整个圆等分之后的 $k$ 份，这就是 $\omega_n^k$，也就是 $n$ 次单位根。因此 $\omega_n^0$ 就是 $0$ 号点，$\omega_n^1$ 就是 $1$ 份，假设单位圆被分成 $x$ 份，那么以此类推会到 $\omega_n^{x-1}$。为了方便我们将 $n$ 变成 $2$ 的次幂形式。

然后观察 $n$ 次单位根的一些性质，第一个性质就是对于任意的 $i\neq j$，$\omega_n^i$ 和 $\omega_n^j$ 一定不同，因为角度不同，所以位置一定。并且由于 $k$ 份的夹角为 $\frac{2k\pi}{n}$，以及单位圆的半径是 $1$，因此第二个性质就是 $\omega_n^k = cos{\frac{2k\pi}{n}}+i\sin{\frac{2k\pi}{n}}$。第三个性质就是 $\omega_n^0=\omega_n^n=1$。

第四个性质就是 $\omega_{2n}^{2k}=\omega_n^k$，这从图像上很容易可以看出。第五个性质是 $\omega_n^{k+\frac{n}{2}}=-\omega_n^k$

通过复数单位根的这五个性质，我们就能实现快速傅里叶变换，这就是我们选择用复数单位根作为 $n$ 个点的原因，同理，如果能找到 $n$ 个不同点满足这五个性质，我们同样可以用它们来实现快速傅里叶变换。

然后我们就可以推导快速傅里叶变换了，设 $A(x)=a_0+a_1x+…+a_{n-1}x^{n-1}$，现在我们要快速的将系数表示法转换成点表示法，也就是快速找到 $n$ 个点，我们这里取的点就是 $w_n^k$，其中 $k=0\sim n-1$，这样就能求出 $n$ 个不同点，从而得到点表示法。我们先将多项式中所有项按照项数的奇偶性分成两类，一类的项数全是奇数，一类的项数全是偶数，如 $A(x)=(a_0+a_2x^2+…+a_{n-2}x^{n-2})+(a_1x+a_3x^3+…+a_{n-1}x^{n-1})$，我们将前半部分用 $A_1(x)$ 来表示，即 $A_1(x)=a_0+a_2x^2+…+a_{n-2}x^{n-2}$，然后我们将 $x^2$ 换元成 $x$，那么就能得到 $A_1(x^2)=a_0+a_2x+a_4x^2…+a_{n-2}x^{\frac{n}{2}-1}$，同理将后半部分用 $A_2(x)$ 来表示，$A_2(x)=a_1x+a_3x^3+…+a_{n-1}x^{n-1}$，然后我们提出一个 $x$，再将 $x^2$ 换元成 $x$，得到 $A_2(x^2)=a_1+a_3x+a_5x^2+…+a_{n-1}x^{\frac{n}{2}-1}$，这里还有一个单独的 $x$ 在 $A_2(x^2)$ 的外面不能忘记。

此时可以发现 $A(x)=A_1(x^2)+xA_2(x^2)$，然后我们将 $\omega_n^k$ 代入，在代入之前，我们先将 $k$ 分成两段，第一段是 $0 \sim \frac{n}{2}-1$，第二段是 $\frac{n}{2} \sim n-1$，我们先将第一段的 $\omega_n^k$ 代入，得到 $A(\omega_n^k)=A_1(\omega_n^{2k})+\omega_n^kA_2(\omega_n^{2k})$，然后通过性质四，可以得到 $A(\omega_n^k)=A_1(\omega_{\frac{n}{2}}^{k})+\omega_n^kA_2(\omega_{\frac{n}{2}}^{k})$。

然后再将第二段的 $\omega_n^k$ 代入，此时我们还是将 $k$ 看作 $0 \sim \frac{n}{2}-1$，那么 $k+\frac{n}{2}$ 就在 $\frac{n}{2} \sim n-1$ 了，因此此时我们代入的是 $w_n^{k+\frac{n}{2}}$，得到 $A(\omega_n^{k+\frac{n}{2}})=A_1(\omega_n^{2k+n})+\omega_n^{k+\frac{n}{2}}A_2(\omega_n^{2k+n})$，由于 $\omega_n^k=\omega_n^{k+n}$，因此得到 $A(\omega_n^{k+\frac{n}{2}})=A_1(\omega_n^{2k})+\omega_n^{k+\frac{n}{2}}A_2(\omega_n^{2k})$，然后通过性质四、五，得到 $A(\omega_n^{k+\frac{n}{2}})=A_1(\omega_{\frac{n}{2}}^k)-\omega_n^k A_2(\omega_{\frac{n}{2}}^k)$，这样我们就得到了左、右两段的方程。

$$ \begin{cases} A(\omega_n^k)=A_1(\omega_{\frac{n}{2}}^{k})+\omega_n^kA_2(\omega_{\frac{n}{2}}^{k}) \newline A(\omega_n^{k+\frac{n}{2}})=A_1(\omega_{\frac{n}{2}}^k)-\omega_n^k A_2(\omega_{\frac{n}{2}}^k) \end{cases} $$

此时可以发现，对于原区间 $0 \sim n-1$ 中的任意一个数 $\omega_n^k$，都能通过 $0 \sim \frac{n}{2}-1$ 和 $\frac{n}{2} \sim n-1$ 这两段中的数计算得到，如果我们要求的 $\omega_n^k$ 在原区间的左半部分，那么他就是 $A_1(\omega_{\frac{n}{2}}^k)+w_n^kA_2(\omega_{\frac{n}{2}}^{k})$，如果在区间的右半部分，那么就是 $A(\omega_n^{k+\frac{n}{2}})=A_1(\omega_{\frac{n}{2}}^k)-\omega_n^k A_2(\omega_{\frac{n}{2}}^k)$，因此可以发现，我们能通过 $A_1(\omega_{\frac{n}{2}}^k)$ 和 $A_2(\omega_{\frac{n}{2}}^k)$ 推出原区间的所有项，并且 $A_1(\omega_{\frac{n}{2}}^k)$ 是原式的所有奇数项，$A_2(\omega_{\frac{n}{2}}^k)$ 是原式的所有偶数项，因此如果我们想求一下长度是 $n$ 的区间，只需要递归求一下两个长度是 $\frac{n}{2}$ 的区间。每次递归下去后总长度除以 $2$，因此我们一共会递归 $logn$ 层，每一层都是 $O(n)$ 的时间复杂度，因此整个的时间复杂度就是 $O(nlogn)$。这样我们进行完了快速傅里叶变换的正变换，也就是将系数表示法转换成点表示法，我们还需要能将点表示法转换成系数表示法。

接下来我们需要进行逆变换，就是将 $n$ 个点转换成原多项式的 $n$ 个系数。假设当前 $n$ 个点是 $(\omega_n^k,A(\omega_n^k))$，设 $x_k=\omega_n^k,y_k=A(\omega_n^k)$，则这 $n$ 个点就是 $(x_k,y_k)$。

我们现在要求的就是 $A(x)=C_0+C_1x+…+C_{n-1}x^{n-1}$ 中的所有系数 $a$，这里有一个公式 $C_k=\sum_{i=0}^{n-1}y_i(\omega_n^{-k})^i$，有了这个公式，那么我们就可以去求所有的 $C_k$，而所有的系数 $a+k$ 就是 $\frac{C_k}{n}$，这样我们就能求出所有的系数，现在我们要求出 $C_k$，我们将这个公式看成一个多项式，设 $B(x)=y_0+y_1x+…+y_{n-1}x^{n-1}$，则 $C_k=B(\omega_n^{-k})$，也就是说 $C_k$ 其实就是 $B(x)$ 这个多项式在 $n$ 个 $n$ 次单位根上的值，其实就是求 $B(x)$ 的点表示法，而这个多项式 $B(x)$ 就和我们前面求的多项式 $A(x)$ 只差了一个负号，如果我们将这个多项式按照上面的步骤重新求一遍，可以发现加了这个负号后结果还是不变的，因此要想求 $B(x)$ 我们只需要再做一遍正向的傅里叶变换即可，就能得到所有的 $C_k$，从而就能求出所有的系数 $a_k$。

以上就实现了傅里叶变换的逆变换，这里我们可以再证明一下为什么 $C_k=\sum_{i=0}^{n-1}y_i(\omega_n^{-k})^i$ 是成立的。

我们将 $\sum_{i=0}^{n-1}y_i(\omega_n^{-k})^i$ 展开，得到 $\sum_{i=0}^{n-1}(\sum_{j=0}^{n-1}a_j(\omega_n^i)^j)(\omega_n^{-k})^i$，然后由于 $(\omega_n^i)=\omega_n^{ij}=(\omega_n^j)^i$，所以得到 $\sum_{i=0}^{n-1}(\sum_{j=0}^{n-1}a_j(\omega_n^j)^i(\omega_n^{-k})^i)$，然后再将 $i$ 次方的两项合在一起，得到 $\sum_{i=0}^{n-1}(\sum_{j=0}^{n-1}a_j(\omega_n^{j-k})^i)$，然后由于两个求和的范围是相同的，因此可以直接交换求和次序，得到 $\sum_{j=0}^{n-1}a_j(\sum_{i=0}^{n-1}(\omega_n^{j-k})^i)$，然后可以发现在第二个求和中 $\omega_{j-k}$ 是一个定值，因此我们可以将他看成 $x$，则此时我们将 $\sum_{i=0}^{n-1}x^i$ 看成 $S(x)$，显然 $S(x)=1+x+x^2+…+x^{n-1}$，然后这里需要分情况讨论，当 $k \neq 0$ 时，此时 $S(\omega_n^k)=\omega_n^0+\omega_n^k+\omega_n^{2k}+…+\omega_n^{(n-1)k}$，而 $\omega_n^kS(\omega_n^k)=\omega_n^k+\omega_n^{2k}+…+\omega_n^{nk}$，由于 $\omega_n^0=\omega_n^{nk}=1$，因此 $S(\omega_n^k)=\omega_n^kS(\omega_n^k)$，所以 $(1-\omega_n^k)S(\omega_n^k)=0$，而 $k \neq 0$，所以 $\omega_n^k \neq 1$，所以 $S(\omega_n^k)=0$。而当 $k=0$ 时，那么 $S(\omega_n^k)=S(1)=n$。

此时我们再回到原式 $\sum_{j=0}^{n-1}a_j(\sum_{i=0}^{n-1}(\omega_n^{j-k})^i)$，根据上面的分析，如果 $j-k=0$，则 $\omega_n^{j-k})^i=n$，如果 $j=k\neq 0$，则 $\omega_n^{j-k})^i=0$，因此整个式子只有在 $j=k$ 时后面那项才不是 $0$，所以上式就等于 $na^k$，由此得出 $C_k=na_k$，所以原多项式的系数 $a_k=\frac{C_k}{n}$，证毕。

以上就是快速傅里叶变换的大致思路，而在实现的时候，如果用递归来实现会发现常数非常大，因此一般都是用迭代的方式来实现快速傅里叶变换。如果要用迭代来实现，那么我们还需要进一步分析。

我们考虑快速傅里叶变换的过程，我们是从原序列每次将区间分成两段，对于每一段递归下去，通过两个子区间来得到当前区间的值，而每次分段都是要将所有偶数项分成一段，所有奇数项分成一段，我们现在不用递归来实现，要用迭代来实现，那么我们就需要知道最底下一层的序列顺序，这样我们才能用最底层的值推出上一层的值，才能一步一步的往上推直到求出原序列的值。所以我们需要将最底层的序列顺序求出来。

举个例子，假设现在有一个长度为 $8$ 的序列 $a_0 \sim a_7$，每一层都按照奇数项和偶数项分成两段，最终会得到的序列是 $a_0,a_4,a_2,a_6,a_1,a_5,a_3,a_7$，我们尝试找规律，最顶层的所有数的二进制表示分别是 $000,001,010,011,100,101,110,111$，而最后一层的所有数的二进制表示分别是 $000,100,010,110,001,101,011,111$，可以发现，最底层的每一项的数刚好是最顶层的每一项翻转过来的值，这个要想证明也比较简单，最开始，如果一个数的二进制最低位是 $1$，那么就是奇数，会被放到后面一半，如果最低位是 $0$，那么就是偶数，会被放到前面一半，而如果被放到后面一半，那么最终它如何移动也只会在最后一半，说明它一定是最高的一半，最高的一半说明它的最高位就应该是 $1$，同理最低位是 $0$ 的数会被放到前面一半，那么它就一定是低的一半，低的一半说明它的最高位就应该是 $0$，然后分完一次后再将每一半区间分成两段，将倒数第二位是 $1$ 的放后面，是 $0$ 的放前面，以此类推，就能发现一个数的编号如果最低位是 $1$，那么最终它的编号的最高位就是 $1$，如果最低位是 $0$，那么最终的最高位就是 $0$，如果次高位是 $1$，那么最终的次高位就是 $1$，以此类推，可以发现最终的编号刚好就是初始编号的二进制翻转。

那么我们如何求出每个数的二进制翻转后的数呢，可以用递推来求，设 $rev_i$ 表示 $i$ 的二进制翻转后的数，那么我们先将 $i$ 的最低位去掉，此时应该是 $i\gg 1$，那么 $rev_{i\gg 1}$ 应该是 $i \gg 1$ 翻转后的数，那么 $i$ 翻转后的值应该是将 $rev_{i \gg 1}$ 右移一位后再将 $i$ 的最低位补到 $rev_{i\gg 1}$ 右移一位后留出来的最高位上，也就是 $rev_i=(rev_{i\gg 1})\vert \big((i \& 1) \ll (bit-1)\big)$，其中 $bit$ 表示所有编号的二进制表示下的最高位数，这样我们就能递推出所有数的二进制翻转后的数，然后就能求出每个数在最底层是的对应编号，这样我们从最底层一直向上计算，就能得到原序列的值