読者です 読者をやめる 読者になる 読者になる

pekempeyのブログ

競技プログラミングに関する話題を書いていきます。

高速多項式除算

自分はそこまで多項式の扱いに詳しい訳ではないので、もっとシンプルで高速なアルゴリズムがあるかもしれません。でも、ひとつもアルゴリズムを知らないよりは知ってた方がいいかなと思って書きました。

概要

高速多項式乗算は競プロではよく使われる。FFT, NTT, Karatsuba あたりが有名だろう。

この記事では高速多項式除算を紹介する。高速多項式除算というのは $f(x)=q(x)g(x)+r(x)$ となる商 $q(x)$ を高速に求めるアルゴリズムのことである。 商を求められればあまりr(x)も高速に求められる。

応用としては高速きたまさ法が有名かと思う。

除算を乗算に変換

除算は乗算を用いて計算できる。$n$ 次多項式 $f(x)=a_0 + a_1x + \cdots + a_n x^n $ に対して、

$$f^{\ast}(x)=a_n+a_{n-1}x+\cdots+a_0x^n$$

と定義する。これは $f(x)=x^nf(1/x)$ とも書ける。$f(x)$を$ m $ 次多項式 $g(x)$ で割った商を $q(x)$、余りを $r(x)$ とすると、$f(x)=q(x)g(x)+r(x)$ が成り立つ。ここで両辺に $\ast$ 変換を掛けると、

$$ f^{\ast}(x)=q^{\ast}(x)g^{\ast}(x)+x^{n-\mathrm{deg}(r)}r^{\ast}(x) $$

となる。いま、$\mathrm{deg}(r) \lt m $ なので、$n-\mathrm{deg}(r) \gt n - m $ である。ゆえに、

$$ f^{\ast}(x)\equiv q^{\ast}(x)g^{\ast}(x) \;(\mathrm{mod}\; x^{n - m + 1}) $$

が成り立つ。$g^{\ast}(x) t \equiv 1 \;(\mathrm{mod}\; x^{n - m + 1})$ となる多項式 $t$ を両辺に掛けることにより、

$$ tf^{\ast}(x)\equiv q^{\ast}(x) \;(\mathrm{mod}\; x^{n - m + 1}) $$

が得られる。いま $\mathrm{deg}(q) \le n-m $ であるから、

$$q^{\ast}(x) = tf^{\ast}(x) \bmod x^{ n+ m - 1 } $$

が成り立つ。$q(x)$ を求めるために多項式乗算と $x^k$ での除算を用いるが、どちらも高速に処理できるものである。 問題は $g^{\ast}$ の逆元の求め方であるが、どのようにすればいいだろうか。実は $ n - m + 1 $ が 2 の冪であるときは簡単な求め方が存在する。そのため、実際には $ n - m + 1 $ が 2 冪になるように調整する。

$g^{\ast}$ の逆元の求め方

$g(x)=a_0+a_1x+\cdots a_{n}x^{n}$ とする。 $\mathrm{mod}\; x$ における $g$ の逆元は $1/a_0$ である。$\mathrm{mod}\; x^k$ における $g$ の逆元 $t$ が得られていたとしよう。このとき $\mathrm{mod}\; x^{2k}$ における逆元は $2t-t^2 g$ である。このことを用いると法が $x,x^2,x^4,x^8,\ldots$ のときの逆元は帰納的に求められる。

証明:$\mathrm{mod} x^k$ における $g$ の逆元を $t$ とする。このとき $gt=qx^k+1$ と書ける。 $$g(2t-t^2g)= 2(qx^k+1)-(q^2x^{2k}+2qx^k+1)=-q^2x^{2k}+1$$ であるから、$2t-t^2g$ は $x^{2k}$ における $g$ の逆元である。

多項式乗算の計算量を $M(n)$ とすると、逆元の計算量は $M(1)+M(2)+M(4)+\cdots+M(n)$ である。たとえば $M(n)=O(n\log n)$ であれば $O(n \log n)$ で逆元も求められる。

高速多項式剰余

$r=f-qg$ であるから、$q$ を求めれば求められる。商を求めるのには乗算 1 回が必要で、さらに $qg$ を計算するので、余りを求めるのには乗算が 2 回必要になる。これに加え $g^{\ast}$ の逆元の計算も加わるが、これは場合によっては使いまわせる。

応用

複数個の点での多項式評価という問題に使われる。これは $n$ 次多項式$f$ が与えられるので、$f(x_0),\ldots,f(x_{n-1})$ での値を求めよ、というものである。Horner 法を用いたとしても $O(n^2)$ となってしまうが、高速多項式剰余を用いることで $O(n\log ^2 n)$ で行える。

どのようにするのだろうか。

$P_{ij}=\prod _{k=i}^{j}(x-x_k)$ と $Q_{ij}(x)=f(x) \bmod P_{ij}(x)$ というものを定義する。

まず以下のことが容易にわかる。

$$ Q_{kk}(x) = f(x_k) $$ $$ Q_{0,n-1}(x)=f(x) $$

$i \le k \le j$ では以下のことが成り立つ。

$$ Q_{ik}(x) = Q_{ij}(x) \bmod P_{ik}(x) $$ $$ Q_{kj}(x) = Q_{ij}(x) \bmod P_{kj}(x) $$

以上のことを踏まえると次のような処理で計算できる。

$Q_{0,n-1}$ は $f(x) $ である。

$$Q_{0,(n-1)/2}=Q_{0,n-1} \bmod P_{0,(n-1)/2}$$ $$Q_{n/2,n-1}=Q_{0,n-1} \bmod P_{n/2,n-1}$$

である。これを再帰的に繰り返せばよい。$n$ 次の計算量を $T(n)$ と書くことにすると、多項式剰余は $O(n\log n)$ でできるから $T(n)=2T(n/2)+O(n\log n)=O(n\log^2 n)$ となる。

オーダーは確かに良いけど、実用的かどうかは分からない。かなり高速な多項式剰余を持ってるなら速くなるのかな…?

参考

http://www.math.kobe-u.ac.jp/Asir/ca.pdf