https://noshi91.hatenablog.com/entry/2023/12/10/163348

載ってないテクニックを募集しています。

前提

$\omega _ n := \exp \left ( 2 \pi i / n \right )$ と置く。以降 $i$ と書いたらそれは虚数単位ではなく添え字である。

数列 $\left ( a _ i \right ) _ { i = 0 } ^ { n - 1 }$ に対し、その DFT $\left ( \hat a _ j \right ) _ { j = 0 } ^ { n - 1 }$ を $\hat a _ j = \sum \limits _ { i = 0 } ^ { n - 1 } \omega _ n ^ { i j } a _ i$ で定義する。

基本事項

以下の事実は式変形で確認できる。

関係式 $\left ( 1 \right )$

$\left ( a _ i \right ) _ { i = 0 } ^ { 2 n - 1 }$ に対し、 $\left ( e _ i \right ) _ { i = 0 } ^ { n - 1 } , \left ( o _ i \right ) _ { i = 0 } ^ { n - 1 }$ を $e _ i = a _ i + a _ { n + i } , o _ i = \omega _ { 2 n } ^ i \left ( a _ i - a _ { n + i } \right )$ で定義する。このとき

\begin{align} \hat a _ { 2 j } &= \hat e _ j \cr \hat a _ { 2 j + 1 } &= \hat o _ j \end{align}

が成り立つ。すなわち、長さ $2 n$ の DFT の偶数部分と奇数部分は、それぞれが長さ $n$ のある列の DFT として与えられる。

参考: Cooley-Tukey Algorithm

関係式 $\left ( 2 \right )$

$\left ( a _ i \right ) _ { i = 0 } ^ { 2 n - 1 }$ に対し、 $\left ( e _ i \right ) _ { i = 0 } ^ { n - 1 } , \left ( o _ i \right ) _ { i = 0 } ^ { n - 1 }$ を $e _ i = a _ { 2 i } , o _ i = a _ { 2 i + 1 }$ で定義する。このとき

\begin{align} \hat a _ j &= \hat e _ { j \bmod n } + \omega _ { 2 n } ^ j \hat o _ { j \bmod n } \end{align}

が成り立つ。すなわち、長さ $2 n$ の DFT は、元の列の偶数部分と奇数部分それぞれの DFT を基に計算できる。

参考: https://hcpc-hokudai.github.io/archive/math_fft_001.pdf

DFT の使いまわし

もっとも基本的なテクニック。同じ列の DFT を複数回使うときは $1$ 度計算すれば十分である (それはそう)。やや気付きにくいケースとして、IDFT を行った列を後々もう一度 DFT し直している場合がある。多項式の $2$ 乗を求めるときや、有理式を通分して足し合わせるときなど、あらゆる場面で使う。

middle product

$\left ( a _ i \right ) _ { i = 0 } ^ { n + m - 1 }$ と $\left ( b _ i \right ) _ { i = 0 } ^ { n }$ を畳み込んだ列を $c$ として、 $\left ( c _ i \right ) _ { i = n } ^ { n + m - 1 }$ だけが必要となることがある。 DFT を用いた畳み込みは巡回畳み込みであるため普通は巡回しないよう十分な長さを取るが、この場合はいくらか巡回しても問題なく、長さが $n + m$ 以上の巡回畳み込みで足りる。特に $n = m = 2 ^ k$ の場合が頻出で、計算量が半分になる。ワイルドカードマッチング、転置原理を用いて導出されたアルゴリズム、オンライン畳み込みなどで用いる。

$2$ 倍の長さの DFT の計算

$\left ( a _ i \right ) _ { i = 0 } ^ { n - 1 }$ について $\hat a$ が求まっているとして、 $a$ を長さ $2 n$ に伸ばして $0$ 埋めした列 $\displaystyle \left ( b _ i \right ) _ { i = 0 } ^ { 2 n - 1 } , b _ i = \begin{cases} a _ i \quad \left ( i \lt n \right ) \cr 0 \quad \left ( i \geq n \right ) \end{cases}$ に対する $\hat b$ を計算したい場面がある。関係式 $\left ( 1 \right )$ より $\hat b _ { 2 j } = \hat a _ j$ であり、 $\left ( c _ i \right ) _ { i = 0 } ^ { n - 1 } , c _ i = \omega _ { 2 n } ^ i a _ i$ に対して $\hat b _ { 2 j + 1 } = \hat c _ j$ である。よって $\hat c$ を計算するだけでよい。 subproduct tree を作るときなどに用いる。

$f \left ( - x \right )$ の DFT の計算

$\left ( a _ i \right ) _ { i = 0 } ^ { 2 n - 1 }$ について $\hat a$ が求まっているとして、 $\displaystyle \left ( b _ i \right ) _ { i = 0 } ^ { 2 n - 1 } , b _ i = \left ( - 1 \right ) ^ i a _ i$ に対する $\hat b$ を計算したい場面がある。 $\hat b _ j = \sum \limits _ { i = 0 } ^ { 2 n - 1 } \omega _ { 2 n } ^ { i j } \left ( - 1 \right ) ^ i a _ i = \sum \limits _ { i = 0 } ^ { 2 n - 1 } \omega _ { 2 n } ^ { i \left ( j + n \right ) } a _ i = \hat a _ { \left ( j + n \right ) \bmod 2 n }$ であるから、 $\hat a$ から $\Theta \left ( n \right )$ で計算できる。 Bostan-Mori 法で用いる。

$f \left ( x ^ 2 \right )$ の DFT の計算

$\left ( a _ i \right ) _ { i = 0 } ^ { n - 1 }$ について $\hat a$ が求まっているとして、 $\displaystyle \left ( b _ i \right ) _ { i = 0 } ^ { 2 n - 1 } , b _ i = \begin{cases} a _ { i / 2 } \quad \left ( i \equiv 0 \pmod 2 \right ) \cr 0 \quad \left ( i \equiv 1 \pmod 2 \right ) \end{cases}$ に対する $\hat b$ を計算したい場面がある。関係式 $\left ( 2 \right )$ より $\hat b _ j = \hat a _ { j \bmod n }$ であるから、 $\hat a$ から $\Theta \left ( n \right )$ で計算できる。

前節の内容と組み合わせると以下の事実が分かる。次数 $2 n - 1$ の多項式 $g$ について $g \left ( x \right ) g \left ( - x \right ) = f \left ( x ^ 2 \right )$ が成り立っているとき、 $g$ の DFT の前半と後半を掛け合わせると $f$ の DFT が得られる。

$x f \left ( x \right )$ の DFT の計算

$\left ( a _ i \right ) _ { i = 0 } ^ { n - 1 }$ について $\hat a$ が求まっているとして、 $\displaystyle \left ( b _ i \right ) _ { i = 0 } ^ { n } , b _ i = a _ { \left ( i - 1 \right ) \bmod n }$ に対する $\hat b$ を計算したい場面がある。 $\hat b _ j = \sum \limits _ { i = 0 } ^ { n } \omega _ { n } ^ { i j } a _ { \left ( i - 1 \right ) \bmod n } = \sum \limits _ { i = 0 } ^ { n } \omega _ { n } ^ { \left ( i + 1 \right ) j } a _ i = \omega _ { n } ^ j \hat a _ j$ であるから、 $\hat a$ から $\Theta \left ( n \right )$ で計算できる。

偶数部分 (奇数部分) だけ IDFT

未知の $\left ( a _ i \right ) _ { i = 0 } ^ { 2 n - 1 }$ について $\hat a$ が分かっているとして、 $\left ( e _ i \right ) _ { i = 0 } ^ { n - 1 } , e _ i = a _ { 2 i }$ や $\hat e$ が求めたい場面がある。関係式 $\left ( 2 \right )$ と同様に、 $a$ の偶数部分と奇数部分をそれぞれ $e , o$ と定義する。関係式 $\left ( 2 \right )$ より、 $0 \leq j \lt n$ について $\hat a _ j = \hat e _ j + \omega _ { 2 n } ^ j \hat o _ j , \hat a _ { n + j } = \hat e _ j + \omega _ { 2 n } ^ { n + j } \hat o _ j = \hat e _ j - \omega _ { 2 n } ^ j \hat o _ j$ であるから、 $\hat e _ j = \left ( \hat a _ j + \hat a _ { n + j } \right ) / 2$ である。よって $\hat e$ は $\Theta \left ( n \right )$ で求まり、それを IDFT することで $e$ も求まる。 $o$ に関するものも同様に計算できる。

低次部分 (高次部分) の DFT の取り出し

未知の $\left ( a _ i \right ) _ { i = 0 } ^ { 2 n - 1 }$ について $\hat a$ が分かっているとして、 $\left ( \mathrm { lo } _ i \right ) _ { i = 0 } ^ { n - 1 } , \mathrm { lo } _ i = a _ i$ に対する $\hat { \mathrm { lo } }$ が求めたい場面がある。 $\left ( \mathrm { hi } _ i \right ) _ { i = 0 } ^ { n - 1 } , \mathrm { hi } _ i = a _ { n + i }$ と定義する。関係式 $\left ( 1 \right )$ より、 $\hat a$ の偶数部分は $\mathrm { lo } + \mathrm { hi }$ の DFT である。また、 $\hat a$ の奇数部分を IDFT することで $\omega _ { 2 n } ^ i \left ( \mathrm { lo } _ i - \mathrm { hi } _ i \right )$ が求まる。ここから $\mathrm { lo } - \mathrm { hi }$ が求まるため、これを DFT して $\hat a$ の偶数部分と足し合わせることで $2 \mathrm { lo }$ の DFT が求まる。 $\hat { \mathrm { hi } }$ も同様に計算できる。

長さ $2 n$ の IDFT ののち長さ $n$ の DFT を行うことと比べて $2 / 3$ 倍の時間計算量になる。

$f \left ( x ^ { - 1 } \right )$ や反転した列の DFT

$\left ( a _ i \right ) _ { i = 0 } ^ { n - 1 }$ について $\hat a$ が求まっているとして、 $\displaystyle \left ( b _ i \right ) _ { i = 0 } ^ { n } , b _ i = a _ { \left ( - i \right ) \bmod n }$ に対する $\hat b$ を計算したい場面がある。 $\hat b _ j = \sum \limits _ { i = 0 } ^ { n } \omega _ { n } ^ { i j } a _ { \left ( - i \right ) \bmod n } = \sum \limits _ { i = 0 } ^ { n } \omega _ { n } ^ { - i j } a _ i = \hat a _ { \left ( - j \right ) \bmod n }$ であるから、 $\hat a$ から $\Theta \left ( n \right )$ で計算できる。

IDFT の計算において回転子を $\omega _ n$ から $\omega _ n ^ { - 1 }$ に取り換えた DFT を計算するがこれも同じ式になるので、FFT の実装において IDFT 用に回転子の切り替えを可能にしておく必要はないことが分かる。もっとも、FFT に凝っている人は DFT の結果を bit 反転された状態で管理すると思われ、その場合 IDFT は別で実装することになるので回転子を一致させる旨味は少ないだろう。

$x f \left ( x \right )$ の DFT の計算と組み合わせると $a$ を反転した列の DFT も $\Theta \left ( n \right )$ で求まる。