FFT と多項式乗算

任意 mod の多項式乗算を 4 回の FFT でやる方法。ここでは、mod は $2^{30}$ 以下を想定。

多項式乗算と FFT

各係数 $c_i$ が $0 \le c_i < \mathrm{mod}$ を満たす多項式 $f(x), g(x)$ が与えられたときに、これらを普通に FFT で畳み込みをしようとすると、最大の係数の計算に必要な bit 数 $$ 2 \log_2(\mathrm{mod} - 1) + \log_2 (\min\{\deg(f), \deg(g)\} + 1)$$ が double の精度である 53 bit を超えてしまうため、計算することができない。

そこで、$f(x), g(x)$ の係数を下位 15 bit と上位 15 bit とで分け、$$ \begin{aligned} f(x) &= f_L(x) + 2^{15} \cdot f_H(x) \\ g(x) &= g_L(x) + 2^{15} \cdot g_H(x) \end{aligned} $$ のように書き直し、$f_L(x) g_L(x)$, $f_L(x) g_H(x)$, $f_H(x) g_L(x)$, $f_H(x) g_H(x)$ の 4 種の畳み込みの結果を用いて、$f(x), g(x)$ の畳み込みの結果を復元する。このとき、畳込みに必要となる bit 数は $30 + \log_2 (\min\{\deg(f), \deg(g)\} + 1)$ 程度で済み、$\deg(f) < 2^{18} \approx 250000$ 程度であれば、誤差が生じずに計算できる場合が多い。

これを普通に計算すると 12 回(使い回しで 8 回とか 6 回 ?)FFT が必要になるのだけど、うまくやると 4 回の FFT で計算することができる。

準備

$n$ 個の複素数の列 $\boldsymbol{z} = (z_0, \ldots, z_{n-1})$ に対する離散フーリエ変換は、フーリエ行列 $F_n$ を $$ (F_n)_{i,j} = \zeta^{(i-1)(j-1)}_n $$ で定義するとき、$F_n \boldsymbol{z}$ を計算することと同じとなる。ここで、$\zeta_n = e^{2 \pi i / n}$ とする(本当は $e^{-2 \pi i / n}$ なのだろうけど、畳込みを考える場合ではどちらでも計算できる)。

フーリエ行列の性質としては、$$ \begin{aligned} \overline{F_n} \boldsymbol{z} & = (F_n \boldsymbol{z})^{'} = F_n \boldsymbol{z}'\\ F_n^{-1} & = \frac{\overline{F_n}}{n} \end{aligned} $$ のようなものがあり、この性質が役に立つ。ここでは、$\overline{F}$ によって $F$ の共役を、 $\boldsymbol{z}'$ によって $\boldsymbol{z}$ の第 2 要素以降を逆順にしたものを表す: $\boldsymbol{z}' = (z_0, z_{n-1}, \ldots, z_1)$。

計算方法

$\boldsymbol{f}_L$, $\boldsymbol{f}_H$ で $f_L(x)$, $f_H(x)$ の係数の列を表すものとする。

$\boldsymbol{f} = \boldsymbol{f}_L + i \boldsymbol{f}_H$ とおくと、$F \boldsymbol{f} = F \boldsymbol{f}_L + i F \boldsymbol{f}_H$ となる(以後 $F_n$ を $F$ で表記する)。$F \boldsymbol{f}$ から $F\boldsymbol{f}_L$ と $F \boldsymbol{f}_H$ の値を別々で取得できればよいのだけど、$F \boldsymbol{f}$ に対して共役をとると、上での性質により $$ \begin{aligned} \overline{F \boldsymbol{f}} &= \overline{F} \boldsymbol{f}_L - i \overline{F} \boldsymbol{f}_H \\ &= (F \boldsymbol{f}_L)' -i (F \boldsymbol{f}_H)' \end{aligned} $$ が成り立ち、つまり、$$ \begin{aligned} F \boldsymbol{f} &= F \boldsymbol{f}_L + i F \boldsymbol{f}_H \\ (\overline{F \boldsymbol{f}})' &= F \boldsymbol{f}_L - i F \boldsymbol{f}_H \end{aligned} $$ が成り立つ。よって、$F \boldsymbol{f_L}, F \boldsymbol{f_H}$ の値は $F \boldsymbol{f}$ から個別に取得することができ、以下のフーリエ逆変換 $$ \begin{aligned} \boldsymbol{f} & \leftarrow F^{-1} ( (F \boldsymbol{f}_L) \circ (F \boldsymbol{g}_L) + i (F \boldsymbol{f}_L) \circ (F \boldsymbol{g}_H) ) \\ \boldsymbol{g} & \leftarrow F^{-1} ( (F \boldsymbol{f}_H) \circ (F \boldsymbol{g}_L) + i (F \boldsymbol{f}_H) \circ (F \boldsymbol{g}_H) ) \end{aligned} $$ を行うことができて、4 回分の畳込みの結果を 4 回の離散フーリエ変換で得られる。ここで、$\boldsymbol{f} \circ \boldsymbol{g}$ で $\boldsymbol{f}$ と $\boldsymbol{g}$ との要素同士の積(アダマール積)を表す。

まとめられる部分は先にまとめてしまってもよく、$$ \begin{aligned} \boldsymbol{f} & \leftarrow F^{-1} ( (F \boldsymbol{f}_L) \circ (F \boldsymbol{g}_L) + i (F \boldsymbol{f}_H) \circ (F \boldsymbol{g}_H) ) \\ \boldsymbol{g} & \leftarrow F^{-1} ( (F \boldsymbol{f}_L) \circ (F \boldsymbol{g}_H) + (F \boldsymbol{f}_H) \circ (F \boldsymbol{g}_L) ) \end{aligned} $$ のように計算することもできる。この方法だと $10$-bit ずつに分解しても 6 回の FFT で済み、より次数の大きい多項式乗算にも対応できる。

実装

資料

その他

  • 任意 mod の場合の多項式乗算の場合、 NTT を使うより この FFT の方法を用いる方が実装も軽く、おそらく速度も早い。
  • NTT の方が適しているのは、専用の mod で問われている場合だと思う。ただ、FFT は誤差が生じる可能性があるので、NTT はそういった面で安心。
  • FFT を覚えていれば)少し付け足すぐらいなので、空で書くのも比較的楽な気がする。