高速フーリエ変換の計算量解析

高速フーリエ変換とは、離散フーリエ変換を高速に行うアルゴリズムです。 具体的にどのくらい高速なのか気になったので解析してみます。

まず、離散フーリエ変換について
素数nの入力ベクトルをx、同じく出力の周波数ベクトルをfとすると、離散フーリエ変換は次のような行列積で表されます。

w_n=e^{-\frac{2 \pi}{n}i}とし、
$$ \begin{pmatrix} f_0 \\ f_1 \\ f_2 \\ \vdots \\ f_{n-1} \end{pmatrix} = \begin{pmatrix} w_n^0 & w_n^0 & w_n^0 & \cdots & w_n^0 \\ w_n^0 & w_n^1 & w_n^2 & \cdots & w_n^{n-1} \\ w_n^0 & w_n^2 & w_n^4 & \cdots & w_n^{2(n-1)} \\ \vdots & \vdots & \vdots & \ddots & \vdots \\ w_n^0 & w_n^{n-1} & w_n^{2(n-1)} & \cdots & w_n^{(n-1)^2} \end{pmatrix} \begin{pmatrix} x_0 \\ x_1 \\ x_2 \\ \vdots \\ x_{n-1} \end{pmatrix} $$
この行列式は、n^2回の掛け算とn(n-1)回の足し算を必要とするため、計算量は明らかに\Theta(n^2)です。
後々のことを考えて、この真ん中の正方行列に名前をつけておきます $$ M_n= \begin{pmatrix} w_n^0 & w_n^0 & w_n^0 & \cdots & w_n^0 \\ w_n^0 & w_n^1 & w_n^2 & \cdots & w_n^{n-1} \\ w_n^0 & w_n^2 & w_n^4 & \cdots & w_n^{2(n-1)} \\ \vdots & \vdots & \vdots & \ddots & \vdots \\ w_n^0 & w_n^{n-1} & w_n^{2(n-1)} & \cdots & w_n^{(n-1)^2} \end{pmatrix} $$ そして、nが2の倍数のとき(仮にn=8の場合)、離散フーリエ変換の式はこのようになります。

$$ \begin{pmatrix} f_0 \\ f_1 \\ f_2 \\ f_3 \\ f_4 \\ f_5 \\ f_6 \\ f_7 \end{pmatrix} = \begin{pmatrix} w_8^0 & w_8^0 & w_8^0 & w_8^0 & w_8^0 & w_8^0 & w_8^0 & w_8^0 \\ w_8^0 & w_8^1 & w_8^2 & w_8^3 & w_8^4 & w_8^5 & w_8^6 & w_8^7 \\ w_8^0 & w_8^2 & w_8^4 & w_8^6 & w_8^8 & w_8^{10} & w_8^{12} & w_8^{14} \\ w_8^0 & w_8^3 & w_8^6 & w_8^9 & w_8^{12} & w_8^{15} & w_8^{18} & w_8^{21} \\ w_8^0 & w_8^4 & w_8^8 & w_8^{12} & w_8^{16} & w_8^{20} & w_8^{24} & w_8^{28} \\ w_8^0 & w_8^5 & w_8^{10} & w_8^{15} & w_8^{20} & w_8^{25} & w_8^{30} & w_8^{35} \\ w_8^0 & w_8^6 & w_8^{12} & w_8^{18} & w_8^{24} & w_8^{30} & w_8^{36} & w_8^{42} \\ w_8^0 & w_8^7 & w_8^{14} & w_8^{21} & w_8^{28} & w_8^{35} & w_8^{42} & w_8^{49} \\ \end{pmatrix} \begin{pmatrix} x_0 \\ x_1 \\ x_2 \\ x_3 \\ x_4 \\ x_5 \\ x_6 \\ x_7 \end{pmatrix} $$ w_8^8=e^{-8\frac{2 \pi}{8}i}=e^{-2 \pi i}=1よりこの式は次の式と等価です。 $$ \begin{pmatrix} f_0 \\ f_1 \\ f_2 \\ f_3 \\ f_4 \\ f_5 \\ f_6 \\ f_7 \end{pmatrix} = \begin{pmatrix} w_8^0 & w_8^0 & w_8^0 & w_8^0 & w_8^0 & w_8^0 & w_8^0 & w_8^0 \\ w_8^0 & w_8^1 & w_8^2 & w_8^3 & w_8^4 & w_8^5 & w_8^6 & w_8^7 \\ w_8^0 & w_8^2 & w_8^4 & w_8^6 & w_8^0 & w_8^2 & w_8^4 & w_8^6 \\ w_8^0 & w_8^3 & w_8^6 & w_8^1 & w_8^4 & w_8^7 & w_8^2 & w_8^5 \\ w_8^0 & w_8^4 & w_8^0 & w_8^4 & w_8^0 & w_8^4 & w_8^0 & w_8^4 \\ w_8^0 & w_8^5 & w_8^2 & w_8^7 & w_8^4 & w_8^1 & w_8^6 & w_8^3 \\ w_8^0 & w_8^6 & w_8^4 & w_8^2 & w_8^0 & w_8^6 & w_8^4 & w_8^2 \\ w_8^0 & w_8^7 & w_8^6 & w_8^5 & w_8^4 & w_8^3 & w_8^2 & w_8^1 \\ \end{pmatrix} \begin{pmatrix} x_0 \\ x_1 \\ x_2 \\ x_3 \\ x_4 \\ x_5 \\ x_6 \\ x_7 \end{pmatrix} $$ そしてこれをこのように変形します(添え字の順番が入れ替わってることに注意) $$ \begin{pmatrix} f_0 \\ f_2 \\ f_4 \\ f_6 \\ f_1 \\ f_3 \\ f_5 \\ f_7 \end{pmatrix} = \begin{pmatrix} w_8^0 & w_8^0 & w_8^0 & w_8^0 & w_8^0 & w_8^0 & w_8^0 & w_8^0 \\ w_8^0 & w_8^2 & w_8^4 & w_8^6 & w_8^1 & w_8^3 & w_8^5 & w_8^7 \\ w_8^0 & w_8^4 & w_8^0 & w_8^4 & w_8^2 & w_8^6 & w_8^2 & w_8^6 \\ w_8^0 & w_8^6 & w_8^4 & w_8^2 & w_8^3 & w_8^1 & w_8^7 & w_8^5 \\ w_8^0 & w_8^0 & w_8^0 & w_8^0 & w_8^4 & w_8^4 & w_8^4 & w_8^4 \\ w_8^0 & w_8^2 & w_8^4 & w_8^6 & w_8^5 & w_8^7 & w_8^1 & w_8^3 \\ w_8^0 & w_8^4 & w_8^0 & w_8^4 & w_8^6 & w_8^2 & w_8^6 & w_8^2 \\ w_8^0 & w_8^6 & w_8^4 & w_8^2 & w_8^7 & w_8^5 & w_8^3 & w_8^1 \\ \end{pmatrix} \begin{pmatrix} x_0 \\ x_2 \\ x_4 \\ x_6 \\ x_1 \\ x_3 \\ x_5 \\ x_7 \end{pmatrix} $$ ここで、 $$ A= \begin{pmatrix} x_0 \\ x_2 \\ x_4 \\ x_6 \end{pmatrix} , B= \begin{pmatrix} x_1 \\ x_3 \\ x_5 \\ x_7 \end{pmatrix} , D_n = \begin{pmatrix} w_{2n}^0 & 0 & 0 & \cdots & 0 \\ 0 & w_{2n}^1 & 0 & \cdots & 0 \\ 0 & 0 & w_{2n}^2 & \cdots & 0 \\ \vdots & \vdots & \vdots & \ddots & \vdots \\ 0 & 0 & 0 & \cdots & w_{2n}^{n-1} \end{pmatrix} $$ とおくと、w_8^4=-1なので $$ \begin{pmatrix} f_0 \\ f_2 \\ f_4 \\ f_6 \\ f_1 \\ f_3 \\ f_5 \\ f_7 \end{pmatrix} = \begin{pmatrix} M_4 A + D_4 M_4 B \\ M_4 A - D_4 M_4 B \end{pmatrix} $$ これを一般に拡張すると、入力ベクトルをX_n、その半分のベクトルをA_n,B_n (\frac{n}{2}次元のベクトル)、出力ベクトルをF_nとして $$ F_n=M_n X_n = \begin{pmatrix} M_{\frac{n}{2}} A_n + D_{\frac{n}{2}} M_{\frac{n}{2}} B_n \\ M_{\frac{n}{2}} A_n - D_{\frac{n}{2}} M_{\frac{n}{2}} B_n \end{pmatrix} $$ ここで、M_{\frac{n}{2}} A_nM_{\frac{n}{2}} B_nはそれぞれ要素数\frac{n}{2}の離散フーリエ変換であり、1度計算すれば2度使いまわすことができます。
よって、要素数nフーリエ変換の計算量をf(n)で表すとき、 $$ f(n)=f(\frac{n}{2})+f(\frac{n}{2})+\frac{n}{2}+n $$ 回の計算が必要なので、 $$ f(n)=2 f(\frac{n}{2})+\frac{3}{2}n $$ 初期値f(1)=1 (掛け算1回で変換終わり)を考えると $$ f(n)=\frac{3}{2}n \log_2 n+n $$ となります。
よって、高速フーリエ変換の計算量は\Theta(n \log n)です。