多次元化FFT
DFT で計算対象になる数値は一般に長さ $n$ の 1 次元配列でもって表される。 FFT の演算ではほぼランダムアクセスに近いメモリアクセスをするため、 データサイズがキャッシュメモリより大きくなると極端に演算速度が遅くなる。
これを回避するため、$n=pq$ と因数分解できると仮定し、 DFT の計算式を次のように書き換える。
\[ k = k_1p + k_2,\ j = j_1q + j_2\quad (0\leq k_1, j_2 \lt q,\ 0\leq k_2, j_1 \lt p) \] \[ \begin{eqnarray} A_{k_1p+k_2} &=& \sum_{j_2=0}^{q-1}\sum_{j_1=0}^{p-1} a_{j_1q+j_2}\omega^{(j_1q+j_2)(k_1p+k_2)}\\ &=& \sum_{j_2=0}^{q-1}\sum_{j_1=0}^{p-1} a_{j_1q+j_2}\omega^{j_1k_2q+j_2k_1p+j_2k_2}\\ &=& \sum_{j_2=0}^{q-1} \left( \left(\sum_{j_1=0}^{p-1} a_{j_1q+j_2}\omega_p^{j_1k_2} \right) \omega^{j_2k_2} \right) \omega_q^{j_2k_1}\\ \end{eqnarray} \]すると内側の括弧内の計算は長さ $p$ の DFT が $q$ 組、 外側の括弧の外は長さ $q$ の DFT が $p$ 組、という形になるので $p$ や $q$ 次第ではキャッシュメモリに収まり演算全体の高速化が期待できる。 また、それぞれの DFT は独立に計算することができるので共有メモリ型の並列計算 も容易にできる。 一般的には $\max (p,q)$ を最小化するため $p\simeq q \simeq \sqrt{n}$ にする。
2 次元化 （Six-step FFT）
上記の式だけでは若干わかりにくいので、入力となる 1 次元データ列を 2 次元に配置し見てみよう。 大きく 6 つの計算ステップを踏むことから Six-step FFT という呼ばれ方で広く知られている。 昔は最初と最後の転置作業を無視して Four-step FFT という呼び方もしていたが、現在そう呼ばれることはあまりない。
具体的な 6 つのステップを以下に示す。中で $[a,b]$ とあるのは $b$ に関して連続に配置される 2 次元配列のインデックスである。 また、$0\leq j_1, k_2 \lt p$、$0\leq j_2, k_1 \lt q$ であり 各配列の長さは該当インデックスの範囲と同じであるとする。
| $a_0[j_1,j_2] = a_{j_1q+j_2}$ | (配列の 2 次元再解釈) | |
| Step 1. | $a_1[j_2,j_1] = a_0[j_1,j_2]$ | (行列の転置) |
|---|---|---|
| Step 2. | $a_2[j_2,k_2] = {\displaystyle \sum_{j_1=0}^{p-1}} a_1[j_2,j_1]\omega_p^{j_1k_2}$ | (長さ $p$ の FFT を $q$ 組) |
| Step 3. | $a_3[j_2, k_2] = a_1[j_2, k_2] \omega^{j_2k_2}$ | (ひねり係数の乗算) |
| Step 4. | $a_4[k_2, j_2] = a_3[j_2, k_2]$ | (行列の転置) |
| Step 5. | $a_5[k_2,k_1] = {\displaystyle \sum_{j_2=0}^{q-1}}a_4[k_2,j_2] \omega_q^{j_2k_1}$ | (長さ $q$ の FFT を $p$ 組) |
| Step 6. | $a_6[k_1,k_2] = a_5[k_2,k_1]$ | (行列の転置) |
| $A_{k_1p+k_2} = a_6[k_1,k_2]$ | (配列の 1 次元再解釈) |
3 回行われる行列の転置作業は素直に考えれば元データと同容量の作業領域を要求することになるが 例えば Step 3 の乗算と Step 4 の転置を同時に行えばキャッシュアウトする回数を減らすことができる。 また、Step 1 から Step 3-4 までについて全データを対象にするのではなく、 $p$ の数倍だけ転置させることで (Step 2 を $q$ 組同時に行うのでなく $t$ 組ずつ行うことになる) オンキャッシュのデータで計算させる手法 (Block 分割 FFT) もある。
3 次元化 FFT (Nine-step FFT)
同様に 3 次元化も可能である。 データ量が非常に大きい場合、2 次元化しても分割後の 1 次元 FFT がキャッシュサイズに載らないこともあり、 3 次元化することで効果が上がるかも知れない。 また分散メモリ型並列計算にも応用させやすい。 式変形は面倒であるが 1 度やっておくと理解が深まる。 変換前後のデータの 3 軸対応をどう取るかについては、 変換後の式を考えて $\omega$ の項が少なくなるように取ると良い。 1 つの分割、計算例を以下に書く。
\[ n = pqr, k = (k_1q + k_2)p + k_3, j = (j_1q + j_2)r + j_3 \] \[ (0\leq k_3,j_1\lt p,\ 0\leq k_2,j_2 \lt q,\ 0\leq k_1,j_3\lt r) \] \[ \begin{eqnarray} A_{k_1pq+k_2p+k_3} &=& \sum_{j_1=0}^{p-1}\sum_{j_2=0}^{q-1}\sum_{j_3=0}^{r-1} a_{j_1qr+j_2r+j_3} \omega^{(j_1qr+j_2r+j_3)(k_1pq+k_2p+k_3)}\\ &=& \sum_{j_1=0}^{p-1}\sum_{j_2=0}^{q-1}\sum_{j_3=0}^{r-1} a_{j_1qr+j_2r+j_3} \omega^{j_1k_3qr+j_2k_2pr+j_3k_1pq+j_2k_3r+j_3k_2p+j_3k_3}\\ &=& \sum_{j_3=0}^{r-1}\left(\left(\sum_{j_2=0}^{q-1}\left(\left(\sum_{j_1=0}^{p-1} a_{j_1qr+j_2r+j_3} \omega_p^{j_1k_3} \right) \omega^{(j_2r+j_3)k_3}\right) \omega_q^{j_2k_2}\right)\omega_{qr}^{j_3k_2}\right)\omega_r^{j_3k_1} \end{eqnarray} \]行列転置を含めて再解釈
3 次元化 FFT では最初と最後を含めて 4 回の行列転置が必要になる。 また、分散メモリ型計算機を用いてこの計算を行う場合はそのうち 何回かはデータ通信を用いての転置となるため、通信速度がネックになりやすい。
| $a_0[j_1,j_2,j_3]=a_{j_1qr+j_2r+j_3}$ | (配列の 3 次元再解釈) | |
| (1) | $a_1[j_3,j_2,j_1] = a_0[j_1,j_2,j_3]$ | (行列の転置) |
|---|---|---|
| (2) | $a_2[j_3,j_2,k_3]={\displaystyle \sum_{j_1=0}^{p-1}} a_1[j_3,j_2,j_1]\omega_p^{j_1k_3}$ | (長さ $p$ の FFT を $qr$ 組) |
| (3) | $a_3[j_3,j_2,k_3] = a_2[j_3,j_2,k_3]\omega^{(j_2r+j_3)k_3}$ | (ひねり係数の乗算) |
| (4) | $a_4[j_3,k_3,j_2] = a_3[j_3,j_2,k_3]$ | (行列の転置) |
| (5) | $a_5[j_3,k_3,k_2] = {\displaystyle \sum_{j_2=0}^{q-1}} a_4[j_3,k_3,j_2] \omega_q^{j_2k_2}$ | (長さ $q$ の FFT を $pr$ 組) |
| (6) | $a_6[j_3,k_3,k_2] = a_5[j_3,k_3,k_2] \omega_{qr}^{j_3k_2}$ | (ひねり係数の乗算) |
| (7) | $a_7[k_3,k_2,j_3] = a_6[j_3,k_3,k_2]$ | (行列の転置) |
| (8) | $a_8[k_3,k_2,k_1] = {\displaystyle \sum_{j_3=0}^{r-1}} a_7[k_3,k_2,j_3] \omega_r^{j_3k_1}$ | (長さ $r$ の FFT を $pq$ 組) |
| (9) | $a_9[k_1,k_2,k_3] = a_8[k_3,k_2,k_1]$ | (行列の転置) |
| $A_{k_1pq+k_2p+k_3} = a_9[k_1,k_2,k_3]$ | （配列の 1 次元再解釈） |