44
Help us understand the problem. What are the problem?

More than 1 year has passed since last update.

posted at

updated at

【競技プログラミング】形式的冪級数の応用テクニック(前編)

0.はじめに

形式的冪級数は、無限に続く収束性を考えない関数$f(x)$を考えて、数え上げなどの問題を解くテクニックです。

最近はyukicoderなどでの出題が増えており、atcoderで出されることがある重要なテクニックになりつつあります。

今回はそんな形式的冪級数の応用的なテクニックをまとめていきます。

以下$[x^k]f(x)$を$f(x)$の$x^k$の係数の意味で使います

わからないところや間違っていることがあればtwitterまでお願いします

0.1. 前提知識

フーリエ変換
形式的冪級数の基礎

0.2. 前編で解説する内容

$1/f(x)$verify
$e^{f(x)}$verify
$\sqrt{f(x)}$verify
微分積分
$\log(f(x))$verify
$f(x)^k$verify
商と余り
$f(x+c)$verify
$f(g(x))$verify
部分和問題verify

0.3. 定義

形式的冪級数においてはf(x)について、f(c)や、g(x)の定数項が0でない時のf(g(x))が定義されないです(この制約のおかげで無限項の計算を有限項に丸めても問題がなくなります)
$[x^k]f(x)$:$f(x)$の$x^k$の項
$1/f(x)$ :$g(x)\times f(x)=1$となるような$g(x)$
$e^{f(x)}:=\sum_{k=0}^{\infty}\frac{f(x)^k}{k!}$($[x^0]f(x)=0$の時のみ定義されます)
$\sqrt{f(x)}$:$g(x)\times g(x)=f(x)$となる$g(x)$(一般には定数項が1のf(x)以外には定義されないらしいのですが、今回はこのように定義します)
$\log(1+f(x)):=-\sum_{k=1}^{\infty}\frac{(-1)^k}{k}f(x)^k$($[x^0]f(x)=0$の時のみ定義されます)
$f(x+c):=\sum_{k=0}^{\infty}\frac{[x^k]f(x)}{k!}(x+c)^k$(形式的冪級数においてではなく多項式において定義されます)
$f(g(x)):=\sum_{k=0}^{\infty}\frac{[x^k]f(x)}{k!}g(x)^k$(形式的冪級数においては$[x^0]g(x)=0$を要請します(多項式においては要請されないです))
$f'(x):=\sum_{k=0}^{\infty}k([x^k]f(x)) x^{k-1}$
$\int f(x)dx:=\sum_{k=0}^{\infty}\frac{x^{k+1}}{k+1}$

によって定義されます

形式的冪級数における定義は基本的に$f(x)$に対する操作について操作後の$x^0 \sim x^N$項目までが、操作前の$f(x)$の$x^0 \sim x^N$項目までによって一意に定まるように定義されています

1.記号的ニュートン法

1.0.概要

形式的冪級数におけるニュートン法は
先頭$1$項を求める
先頭$2$項を求める
先頭$4$項を求める
$\vdots$
先頭$2^n$項を求める

を繰り返すことで
$O(1\log1+2\log2+4\log4+\dots+(N/2)\log(N/2)+N\log N)=O(N\log N)$で
$N$項目まで正確な値を求めるテクニックです

1.1.逆数の場合

$g_n(x)$を$0\sim 2^n-1$項目までの$1/f(x)$の値とします

$g_{n+1}(x)=g_{n}(x)(2-g_{n}(x)\times f(x) \bmod x^{2^{n+1}}) \bmod x^{2^{n+1}}$
が成り立ちます

  • 証明

定義より$g_{n}(x)\times f(x) \bmod x^{2^n}=1$なので、
$1+h_{n}(x)x^{2^n}=g_{n}(x)\times f(x) \bmod x^{2^n}$と置くことができ

$g_{n+1}(x)\times f(x)\bmod x^{2^{n+1}}$
$=g_{n}(x)\times f(x)(2-g_{n}(x)\times f(x) \bmod x^{2^{n+1}})\bmod x^{2^{n+1}}$
$=(1+h(x)x^{2^n})(2-(1+h(x)x^{2^n})) \bmod x^{2^{n+1}}$
$=1-h(x)^2x^{2^{n+1}} \bmod x^{2^{n+1}}$
$=1$
よって数学的帰納法により示せました

1.2 色々な場合

同様の理屈によって

求めたいもの 初項 漸化式
$1/f(x)$ $g_0(x)=1/([x^0]f(x))$ $g_{n+1}(x)=g_n(x)(2-g_n(x)*f(x)) \bmod x^{2^{n+1}}$
$e^{f(x)}$ $g_0(x)=1$ $g_{n+1}(x)=g_n(x)(f(x)+1-\log(g_n(x)))\bmod x^{2^{n+1}}$
$\sqrt{f(x)}$ $g_0(x)=\sqrt{[x^0]f(x)}$ $g_{n+1}(x)=(g_n(x)+f(x)/g_n(x))/2\bmod x^{2^{n+1}}$

が成り立ちます
これは所謂普通のニュートン法と同様の式になっています
注意点として
$e^{f(x)}=\sum_{k=0}^{\infty}f(x)^k/k!$
より、$[x^0]f(x)$が$0$でなければ任意の自然数$N$に対し$[x^0]e^{f(x)}$が$f(x)$の先頭$N$項から決まらないので$e^{f(x)}$は定義されないです
また、$[x^0]f(x)=0$の時、この方法では$\sqrt{f(x)}$が計算出来ないので、$\sqrt{x^2f(x)}=x\sqrt{f(x)}$を使って定数項を0以外にしてください

$e^{f(x)}$の漸化式に出てくる$\log$の計算方法は後述します

2.微分積分

難しそうに見える両者ですが、実は簡単です。
微分は$\frac{d(\sum_{k=0}^{\infty}x^k)}{dx}=\sum_{k=0}^{\infty}k\times x^{k-1}$
積分は$\int (\sum_{k=0}^{\infty}x^k)dx=\sum_{k=0}^{\infty}\frac{x^{k+1}}{k+1}$
によって定義されるため、そのとおりに配列を計算するだけで、$O(N)$です

3.log

今までの知識でできます
$\log(f(x))=\int \frac{f'(x)}{f(x)}dx$
は形式的冪級数においても成り立つので、これを計算すれば$O(N\log N)$です

  • 証明

注意点としては
定義より、$[x^0]f(x)=1$でないといけません

4.(多項式としての)商/余り

$(1/f(x))\times g(x)$との違いは切り捨てる方向です
$(1/f(x))\times g(x)$は$x^N$項目より大きい項を切り捨てる割り算と考えられますが、
$g(x)$を$f(x)$で割った商は$x^0$項目より小さい項を切り捨てる割り算と考えられます
よって$x^{-1}$を代入して$(1/f(x^{-1})) x^{deg(f)-deg(g)}\times g(x^{-1})$を計算することによって商は手に入ります。
実装上は$x^{-k}$を扱えないため、適当に下駄を履かせます
余りは $g(x)-(商)\times f(x)$です

$f(x)\div g(x)=h(x)+(A_1x^{-1}+A_2x^{-2}+A_3x^{-3}+\cdots)$
と表せると考えて
$f(x)=h(x)g(x)+g(x)(A_1x^{-1}+A_2x^{-2}+A_3x^{-3}+\cdots)$
とした時に、
$g(x)(A_1x^{-1}+A_2x^{-2}+A_3x^{-3}+\cdots)$の$deg$は$deg(g)$未満なので余りの要件を満たすと考えると
$f(1/x)/g(1/x)=h(1/x)+(A_1*x+A_2*x^2+A_3*x^3+)$
となるため正当性がわかりやすいと思います

5.累乗

繰り返し二乗法により$O(N\log^2N)$ですが、$O(N\log N)$でもできます
$e^{k\times \log(f(x))}$
を使います
$\exp$と$\log$の適用条件に気をつけてください

6. f(x+c)

$f(x)=\sum_{k=0}^{\infty}\frac{A_k}{k!}x^k$とすると
$f(x+c)$
$=\sum_{k=0}^{\infty}\frac{A_k}{k!}(x+c)^k$
$=\sum_{k=0}^{\infty}\sum_{i=0}^{k}\frac{A_k}{k!}\binom{k}{i}c^{k-i}x^{i}$
$=\sum_{k=0}^{\infty}\sum_{i=0}^{k}A_k\frac{c^{k-i}}{(k-i)!}\frac{x^i}{i!}$
$=\sum_{i=0}^{\infty}\sum_{k'=0}^{\infty}(A_{k'+i}\frac{c^{k'}}{k'!}) \frac{x^i}{i!}$ ($k'=k-i$とおいた)
ここで$(\sum_{k=0}^{\infty}A_kx^k) \times (\sum_{k=0}^{\infty} \frac{c^k}{k!} x^{-k})=\sum_{i=0}^{\infty}\sum_{k'=0}^{\infty}(A_{k'+i}\frac{c^{k'}}{k'!})x^i$
より、求める事ができます
$x^{-k}$の項に関しては実装上は$x^N$などを掛けて下駄を履かせてあげるといいです

7. f(g(x))

7.1.テイラー展開

テイラーの定理より多項式f,定数aに対して、
$f(x+a)=f(a)+f'(a)x+\frac{f''(a)}{2!}x^2+...$
であるが、実は$x$と$a$に形式的冪級数を代入することで、形式的冪級数に拡張することが出来、
$g(x)=p(x)+q(x)x^k$として、$a$に$p(x)$、$x$に$q(x)x^k$を代入すると、$f(p(x)+q(x)x^k)=f(p(x))+f'(p(x))q(x)x^k+\frac{f''(p(x))}{2!}q^2(x)x^{2k}+...$
式の形より$\lceil \frac{N}{k}\rceil$項目まで計算すれば十分である事がわかります
ここでは$k=\sqrt{N/\log N}$とします

7.2.分割統治法(f(p(x))の計算)

$f(x)=s(x)+t(x)x^{\frac{|f(x)|}{2}}$とすると
$f(p(x))=s(p(x))+p(x)^{\frac{|f(x)|}{2}}t(p(x))$
となり、これを再帰的に計算することで$O((N\log N)^\frac{3}{2})$で計算できます

7.3. f(p(x))からf'(p(x))を計算する

二項目以降も同様に計算していたらオーダーが大きくなってしまうのですが、
そこで、$(f(p(x)))'=f'(p(x))p'(x)$を利用して$f'(p(x))=(1/p'(x))*(f(p(x)))'$
とすることで高速化ができます
$O(\sqrt{N\log N})$回$O(N\log N)$の操作を行うので計算量は$O((N\log N)^\frac{3}{2})$です
よってf(g(x))は$O((N\log N)^{3/2})$で計算できました

8. 部分和問題

集合Sに含まれる要素を1つずつ選んでMを作る方法の場合の数を$O(|S|+M\log M)$
で求められます
$[x^M]\prod_{e \in S}(1+x^e)$が求めたい数ですが、
これは$[x^M]e^{\sum_{e \in S}\log(1+x^e)}$です
ここで、logの和の計算が大変そうに見えますが、$\log(1+x^e)$は$x^{ke}$の項以外は$0$なので重複要素を同一視すれば調和級数で抑えられます
よって求められました

9. おわりに

今回は計算に焦点を当てて解説しました
形式的冪級数はとても奥が深くて楽しいので是非、色々研究してみてほしいです。
後編はスターリング数などの数え上げに関する数や、多項式補完やmultipoint evalutionなどの形式的冪級数から少し離れた重要な操作、Subset Convolutionや、包除原理などバラエティに富んだ話題を書きたいと考えています

10. 参考資料

形式的冪級数(Formal-Power-Series) | Luzhiled's memo
僕の形式的冪級数はこれをパク参考にしています
Operations on Formal Power Series - Codeforces

11. 謝辞

記事の公開後、maspyさん/えびちゃんさんの指摘をもとに記事の怪しいところ等を加筆修正しました
ありがとうございます

Why not register and get more from Qiita?
  1. We will deliver articles that match you
    By following users and tags, you can catch up information on technical fields that you are interested in as a whole
  2. you can read useful information later efficiently
    By "stocking" the articles you like, you can search right away
Sign upLogin
44
Help us understand the problem. What are the problem?