FFT と多項式乗算

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

多項式乗算と FFT

$\renewcommand{\vec}[1]{\boldsymbol{#1}}$各係数 $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{align*}
f(x) &= f_L(x) + 2^{15} \cdot f_H(x) \\
g(x) &= g_L(x) + 2^{15} \cdot g_H(x)
\end{align*}
$$
のように書き直し、$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$ 個の複素数の列 $\vec{z} = (z_0, \ldots, z_{n-1})$ に対する離散フーリエ変換は、フーリエ行列 $F_n$ を
$$
(F_n)_{i,j} = \zeta^{(i-1)(j-1)}_n
$$
で定義するとき、$F_n \vec{z}$ を計算することと同じとなる。ここで、$\zeta_n = e^{2 \pi i / n}$ とする(本当は $e^{-2 \pi i / n}$ なのだろうけど、畳込みを考える場合ではどちらでも計算できる)。

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

計算方法

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

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

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

実装

資料

その他

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